基于多源遙感技術的紅層滑坡識別與監測研究

王 慶 芳，鄭 志 軍，董 繼 紅，余 天 彬，劉 文，黃 細 超

(1.成都市地質環境監測站，四川 成都 610041； 2.自然資源部成都地質災害野外科學觀測研究站，四川 成都 610041； 3.四川省地質調查院 稀有稀土戰略資源評價與利用四川省重點實驗室，四川 成都 610081；4.四川省智慧地質大數據有限公司，四川 成都 610081)

0 引 言

滑坡作為一種頻繁發生，破壞最為嚴重的自然災害之一，具有隱蔽性強、破壞性大、突發性高等特點，嚴重危害著人民群眾的生命財產安全，同時也制約著區域經濟的發展[1-3]。紅層是一種外觀以紅色為主色調的陸相碎屑巖沉積地層，龍泉山是典型的低山丘陵紅層區。地層主要是侏羅系、白堊系，巖性以紅色、紫紅色、磚紅色砂巖、粉砂巖、泥巖、礫巖為主[4-5]，該區域的植被覆蓋率較高，人類工程活動強烈。地層具有強度低、親水性強、易崩解、軟化、膨脹等特性，易發生規模較大的順層巖性滑坡[6-8]。據調查顯示，自2008年以來，該區域發生多起大規模的滑坡災害，如雷打石滑坡、五家墳滑坡、油榨房滑坡等。區內滑坡具有分布廣、隱蔽性強、突發性高、識別難度大等特點。因此很有必要利用遙感技術對該區域的滑坡隱患開展大范圍的早期識別研究，對典型的滑坡開展時序監測，分析變形規律及特征。

光學遙感技術因結果準確率高、解譯精度高、獲取信息豐富、可以回溯歷史形變等特點廣泛用于滑坡識別、解譯分析等工作[9-10]，但是該技術自動化程度低、工作量大。而合成孔徑雷達干涉測量(Interferometric Synthetic Aperture Radar，InSAR)技術因其非接觸、不受云霧天氣影響、空間覆蓋范圍廣、監測精度高等優勢，廣泛用于地表形變監測[11-13]。張路等[14]以四川省丹巴縣為例，利用InSAR技術獲取的形變速率圖及差分干涉圖成功識別出17處滑坡隱患，同時總結了影響時序InSAR技術監測效果的主要因素。戴可人等[15]通過SBAS技術獲取的形變速率圖進行隱患的識別圈定，并討論了研究區域內Sentinel-1數據可視性分布情況。朱賽楠等[16]利用光學遙感技術、InSAR技術結合物探及地表位移監測等多個手段對金沙江流域色拉滑坡的特征、變形過程、形變機理進行研究，并對其變形趨勢進行了預判。李媛茜等[17]基于InSAR技術結合野外資料成功圈定了白龍江流域114處滑坡隱患，分析了滑坡發育特征，并通過對牙豁口滑坡結合無人機數據進行分區研究，證實了InSAR技術在滑坡變形監測的有效性和可靠性。同時一些學者也利用Stacking InSAR技術及SBAS InSAR技術進行滑坡隱患的大范圍識別及監測，取得了不錯的結果[18-19]。

基于此，本文以四川省成都市龍泉山區域為實驗區域，收集了覆蓋該區域的ALOS PALSAR-2數據、Sentinel-1數據及高分光學數據，開展多源遙感技術在低山丘陵區紅層滑坡隱患的早期識別，并對比不同波段SAR數據的隱患識別情況；基于SBAS技術結合無人機數據及現場調查結果分析了典型滑坡的變形規律及特征，并給出了相應建議。

1 研究區概況及數據來源

1.1 研究區概況

龍泉山位于四川省成都市東側，西接成都平原區，東部接川東丘陵區，地勢中部較高沿龍泉山山脊展布，東西兩側低，本次研究區主要分布在龍泉山中段(見圖1)，研究區地貌類型主要為構造侵蝕低山地貌，屬于川東低山丘陵區。區內屬中亞熱帶濕潤氣候區，四季分明、濕潤溫暖，降雨豐沛，雨季多集中于6～9月，降雨量約占全年的75%。調查分析發現降雨條件對區域的滑坡災害具有較強的控制作用，滑坡發生的時間主要集中在強降雨期。

圖1 研究區位置Fig.1 Geographical location of the study area

1.2 數據收集與處理

本次研究采用多源遙感技術綜合識別，利用多源SAR數據和高分光學衛星影像開展區域滑坡隱患早期識別。收集了覆蓋該區域2017～2020年的ALOS PALSAR-2降軌數據。同時也收集了可以免費獲取的Sentinel-1數據進行研究區內典型滑坡隱患的形變規律監測分析(見表1)。利用AW3D30 DSM數據來消除InSAR干涉處理中的地形相位及輔助SAR影像進行地理編碼。Sentinel-1衛星的POD精密軌道星歷(POD Precise Orbit Ephemerides)數據被用來輔助Sentinel-1數據的預處理和基線誤差改正。

表1 選用的遙感數據信息

同時在數據處理中對與高程不相干的隨機大氣擾動誤差，利用李振洪教授團隊研制的通用型衛星雷達在線大氣改正系統(Generic Atmospheric Correction Online Service for InSAR，GACOS)來去除[20]。

光學衛星數據主要包括國產高分二號和北京二號數據，用來進行滑坡隱患的光學遙感解譯和識別驗證。

2 低山丘陵區滑坡識別監測方法

結合工作區滑坡發生的特殊性質，如分布廣泛、隱蔽性強、空間尺度較小、受降雨影響大以及突發性強等特點，在基于InSAR技術用于滑坡隱患識別的相關研究和前期數據處理結果分析基礎上，本文結合收集到的SAR數據及其他相關資料，針對低山紅層區突發性滑坡隱患識別、監測研究，給出了本文的滑坡隱患識別與監測處理流程(見圖2)。該流程的關鍵是對基于地表形變速率解譯的滑坡隱患從光學影像中尋找依據，通過SBAS技術判定其變形量及變形規律，通過野外調查、核查識別結果的準確性，并對災害庫進行更新，同時積累該區域的解譯標識庫。

圖2 滑坡隱患識別與監測流程Fig.2 Flow chart of landslide hidden danger identification and monitoring

2.1 Stacking InSAR技術用于滑坡隱患編目

目前開展大范圍滑坡隱患早期識別工作，普遍采用Stacking InSAR技術，該技術具有使用影像數量少、結果獲取效率高等優勢。因此本文在對比不同波段數據獲取結果解譯情況與已有滑坡隱患庫基礎上，建議在進行Stacking計算中采用ALOS-2數據與Sentinel-1數據相結合的形式進行。

其中Stacking技術[21]是由Sandwell等在1998年提出，該技術原理是對解纏相位進行加權平均解算，以達到削弱空間上不相關的噪聲的影響，從而獲取研究區域的平均形變速率[22-23]。主要處理步驟如圖3所示。

圖3 STACKING-INSAR技術流程Fig.3 Stacking-InSAR technology process

2.2 基于光學遙感技術的滑坡隱患識別

結合工作區地質災害特征，對滑坡隱患點的識別，可采用對比分析方法和直接識別方法進行判定。對比分析方法首先利用遙感影像對斜坡失穩前地質環境條件和變形跡象進行解譯，然后對工作區內地質環境條件相似或者變形跡象明顯的地質災害進行遙感解譯。直接識別方法的關鍵是建立遙感影像識別標志，遙感解譯標志包括了地質災害發育的地質環境條件和斜坡變形特征解譯兩個方面，如表2所列。

表2 滑坡光學遙感解譯識別標志

2.3 SBAS InSAR技術用于滑坡監測分析

對于典型滑坡隱患采用SBAS InSAR技術進行長時間序列監測，獲取其變形趨勢及變形量等信息，因為Sentinel-1數據時間分辨率高(時間間隔為6/12 d)、免費獲取，所以本文利用Sentinel-1數據基于SBAS技術開展典型滑坡隱患的時間序列監測分析。

SBAS技術是由Berardino等于2002年提出的，基本原理是假設在監測時間t0～ts內獲取了S+1幅SAR影像，根據干涉組合條件可形成M幅干涉組合，符合式(1)：

(1)

一般根據實際情況設置一定的垂直基線閾值和時間基線閾值對M的值域進行限制，通過這種形式可以降低因垂直基線和時間基線過大引起的失相干現象[23-24]。利用差分干涉計算獲得差分干涉相位，用于第j幅干涉相位是根據tA和tB(tA

3 結果分析

3.1 基于Stacking InSAR技術的滑坡隱患早期識別與效果分析

分別利用ALOS-2數據和Sentinel-1數據進行Stacking InSAR處理。對ALOS-2數據進行距離向和方位向采用2×4的多視處理來提高信噪比，由于ALOS-2數據量較少，因此采用全組合的形式進行差分干涉處理，并進行Stacking InSAR解算獲得沿雷達衛星視線向(Line of Sight，LOS)年平均形變速率結果(見圖4)。對于Sentinel-1數據，在處理中采用8×2進行距離向和方位向多視處理，設置空間基線不大于200 m，時間基線不大于48 d進行干涉對組合，進而提高相干性和運算效率，并通過Stacking InSAR技術獲取了形變速率結果。

圖4 利用ALOS-2數據獲取沿LOS向地表形變速率Fig.4 Surface deformation rate along the LOS direction obtained by ALOS-2 data

根據形變速率圖中的條紋形態判識疑似的地災隱患[29]，結合地形條件及光學影像，對基于不同數據源獲取的地表形變監測結果進行滑坡隱患的識別解譯，其中ALOS-2數據共識別出9處隱患點(見圖4)，野外核查發現8處有明顯變形，1處無明顯變形。Sentinel-1升軌數據識別出8處，經野外驗證5處為滑坡隱患，Sentinel-1降軌數據識別出6處，經野外驗證2處為滑坡隱患，Sentinel-1數據經驗證正確的隱患點與ALOS-2數據識別的隱患點位置一致，具體滑坡隱患信息如表3所列。

表3 滑坡隱患信息

選擇了3處經野外驗證為滑坡隱患進行多源遙感數據識別對比分析(見圖5)。其中獅子山滑坡和宋家灣滑坡為新增隱患點。通過基于ALOS-2數據獲取了這3處滑坡隱患的解譯信息，可以看到明顯的變形信息。相應地，在Sentinel-1升軌數據中可以看到較為明顯的形變信息。其中獅子山滑坡和鄢家溝滑坡的變形信息與ALOS-2數據獲取的變形信息具有較好的一致性，只是Sentinel-1是升軌數據，ALOS-2是降軌數據，所以在速率圖中圖斑顏色正好相反。這3處滑坡隱患在Sentinel-1降軌數據中沒有變形信息，故沒有展示。

宋家灣滑坡隱患是基于ALOS-2數據獲取的新增隱患點，該滑坡在Sentinel-1升軌數據中變形信息并不明顯，為此依據圖2流程進行驗證分析。圖6為該滑坡隱患對應的光學影像，從圖6中可以看到該滑坡隱患具備地質災害發育的地質環境條件，但是斜坡變形特征并不明顯。為此，提取基于InSAR監測變形強烈的部分，再基于SBAS InSAR技術提取其時間序列變形信息(見圖7)。從圖7可以看到：該滑坡一直處于變形狀態，變形具有明顯的周期性，每年7月前后處于加速變形狀態，變形時間為工作區雨季，同時2018年1月至2020年12月累積變形量達到了132 mm，佐證了該點為滑坡隱患點，后期通過野外調查驗證，確定該點為一新增滑坡隱患點。

圖6 宋家灣滑坡光學衛星影像解譯信息(北京2號)Fig.6 Optical satellite image interpretation information of Songjiawan landslide(satellite of Beijing 2)

圖7 宋家灣滑坡時序點變形信息Fig.7 Deformation information of Songjiawan landslide sequence point

3.2 典型滑坡隱患監測

選擇利用時序InSAR技術識別的一處新增隱患點——鄭家房子滑坡為例進行隱患監測分析。該滑坡地理位置為104°13′21.93301″E，30°23′0.09891″N；滑坡形態整體呈簸箕形，滑坡前緣高程685 m，后緣高程770 m，高差85 m，縱長約260 m，橫寬約150 m，坡向約245°，坡度15°～20°，滑坡平面面積為6.25萬 m2，按10 m厚度估算，估算方量為6.25×105m3，為一大型緩傾滑坡。地層巖性為侏羅系中統沙溪廟組(J2s)淺紫紅、黃灰色厚-塊狀中細粒巖屑砂巖、巖屑長石砂巖與紫紅色、灰紫色粉砂質泥巖、粉砂巖互層，滑坡結構為一順向坡，屬于典型的低山紅層區順層滑坡。

同時獲取了該區域的無人機正射影像圖(見圖8)，從光學影像圖中可以看到滑坡后緣呈圈椅狀，可見早期滑動形成的洼地緩坡平臺；滑坡兩側以負地形凹地為界，可見溝道微地貌；滑坡前緣以緩坡平臺為界，堆積體外凸擠向臨空面，陡緩交界線明顯。

圖8 鄭家房子滑坡光學遙感解譯(據2021年無人航測數據)Fig.8 Interpretation of Zhengjiafangzi landslide by optical remote sensing(according to unmanned aerial survey data in 2021)

通過對此滑坡進行現場核查，發現滑坡變形區域與監測識別吻合性較好。滑坡變形區主要集中在左側中部居民區處，房屋局部發生拉裂變形，坡體局部發生溜滑。

據訪問坡體上部右側邊界處的一戶居民，可知該區域斜坡在2018年至2021年底，每年降雨期坡表均有不同程度的溜滑，該居民房屋無變形，屋后堡坎局部發生溜滑。現場調查坡體中部發育裂縫，裂縫走向NE35°，裂縫延伸長度10～25 m，寬10～20 cm，錯落高度5～15 cm。該滑坡為一土質滑坡，斜坡坡向與滑動方向基本一致，屬于龍泉山典型的易滑斜坡區。

為了對鄭家房子滑坡的變形特征進行分析，采用SBAS InSAR技術對Sentinel-1升軌數據進行處理，獲得了該滑坡沿LOS向地表形變速率結果(見圖9)，正值表示滑坡位移靠近衛星運動方向變形，負值表示遠離衛星飛行方向。

圖9 鄭家房子滑坡形變速率Fig.9 Deformation rate of Zhengjiafangzi landslide

將圖9與ALOS-2數據結果對比后發現，兩者在強變形區域具有較好的一致性，圖9中形變速率最大達到了-40 mm/a。在滑坡體上選擇特征點提取其在時間域的累積形變量結果(見圖10)。可以看到特征點的變形信息具有明顯的規律性，每年的5～7月之間發生變形，變形速率最大，其他時間段基本趨于穩定，其變形時間與當地的雨季吻合。其中B、C點形變速率最大，累積變形量也最大，其在2018年1月至2020年12月累積變形量達到134 mm。結合現場調查情況可知該滑坡變形主要集中在每年的汛期，在強降雨作用下，容易發生大規模滑坡，威脅坡體居民區、村道，因此需要在雨季加強持續性監測。

圖10 鄭家房子滑坡體監測點累計形變曲線Fig.10 Cumulative deformation curve of monitoring points of Zhengjiafangzi landslide mass

4 討 論

4.1 低山丘陵區雷達數據適用性及SAR數據選取分析

3.1節中提到Sentinel-1數據和ALOS-2數據識別滑坡隱患的數目及準確率不一致，經分析主要是受差分干涉相干性的影響。ALOS-2數據搭載L波段的傳感器，波長較長，對地表植被具有較好的穿透性。Sentinel-1數據搭載C波段，使得在植被茂密區域相干性差，尤其是夏季更甚。分別在5月與8月，每個月選取2景覆蓋研究區域時間間隔12 d的Sentinel-1數據進行差分干涉處理，獲取了相干系數對比圖(見圖11)。從圖11中可以看到時間間隔12 d之內，相關系數圖在山區相干性較差，大部分區域相關系數低于0.2，因此在進行相位解纏的時候，絕大部分會被掩膜掉，同時滑坡易發區域位于山區，因此會造成監測區域缺少有效信息。結合圖7可知該區域的滑坡隱患變形主要發生在5～7月，但是這段時間因為植被茂密使得基于Sentinel-1數據進行差分干涉處理，獲取到的監測信息效果較差。

圖11 不同時間段相關系數(間隔12 d)Fig.11 Correlation coefficient diagram for different time periods(interval of 12 days)

為了提高該區域的隱患識別能力和監測效果，經過上面的對比和結合文獻[30]，認為在該區域進行滑坡隱患的早期識別，應選擇分辨率更高或者波長更長的SAR數據，比如日本的ALOS-2衛星數據，阿根廷的SAOCOM-1A衛星數據，以及即將發射的國產SAR衛星。而Sentinel-1數據因高時間分辨率、免費獲取可以用來對典型滑坡隱患開展時序監測分析。

4.2 Stacking InSAR技術適用性分析

在本文中用Stacking 技術進行隱患識別的前提是假設識別對象變形為線性形變，進而獲取了地表形變結果。但是對該區域的滑坡通過SBAS InSAR技術進行時序監測，發現其變形具有一定的周期性，這與假設相矛盾，仍然可以識別出滑坡隱患的原因有：① 收集到的ALOS-2數據集，在夏季有影像數據，使得可以很好獲取降雨期間的變形信息。② Sentinel-1數據在5～8月雖然相干性較差，但是時間分辨率高，差分干涉組合多，因此在加權平均過程中可以獲取到有效變形信息。

需要指出來的是Stacking技術在該區域進行滑坡隱患識別時，可能因為平均計算使得獲取的地表形變形變速率結果值偏小，部分變形速率較小的滑坡隱患存在漏判現象，如文中的宋家灣滑坡，在基于Sentinel-1數據的Stacking結果中該滑坡隱患變形效果并不明顯，但是在累積形變量中均有較強烈的表現。

5 結 論

本文以龍泉山區域為例，采用時序InSAR技術利用多源SAR數據，結合光學遙感技術利用高分光學數據對該區域滑坡隱患開展識別研究，基于SBAS技術對典型滑坡隱患開展時序監測分析，并對識別的隱患進行野外驗證。主要有以下結論：

(1) 利用Stacking技術、結合光學遙感技術成功識別出8處滑坡隱患，并通過野外驗證驗證了其可靠性，發現ALOS-2數據的識別效果優于Sentinel-1數據，也證明了時序InSAR技術用于該區域的滑坡隱患識別具有可行性。

(2) 對基于ALOS-2數據識別的鄭家房子滑坡，通過時序監測發現其最大形變速率為-40 mm/a，最大累積形變量為134 mm，并且其變形規律與汛期相吻合。通過野外調查發現，該滑坡在強降雨作用下，容易發生大規模失移，威脅坡體居民區、村道，因此需要在雨季加強持續性監測，同時結果表明多種時序InSAR技術相結合用于滑坡隱患識別監測具有較強可行性和可靠性。

(3) 通過對比不同波段SAR數據的識別結果及影響因素分析，發現影響該區域識別效果主要是干涉相干性，因此建議該區域采用分辨率更高或者波長較長的SAR數據開展滑坡隱患識別工作。

(4) 通過討論Stacking InSAR技術的適用性，分析了部分滑坡隱患漏判的原因，建議采用多種時序InSAR技術相結合形式開展滑坡隱患識別與監測工作。