一种地表形变监测方法及系统与流程

1.本技术属于地质灾害监测技术领域，特别涉及一种地表形变监测方法及系统。


背景技术：

2.为满足经济和社会的高速发展，许多城市进行了大规模的市政工程建设，很多城市由于沿海区填海工程、抽取地下水、铁路、地铁建设等原因，地表形变已成为主要地质灾害之一，严重危及人民的生命安全。目前，对城市地表形变进行监测的技术主要分为实地测量(水准仪、全站仪等)和遥感监测手段。实地测量由于人工成本高、布设风险大、监测点稀疏、覆盖范围小等原因，不适用于大范围、高频次的动态监测。遥感监测手段中，光学传感器受云的影响，在华南沿海城市多云多雨的气候条件下，难以获取长序列的有效光学影像进行监测。
3.星载合成孔径雷达(synthetic aperture radar,sar)具有全天时、全天候的工作能力、覆盖范围广、穿透能力强(能够透过云层获取有效数据)、卫星重访周期短等特点，且很多sar数据可供全球免费下载，因此可实现大范围、高频次动态监测，较传统的地表形变监测具有不可替代的优势。合成空间雷达干涉测量(insar)技术，已经成为地表形变测量的最有效和可行的技术。其中，时间序列合成孔径雷达干涉测量技术在地表形变监测方面取得了丰硕的成果。该类技术大体上可以分为两个类别：永久散射体合成孔径雷达干涉(persistent scatterer insar,ps-insar)技术和小基线合成孔径雷达干涉(small baseline subset insar,sbas-insar)技术。英国学者hooper于2008年提出了一种新的时间序列insar算法——stamps(stanford method for persistent scatterers)，根据应用需求，单独采用或有效结合ps-insar和sbas-insar两种优化模型，是目前使用最广泛的地表形变监测技术之一。在城市建筑密集区域，通常选择使用stamps中的ps-insar技术进行地表形变监测。现有的stamps ps-insar技术存在的问题主要包括：
4.1)为保持sar影像的原始分辨率从而保持较高的ps(persistent scatter,永久散射体)候选点数量，没有对干涉相位进行多视处理，没有在ps-insar处理之前对空间非相关噪声进行抑制。
5.2)ps点二次筛选时，参数敏感性较高，参数在临界点的微调容易导致ps点过少或过多，ps点过少则无法完整呈现形变的空间变化趋势，ps点过多则往往错误保留了水体和植被表面的噪点。
6.3)在大气相位改正环节，stamps现有的方法或单独采用时空滤波方法或单独采用外部大气产品对大气延迟相位进行估计，并不能处理好华南沿海水汽时空变化大的情况。


技术实现要素：

7.本技术提供了一种地表形变监测方法及系统，旨在至少在一定程度上解决现有技术中的上述技术问题之一。
8.为了解决上述问题，本技术提供了如下技术方案：
9.一种地表形变监测方法，包括：
10.基于stamps算法，采用gamma置信区间判别法对监测区域的sar序列影像数据进行同质点提取，生成所述sar序列影像数据的同质点集；
11.基于所述同质点集计算sar序列影像数据中副影像与主影像之间的干涉相位相干系数；
12.按照所述干涉相干系数计算结果对所述sar序列影像数据的干涉图进行自适应多视处理，获取多视干涉相位；
13.计算所述sar序列影像的振幅离散度，根据所述振幅离散度进行ps候选点筛选，得到ps候选点；
14.基于所述同质点集计算干涉相干系数的自适应阈值，采用所述自适应阈值对ps候选点进行筛选，得到有效ps点掩膜；
15.基于所述有效ps点掩膜以及所述多视干涉相位，采用gacos模型和insar时空滤波算法对所述sar序列影像数据进行大气改正；
16.基于stamps算法对所述大气改正后的sar序列影像数据的形变速率及累积形变量进行时间序列反演，得到所述监测区域的地表形变监测结果。
17.本技术实施例采取的技术方案还包括：所述采用gamma置信区间判别法对监测区域的sar序列影像数据进行同质点提取包括：
18.对所述sar序列影像数据中的每一个像素，分别确定一个搜索窗口，以该像素为参考像素，以该像素的时间维采样计算参考样本均值以和设定的α值计算gamma分布置信区间；所述gamma分布置信区间计算公式为：
[0019][0020]
其中1-α为置信水平，表示标准gammag(nl,1)分布的分位点，n为sar序列影像数据的影像数目，l为sar序列影像数据的视数，a＝l；σ2为sar序列影像数据的后向散射强度；
[0021]
判断所述搜索窗口内的邻域像素的时间维样本均值是否落入gamma分布置信区间，如果是，则判定所述邻域像素与参考像素属于同一个统计分布，得到该像素的初始同质像元子集ω
init
；
[0022]
采用所述初始同质像元子集ω
init
中k个像素的样本均值的平均值对所述参考样本均值进行更新，并根据更新后的重新计算gamma分布置信区间；
[0023]
根据新的gamma分布置信区间对所述初始同质像元子集ω
init
中的像素进行二次判断，保留落入新的gamma分布置信区间并与所述参考像素空间相连的邻域像素，生成该像素的精选同质像元子集ω。
[0024]
本技术实施例采取的技术方案还包括：所述基于所述同质点集计算sar序列影像数据中副影像与主影像之间的干涉相干系数包括：
[0025]
对所述sar序列影像数据中的每个像素，取其同质像元子集ω进行干涉相干系数估计，得到像素的干涉相干系数有偏估计；所述干涉相干系数计算公式为：
[0026][0027]
其中z1(l)，z2(l)代表用于干涉的两幅sar影像z1，z2中像素l对应的复数信号，k为该像素所在的同质像元子集ω内部的像素个数；
[0028]
采用双重bootstrapping方法对所述干涉相干系数有偏估计进行偏差纠正，获取所述sar序列影像数据无偏干涉相干系数；所述近似无偏估计的干涉相干系数为：
[0029][0030]
其中，将同质像元子集ω中所有的点表示为x＝(x1,x2,
…
,xk)，k为同质点个数，r为bootstrapping对x反复抽样的次数，为bootstrapping复制的样本，γ为干涉相干系数有偏估计。
[0031]
本技术实施例采取的技术方案还包括：所述按照干涉相干系数计算结果对sar序列影像数据的干涉图进行自适应多视处理包括：
[0032]
对所述sar序列影像数据进行预处理，生成sar序列影像数据的干涉图；所述预处理包括：
[0033]
对所述sar序列影像数据进行时空基线分析，选取所述sar序列影像数据中空间、时间基线最短的一期影像为主影像，其他影像为副影像；
[0034]
按照ps-insar算法的配对准则，利用轨道参数以及dem模型将所述副影像与主影像进行配准，将副影像与主影像组成干涉对，生成干涉图。
[0035]
本技术实施例采取的技术方案还包括：所述按照干涉相干系数计算结果对sar序列影像数据的干涉图进行自适应多视处理包括：
[0036]
对于所述sar序列影像数据中的每个像素，以所述像素所在的同质像元子集为单位，在所述同质像元子集内部对干涉图复信号按照干涉相干系数进行加权，取其相位，获取多视干涉相位；
[0037]
对于任一像素p，所述多视干涉相位的计算公式为：
[0038][0039]
其中k为像素p所在的同质像元子集中的同质像素个数，cohi为第i个像素的干涉相干系数，ifgi为第i个像素的干涉图复信号。
[0040]
本技术实施例采取的技术方案还包括：所述计算所述sar序列影像数据的振幅离散度具体为：
[0041][0042]
其中，da代表振幅离散度，σa和μa分别为sar序列影像数据沿时间维计算的标准差和均值。
[0043]
本技术实施例采取的技术方案还包括：所述基于所述同质点集计算干涉相干系数的自适应阈值，采用所述自适应阈值对ps候选点进行筛选包括：
[0044]
根据设置的真实阈值以及同质像元子集中的同质像元个数，计算出自适应阈值：
[0045][0046]
其中σ
co
h为无偏干涉相干系数的cramer-rao下界标准差；
[0047]
对于所述sar序列影像数据中的每一个像素，根据该像素的自适应阈值，对所述无偏干涉相干系数在时间维形成的序列进行统计，如果时序中干涉相干系数的精估值的个数大于设定比例，则确定该点为有效ps点；
[0048]
将所有有效ps点与所述ps候选点求交集，得到有效ps点掩膜。
[0049]
本技术实施例采取的技术方案还包括：所述基于所述有效ps点掩膜以及多视干涉相位，采用gacos模型和insar时空滤波算法对所述sar序列影像数据进行大气改正包括：
[0050]
采用stamps算法对所述sar序列影像数据进行空间相关相位分量估计，对所述有效ps点掩膜进行相位噪声估计；
[0051]
结合所述振幅离散度和相位噪声估计对所述有效ps点掩膜进行筛选，并对筛选后的ps点进行三维相位解缠：
[0052]
基于所述ps点的相位解缠结果，采用gacos模型和insar时空滤波算法从干涉相位中滤除大气延迟、轨道误差相位，对所述sar序列影像数据进行大气改正。
[0053]
本技术实施例采取的技术方案还包括：所述采用gacos模型和insar时空滤波算法对所述sar序列影像数据进行大气改正包括：
[0054]
采用gacos模型对每一组干涉对的大气延迟相位进行估计，并从对应的干涉相位中去除与地形相关的大气相位分量；
[0055]
采用insar时空滤波算法对所述sar序列影像数据进行时间维低通滤波，得到所述主影像的大气延迟以及轨道误差；
[0056]
采用insar时空滤波算法对所述sar数据进行时间维高通滤波嵌套空间维低通滤波，得到所述副影像的大气延迟、轨道误差以及dem误差；
[0057]
从所述干涉相位中去除主影像的大气延迟以及轨道误差、副影像的大气延迟、轨道误差以及dem误差，对所述sar序列影像数据进行精细化大气改正。
[0058]
本技术实施例采取的另一技术方案为：一种地表形变监测系统，包括：
[0059]
同质点提取模块：用于基于stamps算法，采用gamma置信区间判别法对监测区域的sar序列影像数据进行同质点提取，生成sar序列影像数据的同质点集；
[0060]
多视处理模块：用于基于所述同质点集计算sar序列影像数据中副影像与主影像之间的干涉相干系数，并按照干涉相干系数计算结果对sar序列影像数据的干涉图进行自适应多视处理，获取多视干涉相位；
[0061]
初始ps点筛选模块：用于计算所述sar序列影像的振幅离散度，根据所述振幅离散度进行ps候选点筛选，得到ps候选点；
[0062]
有效ps点筛选模块：用于基于所述同质点集计算干涉相干系数的自适应阈值，采用所述自适应阈值对ps候选点进行筛选，得到有效ps点掩膜；
[0063]
大气改正模块：用于基于所述有效ps点掩膜以及多视干涉相位，采用gacos模型和i nsar时空滤波算法对所述sar序列影像数据进行大气改正；
[0064]
序列反演模块：用于基于stamps算法对所述大气改正后的sar序列影像数据的形变速率及累积形变量进行时间序列反演，得到所述监测区域的地表形变监测结果。
[0065]
相对于现有技术，本技术实施例产生的有益效果在于：本技术实施例的地表形变监测方法及系统基于同质点提取对经典stamps算法进行了改进，在ps点筛选之前，对干涉相位进行基于同质点集的自适应多视处理，通过区分永久散射点ps和分布式散射点ds两种不同散射机制的目标，既有效抑制了相位噪声，又避免了ps点目标的高相干信号被其他点的相位噪声混淆，提高了后续解算精度。采用偏差纠正后的干涉相干系数和自适应阈值辅助ps点筛选，在有效去除水体和植被表面的低相干点的同时，尽量保持了空间形变趋势的完整性，避免了对ps点的过度筛除。通过结合gacos大气模型和insar时空滤波对大气延迟相位进行精细改正，提高了大气改正的精度，从而提高地表形变的监测精度。
附图说明
[0066]
图1是本技术实施例的地表形变监测方法的流程图；
[0067]
图2为本技术实施例结合gacos大气模型和时空滤波的大气精细改正流程图；
[0068]
图3为本技术实施例的地表形变监测系统结构示意图。
具体实施方式
[0069]
为了使本技术的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本技术进行进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本技术，并不用于限定本技术。
[0070]
请参阅图1，是本技术实施例的地表形变监测方法的流程图。本技术实施例的地表形变监测方法包括以下步骤：
[0071]
s10：获取监测区域的sar序列影像数据；
[0072]
本步骤中，获取的sar序列影像数据为sentinel-1iw tops模式slc时间序列数据，slc时间序列数据参数包括时间范围、轨道方向、重访周期、极化方式、入射角、波长、斜距采样间隔以及方位向采样间隔等，具体如表1所示：
[0073]
表1.监测区域的sentinel-1iw slc时间序列数据参数
[0074][0075][0076]
s20：对sar序列影像数据进行预处理，生成sar序列影像数据的干涉图；
[0077]
本步骤中，sar序列影像数据的预处理具体包括：
[0078]
s21：对sar序列影像数据进行时空基线分析，选取sar序列影像数据中空间、时间基线最短的一期影像为主影像，其他影像为副影像；
[0079]
s22：按照ps-insar算法的配对准则，利用轨道参数以及dem(digital elevation model，数字高程模型)模型将副影像与主影像进行配准，将副影像与主影像组成干涉对，生成干涉图；
[0080]
s23：基于dem模型去除干涉图的平地相位和地形相位后，导出配准后的干涉相位文件。
[0081]
s30：在stamps算法中，采用gamma置信区间判别法对sar序列影像数据进行同质点(statistically homogeneous pixel,shp)提取，生成sar序列影像数据的同质点集；
[0082]
本步骤中，sar影像中一个像素(即分辨单元)的后向散射强度，是同属于该像素的大量散射体的后向散射强度的叠加。如果一个像素中全都是分布式散射体，则该像素为ds(distributed scatterer,分布散射点)点目标，可以看作是零均值复高斯随机变量。给定一个具有n期slc(single look complex，单视复数据)格式的sar序列影像数据，对于一个ds点目标来说，其时间维采样{z1,z2,...,zn}
t
可以看作是零均值复高斯随机变量,其协方差矩阵可定义为：
[0083][0084]
其中为第i幅影像的振幅，矩阵的对角线元素是第i幅sar影像的后向散射强度。在多视情况下，假设视数为l，则后向散射强度表示式为服从gamma分布，分布参数为和l：
[0085][0086]
对于一个ds点目标，其后向散射强度的时间维矢量表示为对于一个ds点目标，其后向散射强度的时间维矢量表示为本技术的目的在于找到关于σ的一个置信区间。根据sar统计理论，σ服从gamma分布σ～g(a,b)，其中a＝l,其均值为μ＝σ2，方差为根据中心极限定理，样本均值服从正态分布则对应的的置信区间为：
[0087][0088]
其中是标准整体分布概率密度函数的分位点，1-α为置信水平(例如置信水平为95％，则α＝5％)。对该像素的邻域像素逐个进行判别，如果邻域像素的样本均值落入该置信区间，则认为邻域像素与该像素属于同一个统计分布，将邻域像素与该像素归入一个
同质点子集。
[0089]
当sar影像数据的时序较短，即影像的期数n小于一定值时，gamma分布的样本均值并不服从正态分布此时置信区间是不对称的，则无法根据中心极限定理进行推断。在此情况下，如果一个随机变量服从gamma分布σ～g(a,b)，其中a＝l,则令y服从标准gamma分布y～g(na,1)，则对应的gamma分布置信区间为：
[0090][0091]
其中表示标准gamma分布的分位点。因为b＝σ2/l，则a＝l，代入(4)，经过转换，可以得到的1-α置信区间为：
[0092][0093]
其中1-α为置信水平(例如置信水平为95％，则α＝5％)，表示标准gamma分布g(nl,1)的分位点，n为sar序列影像数据的影像数目，l为sar序列影像数据的视数，σ2为sar序列影像数据的后向散射强度,是σ的样本均值。因此，无论sar序列影像数据的时序长短，可统一采用式(7)作为置信区间进行同质点判别。
[0094]
基于上述的基于gamma分布置信区间的同质点判别原理，本技术实施例中，采用gamma置信区间判别法对sar序列影像数据进行同质点提取的过程具体包括：
[0095]
s31：采用调整型箱线图(adjusted boxplot)在考虑数据分布偏斜度(shrewdness)指标medcouple(mc)的前提下，对sar序列影像数据中的奇异值进行检测和剔除；
[0096]
s32：采用单边anderson
–
darling(ad)检验法判断去除奇异值之后的sar序列影像数据中每个像素的时间维矢量是否服从gamma分布，如果服从gamma分布，执行s33；
[0097]
s33：采用式(7)作为置信区间对sar序列影像数据进行同质点判别，提取出sar序列影像数据的同质点集；
[0098]
其中，具体提取过程为：对sar序列影像数据中的每一个像素，分别确定一个搜索窗口，以该像素为参考像素，以该像素的时间维采样计算参考样本均值以和设定的α值计算gamma分布置信区间；对搜索窗口内的所有邻域像素逐一进行判断，如果邻域像素的时间维样本均值落入gamma分布置信区间，则判定该邻域像素与参考像素属于同一个统计分布，得到该像素的初始同质像元子集ω
init
。由于时间序列的长度n是有限的，样本均值的方差较大。因此，本技术实施例采用初始同质像元子集ω
init
中k个像素的样本均值的平均值对参考样本均值进行更新，然后根据更新后的重新计算gamma分布置
信区间，根据新的gamma分布置信区间对初始同质像元子集ω
init
中的像素进行二次判断，保留落入新的gamma分布置信区间并且与参考像素空间相连的邻域像素，从而生成该像素的精选同质像元子集ω。对sar序列影像数据中的每一个像素分别进行上述处理，生成所有像素的同质点集。本技术实施例中，生成初始同质像元子集ω
init
时，设定α＝50％；生成精选同质像元子集ω时，设定α＝5％，具体α值可根据实际应用进行设定。
[0099]
s40：基于同质点集对sar序列影像数据中副影像与主影像之间的干涉相位相干性进行估计，并采用双重bootstrapping(有放回抽样)方法对干涉相位相干性进行偏差纠正，获取sar序列影像数据无偏干涉相干系数；
[0100]
本步骤中，干涉相位相干性(coherence)是两幅sar影像信号相似性的度量，用于衡量insar相位估计精度以及指导相位滤波和相位解缠。式(3)中的协方差矩阵中，非对角线元素即为相干性，视数为l时，使用滑动窗口的相干性计算公式为：
[0101][0102]
其中幅度值为干涉相干系数，zi(l)为像素l的复信号。使用滑动窗口进行相干性计算时，其中的隐含假设是一个窗口内部的相邻像素都属于同一个统计分布。然而，在纹理比较丰富的场景下，一个窗口内部的像素很可能具有不同的统计分布。对具有不同散射特性的邻域像素做无差别平均，会导致对干涉相位相干性的估计存在偏差，并导致图像分辨率下降。因此，本发明实施例采用同质点集代替传统滑窗对干涉相干系数进行估计，并采用双重bootstrapping(有放回抽样)方法对偏差进行进一步纠正，提高干涉相干系数的估计精度。
[0103]
具体地，基于同质点集对sar序列影像数据进行干涉相位相干性估计和偏差纠正包括：
[0104]
s41：对每个像素，取其同质像元子集ω进行干涉相干系数估计，得到像素的干涉相干系数有偏估计；
[0105]
其中，将式(8)中的滑动窗口替换为同质点集，则干涉相干系数计算公式为：
[0106][0107]
其中z1(l)，z2(l)代表用于干涉的两幅sar影像z1，z2中像素l对应的复数信号，k为该像素所在的同质像元子集ω内部的像素个数。本步骤可以避免sar影像各处质地差异引起的估计偏差。
[0108]
s42：采用非参数化的双重bootstrapping方法对干涉相干系数有偏估计进行偏差纠正，得到无偏干涉相干系数；
[0109]
其中，对于参与干涉的两幅sar影像，将同质点集ω中的每个像素所对应的复数信号可以看作是一对复数观测值xi＝(z
1i
,z
2i
)；则同质点集ω中所有的点可以表示为x＝(x1,x2,
…
,xk)，k为同质点个数，可以看作是一个bootstraping样本即从原始观测x中有放回的抽样k次，此为第一重抽样。对x反复抽样r次后，可以得到一组双重
bootstrapping样本对每个样本采用式(9)计算干涉相干系数计算之后，可以得到r个bootstrapping复制的样本采用这些bootstrapping复制的样本对干涉相干系数有偏估计进行偏差纠正，获得近似无偏估计的干涉相干系数
[0110][0111]
s50：采用同质点集，按照干涉相干系数计算结果对sar序列影像数据的干涉图进行自适应多视处理，获取多视干涉相位；
[0112]
本步骤中，由于研究区域中都存在着植被、水体、人工建筑物等不同散射特性的表面，在进行多视处理时，为了不混淆不同特性的散射体，避免相干性低的散射体的相位噪声对相干性高的散射体造成干扰，本技术实施例采用同质点集代替规则滑窗，对时间序列中的所有干涉图进行自适应多视处理，通过区分永久散射点ps和分布式散射点ds两种不同散射机制的目标，既有效抑制了相位噪声，又避免了ps点目标的高相干信号被其他点的相位噪声混淆，提高了后续解算精度。
[0113]
具体地，对于每个像素，以其所在的同质像元子集为单位，在同质像元子集内部对干涉图复信号按照精确估计的干涉相干系数进行加权，取其相位，并对sar序列影像数据中的每一个像素重复此操作，获取自适应多视后的干涉相位。对于任一像素p，多视干涉相位的计算公式如下：
[0114][0115]
其中k为该像素所在的同质像元子集中的同质像素个数，cohi为第i个像素的干涉相干系数，ifgi为第i个像素的干涉图复信号。
[0116]
基于上述，对于每个像素，以其所在的同质像元子集对干涉相位进行加权，权重为该像素的相干系数，ps点位于独属于自己的同质像元子集内部而没有受到其他点目标相位噪声的干扰；而ds点在所属的同质像元子集内部进行相干系数加权的相位平均，显著降低了相位噪声。与现有技术相比，本技术实施例在抑制相位噪声的同时，保留了局部细节，避免减少ps点的数量。
[0117]
s60：计算sar序列影像数据的振幅离散度(amplitude dispersion)，并根据振幅离散度对sar序列影像数据进行ps候选点筛选，得到ps候选点；
[0118]
本步骤中，振幅离散度计算公式如下：
[0119][0120]
其中，da代表振幅离散度，σa和μa分别为sar序列影像数据沿时间维计算的标准差和均值。本技术实施例中，da的阈值设定为0.35～0.4。
[0121]
s70：基于同质点集计算干涉相干系数的自适应阈值，采用自适应阈值对ps候选点进行筛选，得到有效ps点掩膜；
[0122]
本步骤中，首先根据设置的真实阈值以及同质点集中的同质像元个数，计
算出一个自适应阈值：
[0123][0124]
其中σ
coh
为无偏干涉相干系数的cramer-rao下界标准差，计算公式为：
[0125][0126]
对于sar序列影像数据中的每一个像素，根据式(13)计算出该像素的自适应阈值，对基于同质点集计算并经过偏差纠正的干涉相干系数在时间维形成的序列进行统计，如果时序中干涉相干系数的精估值的个数大于设定比例(可根据情况自行调整，例如85％)，则确定该点为有效ps点。对所有像素进行逐个分析之后，与ps候选点求交集，得到有效ps点掩膜。
[0127]
基于上述，本技术实施例采用偏差纠正后的干涉相干系数和自适应阈值辅助ps点筛选，在有效去除水体和植被表面的低相干点的同时，尽量保持了空间形变趋势的完整性，避免了对ps点的过度筛除。
[0128]
s80：采用stamps算法对sar序列影像数据进行空间相关相位分量估计，对有效ps点掩膜进行相位噪声估计，并获取地形残余相位估计；
[0129]
s90：结合振幅离散度和相位噪声估计对有效ps点掩膜进行筛选，并对筛选后的ps点进行三维(空间二维 时间维)相位解缠；
[0130]
s100：基于ps点的相位解缠结果，采用gacos模型和insar时空滤波算法从干涉相位中滤除大气延迟、轨道误差相位等，对sar序列影像数据进行精细化大气改正；
[0131]
本步骤中，首先采用gacos模型对每一组干涉对的大气延迟相位进行估计，去除与地形相关的大气相位分量，采用stamps中的时空滤波方法对残余大气相位进行滤除，提高了大气改正的精度，从而提高地表形变的监测精度。具体如图2所示，为本技术实施例的结合gacos大气模型和insar时空滤波算法的精细化大气改正流程图，具体包括以下步骤：
[0132]
s101：首先基于gacos大气模型，对每景sar数据获取日期的大气延迟相位进行估计，并从对应的干涉相位中去除与地形相关的大气相位分量，完成sar序列影像数据的初步大气改正；
[0133]
s102：采用insar时空滤波算法对sar数据进行时间维低通滤波，得到主影像的大气延迟以及轨道误差；
[0134]
s103：采用insar时空滤波算法对sar数据进行时间维高通滤波嵌套空间维低通滤波，得到其他副影像的大气延迟、轨道误差以及dem误差；
[0135]
s104：从干涉相位中去除主影像的大气延迟以及轨道误差、副影像的大气延迟、轨道误差以及dem误差，完成sar序列影像数据的精细化大气改正。
[0136]
s110：基于stamps算法对精细化大气改正后的sar序列影像数据的年形变速率及累积形变量进行时间序列反演，得到监测区域的地表形变监测结果。
[0137]
基于上述，本技术实施例的地表形变监测方法基于同质点提取提出了shp stamps算法，对经典stamps算法进行了改进，在ps点筛选之前，对干涉相位进行基于同质点集的自适应多视处理，通过区分永久散射点ps和分布式散射点ds两种不同散射机制的目标，既有
效抑制了相位噪声，又避免了ps点目标的高相干信号被其他点的相位噪声混淆，提高了后续解算精度。采用偏差纠正后的干涉相干系数和自适应阈值辅助ps点筛选，在有效去除水体和植被表面的低相干点的同时，尽量保持了空间形变趋势的完整性，避免了对ps点的过度筛除。通过结合gacos大气模型和insar时空滤波对大气延迟相位进行精细改正，提高了大气改正的精度，从而提高地表形变的监测精度。
[0138]
请参阅图3，为本技术实施例的地表形变监测系统结构示意图。本技术实施例的地表形变监测系统40包括：
[0139]
同质点提取模块41：用于基于stamps算法，采用gamma置信区间判别法对监测区域的sar序列影像数据进行同质点提取，生成sar序列影像数据的同质点集；
[0140]
多视处理模块42：用于基于同质点集计算sar序列影像数据的干涉相干系数，并按照干涉相干系数计算结果对sar序列影像数据的干涉图进行自适应多视处理，获取多视干涉相位；
[0141]
初始ps点筛选模块43：用于计算sar序列影像数据的振幅离散度，根据振幅离散度对sar序列影像数据进行ps候选点筛选，得到ps候选点；
[0142]
有效ps点筛选模块44：用于基于所述同质点集计算干涉相干系数的自适应阈值，采用自适应阈值对ps候选点进行筛选，得到有效ps点掩膜；
[0143]
大气改正模块45：用于基于所述有效ps点掩膜以及所述多视干涉相位，采用gacos模型和insar时空滤波算法对所述sar序列影像数据进行大气改正；
[0144]
序列反演模块46：用于基于stamps算法对所述大气改正后的sar序列影像数据的形变速率及累积形变量进行时间序列反演，得到所述监测区域的地表形变监测结果。
[0145]
对所公开的实施例的上述说明，使本领域专业技术人员能够实现或使用本技术。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的，本技术中所定义的一般原理可以在不脱离本技术的精神或范围的情况下，在其它实施例中实现。因此，本技术将不会被限制于本技术所示的这些实施例，而是要符合与本技术所公开的原理和新颖特点相一致的最宽的范围。