金沙江结合带巴塘段滑坡群InSAR探测识别与形变特征

杨成生，董继红，朱赛楠，熊国华

(1. 长安大学 地质工程与测绘学院，陕西 西安 710054； 2. 中国地质环境监测院，北京 100081)

0 引 言

滑坡作为一种极具破坏性的地质灾害，严重制约着区域经济的发展，并对社会公共安全构成严重威胁[1]。2019年中国发生滑坡灾害共4 220起，占全国地质灾害的68.27%，造成经济损失27.7亿元[2]。因此，开展滑坡灾害大范围调查与监测，掌握滑坡的分布及活动规律，对防灾减灾十分必要。合成孔径雷达干涉测量(Interferometric Synthetic Aperture Radar，InSAR)作为一种大范围地表形变监测新技术，已经被广泛用于地震[3]、地面沉降[4]、地裂缝[5]等地质灾害调查与监测，并取得了很好的效果。自Achache等利用干涉图对法国阿尔卑斯地区La Clapiere滑坡位移探测[6]以来，InSAR技术被广泛用于滑坡监测，大量学者开展了相关应用研究。Delacourt等利用D-InSAR技术和光学影像相关法，对法国Hellbourg附近的数个滑坡进行了监测[7]；Zhao等利用L波段的ALOS/PALSAR对美国加利福利亚州北部和俄勒冈州南部区域的滑坡进行了大范围的识别和监测，发现降雨是影响当地滑坡活动的主要因素之一[8]。刘晓杰等运用点目标分析(IPTA)技术获取了2017年中国四川茂县滑坡发生前的运动特征，并指出持续性降雨是诱发该滑坡的主要因素[9]。尽管InSAR技术用于滑坡识别和监测越来越广泛，但是该技术在实际应用中仍受到地表植被覆盖、地形起伏、大气延迟等因素的影响[10-19]。

金沙江结合带是川滇块体的西部边界断裂，断裂带走向总体近SN向，北段呈 NNW向，南段略向 NNE向延伸，是由多条主干断裂组成的一条复杂构造带[20]。金沙江结合带是中国滑坡灾害的高发区之一。由于流域内沿金沙江两侧地势陡峭，滑坡灾害的发生极易造成堵江事件，并诱发二次灾害，造成经济损失。例如，2018年10月10日和11月3日，西藏自治区江达县波罗乡白格村附近先后两次发生山体滑坡(简称“白格滑坡”)，造成堵江断流，形成堰塞湖，虽然该滑坡未造成人员伤亡，但造成的经济损失超过42亿元[21]。因此，开展金沙江结合带潜在滑坡的早期识别，对保护金沙江沿岸居民生命和财产安全，以及指导当地政府防灾减灾都具有重要意义。本文以金沙江结合带巴塘段为试验区，开展基于InSAR技术的潜在滑坡识别与监测研究，为金沙江结合带其他区域的滑坡调查与监测提供参考。

1 研究区概况及数据来源

1.1 研究区概况

金沙江位于长江上游，流经中国西藏、云南及四川等地，全长3 479 km。流域内坡陡水急，落差高达5 100 m，形成了切割强烈的深V峡谷。同时，金沙江蕴含丰富的水能资源，是中国重要的水电供给基地。然而，流域内岩体结构复杂破碎，软弱岩层发育，流域性特大高位地质灾害频繁发生[22]。研究区位于四川省巴塘县中心绒乡(图1)，属于金沙江的中游，区域内最大高差超过1 000 m。该地区是典型的金沙江干旱河谷气候，光热资源丰富，炎热少雨，植被不甚发育，水土流失严重，生态十分脆弱。

图1 金沙江结合带巴塘段滑坡隐患分布Fig.1 Distribution of Possible Landslides in Batang Section of Jinsha River Convergence Zone

先前研究表明，区域内受地震影响的滑坡堆积体广泛发育[23]，巴塘断裂和河流对金沙江两岸的地质灾害具有明显的控制作用，滑坡灾害呈带状密集分布。滑坡点主要发育在地形坡度15°～40°，相对高差大于500 m的地区。同时该区域也受金沙江断裂带影响，构造活动强烈，冲沟侵蚀下切严重，使得该区域成为滑坡、崩塌、泥石流灾害高发区[24]。

1.2 数据来源

为了探测和监测研究区内滑坡隐患点的分布及形变规律，本次研究共收集了覆盖研究区的欧洲航天局Sentinel-1A卫星升轨SAR影像64景，时间从2018年1月至2020年2月，同时还收集了69景获取于2017年11月至2020年3月的Sentinel-1A卫星降轨SAR影像。表1列出了研究中所使用的SAR数据集基本参数。研究中使用AW3D30(ALOS Global Digital Surface Model “ALOS World 3D-30m”) 数字高程模型(DEM)来消除In-SAR干涉处理中的地形相位及辅助SAR影像进行地理编码。AW3D30数字高程模型数据水平分辨率为 30 m(1″)，高程精度为5 m，是目前精度较高的地形数据之一。同时，该数字高程模型数据在本文中还被用来计算SAR数据的叠掩与阴影区域。Sentinel-1A卫星的POD精密轨道星历(POD Precise Orbit Ephemerides)数据被用来辅助Sentinel-1A数据的预处理和基线误差改正。

表1 SAR数据集的基本参数

2 研究方法

2.1 基于堆叠InSAR技术滑坡隐患识别

本次研究使用堆叠InSAR技术(Stacking-InSAR)进行大范围地区滑坡隐患探测。堆叠InSAR技术是由Sandwell等在1998年提出的[25]，其原理是对多幅差分干涉图的解缠相位进行加权平均，从而获取研究区域的形变速率。在解缠相位进行堆叠处理之前，需将所有SAR影像采样至同一坐标系下。同D-InSAR技术相比，该技术可以有效减弱大气延迟误差和数字高程模型误差的影响。堆叠InSAR技术的数学模型[26-27]可以写为

(1)

式中：v为年均相位形变速率；Δti为第i个差分干涉图时间基线；φi为第i个差分干涉图解缠相位值。

本文设置空间基线不大于150 m和时间基线不超过50 d，对升轨、降轨影像组合进行差分干涉处理。顾及到滑坡范围一般较小的特点，数据处理中采用了32×32像素窗口对干涉图进行自适应滤波来削弱噪声的影响。使用狄洛尼三角剖分(Delaunay Triangulation)最小费用流方法进行干涉图相位解缠，并采用GACOS在线大气改正系统进行残余大气延迟改正[28]。

由于SAR影像是采用侧视雷达成像，故在地形起伏较大的区域存在叠掩、阴影等几何畸变现象[图2(a)]。例如，当滑坡体位于叠掩区域时，会造成滑坡形变体透视收缩或监测点位稀疏，而阴影区是雷达视线(Light of Sight,LOS)向无法探测的区域，容易形成滑坡隐患判识的盲区[29]。因此，升轨、降轨组合可以提高大范围滑坡判识的准确度[图2(b)]，减少滑坡漏判率。为了避免在数据处理过程中对阴影区引入处理噪声，根据升轨、降轨卫星的几何参数计算了干涉图的阴影区域(图3)，并对每幅干涉对进行阴影区掩膜。以干涉图相干性高和解缠相位连续为条件，分别选取了183幅升轨干涉对和210幅降轨干涉对用于堆叠InSAR技术处理，从而获取研究区域的形变速率。通过对形变速率设置合适的阈值，初步圈定疑似滑坡区域，将初选结果叠加至光学影像进一步判识，获取最终潜在滑坡区域。

图2 SAR成像几何关系及升轨、降轨联合监测模式示意图Fig.2 Schematic Views of SAR Imaging Geometrical Relationship and Joint Monitoring Mode of Ascending and Descending Orbits

图3 SAR影像几何畸变现象分布Fig.3 Distributions of SAR Image Geometric Distortion

2.2 二维形变解算

由于SAR卫星侧视成像的特点，导致单一轨道InSAR结果很难准确反映真实形变特征。升轨、降轨SAR影像组合则为坡体形变监测提供了不同观测视角。因此，本次研究使用多维小基线子集(Multidimensional Small Baseline Subset，MSBAS)技术对升轨、降轨SAR影像进行联合处理，从而获取滑坡的二维形变。

多维小基线子集技术是由Samsonov等在2013年所提出的[30]。相比常规的一维小基线子集(Small Baseline Subset,SBAS)技术[31]，该技术可以计算地表二维(水平东西向和垂直向)变形速率和时间序列，称之为MSBAS-2D技术。

MSBAS-2D技术是根据多个传感器在相同区域和时间段内采集到的升轨、降轨SAR影像联合处理来获取地表二维变形时间序列。MSBAS-2D矩阵形式可以表示为

(2)

(3)

s={SE,SU}={-cosθsinφ,cosφ}

(4)

将Tikhonov正则化矩阵L乘以正则化参数λ，可实现时间序列的正则化，其效果类似于低通滤波，消除干涉对中高频信号噪声。当升轨和降轨数据获取时间不同时，需要进行正则化处理。未知量水平东西向速率和垂直向速率则通过应用奇异值分解来求取，形变时间序列则通过对形变速率的数值积分来重建。

本次研究数据处理中，首先按照堆叠InSAR技术中的设置处理，分别选取了183幅升轨和210幅降轨的高质量干涉图进行相位解缠。同时，将升轨、降轨解缠相位重采样至公共区域，再按照式(2)进行MSBAS-2D技术求解，进而获取研究区域二维形变速率和时间序列。试验技术流程如图4所示。

图4 技术流程示意图Fig.4 Schematic Chart of Technique Flow

3 结果分析

3.1 中心绒乡滑坡隐患点大范围判识

通过挑选高质量的相位解缠图，利用堆叠InSAR技术获取了金沙江结合带巴塘段中心绒乡沿卫星视线向年平均形变速率，结果如图5所示。其中，负值表示滑坡位移远离卫星方向；正值表示滑坡位移靠近卫星方向。通过对形变速率设置合适的形变阈值，可实现潜在滑坡的自动识别。就本次研究区域而言，考虑到监测误差的影响及区域内滑坡的活动性，依据经验选定形变速率20 mm·年-1为阈值，避免误差对结果的影响。当形变速率绝对值大于20 mm·年-1时，将被确定为疑似滑坡。将升轨、降轨数据自动获取的潜在滑坡点叠加至光学影像，并通过将InSAR形变监测结果与光学影像对比分析，最终获取潜在滑坡的位置与边界。本次确定的疑似滑坡灾害点位置和范围如图1、5所示。

图5 中心绒乡年平均形变速率分布Fig.5 Distributions of Annual Average Deformation Rates in Zhongxinrong Township

图6 典型滑坡遥感影像及InSAR形变监测结果Fig.6 Remote Sensing Images and InSAR Deformation Monitoring Results of Typical Landslides

选取了4处滑坡的InSAR形变监测结果与光学影像进行叠加，结果如图6所示。其中，图6(a)为旺各滑坡，从光学影像可以看到滑坡体形态明显，且InSAR形变监测结果与光学影像判识结果具有较好的一致性；图6(b)为霍荣村滑坡，该滑坡体呈长舌状，且在光学影像中可以看到明显的滑坡侧壁，InSAR形变监测结果可以明显地展示出滑坡的变形区域；图6(c)为那那贡村滑坡在光学影像中的形态，主要变形区域位于滑坡体后缘；图6(d)为锐哇村滑坡，该滑坡存在两个主要变形区域，变形较大区域位于滑坡体后缘，光学影像中可以看到明显的拉张裂缝。

3.2 中心绒乡区域滑坡群二维时序形变监测结果

从图5可以发现，中心绒乡附近升轨、降轨数据均存在明显的滑坡形变信息，其中包括贡伙村滑坡1和2、安里克米滑坡、仁娘村滑坡等多处。为此，选取该区域对升轨、降轨数据进行MSBAS-2D技术处理，获取了该区域2018年1月至2020年3月二维形变速率图(图7)。从图7(a)可以发现，所有滑坡在垂直向均表现为负值，表明滑坡体在沿近似坡向向下的主滑方向运动；从图7(b)可以发现，安里克米滑坡和贡伙村滑坡1在水平东西向表现为正值，也就是说，这两个滑坡体在水平方向上向东运动，结合实际地形可以发现这两个滑坡位于沿山脊线的东侧，即这两个滑坡体沿滑坡方向向下移动。

图7 重点区域垂直向与水平东西向年均形变速率分布Fig.7 Distributions of Vertical and East-west Annual Average Deformation Rates in the Key Areas

图8 中心绒乡滑坡群垂直向与水平东西向形变时间序列Fig.8 Time Series of Vertical and East-west Deformation of Landslide Group in Zhongxinrong Township

针对中心绒乡滑坡群，在滑坡体上选取了4个特征点(A、B、C、D)[图7(a)]，并提取了它们的形变时间序列结果(图8)。结果显示：安里克米滑坡在水平东西向为正值，结合地形可知该滑坡体近似向东运动，在两年时间段内最大累积量达到了44 mm；仁娘村位于仁娘村滑坡体下端，两年时间段内在垂直向最大累积形变量达到了88 mm，且仍然呈现出加速变形的趋势，后续仍需加强对该区域的监测；贡伙村滑坡是中心绒乡滑坡群比较大的滑坡，依据光学影像和InSAR形变监测结果显示该滑坡体面积约为1.9×106m2，滑坡形态为狭长带状。据此，选取了两个特征点进行时间序列分析。贡伙村滑坡1上特征点位于中心绒乡贡伙村小学附近，从时间序列结果可以发现该点无论是垂直向还是水平东西向，变形状态仍然呈现加速变形，其中垂直向最大累积形变量为80 mm。仁娘村滑坡在水平方向变形规律呈现先负后正再负的特征，结合地形发现该滑坡体呈近SN向，因此，沿着水平方向(东西向)运动比较弱小，同时受地表微地貌的影响，其水平运动表现为向西→向东→向西的位移过程。由于该区域没有实地监测数据，为了评价InSAR形变监测结果的精度，选取并计算了形变参考区[图5(a)]内所有点的形变速率标准差(表2)。从表2可以看出，对于升轨、降轨数据，形变参考区内的形变速率在0.634 mm·年-1以内，属于稳定区域，同时标准差未超过1 mm·年-1，远远小于滑坡隐患判识的阈值，证明了InSAR形变监测结果的可靠性。同时，提取了形变参考区形变值的二维形变时间序列(图8)，结果显示形变参考区内的时间序列形变在垂直向和水平方向上均位于0附近，最大变形量不超过4 mm，证明InSAR形变监测结果是可靠的。

表2 形变参考区内的形变速率标准差

3.3 典型滑坡活动特征

图9 贡伙村滑坡2的InSAR形变监测二维时间序列与降雨关系Fig.9 Relationship Between InSAR Deformation Monitoring Two-dimensional Time Series of Gonghuocun Landslide 2 and Rainfall

贡伙村滑坡2经纬度为(29°11′26″N，99°8′18″E)，面积约为3.1×106m2，长度为2 045 m，宽度为2 347 m，高差为924 m。2018年1月至2020年3月的升轨、降轨Sentinel-1A卫星监测结果均显示该滑坡形变特征明显。据此，在该滑坡体上选择了3个特征点(P1、P2、P3)[图7(a)]进行时间序列分析。其中，P1点位于滑坡体后缘；P2点位于滑坡体中部，也是InSAR形变量最大区域；P3点位于滑坡体下缘。从时间序列可以发现，这3个特征点二维形变具有较好的相关性(图9)。其中，P2点垂直向累积形变量最大，达到了77 mm，水平东西向该点仍然在加速变形；P3点水平东西向的累积位移量最大，其次是P2点，分析其原因可能为P3点位于滑坡体下缘，滑坡体的物质在此处堆积。

图10 贡伙村滑坡2的InSAR形变监测与光学影像解译结果Fig.10 Results of InSAR Deformation Monitoring and Optical Image Interpretation for Gonghuocun Landslide 2

从光学影像不难发现，贡伙村滑坡2滑坡体上存在数条冲沟及拉裂缝，这些地形特征为滑坡形变提供了有利条件。为了探讨降雨对该滑坡的影响，收集了与InSAR数据覆盖时间段一致的降雨数据，如图9所示。2018年，由于没有明显降水峰值，与之相对应的滑坡累积位移变化不明显，但是2019年6月至9月出现了较为密集的降雨(图9阴影区域)，出现了3个主要降雨峰值(H1、H2、H3所处位置)，并且在8月6日出现了最大降雨量(9.225 mm)，以此为区间中心，以P2点垂直向位移为主要分析对象(P2点位于滑坡体中部，具有较高的代表性)，分析连续降雨时间段内降雨对累积位移变化的影响。结果发现，P2点在峰值H1之前呈现线性加速变形状态，在出现第一个降雨峰值(H1)时，P2点处于稳定状态，在第二个峰值(H2)出现后，滑坡体在垂直向出现短暂平稳趋势，在第三个峰值(H3)之后，P2点又恢复了线性变化趋势，其原因是进入9月之后，该区域降雨减少。结合图9和图10可知，由于冲沟和拉裂缝的存在，为雨水入渗提供了便利，增大岩土体自重，使得岩土体在饱水状态下易于静态液化，降低岩土体的强度。

通过以上分析不难发现，强降雨对滑坡活动有一定的短暂影响。同时，先前的调查发现，该区域冻融作用明显，冻和融的往复降低了岩土体完整性与强度[32-33]，同时新构造活动、历史地震和人类活动(主要是农田耕种和工程削坡扰动)[34-35]也为滑坡发生提供了良好的地质条件，因此，该区域滑坡极易发生。

从滑坡的时间序列形变结果来看，该滑坡形变速率稳定，并在两年监测期间内保持持续位移，处于长期蠕滑状态。由于贡伙村滑坡2的滑坡体已发育有多条裂缝，加之滑坡所处的区域地质条件，使得在强降雨状态下该滑坡体极易失稳复活。因此，建议当地加强对该滑坡监测预警工作，防止发生重大灾害性滑坡事件。

4 结 语

(1)以金沙江结合带巴塘段为研究区，采用堆叠InSAR技术进行金沙江流域滑坡调查与监测，证明升轨、降轨Sentinel-1A卫星联合监测可以有效避免因SAR影像几何畸变造成滑坡隐患的漏判、错判。

(2)选取典型滑坡体利用MSBAS-2D技术获取了二维(水平东西向和垂直向)形变监测结果，克服了传统只能获取雷达视线向形变结果的局限性，发现中心绒乡滑坡群利用InSAR技术探测到的滑坡隐患目前处于缓慢蠕滑阶段，同时结合时间序列形变结果可知，贡伙村滑坡1和2处于加速变形。

(3)监测获取的中心绒乡滑坡群时间序列形变结果表明，该处滑坡群隐患目前仍处于蠕滑状态。由于滑坡群所处的地质条件脆弱，在强降雨条件下极易导致滑坡失稳。因此，对该滑坡群开展持续监测非常必要。

(4)本文仅从InSAR技术角度进行该区域滑坡监测，下一步将结合野外调查资料和现场监测结果多角度分析该区域滑坡变形特征。

2001年，我高考第一志愿报考了长安大学，并如愿以偿来到文明古都西安。学校报到的第一天，我用非常不标准的普通话办完了入学手续，从此长安大学便成了我的家。记忆深刻的是2002年夏天，中国队第一次进入足球世界杯，每天晚上宿管阿姨极配合地早早将电视机搬到宿舍楼前，黑压压的一群人搬着凳子围在电视前为中国队加油。校园里欢呼声、呐喊声、长吁短叹声，如音律一致的擂鼓在空中回荡。头脑灵光的舍友还包场了报告厅，我负责卖门票，最终的收益只够我们八个舍友一次欢快的聚餐。本科毕业后，我继续在母校攻读硕士、博士学位，并留校工作，在母校怀里一待就是二十年。这二十年里，陪伴着母校四季更迭、岁月流转，有幸目睹了母校的发展与辉煌。春风绿树，桃李满枝，值此母校七十周年华诞，祝愿母校长风破浪、风雨化虹、一往无前，再铸璀璨明天！