基于精細DEM 的崩塌滑坡災害識別及主控因素分析
——以雅魯藏布江縫合帶加查—朗縣段為例*

王瑞琪 王學良 劉海洋 孫娟娟 王新輝 張 蘇

( ①中國科學院地質(zhì)與地球物理研究所，中國科學院頁巖氣與地質(zhì)工程重點實驗室 北京 100029)

( ②中國科學院大學 北京 100049)

( ③中國科學院地球科學研究院 北京 100029)

( ④北京洛斯達數(shù)字遙感技術有限公司 北京 100120)

( ⑤國網(wǎng)經(jīng)濟技術研究院有限公司 北京 102209)

0 引 言

雅魯藏布江縫合帶加查—朗縣段位于青藏高原東南部，崩塌、滑坡等地質(zhì)災害十分發(fā)育，伴隨著青藏高原地區(qū)人類工程的日益增多，該區(qū)域的地質(zhì)災害識別、發(fā)育特征與成因機理研究迫在眉睫，以求為后續(xù)的風險評估工作打下基礎。

近年來，針對雅魯藏布江縫合帶的工程地質(zhì)條件及地質(zhì)災害分布情況，有許多學者展開了調(diào)查研究。如張路青等( 2004) 對川藏公路南線八宿—林芝段滾石災害進行了工程地質(zhì)調(diào)查，給出了初步的防護建議; 袁廣祥等( 2010) 對雅魯藏布江大拐彎北部川藏公路地質(zhì)災害發(fā)育與分布開展了研究; 陳洪凱等( 2011) 運用層次分析法分析了川藏公路沿線的地質(zhì)災害的危險性，并建立了相應的評價模型和評價數(shù)據(jù)庫; 吳瑞安等( 2017) 對川藏鐵路加查—朗縣段幾處重大滑坡災害進行了詳細闡述; 李秀珍等( 2017) 結(jié)合現(xiàn)場調(diào)查和遙感解譯等手段，研究了川藏交通廊道康定至林芝段滑坡崩塌的空間分布特征和危害方式。限于數(shù)據(jù)精度和工作程度等方面的原因，如何在研究區(qū)開展災害面線上的精細識別和主控因素分析尚顯不足。本文利用機載雷達獲取了精細DEM( 10 m 分辨率) ，UAV 獲取高清航攝數(shù)據(jù)，在已有的崩塌滑坡災害遙感目視解譯識別方法的基礎上，提出了區(qū)域尺度上基于統(tǒng)計學方法的崩塌滑坡災害的精細識別方法; 開展了研究區(qū)崩塌滑坡主控因素分析，尤其是各地形地貌因子對崩塌滑坡的影響程度，并進行了可靠性驗證。

1 地質(zhì)背景條件

加查—朗縣段位于藏東南地區(qū)，雅魯藏布江中上游地帶，區(qū)內(nèi)河流下切與側(cè)蝕作用強烈，山高谷深，斜坡陡立，地形起伏大，坡體植被發(fā)育較少，巖體裸露。結(jié)合區(qū)域地質(zhì)資料，區(qū)內(nèi)分布有3 條與雅魯藏布江縫合帶大致平行的大型伴生斷裂( 圖1) ，新構造運動強烈，斷裂帶附近巖體破碎。區(qū)內(nèi)巖性組成較為簡單，主要出露千枚巖、礫巖和花崗巖。研究區(qū)主要位于雅魯藏布江兩岸，高程在3060～5330 m 之間，主要集中在3100～4400 m 的高程范圍。在內(nèi)外動力地質(zhì)作用下，研究區(qū)容易發(fā)生崩塌、滑坡等地質(zhì)災害。

2 研究區(qū)地質(zhì)災害高精度識別

2.1 崩塌、滑坡點識別

研究區(qū)的崩塌、滑坡災害的識別技術主要采取遙感技術，包括衛(wèi)星遙感和無人機遙感技術( Mckean et al.，1991; Mantovani et al.，1996; Akgun et al.，2012; Scaioni et al.，2014; 黃海峰等，2017; 許強等，2018; 唐堯等，2019) 。數(shù)據(jù)來源主要為機載雷達獲取的精細DEM 數(shù)據(jù)( 10 m 精度) 、高清航拍數(shù)據(jù)，以結(jié)合衛(wèi)星遙感影像作為輔助數(shù)據(jù)。根據(jù)王學良( 2018) ，本文采用的歷史崩塌滑坡識別方法步驟如下:( 1) 建立研究區(qū)的崩塌滑坡解譯標志。主要依據(jù)高精度DEM 顯示的歷史崩塌、滑坡形態(tài)，輔以遙感影像的災害特征，建立相應的解譯標志。( 2) 利用ArcGIS 平臺，遵循建立的解譯標志，進行研究區(qū)的歷史崩塌滑坡識別。( 3) 現(xiàn)場驗證與成果完善，采用無人機和傳統(tǒng)工程地質(zhì)調(diào)查手段，進行現(xiàn)場驗證，驗證數(shù)為26 處，準確數(shù)為24 處，遙感解譯準確率可達92.3%。( 4) 對解譯的歷史崩塌滑坡進行細化圈定，基于10 m 精度DEM，利用統(tǒng)計學方法，分析崩塌滑坡與各因素的相關性，從而形成精細地形數(shù)據(jù)的崩塌滑坡識別方法。

對研究區(qū)滑坡、崩塌災害進行遙感解譯，共計解譯崩塌41 處，滑坡92 處，分布情況如圖1 所示。

2.2 崩塌、滑坡規(guī)模識別

基于ArcGIS 平臺，圈定崩塌滑坡的范圍之后，可以在屬性表中自動計算各災害點的面積，為獲取各災害點的規(guī)模，還需獲取崩塌滑坡的厚度。利用高精度DEM 數(shù)據(jù)，對這些崩塌滑坡截取剖面。滑坡深度的估算主要基于現(xiàn)場調(diào)查，依據(jù)滑坡前緣高差和后緣陡坎的高差，通過在GIS 中截取滑坡剖面并獲取其陡坎高差實現(xiàn)，而崩塌堆積物主要呈錐形，其平均厚度取最大厚度值的二分之一。

圖1 研究區(qū)崩塌滑坡識別結(jié)果Fig. 1 Identification results of collapse and landslides in the study area

研究區(qū)內(nèi)部分滑坡由于風化作用，后緣陡坎已經(jīng)漸趨平緩，后緣不存在明顯陡坎，而能夠在坡腳形成堆積體的大型崩塌較少，因此選取其中50 處具有明顯陡坎的滑坡以及6 處具有錐形堆積體的大型崩塌，共計56 處崩塌滑坡災害點進行統(tǒng)計。分析發(fā)現(xiàn)災害面積主要集中在105～106im2，占災害總數(shù)的66.04%。將面積與厚度進行乘法計算后，發(fā)現(xiàn)災害體積主要集中在105～107im3，占災害總數(shù)的84.91%。

對這56 處崩塌滑坡的面積與規(guī)模進行分析，發(fā)現(xiàn)地質(zhì)災害體積與面積大致呈現(xiàn)指數(shù)關系，Larsen et al. ( 2010) 認為，滑坡( 包括巖質(zhì)滑坡和土質(zhì)滑坡) 的體積和面積滿足V = αAλ，對統(tǒng)計結(jié)果進行擬合( 圖2) ，根據(jù)擬合結(jié)果，α = 5.6367，λ = 1.154，R2= 0.7558，即V = 5.636 7A1.154，擬合效果較好，反映了數(shù)據(jù)的準確性和規(guī)律性，并可以通過此公式估算區(qū)內(nèi)崩塌滑坡體積。

圖2 研究區(qū)地質(zhì)災害體積與面積關系擬合Fig. 2 Fitting curve of geological hazard volume and area

3 研究區(qū)崩塌滑坡發(fā)育主控因素分析

3.1 地層巖性

崩塌主要集中發(fā)育在花崗巖、礫巖等硬脆巖體之中( Wang et al.，2012，2014) 研究區(qū)的41 處崩塌均發(fā)生在硬脆的礫巖和花崗巖區(qū)域，花崗巖位于高海拔位置，坡度通常較大，受到內(nèi)外動力地質(zhì)作用，極易發(fā)生崩塌，礫巖巖體由于坡度較陡也常伴有崩塌發(fā)生，多為中型崩塌。千枚巖具有千枚構造，在構造作用下，巖體破碎，且坡度通常較緩，因而92 處滑坡均在千枚巖中，以中型滑坡為主，邊界通常為硬巖和軟巖的接觸界面( 圖3) 。

圖3 研究區(qū)特征剖面( 扎西林村)Fig. 3 Typical section of the study area( Zhaxilin Village)

3.2 地質(zhì)構造

區(qū)內(nèi)伴生雅魯藏布江縫合帶的斷裂，走向呈現(xiàn)近東西向，其傾向大致為180°～220°，傾角在65°左右，斷裂帶對區(qū)內(nèi)地質(zhì)災害發(fā)育具有控制作用。

3.2.1 斷層

結(jié)合現(xiàn)場對斷層位置的進一步確認和修正，對92 處滑坡與斷裂的相對位置進行統(tǒng)計分析發(fā)現(xiàn)，以斷層為前邊界的牽引式滑坡數(shù)量最多，有46 處，占滑坡總數(shù)的50%; 以斷層為后邊界的推動式滑坡最少，為13 處，約占滑坡總數(shù)的14.1%; 滑坡體前部和后部均有斷層穿過的混合式滑坡為33 處，占滑坡總數(shù)的35.9%。對滑坡的滑坡方向與斷裂的走向進行統(tǒng)計分析，發(fā)現(xiàn)有78 處滑坡的滑動方向垂直于斷裂走向，約占滑坡總數(shù)的84.8%，如扎西林村滑坡、拉崗村滑坡等。在斷裂帶的影響下，斷裂附近崩塌分布較為密集，且普遍規(guī)模較大。

3.2.2 巖體結(jié)構

研究區(qū)內(nèi)主要發(fā)育滑移式崩塌、傾倒式崩塌和拉裂式崩塌。以加查鎮(zhèn)某崩塌點為例，分析該類崩塌的破壞模式和控制特征。崩塌點的破壞模式分析主要依靠無人機航攝技術，利用小型多旋翼無人機獲取巖體的結(jié)構面信息( 劉海洋等，2017) 。工作區(qū)內(nèi)高差大，為保證地面分辨率能達到設計要求，采取旁向65%重疊度，航向75%重疊度，航高大致保持在300 m。航空影像數(shù)據(jù)主要包括正射影像數(shù)據(jù)和三維模型影像數(shù)據(jù)。

三維模型的生產(chǎn)主要通過影像密集匹配，三角網(wǎng)構建、自動賦予紋理等步驟，最終得到三維模型。應用AgisoftPhotoscan 建模軟件對三維模型進行制作加工，處理成果影像，完成三維實景模型制作( 圖4a) 。選取巖體結(jié)構明顯的部分坡體進行結(jié)構面提取分析，將點云數(shù)據(jù)以文本文件的形式，從Photoscan 中導入Coltop3D( 圖4b) ，利用Coltop3D軟件提取巖體結(jié)構產(chǎn)狀信息。

圖4 加查鎮(zhèn)某崩塌三維點云模型Fig. 4 A 3D point cloud model of a collapse in Jiacha Town

此崩塌點以礫巖為主，礫石成分為花崗質(zhì)，主要發(fā)育3 組結(jié)構面，其中J1為礫巖層面。利用Coltop軟件進行產(chǎn)狀統(tǒng)計，共統(tǒng)計57 條結(jié)構面。J1的產(chǎn)狀統(tǒng)計中，最大值為199°∠46°，最小值為180°∠35°，均值為190°∠40°; J2的產(chǎn)狀統(tǒng)計中，最大值為231°∠80°，最小值為208°∠59°，均值為221°∠74°; J3的產(chǎn)狀統(tǒng)計中，最大值為289°∠90°，最小值為263°∠70°，均值為271°∠87°; 取均值作為這3 組結(jié)構面產(chǎn)狀，即:J1=190°∠40°，J2=221°∠74°，J3=271°∠87°，坡面產(chǎn)狀為203°∠50°。

統(tǒng)計3 組結(jié)構面間距，其中J1結(jié)構面間距在3.2～20.1 m 之間，在5 ～6 m 間距范圍內(nèi)的結(jié)構面較多，此外，間距范圍在3 ～9 m 各間距范圍內(nèi)的結(jié)構面間距分布較均勻，均值為5.5 m，間距較小; J2結(jié)構面間距在1.3 ～10.7 m 之間，主要集中在4 ～5 m 之間，均值為8.5 m，間距較大; J3結(jié)構面間距在1.7～72.1 m 之間，主要集中在2～3 m 范圍之間，均值為18.0 m，間距大。就斜坡整體而言，各組結(jié)構面間距較大，切割較為強烈，巖體較破碎，容易發(fā)生大塊石的崩落。

圖5 主要結(jié)構面的穩(wěn)定性分析Fig. 5 Stability analysis of main structural planes

利用這57 條結(jié)構面的產(chǎn)狀和坡面產(chǎn)狀繪制赤平極射投影圖( 圖5) 。由于上半球投影較為常用，因而采用此投影方法。其中，綠色線段DO，EO，F(xiàn)O的方向分別為J1、J2、J3傾向，坡面傾向為MO，與J1產(chǎn)狀相近。J1、J2、J3分別相交于A、B、C 3 點，連接3點與圓心O 點，形成的紅色線段即為交線AO，BO和CO。根據(jù)李智毅等( 1994) ，J2與J3的交線BO傾向與坡面相近，但傾角大于坡角，因而不易沿J2、J3交線滑動; J1、J2交線AO 和J1、J3的交線CO 傾角均小于坡角，AO 位于J1、J2傾向線DO 和EO 的夾角之外，因而楔形體沿著J1發(fā)生滑動，J2僅起到切割作用; CO 在J1、J3傾向線DO 和FO 的夾角之外，因而楔形體沿J1發(fā)生滑動，J3同樣僅起到切割作用。綜合上述分析，結(jié)構面切割形成的楔形體主要沿J1發(fā)生滑移式崩塌。

3.3 地形地貌

分析崩塌滑坡識別結(jié)果可以發(fā)現(xiàn)，滑坡主要發(fā)育在雅魯藏布江南岸以及北岸坡體的中下部，而崩塌主要發(fā)生在北岸坡體的中上部?？刂票浪?、滑坡發(fā)生的地形地貌因素主要為坡度、坡向和高程，為準確確定崩塌、滑坡受坡度、坡向和高程的控制情況，對于崩塌物源區(qū)的2390 個柵格和滑坡面內(nèi)的311 792 個柵格進行統(tǒng)計分析。利用ArcGIS 軟件的轉(zhuǎn)換工具和提取工具，將這些柵格的坡度、坡向以及高程信息提取至點文件。將坡度數(shù)據(jù)按5°間隔，劃分為18 個等級; 將坡向按照20°間隔，劃分為18 個等級; 將高程按照100 m 間隔，劃分為24 個等級。

按照這個分類對于全區(qū)內(nèi)坡度、坡向和高程數(shù)據(jù)進行重分類，并統(tǒng)計各數(shù)據(jù)范圍內(nèi)的面積，以崩塌和滑坡面范圍內(nèi)各坡度、坡向、高程內(nèi)的柵格面積與區(qū)內(nèi)各坡度、坡向、高程范圍內(nèi)的柵格面積作比，作為崩塌和滑坡的面密度，由于每個柵格面積相同，即可轉(zhuǎn)化為柵格個數(shù)的比值。

分析發(fā)現(xiàn)，在0°～80°的坡度范圍內(nèi)，崩塌面密度隨坡度上升而大致呈現(xiàn)增大趨勢，而滑坡面密度在20°～25°坡度范圍內(nèi)之間呈現(xiàn)一個峰值，而后隨坡度上升呈現(xiàn)下降趨勢( 圖6) 。不同坡向的崩塌面密度變化趨勢與崩塌數(shù)量一致，主要分布在140°～260°的坡向范圍內(nèi); 滑坡面密度與坡向關系較小，不同坡度范圍內(nèi)滑坡分布較為離散( 圖7) 。崩塌面密度隨高程上升而大致呈現(xiàn)增大趨勢，主要分布在3800～4900 m 高程范圍內(nèi); 滑坡在3200～4000 m 高程范圍內(nèi)的分布情況較為分散，在3800 m 以下和4000 m 以上坡度范圍內(nèi)滑坡面密度較低( 圖8) 。

圖6 災害與坡度關系圖Fig. 6 Relationship between hazards and slope

圖7 災害與坡向關系圖Fig. 7 Relationship between hazards and slope direction

圖8 災害與高程關系圖Fig. 8 Relationship between hazards and elevation

為研究各因子對崩塌滑坡發(fā)生的影響程度，應用曲線下面積法( AUC，Area Under the Curve) 。以各因子不同等級柵格累計比值為橫軸，災害( 崩塌或滑坡) 柵格累計比值為縱軸，面密度從高到低排列，用每組的面密度值高低表征發(fā)生崩塌或滑坡可能性的大小。利用Origin2019 專業(yè)繪圖軟件，對兩組數(shù)據(jù)繪制折線圖，即為受試者工作特征曲線( ROC曲線，Receiver-operating Characteristic Curve) 。

對于各ROC 曲線求其AUC 值，分析發(fā)現(xiàn)，針對坡度而言，崩塌的AUC 值為0.914，略大于滑坡的AUC 值( 0.835) ，表明坡度對于崩塌的控制作用明顯大于滑坡; 就坡向和高程而言，崩塌的AUC( 0.787 和0.848) 值均遠大于滑坡AUC 值( 0.528和0.682) ，表明坡向和高程對于崩塌具有一定的控制作用，而對滑坡的影響更小。因此，可以認為研究區(qū)內(nèi)的崩塌由坡度、坡向和高程共同控制，其中坡度為主控因素; 而滑坡與坡度的相關性較大，高程和坡向?qū)碌陌l(fā)育具有一定的影響，但影響較小。

4 結(jié) 論

(1) 利用遙感影像數(shù)據(jù)與10 m 精度無人機航攝數(shù)據(jù)，對崩塌滑坡進行識別，共計解譯崩塌41 處，滑坡92 處，在構造作用和巖性控制下，呈現(xiàn)條帶狀分布。崩塌主要集中在斷層帶附近的礫巖、花崗巖巖體，以滑移式崩塌為主，而滑坡主要發(fā)育在千枚巖區(qū)域，以牽引型為主，84.8%的滑坡滑動方向與斷裂走向大致垂直。

(2) 基于10 m 精度的DEM 數(shù)據(jù)，通過截取崩塌滑坡剖面獲取其厚度( 或深度) ，發(fā)現(xiàn)區(qū)內(nèi)崩塌滑坡的體積和面積滿足公式V=5.636 7A1.154，擬合效果較好，可為崩塌滑坡災害的規(guī)律識別提供參考。

(3) 從地層巖性、地質(zhì)構造、地形地貌3 個方面分析崩塌的主控因素，表明崩塌的主控因素較多，地層巖性決定了研究區(qū)崩塌和滑坡的發(fā)育類型，地形中的坡度、坡向和高程因素均對崩塌的形成重要影響作用，其中，坡度與崩塌發(fā)生的相關性最大，滑坡的發(fā)育受斷層的控制作用明顯，地形地貌中的坡度因子對滑坡影響作用較大，而坡向和高程對于滑坡的貢獻較小。