CN111998766A - 一种基于时序InSAR技术的地表形变反演方法 - Google Patents

Abstract

本发明涉及一种基于时序InSAR技术的地表形变反演方法，包括以下步骤：地理编码及影像配准，生成差分干涉对，选取PS点，时空解缠，误差分离。具体的，在影像配准方式上，采用粗配准加精配准，保证TOPS成像模式下的配准精度需要；对于大气延迟误差，分别对分层大气和浮动大气利用不同的算法进行处理，最大程度上削弱大气误差的影响，提高形变反演的精度；若提供GNSS数据，可以利用该观测数据对大气延迟进行改正，同时在InSAR观测方程中作为约束条件进行联合平差，保证了形变结果的精度和可靠性，实现全天候、大面积、连续、高空间分辨率的地表形变监测。

Description

一种基于时序InSAR技术的地表形变反演方法

技术领域

本发明涉及地表形变监测领域，具体为一种基于时序InSAR技术的地表形变反演方法。

背景技术

近年来，由地表形变引起的地质灾害受到广泛关注，其成因包括地壳运动、岩性、褶皱和断裂构造等外力作用的自然因素，以及由于地下水、矿物过量开采等人为因素。严重的地表形变会对社会、经济造成严重的损失，因此开展形变监测和反演具有重要的意义。

传统地表形变的监测方法，尤其是地震同震形变场、城市地面沉降和矿区地表沉降多采用一、二等精密水准或者GNSS监测。精密水准测量适用于面积较小的工程测量范围，在早期的城市沉陷监测中得到了广泛应用，随着技术的发展，精度高、布网迅速和周期短的GNSS技术也被成功应用在形变监测中。但是无论水准测量技术还是GNSS技术，传统监测方法仍然存在着明显的缺陷，主要表现在(1)观测点容易受到损坏，造成后续观测成果缺失；(2)仅对点线进行观测，不能满足大范围监测要求。采用点线监测，通常需要进行内插以估算整体形变趋势，因此监测点的分布以及密度极大的影响了监测结果的精度和可靠性，若监测点过于稀疏，内插获取的形变面不能有效反映整体趋势；(3)传统观测方法劳动强度大、耗时长、工作效率低；(4)对于人力不可到达的复杂地形区域，可利用的监测点较少。因此，如何提供一种地表形变反演方法，在弥补传统监测手段的不足，大范围、高效率地进行地表形变监测，具有重要意义。

合成孔径雷达干涉测量技术(InSAR)是一种新兴的空间大地测量技术，相比于传统地表形变监测手段，其在监测范围、监测效率等方面有了很大的进步。结合了高分辨率成像技术、合成孔径雷达技术及干涉测量技术，InSAR目前已被广泛用于如地震，火山，滑坡，水土流失等民用、军用监测领域。该技术能够解决传统测量如水准、GNSS监测等面临的监测网络空间密度低、作业周期长等问题，能够全天候、大面积、连续、高空间分辨率地获取地表形变信息，实现全域监测。

传统差分干涉测量技术(DInSAR)由于其监测两成像时刻间的形变信息，在较长时间下两景影像的观测会受到如不同的大气环境、空间基线、散射体失相干等因素制约，限制了该技术的测量精度。时空失相关是对两次成像各自附加的噪声不同或不相关，相位差分无法抵消从而造成信噪比低，干涉图不明显，影响最终的形变反演精度；大气延迟则体现在易变大气条件会导致不同的相位延迟，可能会掩盖形变或地形信息。因此，如何克服上述时空失相关、减弱大气延迟对形变信息提取的影响，最大限度利用SAR影像数据建立长时间序列的差分干涉图是需要解决的问题。

PSInSAR(Permanent Scatter Interferometry)及SBAS(Small BaselineSubset)等时序InSAR技术通过选取高相干点、采用短基线减弱了时空失相关的影响。但是随着星载SAR平台的发展，不同成像方式对于影像配准也提出了更高要求，普通的多项式配准方法并不能达到精度要求。同时，对于大气误差的处理，现有的滤波处理方法并没有考虑分层大气和浮动大气的不同特征分别进行处理。另外，考虑到GNSS和InSAR数据在时间和空间分辨率上的互补特点，如何融合这两种技术得到高精度的形变序列，也是值得解决的问题。

发明内容

有鉴于此，为克服现有技术所存在的缺陷，本申请提供如下技术方案。

一种基于时序InSAR技术的地表形变反演方法，包括以下步骤：

S1、地理编码及影像配准，利用外部DEM进行几何配准，基于影像的幅度及强度信息进行相关配准，增强谱分集的精配准；

S2、生成差分干涉对；

S3、选取PS点，分别利用影像强度和相位信息进行PS点的迭代选取；

S4、时空解缠，在二维平面基础上增加时间维度，形成三维空间搜索路径，确定残差点进行三维解缠；

S5、误差分离，针对分层大气和浮动大气，分别采用参数化和时空滤波的方式进行大气误差的分离。

进一步的，若提供GNSS数据，上述地表形变反演方法还包括GNSS-InSAR融合，主要体现在大气延迟校正和带GNSS约束的InSAR数据平差。

优选地，所述步骤S1中利用外部DEM进行几何配准，具体地：

将外部DEM转换到主影像的影像坐标系下，建立查找表，然后通过卫星轨道参数计算出某一像素的位置，在主影像设置合适的窗口大小，均匀选取像素点，通过距离方程、多普勒方程、地球模型方程进行联立求解出像素点的坐标，反算这些像素点在副影像坐标系下的坐标；

利用Delaunay三角网内插，计算全部像素点的偏移量，利用最终的偏移量文件对副影像的改正。

优选地，所述步骤S1中基于影像的幅度及强度信息进行相关配准，具体地：

在主影像中选择合适的窗口，计算窗口的相干系数；

利用极大似然估计进行主副影像间的偏移量；

根据计算出的偏移量得到转换矩阵，进行副影像的配准，同时，确保达到0.05像元的配准精度。

优选地，所述步骤S1中增强谱分集的精配准具体地：

提取主影像和重采样后的副影像burst间的重叠区域；

主副影像干涉形成重叠区域的前后视图，再次干涉并滤波；

估计ESD相位和配准误差；

利用多项式拟合所有重叠区域偏差，计算方位向偏移量。

优选地，所述步骤S2生成差分干涉对包括：

主影像和副影像形成的时间基线、空间基线以及频率基线是影响干涉效果的重要因素，为挑选最佳的公共主影像，需要使时序干涉图的相关系数之和达到最大：

式中，f(γ_t)，f(γ_B)，f(γ_DC)分别为时间基线系数、垂直基线系数和多普勒质心频率差系数，找到最佳主影像后形成干涉对；

外部DEM转到雷达坐标系下形成地形相位，干涉相位减去地形相位以及平地相位，形成差分干涉对。

优选地，所述步骤S3选取PS点包括：

利用振幅离差系数进行粗选，保留高于阈值的点目标；

通过差分干涉后的差分相位、滤波相位以及地形残余相位来计算相位点在时间上的相干系数，通过选定合适的相干系数阈值，高于该阈值的点则被识别为PS目标点，该步骤需要通过滤波以及迭代运算实现对PS目标点的最终选取。

优选地，所述步骤S4中时空解缠包括：

在去除地形相位后，干涉对中单点相位构成表示为：

W{φ}＝W{φ_d+φ_τ+φ_a+φ_o+φ_n}

式中，W{·}为相位缠绕，φ为单点干涉相位，φ_d、φ_τ、φ_a、φ_o、φ_n分别表示形变、地形误差、大气、轨道、和噪声相位。

其中φ_d、φ_a和φ_o在空间具有相关性，可利用空间低通滤波将其分离，估算出空间相关相位

后，将其去除得到：

式中，上标u代表空间非相关，δ代表未去除的空间不相关的形变、大气、轨道误差的总量；因为这三者具有很强的空间相关性，所以δ值很小；

构造相干系数：

由于δ很小，相干系数γ可以反映单点受噪声的影响大小，噪声均值越大，相干系数值越大，利用空间搜索法可以求算出γ的最大值，并求出对应的

在去除

误差后，进行时空三维解缠。

优选地，步骤S5中误差分离包括：

解缠后的相位值为：

式中，φ_uw为解缠后的相位，

为残余具有空间相关性的地形残差相位，k为单差模糊度整数；如果三维解缠正确，则大部分点的k应该相等。为分离出形变相位，需要在弧段上进行时空域滤波。

式中，

代表两点间做差，然后构建三角网进行时间域低通滤波，将其分离后得到

式中，[·]_{LP_time}表示时间域低通滤波，

和

分别表示副影像大气和轨道分量。

在时间域上表现为较强的相关性，利用时间域低通可以获取形变和主影像大气；其余量在时间域上表现为不相关，要获取副影像大气，

在空间表现为相关性，

表现为随机噪声误差，因此可以通过空间的高通滤波分离出副影像的大气延迟相位和轨道误差；最后通过从弧段到点上的积分可以获取每个点的φ_d及

优选地，若误差分离中高低通滤波后差分相位仍然出现较为明显的长波趋势，即分层大气延迟干扰，则通过采用线性函数进行吸收并去除。

优选地，该地表形变反演方法还包括GNSS-InSAR融合。GNSS具有全天候、高精度、高效率等特点，能够长时间连续观测，提供高精度观测数据。由于SAR影像时间周期较长，很难提供足够时间分辨率的信息，利用GNSS-InSAR融合技术，可以在改正InSAR数据本身难以消除的误差(如大气延迟等)的同时，实现GNSS高时间分辨率、高位置精度与InSAR高空间分辨率的有效统一，从而得到高可靠性、高精度地表形变结果。GNSS和InSAR的技术融合主要体现在：GNSS应用于大气延迟校正；带GNSS约束的InSAR数据平差。

具体为，在获得解缠相位之后，若有GNSS外部数据，通过GNSS数据得到ZWD，反演大气延迟，将其从解缠相位中剔除；

将DEM误差相位加回，并重新进行参数化，利用最小二乘进行形变相位解算，将GNSS观测方程作为约束条件，和InSAR观测方程进行联合平差。

本发明所获得的有益技术效果：

1)本发明解决了传统的水准及GNSS形变监测手段监测周期长、空间点位稀疏、耗费人力物力，以及传统的DInSAR监测方法受时空失相干及大气延迟影响严重的问题。本发明分别在配准方法、大气延迟处理算法以及GNSS-InSAR融合技术方面具有较强优势，能够全天候、大面积、连续、高空间分辨率地监测地表形变，减弱误差的影响，提高监测效率；

2)本发明在影像配准方式上，采用粗配准加精配准的方法。在粗配准基础上，利用增强谱分解的配准方法，保证TOPS成像模式下的配准精度需要；

3)本发明对于大气延迟误差，分别对分层大气和浮动大气利用不同的算法进行处理，最大程度上削弱大气误差的影响，提高形变反演的精度；

4)若提供GNSS数据，可以利用GNSS数据对大气延迟进行改正，同时在InSAR观测方程中作为约束条件进行联合平差，保证了形变结果的精度和可靠性。

上述说明仅是本申请技术方案的概述，为了能够更清楚了解本申请的技术手段，从而可依照说明书的内容予以实施，并且为了让本申请的上述和其他目的、特征和优点能够更明显易懂，以下以本申请的较佳实施例并配合附图详细说明如后。

根据下文结合附图对本申请具体实施例的详细描述，本领域技术人员将会更加明了本申请的上述及其他目的、优点和特征。

附图说明

为了更清楚地说明本申请实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图是本申请的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。在所有附图中，类似的元件或部分一般由类似的附图标记标识。附图中，各元件或部分并不一定按照实际的比例绘制。

图1是本公开一种实施例中基于时序InSAR技术的地表形变反演方法的流程图；

图2是本公开一种实施例中影像配准流程图。

具体实施方式

为使本申请实施例的目的、技术方案和优点更加清楚，下面将结合本申请实施例中的附图，对本申请实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本申请一部分实施例，而不是全部的实施例。

此外，本申请可以在不同例子中重复参考数字和/或字母。这种重复是为了简化和清楚的目的，其本身并不指示所讨论各种实施例和/或设置之间的关系。

还需要说明的是，在本文中，诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含。

实施例1

如附图1所示，一种基于时序InSAR技术的地表形变反演方法，包括以下步骤：

S1、地理编码及影像配准，利用外部DEM进行几何配准，基于影像的幅度及强度信息进行相关配准，利用增强谱分集进行精配准，如附图2所示。

针对Sentinel-1A IW模式下的TOPS扫描方式，每个burst在进行扫描时，天线波束沿方位向从后向前旋转扫描，压缩了方位向波束方向图并且抑制了成像时所造成的扇贝效应。由于TOPS模式特殊的性导致多普勒质心偏移，使得数据在方位向配准精度要求为1/1000，因此采用粗配准和精配准联合的方式处理影像。对于粗配准选择利用SAR成像几何、轨道数据及外部DEM的几何配准方式，从而获得主副影像之间的偏移量，为避免在使用ESD进行精配准过程中出现的相位缠绕问题，应尽可能使粗配准精度达到0.05个像元。

精配准主要是为了改正由于轨道不平行引起的系统误差及TOPS模式方位向特殊性的问题。采用增强谱分集的方式，该方法主要在频率域处理一维信号，利用响应脉冲函数，对主辅影像方位向上两个相邻的burst的重叠区域分别做干涉处理，再对该干涉处理的结果再次进行干涉处理，到重叠区域的干涉相位差，最终获得方位向的配准误差。

所述利用外部DEM进行几何配准，具体地：

将外部DEM转换到主影像的影像坐标系下，建立查找表，然后通过卫星轨道参数计算出某一像素的位置，在主影像设置合适的窗口大小，均匀选取像素点，通过距离方程、多普勒方程、地球模型方程进行联立求解出像素点的坐标，反算这些像素点在副影像坐标系下的坐标；

利用Delaunay三角网内插，计算全部像素点的偏移量，利用最终的偏移量文件对副影像的改正。

所述基于影像的幅度及强度信息进行相关配准，具体地：

在主影像中选择合适的窗口，计算窗口的相干系数；

利用极大似然估计进行主副影像间的偏移量；

根据计算出的偏移量得到转换矩阵，进行副影像的配准，同时，确保达到0.05像元的配准精度。

所述增强谱分集的精配准，具体地：

提取主影像和重采样后的副影像burst间的重叠区域；

主副影像干涉形成重叠区域的前后视图，再次干涉并滤波；

估计ESD相位和配准误差；

利用多项式拟合所有重叠区域偏差，计算方位向偏移量。

S2、生成差分干涉对。

具体地，主影像和副影像形成的时间基线、空间基线以及频率基线是影响干涉效果的重要因素，为挑选最佳的公共主影像，需要使时序干涉图的相关系数之和达到最大：

式中，γ为相干系数，f(γ_t)，f(γ_B)，f(γ_DC)分别为时间基线系数、垂直基线系数和多普勒质心频率差系数，找到最佳主影像后形成干涉对。

外部DEM转到雷达坐标系下形成地形相位，干涉相位减去地形相位以及平地相位，形成差分干涉对。

S3、选取PS点，分别利用影像强度和相位信息进行PS点的迭代选取。

利用振幅离差系数进行粗选，保留高于阈值的点目标。

通过差分干涉后的差分相位、滤波相位以及地形残余相位来计算相位点在时间上的相干系数，通过选定合适的相干系数阈值，高于该阈值的点则被识别为PS目标点，该步骤需要通过滤波以及迭代运算实现对PS目标点的最终选取。具体的，对于振幅离差选点，当像元存在主散射体时，相位标准差与振幅存在如下关系：

式中，

是相位标准差；σ_A是幅度A的标准差；m_A为时序SAR影像的幅度均值；D_A为离差指数，当D_A小于阈值时，选取为PS待定点。

振幅离差在挑选振幅较强的PS点时具有很高的成功率，但是对于低信噪比的散射体，振幅离差和相位稳定之间的关系较低，因此还需要进一步选取。

相位分析法主要利用不同相位分量时间及空间相关性特征，采用低通滤波分离出噪声相位，选择受噪声影响较小的像元为PS点。

式中，

为滤波得到的相位；

为由基线引起的误差；N为SAR影像数量。γ反映了目标点在时间序列上的稳定性，γ越大表明待定PS点的噪声越小，因此选取高于阈值的点PS点。该阈值的确定通常需要迭代及滤波运算实现。

S4、时空解缠，在二维平面基础上增加时间维度，形成三维空间搜索路径，确定残差点进行三维解缠。具体包括：

在去除地形相位后，干涉对中单点相位构成表示为：

W{φ}＝W{φ_d+φ_τ+φ_a+φ_o+φ_n} (4)

式中，W{·}为相位缠绕，φ为单点干涉相位，φ_d、φ_τ、φ_a、φ_o、φ_n分别表示形变、地形误差、大气、轨道、和噪声相位。其中φ_d、φ_a和φ_o在空间具有相关性，可利用空间低通滤波将其分离，估算出空间相关相位

后，将其去除得到：

式中，上标u代表空间非相关，δ代表未去除的空间不相关的形变、大气、轨道误差的总量；因为这三者具有很强的空间相关性，所以δ值很小；

构造相干系数：

由于δ很小，相干系数γ可以反映单点受噪声的影响大小，噪声均值越大，相干系数值越大，利用空间搜索法可以求算出γ的最大值，并求出对应的

在去除

误差后，进行时空三维解缠。

传统的时序InSAR方法多单独采用时间维(1D)、空间维解缠(2D)算法，或者进行1D+2D伪三维解缠，无法同时考虑缠绕相位的时空特征。本实施例采用真三维(3D)解缠算法，即残差点搜索路径从平面扩展到空间。该方法在三维空间中搜索残差路径，生成残差点，相比于平面枝切法，该方法生成的三维积分路径更能够保证解缠的准确性。

三维解缠分别体现在空间和时间维度上，基于最小化L-P范框架，将相位解缠看做最优化问题，找到目标函数的最小值，从而获取解缠后的相位值。目标函数为：

式中，Δφ为解缠相位梯度，Δψ为缠绕相位梯度，x和y分别表示二维的两个方向，z表示第三维方向，i和j表示点坐标，w表示权重。

和传统的L-P范分布不同之处在于，这里加入了第三维度的相位梯度信息，通常为时间序列信息。这里，参数P决定如何处理Δφ和Δψ的差值关系，当P＝0时，若仅考虑二维信息，上式为类似枝切解缠算法的L-0范数目标函数；当P＝2时，即为最小二乘解缠算法。

针对非线性形变，在先验信息不充足的情况下，形变解算中通常假设形变为线性模型，然而实际的形变过程很难保证该假设成立，从而影响高程和形变参数的解算。在周期形变和复杂形变场景下，函数模型对DEM误差估计的影响程度较大。针对模型误差问题，采用迭代处理，利用地形、大气、形变不同的时空特征，通过时空域高低通滤波的方式，分离出形变信号，减弱了模型对于参数估计的影响，从而确保形变提取的可靠性。另外，形变信号具有空间相关性特征，相邻点形变相似。为避免先验形变模型引入的误差，也可以加入空间形变约束，直接反演形变量。

S5、误差分离，针对分层大气和浮动大气，分别采用参数化和时空滤波的方式进行大气误差的分离。

在影响InSAR地面沉降监测精度的诸多因素中，大气延迟是主要误差之一。考虑到测区气象条件的复杂多变，可能在InSAR干涉图内存在较为严重的大气延迟误差。大气延迟主要是对流层影响，对于对流层的分层大气，因其具有很强的规律性，可采用一个高程相关的参数加以吸收；浮动大气受水汽变化影响，模型估计精度较低。本实施例针对两部分延迟分别提供了处理方案。

分层大气：采用高程相位相关分析可以有效地分离垂直分层效应造成的地形相关的大气延迟相位。根据垂直分层效应的原理，一般认为大气延迟相位与高程具有线性关系，因此使用线性回归建立高程和相位的关系模型反演大气延迟相位，也可以使用二次或高次多项式进行拟合。传统的方法对该分层大气进行单独估计，但是由于形变和地形误差相位的存在会使得延迟参数K的估计存在误差；另外，该分层大气在空间的变化并不一致，因此单一参数无法应用于整幅干涉对的延迟量估计。针对于以上两个问题，采用分层大气延迟和形变、地形参数联合估计的方法，在模型中加入了大气参数进行吸收。

Δφ＝Δφ_d+τ+K×Δh

式中，Δ代表双差观测，Δφ_d+τ为形变和地形相位的和，Δh为高程差，为确保参数K估计的准确性同时减弱浮动大气的影响，采用四叉树算法，对整幅干涉对分区解算，即减少了运算负担也保证了参数估值的准确性。

浮动大气：该部分大气在空间和时间上表现为随机性，很难用确定性的模型进行描述。尤其低纬地区在夏季受大气影像严重，且多表现为浮动大气影响。因此在保证足够影像数量的基础上，可以适量删除个别严重失相关的干涉对以减少大气延迟的影响。对残余大气利用其统计特性进行分离，即根据大气延迟和形变相位在时空维的相关性差异，主要是通过时空滤波的方法来分离。

所述误差分离具体包括：

解缠后的相位值为：

式中，φ_uw为解缠后的相位，

为残余具有空间相关性的地形残差相位，k为单差模糊度整数。如果三维解缠正确，则大部分点的k应该相等。为分离出形变相位，需要在弧段上进行时空域滤波。

式中，

代表两点间做差，然后构建三角网进行时间域低通滤波，将其分离后得到

式中，[·]_{LP_time}表示时间域低通滤波，

和

分别表示副影像大气和轨道分量；

在时间域上表现为较强的相关性，利用时间域低通可以获取形变和主影像大气；其余量在时间域上表现为不相关，要获取副影像大气，

在空间表现为相关性，

表现为随机噪声误差，因此可以通过空间的高通滤波分离出副影像的大气延迟相位和轨道误差；最后通过从弧段到点上的积分可以获取每个点的φ_d及

进一步的，该地表形变反演方法还包括GNSS-InSAR融合。在获得解缠相位之后，若有GNSS外部数据，通过GNSS数据得到ZWD，反演大气延迟，将其从解缠相位中剔除。

GNSS具有时间采样率高的特征，能够提供更高时间分辨率和更高精度的ZWD估计，减轻甚至避免水汽产品与雷达成像时间不一致带来的不确定性。但是GNSS反映的是单点的ZTD估计，将其应用于InSAR大气相位校正需要进行空间插值，因此对GNSS站点的密度和分布有较高的要求。

将DEM误差相位加回，并重新进行参数化，利用最小二乘进行形变相位解算，将GNSS观测方程作为约束条件，和InSAR观测方程进行联合平差。具体地，以GNSS先验监测结果作为强约束，构建形变解算的约束条件，能够提高最小二乘的解算精度。该约束与InSAR数据的融合可以有三种方案：(1)以GNSS观测值作为条件约束，进行函数模型约束；(2)利用GNSS观测量的协方差矩阵进行附加随机模型约束的解算；(3)同时进行函数和随机模型双约束的融合解算。约束误差方程可以写为：

式中，L_GNSS和L_InSAR分别为GNSS和InSAR观测值，ε_GNSS和ε_InSAR为观测误差，X为待估参数，包含形变和地形参数，利用最小二乘解算得到最终的形变信息。

上述基于时序InSAR技术的地表形变反演方法，能够全天候、大面积、连续、高空间分辨率地监测地表形变。具体的，在影像配准方式上，采用粗配准加精配准，保证不同成像模式下的配准精度需要；对于大气延迟误差，分别对分层大气和浮动大气利用不同的算法进行处理，最大程度上削弱大气误差的影响，提高形变反演的精度；若提供GNSS数据，可以利用GNSS数据对大气延迟进行改正，同时在InSAR观测方程中作为约束条件进行联合平差，保证了形变结果的精度和可靠性。与传统的时序技术相比，本发明分别在配准方法、大气延迟处理算法以及GNSS-InSAR融合技术方面具有较强优势。

以上所述仅为本发明的优选实施例而已，其并非因此限制本发明的保护范围，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，通过常规的替代或者能够实现相同的功能在不脱离本发明的原理和精神的情况下对这些实施例进行变化、修改、替换、整合和参数变更均落入本发明的保护范围内。

Claims (10)

1.一种基于时序InSAR技术的地表形变反演方法，其特征在于，包括以下步骤：

S1、地理编码及影像配准，利用外部DEM进行几何配准，基于影像的幅度及强度信息进行相关配准，利用增强谱分集的方法进行精配准；

S2、生成差分干涉对；

S3、选取PS点，分别利用影像强度和相位信息进行PS点的迭代选取；

S4、时空解缠，在二维平面基础上增加时间维度，形成三维空间搜索路径，确定残差点进行三维解缠；

S5、误差分离，针对分层大气和浮动大气，分别采用参数化和时空滤波的方式进行大气误差的分离。

2.根据权利要求1所述的基于时序InSAR技术的地表形变反演方法，其特征在于，所述步骤S1中利用外部DEM进行几何配准，具体地：

将外部DEM转换到主影像的影像坐标系下，建立查找表，然后通过卫星轨道参数计算出某一像素的位置，在主影像设置合适的窗口大小，均匀选取像素点，通过距离方程、多普勒方程、地球模型方程进行联立求解出像素点的坐标，反算这些像素点在副影像坐标系下的坐标；

利用Delaunay三角网内插，计算全部像素点的偏移量，利用最终的偏移量文件对副影像改正。

3.根据权利要求1所述的基于时序InSAR技术的地表形变反演方法，其特征在于，所述步骤S1中基于影像的幅度及强度信息进行相关配准，具体地：

在主影像中选择合适的窗口，计算窗口的相干系数；

利用极大似然估计进行主副影像间的偏移量；

根据计算出的偏移量得到转换矩阵，进行副影像的配准，同时，确保达到0.05像元的配准精度。

4.根据权利要求1所述的基于时序InSAR技术的地表形变反演方法，其特征在于，所述步骤S1中增强谱分集的精配准，具体地：

提取主影像和重采样后的副影像burst间的重叠区域；

主副影像干涉形成重叠区域的前后视图，再次干涉并滤波；

估计ESD相位和配准误差；

利用多项式拟合所有重叠区域偏差，计算方位向偏移量。

5.根据权利要求1所述的基于时序InSAR技术的地表形变反演方法，其特征在于，所述步骤S2生成差分干涉对包括：

主影像和副影像形成的时间基线、空间基线以及频率基线是影响干涉效果的重要因素，为挑选最佳的公共主影像，需要使时序干涉图的相关系数之和达到最大：

式中，γ为相干系数，f(γ_t)，f(γ_B)，f(γ_DC)分别为时间基线系数、垂直基线系数和多普勒质心频率差系数，找到最佳主影像后形成干涉对；

外部DEM转到雷达坐标系下形成地形相位，干涉相位减去地形相位以及平地相位，形成差分干涉对。

6.根据权利要求1所述的基于时序InSAR技术的地表形变反演方法，其特征在于，所述步骤S3选取PS点包括：

利用振幅离差系数进行粗选，保留高于阈值的点目标；

通过差分干涉后的差分相位、滤波相位以及地形残余相位来计算相位点在时间上的相干系数，通过选定合适的相干系数阈值，高于该阈值的点则被识别为PS目标点，该步骤需要通过滤波以及迭代运算实现对PS目标点的最终选取。

7.根据权利要求1所述的基于时序InSAR技术的地表形变反演方法，其特征在于，所述步骤S4中时空解缠包括：

在去除地形相位后，干涉对中单点相位构成表示为：

W{φ}＝W{φ_d+φ_τ+φ_a+φ_o+φ_n}

式中，W{·}为相位缠绕，φ为单点干涉相位，φ_d、φ_τ、φ_a、φ_o、φ_n分别表示形变、地形误差、大气、轨道、和噪声相位；其中φ_d、φ_a和φ_o在空间具有相关性，可利用空间低通滤波将其分离，估算出空间相关相位

后，将其去除得到：

式中，上标u代表空间非相关，δ代表未去除的空间不相关的形变、大气、轨道误差的总量；因为这三者具有很强的空间相关性，所以δ值很小；

构造相干系数：

由于δ很小，相干系数γ可以反映单点受噪声的影响大小，噪声均值越大，相干系数值越大，利用空间搜索法可以求算出γ的最大值，并求出对应的

在去除

误差后，进行时空三维解缠。

8.根据权利要求1所述的基于时序InSAR技术的地表形变反演方法，其特征在于，所述步骤S5中误差分离包括：

解缠后的相位值为：

式中，φ_uw为解缠后的相位，

为残余具有空间相关性的地形残差相位，k为单差模糊度整数；如果三维解缠正确，则大部分点的k应该相等；为分离出形变相位，需要在弧段上进行时空域滤波；

式中，

代表两点间做差，然后构建三角网进行时间域低通滤波，将其分离后得到

式中，[·]_{LP_time}表示时间域低通滤波，

和

分别表示副影像大气和轨道分量；

在时间域上表现为较强的相关性，利用时间域低通可以获取形变和主影像大气；其余量在时间域上表现为不相关，要获取副影像大气，

在空间表现为相关性，

表现为随机噪声误差，因此可以通过空间的高通滤波分离出副影像的大气延迟相位和轨道误差；最后通过从弧段到点上的积分可以获取每个点的φ_d及

9.根据权利要求1所述的基于时序InSAR技术的地表形变反演方法，其特征在于，若误差分离中高低通滤波后差分相位仍然出现较为明显的长波趋势，即分层大气延迟干扰，则通过采用线性函数进行吸收并去除。

10.根据权利要求1所述的基于时序InSAR技术的地表形变反演方法，其特征在于，该地表形变反演方法还包括GNSS-InSAR融合，在获得解缠相位之后，若有GNSS外部数据，可以优选通过GNSS数据得到ZWD，反演大气延迟，将其从解缠相位中剔除；

将DEM误差相位加回，并重新进行参数化，利用最小二乘进行形变相位解算，将GNSS观测方程作为约束条件，和InSAR观测方程进行联合平差，以获取最终的形变信息。

Images