2.State Key Laboratory of Information Engineering in Surveying, Mapping and Remote Sensing, Wuhan University, Wuhan 430079, China 3.Beijing North-Star Digital Remote Sensing Technology Co. Ltd, Beijing 100120, China

In recent years, geological disasters like landslides have often occurred within the Markam area, which is located in the southeastern part of Tibet Autonomous Region. The stability of these landslides is affected by factors such as physical geographical/meteorological conditions and anthropogenic activities. These factors pose great threats to power grid construction/operation, major transportations, and public security. Hence, effective technologies are needed to detect potential landslide hazards within this area and provide key supports for decision making regarding disaster prevention and reduction. This study employed the small baseline subset InSAR technique to process archived SAR data stacks acquired by ALOS PALSAR and ENVISAT ASAR covering the Markam area. A few suspected landslide hazard sites are identified along National Highway 318 and the Valley of Jinsha River. Spatial distribution and temporal evolution patterns of the surface displacements upon these unstable slopes are characterized.

remote sensing , landslide detection , InSAR , SBAS , time series analysis , displacement characters

滑坡是指边坡上的岩土体在自然环境和人类活动等因素的作用下，沿着贯通的破坏面整体或分散向下滑动的现象( 黄润秋和许强，2008 )。滑坡作为一种山区常见的地质灾害，具有分布广、危害大的特点( 李媛 等，2004 )，而且在演化成灾前通常会产生表面变形现象。因此，开展以形变探测为主的滑坡隐患识别监测，可以为滑坡灾害的防范提供科学的决策依据，对于有效规避灾害风险，保障人民生命财产安全，具有十分重要的现实意义。

目前滑坡形变监测通常采用多种监测手段相结合的方法，例如GPS测量、水准测量、遥感影像解译、伸缩计深部位移监测等( Wang和Li，2009 )。但这些方法通常受到地形条件复杂、人力物力成本高等因素制约，仅适用于小范围区域，且大多只能获取离散点观测数据，无法进行大范围高密度的形变监测。作为一种新型对地观测手段，雷达干涉测量InSAR(Synthetic Aperture Radar Interferometry)采用星载主动式微波遥感成像传感器，能够全天时、全天候成像，大范围、高精度获取地表形变信息，可在滑坡监测应用中对传统方法构成很好的补充( 廖明生和王腾，2014 ； 史绪国，2016 )。

藏东地区位于西藏自治区东南部，是从川西和滇西北进入西藏的咽喉门户，同时也是川藏高速公路、川藏铁路、川藏电网联网工程等规划建设中的国家重大基础设施的必经之地，地理位置十分重要。然而，该区域地形复杂，气候多变，地质构造运动活跃，地质环境高度破碎，被公认为是中国滑坡、崩塌、泥石流等地质灾害最严重的地区之一。采用有效的技术手段，对该地区分布的滑坡等地质灾害隐患进行识别、编目、监测和预警，对于开展防灾减灾，主动规避灾害风险，保障人民群众生命财产安全和重大基础设施建设运行安全具有重要意义。

目前对藏东地区地质灾害的研究主要集中在对已知的灾害事件进行记录，对已有资料的补充完善，以及探讨环境影响因素与滑坡形变和发育特征之间的相互关联。近些年来，随着遥感技术的普及，有学者采用遥感影像对该地区滑坡和泥石流灾害进行研究，综合利用数字正射影像(DOM)和数字高程模型(DEM)数据，发现变化区域，对感兴趣地区进行定性定量分析( 陈楚江 等，2004 ； 张明华，2004 )。但上述研究主要靠人工解译识别滑坡，得到的信息存在较大的不确定性，对于研究滑坡灾害特征还缺乏普遍适用性和可靠性。

本文采用StaMPS软件包中提供的小基线集SBAS(Small BAseline Subsets)时序InSAR方法，处理分析了2007年—2011年获取的覆盖藏东芒康县海通沟和金沙江河谷地区的18景L波段ALOS PALSAR影像数据，探测到了该地区几处高度疑似的滑坡体，并利用同时期获取的17景C波段ENVISAT ASAR数据，开展了滑坡形变探测结果的交叉验证，分析了这些不稳定坡体在时间序列上的形变特征。

小基线集InSAR是在差分干涉测量(DInSAR)和永久散射体干涉测量PS-InSAR(Persistent Scatterers InSAR)技术的基础上发展起来的一种时间序列InSAR分析方法。传统的DInSAR容易受到时空去相关及大气效应的影响，而PS-InSAR技术难以适用于相干目标点分布稀疏的地区。针对上述问题， Berardino等(2002) 提出了基于小基线集的时序InSAR方法。

小基线集方法的基本原理如下：假设雷达卫星重复轨道观测获取了覆盖同一地区的 N +1幅SAR影像，选取其中的一幅作为主影像，将其余的 N 幅影像配准到主影像的像空间上。通过设置时间和空间基线的合理阈值，对序列SAR影像进行组合，处理生成 M 幅相干性较好的差分干涉图。

这里假设输入为解缠后的序列干涉图，对于由 t _A 和 t _B 时刻获取的主从影像生成的第 j 幅差分干涉图，其中任意一点的干涉相位可以表示为 $ \begin{split} \delta {\varphi _j} =& {\varphi _{{t_{\rm{A}}}}} - {\varphi _{{t_{\rm{B}}}}} = \frac{{4{\text{π}}}}{\textit{λ}}\left( {{d_{{t_{\rm{A}}}}} - {d_{{t_{\rm{B}}}}}} \right) + \\ & {\rm{\Delta }}\phi _j^{{\rm{topo}}} + {\rm{\Delta }}\phi _j^{{\rm{atm}}} + {\rm{\Delta }}{n_j} \end{split} $

式中，1≤ j ≤ M ， φ _{t A} 和 φ _{t B} 是 t _A 和 t _B 时刻相对于 t ₀ 时刻的累积形变相位。Δ ϕ _j ^topo 表示高程误差相位，Δ ϕ _j ^atm 表示大气相位，二者可以按照永久散射体干涉测量中的方法进行估计去除( Ferretti 等，2000 ； Hooper 等，2004 )。Δ n _j 为处理过程引入的随机噪声误差。去掉差分相位中的高程误差和大气效应分量后，式(1)可以简化为 $ \delta {\varphi _j} = \frac{{4{\text{π}} }}{\textit{λ}}\left( {{d_{{t_{\rm{A}}}}} - {d_{{t_{\rm{B}}}}}} \right) + {\rm{\Delta }}{n_j} $

假设观测序列中两个相邻时间点之间的形变是线性的，即整个时间段内是分段线性的。那么第 j 幅干涉图对应的形变相位观测值可以写成 $\delta {\varphi _j} = \mathop \sum \limits_{k = {t_{{\rm{A}} + 1}}}^{{t_{\rm{B}}}} \left( {{t_k} - {t_{k - 1}}} \right){\nu _k} + {\rm{\Delta }}{n_j} $

将所有干涉图的差分干涉相位组合成矩阵形式,可以得到

式中， B 是一个大小为 M × N 的矩阵，记录了所有干涉对的时间基线长度。 ν 为以相位表示的待估计求解的各时间段上的形变速率， δϕ 为所有干涉图的差分干涉相位， ν 和 δϕ 均为一维向量。 Δ 为待最小化的随机误差向量。由于在选取干涉组合时，设置了空间及时间基线的阈值，因此有可能出现孤立不连通的干涉图子集，此时矩阵 B 就会成为秩亏矩阵。采用奇异值分解方法SVD(Singular Value Decomposition,)，可以在矩阵秩亏的情况下求出 B 的广义逆矩阵，从而求得最终的形变速率向量。

SBAS方法采用多主影像策略，能够组合出足够多的干涉像对，增加有效的观测数据量，而且由于选择的干涉对时空基线相对较短，在一定程度上可以减小由长基线所导致的失相干。因此，SBAS方法在SAR影像数量较少的情况下仍能够较准确地提取形变速率。

芒康县位于西藏自治区东南部，昌都市最东部，海拔高度在4000—4500 m，地处川、滇、藏3省交汇处。位于98°00′—99°05′E，28°37′—30°20′N之间。东与四川省巴塘县隔金沙江相望，南与云南省德钦县毗邻，西与左贡县接壤，北与贡觉、察雅两县相接。芒康县东部海通沟地区为高山峡谷地貌，河流下蚀作用强烈，新构造活动频繁，是中国地质灾害最严重的地区之一，被誉为川藏线上的“鬼门关”( 王培清和黎普明，2002 ； 杨志法 等，2006 )。

据相关新闻报道，芒康县境内318国道海通沟至朱巴龙段近年来多次发生滑坡、泥石流、塌方等灾害，造成该路段多处路基被冲毁，上百米路面被泥石流堆积物覆盖，交通中断，严重威胁当地人民生命财产和过往车辆安全。

为探测实验区内分布的滑坡灾害隐患，收集了两组不同的历史存档星载SAR数据，包括日本宇航局(JAXA)的ALOS卫星搭载的L波段PALSAR系统和欧洲空间局(ESA)的ENVISAT卫星搭载的C波段ASAR系统观测获取的数据。

图1 所示为两种SAR数据的大致覆盖范围。其中两个方框分别表示升轨获取的PALSAR数据和降轨获取的ASAR数据集的覆盖区，从东北折而向西的曲线表示318国道，黑色圆点标记出了研究的重点区域海通沟地区，位于海通沟东侧从北向南纵贯实验区的河流为金沙江。选取的影像数据和干涉对基本信息参见 表1 和 表2 ，零基线对应的为主影像。本实验使用的PALSAR数据均为零级产品，其中双极化FBD模式数据在聚焦成像前需进行距离向两倍过采样，以获得与单极化FBS模式相同的分辨率，从而使得两种观测模式数据混合干涉处理成为可能。本实验中PALSAR数据成像处理采用的是由JPL开发的开源软件ROI_PAC，时序干涉处理分析则采用了英国利兹大学Hooper教授开发的StaMPS软件包( Hooper 等，2004 )，其中单对数据差分干涉处理调用了荷兰Delft理工大学开发的开源软件DORIS。

在山区滑坡形变探测中采用PS-InSAR方法往往受植被覆盖和时空失相干效应的影响，提取到的高相干点目标稀少，相比之下采用SBAS方法效果较好。 图2 所示为本实验所用两个数据集分别采用小基线集策略生成的干涉对组合及其时空基线分布。图中圆圈代表单次SAR观测，连线表示干涉组合。

L波段PALSAR数据集临界基线较长，设置阈值为时间基线1000天，空间基线1800 m，共生成24个干涉对。C波段ASAR数据集的临界基线较短，采用阈值为时间基线400天，空间基线600 m，共获得42个干涉对。对这些干涉对进行差分处理生成差分干涉图，然后建立高相干点目标观测方程，最后利用奇异值分解估算出高相干点形变速率并分解残余相位。

已有研究表明，去除DEM误差和线性形变速率后，残余相位表现为非线性形变和大气相位( Hooper 等，2004 )。大气相位在空间上具有较强的自相关性，表现为低频信号，而在时间上不相关，表现为高频信号，因此可采用空间低通、时间高通的滤波器，来估算分离大气延迟相位，将其从残余相位中减去，可得到非线性形变相位。线性与非线性形变分量相加即为完整的形变相位。

采用SBAS分析提取到高相干点目标上的视线向形变速率分布如 图3 所示，背景为DEM山影。 图3 中红色代表形变方向为沿视线方向远离卫星，而蓝色表示沿视线方向朝向卫星。PALSAR数据为升轨右视观测获取，而ASAR数据采用的是降轨右视观测，二者的视线方向存在很大差异，因此两种数据集探测到的形变速率无法进行定量化的直接比较，而只能通过定性对比两个形变探测结果的空间分布来判断它们之间是否具有一致性。

另一方面，由于时序InSAR处理实现的是一种相对测量，需要依赖于独立的参考基准才能将形变测量结果绝对化。StaMPS软件在数据处理时默认采用将整个场景中所有测量点的平均形变速率作为相对参考基准的处理策略，因此在对其解算结果进行解译时必须考虑到这一点，通过分析一组空间聚集的测量点上的形变速率是否具有相似性并表现出相对于周边背景的显著差异，来判断疑似滑坡的存在。

从PALSAR数据集处理结果中，我们判读识别出5处主要的疑似滑坡分布区域，如 表3 所示。这一结果的有效性，将在后文中选取白江岗和草地贡两处区域为例，采用高分辨率光学影像目视判读解译来加以验证。

从 图3(a) 中可以看出，PALSAR数据集提取出的高相干点目标分布较均匀，雷达视线向(LOS)形变速率在－40—40 mm/a。在该区域约700 km ² 的范围内，总共提取出了844464个相干点目标，其分布密度约为每平方公里1200个，能够满足时序InSAR技术在山区应用的需求。

目前InSAR形变探测结果的验证一般采用水准或GPS观测数据进行对比分析。由于本研究没有此类地面观测数据可供参考，因此采用了PALSAR与ASAR两种不同数据InSAR探测结果相互比对进行交叉验证的途径来间接检验PALSAR数据集探测结果的有效性。

图3(b) 为采用17景ENVISAT ASAR数据处理得到的形变速率分布图。对比 图3(a) 和 图3(b) ，可以发现两个数据集探测结果在几处主要变形区的空间分布上大致能够相互对应，如 图3(b) 中标注的达拉仲村、草地贡村和白江岗村等区域都发现了较明显的变形趋势，但形变量相对于PALSAR探测结果存在较大差异。产生差异的原因是多方面的，包括影像观测获取时间不一致，雷达观测方向和视角不同，且PALSAR数据集干涉成像质量明显优于ASAR数据。在降轨观测获取的ASAR数据集中，角龙桥村位于叠掩区，无法有效提取高相干点目标，因此未能探测到该疑似滑坡。

对比分析两种数据集的探测结果，可以得出初步结论：从初步定性分析结果来看，利用时序InSAR技术可以探测到一些大型的滑坡体，在海通沟地区探测出的这几个地点为高度疑似的潜在滑坡体。PALSAR数据集相比ASAR数据集能够提取更高密度的点状目标，从而可获得细节信息更丰富的形变监测结果。

如前所述，由于本研究实验区缺少水准或GPS测量数据，因此无法利用地面测量资料对InSAR形变探测结果进行定量化验证。另一方面，由于两个数据集在成像几何等方面存在巨大差异，从中提取到的视线向形变速率很难直接进行定量化比较。

同时，实验区地理位置偏远，地形复杂，交通条件较差，开展实地勘测非常困难，而且由于所用的SAR数据是7年以前观测获取的，即使现在能到达现场考察，时效性也稍显不足。为此，只能尝试收集与SAR数据同时期获取的、覆盖实验区的光学遥感影像，通过目视判读解译，来对InSAR探测结果进行有效性验证分析。

首先选取了Google Earth使用的由DigitalGlobe提供的唯一一幅覆盖白江岗区域的高分辨率光学卫星影像( 图6 )。该影像观测获取时间为2011年3月2日。 图6 中红色虚线椭圆圈出的不规则白色线条，根据我们以往研究获得的解译经验，实际上是滑坡体变形错动后形成的后缘拉裂缝。与 图4(c) 所示形变场对照，不难发现拉裂缝正好位于形变速率大的区域边缘，二者吻合程度很高，充分说明PALSAR数据时序干涉处理提取的形变信息是有效的。

除白江岗外，对草地贡片区也进行了验证分析，选取的影像包括一幅Google Earth使用的米级分辨率DigitalGlobe光学卫星影像和一幅亚米级高分辨率真彩色航空相片，分别如 图7(a) 和 (b) 所示。前者观测成像时间为2009年10月26日，后者则是由北京洛斯达数字遥感技术有限公司于2013年7月飞行获取的。 图7(a) 中红色虚线椭圆圈出的很明显是一条规模较大的冲沟，在 图7(b) 的航片中能够看到更多细节信息。综合色彩、形状、纹理等特征，判断这条沟两侧地表存在强烈的风化剥蚀作用过程，很容易产生变形破坏。经与 图4(d) 对比发现，形变速率较大的相干点目标主要就集中在这条沟附近，二者在空间分布上基本吻合，再次说明了L波段InSAR形变探测结果的有效性。

已有研究表明，滑坡灾害的发生主要受到降雨等外部诱发因素的影响。针对前述从雷达数据中探测到的海通沟地区疑似滑坡，利用该地区历史气象资料，探讨降雨与滑坡形变之间在时间序列上的相互关联。

降水数据来源于芒康气象站(站点编号56227)记录存档的2007年—2010年累月降雨量数据集，下载于中国气象数据网( www.cma.gov.cn [2018-02-27])。 图8 为芒康县降水量统计分布图，从 图8 中可以看出每年的6月—9月是当地降水强度最大的季节，这4个月份月平均降水量达1000 mm。这种降水量分布的季节性差异主要是因为芒康属于高原温带半湿润季风气候区，夏季气候温和湿润，冬季寒冷干燥。

具体做法是从PALSAR数据探测到的疑似滑坡区域各选取一个高相干点目标作为代表，将提取到的达拉仲、角龙桥和白江岗3个滑坡体上的雷达视线向形变序列绘制成折线图( 图9 )。 图4(a) — (c) 中的白色十字星标出了这3个点的位置。

由于PALSAR数据观测获取频度较低，难以开展降雨量与形变量之间的定量关联分析，只能从时间演化趋势上进行定性解读判断。从 图9 中黑色虚线框标出的折线下滑区段来看，2008年—2010年期间每年6月—9月之间形变量明显增大，而当地的降雨也多集中在这几个月份，在时间上是吻合的。

滑坡的发生和发育与滑坡体的地质地貌和降雨的作用存在密切的关系，一般滑坡灾害发生具有随机性、季节性和地区性。芒康地区的地质背景偏软质裸土、植被覆盖较少，地形坡度大，受降雨影响以及河流下蚀共同作用，该地区容易遭受季节性降雨引起的滑坡灾害。在雨季，软质裸土受到冲蚀，降水在径流过程中下渗，会显著增加坡体含水量，从而加剧坡体失稳下滑。一旦滑坡体被部分破坏，坡面的水分继续下渗到裂缝和空洞中，更易导致二次滑坡的发生。

本文采用SBAS技术，从时间序列PALSAR数据集中探测藏东芒康县海通沟地区典型滑坡灾害，成功探测到了该地区318国道沿线的几处高度疑似滑坡体，提取了滑坡体表面形变特征，与同期ASAR数据集的形变探测结果进行了交叉比对，并通过对同时期高分辨率光学遥感影像的目视判读解译初步验证了InSAR探测滑坡形变结果的有效性。实验结果表明，长波长的L波段PALSAR数据在藏东复杂地区的滑坡探测中具有很好的应用潜力。

滑坡形变探测结果表明，在芒康海通沟地区存在几处形变异常的区域，如达拉仲、角龙桥等滑坡。据已有的资料表明，该地区曾经发生过滑坡事件，严重影响了318国道线的安全畅通。本论文发现降雨与滑坡形变演化趋势存在较大的关联性，为进一步研究滑坡预警和防范提供重要依据。

从本文实验研究的结果来看，在高频发地质灾害的藏东地区开展以时序InSAR方法为主要手段的潜在滑坡体探测是可行有效的，可为该地区电网建设运行、公路的稳定性维护及地质灾害防治部门提供重要支撑信息。相比传统方法，时序InSAR技术可以实现大范围探测并提取高密度的点目标，用于更好地发现和识别潜在滑坡体，在滑坡灾害早期识别和监测方面具有重要的应用价值。

Berardino P, Fornaro G, Lanari R and Sansosti E. 2002. A new algorithm for surface deformation monitoring based on small baseline differential SAR interferograms. IEEE Transactions on Geoscience and Remote Sensing, 40 (11): 2375–2383. [DOI: 10.1109/TGRS.2002.803792] Chen C J, Xue C S and Yu S H. 2004. Identification of geological hazards with remote sensing in the Motuo, Tibet. Journal of Engineering Geology, 12 (1): 57–62. [DOI: 10.3969/j.issn.1004-9665.2004.01.011] ( 陈楚江, 薛重生, 余邵准. 2004. 西藏墨脱公路的灾害地质遥感识别. 工程地质学报, 12 (1): 57–62. [DOI: 10.3969/j.issn.1004-9665.2004.01.011] ) Ferretti A, Prati C and Rocca F. 2000. Nonlinear subsidence rate estimation using permanent scatterers in differential SAR interferometry. IEEE Transactions on Geoscience and Remote Sensing, 38 (5): 2202–2212. [DOI: 10.1109/36.868878] Hooper A, Zebker H, Segall P and Kampes B. 2004. A new method for measuring deformation on volcanoes and other natural terrains using InSAR persistent scatterers. Geophysical Research Letters, 31 (23): L23611 [DOI: 10.1029/2004GL021737] Li Y, Meng H, Dong Y and Hu S E. 2004. Main types and characterisitics of geo-hazard in China-based on the results of geo-hazard survey in 290 counties. The Chinese Journal of Geological Hazard and Control, 15 (2): 29–34. [DOI: 10.3969/j.issn.1003-8035.2004.02.005] ( 李媛, 孟晖, 董颖, 胡树娥. 2004. 中国地质灾害类型及其特征——基于全国县市地质灾害调查成果分析. 中国地质灾害与防治学报, 15 (2): 29–34. [DOI: 10.3969/j.issn.1003-8035.2004.02.005] ) Liao M S and Wang T. 2014. Time Series InSAR Technology and Application. Beijing: Science Press (廖明生, 王腾. 2014. 时间序列InSAR技术与应用. 北京: 科学出版社) Shi X G. 2016. Monitoring Method and Application of Landslide Surface Deformation Based on Spaceborne Radar Remote Sensing. Wuhan: Wuhan University (史绪国. 2016. 基于星载雷达遥感的滑坡表面形变监测方法与应用. 武汉: 武汉大学) Wang P Q and Li P M. 2002. Analysis of geological disasters in southeast Tibet. Advances in Science and Technology of Water Resources, 22 (4): 21–22, 62. [DOI: 10.3880/j.issn.1006-7647.2002.04.008] ( 王培清, 黎普明. 2002. 藏东南地区地质灾害浅析. 水利水电科技进展, 22 (4): 21–22, 62. [DOI: 10.3880/j.issn.1006-7647.2002.04.008] ) Yang Z F, Shang Y G, Zhang L Q and Xu B. 2006. Study on the Geological Disaster of Sichuan-Tibet Highway and Its Prevention and Control Countermeasures: A Case Study of the Basu-Nyingchi Section. Beijing: Science Press (杨志法, 尚彦军, 张路青, 许兵. 2006. 川藏公路地质灾害及其防治对策研究: 以八宿至林芝路段为例. 北京: 科学出版社) Zhang M H. 2004. A remote sensing survey of geological disasters of Motuo County highway in Tibet. Highway (5): 91–96. [DOI: 10.3969/j.issn.0451-0712.2004.05.025] ( 张明华. 2004. 西藏墨脱公路地质灾害遥感勘察. 公路 (5): 91–96. [DOI: 10.3969/j.issn.0451-0712.2004.05.025] ) 基金项目: 国家重点研发计划(编号：2017YFB0502700)；国家自然科学基金(编号：61331016，41774006)；国家电网公司科技项目资助(编号：SG XN 0000 AQ JS 16 00002)；ALOS RA项目(编号：1247，1440，3248)；“龙”计划项目(编号：32778) Supported by: National Key R&D Program of China (No. 2017YFB0502700);National Natural Science Foundation of China (No. 61331016, 41774006); State Grid Science & Technology Research Program (No. SG XN 0000 AQ JS 16 00002) 张亚迪, 李煜东, 董杰, 范强, 车彬, 张路, 周杨, 廖明生. 2019. 时序InSAR技术探测芒康地区滑坡灾害隐患.遥感学报, 23(5): 987-996
Zhang Y D, Li Y D, Dong J, Fan Q, Che B, Zhang L, Zhou Y and Liao M S. 2019. Landslide hazard detection in Markam with time-series InSAR analyses.Journal of Remote sensing, 23(5): 987-996