基于數字高程模型的古滑坡區提取

張路路，劉召芹，劉峰，李巍

（1.山東科技大學，山東青島266590；2.中國科學院遙感與數字地球研究所遙感科學國家重點實驗室，北京100101）

基于數字高程模型的古滑坡區提取

張路路1，劉召芹2，劉峰1，李巍2

（1.山東科技大學，山東青島266590；2.中國科學院遙感與數字地球研究所遙感科學國家重點實驗室，北京100101）

針對古滑坡的滑前影像無法獲得，植被、紋理信息都已恢復，無法通過對比分析滑坡滑動前后的植被、紋理等信息的變化來提取滑坡區域的問題，提出了一種新的基于數字高程模型的滑坡區域范圍提取方法。該方法基于簡化的滑坡體模型及特征分析，對滑坡區進行水流方向、坡度、山脊山谷線提取，通過流域分析獲取滑坡區域范圍；利用坡度圖實現滑坡壁與滑坡體提取。實驗利用全球30m分辨率ASTER GDEM數據，提取了四川理縣3個古滑坡體區域范圍，驗證了該方法的有效性。

古滑坡；滑坡壁；滑坡體；數字高程模型；流域分析

0 引 言

據初步統計，全國至少有400多個市、縣、區、鎮，10000多個村莊受到滑坡災害嚴重侵害［1］，有證可查的滑坡災害點約為41萬多處，總面積為1.7352×104km2，占國土總面積的18.10%［2］（截至2000年）。上世紀八十年代以來，我國造成10人以上死亡和經濟損失超過500萬元的滑坡就超過130起。這些大型災難性滑坡的發生，不僅帶來重大的人員傷亡和財產損失，而且引發了嚴重的社會和公共安全問題。除此之外，古滑坡的潛在危害也是不容忽視的，古滑坡在強烈的外部作用下能夠重新被激活，并演化成巨型滑坡［3］。公路工程建設中，大量的邊坡開挖，古滑坡體的平衡與穩定遭到嚴重改變和破壞，將會形成規模更大，破壞更嚴重，難于治理的新滑坡，嚴重影響工程建設質量和安全［4］。而我國西部山區有很多村莊都位于古滑坡體上，并且這些滑坡現在仍在變動中。因此，為減少和避免滑坡災害造成的損失，我們首先要勘查圈定滑坡邊界，以估算滑坡致災面積，推測滑坡的規模，為救災快速應急響應提供有效數據支持。

目前已有不少學者開展了滑坡邊界的提取研究，并取得了一定的成果。葉潤青等［5］以遙感影像光譜信息為基礎，利用遙感圖像分類方法將滑坡區分為植被、土體和水體來圈定滑坡邊界；陳瑩［6］以遙感影像為基礎，通過災前和災后水體的提取，檢測河流堵塞區域，獲得準滑坡范圍，再利用簡單的數學形態學處理方法與人機交互判斷，最終識別出滑坡區域；Mondini［7］利用超高分辨率全色影像和高分辨率多光譜影像，通過對比分析滑前滑后影像的變化，圈定了滑坡區范圍；魏冬梅［8］以無人機彩色圖像為基礎，進行滑坡邊緣的色度學差異研究，利用圖像分割技術對滑坡區域進行提??；葉樹林［9］依據放射性核素鈾、鐳、氡等在地質體中普通存在并有一定差異的基本事實提出了利用放射性勘查方法來圈定滑坡邊界；徐金明［10］針對巖體滑坡提出了利用流體包裹體來確定滑坡邊界的方法，主要就是對滑坡滑動過程中流體包裹體跡面特征參數和包裹體熱動力學參數進行了類別分析與對應分析，進而確定滑坡邊界。

綜合國內外文獻資料來看，勘查滑坡邊界的方法主要有兩種：①通過遙感影像來監測滑坡區河流、植被、陰影、建筑物等信息的變化來確定滑坡邊界；②通過地質調查鉆孔分析其巖性來確定滑坡邊界。但對于古滑坡而言，無法得到滑前遙感影像，對于一些光禿滑坡區而言，其上河流、植被、陰影、建筑物等信息較缺乏，無法利用影像分析確定邊界，而實地踏勘將耗費大量人力物力，且作業艱苦，周期較長。針對這一問題，本文從古滑坡區地形特征出發，結合DEM水文分析，提出了一種利用DEM水文分析進行滑坡區域范圍提取的簡便方法。

1 理論與方法

1.1 古滑坡特征分析

在自然地質作用和人類活動等因素的影響下，斜坡上的巖土體在重力作用下沿一定的軟弱面“整體”或局部保持巖土結構完整而向下滑動的過程和現象及其形成的地貌形態，稱之為滑坡［8］，其中發生在一級河流階地侵蝕時期或更早階地侵蝕或堆積期的滑坡稱為古滑坡［11］。因古滑坡發生時間較長，受人類活動影響較大，除滑坡壁、滑坡體、滑坡周界外，其他要素難以識別，為方便研究，本文將古滑坡簡化為如圖1所示模型，并對模型中的三大要素進行提取。

圖1 古滑坡簡化模型

具有與山體分離的部分斜坡可能向前運動的臨空面是滑坡發育的三個基本地質環境條件之一［11］，故山體滑坡大多發生在溝谷兩側，周圍由山脊高陡的巖壁圍繞。受雨水沖刷，人類活動等影響，圖1所示的古滑坡壁坡度已不似新滑坡陡峭，且上沿模糊不清，難以識別；滑坡末端靠近溝谷區域，受人類活動或河水長期沖刷侵蝕，難以辨識。為提取滑坡范圍邊界，此處將滑坡壁上沿延伸至山脊，以山體與溝谷交界線作為滑坡前緣邊界，二者共同構成滑坡區域范圍。因此滑坡區域提取分為兩步：滑坡壁上沿提取和滑坡前緣邊界提取。其提取流程如圖2所示。

圖2 滑坡區域提取流程圖

1.2 提取方法

（1）滑坡壁上沿提取

對于古山體滑坡而言，將滑坡壁上沿延伸至山脊線，即為分水嶺，因此，本文采用流域提取的方法進行滑坡壁上沿的提取。目前應用最廣泛的于DEM數據上提取水文特征的方法是基于地表徑流模擬的方法［12］。地表徑流模擬的方法主要包括三個步驟：一是模擬確定地表徑流方向；二是計算匯水面積累加矩陣；三是確定匯水面積閾值，定義匯水面積超過閾值的像素為地表河網。

①確定地表徑流方向。目前確定水流方向的方法有兩種：單向流法和多向流法［13］。本文采用前者，即D8算法。D8法假設單個網格中的水流方向只有8種可能，即流入它的8個鄰域方向。使用D8法確定地表徑流方向有很多實現方法，本文使用Jenson的方法［14］，通過計算中心網格與各相鄰網格間的距離權落差來確定最陡坡度流向?？珊唵蚊枋鰹椋涸?×3掃描窗口中，設中心柵格為C，鄰域8個水流可能流出的方向分別為E、SE、S、SW、W、NW、N、NE（東、東南、南、西南、西、西北、北、東北），為方便說明和計算，每個流向用一個特定編碼來表示：i（i＝1，2，4，8，16，32，64，128；東，東南，南，西南，西，西北，北，東北）［15］。

圖3 流向及其編碼

②計算匯水面積累加矩陣。確定格網水流方向后，便可進行河網提取。匯水面積累加矩陣值直觀上表示區域地形給水點的水流累積量。其基本思想是：以規則格網表示的數字地面高程模型每點處有一個單位的水量，按照自然水流從高處往低處的自然規律，根據區域地形的水流方向數據計算每點處所流過的水量數值，便得到了該區域的匯流累積量。

③河網的提取。設定合適的集流閾值，從水流累積矩陣圖中提取河流柵格圖。當柵格特征值大于此閾值時，柵格賦值為1，表示河道，小于該閾值的柵格賦值為0，表示產流區。各水道按有效的水流方向連接，便可生成連續的河網［16］。

以上步驟在ArcGIS中即可完成，由于要對5個實驗區進行操作，工作量較大，研究采用ArcGIS10.1中的空間分析建模功能，其模型生成器為設計和實現空間處理模型提供了一個圖形化的建?？蚣?，能將一系列的工具和數據串起來以創建高級的功能和流程［12］。河網提取模型如圖4所示。

圖4 河網提取模型

圖5 實驗區1河網提取結果

由于這里是對反地形進行提取的，所以開始輸入的DEM應該為進行反地形計算后的DEM，柵格計算要根據需求設定合適的閾值。圖5是河網提取結果圖。

圖中水流方向以及河網都是基于反地形數據進行提取的，由圖5可以看出，河網支流較多，滑坡壁位于河網干流，因此首先進行干流提取。由圖5可知，滑坡壁由左右兩條河流主干組成，它們分別為兩邊最長的那根，即由m點和n點出發的河流，它們最后匯聚于p點，一起流向q點，二者在pq處有一小部分的重疊。以此為基礎，本文對各個節點進行遍歷分析，首先找出一條最長的河流線mq（nq），然后找出一條與該河流重疊長度不大于其本身長度10%的另一條最長河流nq（mq）。刪除pq段，即可得到河流主干mpn。

由于ASTER GDEM（V1）數據本身存在一些缺陷，邊界堆疊導致數據顯示為直線、坑、隆起或其他異常的幾何形狀等［17］，導致提取的河網信息有時發生一些錯誤，因此還要結合坡度以及正地形水流方向信息，對其做進一步的判識，對明顯錯誤區進行編輯，從而得到滑坡壁上沿，如圖6所示。

圖6 實驗區1滑坡壁邊界

（2）滑坡前緣邊界提取

滑坡前緣位于山體與溝谷的交界處，滑坡前緣以上的坡面，其坡度大于滑坡前緣以下坡面的坡度?；谶@一特征，本文利用坡度這一基本地形因子對滑坡前緣進行提取，但坡度是一種微觀的地形因子，表示的僅僅是局部地面傾斜，與滑坡前緣具有的宏觀特性不相匹配，所以需要對其作用形式稍作修正［18］?？紤]到實際地形復雜多變，某些地區存在局部突變，同時也為了從滑坡體與溝谷的宏觀地貌差異入手尋找滑坡前緣，需要一種表征地面坡度整體變化的統計指標，因此選擇計算平均坡度，來考察不同區域范圍坡度陡緩情況。具體步驟如下：

①用圖4河網提取模型所示的河網提取模型對滑坡區DEM進行河網提取，得到河網主干線；

②以河網線為基線，向坡體方向做平行線，平行線間距同DEM格網大小，由河網線到坡體方向依次標記為a，b，c，…h；

③計算各平行線處DEM的平均坡度；

④通過分析圖7所示的各個面內DEM數據的平均坡度圖，橫軸為分割面序號，縱軸為坡度值。比較各相鄰面內平均坡度差異，找出變化最陡處，即差值最大者對應的相鄰面（面d和面e），它們之間的分界線即可視為滑坡前緣邊界。

由以上步驟得到滑坡前緣邊界后，與滑坡壁上沿進行疊加，便可得到滑坡區范圍。

（3）滑坡壁與滑坡體提取

滑坡壁是由滑坡體脫離未移動山體而形成的陡壁，坡度較高；滑坡體脫離山體向下滑動，其坡度較緩，因此在滑坡壁與滑坡體交界處，存在坡度陡變。據此，可利用滑坡坡度圖實現滑坡壁的提取，提取算法流程如圖8所示。

圖7 DEM平均坡度圖

圖8 滑坡壁提取流程

提取過程中，首先要對坡度圖像進行二值化，滑坡壁提取要選擇坡度大于某一閾值的數據，滑坡體則相反。由于圖像二值化會造成邊緣損失，這里將二值化后的圖像與原圖像邊緣進行疊加進行邊緣填充，以保證滑坡壁邊緣完整性。而滑坡體本身就處于滑坡內部，不涉及邊緣，故進行滑坡體提取時不需要進行邊緣填充。

2 實驗與討論

2.1 實驗區地理位置

所選的3個滑坡實驗區位于四川省阿壩藏族羌族自治州理縣，地質結構屬龍門山斷裂帶中斷，地形呈蜿蜒起伏的立體單元，地表由西北向東南傾斜。地貌類型為低中山－中山－高山－極高山，是典型的中高山峽谷區。其地理位置坐標如表1所示，圖9展示了其位置分布。

表1 滑坡實驗區地理位置

圖9 研究區地理分布

2.2 數據

本文采用中國科學院計算機網絡信息中心國際科學數據鏡像網站（http：／／www.gscloud.cn）2009年發布的全球空間分辨率為30m的數字高程數據產品，此數據集利用ASTER GDEM第一版本（V1）的數據進行加工得來，ASTER GDEM是使用全自動化的方法對所有從Terra衛星發射后到2008年8月獲取的150萬景ASTER近紅外影像進行處理得到的［19］，覆蓋范圍為83°N～83°S的所有陸地區域，達到了地球陸地表面的99%，其數據精度估計在95%誤差置信水平下，平面誤差為30m，高程誤差為15m～30m［20］。

2.3 實驗結果討論

通過對滑坡區DEM的分析，獲得了滑坡實驗區范圍邊界、滑坡壁以及滑坡體，顯示結果如圖10所示。

由圖7中顯示的3個實驗區的結果圖可以清楚的看出，利用本文所介紹的邊界提取方法獲得的結果與實際邊界基本吻合，可以有效地表達滑坡區域。為了進一步評價提取的邊界精度，利用人工目視解譯出滑坡范圍邊界，由于目視解譯存在人為誤差，本文以目視解譯邊界的3×3鄰域進行判定，如果上述方法提取出來的邊界點落入目視解譯邊界的鄰域中，即可認定該邊界點為正確的邊界點，計算出正確的百分比。經計算得出精度為表2所示。

圖10 滑坡實驗區提取結果

表2 提取精度評估

對邊緣定位的精度評估只是一種定性的評估，從精度的結果來看，本文方法實現的邊界結果與目視解譯相差不大，基本滿足解譯要求。

3 結束語

滑坡壁和滑坡體是滑坡的重要地質參數，它們的確定對滑坡災害的防治有著重要的意義。以上結果表明利用DEM圈定滑坡范圍邊界是有一定理論依據且方法是可行的、有效的。利用DEM水文分析方法在滑坡地區，通過河網提取、水流方向分析、坡度分析對滑坡范圍邊界進行圈定；利用DEM坡度圖對滑坡壁和滑坡體進行提取，此方法無需依靠滑前影像特征，也無需分析紋理特征，克服了傳統方法數據較難獲取，紋理特征不易分析等困難，為滑坡體調查和滑坡制圖提供了一種有效的方法。

［1］黃潤秋.20世紀以來中國的大型滑坡及其發生機制［J］.巖石力學與工程學報，2007，26（3）：433－454.

［2］段永侯.我國地質災害的基本特征與發展趨勢［J］.第四紀研究，1999（3）：208－216.

［3］LIN C W，TSENG C M，TSENG Y H，et al.Recognition of large scale deep－seated landslides in forest areas of Taiwan using high resolution topography［J］.Journal of Asian Earth Sciences，2013，62：389－400.

［4］李勇飛，嵇其偉.淺談古滑坡的判別［J］.西部探礦工程，2006（12）：286－287.

［5］葉潤清，鄧清祿，王海慶.基于圖像分類方法滑坡識別與特征提?。跩］.工程地球物理學報，2007，4（6）：574－577.

［6］陳瑩.基于遙感影像的變化檢測方法在滑坡體提取中的應用［D］.重慶：西南大學，2011.

［7］MONDINI A，GUZZETTI F，REICHENBACH P，et al.Semi－automatic recognition and mapping of rainfall induced shallow landslides using optical satellite images［J］.Remote Sensing of Environment，2011，115（7）：1743－1757.

［8］魏冬梅.基于遙感影像滑坡邊界自動提取方法的研究［D］.成都：西南交通大學，2013.

［9］葉樹林，黃國夫.放射性勘查方法圈定滑坡邊界的機理及應用［J］.中國地質災害與防治學報，2000，11（3）：86－88.

［10］徐金明.流體包裹體法在巖體滑坡周界預測中的應用［J］.土木工程學報，2007，40（1）：69－73.

［11］王治華.滑坡遙感［M］.北京：科學出版社，2012.

［12］李麗，郝振純.基于DEM的流域特征提取綜述［J］.地球科學進展，2003，18（2）：251－256.

［13］黃娜娜，寧芊.基于DEM的數字河網提取方法及應用研究［J］.人民長江，2011，42（24）：50－53.

［14］JENSON S，DOMINGUE J.Extracting topographic structure from digital elevation data for geographic information system analysis［J］.Photogrammetric Engineering and Remote Sensing，1988，54（11）：1593－1600.

［15］張宏鳴.流域分布式土壤侵蝕學坡長提取與分析［J］.咸陽：西北農林科技大學，2012.

［16］沈中原，李占斌，李鵬，等.基于DEM的流域數字河網提取算法研究［J］.水資源與水工程學報，2009，20（1）：20－23.

［17］張朝忙.中國地區SRTM3DEM與ASTER GDEM高程精度質量評價［D］.武漢：華中農業大學，2013.

［18］蘇婧.地形特征線的提取方法研究［D］.西安：西安建筑科技大學，2012.

［19］趙海濤，張兵，左正立，等.中國及周邊區域ASTER GDEM與STRM DEM高程對比分析及互補修復［J］.測繪科學，2012，37（1）：8－11.

［20］FOURNIADIS I，LIU J，MASON P.Landslide hazard assessment in the Three Gorges area，China，using ASTER imagery：Wushan－Badong［J］.Geomorphology，2007，84（1）：126－144.

Extraction of Ancient Landslide Based on DEM

ZHANG Lu－lu1，LIU Zhao－qin2，LIU Feng1，LI Wei2
（1.Geometric College，Shandong University of Science and Technology，Qingdao 266590；2.State Key Laboratory of Remote Sensing Science，Institute of Remote Sensing and Digital Earth，Chinese Academy of Sciences，Beijing100101）

This paper presents a new method of extracting landslide area based on DEM.Flow direction，slope，ridge line and valley line are extracted based on the landslide model and feature analysis，then the landslide area is obtained through watershed analysis；finally the landslide scarp and the displaced material are acquired based on slope map.The ranges of three ancient landslides of Li county from Sichuan province are extracted from ASTER GDEM data with resolution of 30m，verifying the reliability of this method.

ancient landslides；landslide scarp；displaced material；DEM；watershed analysis

10.3969／j.issn.1000－3177.2015.06.007

TP79

A

1000－3177（2015）142－0037－05

2014－09－29

2014－11－17

973計劃課題（2013CB733202、2013CB733201）。

張路路（1991—），女，碩士研究生，主要從事遙感信息提取研究。

E－mail：skzll08＠163.com

劉召芹（1973—），男，副研究員，主要從事遙感制圖與導航定位研究。

E－mail：liuzq＠radi.ac.cn