滑坡是在恶劣环境条件(如地震、强降雨等)下，坡体上不稳定岩土受自身重力作用沿滑动面向下运动的自然灾害现象^[1]。我国西南地区地势起伏较大、地质活动强烈、地层断裂带和高山峡谷地貌普遍发育，存在滑坡地质灾害隐患^[2]。

阿海电站自2012年开始在重点区域布设GPS监测网，对包括新建滑坡体在内的水库影响待观区滑坡体进行监测^[3]。自监测以来，新建滑坡体持续出现较为明显的形变且滑坡体下游靠近革囊渡大桥部位，公路挡墙和抗滑桩已严重变形。GPS观测周期长、监测成本高、效率和空间分辨率低，而覆盖范围广、空间分辨率高的雷达遥感观测数据能够较好地弥补该不足^[4]。本文尝试利用InSAR观测数据监测阿海电站滑坡形变及其演化过程，首先收集研究区2017~2019年共54景Sentinel-1A观测数据并进行处理，提取得到滑坡体变形速率和位移时间序列；然后利用GPS观测数据从均方根误差(RMSE)和相关系数对InSAR观测结果的可靠性进行验证^[5]；最后对滑坡变形机制及电站安全性进行分析。

1 研究区概况及实验数据 1.1 研究区概况

阿海电站位于云南省丽江市玉龙县(右岸)与宁蒗县(左岸)交界处，属金沙江中游河段。相关研究资料表明，阿海电站周围已发生多起山体滑坡等地质灾害事件^[6]，对阿海电站大坝安全构成严重威胁。本文选择新建滑坡为主要研究对象，并对InSAR监测结果进行可靠性分析，选择新建滑坡的主要原因为：1)新建滑坡坡向为EW向，形变主要表现为EW向和垂直向，能够在一定程度上克服InSAR对NS向位移不敏感的问题；2)新建滑坡体表面植被稀疏，具有大量的永久散射体，能够削弱时间失相干误差的影响；3)新建滑坡体GPS监测站布设较为完整，有利于后续InSAR结果的验证。

1.2 实验数据

本文实验数据包括Sentinel-1A卫星影像、谷歌地球光学影像、外部数字高程模型、精密轨道数据、GPS监测数据、地质资料和实地勘察资料等。Sentinel-1A卫星使用54景升轨SLC数据，时间范围为2017-07-04~2019-04-25；DEM数据采用美国太空总署(NASA)发布的SRTM-30数据，分辨率为30 m，用于去除地形相位；精密轨道数据使用哨兵卫星对应的轨道数据。本文收集2017-07-27~2019-04-04时段16个布置在新建滑坡体上的GPS站点数据，共获取10期静态观测数据。

2 时序InSAR数据处理 2.1 实验数据处理

本文收集研究区2017-07~2019-04共54景Sentinel-1A SLC数据，成像模式为IW，极化类型为VV。SAR影像数据处理采用SBAS-InSAR技术，对输入的54景SLC影像进行干涉配对，选择时间为2017-10-08的SAR影像作为超级主影像，剩余53景影像作为辅助影像，与超级主影像进行配准。时间基线阈值设置为60 d，时空基线如图 1所示，图中黄色点为超级主影像。采用Goldstein滤波去除平地和地形相位，相干系数设置为0.2，利用连接图编辑工具去除相干性较低的像对。最后经过2次SBAS反演及地理编码，获得完整的研究区形变值。

2.2 SBAS-InSAR技术原理

SBAS-InSAR技术是一种基于时间序列的InSAR分析技术，具体原理参考文献[7]，其技术路线如图 2所示。

对研究区内N+1幅影像按t₀, …, t_n顺序进行排序，选取一幅影像作为主影像，将剩余SAR影像与主影像进行配准，生成M幅差分干涉图，且满足：

$ \frac{N+1}{2} \leqslant M \leqslant \frac{N(N+1)}{2} $

(1)

将经过2个不同时刻的辅影像t_A和主影像t_B(t_B>t_A)获取的SAR影像生成第j景差分干涉图，其x方向坐标和r距离向坐标的干涉相位为：

$ \begin{gathered} \delta \varphi_{j}(x, r)=\varphi_{B}(x, r)-\varphi_{A}(x, r) \approx \\ \frac{4 \pi}{\lambda}\left[d\left(t_{B}, x, r\right)-d\left(t_{A}, x, r\right)\right]+\Delta \varphi_{\mathrm{topo}}^{j}(x, r)+ \\ \Delta \varphi_{\mathrm{atmo}}^{j}(x, r)+\Delta \varphi_{\text {noise }}^{j}(x, r) \end{gathered} $

(2)

式中，j∈(1，…，M)，λ为雷达信号中心波长，d(t_B, x, r)和d(t_A, x, r)为t_B和t_A时刻相对于t₀时刻雷达方向的形变量，Δφ_topo^j(x, r)为差分过程中的地形相位，Δφ_atmo^j(x, r)为大气延迟相位，Δφ_noise^j(x, r)为噪声相位。当忽略地形相位、噪声相位及大气延迟相位等因素时，可将式(2)简化为：

$ \begin{gathered} \delta \varphi_{j}(x, r)=\varphi_{B}(x, r)-\varphi_{A}(x, r) \approx \\ \frac{4 \pi}{\lambda}\left[d\left(t_{B}, x, r\right)-d\left(t_{A}, x, r\right)\right] \end{gathered} $

(3)

处理后的相位可用矩阵形式表示为：

$ \delta \boldsymbol{\varphi}=\boldsymbol{A} \boldsymbol{\varphi} $

(4)

式中，A为M×N矩阵。采用奇异值分解(SVD)可以求出最小二乘解，从而有效解决秩亏问题，以得到最终的形变量信息。

3 InSAR形变监测可靠性分析

由于InSAR地面沉降测量结果为形变场，其测量值大小是相对于基准点而言，而该基准点就是研究区内所有形变速率测量的参考基准，其自身必须为稳定或变化状况己知的点位。在形变量计算时选择GPS基准站作为形变基准点。基线处理采用GPS数据处理软件GAMIT和精密卫星星历，GPS点位分布如图 3所示，计算完成后将GPS监测结果视为真值，并将其3个方向的结果投影到LOS向，投影公式为：

$ \begin{aligned} & \text { LOS }=\left[-\sin \theta \cos \left(a-\frac{3}{2} \pi\right)-\right. \\ & \left.\sin \theta \sin \left(a-\frac{3}{2} \pi\right) \cos \theta\right]\left[\begin{array}{l} \Delta x \\ \Delta y \\ \Delta h \end{array}\right] \end{aligned} $

(5)

式中，θ为卫星入射角，α为卫星航向角，Δx、Δy、Δh分别为GPS沿NS向、EW向和垂向的形变量。

精度评定方法采用最近距离法，以GPS监测点为中心，建立半径为100 m的缓冲区，在该缓冲区内搜索高相干度的目标点，并将查找到的所有高相干点的累积沉降时序平均值与同期GPS监测结果进行对比分析，采用均方根误差(RMSE)和相关系数(R²)对结果进行精度评定：

$ \mathrm{RMSE}= \pm \sqrt{\frac{1}{N} \sum\limits_{i=1}^{n}\left(\Delta_{i}\right)^{2}} $

(6)

式中，Δ_i为GPS观测值与InSAR观测值的较差，N为较差值个数。

$ \begin{gathered} R^{2}(X, Y)=\frac{\operatorname{cov}(X, Y)}{\sigma_{x} \sigma_{y}}= \\ \frac{E[(X-E(X))(Y-E(Y))]}{\sigma_{x} \sigma_{y}} \end{gathered} $

(7)

式中，cov(X, Y)为X和Y之间的协方差。相关系数取值范围为[－1, 1]，若相关系数为0，则2个变量之间不具有相关性。

3.1 平均形变速率结果验证

利用16个布置在阿海水电站的GPS站点数据对InSAR平均形变速率结果进行验证，比较两者在重叠时段的年均形变速率，对比结果见表 1和图 4。

由表 1可以看出，新建滑坡体的GPS与InSAR形变速率RMSE为22.9 mm/a，其中单点最大偏差为42.8 mm/a，单点最小偏差为0.6 mm/a，且偏差值小于20 mm/a的测站占比50%，有3个测站偏差值大于30 mm/a。图 4为GPS与InSAR数据年均形变速率对比结果，可以发现，新建滑坡体年均形变速率在120 mm/a以内，当GPS年均形变速率超过80 mm/a时，InSAR与GPS监测结果会出现较大偏差(如XJ06站、XJ17站及XJ20站)，其余站点年均形变速率结果较为吻合且趋势接近。线性回归分析结果(图 5)表明，InSAR与GPS实测数据相关系数为0.898 4，表现出很好的一致性。

3.2 时序位移量结果对比验证

利用布设在新建滑坡体上16个GPS站点获取的位移序列，采用最优邻近法对InSAR数据进行提取验证，对比结果如图 6所示。由图可知，最大RMSE为58.5 mm，最小RMSE为4.8 mm，平均RMSE为25.3 mm。对各站点RMSE进行详细分析可知，XJ14站、GND2站、GND1站、XJ07站、RZLB站、XJ04站、XJ18站及RZLA站InSAR数据与GPS站点一致性较好，RMSE在0~20 mm之间，且在形变趋势上与GPS站点十分吻合。XJ03站、XJ05站及XJ09站RMSE分别为22.6 mm、28.7 mm、42.5 mm，出现偏差的原因为InSAR相干点与GPS站点距离过大，相干点无法有效反映2个GPS站点的真实形变状况。GND3站和XJ19站分别在2018-06和2018-09之后出现与InSAR监测结果完全相反的形变趋势，原因为GPS站点附近出现滑坡、崩塌等突发性形变，而选取的InSAR相干点未落入最优邻近方法范围内。在形变梯度较大区域，XJ06站、XJ17站及XJ20站RMSE出现较大偏差，分别为51.2 mm、53.6 mm、58.5 mm，其原因主要为InSAR技术自身具有局限性，相关站点处出现大于两景影像可探测最大形变阈值的形变，导致无法监测到正确的形变量。

为进一步评价InSAR形变结果的可靠性，图 7为各站点GPS数据与InSAR数据的相关性。由图可知，除GND3站和XJ19站外，其余站点GPS数据与InSAR数据具有较好的相关性，相关系数均大于0.8，其中XJ06站相关系数最高，为0.985；而XJ19站和GND3站相关系数出现明显偏差，分别为0.110和0.003，可能为InSAR技术未探测到滑坡体表面的突发性形变。

4 滑坡运动监测结果及电站安全性分析 4.1 滑坡运动监测结果分析

采用SBAS-InSAR技术进行数据处理，得到2017-07-04~2019-04-25阿海地区新建滑坡体形变速率分布，结果如图 8所示，图中形变速率为负值表示该区域远离卫星运动，形变速率为正值表示靠近卫星运动。从图中可以看出，研究区整体形变速率在－29.8~136.7 mm/a之间。为进一步研究新建滑坡体的构成机制，将其分为A、B两个滑坡体，并对2个滑坡体的空间分布特征进行详细分析。

4.1.1 新建滑坡体A

根据变形幅度，将新建滑坡体A分成Ⅰ、Ⅱ、Ⅲ三个子区域，形变量由滑坡体后缘向前缘逐渐递增。区域Ⅰ位于整个坡体前缘冲沟裂缝以北，形变最为严重，形变速率为85.8~136.7 mm/a，整体稳定性较差；区域Ⅱ位于滑坡体中部新建村左侧公路下段，与区域Ⅰ相比形变值较小，但仍处于相对不稳定状态，形变速率为30.6~52.1 mm/a；区域Ⅲ整体处于相对稳定状态，形变速率为15.1~30.5 mm/a。在3个子区域内选择3个特征点，并绘制位移时间序列曲线，由图 9可以看出，P₁和P₂点具有显著的变形趋势，累积形变量分别达216 mm和80 mm，P₃点则相对稳定，累积形变量仅为21 mm。其原因可能为2017~2019年水库蓄水且降雨增多，导致金沙江水位抬升，滑坡体前缘处于饱水状态，当金沙江水位骤降，滑坡体前缘受力性质发生改变，在水位突降情况下由于拖拽作用发生形变，而滑坡体中部(区域Ⅱ)由于前缘(区域Ⅰ)失稳产生局部坍塌，也发生相应形变。

4.1.2 新建滑坡体B

新建滑坡体B主要沿滑坡体中部冲沟裂缝分布，形变量明显高于周围区域，选取形变特征点绘制时间序列曲线，由图 10可以看出，P₄点最大形变速率为54.74 mm/a，累积形变量超过60 mm。光学影像及地质资料显示，该滑坡体主要由玄武岩构成，基岩裸露且风化现象强烈，岩体表面疏松，滑坡体B发生形变可能是由于降雨或水库蓄水诱发滑坡体抗剪强度下降所致。

4.2 阿海电站安全性分析

新建滑坡体位于阿海电站库区左岸，距离电站约37 km，顺河长约1.5 km，宽约600 m，堆积体总方量约6 970×10⁴ m³，属特大型堆积体。结合光学遥感及地质数据，对新建滑坡体形变现状进行野外实地勘察，现场照片见图 11。图 11(a)和11(b)为新建滑坡体B实地勘察照片，可以看出，滑坡体前缘靠近公路部分已发生明显滑动，XJ14站附近区域出现宽约3 cm的裂缝；图 11(c)和11(d)为新建滑坡体A实地勘察照片，可以发现，滑坡体前缘下方有许多坍塌碎石，而靠近革囊渡大桥附近路面出现较大裂缝，且公路挡墙已经变形。截至2019-04，新建滑坡体各监测点的形变量呈现随时间增大的趋势，说明滑坡体形变仍在持续，若金沙江上游出现特大暴雨，新建滑坡体可能发生滑坡，堵塞金沙江并形成堰塞湖，对上游居民及下游水电站造成巨大威胁。

5 结语

本文利用SBAS-InSAR技术对2017-07-04~2019-04-25阿海电站库区54景Sentinel-1A影像进行处理，获得研究区新建滑坡体及周边区域地表形变信息，通过对滑坡体坡面形变情况及InSAR技术可靠性进行验证分析，得到以下结论：

1) 对滑坡体表面GPS站点数据和InSAR数据从平均速率及时序位移量2个方面进行对比分析，结果表明，监测时段内两者的线性平均速率RMSE为22.9 mm/a，最小偏差为0.6 mm/a，相关系数达0.898 4；16个GPS站点平均时序位移RMSE为25.3 mm，在形变量较小区域InSAR数据与GPS监测结果基本相符，监测精度较高；在形变较为严重区域，由于受SBAS-InSAR自身技术限制，其与GPS监测结果存在一定偏差。

2) 根据新建滑坡体特征将其分为A、B两个部分，对2017-07-04~2019-04-25滑坡体形变速率及累积形变量进行分析发现，新建滑坡体A中区域Ⅰ形变最严重，最大形变速率达136.7 mm/a，整体形变量由滑坡体后缘向前缘逐渐递增；新建滑坡体B整体较稳定，主要形变位置发生在滑坡体中部，最大形变速率为54.74 mm/a。

3) 新建滑坡体形变呈线性沉降趋势，监测时段内出现显著的沉降加速现象。结合实地勘察发现，各GPS站点的形变量呈现随时间增大的趋势，说明滑坡体形变仍在持续，可能会影响下游阿海电站安全，建议加强巡查及相关形变监测。