CN110146885B - 基于多卫星平台MT-InSAR融合的新成陆区超长形变时序提取方法 - Google Patents

Abstract

本发明提供了一种基于多卫星平台MT‑InSAR融合的新成陆区超长形变时序提取方法，包括以下步骤：第一步：各个单一卫星平台MT‑InSAR形变时间序列反演；第二步：提取各个单一卫星平台MT‑InSAR形变时间序列的公共时间跨度和公共空间覆盖范围；第三步：选取各个卫星平台MT‑InSAR公共高相干点；第四步：新成陆区有时间重叠的两个卫星平台MT‑InSAR形变时序融合；第五步：新成陆区有时间间隙的各个卫星平台MT‑InSAR形变时序融合方法。本发明能够有效解决新成陆区持续时间跨度长达十几年甚至几十年的大面积不均匀缓慢地面沉降的完整连续监测问题，以及预测新成陆区累积固结沉降量和固结沉降结束时间问题，突破了基于单一卫星平台监测技术的局限，极大提高了高精度地表形变连续量测的时间长度。

Description

基于多卫星平台MT-InSAR融合的新成陆区超长形变时序提取 方法

技术领域

本发明涉及多时相遥感和大地测量技术领域，具体涉及一种针对新成陆区基于多个卫星平台MT-InSAR融合以生成超长(≥10年)形变时序的方法。本发明能够有效解决新成陆区持续时间跨度长达十几年甚至几十年的大面积不均匀缓慢地面沉降的完整连续监测问题，以及预测新成陆区累积固结沉降量和固结沉降结束时间问题，突破了基于单一卫星平台监测技术的局限，极大提高了高精度地表形变连续量测的时间长度。

背景技术

多时域合成孔径雷达干涉测量技术(Multi-Temporal Synthetic ApertureRadar Interferometry，MT-InSAR)是新近发展的一种能够大面积、精确监测地表形变的技术方法。相较于传统的差分干涉测量技术，MT-InSAR技术提高了地表形变的量测精度，较好地抑制了干涉图中时空失相关和大气噪声的影响。然而,MT-InSAR技术通常是基于单一平台多景雷达影像进行联合分析，对于持续十年甚至数十年之久的缓变型沉降等现象，单一平台MT-InSAR难以量测其长期连贯形变过程和掌握相关的动态变化规律。例如，河口海岸疏浚的淤泥和泥沙通过水力吹填在岸边堆积形成的新成陆区,新成陆区内新近形成的冲填土压缩性高，易发生大面积不均匀地面沉降，且沉降持续时间跨度较长(十几年至几十年)，形变速率缓慢，准确掌握新成陆区大面积、长时期冲填土固结沉降时空演变过程和特征，是新成陆区急需解决的突出环境地质安全问题。

随着越来越多的合成孔径雷达(Synthetic Aperture Radar，SAR)卫星发射升空，如欧空局的Sentinel-1,日本宇航局的ALOS-2，意大利空间局的COSMOS-SkyMed,以及德国空间局的TerraSAR-X等等，为同一地区形变监测提供了多平台丰富的SAR数据。然而不同卫星提供的SAR数据的成像几何、时空分辨率以及雷达波长等都不尽相同，不能够作为一个干涉数据集进行解算。对于同一区域，不同卫星数据各自得到的MT-InSAR量测结果则在时空尺度上存在不一致性。如何实现多平台MT-InSAR形变时间序列数据的融合，生成连贯一致的多平台MT-InSAR形变超长时间序列数据集，是当前亟待解决的难题之一。

发明内容

本发明的目的就是提供一种基于多卫星平台MT-InSAR融合的新成陆区超长形变时序提取方法。它能够有效解决新成陆区持续时间跨度长达十几年甚至几十年的大面积不均匀缓慢地面沉降的完整连续监测问题，以及预测新成陆区累积固结沉降量和固结沉降结束时间问题，突破了基于单一卫星平台监测技术的局限，极大提高了高精度地表形变连续量测的时间长度。

实现本发明目的的具体技术方案是：

一种基于多卫星平台MT-InSAR融合的新成陆区超长形变时序提取方法，能够针对新成陆区融合有时间重叠和有时间间隙的各个卫星平台MT-InSAR形变时序，生成时间跨度长达十几年甚至几十年的形变时间序列，该方法包括以下步骤：

步骤1：各个单一卫星平台MT-InSAR形变时间序列反演

利用高级时序合成孔径雷达干涉测量方法即MT-InSAR和各个卫星平台获取的新成陆区合成孔径雷达即SAR影像时间序列，分别反演新成陆区地表形变场和形变时间序列；

步骤2：提取各个单一卫星平台MT-InSAR形变时间序列的公共时间跨度和公共空间覆盖范围

在时间上，从步骤1得到的各个卫星平台MT-InSAR形变时间序列中提取出重叠时间范围内的形变时间序列子集；在空间上，提取各个卫星平台MT-InSAR地表形变场的公共空间覆盖范围，以及在此公共范围内的各个卫星平台MT-InSAR选取的高相干点，即时间相干系数阈值大于给定阈值；

步骤3：选取各个卫星平台MT-InSAR公共高相干点

各个卫星平台MT-InSAR选取的公共高相干点是基于像元的点目标，将各个单一平台MT-InSAR提取的高相干点进行空间匹配，将不同分辨率MT-InSAR形变场重采样到统一空间分辨率，统计和保存各平台MT-InSAR公共高相干点的像元位置并保留像元位置数据；

步骤4：新成陆区有时间重叠的两个卫星平台MT-InSAR形变时序融合

基于重叠时间段内两个卫星平台的形变时间序列子集的累积分布函数差异建立其偏差的统计分布特征，估算和去除后一个卫星平台MT-InSAR形变时序相对于前一个卫星平台MT-InSAR形变时序的偏差，在每一个公共高相干点上融合生成时间上连贯一致的长时间MT-InSAR形变序列数据集；详细步骤包括：

1)计算重叠时间段内两个卫星平台MT-InSAR形变时间序列子集的累积分布函数

设有N个卫星平台合成孔径雷达数据集，通过MT-InSAR方法，基于每个合成孔径雷达数据集反演的形变时间序列表示为：

d^(x)(P)为第x个合成孔径雷达卫星平台数据集反演得到的P点的形变时间序列；

第x个合成孔径雷达卫星平台数据集反演得到的P点的形变时间序列中第1,…Q_x个形变值；其中Q_x是第x个卫星平台合成孔径雷达影像的数量，x为合成孔径雷达卫星平台的序号，N是卫星平台的总数；S是步骤3中选取的各个卫星平台MT-InSAR公共高相干点的集合，P是该集合中任意一个高相干点。

若要融合两个在时间上有重叠的MT-InSAR形变时序，则x＝1,2,其中1代表在时间上先获取的合成孔径雷达数据集，2代表时间上后获取的合成孔径雷达数据集；第x个形变数据d^(x)是相对于时间参考点

的形变量，时间参考点

的形变量为

要融合时间重叠的两个卫星平台MT-InSAR形变时序，需要确定第二个平台的MT-InSAR形变时序中每一个时间点的形变量相对于第一个平台MT-InSAR形变时序的时间参考点的形变量；

根据两个卫星平台的形变时间序列子集在重叠时间段内的配对形变数据，计算各自的累积分布函数；

2)基于累积分布函数匹配建立偏差分布特征，估算和去除MT-InSAR形变时序之间的偏差；

是基于第二个合成孔径雷达数据集反演的MT-InSAR形变时间序列中第j个形变值，需要将其校正到基于第一个平台合成孔径雷达数据集反演的MT-InSAR形变水准，去除相对偏差；j＝1,…,K，K是第二个合成孔径雷达数据集中影像的数量；改正后的第二个合成孔径雷达数据集反演的MT-InSAR形变时间序列中第j个形变值表示为

这一过程表示为：

Δdef_j是基于第二个合成孔径雷达数据集反演的MT-InSAR形变时序中第j个形变值相对于基于第一个平台合成孔径雷达数据集反演的MT-InSAR形变时序中第j个形变值的偏差；

要确定两个时序之间的相对偏差，首先，在两个形变时序的重叠时间跨度内，建立两个形变时序之间的相关的控制偏差量；公共时间跨度内的基于第一个平台合成孔径雷达数据集反演的MT-InSAR形变时序和基于第二个合成孔径雷达数据集反演的MT-InSAR形变时序分别定义为C_master和C_slave；C_slave相对于C_master的相对偏差Δdef计算过程如下：

其中

IQR_slave＝C_slave|_p＝75％-C_slave|_p＝25％ (10)

IQR_master＝C_master|_p＝75％-C_master|_p＝25％ (11)

在公式(3)到(11)中，g′和f′是调整校正后的时间序列分布的两个调制参数；

和

是两个控制观测数据即时间重叠时间段内观测值的中值；

是两个控制观测数据即时间重叠时间段内观测值的中值的差值；

公共时间跨度内的基于第二个合成孔径雷达数据集反演的MT-InSAR形变时序中第k个形变值；

公共时间跨度内的基于第一个合成孔径雷达数据集反演的MT-InSAR形变时序中第k个形变值；Δ_kk是公共时间跨度内的第二个合成孔径雷达数据集反演的MT-InSAR形变时序中第k个形变值与基于第一个合成孔径雷达数据集反演的MT-InSAR形变时序中第k个形变值的差值；Δ′_kk是公共时间跨度内两个平台数据集反演的MT-InSAR形变时序中第k个形变值的差值与时间重叠时间段内观测值的中值的差值之差；kk是p中对应的第k个百分位数，是通过找到C_slave的数据值相等或接近给定的观测数据值

来计算；p是C_slave对应的累积分布函数；IQR_master和IQR_slave分别是参考形变时间序列和待校正形变时间序列的形变控制观测值的四分位距，分别为75％和25％。

经过以上估算，确定和去除后一个卫星平台MT-InSAR形变时序相对于前一个卫星平台MT-InSAR形变时序的偏差，将后一个卫星平台MT-InSAR形变时序校正到前一个卫星平台MT-InSAR形变时序的形变水准，融合生成时间上连贯一致的长时间序列MT-InSAR形变数据集；步骤5：新成陆区有时间间隙的各个卫星平台MT-InSAR形变时序融合

融合有时间间隔的各个卫星平台MT-InSAR形变时序，将各个卫星平台的MT-InSAR形变时序与外部沉降经验模型进行非线性拟合；新成陆区冲填土长期缓慢的固结过程会引起大范围不均匀地面沉降，通过土力学实验模拟建立的新成陆区长期缓慢地表形变经验模型(12)：

式中，t为新成陆区冲填土固结时间；S_CMT(t)为t时刻的累积沉降量；S_m是累积形变总量；k和λ为模型参数，决定沉降模型的曲率和曲线的形状；δ为时间延迟，对应于新成陆区地面沉降开始的时间；

在拟合时，首先对公式(12)中的模型改进为，

d为有时间间隙的形变时间序列间的形变差值；m为参数，先获取的形变时间序列对应的m为0，后获取的形变时间序列对应的m值为1；

进而，设置模型参数λ的取值范围为0～2，在进行模型拟合时λ的步长设置为0.1；针对不同的λ值，采用迭代的莱文贝格-马夸特最小二乘最优化算法(Levenberg-Marquardt,LM)对有时间间隙的各个平台MT-InSAR形变时间序列和外部沉降模型进行拟合；对于每一个给定的λ值，拟合得到一个由各个平台MT-InSAR形变时序约束的形变模型，从所有λ值对应的模型中选取最优形变模型，拟合后的各个平台MT-InSAR形变时序即为融合得到的多平台MT-InSAR形变时序；其中最优形变模型选取依据是模型拟合得到的各个平台形变时间序列与根据公式(13)中参数m得到的MT-InSAR形变时间序列之间的均方根误差值，即均方根误差值越小、模型非线性拟合效果越好，所得模型最优；

将各个平台的MT-InSAR形变时序与新成陆区地表形变经验模型进行非线性拟合后，实现新成陆区有时间间隙的多卫星平台MT-InSAR形变时序融合，生成时间上连贯一致的长时间MT-InSAR形变序列数据集。

本发明提出的一种基于多卫星平台MT-InSAR融合的新成陆区超长形变时序提取方法，能够有效解决新成陆区持续时间跨度长达十几年甚至几十年的大面积不均匀缓慢地面沉降的完整连续监测问题，以及预测新成陆区累积固结沉降量和固结沉降结束时间问题。本发明的优点在于：(i)突破了基于单一卫星平台监测技术的局限，极大提高了高精度地表形变连续量测的时间长度，实现了新成陆区连贯一致的超长时间序列(≥10年)形变数据集的生成，可以为新成陆区缓变型沉降长期量测提供新方法和可靠数据；(ⅱ)提出的时间维度上多平台MT-InSAR地表形变时间序列的融合方法，将时间维度上多卫星平台MT-InSAR融合分为两种情形，即有时间重叠的多平台MT-InSAR形变时序融合和存在时间间隙的多平台MT-InSAR形变时序融合，无需同时处理大量干涉纹图，本发明亦适用于持续十年甚至数十年之久的缓变型沉降等其他长时间动态变化地学现象的观测与分析。

附图说明

图1本发明流程图；

图2累积分布函数计算示意图；

图3本发明反演地表形变的流程图；

图4本发明实施例效果图。

具体实施方式

实施例

以上海市在东部沿海滩涂湿地冲填形成的新成陆区长期缓慢形变时序融合为例，说明本发明在实施例中的具体实施步骤。图1是针对该实例的多卫星平台MT-InSAR数据获取超长形变时序融合方法流程图。图2是累积分布函数计算示意图。

步骤1：利用多时相合成孔径雷达干涉测量方法(MT-InSAR)反演各个单一卫星平台MT-InSAR形变时间序列；

长期、连续覆盖上海东部新成陆区的SAR卫星影像有：欧空局ENVISAT ASAR(ENV)数据集，时间跨度为2007年2月26日至2010年9月13日的35景IM模式影像(升轨、VV极化)，分辨率约为30m，影像幅宽约为100km，卫星视角和航向角分别为23°和348°。意大利空间局Cosmo-SkyMed(CSK)卫星数据集，时间跨度为2013年12月7日至2016年3月18日的61景条带模式数据(降轨、HH极化)，分辨率为3m，卫星视角和航向角分别为29°和198°。欧空局Sentinel-1A(S1A)卫星，时间跨度为2015年2月26日至2017年4月4日的33景干涉宽幅模式数据(升轨、VV极化)，分辨率约为20m，卫星视角和航向角分别为39°和348°。

利用多时相合成孔径雷达干涉测量方法(MT-InSAR)，地表形变反演流程和步骤如图3所示，对三个卫星平台的SAR数据集分别反演地表形变时间序列。针对ENV数据集和MT-InSAR方法反演形变时序。为了减小噪声因素的影响，对各个干涉图进行了多视处理(方位向20视、距离向4视)，多视后的影像分辨率约为90m×90m左右。利用90m×90m分辨率的SRTM(Shuttle Radar Topography Mission)高程数据对各个干涉影像对进行了差分干涉处理，以消除地形相位带来的影响。选取的高相干点的时间相干系数大于0.65。针对CSK数据集和MT-InSAR方法反演形变时序，设置方位向视数为30视、距离向视数为30视，对干涉图进行了多视处理，多视后干涉图的分辨率约为90m×90m。在利用90m×90m分辨率的SRTM DEM高程数据去除高程相位，选取的高相干点的时间相干系数大于0.65。针对S1A数据集和MT-InSAR方法反演形变时序，设置了方位向视数为4视、距离向视数为16视对干涉图进行了多视处理，多视后干涉图的分辨率约为80m×80m。同样利用SRTM DEM高程数据去除高程相位，选取的高相干点的时间相干系数大于0.65。

步骤2：提取各个单一卫星平台MT-InSAR形变时间序列的公共时间跨度和公共空间覆盖范围

步骤3：选取多卫星平台MT-InSAR公共高相干点

步骤4：针对有时间重叠的多卫星平台MT-InSAR形变时序，利用基于累积分布函数的时序融合方法融合各个平台MT-InSAR形变时序在每个公共高相干点上的形变时序。以融合CSK和S1A的MT-InSAR形变时序为实例，CSK数据集获取时间跨度为2013-2016，S1A数据集获取时间跨度为2015-2017。两个数据集重叠时间跨度为2015-2016。融合CSK MT-InSAR和S1A MT-InSAR两个形变时间序列数据集，令N＝1,2，将CSK MT-InSAR形变时序记为第一个形变时序

S1A MT-InSAR形变时序记为第二个形变时序集记为

两个形变时序有各自的形变时间参考点，在各自的形变时间参考点上的形变值都设为0，即

要融合两个卫星平台MT-InSAR形变时序，关键是要确定第二个形变时序相对于第一个形变时序时间参考点的形变偏移量，用Δdef_i表示；将第二个形变时序改正为以第一个形变时序时间参考点的形变时序，表示为

则改正后的第二个形变时序和未改正的第二个形变时序之间的关系为：

为了确定将第二个形变时序的每一个时间点的形变量相对于第一个形变时序时间参考点的形变量，在两个数据集的重叠时间跨度2015-2016内，建立两个形变时序集各自的累积分布函数，通过匹配两个形变时序在累积分布函数曲线的不同分位数上的形变数值，作为两个形变时序之间偏差控制量，并确定两者的非线性关系。在每一个公共高相干点上，通过确定和去除第二个形变时序相对于第一个形变时序的相对偏差，将第二个形变时序的每一个形变观测量统一到与第一个形变时序相同的形变时间参考点，实现时间重叠的两个卫星平台形变时序的融合。S1A MT-InSAR形变时序相对与CSK MT-InSAR时序的相对偏差Δdef在重叠时间跨度的计算过程如下：

其中

IQR_S1A＝C_S1A|_p＝75％-C_S1A|_p＝25％ (10)

IQR_CSK＝C_CSK|_p＝75％-C_CSK|_p＝25％ (11)

p是C_S1A对应的累积分布函数；kk是p中的第k个对应的百分位数，其是通过找到C_S1A的数据值相等(或接近)给定的观测数据值

来计算；

和

是两个控制观测数据(即时间重叠时间段内观测值)的中值；g′和f′是调整校正后的时间序列分布的两个调制参数；IQR_CSK和IQR_S1A分别是参考形变时间序列和待校正形变时间序列的形变控制观测值的四分位距，分别为75％和25％。

步骤5：对有时间间隙的各个卫星平台MT-InSAR形变时序，以基于ENVISAT ASAR(ENV)数据集反演的MT-InSAR形变时序和步骤4融合好的CSK+S1A MT-InSAR形变时序为实例。其中ENV数据集获取的时间段为2007-2010，与融合好的CSK+S1A MT-InSAR形变时序存在三年的时间间隙。融合有时间间隔的各个卫星平台MT-InSAR形变时序，需要将各个平台的MT-InSAR形变时序与外部沉降经验模型(12)进行非线性拟合。在拟合时，采用了如下策略进行模型最优非线性拟合：首先对公式(12)中的模型进行改进，

d为有时间间隔的两个时间序列之间的差值，在该实例中以ENV时间序列的初始时刻为统一时间参考点，即为融合好的CSK+S1A时间序列相对初始时刻对应的相对形变量。m为已知参数，时间上先获取的形变时间序列(ENV时间序列)对应的m为值为0，时间上后获取的时间序列(CSK+S1A时间序列)对应的m值为1。

设置模型参数λ的取值范围为0～2，在进行模型拟合时λ的步长可设置为0.1。针对不同的λ值，采用迭代的莱文贝格-马夸特最小二乘最优化算法(Levenberg-Marquardt,LM)对有时间间隙的三个卫星平台MT-InSAR形变时间序列和外部沉降模型进行非线性拟合。对于每一个给定的λ值，都可以拟合得到一个由三个卫星平台MT-InSAR形变时序约束的形变模型，从所有λ值对应的模型中选择一个模型作为最优形变模型，拟合得到的三卫星平台MT-InSAR形变时序即为融合后的三卫星平台MT-InSAR形变时序。最优形变模型的评判依据是最优模型拟合得到的三卫星平台形变时间序列与根据公式(13)中参数m得到的MT-InSAR形变时间序列之间的均方根误差值，均方根误差值越小、模型非线性拟合效果越好，所得模型最优。

图4左图为实例上海东部新成陆区三卫星平台MT-InSAR融合后得到的累积沉降图，右图为三卫星平台MT-InSAR融合后的四个高相干点的形变时间序列。

本说明书未作详细描述的内容属于本领域专业技术人员公知的现有技术。

Claims (1)

1.一种基于多卫星平台MT-InSAR融合的新成陆区超长形变时序提取方法，其特征在于，该方法包括以下步骤：

步骤1：各个单一卫星平台MT-InSAR形变时间序列反演

利用高级时序合成孔径雷达干涉测量方法即MT-InSAR和各个卫星平台获取的新成陆区合成孔径雷达即SAR影像时间序列，分别反演新成陆区地表形变场和形变时间序列；

步骤2：提取各个单一卫星平台MT-InSAR形变时间序列的公共时间跨度和公共空间覆盖范围

在时间上，从步骤1得到的各个卫星平台MT-InSAR形变时间序列中提取出重叠时间范围内的形变时间序列子集；在空间上，提取各个卫星平台MT-InSAR地表形变场的公共空间覆盖范围，以及在此公共范围内的各个卫星平台MT-InSAR选取的高相干点，即时间相干系数阈值大于给定阈值；

步骤3：选取各个卫星平台MT-InSAR公共高相干点

各个卫星平台MT-InSAR选取的公共高相干点是基于像元的点目标，将各个单一平台MT-InSAR提取的高相干点进行空间匹配，将不同分辨率MT-InSAR形变场重采样到统一空间分辨率，统计和保存各平台MT-InSAR公共高相干点的像元位置并保留像元位置数据；

步骤4：新成陆区有时间重叠的两个卫星平台MT-InSAR形变时序融合

基于重叠时间段内两个卫星平台的形变时间序列子集的累积分布函数差异建立其偏差的统计分布特征，估算和去除后一个卫星平台MT-InSAR形变时序相对于前一个卫星平台MT-InSAR形变时序的偏差，在每一个公共高相干点上融合生成时间上连贯一致的长时间MT-InSAR形变序列数据集；详细步骤包括：

1)计算重叠时间段内两个卫星平台MT-InSAR形变时间序列子集的累积分布函数

设有N个卫星平台合成孔径雷达数据集，通过MT-InSAR方法，基于每个合成孔径雷达数据集反演的形变时间序列表示为：

d^(x)(P)为第x个合成孔径雷达卫星平台数据集反演得到的P点的形变时间序列；

第x个合成孔径雷达卫星平台数据集反演得到的P点的形变时间序列中第1,…Q_x个形变值；其中Q_x是第x个卫星平台合成孔径雷达影像的数量，x为合成孔径雷达卫星平台的序号，N是卫星平台的总数；S是步骤3中选取的各个卫星平台MT-InSAR公共高相干点的集合，P是该集合中任意一个高相干点；

若要融合两个在时间上有重叠的MT-InSAR形变时序，则x＝1,2,其中1代表在时间上先获取的合成孔径雷达数据集，2代表时间上后获取的合成孔径雷达数据集；第x个形变数据d(x)是相对于时间参考点

的形变量，时间参考点

的形变量为

要融合时间重叠的两个卫星平台MT-InSAR形变时序，需要确定第二个平台的MT-InSAR形变时序中每一个时间点的形变量相对于第一个平台MT-InSAR形变时序的时间参考点的形变量；

根据两个卫星平台的形变时间序列子集在重叠时间段内的配对形变数据，计算各自的累积分布函数；

2)基于累积分布函数匹配建立偏差分布特征，估算和去除MT-InSAR形变时序之间的偏差；

是基于第二个合成孔径雷达数据集反演的MT-InSAR形变时间序列中第j个形变值，需要将其校正到基于第一个平台合成孔径雷达数据集反演的MT-InSAR形变水准，去除相对偏差；j＝1,…,K，K是第二个合成孔径雷达数据集中影像的数量；改正后的第二个合成孔径雷达数据集反演的MT-InSAR形变时间序列中第j个形变值表示为

这一过程表示为：

Δdef_j是基于第二个合成孔径雷达数据集反演的MT-InSAR形变时序中第j个形变值相对于基于第一个平台合成孔径雷达数据集反演的MT-InSAR形变时序中第j个形变值的偏差；

要确定两个时序之间的相对偏差，首先，在两个形变时序的重叠时间跨度内，建立两个形变时序之间的相关的控制偏差量；公共时间跨度内的基于第一个平台合成孔径雷达数据集反演的MT-InSAR形变时序和基于第二个合成孔径雷达数据集反演的MT-InSAR形变时序分别定义为C_master和C_slave；C_slave相对于C_master的相对偏差Δdef的计算过程如下：

其中

IQR_slave＝C_slave|_p＝75％-C_slave|_p＝25％ (10)

IQR_master＝C_master|_p＝75％-C_master|_p＝25％ (11)

在公式(3)到(11)中，g′和f′是调整校正后的时间序列分布的两个调制参数；

和

是两个控制观测数据即时间重叠时间段内观测值的中值；

是两个控制观测数据即时间重叠时间段内观测值的中值的差值；

公共时间跨度内的基于第二个合成孔径雷达数据集反演的MT-InSAR形变时序中第k个形变值；

公共时间跨度内的基于第一个合成孔径雷达数据集反演的MT-InSAR形变时序中第k个形变值；Δ_kk是公共时间跨度内的第二个合成孔径雷达数据集反演的MT-InSAR形变时序中第k个形变值与基于第一个合成孔径雷达数据集反演的MT-InSAR形变时序中第k个形变值的差值；Δ′_kk是公共时间跨度内两个平台数据集反演的MT-InSAR形变时序中第k个形变值的差值与时间重叠时间段内观测值的中值的差值之差；kk是p中对应的第k个百分位数，是通过找到C_slave的数据值相等或接近给定的观测数据值

来计算；p是C_slave对应的累积分布函数；IQR_master和IQR_slave分别是参考形变时间序列和待校正形变时间序列的形变控制观测值的四分位距，分别为75％和25％；

经过以上估算，确定和去除后一个卫星平台MT-InSAR形变时序相对于前一个卫星平台MT-InSAR形变时序的偏差，将后一个卫星平台MT-InSAR形变时序校正到前一个卫星平台MT-InSAR形变时序的形变水准，融合生成时间上连贯一致的长时间序列MT-InSAR形变数据集；步骤5：新成陆区有时间间隙的各个卫星平台MT-InSAR形变时序融合

融合有时间间隔的各个卫星平台MT-InSAR形变时序，将各个卫星平台的MT-InSAR形变时序与外部沉降经验模型进行非线性拟合；新成陆区冲填土长期缓慢的固结过程会引起大范围不均匀地面沉降，通过土力学实验模拟建立的新成陆区长期缓慢地表形变经验模型(12)：

式中，t为新成陆区冲填土固结时间；S_CMT(t)为t时刻的累积沉降量；S_m是累积形变总量；k和λ为模型参数，决定沉降模型的曲率和曲线的形状；δ为时间延迟，对应于新成陆区地面沉降开始的时间；

在拟合时，首先对公式(12)中的模型改进为，

d为有时间间隙的形变时间序列间的形变差值；m为参数，先获取的形变时间序列对应的m为0，后获取的形变时间序列对应的m值为1；

进而，设置模型参数λ的取值范围为0～2，在进行模型拟合时λ的步长设置为0.1；针对不同的λ值，采用迭代的莱文贝格-马夸特最小二乘最优化算法对有时间间隙的各个平台MT-InSAR形变时间序列和外部沉降模型进行拟合；对于每一个给定的λ值，拟合得到一个由各个平台MT-InSAR形变时序约束的形变模型，从所有λ值对应的模型中选取最优形变模型，拟合后的各个平台MT-InSAR形变时序即为融合得到的多平台MT-InSAR形变时序；其中最优形变模型选取依据是模型拟合得到的各个平台形变时间序列与根据公式(13)中参数m得到的MT-InSAR形变时间序列之间的均方根误差值，即均方根误差值越小、模型非线性拟合效果越好，所得模型最优；

将各个平台的MT-InSAR形变时序与新成陆区地表形变经验模型进行非线性拟合后，实现新成陆区有时间间隙的多卫星平台MT-InSAR形变时序融合，生成时间上连贯一致的长时间MT-InSAR形变序列数据集。

Images