基于时序InSAR技术的南水北调禹州段沿线稳定性研究

2022-10-18刘蔚

水利科技与经济 2022年10期

刘蔚

(山东科技大学测绘与空间信息学院，山东青岛 266590)

1 概述

南水北调中线工程是我国水资源调动的重大战略性工程，它的正式通水缓解了我国华北平原水资源短缺的问题，极大地改善受水区域的生态环境和投资环境，优化了我国水资源配置的格局，推动了社会经济的可持续发展。由于干渠距离长、流经范围广、流经区域地质环境复杂，使得干渠及干渠两侧容易出现不稳定形变。为了应对和预防不必要的损失，对干渠及沿线一定区域进行稳定性监测对于保障其安全运营具有重要意义[1]。

传统的水利工程形变监测多是利用水准测量、GPS测量等地面监测手段，观测周期长[2]，且只能获取到离散点的形变结果，难以反映水利工程沿线的整体形变规律，无法从宏观上进行掌握。而合成孔径雷达干涉测量技术自提出以来，就凭借其监测范围大、成本低、精度高以及全天候观测等优势[3]，被广泛应用于各个领域。众多专家学者更是依据InSAR技术探讨了大型线性工程沿线的地表沉降情况。刘朋俊等利用时序InSAR技术提取了南水北调辉县段湿润性黄土地质区的形变规律，并结合土壤湿度和降雨量进行了相关性分析[4]。徐东彪等利用DS-InSAR技术监测到小浪底大坝区域的形变监测结果[5]。郭在洁等利用SBAS-InSAR技术获取了地铁沿线的地表形变信息，并对地铁沿线的稳定性进行了分析[6-7]。

南水北调中线工程禹州段干渠需要穿越湿润性黄土、膨胀岩土以及煤矿采空区，区域内地质构造复杂，容易发生变形问题[8]。因此，本文利用SBAS-InSAR技术对覆盖禹州段干渠的27景Sentinel-1A影像进行了处理，获得干渠沿线的形变情况，并结合形变成因及形变趋势等对重点形变区域进行详细分析。

2 研究区概况及数据来源

2.1 研究区概况

南水北调禹州段干渠全长42.8 km，干渠自平顶山市郏县流入禹州市，途径鸿畅镇、张得乡、梁北镇、火龙镇、韩城镇、朱阁镇、郭连镇、古城镇等地后流入长葛市，大致呈西南-东北走向[9]。

2.2 数据来源

本次实验采用的数据为欧空局提供的27景Sentinel-1A SLC升轨影像，极化方式为VV极化，数据的时间跨度为2019年1月20日至2020年10月15日，选取影像成像的时间间隔为24 d，影像的入射角为38.96°。利用美国航空航天局提供的90 m分辨率的SRTM DEM进行轨道校正和地形相位去除，利用POD精密轨道文件来消除轨道误差带来的噪声影响。

3 SBAS-InSAR数据处理

3.1 技术原理

SBAS-InSAR技术将所有获得的SAR影像组合成若干个集合，再利用奇异值分解(SVD)方法将多个小基线集联合起来求解，以获得形变速率和时间序列形变的方法。该方法有效地削弱了空间去相干的影响，同时减小地形对差分的影响。该项技术的技术原理[6]如下：

假设有N+1个按照时间先后顺序且覆盖相同区域的SLC影像，对应的影像获取时刻记为t0，…，ti…，tN，通过设置相应的空间基线阈值和时间基线阈值，可以将N+1个SAR影像划分为L个短基线集。通过对各个短基线集内的影像进行差分干涉处理，可获得M幅差分干涉图，其中M满足的条件如下：

(1)

假设第j景干涉图是由tA和tB时刻获得的两景SAR影像干涉生成的，且已去除了地形相位分量，假设tB>tA，则j在原始的方位向-距离向的图像坐标系(x,r)处的差分干涉相位为：

δφj(x,r)=φ(tB,x,r)-φ(tA,x,r)

(2)

式中：λ为雷达波长；d(tB,x,r)和d(tA,x,r)分别为tB和tA时刻相对于初始时刻t0的视线像(LOS)累积形变量，因而有d(t0,x,r)≡0。

为了便于计算，式(2)中既没有考虑地形相位的影响，也没有大气相位的变化和失相干现象。

图6是不同粘结剂LFP电极的CV和EIS曲线。从图6a中可以看出，相比于PVDF和PVA-g-PAA，交联的PVA-g-PAA-c-5% PER的氧化还原电位差更小，对应的电流更大，说明对应电极的极化更小，可逆性更好，有更好的电化学动力学性能。图6b为使用不同粘结剂的LFP电极在0.2 C倍率下循环100圈后的EIS图，可以看出，PVA-g-PAA-c-5%PER的电阻要比PVDF和PVA-g-PAA的电阻要小很多。说明交联后的 PVA-g-PAA-c-5%PER拥有更低的锂离子电阻，可以减少LFP电极的极化，提供一个良好的动力学过程。

根据式(2)，M幅差分干涉图可以得到M个方程，该方程组的矩阵形式可表示为：

δφ(x,r)=Aφ(x,r)

(3)

式中：A是一个M×N阶矩阵；δφ(x,r)和φ(x,r)为N阶向量。

矩阵A是一个近似关联矩阵，它可以直接由已知的差分干涉图获取得到，若所有的数据都属于一个基线集，则L=1，M≥N,此时的A为N阶矩阵，当M=N时，方程有唯一解；当M>N时，方程组是超定解的，可用最小二乘法求解φ(x,r)的估计值：

φ(x,r)=(AAT)-1ATδφ(x,r)

(4)

当L>1，此时的ATA是一个奇异矩阵，则A的秩为N-L+1，此时方程有无穷多解。为得到唯一解，SBAS技术利用矩阵的奇异值分解法，求出最小范数意义上的最小二乘解，进而得到累积形变量。

3.2 数据处理

由于单景Sentinel-1A影像覆盖区域较大，为了提高数据的处理效率,首先将原始哨兵数据进行裁剪得到实验研究区域，然后对裁剪后得到的27景Sentinel-1A影像进行SBAS-InSAR处理。实验过程中，根据设置的时间基线阈值和空间基线阈值，自动选取2019年7月7日的影像为超级主影像，其他日期的影像为从影像。通过对获得的干涉对进行筛选剔除，共获得161个差分干涉对。为了提高辐射分辨率降低空间分辨率，本次实验设置方位向和距离向为4∶1的多视处理来抑制噪声；利用最小费用流法进行相位解缠，解缠的相关系数阈值设置为0.2，有利于相干性较低的区域更好地进行相位解缠；为了获得更清晰的干涉图，在滤波时选择Goldstein方法，在质量较好的干涉图上人工选择GCP点进行轨道精炼和重去平，GCP点要选在没有相位跃变且远离形变区域的位置，以免影响解缠效果；对解缠后得到的形变结果进行两次反演，利用第一次反演来估算形变速率和残余地形，在第一次反演的基础上，利用高通滤波和低通滤波来减小大气相位的影响，最终得到更为准确的研究区域年均形变速率和累积沉降量变化；对所得结果进行地理编码，并将其投影至WGS84坐标系下。

4 监测结果和分析

4.1 禹州市整体形变结果

图1 禹州市形变速率

4.2 禹州段干渠及两侧形变分析

为了分析禹州段干渠的稳定性，在干渠两侧沿线建立1 km的缓冲区，见图2(a)。可以看到，干渠沿线有不同程度的形变发生，大部分区域的形变速率在-20～19 mm/a，形变最严重的区域出现在古城镇，该区域的最大沉降速率可达-68.17 mm/a。为了全面了解干渠渠堤及两侧1 km范围内的形变情况，提取缓冲区内的全部特征沉降点，共计187 238个，绘制年均形变速率的频率分布统计图，见图2(b)。可以看出，有99.459%的特征点的年均形变速率集中在-20～19 mm/a，有0.541%的特征点的年均形变速率大于-20 mm/a。

图2 干渠形变分析图

由于禹州段干渠在建设时除了流经颍河等较大河流河谷外，沿线地表基本上均分布有湿陷性黄土，而黄土遇水后会破坏土体结构使其强度降低，容易使渠道地基产生湿陷变形。此外，禹州段干渠是典型的强膨胀岩土渠段，土体的含水量与孔隙之比较高，膨胀土在吸水膨胀、失水收缩的过程中，也会使渠道地基出现变形。

在干渠设计时，南水北调禹州段干渠经过原新峰矿务局二矿、禹州市梁北镇郭村煤矿、梁北镇工贸公司煤矿和梁北镇福利煤矿等4处采空区，累计长度3.11 km，在干渠选线之初，就基本选在采空区的稳定区内，且对采空区进行了注浆处理等[10]。在监测时段内，除古城镇区域外，其他区域也未监测到明显沉降。由此可见，建立在这4处采空区上的干渠还是比较稳定的，这也与《南水北调中线工程不良渠段风险评估》中所述一致，侧面反映了时序InSAR技术监测的可靠性。

4.3 古城镇重点形变区域分析

由图2(a)可知，南水北调禹州段经过古城镇时，其干渠北侧存在形变严重的区域，该沉降区域距离干渠较近，容易对干渠的稳定性产生影响，因此了解其形成原因及形变趋势对于保护干渠稳定性具有重要意义。由于沉降区域出现在古城镇干渠北侧，因此仅在干渠北侧建立3 km的缓冲区，见图3(a)。可以看到，在3 km的缓冲区范围内有两个比较明显的沉降区域，将其分别记为A区域和B区域。通过结合实际地理位置可知，A区域位于古城镇张堂村，累积沉降量最大值为-102.83 mm，位于E 113°33′2″，N34°12′9″处；B区域位于古城镇龙屯村，B区域累积沉降量最大值为-118.63 mm，位于E 113°32′10″，N 34°12′19″处。

图3 大沉降量区域分析图

据资料显示，古城镇地下蕴含上千万吨的煤炭资源，储量极其丰富，煤炭资源的长期开发会导致地面塌陷，会使地质结构发生变化，且A、B两区域内并无其它可导致大面积沉降的因素，据此推断A、B两区域均是由于煤炭开发形成采空区所导致的地表沉降。对A、B区域的累积最大沉降点进行时间序列分析，见图3(b)。在2019年初，沉降最大值点是有略微抬升的，这可能与禹州市2018年以来对矿山进行生态修复有关，但持续的采矿活动还是导致了地面沉降的进一步发展。

为了进一步探究研究时段内最大沉降点的变化有无减缓趋势，将累积沉降量与时间进行线性拟合分析。因为在选取Sentinel-1A影像时，数据的时间间隔为均匀的24 d，为了便于显示，选取每48 d为一个时间间隔，共计14个时间点。A区域最大沉降点随时间变化的线性拟合曲线见图4(a)，B区域最大沉降点随时间变化的线性拟合曲线见图4(b)。可以看出，无论是A区域还是B区域，随着时间的变化，最大沉降点的累计沉降量随着时间变化在拟合直线左右微小幅度地偏离，可见最大沉降点在研究时段内其沉降并无减缓趋势。因此，应加强对A、B两沉降区域的生态修复工作，通过回填、植树造林等方法减缓其沉降的发生，尤其是对于A区域西南侧，其沉降位置距离干渠直线距离小于500 m，更要注意加强生态修复，以免对干渠的稳定性产生不良影响。

4.4 干渠两侧危险性评价

为了直观地展示干渠及沿线一定范围内的风险性，研究干渠周边沉降对干渠的影响，依据《地面沉降雷达数据技术规程》(DD 2014-11)中指定的地面沉降严重等级分类方法，小于10 mm/a为低沉降区，10～30 mm/a为较低沉降区，30～50 mm/a为中等沉降区，50～80 mm/a为较高沉降区，大于80 mm/a为高风险区。依据SBAS-InSAR技术监测得到的形变速率，对干渠两侧1 km缓冲区范围进行危险性评价分区，结果见图5。由图5可知，大部分位于低风险区，张得乡和畅鸿镇区域零星分布较低风险区，且距离干渠有一定距离，古城镇区域存在较高风险区，虽然沉降的发展趋势距离干渠越来越近，但在研究时段内仍未威胁到干渠渠堤的稳定。由此可知，在研究时段内南水北调禹州段干渠渠堤整体稳定性良好。