差分干涉技术在变形观测中的应用

差分干涉技术在变形观测中的应用地面沉降监测是对地观测的方法之一，经典的方法是采纳精密水准测量，观测一、二等水准网，然后计算平差，得到微小地面沉降的变化值。近代GPS技术开头融入到地面沉降监测工作中。然而，受到若干条件的限制，GPS测量的数据密度往往不够，而且目前其空间域的辨别率也难以满意高空间辨别率形变监测的要求。

星载雷达差分干涉测量（D-InSAR）技术是20世纪90年月进展起来的1项对地观测新技术，其是利用雷达图像的相位差信息，精确测量出图像上每个点的三维位置，用于提取地面目标点微小地形变化的技术，精度能达到厘米级甚至毫米级，干涉图质量较好时，可以从干涉图识别出毫米级的形变。

D-InSAR技术的特点是面积大、快速、精确     、成本低，其有用性大大优于常规地面沉降监测技术。目前，SAR数据已掩盖了全球大部分区域，在无常规测量数据区域的地震监测中发挥着显著优势。笔者在深化讨论D-InSAR技术的基础上，对在地面沉降监测应用中的讨论，解决了D-InSAR技术中影像处理面临的一些关键问题，得到最终的地表形变结果。

1、D-InSAR的三轨法

通常，进行地面沉降差分干涉，需要猎取被观测地区沉降前后的2幅影像，通过比较得出形变量。

在2幅影像干涉图中，同时包含有相位信息和地形信息。欲得到形变量，需除去其中的地形信息，因此，必需先太清晰地面的DEM值，被称为“二轨法”。

在难以获得地面DEM值时，可采纳“三轨法”。三轨法需要用到同一区域的3幅影像，分别是形变前的2幅和形变后的1幅。然后，利用这3幅影像可生成2幅干涉纹图，其中，1幅呈现形变前的地形信息，另1幅呈现形变后的地形信息。同样的，对2幅干涉图进行匹配、重采样、去除平地效应和相位解缠，最终将2幅影像进行差分，便可得到地表沿视线方向的地形变化。

这方法的优点是无需地面信息，如DEM，对于一些无地形数据的变化监测尤为重要;缺点是相位解缠较简单，其质量将影响到最终的形变结果。

2、D-InSAR技术在巴姆地区的应用

2.1被观测数据的选取

针对不同的干涉应用需选取合适的SAR干涉影像对。伊朗东南部城市巴姆地区曾于2023年发生了强度为6.6级的大地震，为了利用D-InSAR技术猎取巴姆地震引起的地面形变，欧空局（ESA）供应了7幅ENVISATASAR数据。其中，4幅为降轨数据,巴姆古城位于图像中部，Track号为120;3幅为升轨数据,巴姆古城位于图像的近距端，Track号为385。

进行干涉处理，要求像对必需相干，因此，选取像对时必需考虑2个因素：a）临界基线距的限制。过短的基线将无法监测出形变变化量，过长的基线将导致2幅影像间的相关系数过低而无法形成干涉;b）存在时间间隔的影像。对于地形变监测为目的的干涉测量，其时间基线的选取要尽可能地短，以减小时间的去相关影响。

利用Eoli-sa软件进行数据选取时，依据被观测区域的经纬度和时间，查找合适的干涉数据。其中，Orbit确定轨道，景号（Frame）确定雷达照耀区域位置，航迹号（track）确定成像时间。用于干涉测量的SAR数据必需具有相同的track和frame。

从中看出，形变前的影像9192和6687以及形变后的影像9693具有相同的Track号。6687的Frame为3020，即中心结点号为3020，该幅影像的结点跨度3011～3029。由于，ENVISAT卫星的轨道被分割成7200个结点，通过与影像中心相近的400个结点来标识SAR的Frame号。因此，2个相邻的Frame用18个结点来标识。同理，9192影像的结点跨度3004～3021，即这2幅影像有共同的结点，即有共同区域，进而得出9693与9192也有共同的区域。

2.2试验过程

笔者实行三轨法进行差分处理，分别利用9192和6687构成的像对，9192和9693构成的像对，得出地震形变前后的地形，然后，再对这2个像对影像进行差分，获得地面点高程的形变信息。差分的关键步骤包括干涉影像的配准、干涉图的滤波、去除平地效应、相位解缠，最终把相位信息转换到高程。

2.2.1SAR复数影像配准

在雷达干涉处理中，配准的目的是将2幅影像上的对应点处理成地面上同名地物点。先采纳匹配方法找到2幅影像上的同名点，继而确定2点的偏移量，再对其中1幅影像进行重采样，使2幅影像上的同名点一一对应。

干涉雷达复图像的匹配，分粗匹配和精匹配两步进行[2]。确定图像搜寻窗内按行列不同的整像元偏移量为粗匹配，计算复图像对之间的匹配质量评价指标,得到精度在像元级的匹配;精匹配通过对复图像作亚像元插值，再计算匹配指标。理论上，匹配度需要达到子像素级（1/10像素）[3]。匹配策略如下：

a）粗匹配采纳基于卫星轨道参数的方法。由于卫星飞行的轨道参数是可确定的，且较精确     和稳定，因此是最基本的粗匹配方法。其优点是能自动完成偏移量的估算，缺点是匹配精度较低。

式中，Pm（x,y）和Ps（x,y）对应主影像和辅影像的行列值。粗定位就是利用卫星轨道参数计算同名点对的偏移量offset（x,y）;

b）像元级匹配。干涉雷达测量中复图像的像元级配准通常采纳基于窗口的配准技术,利用主图像和辅图像在空间域或频率域进行配准。笔者采纳的匹配相像性评价参数是较基本相关系数法，它也是其他许多匹配算法的基础，且计算起较简洁，不会消失太大误差;

c）亚像元级匹配。亚像元级匹配是在像元级匹配的基础上进行的。影像对（9192和9693）配准的统计分析见表1。

随便选取10个精配准点，检查其偏移量及精匹配的结果。匹配完成后，2幅影像上的同名点仅存在较小差异，用这些点作为掌握点对，用二次多项式拟合确定2幅影像中任何1个匹配位置的相对偏移量;

2.2.2干涉图的生成和处理

计算干涉系数由2幅影像上对应点共轭相乘得到，是作为评价干涉图质量的依据。干涉系数产生的相干图能直观地显示干涉图的质量，因此，相干图也可作为解缠质量的依据。通过处理得到，相干图上相干系数较大的区域是干涉效果较好的区域，干涉条纹较为明显，解缠也相对简单，反之，相干系数小或者干涉条纹不清晰的地方，一般不易解缠胜利。

从中看出，震前震后的2幅干涉图中垂直基线都很大，在干涉图相位中的山区部分即右下角，基线去相关现象较严峻，致使相干性很差，表现为黑色。

解缠前还要进行滤波，滤波的目的不仅为提高信噪比，且能为解缠计算供应较好的数据，在去除干涉

图噪声的基础上能保持边缘信息，防止产生误差。笔者采纳Goldstein滤波算法，去除噪声后条纹清楚，很好地保留了地形信息。图4为主、辅影像干涉图去除平地效应的结果，为进行dstein滤波的结果。

将震前后去除平地效应的地形相减，得到形变的地形。这时，显示的还是周期性的条纹，表示地形的相位还是缠绕的，没有得到真实地形，需要通过解缠来得到地形的真实高程。

2.2.3形变地形的相位解缠

相位解缠的结果直接影像地面沉降量的精度,是InSAR数据处理过程中的关键。相位解缠的基本思想从数学方面讲，假设给定1个二维相位矩阵Ψi,j，为了获得该矩阵的真实值，需要对每个点（i,j）加上或减去2π的整数倍，得到1个连续的函数Φi,j[6]，即

其中，ki,j为整数，-π≤Ψi,jπ，i∈[0,M-1]，j∈[0,N-1]。图6给出的是形变干涉的结果。

进行相位解缠计算时必需兼顾到全都性和精确性。全都性是指解缠后任意2点之间的相位差与2点之间路径无关，精确性则是指解缠后的相位要能忠实地恢复原始相位函数。将解缠相位梯度和缠绕相梯度之间的差异最小化，即最小成本网络流法，见图7。

3、试验结果

从干涉图可以看出,相隔175天和相隔35天的干涉图其相关性依旧保持很好,这是由于巴姆地区气候干旱,植被稀有，保持了很好的相干性。

从中可以看到，条纹变化较大的为形变区。依据条纹的分布也可得出形变结果。ASAR雷达采纳的是C波段微波成像，波长λ为5.666cm，可利用式（3）计算出，

每个周期的条纹变化代表了2.8cm的形变量。从图6中还可以看到，干涉条纹的颜色周期是不同的。依据GMT画出的图8、图9表示形变的三维状况。其中，北边区域为下