基于改进PSInSAR技术的非城区地表形变监测

张红峰 刘 瀛,2

1 山西省第五地质工程勘察院，山西省临汾市广宣街26号，041000 2 华北理工大学矿业工程学院，河北省唐山市建设南路80号，063210

PSInSAR技术通过分析SAR影像序列中的时间稳定永久散射体(permanent scatterer,PS)进行地表形变监测，广泛应用于城区、滑坡等形变监测及高程修正[1-2]领域。但非城区中稀疏的稳定散射体导致在地表监测过程中PS点分布密度不足，难以精确估计速率及残差等信息，从而使PSInSAR技术在非城区地形监测应用中受限[3]。为提高非城区PS点的分布密度，相邻影像连续干涉对[4-5]、小尺寸PS点结合角反射器优化PS点分布[6]、残差估计与解缠绕模型优化提高空间相干PS点识别率[7]、多尺度频率自适应相干估计提高PS点识别[8]等提高PS点分布密度的改进算法被相继提出。

在已有研究基础上，为解决非城区PS点分布不足而导致PSInSAR技术在非城区地表形变监测误差较大的问题，提出基于分时散射体(partial time scatterer, PTS)的改进PSInSAR算法。首先对SAR影像进行基于改进经验模态分解(empirical mode decomposition, EMD)算法的滤波降噪，然后采用双阈值与可信概率估计进行PTS提取，最后通过参数差分迭代估计法对PTS目标进行相位分离和形变速率估计，从而得到非城区监测区的地表形变。

1 边缘保持平滑滤波

在经典PSInSAR算法的SAR影像干涉对(i,k)中，像元P的邻域像素具有平稳随机过程特性[9]，因此其相干系数可表示为：

(1)

SAR影像中通常含有较多的干扰噪声，因此在PTS目标提取前需对SAR影像进行滤波降噪。将传统EMD算法直接应用于SAR影像时，易造成对SAR影像边缘区域像素值的过度滤波[5]，并且其模式混叠及以经验为依据的分量识别会使重建信号的准确度较低。为此，在传统算法行列分解的基础上，引入45°和135°两个分解方向，同时为优化IMF分量分界时相关系数法的复杂性和单一指标的不确定性，本文联合SAR影像的细节信息与逼近信息，采用基于均方根R和平滑度r估计的IMF分界K值自适应定位方法。假设干涉图经过改进EMD处理后的分量为dj(x)，等权滤波并重构后可得：

(2)

式中，z和ω(x)分别为平滑滤波器的系数及窗口值。R和r的表达式为：

(3)

(4)

式中，x(t)为含噪SAR干涉对影像序列。归一化R和r，并由变异系数法进行赋权处理，即

(5)

式中，μ和σ为相应的均值与标准差。由此可得，IMF分界的复合指标为Fo=WR×PR+Wr×Pr，相应的分界阈值为：

Tk-1=|Fo(k-1)/Fo(k)|

(6)

则滤波时的约束权系数ω(x)为：

(7)

2 PTS目标提取

对SAR影像干涉图进行边缘保持滤波后，首先由离差指数和相干系数对PTS目标点进行联合初识，获得PTS的候选目标[9]，以缩小后续目标点的搜索范围，然后由目标点的相位信息推导其可信概率作为PTS可靠筛选的依据。

(8)

(9)

(10)

式中，Δω为观测值与拟合值的相位差，即

(11)

进而可得候选目标为可靠PTS的概率为：

Pg=1-[(1-α)pR(rg)/p(rg)]

(12)

式中，p(·)与pR(·)分别为PTS像元与非PTS像元的概率函数，α为调节系数。

3 基于PTS的PSInSAR监测

根据大气自相关特性对由提取的PTS目标构建的Delaunay三角网进行二次差分，并建立三角网连接边上像元的观测方程[8]。在提取PTS目标点时，由于PTS目标的季节、气象等的分时稳定性是主要参考因素，传统的基于共同主影像PSInSAR算法的参数解算过程在解析PTS贡献时较为困难，因此借助连续干涉对思想[4]，将干涉图聚类为N-1个具有不同时相特性的组合(图1)，并以相干度作为组内各影像贡献度的衡量权重，以削弱低相干影像的贡献。

图1 影像干涉对时间基线及其分组Fig.1 Image interference pair time baseline and its grouping

基于以上分析，以提取的PTS目标对时序εp进行解算，其解算模型为：

(13)

(14)

采用经典PSInSAR算法进行初始值解缠，通过Delaunay迭代求解高程修正与形变速率，即可解算出非线性形变相位[2]。

4 实验分析

为验证本文算法在非城区形变监测中的有效性，以陕西靖边县周边区域作为测试实验区(图2)，实验数据为24景C波段干涉宽模式的SLC影像数据，采用ESA精密轨道数据处理配准及相关相位去除等操作。

图2 实验区光学影像Fig.2 Optical image of the test area

4.1 PTS目标点提取实验

采用本文算法和传统PSInSAR算法分别提取PTS点和PS点，实验结果如图3所示。从图中可以看出，在A区域，建筑物等可被较好地识别为PTS点和PS点；在B区域，植被较为稀少，PS点因失相干性而较少，密度较低，但本文算法仍可识别出较多的PTS点。对PTS点和PS点进行编号和经纬度重叠分析后发现，大部分传统PS点均被标记为PTS点，即这些同名的PTS点与PS点具有相同的空间分布特性，说明本文的提取算法具有有效性。

4.2 实测数据对比分析

为验证本文算法监测沉降的精度和可靠性，以实验区内2016年二等水准资料为基础，复测2016年布设的水准路线，以复测水准实际值作为精度和可靠性实验数据的基准参数，并将水准监测的垂线形变方向投影到SAR卫星LOS方向，选取与实验时间相同的水准历史资料对复测数据进行精度评价。

DSInSAR地形形变监测技术引入了雷达最新的分布式识别和滤波技术，可将裸土、沥青路面，特别是新围垦区等区域识别为DS点，因此将DSInSAR技术引入到实验中作为一种性能比较算法。提取监测区内水准点对应的形变数据点，选用平均误差ϑ和中误差m作为评价指标[10]。

表1(单位mm/a)为实验结果，从表中可以看出，在A区域，由于城区内永久相干目标较多，3种算法均获得了较好的误差精度；在B区域，PS点较少，形变监测精度的提高需依赖PTS点或DS点，因此本文算法和DSInSAR技术算法具有明显的优势，但本文算法更优，原因为本文算法的PTS点可充分利用测量点在不同时刻的贡献优势；在C区域，PTS点极少，沥青路面、新围垦区等具有雷达分布特性的DS点也较少，因此3种算法的监测误差均较大，而本文算法通过分时分析，在测量点较少的情况下充分利用各测量点在特征最优时刻的监测贡献度，避免监测点在特征较差时段对精度的不利影响，从而有效提高了形变监测的精度。

图3 PTS与PS提取实验对比Fig.3 Comparison of PTS and PS extraction experiments

表1 实验区内各算法的监测精度

图4为B区域水准测量结果与各算法监测结果的互差直方图，进一步验证了本文算法的监测精度较好，结果具有可靠性。

图4 实测值与各算法监测结果互差直方图Fig.4 Histogram of mutual difference between measurement and monitoring results of each algorithm

5 结 语

为解决PSInSAR技术在非城区监测误差较大的问题，提出基于PTS的改进PSInSAR算法。该算法通过改进的EMD对SAR影像进行边缘保持平滑滤波降噪，通过可信概率联合提取PTS目标；然后通过迭代估计PTS相位和形变速率，从而实现非城区地表形变的准确监测。实验结果表明，本文算法在非城区可获得比PSInSAR技术和DSInSAR技术更优的监测精度和可靠性，算法具有有效性。