許東麗， 王 濤， 達瓦桑杰， 楊愛平， 李金鎖， 朱建東

(1.四川地質礦產(chǎn)勘查開發(fā)局207地質隊， 樂山 614000； 2.四川嘉源蓉創(chuàng)地質科技有限公司， 成都 610000；3.西藏自治區(qū)地質礦產(chǎn)勘查開發(fā)局第二地質大隊, 拉薩 850003)

中國是滑坡災害頻發(fā)國家，尤其是中國西部地區(qū)的滑坡，往往規(guī)模大、機制復雜，危害大[1]?；略谛纬砂l(fā)展過程中，往往伴隨著地表形變，通過對這些形變的監(jiān)測是識別及預測滑坡的重要手段。目前，滑坡監(jiān)測的主要方法有精密水準測量、全球定位系統(tǒng)(global positioning system，GPS)、位移監(jiān)測等，但這些方法監(jiān)測覆蓋范圍小、工作量大，不適用于大型滑坡監(jiān)測。近年來，隨著合成孔徑雷達(interferometric synthetic aperture radar,InSAR)技術的發(fā)展，研究者們拓展了其應用范圍，將原本主要用于地面沉降監(jiān)測的技術——星載InSAR用于滑坡監(jiān)測，取得了一系列良好的結果。

在InSAR變形監(jiān)測中，早期應用的是D-InSAR[2-5](differential InSAR)技術，在D-InSAR基礎上發(fā)展了長時間序列InSAR技術并運用在滑坡監(jiān)測中，主要包括永久散射技術PS-InSAR(permanent scatterer InSAR)和小基線集技術(SBAS-InSAR)。其中，PS-InSAR適用于城市等干涉條件和輻射比較穩(wěn)定的區(qū)域[6]，牛全福等[7]使用PS-InSAR技術對蘭州市地面形變進行監(jiān)測；而SBAS-InSAR是一種基于多主影像的InSAR時序方法，只利用時空基線較短的干涉對來提取地表形變信息，對數(shù)據(jù)的數(shù)量要求較低，適用于山區(qū)等干涉條件不穩(wěn)定的地區(qū)。因此，目前主要將其運用于非城區(qū)的線狀工程以及地質災害(滑坡為主)的形變監(jiān)測。文獻[8-11]使用SBAS-InSAR技術對川東地區(qū)、金沙江流域滑坡和鮮水河地質災害進行識別和監(jiān)測;文獻[12-13]分別使用SBAS-InSAR技術對地鐵及高速公路進行沉降監(jiān)測，以上研究均證明了SBAS-InSAR在低相干區(qū)也具有準確估計形變的能力。

綜上所述，中外眾多學者采用InSAR技術對不同類型工程和災害的形變特征進行研究，取得了較滿意的研究成果。但將InSAR技術應用于青藏高原特大型滑坡形變特征的研究相對較少。因此以規(guī)模大、植被覆蓋率低、處于變形加速階段的西藏自治區(qū)山南市龐村特大型滑坡為研究對象，采用SBAS-InSAR技術，選取2018年1月—2020年7月共75景Sentinel-1A數(shù)據(jù)，對滑坡區(qū)域的形變速率、累計形變量、時序形變特征及空間形變特征進行分析，并與地面監(jiān)測站實測的位移數(shù)據(jù)進行精度對比，驗證了SBAS-InSAR技術應用于青藏高原特大型滑坡形變監(jiān)測的可靠性，以期為青藏高原特大型滑坡的監(jiān)測工作提供有益參考。

1 滑坡特征與數(shù)據(jù)源

1.1 滑坡基本特征

龐村特大型滑坡位于西藏自治區(qū)山南市隆子縣的加玉鄉(xiāng)東側，加玉河右岸(圖1)。該滑坡為一古滑坡，滑坡整體地形呈南高北低，上陡下緩，坡度20°～40°，局部地段可達50°以上，平均坡度約23°，整體面積約0.3 km2，整體體積約1.1×107m2，屬于特大型滑坡?；缕矫嫘螒B(tài)呈圈椅狀，剖面形態(tài)呈“階梯”形，平面形態(tài)呈舌形，主滑方向約14°?；麦w后緣拉張裂縫已基本貫通并形成下錯陡坎，基巖出露，最大下錯達5.2 m，滑坡后壁可見明顯擦痕；滑坡前緣剪出口因G219改道開挖，出現(xiàn)不同程度垮塌(圖2)。當?shù)卮迕裼?019年4月發(fā)現(xiàn)滑坡裂縫，隨著雨季和地震的影響，裂縫持續(xù)增大，滑坡整體處于不穩(wěn)定狀態(tài)。

圖1 研究區(qū)域衛(wèi)星影像圖Fig.1 Satellite imagery of the study area

1.2 數(shù)據(jù)源

1.2.1 雷達數(shù)據(jù)

本文中數(shù)據(jù)源選用Sentinel-1A衛(wèi)星數(shù)據(jù)，Sentinel-1A衛(wèi)星是由歐空局發(fā)射的免費對地觀測雷達衛(wèi)星，軌道高度693 km，多極化，影像采用TOPS(terrain observation with progressive scans in azimuth)模式成像，重訪周期為12 d，雷達波段為C波段，波長5.6 cm。數(shù)據(jù)空間分辨率30 m，本文中所采用的數(shù)據(jù)時間跨度為2018年1月3日—2020年7月3日，總計75景影像，數(shù)據(jù)來源(https://qc.sentinel1.eo.esa.int/)。數(shù)據(jù)的具體時序信息如表1所示，覆蓋范圍如圖1所示。

表1 數(shù)據(jù)時間序列列表Table 1 Time series list of data

同時，由于研究區(qū)域所處區(qū)域地形起伏較大，為消除地形相位對形變相位的影響，本文中采用30 m分辨率的SRTM(shuttle radar topography mission)數(shù)據(jù)作為輔助數(shù)據(jù)，計算地形相位分量。

1.2.2 地面監(jiān)測數(shù)據(jù)

相關管理部門于2019年8月初在滑坡體上安裝4臺一體化GNSS監(jiān)測站以及3臺拉線式裂縫計，對滑坡的位移狀況進行監(jiān)測。本次搜集該地面監(jiān)測儀器從2019年8月1日—2020年7月3日的地面監(jiān)測數(shù)據(jù)，覆蓋InSAR數(shù)據(jù)的后半段，對SBAS-InSAR解譯結果進行精度驗證。地面監(jiān)測儀器安裝位置如圖3所示。

1.2.3 降雨量數(shù)據(jù)

本次實驗還搜集了滑坡周邊日降雨量數(shù)據(jù)，該數(shù)據(jù)的采集地點位于滑坡威脅區(qū)外，由雨量計測得。數(shù)據(jù)時間為2019年8月1日—2020年7月3日，同樣覆蓋InSAR數(shù)據(jù)的后半段。

2 數(shù)據(jù)處理

2.1 時序InSAR形變監(jiān)測原理

InSAR技術能夠獲取大面積、高精度的數(shù)字高程模型，可用來監(jiān)測地震、火山、地面沉降等造成的微小地表形變信息。用于形變分析時，InSAR技術可獲得監(jiān)測區(qū)域厘米級甚至毫米級的形變量，因此許多研究者將其用于滑坡的緩慢形變監(jiān)測，分析滑坡的滑動趨勢，為滑坡有效治理提供可靠的技術依據(jù)。

在滑坡監(jiān)測領域，最早使用的是D-InSAR技術，它獲取同一研究區(qū)域不同時間獲取的三幅或以上的影像進行干涉處理，獲得干涉像對，再對干涉像對進行差分處理，消除地形起伏和幾何畸變影響，最終獲得該時段的滑坡形變信息。其原理如下：設形變前衛(wèi)星在A1、A2處獲取影像，斜距為r1、r2,θ為衛(wèi)星在A1位置成像的側視角；形變發(fā)生后在A3位置獲取影像，得到斜距r3、空間基線B2,δ1、δ2為基線B1、B2的基線傾角，B‖、B⊥和B′‖、B′⊥分別為空間基線B1、B2沿雷達視線和垂直于雷達視線方向上的分量，ΔRd為地面點的形變量(圖4)。則D-InSAR形變監(jiān)測原理可表達為

圖4 D-InSAR原理Fig.4 Principle of D-InSAR

r1r2?B1sin(θ-δ1)?B‖

(1)

由式(1)可得，干涉測量相位差與視線方向基線分量成正比。

將在位置A1和A2處獲得的地表未形變前的影像進行干涉處理，相位僅包含地形分量，相位差為

(2)

式(2)中：λ為波長；P為重復軌道，P=2。

地表發(fā)生形變后在位置A3獲取影像，將其與未發(fā)生形變時在A1處獲取的影像進行干涉處理，得到相位包含地形和形變分量，相位差為

(3)

設視線方向的形變量為ΔRLOS,則由式(1)～式(3)，可得形變量引起的干涉圖相位差為

(4)

式(4)中：ΔφLOS可根據(jù)相位解纏獲得，則根據(jù)式(1)～式(4)可計算沿雷達視線方向的形變量ΔRLOS。

由式(4)可得，當ΔRLOS=λ/2時，ΔφLOS=2π，因此當形變量為半波長時，即可引起一個周期的相位變化，說明差分相位對地表形變較為敏感，可用于探測微小形變。

由D-InSAR形變監(jiān)測原理可知，該技術僅僅使用三期數(shù)據(jù)進行差分處理，適用于分析緩慢形變滑坡較短時間周期內的形變，當滑坡形變量較大和監(jiān)測周期較長時，使用D-InSAR技術會導致空間失相干和時間失相干，降低形變監(jiān)測精度。為了解決以上問題，研究者們提出了長時間序列InSAR技術，主要包括PS-InSAR技術和SBAS-InSAR技術，這兩種方法對時空失相干和大氣擾動兩方面進行限制，顯著提高了形變監(jiān)測能力。其中PS-InSAR利用時間相干原理，對影像數(shù)量要求高，適用于城區(qū)等獲取時間規(guī)律、連續(xù)或干涉條件和輻射比較穩(wěn)定的區(qū)域；SBAS-InSAR技術利用空間相干原理，對數(shù)據(jù)數(shù)量要求低。

由于本文中研究區(qū)域面積較大、地形起伏大、滑坡持續(xù)時間長，干涉條件不穩(wěn)定，因此選用SBAS-InSAR進行形變監(jiān)測和分析。

2.2 SBAS處理流程

SBAS技術是短基線差分干涉方法，將所有的SAR影像依據(jù)空間和時間基線分成不同的短基線子集，各子集內進行差分干涉，再對纏繞的干涉相位進行解纏，根據(jù)各相干像元的相位與觀測時間的關系，利用奇異值分解方法鏈接各差分干涉圖，抑制DEM(digital elevation model)和大氣相位對形變相位的影響，獲得最小二乘解，其處理流程見圖5。

圖5 SBAS-InSAR處理流程Fig.5 Data processing flow of SBAS-InSAR

3 研究結果與分析

3.1 地表形變速率與累計形變量

3.1.1 龐村滑坡及其周邊區(qū)域形變情況

搜集Sentinel-1A降軌數(shù)據(jù)75景，基于GAMMA軟件處理平臺，使用SBAS-InSAR數(shù)據(jù)處理方法，獲得了龐村滑坡及其周邊區(qū)域的雷達視線(line of sight, LOS)方向形變速率圖(圖6)，其中負值表示形變方向為沿雷達視線并遠離衛(wèi)星方向。在 2018年1月3日—2020年7月3日期間，龐村及其周邊區(qū)域最大形變速率約為-128.09 mm/年。

圖6 龐村滑坡及其周邊地區(qū)地表形變速率圖Fig.6 The surface deformation rate of Pangcun landslide and the surrounding area

3.1.2 地表形變速率與累計形變量

選取龐村滑坡處的解譯結果進行詳細分析，結果如圖7所示，在2018年1月3日—2020年7月3日期間，龐村滑坡區(qū)域的最大形變速率為-38.182 mm/年。從解譯結果中可以看出滑坡形變主要發(fā)生在右側；在形變區(qū)域存在明顯的形變分區(qū)，形變速率較大的地區(qū)主要位于滑坡右側的中下部，形變速率大于23.55 mm/年，其中形變速率最大的區(qū)域主要有三處，形變速率均大于33.340 mm/年(圖7)。

依據(jù)形變速率，對滑坡區(qū)域2018年1月3日—2020年7月3日的累計形變量進行統(tǒng)計分析，可得滑坡區(qū)域內最大累計形變量為110.393 mm(圖8)。

圖8 龐材滑坡地表累計形變量圖Fig.8 The cumulative displacement of Pangcun landslide

3.2 滑坡形變空間分布特征

根據(jù)形變速率以及累計形變量圖，可將滑坡體劃分為3個不同形變程度的分區(qū)，如圖9所示。

圖9 龐材滑坡地表滑坡形變分區(qū)Fig.9 The deformation zone of Pangcun landslide

為進一步研究滑坡體形變的空間分布特點，根據(jù)初步劃定的形變分區(qū)，沿滑坡主滑方向取兩條剖面線A—A′和B—B′(圖10),長度分別為895 m和876 m，從起點開始，每隔25 m取剖面點進行形變速率和累計形變量分析，結果如圖11所示。

圖10 剖面線分布圖Fig.10 Section line distribution

圖11 剖面各點形變速率與累計形變量Fig.11 Deformation rate and cumulative displacement at each point of the section line

在剖面線A—A′上，形變速率和累計形變量最大的點距離起點600 m，形變速率達26.036 mm/年，累計形變量73.54 mm；剖面線B—B′上，形變速率和累計形變量最大的點位于距起點500 m處，形變速率和累計形變量分別為28.786 mm/年，累計形變量80.086 mm。

由剖面點形變速率和累計形變量分析結果可知，剖面上形變最大的點位于滑坡前緣區(qū)域，與實地調查結果一致，該區(qū)域主要為G219國道改造，人工修路開挖造成邊坡失穩(wěn)，形成陡坎并下滑。

3.3 滑坡形變時間分布特征

根據(jù)滑坡形變分區(qū)結果，在重點形變區(qū)(Ⅳ、Ⅴ、Ⅵ表示)內選取兩個特征點P1、P2(圖12)進行時序分析，研究滑坡形變的時間變化特點，結果如圖13所示。

圖12 重點形變區(qū)及特征點分布Fig.12 Distribution of the main deformation areas and feature points

由圖13可以看出，滑坡在兩個時段存在明顯的形變加速，分別為2018年11月29日—2019年1月22日(Ⅰ期)和2019年10月1日—2019年11月18日(Ⅱ期)，其余時期為均勻沉降。

根據(jù)實地調查結果以及降雨量數(shù)據(jù)：Ⅰ期的形變加速主要是受2018年10月對G219國道的改造工作以及隆子縣2018年12月1日地震(28.12°N， 92.62°E，震源深度5 km，據(jù)龐村滑坡23 km)影響，國道G219改造工作對滑坡前緣邊坡開挖，造成邊坡失穩(wěn)，并在地震的作用下加速下滑。Ⅱ期形變加速主要受2019年7月19日錯那5.6級地震(27.67°N， 92.89°E，震源深度10 km，據(jù)龐村滑坡72 km)以及研究區(qū)域連續(xù)降雨(圖13)，降低了巖土體的穩(wěn)定性，滑坡體加速形變。

圖13 特征點累計形變量時序圖Fig.13 Time series diagram of cumulative displacement of feature points

3.4 精度驗證

精度驗證部分，將本文中SBAS-InSAR解譯結果與地面監(jiān)測結果進行對比，得到如下結果。

(1)SBAS與GNSS監(jiān)測結果。如圖14所示，將地面GNSS監(jiān)測站監(jiān)測結果與SBAS在相應區(qū)域的解譯結果進行對比發(fā)現(xiàn)：四個區(qū)域均處于持續(xù)形變階段，并且儀器GNSS-4處的累計形變量比其他三處小，SBAS解譯結果可得到相同結論。

圖14 SBAS-InSAR與GNSS結果對比Fig.14 Comparison of SBAS-InSAR and GNSS results

(2)SBAS與拉線式位移計監(jiān)測結果。如圖15所示，將地面拉線式位移計監(jiān)測結果與SBAS在相應區(qū)域的解譯結果進行對比發(fā)現(xiàn)：LF-1與LF-2兩個區(qū)域處于持續(xù)形變階段，LF-3在前期形變后處于穩(wěn)定階段，SBAS解譯結果可得到相同結論；同時結合監(jiān)測點位分布圖(圖3：LF-1LF-2處于滑坡后緣右側，LF-3處于滑坡后緣左側)，驗證了SBAS解譯結果以及根據(jù)解譯結果劃分的形變分區(qū)與地面拉線式位移監(jiān)測結果的一致性，即滑坡形變主要集中在滑坡右側。

圖15 SBAS-InSAR與拉線式位移計結果對比Fig.15 Comparison of SBAS-InSAR and linear displacement results

由兩種對比結果可以看出，地面監(jiān)測結果在數(shù)值上比SBAS解譯結果大，這是由于地面監(jiān)測的主要為滑坡滑動方向的形變量，而SBAS解譯結果為沿著雷達視線方向的形變量；并且由于Sentinel-1A數(shù)據(jù)的地面分辨率較低，區(qū)域內形變結果代表了該對應像素分辨率內地物的平均形變量值。但二者在總體形變趨勢、形變的空間分布以及時間分布特征上基本一致。

因此運用SBAS-InSAR對免費的Sentinel-1數(shù)據(jù)進行解譯處理，所得的滑坡形變結果可為滑坡形變監(jiān)測與治理提供可靠的技術依據(jù)；并且相對于地面監(jiān)測站來說，成本更低，結果較為準確，可在部分無地面監(jiān)測站點的大型滑坡監(jiān)測中，使用SBAS解譯結果代替地面監(jiān)測結果。

4 結論

以西藏自治區(qū)龐村滑坡為研究對象，采用SBAS-InSAR技術進行特大型滑坡形變監(jiān)測，獲得了研究區(qū)域2018年1月3日—2020年7月3日時間段內的地表形變速率與累計形變量圖，并對結果進行分析，得出了以下結論。

(1)在2018年1月3日—2020年7月3日期間，滑坡區(qū)域的最大形變速率為-38.182 mm/年?；滦巫冎饕l(fā)生在滑坡右側中下部。

(2)在2018年1月3日—2020年7月3日期間，滑坡區(qū)域內最大累計形變量為110.393 mm，滑坡區(qū)域存在明顯的分界線，可分為三個重點形變區(qū)。

(3)空間分布上，剖面線A—A′和B—B′的起點和終點均不存在形變，形變最大的區(qū)域分別位于距離剖面A—A′和B—B′起點600 m和500 m處。

(4)時間分布上，在2018年1月3日—2020年7月3日期間，滑坡表面存在兩個明顯的形變加速時段，分別為2018年12月29日—2019年1月22日和2019年10月1日—11月18日，其余時期為均勻沉降。

(5)精度驗證上，本文中SBAS-InSAR解譯結果與地面監(jiān)測所獲得的形變結果基本一致，驗證了使用免費的Sentinel-1A數(shù)據(jù)，基于SBAS-InSAR技術對大型滑坡進行形變監(jiān)測的有效性。

因此，使用免費的Sentinel-1A數(shù)據(jù)，基于SBAS-InSAR技術對西藏自治區(qū)龐村特大型滑坡進行形變監(jiān)測與分析，充分展示了該滑坡形變的時間和空間分布特征，為探究滑坡的滑動機理以及后續(xù)滑坡治理工程提供了可靠的技術依據(jù)；將形變監(jiān)測結果與地面調查以及地面監(jiān)測結果進行對比，驗證了解譯結果的可靠性，可將監(jiān)測結果進行由單點到面域的轉換；并且在部分無地面監(jiān)測站點的滑坡災害中，可將SBAS-InSAR監(jiān)測結果在一定程度上代替地面監(jiān)測結果，為后續(xù)類似的特大型滑坡監(jiān)測提供新的思路和參考。