基于时序Sentinel-1数据的锦屏水电站左岸边坡形变探测与特征分析

王振林， 廖明生， 张路， 罗恒， 董杰

(武汉大学测绘遥感信息工程国家重点实验室，武汉 430079)

0 引言

我国西部山区由于其复杂的地形、地质及气候条件，已成为我国地质灾害最为频发的地区之一。但在这些区域蕴含着丰富的水力资源，因此建设有大量的水电站。这些水电站多位于深切的峡谷之中，也是地质灾害发生最为集中的区域。而在所有地质灾害中，滑坡是最常见的地质灾害类型，占地质灾害总量的70%以上[1-2]。这些滑坡与大坝的修建之间又有多种模式的相互影响关系。一方面，原有的活动滑坡及降雨或地震触发的滑坡会对大坝的安全运行构成威胁； 另一方面，坝区施工可能会引发坡体失稳，且蓄水造成的水位上涨也会诱发更多滑坡。因此，在水电站修建前的勘查选址阶段需要开展库区滑坡灾害隐患排查，尽量远离高风险灾害点； 同时在水电工程施工完成开始蓄水后，仍需要对库区边坡稳定性进行监测，以保证大坝的安全运行。

在水电站库区边坡的稳定性监测中，边坡形变的探测是滑坡灾害早期识别的关键。常见的形变探测方法包括水准测量、GPS和倾斜仪测量等[3]，但是这些方法普遍存在的问题是观测方式属于点观测，观测成本高且效率较低。而光学遥感容易受到云雾的影响，且不能有效探测缓慢形变，因此也具有一定的局限性。合成孔径雷达干涉测量(interferometric synthetic aperture Radar，InSAR)技术则能够有效克服以上问题，全天候全天时地实现大范围高精度的测量[4]。

目前，InSAR技术在形变探测方面的应用途径主要包括传统差分干涉测量(differential InSAR，D-InSAR)和时间序列InSAR分析[5]。传统的D-InSAR技术已经在一些库区滑坡形变监测中得到应用[6]，但对于植被覆盖较多的区域，该技术往往受到失相干的影响，无法准确获得有效的形变信息。相比之下，时间序列InSAR技术通过对多幅干涉图进行整合分析，能够有效减少失相干的影响，更准确地提取滑坡的形变信息。而常见的时间序列InSAR技术又分为永久散射体干涉测量技术(permanent scatterer interferometry，PSI)[7]和小基线集技术(small baseline subsets，SBAS)[8]。其中PSI技术采用单一主影像，在相干性较好的城区得到了广泛应用[9]，但在相干性较差的山区往往容易受到失相干影响，无法提取到足量的相干点。相比之下，SBAS技术采用多个主影像构造干涉图集，避免了单一主影像造成的时空基线过长而引起的失相干效应，因此该技术被广泛应用于地面沉降[10]和滑坡形变监测等领域[11]。

本文采用SBAS时序InSAR技术，对锦屏水电站库区左岸边坡进行形变探测，获取该边坡沿雷达视线向(line of sight，LOS)的累积形变量及年平均形变速率图，并对滑坡形变特征进行分析，进一步讨论该滑坡对锦屏一级水电站大坝安全构成的潜在威胁。

1 SBAS技术原理

假设有同一区域N+1幅时间顺序为t0，…，tN的合成孔径雷达(synthetic aperture Radar,SAR)影像，SBAS方法按照一定的时空基线及相干性阈值进行组合，最终生成M幅差分干涉图，并得到其解缠干涉图。再假设相邻最短时间间隔内的形变速率是线性的，则第i幅干涉图任意点P的解缠相位可表示为

(1)

式中：φi表示第i幅干涉图P点的解缠相位；a和b分别表示第i幅干涉图主辅影像对应的序号；vk表示影像序列中第k个时间间隔的平均相位速率。将式(1)扩展到所有干涉图，写成矩阵形式为

Φ=B·V，

(2)

式中：Φ表示解缠相位矩阵；B表示系数矩阵；V表示形变相位速率矩阵。

由于SBAS生成的干涉图子集具有多个主影像，若子集之间时间连续，则可采用最小二乘法将所有干涉图子集的形变速率连接起来得到最终的形变速率[12]； 若干涉图子集之间存在不连续的情况，则需要根据最小范数准则，采用奇异值分解的方法进行解算[8]。解算完成后再计算任意点P的总平均形变相位速率，即

(3)

2 研究区及其数据源

2.1 研究区概况

锦屏一级水电站位于四川省盐源县与木里藏族自治县的交界处，其主坝为混凝土双曲拱坝，坝高305 m，为目前世界同类坝型中的第一高坝，是雅砻江干流下游河段的控制性水电站[13]。水电站地处青藏高原向四川盆地过渡的大陆坡降带，属于典型的高山峡谷地貌地形。大坝选址为普罗斯沟与手爬沟之间，坝区海拔主要在1 600～3 000 m之间，坡度在20°～50°之间[14]。该区域雨季集中在每年的6—9月份。

锦屏一级水电站工程于2005年11月12日正式开工建设，枯水期水位约1 650 m。在2012年11月底建成后开始蓄水，分阶段蓄水至1 800 m，并于2013年9月正式投产并网发电，在2014年8月底进一步蓄水至正常水位1 880 m。本研究区位于大坝上游约1.5 km的左岸边坡上，总面积约200万m2，海拔高度在1 800～2 800 m之间。

2.2 数据源

选取了56景Sentinel-1干涉宽幅模式数据。该模式数据幅宽为250 km，分为3个子条带，每个子条带宽约85 km，由9或10个burst构成。本次选取了第2子条带中的第8个burst。同时选取22景ALOS-PALSAR存档数据，用于对比分析。2种数据集的具体参数如表1所示。实验研究范围如图1所示，图1中绿色框为Sentinel-1数据所选burst覆盖范围，黄色框为ALOS-PALSAR数据覆盖范围，红色框表示具体研究区域。参考数字高程模型(digital elevation model，DEM)采用SRTM 90 m空间分辨率数据[15]。

表1 星载SAR数据参数

3 结果与分析

3.1 相干点形变速率分布

采用StaMPS-SBAS[12]方法进行数据处理，从Sentinel-1数据集中组合生成了101幅干涉图，从中提取得到了16 764个相干点的LOS形变速率，如图2(a)所示。图2(a)中红色表示远离卫星的形变，蓝色表示朝向卫星的形变。从图2(a)中可以看出，库区左岸边坡有一块较大的滑坡形变区域。滑坡体顶部形变速率较大，最大值超过200 mm/a，而中部形变速率较小，大约在80～120 mm/a之间，底部形变速率最小，大约在30～50 mm/a之间，其他区域相对稳定。

(a) Sentinel-1数据(b) ALOS-PALSAR数据

对ALOS-PALSAR数据做相同处理，生成了35幅干涉图，共提取到2 762个相干点，形变速率如图2(b)所示。可以看出，该边坡在ALOS-PALSAR数据获取时间段内基本保持稳定。

3.2 形变场时空分布特征

在时序InSAR探测出的不稳定坡体上部取一点P(图2)。选取该点50 m范围内所有的相干点，计算其平均累积形变序列作为P点的累积形变序列。得到2种数据集下P点的LOS形变序列分别如图3所示。从图3(a)中可以发现该点在2015年10月26日—2018年4月13日期间LOS累积形变量超过500 mm。在一般的近岸滑坡中，滑坡的形变速率往往受到库水位变化和降雨的影响，呈现较明显的周期性[16]。但该滑坡近3 a间的形变并未表现出明显周期性，而是近似呈线性滑动趋势，这说明了其形变与水位波动及降雨没有显著的相关性。

(a) Sentinel-1数据 (b) ALOS-PALSAR数据

而从图3(b)可以看出该点在2006年12月11日—2011年2月6日期间基本保持稳定，说明该点在大坝蓄水前未发生形变。

为进一步研究滑坡累积形变量与高程之间的关系，计算Sentinel-1数据结果中各相干点的累积形变值，并采用反距离加权法对相干点的累积形变量进行插值，得到研究区域表面的累积形变分布，如图4所示。并选择图4中2条剖面线详细分析高程及累积形变变化(图5)。从图5(a)A-A′剖面线可以看出，山脊处的形变略大于山谷处，整个形变体最宽处宽度约为1 100 m。而从图5(b)B-B′剖面线可以看出，该滑坡顶部形变较大，中下部较缓，且最大形变集中在坡度转折点附近，整个形变体最长处长度约为1 500 m。

图4 Sentinel-1LOS累积形变

(a) A-A′剖面线 (b) B-B′剖面线

考虑到研究区域位于川滇山区，雷达回波相位容易受到大气传播延迟的影响。因此，利用GACOS(generic atmospheric correction online service for InSAR)方法提供的对流层水汽分布数据，对形变探测结果进行大气校正与分析[17]。绘制上述P点的大气延迟相位和月降雨量如图6所示，可以看出二者随时间的变化趋势有一定的相关性，都表现出随季节变化的周期性。

图6 大气延迟相位与降雨量

进一步比较大气改正前后该点的累积形变曲线如图7所示，可以看出经过大气校正后的形变累积曲线在夏季有一定的改善，但整体趋势保持不变。分析其原因主要是该滑坡体形变速率较大，形变信号相对于大气信号占主导地位，因此受大气影响并不明显。

图7 累积形变大气校正前后对比

3.3 形变驱动因素分析与风险评价

对比图3可以发现，该滑坡体所在边坡在2011年之前并未发生形变，而在2015年10月后的观测数据中表现出明显的形变，但由于2012—2014年间处于观测数据空白期，因此无法对其进行连续跟踪监测，确定具体的形变开始时间。考虑到大坝竣工后库区水位从1 600多m急剧上升到1 800多m，因此推测可能是巨大的水位上涨对两岸山体产生影响[18]，从而诱发了该滑坡复活。

根据形变速率图可以得出滑坡体的大致边界如图2中黑色虚线所示，经过量测得到其表面积约为75万m2。假设该滑坡体的平均厚度为2 m，那么一旦滑坡体发生滑落，将使得约150万m3的山体突然落入水库，从而引发严重的涌浪效应，威胁到大坝的安全运行。因此，为保证大坝的安全运行，需要加强该区域的实地踏勘及现场测量工作，以查明滑坡具体的滑动程度及深部位移情况，并开展长期连续监测，做好必要的防范措施和灾害预警。

4 结论

采用小基线集时序InSAR技术对锦屏一级水电站左岸边坡进行了形变探测与分析，获取了研究区雷达视线方向的平均形变速率及累积形变特征，验证了该技术在川滇山区滑坡形变监测中的有效性。通过对形变探测结果进行分析，可得出以下结论：

1)距离水电站大坝上游约1.5 km处存在一处表面积超过75万m2的大型滑坡体，且最大形变速率超过200 mm/a。

2)该滑坡形变出现在大坝蓄水之后，初步推测可能是由于大坝建成后库区水位的急剧上升对两岸山体产生了影响，从而诱发了该滑坡复活。该滑坡在后续的滑动特征无明显周期性，基本呈匀速线性滑动。

3)该滑坡滑动范围广，形变速率大，对于锦屏一级水电站的安全运行造成一定的威胁，需要加强该区域的现场监测和实地调查。

由于Sentinel-1观测数据的时空基线较短，干涉图相干性较高，小基线集时序InSAR技术能够精确地获取研究区形变信息的时间序列。随着雷达卫星传感器技术的不断发展，SAR影像的时空分辨率越来越高，小基线集时序InSAR技术将能够更好地应用于山区的滑坡形变监测中，为地质灾害防治工作提供更可靠更精准的参考信息。