CN111077525B - 融合sar与光学偏移量技术的地表三维形变计算方法及系统 - Google Patents

融合sar与光学偏移量技术的地表三维形变计算方法及系统 Download PDF

Info

Publication number
CN111077525B
CN111077525B CN201911324493.3A CN201911324493A CN111077525B CN 111077525 B CN111077525 B CN 111077525B CN 201911324493 A CN201911324493 A CN 201911324493A CN 111077525 B CN111077525 B CN 111077525B
Authority
CN
China
Prior art keywords
deformation
sar
east
south
north
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201911324493.3A
Other languages
English (en)
Other versions
CN111077525A (zh
Inventor
赵超英
陈立权
刘晓杰
陈恒祎
王宝行
丁明涛
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Changan University
Original Assignee
Changan University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Changan University filed Critical Changan University
Priority to CN201911324493.3A priority Critical patent/CN111077525B/zh
Publication of CN111077525A publication Critical patent/CN111077525A/zh
Application granted granted Critical
Publication of CN111077525B publication Critical patent/CN111077525B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/9021—SAR image post-processing techniques
    • G01S13/9023—SAR image post-processing techniques combined with interferometric techniques
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
    • G01B11/00—Measuring arrangements characterised by the use of optical techniques
    • G01B11/16—Measuring arrangements characterised by the use of optical techniques for measuring the deformation in a solid, e.g. optical strain gauge
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
    • G01B21/00—Measuring arrangements or details thereof, where the measuring technique is not covered by the other groups of this subclass, unspecified or not relevant
    • G01B21/32—Measuring arrangements or details thereof, where the measuring technique is not covered by the other groups of this subclass, unspecified or not relevant for measuring the deformation in a solid
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/86—Combinations of radar systems with non-radar systems, e.g. sonar, direction finder

Abstract

本发明公开了融合SAR与光学偏移量技术的地表三维形变计算方法及系统，包括：获取待监测区域形变前后的两景SAR数据，利用SAR偏移量跟踪方法根据两景单视复数影像计算待监测区域地表的方位向形变DAZI和视线向形变DLOS；获取待监测区域形变前后的两景光学，根据两景光学数据计算地表的南北向形变DNS和东西向形变DEW；最后通过公式(1)，获得地表的真实垂直向形变Dup、真实南北向形变Dns和真实东西向形变Dew。本发明方法和系统能够获得地表三维真实形变场，为诸如地震、滑坡等地质灾害的研究提供重要的技术支撑。

Description

融合SAR与光学偏移量技术的地表三维形变计算方法及系统
技术领域
本发明属于地表三维形变监测技术领域，具体涉及一种融合SAR与光学偏移量技术的地表三维形变计算方法及系统，可应用于地震、滑坡、火山等地质灾害监测领域的研究中。
背景技术
目前，遥感技术能够有效应用于针对诸如地震、滑坡、火山、地裂缝等地质灾害的调查、监测及预警工作中。其中基于SAR数据和光学数据的两类技术作为遥感领域重要的两个分支，在地质灾害研究中均充当了重要的角色。
针对光学的应用主要通过目视解译或多个波段数据的组合变换进而提取表征地质灾害体特征的因素，此类方法皆是对灾害体展开定性的研究，利用光学偏移量技术可以实现地质灾害所引起地表形变的定量分析，然而传统的方法往往只能获取地表二维形变。
基于SAR数据的地表形变监测技术主要包括合成孔径雷达干涉测量(InSAR)技术及SAR偏移量技术等，目前利用SAR数据的合成孔径雷达干涉测量(InSAR)及偏移量技术能够有效的应用于地表三维形变场的恢复，但为了获取研究区域地表三维形变，传统的方法至少需要三个来自于不同雷达成像几何结构的SAR数据。当可用SAR卫星数据源较少，造成该方法在实际应用中不足以满足三维形变场计算条件时，这些技术便无法开展。
发明内容
针对现有技术的缺陷和不足，本发明给出了一种融合SAR与光学偏移量技术的地表三维形变计算方法及系统，解决现有的地表三维形变解算在SAR数据源不足时使用受限的问题。
为了实现上述目的，本发明采用如下技术方案予以实现：
融合SAR与光学偏移量技术的地表三维形变计算方法，包括以下步骤：
步骤1，获取待监测区域形变前后的两景SAR数据，将两景SAR数据处理为两景单视复数影像，根据两景单视复数影像计算待监测区域地表的方位向形变DAZI和视线向形变DLOS；
步骤2，获取待监测区域形变前后的两景光学数据，两景光学数据与两景SAR数据覆盖区域相同，根据两景光学数据计算地表的南北向形变DNS和东西向形变DEW；
步骤3，对步骤1获得的方位向形变DAZI、视线向形变DLOS和步骤2获得的南北向形变DNS和东西向形变DEW进行重采样，使其具有相同的空间分辨率；
步骤4，通过公式(1)，获得地表的真实垂直向形变Dup、真实南北向形变Dns和真实东西向形变Dew；
Figure GDA0003924271020000021
式中，θ表示卫星入射角，
Figure GDA0003924271020000022
表示卫星飞行方位角。
具体的，所述的步骤1中利用SAR偏移量跟踪方法计算待监测区域地表的方位向形变DAZI和视线向形变DLOS。
具体的，所述的步骤2包括：
步骤2.1，对于获取的两景光学数据，判断两景光学数据的覆盖范围大小是否相同，若相同，进行步骤2.2；若不同，裁剪出两景数据相同的公共区域，对两景公共区域进行步骤2.2；
步骤2.2，使用基于ENVI软件的COSI-Corr光学偏移量计算工具获得地表的南北向形变量DNS和东西向形变量DEW。
优选的，所述的步骤2还包括对南北向形变DNS和东西向形变DEW的去除误差过程：根据得到的南北向形变DNS和东西向形变DEW掩模形变区域，在远离形变区域的位置选取若干像素点，根据这些像素点的形变值及地理位置利用一次拟合多项式拟合整个监测区域的轨道误差趋势面，然后将原始形变图与此误差趋势面相减，得到去除误差后的南北向形变DNS和东西向形变DEW。
优选的，所述的步骤4中，在使用公式(1)计算Dup、Dns和Dew时，采用赫尔默特方差分量估计的方法迭代来最终确定SAR与光学两类形变观测值的权，进而求得地表的真实垂直向形变量Dup、真实南北向形变量Dns和真实东西向形变量Dew。
本发明还公开了融合SAR与光学偏移量技术的地表三维形变计算系统，该系统包括：
SAR偏移量计算模块，用于获取待监测区域形变前后的两景SAR数据，将两景SAR数据处理为两景单视复数影像，根据两景单视复数影像，利用SAR偏移量跟踪方法计算待监测区域地表的方位向形变DAZI和视线向形变DLOS；
光学偏移量计算模块，用于获取待监测区域形变前后的两景光学数据，两景光学数据与SAR偏移量计算模块中两景SAR数据覆盖区域相同，根据两景光学数据计算地表的南北向形变DNS和东西向形变DEW；
重采样模块，用于对SAR偏移量计算模块获得的方位向形变DAZI、视线向形变DLOS和光学偏移量计算模块获得的南北向形变DNS和东西向形变DEW进行重采样，使其具有相同的空间分辨率；
地表三维形变计算模块，用于通过公式(1)，获得地表的真实垂直向形变Dup、真实南北向形变Dns和真实东西向形变Dew；
Figure GDA0003924271020000041
式中，θ表示卫星入射角，
Figure GDA0003924271020000042
表示卫星飞行方位角。
具体的，所述的SAR偏移量计算模块中利用SAR偏移量跟踪方法计算待监测区域地表的方位向形变DAZI和视线向形变DLOS。
具体的，所述的光学偏移量计算模块包括：
对于获取的两景光学数据，判断两景光学数据的覆盖范围大小是否相同，若相同，使用基于ENVI软件的COSI-Corr光学偏移量计算工具获得地表的南北向形变量DNS和东西向形变量DEW；若不同，裁剪出两景数据相同的公共区域，再使用基于ENVI软件的COSI-Corr光学偏移量计算工具生成地表的南北向形变量DNS和东西向形变量DEW。
进一步的，所述的光学偏移量计算模块还包括误差去除模块，用于去除南北向形变DNS和东西向形变DEW中的误差，具体为：根据得到的南北向形变DNS和东西向形变DEW掩模形变区域，在非形变区域的位置选取若干像素点，根据这些像素点的形变值及地理位置利用一次拟合多项式模拟整个待监测区域的轨道误差趋势面，然后将原始形变图与此误差趋势面相减，得到的去除误差后的南北向形变DNS和东西向形变DEW。
优选的，所述的地表三维形变计算模块中，在使用公式(1)计算Dup、Dns和Dew时，采用赫尔默特方差分量估计的方法迭代来最终确定SAR与光学两类形变观测值的权，进而求得地表的真实垂直向形变量Dup、真实南北向形变量Dns和真实东西向形变量Dew。
与现有技术相比，本发明的有益效果是：
本发明通过融合SAR与光学偏移量技术得到一种用于计算地表三维形变的方法，该方法的优点在于利用SAR与光学数据，通过偏移量技术获得地表形变，然后通过建立地表真实三维形变求解模型，且通过赫尔默特方差分量估计的方法最终确定SAR与光学两类形变观测值的权重，最终能够有效地获得地表三维真实形变，解决了现有仅依赖SAR数据的地表三维形变解算在SAR数据源不足时无法恢复地表三维形变的问题，为地表三维形变计算提供了一种新的方法，此工作对于研究地表形变具有重要意义。
附图说明
图1为本发明实施例记载的地表三维形变计算方法流程图。
图2为实施例记载的凯库拉地震沿方位向形变监测结果。
图3为实施例记载的凯库拉地震沿视线向形变监测结果。
图4为实施例记载的凯库拉地震在南北向上的形变监测结果，图中沿北方向形变为正，沿南方向形变为负。
图5为实施例记载的凯库拉地震在东西向上的形变监测结果，图中沿东方向形变为正，沿西方向形变为负。
图6为利用本发明所得凯库拉地震地表三维形变场；图中黑色箭头表示地表水平向形变，箭头方向为水平形变方向，长度表示形变大小，颜色变化表示垂直向形变。
以下结合说明书附图和具体实施方式对本发明做具体说明。
具体实施方式
以下给出本发明的具体实施方式，需要说明的是本发明并不局限于以下具体实施例中，凡在本申请技术方案基础上做的等同变换均落入本发明的保护范围。
如图1所示，本发明实施例提供的一种融合SAR与光学偏移量技术的地表三维形变计算方法包括以下步骤：
步骤1，获取待监测区域的两景SAR数据，如地震前后同一区域的两景数据；将两景SAR数据处理为两景单视复数影像，利用SAR偏移量跟踪方法根据两景单视复数影像计算待监测区域地表的方位向形变DAZI和视线向形变DLOS。
SAR偏移量跟踪方法的具体思路是：选定其中一景作为主影像，将另一景从影像配准重采样到该主影像，使主影像和从影像具有相同的空间几何结构；利用SAR幅度信息对两景SAR影像进行粗配准，计算主、从影像之间的初始偏移量，然后在粗配准的基础上，选取一定大小的搜索窗口并计算主从影像之间的归一化互相关系数ρ(x,y)，如式(2)所示，
Figure GDA0003924271020000061
式中，x和y为主影像坐标，xs和ys为两个方向上的偏移量，m和n为窗口大小，fm和fs分别为主影像幅度和从影像幅度。
当互相关系数ρ(x,y)达到最大时便可获得精确的子像素级的偏移量，最终获得地表的方位向形变DAZI和视线向形变DLOS。
步骤2，获取待监测区域形变前后的两景光学数据，两景光学数据与步骤1中两景SAR数据获取时间点大致相同、覆盖区域相同。根据两景光学数据计算地表的光学偏移量，即南北向形变DNS和东西向形变DEW；
需要注意的是，在利用两景光学数据计算地表的光学偏移量时，首先判断两景光学数据的覆盖范围大小是否相同，若相同，则使用基于ENVI软件的COSI-Corr光学偏移量计算工具生成地表的南北向形变量DNS和东西向形变量DEW；若不同，需要先裁剪出两景数据相同的公共区域，对两景公共区域再使用基于ENVI软件的COSI-Corr光学偏移量计算工具生成地表的南北向形变量DNS和东西向形变量DEW。
对于得到的南北向形变量DNS和东西向形变量DEW，由于一些光学影像数据没有经过严格的正射校正和几何校正，所以在计算所得形变场中存在明显的轨道误差，因此需要去除该误差，具体为：根据得到的南北向形变DNS和东西向形变DEW掩模形变区域，在非形变区域的位置选取若干像素点，根据这些像素点的形变值及地理位置利用一次拟合多项式拟合整个待监测区域的轨道误差趋势面，然后将原始形变图与此误差趋势面相减，得到去除误差后的南北向形变DNS和东西向形变DEW。
步骤3，对步骤1获得的方位向形变DAZI、视线向形变DLOS和步骤2获得的南北向形变DNS和东西向形变DEW进行重采样，使其具有相同的空间分辨率，并裁剪为相同的范围。
另外，若光学偏移量结果和SAR偏移量结果的坐标系不一致时，需要将两者转化为相同的坐标系，如光学偏移量结果通常为UTM坐标系，而经过地理编码后的SAR偏移量结果往往为WGS 84地理坐标系，因此首先需要将光学偏移量结果由UTM坐标系转为与SAR结果相同的地理坐标系。
步骤4，根据SAR卫星成像几何关系建立方位向形变和视线向形变与地表真实三维形变之间的关系式(3)和(4)，
Figure GDA0003924271020000081
Figure GDA0003924271020000082
式中，Dup表示地表的真实垂直向形变，Dns表示真实南北向形变，Dew表示真实东西向形变；
再结合光学所得南北向和东西向形变，可建立SAR与光学偏移量联合求解得到地表三维形变的观测方程，如式(1)所示。
Figure GDA0003924271020000083
式中，θ表示卫星入射角，
Figure GDA0003924271020000084
表示卫星飞行方位角。
最后通过公式(1)，获得地表的真实垂直向形变Dup、真实南北向形变Dns和真实东西向形变Dew。
在使用公式(1)计算Dup、Dns和Dew时，采用赫尔默特方差分量估计的方法通过迭代定权，最终确定SAR与光学形变观测值在计算过程中的权重，并最终求得地表三维形变场Dup、Dns和Dew。具体步骤如下：
步骤1)，首先确定初始方差，通过给定一个遍历窗口，分别计算方位向、视线向、东西向及南北向形变量文件中每个中心像元在该窗口内的方差作为初始方差，利用式(5)权计算公式计算初始权，式中c为任意给定常数，
Figure GDA0003924271020000091
为第i类观测方差，pi即为权值。令计算所得SAR与光学两类形变观测值的权分别为p1和p2；
Figure GDA0003924271020000092
步骤2)，令式(4)中L1＝[DLOS DAZI]T，L2＝[DNS DEW]T，
Figure GDA0003924271020000093
A为系数矩阵，A＝[A1 A2]T，其中，
Figure GDA0003924271020000094
Figure GDA0003924271020000095
将式(4)改写为误差方程，如式(8)所示，
Figure GDA0003924271020000096
由式(8)可得法方程为
Figure GDA0003924271020000097
式中，
Figure GDA0003924271020000098
Figure GDA0003924271020000099
步骤3)，最终建立残差平方和
Figure GDA00039242710200000910
和
Figure GDA00039242710200000911
之间的关系式如式(10)所示，其中，
Figure GDA00039242710200000912
S矩阵如式(11)所示，
Figure GDA00039242710200000913
Figure GDA00039242710200000914
根据公式(10)求得两类观测值迭代后单位权方差，计算公式如式(12)所示，
Figure GDA0003924271020000101
步骤4)，根据所求得两类观测值迭代后单位权方差，利用公式(13)再次定权，式中，k为迭代次数，
Figure GDA0003924271020000102
然后重复步骤2)和步骤3)，进行迭代，直至两类观测值迭代后的单位权方差相等迭代停止，最终求得
Figure GDA0003924271020000103
即求得地表三维形变。
本发明的实施例还提供了一种融合SAR与光学偏移量技术的地表三维形变计算系统，该系统包括：
SAR偏移量计算模块，用于获取待监测区域形变发生前后的两景SAR数据，将两景SAR数据处理为两景单视复数影像，利用SAR偏移量跟踪方法根据两景单视复数影像中计算待监测区域地表的方位向形变DAZI和视线向形变DLOS。SAR偏移量跟踪方法的具体思路同上文记载。
光学偏移量计算模块，用于获取待监测区域形变发生前后的两景光学数据，两景光学数据与SAR偏移量计算模块中两景SAR数据覆盖区域相同，根据两景光学数据，利用光学偏移量技术计算地表的南北向形变DNS和东西向形变DEW；
需要注意的是，在利用两景光学数据计算地表的光学偏移量时，首先判断两景光学数据的覆盖范围大小是否相同，若相同，则使用基于ENVI软件的COSI-Corr光学偏移量计算工具生成地表的南北向形变量DNS和东西向形变量DEW；若不同，需要先裁剪出两景数据的公共区域，对两景公共区域再使用基于ENVI软件的COSI-Corr光学偏移量计算工具生成地表的南北向形变量DNS和东西向形变量DEW。
对于得到的南北向形变量DNS和东西向形变量DEW，由于一些光学影像数据没有经过严格的正射校正和几何校正，所以在计算所得形变场中存在明显的轨道误差，因此还设置了误差去除模块，误差去除模块用于通过掩模得到的南北向形变DNS和东西向形变DEW形变区域，在非形变区域的位置选取若干像素点，根据这些像素点的形变值及地理位置利用一次拟合多项式拟合整个待监测区域的轨道误差趋势面，然后将原始形变图与此误差趋势面相减，得到去除误差后的南北向形变DNS和东西向形变DEW。
重采样模块，用于对SAR偏移量计算模块获得的方位向形变DAZI、视线向形变DLOS和光学偏移量计算模块获得的南北向形变DNS和东西向形变DEW进行重采样，使其具有相同的空间分辨率；
另外，若光学偏移量结果和SAR偏移量结果的坐标系不一致时，需要设置一个坐标系转化模块，该模块用于将两个不同的坐标系转化为相同的坐标系，如光学偏移量结果通常为UTM坐标系，而经过地理编码后的SAR偏移量结果往往为WGS 84地理坐标系，因此首先需要将光学偏移量结果由UTM坐标系转为与SAR结果相同的地理坐标系。
地表三维形变计算模块，用于通过公式(1)，获得地表的真实垂直向形变Dup、真实南北向形变Dns和真实东西向形变Dew；
Figure GDA0003924271020000111
式中，θ表示卫星入射角，
Figure GDA0003924271020000121
表示卫星飞行方位角。
在地表三维形变计算模块中，使用公式(1)计算Dup、Dns和Dew时，采用赫尔默特方差分量估计的方法通过迭代定权，最终确定SAR与光学形变观测值在计算过程中的权重，并最终求得地表三维形变场Dup、Dns和Dew。具体步骤同上文中的步骤1)至步骤4)。
实施例
本实施例选择位于新西兰东北部的凯库拉县地区作为实验区域，该地区于2016年11月14日发生了MW7.8级大地震，地震发生之后，地表破裂明显，且形变量级较大。实验选用的SAR数据为2016年2月20日和2016年12月10日覆盖该地区的两景地震发生前后Sentinel-1数据，光学数据为2016年2月19日和2017年2月23日覆盖该地区的两景地震发生前后Sentinel-2数据，由于本实施例仅研究该地震同震三维形变的计算，且同震形变量级远远大于震间形变，因此震后的SAR与光学数据获取时间节点相差时间内地表发生的震间形变可忽略不计，即可视为数据获取时间大致相同。
实验过程
步骤1，首先对于获取的两景Sentinel-1数据进行预处理生成单视复数影像(SLC)，然后选定2016年2月20日影像作为主影像，将另一景从影像配准重采样到该主影像，使主从影像具有相同的空间几何结构。然后利用SAR偏移量技术，设置搜索窗口为64*64，搜索步长为8*8，计算方位向和视线向二维形变场，分别如图2和图3所示。图中Azimuth为卫星飞行方向，即为方位向，LOS为沿雷达视线方向，即为视线向，远离卫星方向形变量为负，靠近卫星方向形变量为正。
步骤2，获取与SAR数据覆盖相同地区的2016年2月19日和2017年2月23日两景Sentinel-2数据，首先进行数据预处理，包括辐射定标和大气校正，最终提取第8波段进行光学偏移量的计算。由于数据覆盖范围大小不同，因此需裁剪出两景数据的公共区域。然后利用基于ENVI软件的COSI-Corr光学偏移量计算工具计算光学偏移量，设置搜索窗口大小为32*32，移动步长为8*8，掩膜阈值为0.9，迭代次数为2次，计算获得地表在南北向和东西向上的形变。由于一些光学影像数据没有经过严格的正射校正和几何校正，所以在计算所得形变场中存在明显的轨道误差，因此在非形变区域的位置选取若干像素点，根据这些像素点的形变值及地理位置利用一次拟合多项式模拟整个研究区域的轨道误差趋势面，然后将原始形变图与此误差趋势面相减，达到去除轨道误差的目的，最终获得地表在南北向和东西向上的形变如图4和图5所示。
步骤3，对SAR与光学偏移量结果进行重采样，使二者具有相同的空间分辨率和坐标系统。由于Sentinel-2光学偏移量的结果为UTM坐标系，因此首先需要将光学偏移量结果由UTM坐标系转为WGS 84下的地理坐标系，然后将SAR与光学结果重采样为相同的空间分辨率，并且裁剪为相同的范围，最终导出为TIFF格式文件。
步骤4，根据式(1)，通过步骤1)至步骤4)获得地表的真实垂直向形变Dup、真实南北向形变Dns和真实东西向形变Dew，如图6所示为地表的三维形变场。由图6中矢量箭头的方向及垂直方向的形变可以看出地震造成各板块之间相互作用，引起各板块的地表发生不同的形变，从而能够有效地解释地震所造成板块运动的机理。

Claims (6)

1.融合SAR与光学偏移量技术的地表三维形变计算方法，其特征在于，包括以下步骤：
步骤1，获取待监测区域形变前后的两景SAR数据，将两景SAR数据处理为两景单视复数影像，利用SAR偏移量跟踪技术，根据两景单视复数影像计算待监测区域地表的方位向形变DAZI和视线向形变DLOS；
步骤2，获取待监测区域形变前后的两景光学数据，两景光学数据与两景SAR数据覆盖区域相同，利用偏移量技术，根据两景光学数据计算地表的南北向形变DNS和东西向形变DEW；
步骤3，对步骤1获得的方位向形变DAZI、视线向形变DLOS和步骤2获得的南北向形变DNS和东西向形变DEW进行重采样，使其具有相同的空间分辨率；
步骤4，通过公式(1)，获得地表的真实垂直向形变Dup、真实南北向形变Dns和真实东西向形变Dew；
Figure FDA0003933231660000011
式中，θ表示卫星入射角，
Figure FDA0003933231660000012
表示卫星飞行方位角；
在使用公式(1)计算Dup、Dns和Dew时，采用赫尔默特方差分量估计的方法迭代来最终确定SAR与光学两类形变观测值的权，进而求得地表的真实垂直向形变量Dup、真实南北向形变量Dns和真实东西向形变量Dew。
2.如权利要求1所述的融合SAR与光学偏移量技术的地表三维形变计算方法，其特征在于，所述的步骤2包括：
步骤2.1，对于获取的两景光学数据，判断两景光学数据的覆盖范围大小是否相同，若相同，进行步骤2.2；若不同，裁剪出两景数据相同的公共区域，对两景公共区域进行步骤2.2；
步骤2.2，使用基于ENVI软件的COSI-Corr光学偏移量计算工具获得地表的南北向形变量DNS和东西向形变量DEW。
3.如权利要求1所述的融合SAR与光学偏移量技术的地表三维形变计算方法，其特征在于，所述的步骤2还包括对南北向形变DNS和东西向形变DEW的去除误差过程：根据得到的南北向形变DNS和东西向形变DEW掩模形变区域，在远离形变区域的位置选取若干像素点，根据这些像素点的形变值及地理位置利用一次拟合多项式拟合整个监测区域的轨道误差趋势面，然后将原始形变图与此误差趋势面相减，得到去除误差后的南北向形变DNS和东西向形变DEW。
4.融合SAR与光学偏移量技术的地表三维形变计算系统，其特征在于，该系统包括：
SAR偏移量计算模块，用于获取待监测区域形变前后的两景SAR数据，将两景SAR数据处理为两景单视复数影像，利用SAR偏移量跟踪方法，根据两景单视复数影像计算待监测区域地表的方位向形变DAZI和视线向形变DLOS；
光学偏移量计算模块，用于获取待监测区域形变前后的两景光学数据，两景光学数据与SAR偏移量计算模块中两景SAR数据覆盖区域相同，利用偏移量技术，根据两景光学数据计算地表的南北向形变DNS和东西向形变DEW；
重采样模块，用于对SAR偏移量计算模块获得的方位向形变DAZI、视线向形变DLOS和光学偏移量计算模块获得的南北向形变DNS和东西向形变DEW进行重采样，使其具有相同的空间分辨率；
地表三维形变计算模块，用于通过公式(1)，获得地表的真实垂直向形变Dup、真实南北向形变Dns和真实东西向形变Dew；
Figure FDA0003933231660000031
式中，θ表示卫星入射角，
Figure FDA0003933231660000032
表示卫星飞行方位角；
所述的地表三维形变计算模块中，在使用公式(1)计算Dup、Dns和Dew时，采用赫尔默特方差分量估计的方法迭代来最终确定SAR与光学两类形变观测值的权，进而求得地表的真实垂直向形变量Dup、真实南北向形变量Dns和真实东西向形变量Dew。
5.如权利要求4所述的融合SAR与光学偏移量技术的地表三维形变计算系统，其特征在于，所述的光学偏移量计算模块包括：
对于获取的两景光学数据，判断两景光学数据的覆盖范围大小是否相同，若相同，使用基于ENVI软件的COSI-Corr光学偏移量计算工具获得地表的南北向形变量DNS和东西向形变量DEW；若不同，裁剪出两景数据相同的公共区域，再使用基于ENVI软件的COSI-Corr光学偏移量计算工具生成地表的南北向形变量DNS和东西向形变量DEW。
6.如权利要求4所述的融合SAR与光学偏移量技术的地表三维形变计算系统，其特征在于，所述的光学偏移量计算模块还包括误差去除模块，用于去除南北向形变DNS和东西向形变DEW中的误差，具体为：根据得到的南北向形变DNS和东西向形变DEW掩模形变区域，在非形变区域的位置选取若干像素点，根据这些像素点的形变值及地理位置利用一次拟合多项式模拟整个待监测区域的轨道误差趋势面，然后将原始形变图与此误差趋势面相减，得到的去除误差后的南北向形变DNS和东西向形变DEW。
CN201911324493.3A 2019-12-20 2019-12-20 融合sar与光学偏移量技术的地表三维形变计算方法及系统 Active CN111077525B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201911324493.3A CN111077525B (zh) 2019-12-20 2019-12-20 融合sar与光学偏移量技术的地表三维形变计算方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201911324493.3A CN111077525B (zh) 2019-12-20 2019-12-20 融合sar与光学偏移量技术的地表三维形变计算方法及系统

Publications (2)

Publication Number Publication Date
CN111077525A CN111077525A (zh) 2020-04-28
CN111077525B true CN111077525B (zh) 2022-12-27

Family

ID=70316172

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201911324493.3A Active CN111077525B (zh) 2019-12-20 2019-12-20 融合sar与光学偏移量技术的地表三维形变计算方法及系统

Country Status (1)

Country Link
CN (1) CN111077525B (zh)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111736152B (zh) * 2020-08-17 2020-12-22 深圳大学 一种道路边坡稳定性监测方法及车载平台装置
CN112433213B (zh) * 2020-10-21 2023-05-23 中国电子科技集团公司第二十九研究所 一种sar干涉测量结果与光学影像位置偏移综合纠正方法
CN112685819A (zh) * 2020-12-29 2021-04-20 四川大学 一种用于大坝及滑坡变形gb-sar监测的数据后处理方法及系统
CN114964028B (zh) * 2022-05-24 2023-03-10 长安大学 一种综合遥感快速解译地震地表破裂带的方法

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US4975704A (en) * 1990-01-26 1990-12-04 The United States Of America As Represented By The Administrator Of The National Aeronautics And Space Administration Method for detecting surface motions and mapping small terrestrial or planetary surface deformations with synthetic aperture radar
CN104035131A (zh) * 2014-05-23 2014-09-10 中国科学院电子学研究所 一种地壳偏移量估计方法及装置
CN104237887A (zh) * 2014-09-29 2014-12-24 中国测绘科学研究院 一种sar遥感影像匹配方法
CN106338775A (zh) * 2016-09-07 2017-01-18 民政部国家减灾中心（民政部卫星减灾应用中心） 基于干涉合成孔径雷达数据的建筑物倒损程度评估方法
CN109711446A (zh) * 2018-12-18 2019-05-03 中国科学院深圳先进技术研究院 一种基于多光谱影像和sar影像的地物分类方法及装置
CN110058236A (zh) * 2019-05-21 2019-07-26 中南大学 一种面向三维地表形变估计的InSAR和GNSS定权方法
CN110109105A (zh) * 2019-04-23 2019-08-09 中国电力科学研究院有限公司 一种基于时序的InSAR技术监测地表形变的方法
CN110120046A (zh) * 2019-03-27 2019-08-13 长安大学 一种融合dem、光学遥感和形变信息的潜在滑坡识别方法

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US4975704A (en) * 1990-01-26 1990-12-04 The United States Of America As Represented By The Administrator Of The National Aeronautics And Space Administration Method for detecting surface motions and mapping small terrestrial or planetary surface deformations with synthetic aperture radar
CN104035131A (zh) * 2014-05-23 2014-09-10 中国科学院电子学研究所 一种地壳偏移量估计方法及装置
CN104237887A (zh) * 2014-09-29 2014-12-24 中国测绘科学研究院 一种sar遥感影像匹配方法
CN106338775A (zh) * 2016-09-07 2017-01-18 民政部国家减灾中心（民政部卫星减灾应用中心） 基于干涉合成孔径雷达数据的建筑物倒损程度评估方法
CN109711446A (zh) * 2018-12-18 2019-05-03 中国科学院深圳先进技术研究院 一种基于多光谱影像和sar影像的地物分类方法及装置
CN110120046A (zh) * 2019-03-27 2019-08-13 长安大学 一种融合dem、光学遥感和形变信息的潜在滑坡识别方法
CN110109105A (zh) * 2019-04-23 2019-08-09 中国电力科学研究院有限公司 一种基于时序的InSAR技术监测地表形变的方法
CN110058236A (zh) * 2019-05-21 2019-07-26 中南大学 一种面向三维地表形变估计的InSAR和GNSS定权方法

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
基于光学影像匹配技术的地震形变场提取研究;王颖洁 等;《矿山测量》;20120831(第4期);第20-24页 *
基于光学影像匹配技术的地震形变监测研究;陈志阳等;《矿山测量》;20160228(第01期);第77-79、83页 *
彭颖 等.融合光学与SAR遥感影像的同震地表三维形变监测适用性研究.《测绘与空间地理信息》.2019,第42卷(第10期),第86-89页. *
融合光学与SAR遥感影像的同震地表三维形变监测适用性研究;彭颖 等;《测绘与空间地理信息》;20191031;第42卷(第10期);第86-89页 *

Also Published As

Publication number Publication date
CN111077525A (zh) 2020-04-28

Similar Documents

Publication Publication Date Title
CN111077525B (zh) 融合sar与光学偏移量技术的地表三维形变计算方法及系统
CN110058236B (zh) 一种面向三维地表形变估计的InSAR和GNSS定权方法
WO2022213673A1 (zh) 融合无人机dom和星载sar影像的地表三维形变提取方法
Kääb et al. Coseismic displacements of the 14 November 2016 M w 7.8 Kaikoura, New Zealand, earthquake using the Planet optical cubesat constellation
CN108983232B (zh) 一种基于邻轨数据的InSAR二维地表形变监测方法
CN102607534B (zh) 基于运动恢复结构的卫星相对姿态测量方法
Giribabu et al. DEM generation using Cartosat-1 stereo data: issues and complexities in Himalayan terrain
CN107917880B (zh) 一种基于地基云图的云底高度反演方法
Sheng et al. Automated georeferencing and orthorectification of Amazon basin-wide SAR mosaics using SRTM DEM data
Li et al. A new analytical method for estimating Antarctic ice flow in the 1960s from historical optical satellite imagery
KR101804522B1 (ko) 다중 변위 추정 커널을 이용한 위성레이더 오프셋트래킹 장치 및 그 방법
CN110363758B (zh) 一种光学遥感卫星成像质量确定方法及系统
CN115060208A (zh) 基于多源卫星融合的输变电线路地质灾害监测方法及系统
CN111473802A (zh) 一种基于线阵推扫的光学传感器内方位元素定标方法
CN111650579B (zh) 一种岩移参数自适应获取的InSAR矿区三维形变估计方法、装置及介质
Chrysoulakis et al. Validation of ASTER GDEM for the Area of Greece
Ye et al. Improved geometric modeling of 1960s KH-5 ARGON satellite images for regional Antarctica applications
CN115201825A (zh) 一种InSAR震间形变监测中的大气延迟校正方法
CN115469308A (zh) 多轨道InSAR震间形变速率场拼接方法、装置、设备及介质
CN108646244B (zh) 测量建筑五维形变的分析方法及系统
CN111368716B (zh) 一种基于多源时空数据的地质灾害灾毁耕地提取方法
CN112505696A (zh) 一种基于星载sar的露天矿边坡形变动态监测方法
CN110927765A (zh) 激光雷达与卫星导航融合的目标在线定位方法
Restrepo et al. Building Lunar Maps for Terrain Relative Navigation and Hazard Detection Applications
CN107945235B (zh) 一种高轨大面阵静止卫星影像几何定位仿真方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant