ArcGIS環境下基于DEM的地形指標和水文信息提取分析

周林虎，張秉來，曹榮泰，張 建，祁兆鑫，劉宇平，胡榮澤，張雯昕

(中國電建集團青海省電力設計院有限公司，青海 西寧 810008）

0 引言

近年來，“3S技術”不斷發展成熟，被越來越多的工程設計和勘察單位所應用，導致地質勘察與測量鄰域也進入了數字時代，尤其是地理信息系統(geographic information systems，GIS)和數字高程模型(digital elevation model，DEM)在地質勘察方面得到了廣泛應用[1]。地理信息，即與空間地理分布有關的信息，具有區域性、多維性和動態性等特點，是指地表物體和環境固有的數量、質量、分布特征，以及聯系和規律的數字、文字、圖形、圖像等的總稱[2]。DEM即指利用有限的高程數據實現地形三維可視化模擬，是一種實體地面模型，為數字地形模型(digital terrain model，DTM)的一個分支[3]。

基于DEM的地表特征提取和洪水淹沒模擬技術主要是根據DEM提供的高程信息，對研究區的高程、坡度變率、坡向變率、地形起伏度、地面粗糙度等基本地表特征信息，以及山脊線、山谷線和山頂點等地形特征參數進行提取和分析，以專題圖的形式進行矢量或柵格數據的輸出，該技術不僅可以對研究區地形地貌分布情況與變化情況進行描述和分析，還可以對區域水土流失、滑坡、崩塌和泥石流等地質災害進行評估，提出合理的避讓和防治措施。

另外，通過模擬洪災發生時的洪水淹沒范圍和深度，可預先獲知和劃分區域洪水受災程度，為洪災評估、預警和防治提供參考[3]。對洪水淹沒情況進行分析亦是防汛抗洪工作的重要內容，是當前GIS 技術在水利應用鄰域的研究熱點。早期的洪水淹沒分析大多以水動力學模型為主，但水動力學模型具有較強的針對性，且構建復雜、所需參數多，無法進行廣泛推廣。隨著對洪水淹沒分析的深入研究，不少學者發現洪水淹沒程度主要取決于洪水特征以及地形因子，用DEM 能夠較好地反映地表起伏狀態[4]。特別是近些年，隨著暴雨和強降雨的不斷加劇，我國多地洪澇災害數量增多、程度加重，因此對區域洪水進行模擬，科學、準確地預測洪水淹沒范圍和洪水深度分布對防洪減災具有重要的現實意義[5-6]。

基于此，以青海省某一區域為研究區，利用ArcGIS軟件提取和分析了該區域基本地表特征信息，生成了相關專題圖件；并在此基礎上進行了研究區洪水淹沒分析，為該區域的城鎮發展、農業規劃、資源開發、災害預測及治理等提供指導性建議。

1 研究區概況

本研究區海拔為2 071.96～2 344.55 m，高低起伏較大，最大高差達272.59 m。區內地層以黃土和泥巖為主，滑坡、泥石流和崩塌等地質災害發育較多，尤其6～9月份降雨較多，易發生地質災害。該區域地層裸露較多，地表植被相對較少，為水土流失較嚴重的地區。長期而頻繁的水土流失、滑坡和泥石流等地質災害的發生，使得該區域土地退化嚴重，生態環境不斷惡化，嚴重影響和阻礙了該區域工程設施的建設和社會經濟的可持續發展，因此對該地區地形地貌特征進行分析，對地質災害進行調查和評價，并提出合理的防治和治理措施具有重要意義。

2 傾斜攝影測量技術

1)優點和應用前景

在以往的工程地質災害勘察工作中，地質災害調查分析多以勘察員的工程經驗判斷為主，受現場視野的局部性、地形、現場復雜條件、勘察手段和人員安全因素等限制，勘察人員往往無法到達現場，野外實施過程中地質災害點要素描寫欠全面、規范，災害成因分析不足，資料收集不夠完整，地質災害勘察結果只是定性分析得出的結果，而無法得出定量分析結果。利用無人機傾斜攝影測量技術，通過獲得的航空影像和DEM數據，可準確地對研究區坡度、坡向、流向和流量等地表特征進行分析和計算，從而調查地質災害點區域自然環境特征和孕育地質災害的地質背景，查明地質災害的空間分布現狀和強度等發育特征，得出地基穩定性評價，節約成本，提質增效。

2)關鍵技術難點

近幾年，傾斜攝影測量技術得到了迅速發展，傾斜攝影測量技術能夠獲取建筑物、樹木等地理實體的紋理細節，不但豐富了影像數據源信息，同時高冗余度的航攝影像重疊為高精度的影像匹配提供了條件，使得基于人工智能的三維實體重建成為了可能。分層顯示技術、紋理映射技術成為傾斜攝影測量和建模的關鍵支撐點，極大地提升了三維建模的效率，同時也降低了建模生產成本。目前，基于傾斜攝影測量成果的應用還比較少，因此，應用創新點挖掘還需要逐步深入。同時不能局限于傾斜攝影測量技術，還需要研究其與雷達、紅外、多光譜、高光譜等多種傳感器的結合，集成在更小的無人機上拓寬攝影測量技術的應用范圍。未來，基于點云數據計算的大規模三維數據生產使得工程測量、三維建模等工作進行變革，開啟三維遙感的新時代。

3 地形指標提取

地形指標主要包括坡度變率、坡向變率、地形起伏度和地面粗糙度等4個基本地形指標，即最基本的自然地理要素，同時也是對人類生產、生活影響最大的自然要素[7]。地形特征值廣泛應用于諸多研究和應用鄰域，地形指標的提取對水土流失、土地利用、土地資源評價、城市規劃等方面的研究起著重要作用。

3.1 坡度變率

坡度變率，又名剖面曲率，是指坡度在微分空間上的變化率，地表面高程相對于水平面變化的二階導數，即坡度之坡度(slope of slope，SOS)[8]。其能很好地反應出地面高程變化率的大小，從而反應出剖面曲率信息。根據ArcGIS表面分析工具，提取出的坡度和坡度變率如圖1所示。

圖1 研究區坡度和坡度變率

對于坡度和坡度曲率的研究和分析不僅能直觀判斷和分析研究區地形高低起伏情況、微地形地貌，對于土地利用規劃和資源開發具有參考價值，還能作為衡量和評價區域物質能量轉移、水土流失、土壤侵蝕以及洪澇災害的重要基礎數據，對于地質災害的評價和防治具有重要意義[9]。在我國，坡度為25°是水土保持法規定的臨界坡度[10]，即坡度大于25°的地區，其水土保持能力相對較弱，相同條件下發生水土流失、土壤侵蝕、崩塌、滑坡和泥石流等地質災害的概率亦較大。通過分析研究區坡度數據分布情況，得到區內坡度大于25°的區域面積占總面積的49.86%，由此可知該流域內一半左右區域坡度較大，水土保持能力較弱，且結合現場植被分布情況可知，該部分區域植被覆蓋率較低，水土流失情況嚴重，區內存在較多沖溝，加之豐水季節雨水較多，因此判定該區域發生地質災害的概率較大。為防止水土流失進一步加劇，需采取專門的水土保持措施。通過坡度變率圖可知，圖中紅色高值即坡度圖中坡度變化較顯著、落差較大的區域，即地勢急劇變化的區域。通過計算得到，該區域坡度變率為80～90的區域占比82.75%，見表1所列，說明該區域地形起伏較大，坡度變化較顯著，多為山地、丘陵和溝谷地貌，不適宜發展大規模機械化種植業，可發展林業、山地畜牧業以及林下種植業。

表1 不同坡度變率面積占比

3.2 坡向變率

坡向變率，又名平面曲率，是指坡向在微分空間上的變化率，對提取的坡向數據再進行一次坡度提取，即坡向之坡度(slope of aspect，SOA)，其可以很好地反映等高線的彎曲程度[8]。坡向和坡向變率提取方法如下：①通過坡向工具提取坡向數據DEM1；②提取坡向變率SOA1，即對坡向數據DEM1提取坡度，得到坡向變率SOA1；③提取反地形數據FDEM通過柵格計算器工具，輸入計算式為FDEM=2 344.55-DEM，2 344.55 m即本流域最大高程；④提取反地形坡向變率SOA2，即對FDEM提取反地形坡向FDEM1，然后對坡向數據FDEM1提取坡度，得到反地形坡向變率SOA2；⑤計算SOA。使用柵格計算器工具，輸入式為SOA=[(SOA1+SOA2)-ABS(SOA1-SOA2)]/2，即可得到無誤差的DEM坡向變率。

提取出的坡向和坡向變率如圖2～圖4所示，通過SOA可知，其與坡度變率分布情況基本相似，但是通過計算不同坡度變率區域面積占比情況得出，坡向變率為0～10和80～90的區域面積占比為8.29%和88.68%，均高于坡度變率為0～10和80～90的面積占比，而坡向變率為10～80的區域面積總體占比低于坡度變率為10～80的面積，見表2所列。

表2 不同坡向變率面積占比

圖2 研究區DEM1和SOA1

圖3 研究區FDEM1和SOA2

圖4 研究區SOA

3.3 地形起伏度

地形起伏度是描述某一特定區域內高程極差(即最高點和最低點之間高差)的地形指標[11-12]。鄰域范圍可根據研究區具體情況和實際需求進行設置，鄰域范圍越大，求得的起伏度亦越大。本項研究將鄰域范圍設為12×12，而研究區像元柵格總數為200 970 931，即將研究區分為1 395 631個矩形柵格，求解這些矩形柵格像元內的像元值極差，統計結果即該區域內的相對高差(即起伏度)。通過ArcGIS焦點統計工具可直接生成起伏度成本數據，如圖5所示。由該圖可知，地形起伏度范圍為0～118.75 m，其中起伏度為0～13.19 m的區域面積占98.946 3%，而隨著地形起伏度不斷增大，所占面積比例也越來越小，見表3所列。

表3 不同地形起伏度面積占比

圖5 研究區地形起伏度

3.4 地面粗糙度

地面粗糙度是反映區域表面形態的一個宏觀指標，指某一特定區域內地球表面積與其投影面積之比[13]。根據地面粗糙度的定義，求每個柵格單元的表面積與其投影面積之比。如圖6所示，假如△ABC是一個柵格單元的縱剖面，a為此柵格單元的坡度，則AB面的面積(AB×AC)為此柵格的表面積，AC面的面積(AC×AC)為此柵格的投影面積(也既是此柵格分辨率平方)，此柵格單元的地面粗糙度M=“AB面的面積”/“AC面的面積”=(AC×AB)/(AC×AC)=AB/AC=1/cosa。根據ArcGIS表面分析工具計算出研究區地面粗糙度如圖7所示，可知研究區地面粗糙度為1～602.26，但大部分區域介于1～67.81，隨著地面粗糙度不斷增大，所占面積比例亦逐漸減少，這和地形起伏度表現出相同特征，見表4所列。

表4 不同地面粗糙度面積占比

圖6 △ABC

圖7 地面粗糙度

3.5 山脊線和山谷線

特征地形要素主要是指對地形在地表的空間分布特征具有控制作用的點、線或面狀要素[14]。特征地形要素構成地表地形與起伏變化的基本框架，與地形指標的提取主要采用小范圍鄰域分析不同的是，特征地形要素的提取更多地應用較為復雜的技術方法，如山谷線、山脊線等的提取采用了全局分析法，成為柵格數據地學分析中很具特色的數據處理內容[15]。

山脊線和山谷線的提取步驟為：①點擊初始DEM數據，使用焦點統計工具，設置統計類型為平均值MEAN，鄰域的類型為矩形(也可以為圓)，鄰域的大小為12×12(根據實際情況而定)，則可得到一個鄰域為12×12的矩形的平均值數據層，記為MEAN，如圖8(a)所示；②使用柵格計算器計算，輸入公式ZFDX=DEM-MEAN，即可求出正負地形分布區域，如圖8(b)所示；③采用柵格計算器工具，輸入公式shanji=(ZFDX＞0)&(SOA＞80)，即可提取出山脊線，如圖9(a)所示；④同理，使用柵格計算器工具，輸入公式shangu=(ZFDX＜0)&(SOA＞80)，即可求出山谷線如圖9(b)所示。通過圖9可看出研究區地形地貌的山梁、山坡、溝谷和洼地等分布情況，從而直觀分析出區內匯流和分流情況，為后續洪水淹沒分析結果的可靠性和真實性提供參考。

圖8 MEAN和ZFDX

圖9 研究區山脊線和山谷線

4 洪水淹沒分析

判別給定降水條件下某一區域的潛在淹沒區，既可以規避盲目建設所帶來的雨洪淹沒風險，也能為生態保護與修復提供可選的項目突破點，因而在區域規劃設計中具有重要作用。洪水淹沒分析一般分為無源淹沒和有源淹沒2種分析方法，洪水淹沒分析是進行水文預測預報、洪水災害評估的一項重要內容，也是GIS地形三維仿真系統中的一個重要功能[16-17]，前者是指低于某一高程的所有柵格均被劃入淹沒區，相當于區域大面積均勻降水，所有低洼處都有可能積水成災；后者考慮了柵格之間的物理連通性，因而某些柵格受到地形阻擋而不會被淹沒，且除降水外還考慮地表水，如河流、水庫等的外溢現象[18]。與有源淹沒相比，無源淹沒方法更為簡單，且在受地表水流動影響較小的地區能實現良好的模擬效果，因而得到廣泛應用。由于本項研究區周圍未發現能形成洪水的水庫和河流，洪水基本上由當地強降雨引起，因此洪水淹沒模擬采用無源淹沒方法。

由于本項研究區高程差較大，高低起伏較顯著，為直觀顯示出不同淹沒深度條件下洪水淹沒區域和面積占比，本次分析設置了10 m、50 m、100 m、150 m和200 m等5個梯度。在實際分析和模擬中，則應根據實際高程數據和地形情況設置合理的淹沒深度進行分析。表5為不同淹沒深度條件下淹沒面積占比，由該表可知，隨著淹沒深度不斷增大，淹沒面積亦顯著增加。當淹沒模擬深度為10 m時，淹沒面積僅為0.62%，該區域主要集中于研究區南部一小部分區域，隨著淹沒模擬深度增大100 m，淹沒區域逐漸擴散至南部大部分區域和西北部，淹沒面積達到35.38%，由此亦可分析出本項研究區高程由低至高依次為南部、西北部和中東部，與實際地形相符，如圖10所示。

表5 不同淹沒深度條件下淹沒面積占比

圖10 不同淹沒深度條件下淹沒面積示意圖

值得討論的是，本項洪水淹沒分析對于相互聯通的小面積區域具有良好的模擬效果，但對于高程變化較大的大面積區域，由于有許多分水嶺和匯水溝的存在，一般其內部并非完全貫通，因此不能形成完整的水域，與實際情況會略有偏差。

5 結語

基于ArcGIS的3D Analyst、表面分析和柵格計算分析模塊，提取和分析了研究區坡度、坡度變率、坡向、坡向變率、地形起伏度和地面粗糙度等地形變化指標，以及山脊線和山谷線等特征地形要素，從而得到了區內地形高低起伏情況和微地形地貌特征，在此基礎上通過ArcGIS的柵格計算工具，提取了不同深度洪水情況下的淹沒面積和洪水位，結合區域地形情況，驗證了洪水淹沒分析的真實性和可靠性。該方法不需要大量的現場勘察和搜集資料工作即可得到研究區整個地形地貌、高程分布情況和洪水淹沒情況，節省了大量時間和人力。