基于InSAR技術的瀾滄江卡貢鄉(xiāng)—如美鎮(zhèn)段崩滑隱患探測

陳俊伊， 李為樂， 陸會燕， 李鵬飛， 鐵永波

(1.成都理工大學地質(zhì)災害防治與地質(zhì)環(huán)境保護國家重點實驗室，四川 成都 610059；2.中國電建集團貴陽勘測設計研究院有限公司，貴州 貴陽 550081；3.中國地質(zhì)調(diào)查局成都地質(zhì)調(diào)查中心，四川 成都 610081)

0 引言

水力發(fā)電是我國能源供給的重要組成部分，是我國水電工程建設的主戰(zhàn)場，主要集中在西部山區(qū)，特別是青藏高原東南部橫斷山脈的三江地區(qū)。該區(qū)域具有地勢高差大、河谷切割深、地質(zhì)構(gòu)造復雜、地震頻發(fā)等特點，導致該地區(qū)斜坡穩(wěn)定性差，地質(zhì)災害分布廣泛，嚴重影響了水電工程的施工建設和運營安全[1-2]。例如： 2018年10月11日和11月3日金沙江上游白格滑坡發(fā)生兩次失穩(wěn)，形成了滑坡-堰塞湖-潰決洪水災害鏈，不僅威脅了下游居民的生命財產(chǎn)安全，還影響了川藏鐵路及下游水電站的規(guī)劃建設[3]。瀾滄江豐富的水力資源對于我國的水電工程有著非常重要的意義，但復雜的地質(zhì)環(huán)境也為水電工程的建設增添了難度，為了避免出現(xiàn)滑坡等災害威脅工程建設，地質(zhì)災害的早期識別尤為重要。

InSAR(interferometric synthetic aperture rader)技術通過分析干涉信息來獲取地表形變信息，可以在地質(zhì)災害識別工作中常規(guī)地面調(diào)查手段難以快速大面積實施的區(qū)域進行崩滑隱患識別[4-5]。InSAR技術對地表形變的探測方法主要有雷達差分干涉測量(D-InSAR)技術及時序InSAR技術等。時序InSAR技術基于雷達影像序列探測地表形變的時空演化規(guī)律，該技術使InSAR技術對地表形變探測的精度達到亞厘米級[6]。時序InSAR技術包括永久性散射體合成孔徑干涉測量(PS-InSAR)[7]、小基線集干涉測量(SBAS-InSAR)[8]以及干涉堆疊法(Stacking-InSAR)[9]等。PS-InSAR技術需要通過選取永久散射體(permanent scatters,PS點)提高相干性，減少空間失相干對結(jié)果的影響，要求PS點在長時間序列上保持穩(wěn)定散射[10]，在山區(qū)環(huán)境PS點目標主要為裸露的巖土體，其分布并不廣泛，因而利用PS-InSAR技術進行廣域崩滑隱患識別效果并不理想。相對來說，Stacking-InSAR和SBAS-InSAR技術主要利用在山區(qū)環(huán)境廣泛分布的分布式目標(distributed scatters,DS點)進行干涉，在廣域崩滑隱患識別和監(jiān)測應用中效果較好[1,11-13]。因此本文選用SBAS-InSAR技術和Stacking-InSAR技術對覆蓋研究區(qū)的SAR影像數(shù)據(jù)進行處理，對瀾滄江卡貢鄉(xiāng)—如美鎮(zhèn)段進行崩滑隱患識別和重點區(qū)域時序形變監(jiān)測，為下游水電工程的建設提供保障，也為研究具有較強隱蔽性的高位、高危崩滑災害提供參考。

1 研究區(qū)概況

本文主要以瀾滄江干流卡貢鄉(xiāng)—如美鎮(zhèn)段及干流兩側(cè)5 km的范圍內(nèi)為研究區(qū)。瀾滄江發(fā)源于青藏高原唐古拉山脈，位于我國西南部[14]。研究區(qū)卡貢鄉(xiāng)—如美鎮(zhèn)段為瀾滄江中游，位于青藏高原東南部西藏自治區(qū)昌都市境內(nèi)，屬于橫斷山系，此處因金沙江、瀾滄江、怒江自北向南并行流淌，也稱為三江并流區(qū)[15](圖1)。

研究區(qū)屬于高山峽谷區(qū)，河谷深切于橫斷山脈，地形起伏大且岸坡高陡，局部甚至為陡壁[16]。區(qū)域內(nèi)的地形地貌由地質(zhì)構(gòu)造運動等控制，而構(gòu)造運動則會在區(qū)域內(nèi)形成活動斷裂。根據(jù)地震活動斷層探察數(shù)據(jù)中心提供的數(shù)據(jù)可知，研究區(qū)內(nèi)的活動斷層主要為發(fā)育于全新世的巴青—類烏齊斷裂和瀾滄江斷裂，其中瀾滄江斷裂貫通了研究區(qū)，且該斷裂活動較為強烈，曾多次發(fā)生4.7級以上地震，對研究區(qū)山體的穩(wěn)定性造成較大影響[17]。研究區(qū)按照區(qū)域地層劃分，主要屬于華南地層大區(qū)，羌北—昌都—思茅地層區(qū)，唐古拉—昌都地層分區(qū)。研究區(qū)沉積巖分布較廣，出露地層主要有泥盆系、石炭系、二疊系、三疊系、侏羅系、白堊系、古近系及新近系，第四系出露相對較少[18]。區(qū)域內(nèi)屬于典型的高原溫帶半濕潤季風性氣候，干濕季分明，6—9月集中降雨，其余時間降水量極少，降水峰值月的多年平均降水量比降水低值月的多年平均降水量多152.8 mm，降水量分配極為不均[19]。

圖1 研究區(qū)位置和雷達影像覆蓋范圍

2 數(shù)據(jù)與方法

本研究主要收集了覆蓋研究區(qū)的Sentinel-1雷達衛(wèi)星影像數(shù)據(jù)。利用GAMMA軟件通過Stacking-InSAR以及SBAS-InSAR兩種方法對收集的雷達影像數(shù)據(jù)進行處理，探測出研究區(qū)地表形變異常區(qū)，進而將探測出的地表形變區(qū)與2009年—2022年5月的Google Earth三維衛(wèi)星影像平臺(分辨率0.5～1.0 m)和天地圖二維衛(wèi)星影像平臺(分辨率1.0 m)多時相高分辨率光學遙感影像進行對比，利用崩滑隱患的形態(tài)特征和形變跡象，圈定研究區(qū)崩滑隱患點。

2.1 雷達衛(wèi)星數(shù)據(jù)收集

本次研究采用歐空局發(fā)布的Sentinel-1雷達衛(wèi)星影像數(shù)據(jù)[20]。本文選用的雷達衛(wèi)星數(shù)據(jù)為IW工作模式，通過漸進式掃描地形觀測獲取空間分辨率為5 m×20 m的3個子條帶，使其單景覆蓋范圍達250 km×170 km[21]。由于研究區(qū)為高山峽谷地區(qū)，為防止產(chǎn)生的幾何畸變加大形變探測的難度，甚至形成監(jiān)測盲區(qū)，采用聯(lián)合升軌、降軌探測結(jié)果對研究區(qū)進行解譯，提高隱患探測的準確性[22]。本次研究共收集了覆蓋研究區(qū)2018年1月5日—2021年6月30日的103景升軌影像數(shù)據(jù)，以及覆蓋研究區(qū)2018年1月12日—2021年6月25日的104景降軌影像數(shù)據(jù)，影像覆蓋范圍如圖1所示，同時收集了30 m分辨率的SRTM-DEM數(shù)據(jù)以提供地形地貌特征參數(shù)。

2.2 InSAR處理方法

為避免大氣效應等因素對解譯結(jié)果造成影響[13]，同時相對快速地獲取地表形變信息，首先通過Stacking-InSAR方法進行大范圍形變探測，并獲取隱患點的空間形態(tài)特征以及年均形變速率，再通過SBAS-InSAR抑制大氣效應以及DEM誤差的干擾對典型隱患點進行詳細形變探測，獲取隱患點在時間序列上的形變特征。

2.2.1 基于Stacking-InSAR的形變探測方法

由于地表形變相位與時間成線性變化，而大氣效應的影響與時間無關，因此通過對多幅差分干涉圖所獲取的解纏相位進行加權(quán)平均，可以得到更高精度的地表形變相位。具體步驟為： ①通過空間基線與時間基線設置，將配準后的SAR影像進行差分干涉，獲得N個干涉對； ②對干涉對進行濾波、解纏； ③選取解纏相位圖較為理想、解纏結(jié)果在工作區(qū)無明顯跳變的組合，加權(quán)求取平均值，估算線性相位速率，計算線性相位速率的公式為[23]

。

(1)

式中：Vφ為線性形變對應的相位變化速率；N為參與Stacking計算的干涉對數(shù)；ΔTi為第i個干涉對的成像時間間隔，d；pi為第i個差分干涉圖的解纏相位。

2.2.2 基于SBAS-InSAR的形變探測方法

Stacking-InSAR技術的計算過程較為簡單，能夠快速篩查出地災隱患點，但不能反映隱患點在時間序列上的形變特征。因此，為了更好地了解其形變過程，需要對典型形變點進行SBAS-InSAR處理。

。

(2)

利用外部DEM模擬地形相位，去除地形相位和參考橢球面引起的平地相位，得到m幅差分干涉圖，對干涉圖進行空間濾波后選取相干性較高的目標點進行相位解纏。根據(jù)時間基線、空間基線與干涉相位估算線性形變速率以及殘余相位，其中殘余相位中主要包括大氣相位、非線性形變相位以及噪聲； 對殘余相位以空域低通、時域高通的方式進行濾波，分離出大氣相位并將其去除，獲取非線性形變相位及噪聲； 通過多視及濾波等方法去除噪聲，得到非線性形變相位，將其與線性形變相加得到完整的形變相位，獲取累積形變量[11]。

3 Stacking-InSAR形變探測結(jié)果

Stacking-InSAR技術通過設置時空基線閾值形成小基線對，提高探測結(jié)果的準確性。本文設置升軌影像數(shù)據(jù)的最大時間間隔閾值為72 d，最大空間基線閾值為186.2 m，共生成500對干涉對，通過對解纏結(jié)果進行篩選，篩除141對有明顯相位跳躍及大范圍大氣噪聲的干涉對。設置降軌影像數(shù)據(jù)的最大時間間隔閾值為84 d，最大空間基線閾值為199.5 m，共生成505對干涉對，通過對解纏結(jié)果進行篩選，篩除118對有明顯相位跳躍及大氣噪聲的干涉對。分別對得到的升軌、降軌Stacking-InSAR結(jié)果進行解譯，并將解譯結(jié)果與Google Earth衛(wèi)星影像上獲取的高精度多時序光學影像圖進行對比，提高解譯結(jié)果的精確性并判別隱患點類型。滑坡隱患在光學影像中的主要解譯標志為： 后緣出現(xiàn)高亮裂縫，側(cè)緣邊界出現(xiàn)剪切裂縫，且有擴大貫通趨勢，部分坡體上出現(xiàn)新鮮下錯臺坎，局部有垮塌現(xiàn)象，整體形態(tài)多呈圈椅狀、梨狀、舌狀，有時呈不規(guī)則形狀[24-25](圖2)。崩滑隱患一般發(fā)育于地形陡峭的山崖或斜坡上部，巖體節(jié)理裂隙較為發(fā)育[26](圖3)。

(a) Stacking-InSAR年均形變相位 (b) 光學影像

(a) Stacking-InSAR年均形變相位 (b) 光學影像

研究區(qū)在觀測期內(nèi)升軌、降軌影像共探測出明顯隱患點共149處(圖4、圖5)，其中升軌影像共探測出隱患點68處，降軌影像共探測出隱患點95處，升軌、降軌探測出的相同隱患點14處，且隱患點多分布于研究區(qū)東岸。瀾滄江卡貢鄉(xiāng)—如美鎮(zhèn)段沿江共探測出隱患點65處，這些隱患點皆位于瀾滄江干流的陡壁上，范圍廣，形變量較大，具有堵江風險，嚴重威脅了沿線水電站的施工建設和運營安全。

圖4 Sentinel-1升軌影像重點區(qū)Stacking-InSAR年均形變相位

圖5 Sentinel-1降軌影像重點區(qū)Stacking-InSAR年均形變相位

4 典型隱患點形變分析

為了進一步分析研究區(qū)隱患點的形變特征，本文對研究區(qū)內(nèi)的隱患點進行篩選，選取其中3處較為典型的隱患點進行SBAS-InSAR處理，獲取隱患點在時間序列上的形變信息。

4.1 YH1滑坡隱患點

YH1滑坡隱患位于瀾滄江干流左岸，整體呈長舌狀，為古滑坡堆積體復活。隱患點長約3 337 m，寬約1 827 m，面積約2.9 km2，隱患點底部海拔3 197 m，頂部海拔4 302 m。隱患點坡體后緣有明顯下錯跡象，坡體上無植被覆蓋，沖溝遍布，坡體上多處可見崩塌，受風化侵蝕較嚴重，巖體質(zhì)量較低，為泥石流發(fā)育提供了充足的物源(圖6(a))。根據(jù)SBAS-InSAR探測結(jié)果(圖6(b))，YH1滑坡隱患點視線向年均形變速率最大達到30.2 mm/a，如圖6(c)所示。在隱患點內(nèi)選取3處特征點P1、P2及P3繪制時間序列曲線(圖7)，3處特征點的形變趨勢相同。在監(jiān)測期內(nèi)該隱患處于勻速形變階段，最大累積形變量已達到102.6 mm，表明該滑坡隱患形變較為強烈，整體不穩(wěn)定，且該隱患點范圍大，具有極大堵江風險，影響了坡體上的房屋安全，同時對瀾滄江下游水電站工程建設具有較大威脅。

(a) 光學影像(b) Stacking-InSAR年均形變相位 (c) InSAR年均形變速率

圖7 YH1滑坡隱患特征點時間序列曲線

4.2 YH2滑坡隱患點

YH2滑坡隱患點位于瀾滄江右岸，整體呈圈椅狀，坡體巖層屬于反傾巖層，其變形破壞模式為傾倒變形，隱患點長約1 174 m，寬約1 285 m，底部海拔3 113 m，頂部海拔3 744 m，面積約1.3 km2。光學影像上可見隱患點坡體后緣有下錯跡象，坡體上沖溝遍布，無植被覆蓋，坡體受風化侵蝕較嚴重，前緣多處發(fā)生崩塌(圖8(a))。選取隱患點內(nèi)的3處特征點P1、P2及P3繪制時間序列曲線(圖9)，時間序列曲線顯示3處特征點均呈勻速變形階段，在曲線上某些觀測值出現(xiàn)小幅度跳動，可能由大氣效應等因素形成的噪聲未徹底清除引起。根據(jù)SBAS-InSAR結(jié)果可知(圖8(b))，該隱患點形變速率較大，視線向最大年均形變速率達54.8 mm/a(圖8(c))，在觀測時間內(nèi)特征點的最大累積形變量達127.1 mm，可見該隱患點形變較為強烈，且形變范圍較大，有堵江風險，對道路、對岸居民及下游水電站建設具有較大威脅。

(a) 光學影像(b) Stacking-InSAR年均形變相位(c) InSAR年均形變速率

圖9 YH2滑坡隱患特征點時間序列曲線

4.3 YH3滑坡隱患點

YH3滑坡隱患點位于瀾滄江左岸，寬度約813 m，長度約1 114 m，面積約0.8 km2，隱患點底部海拔2 828 m，頂部海拔3 699 m，地形起伏大。隱患點整體呈圈椅狀，坡體上有沖溝，坡體局部有淺層崩塌，坡腳受江水掏蝕較劇烈，致使坡腳較陡(圖10(a))。根據(jù)SBAS-InSAR技術對隱患點進行形變探測結(jié)果(圖10(b))，坡體最大視線向年均形變量達41.7 mm/a(圖10(c))，圖中紅色點代表遠離衛(wèi)星視線方向運動。利用SBAS-InSAR技術在坡體上探測出一處強烈形變區(qū)，與Stacking-InSAR探測結(jié)果一致(圖10(b))，強烈形變區(qū)最寬處約448 m，長度約為532 m，面積約為0.2 km2。選取強烈形變區(qū)中的3處特征點P1、P2及P3繪制時間序列曲線(圖11)，時間序列結(jié)果顯示，3處特征點的形變趨勢相同，整體呈勻速變形，在監(jiān)測時間內(nèi)最大累積形變量達135.0 mm。對比所收集的光學影像，強烈形變區(qū)后緣出現(xiàn)裂縫，區(qū)域內(nèi)有局部淺層崩塌，證明該區(qū)域處于不穩(wěn)定狀態(tài)。結(jié)合光學及InSAR解譯結(jié)果，可以認定該隱患點具有較大的堵江風險。

(a) 光學影像 (b) Stacking-InSAR年均形變相位 (c) InSAR年均形變速率

圖11 YH3滑坡隱患特征點時間序列曲線

5 結(jié)論與現(xiàn)存問題

5.1 結(jié)論

本文利用Stacking-InSAR技術對覆蓋瀾滄江卡貢鄉(xiāng)至如美鎮(zhèn)段的Sentinel-1升軌、降軌影像進行處理，結(jié)合光學影像及Stacking-InSAR結(jié)果解譯出研究區(qū)范圍內(nèi)的崩滑隱患點。通過SBAS-InSAR技術對研究區(qū)范圍內(nèi)的3處典型滑坡隱患點進行處理，獲取其時間序列曲線，進一步分析隱患點的形變特征。綜合Stacking-InSAR技術、SBAS-InSAR技術及光學影像的結(jié)果，主要獲得以下幾點認識:

(1)研究區(qū)內(nèi)升軌、降軌影像共探測出隱患點149處，其中升軌影像共探測出隱患點68處，降軌影像共探測出隱患點95處。

(2)研究區(qū)遙感影像探測出的隱患類型主要為崩塌和滑坡，崩塌隱患共54處，滑坡隱患共95處。

(3)研究區(qū)滑坡隱患點多處于勻速變形階段，坡體多有明顯變形跡象，坡體巖體較為破碎，多有局部崩塌現(xiàn)象，視線向年均形變速率與累積形變量較大，坡體較不穩(wěn)定，且面積較大，有堵江風險，對瀾滄江及下游水電站建設威脅較大。

5.2 現(xiàn)存問題

通過Stacking-InSAR技術及SBAS-InSAR技術可以較為快速準確地探測出研究區(qū)崩滑隱患點的數(shù)量及形態(tài)，但由于衛(wèi)星飛行方向固定及側(cè)視成像的特點，會不可避免地產(chǎn)生幾何畸變，尤其是在高山峽谷地區(qū)。研究區(qū)處于高山峽谷地貌，地形起伏較大，本文所解譯的隱患點結(jié)果可能由于幾何畸變等因素的影響出現(xiàn)遺漏，同時大氣效應等噪聲未能徹底清除也會對探測結(jié)果產(chǎn)生影響。因此，為了能確保探測結(jié)果的準確性，還需要后期通過機載LiDAR探測以及野外勘察等手段進行驗證。