ArcGIS環(huán)境下基于DEM的流域特征提取

王玉富,王翰釗

(鄭州測繪學(xué)校，河南 鄭州 450015)

DEM(Digital Elevation Model)是用一組有序數(shù)值陣列形式表示地面點(diǎn)的平面坐標(biāo)(x，y)和高程z的一種實(shí)體地面模型.它包含了大量的地理信息，其用途十分廣泛，可用于諸如土方量計(jì)算、通視分析、水流流量分析、水系網(wǎng)絡(luò)分析、降雨量分析等，甚至移動通訊基站的布設(shè)也可以用DEM進(jìn)行分析.

目前，利用DEM 進(jìn)行流域分析的工具很多，ArcGIS的水文分析模塊(Hydro logy model)是美國環(huán)境系統(tǒng)研究所公司(ESRI)為ArcGIS推出的一個(gè)水文分析模塊，主要用于地形和河流網(wǎng)系的提取和分析，實(shí)現(xiàn)地形模型可視化,其強(qiáng)大的流域特征分析功能可以滿足各種流域DEM處理的需要.

1 DEM構(gòu)建方法

DEM構(gòu)建的數(shù)學(xué)過程是：根據(jù)區(qū)域內(nèi)所有的高程點(diǎn)數(shù)據(jù)，用傅立葉級數(shù)和多項(xiàng)式函數(shù)擬合地面高程曲面，即用規(guī)則的數(shù)學(xué)曲面擬合不規(guī)則的實(shí)際地面，然后再將其分成正方形規(guī)則區(qū)域或面積大致相等的不規(guī)則區(qū)域進(jìn)行分塊搜索，根據(jù)有限個(gè)點(diǎn)進(jìn)行擬合，形成詳細(xì)高程曲面.主要模型有：規(guī)則格網(wǎng)模型、等高線模型、不規(guī)則三角網(wǎng)模型和層次模型等.

根據(jù)工作目的、信息源情況和數(shù)據(jù)采集設(shè)備情況，目前構(gòu)建DEM主要有以下三種方法：①影像數(shù)據(jù)源；②地形圖數(shù)據(jù)源；③地面實(shí)測記錄[1，2].基于前兩種方法的DEM構(gòu)建，需要將柵格數(shù)據(jù)處理為矢量數(shù)據(jù)，并對矢量數(shù)據(jù)進(jìn)行檢查后方可進(jìn)行DEM模型生成；第三種則可直接進(jìn)行DEM模型構(gòu)建，但需要實(shí)測的數(shù)據(jù)量大，相應(yīng)的工作量也很大.目前多用前兩種方式進(jìn)行DEM構(gòu)建，其流程如圖1所示.

圖1 DEM構(gòu)建流程

2 ArcGIS下DEM流域特征提取

通過DEM特征提取進(jìn)行流域特征分析，主要有兩個(gè)過程：一是地形特征的提取，二是水系特征的提取[3].

2.1 地形特征的提取

地形特征點(diǎn)主要包括山頂點(diǎn)(peak)、脊點(diǎn)(ridge)、谷點(diǎn)(channel)、凹陷點(diǎn)(pit)、鞍點(diǎn)(pass)和平地點(diǎn)(plane)等.利用DEM提取地形特征點(diǎn)，可通過一個(gè)3×3或更大的柵格窗口，通過中心網(wǎng)格與8個(gè)相鄰網(wǎng)格的高程關(guān)系來進(jìn)行判斷后獲取.

對DEM提取的地形特征點(diǎn)進(jìn)行數(shù)學(xué)分析，可以獲得相應(yīng)區(qū)域的坡度、坡向等參數(shù).坡度與坡向是互相聯(lián)系的兩個(gè)參數(shù)，坡度反映斜坡的傾斜程度，坡向反映斜坡所面對的方向.當(dāng)基于DEM計(jì)算坡度、坡向時(shí)，通常采用擬合曲面法求解[4]，ArcGIS的水文模塊也是采用的這種方法.

設(shè)擬合曲面的方程為z=f(x,y)，則其在點(diǎn)(x0,y0,z0)處的切平面方程為：

fx(x0,y0)(x-x0)+fy(x0,y0)(y-y0)-(z-z0)=0

其中fx(x0,y0)，fy(x0,y0)分別為函數(shù)z=f(x,y)對x,y的偏導(dǎo)數(shù).

該點(diǎn)的坡度為：

坡向?yàn)椋?/p>

除坡度、坡向的計(jì)算外，地形特征的提取還有表面積的計(jì)算、網(wǎng)點(diǎn)坡度變化率的計(jì)算以及地表粗糙度的計(jì)算等.其中ArcGIS在進(jìn)行表面積的計(jì)算時(shí)，是將DEM格網(wǎng)分解成三角網(wǎng)，使用海倫公式:

計(jì)算三角形的面積，其中a,b,c是三角形的三條邊的長度，p=(a+b+c)/2.

從而求得地表面積.

在計(jì)算任一網(wǎng)格點(diǎn)坡度變化率時(shí)，是取該網(wǎng)格點(diǎn)相鄰所有網(wǎng)格點(diǎn)坡度變化率中絕對值最大的一個(gè)網(wǎng)格點(diǎn)的值和它的符號.

網(wǎng)格表面的地表粗糙度：網(wǎng)格實(shí)際表面積S表與其水平投影面積S投之比，記為：

CZ=S表/S投.

2.2 水系特征的提取

ArcGIS中使用的DEM水系特征提取方法基于O’Callaghan與Marks(1984)提出的坡面徑流模擬方法，該方法的原理是先計(jì)算出各個(gè)網(wǎng)格單元上的匯水面積，然后結(jié)合匯水面積閾值的設(shè)定來判定河網(wǎng).該方法使用了水文學(xué)中匯流的概念來判斷地表徑流的路徑，而且能夠生成連續(xù)的流路，因此被認(rèn)為是水系提取中較好的一種方法.ArcGIS 將水系特征提取的相關(guān)命令集成在Hydrology 菜單中，使用時(shí)可直接調(diào)用相應(yīng)的函數(shù).

3 流域特征提取過程

3.1 洼地填充處理

洼地一般是由采樣誤差造成的，也有的洼地是地表形態(tài)的真實(shí)反映.因此，在洼地填充前要先計(jì)算洼地深度，由此來判斷洼地是由數(shù)據(jù)誤差造成的還是地表的真實(shí)反映，從而在進(jìn)行洼地填充時(shí)設(shè)置合理的閾值.

洼地填充的原理是：掃描每個(gè)網(wǎng)格時(shí)，比較該網(wǎng)格與相鄰的8個(gè)網(wǎng)格的高程，如果中心網(wǎng)格的高程比相鄰的8個(gè)網(wǎng)格的高程都低，則將相鄰8個(gè)網(wǎng)格中高程最低的那個(gè)網(wǎng)格的值賦予中心網(wǎng)格.

ArcGIS中洼地填充步驟如下：選擇Hydrology工具中的Fill工具，使用可選項(xiàng)Zlimit確定要被填平的洼地的閾值，步驟如下：

1)FlowDirection：以原始DEM數(shù)據(jù)作為輸入，計(jì)算初始流向;

2)Sink：以(1)的結(jié)果作為輸入，識別出所有的洼地(Sinks);

3)Watershed：找出洼地的影響區(qū)域;

4) Zonal Fill(位于ArcToolboxSpatial AnalystToolsonal下)：用于找出洼地最深深度;

5) ZonalStatistics(位于ArcToolboxSpatial Analyst Toolsonal下)：用于找出洼地最淺深度;

6)Minus(位于ArcToolbox3D Analyst ToolsRasterMath)下：求出最深與最淺洼地的差值.

進(jìn)行完以上操作后，即可完成洼地填充的處理，生成的數(shù)據(jù)稱為無凹陷柵格數(shù)據(jù)(Depressionless DEM).

圖2 水流方向編碼

圖3 匯流累積計(jì)算過程

3.2 水流方向計(jì)算

對于每個(gè)網(wǎng)格，水流方向是指水流離開此網(wǎng)格時(shí)的指向.確定流向的算法主要有單流向算法(Simple Flow Direction，SFD)和多流向算法( Multiple Flow Direction，MFD) 兩種， 其中單流向算法因其簡單方便而得到廣泛的應(yīng)用.單流向法假定一個(gè)網(wǎng)格中的水流只從一個(gè)方向流出網(wǎng)格，然后根據(jù)網(wǎng)格高程判斷水流方向.目前應(yīng)用最廣泛的單流向算法是D8算法.D8算法中，對水流離開每一網(wǎng)格的方向簡化為正東、東南、正南、西南、正西、西北、正北、東北8個(gè)方向，并依次用代碼1,2,4,8,16,32,64,128來表示，如圖2所示.無洼地的DEM上每一個(gè)網(wǎng)格的水流方向都是唯一的，一個(gè)網(wǎng)格的水流方向編碼取且僅取這8個(gè)數(shù)字中的1個(gè).

1)在ArcMap中用左鍵單擊ArcToolbox圖標(biāo)，啟動ArcToolbox；

2)打開水文分析模塊.啟動ArcToolbox，展開Analysis Tools工具箱，打開hydrology工具集；

3)雙擊Flow Direction工具，打開水流方向(Flow Direction)計(jì)算對話框.①Input surface data文本框中選擇輸入數(shù)據(jù)dem;②Output flow direction raster文本框中命名計(jì)算出來的水流方向文件名為flowdir，并選擇保存路徑;③在Force all edge cells to flow outward(Optional)前的復(fù)選框前打鉤，所有在DEM數(shù)據(jù)邊緣的網(wǎng)格的水流方向全部是流出DEM數(shù)據(jù)區(qū)域.默認(rèn)為不選擇;④drop raster輸出.drop raster是該網(wǎng)格在其水流方向上與其臨近的網(wǎng)格之間的高程差與距離的比值，以百分比的形式記錄.它反映了在整個(gè)區(qū)域中最大坡降的分布情況;⑤單擊OK，進(jìn)行水流方向計(jì)算過程.

3.3 匯集累積量計(jì)算

在地表徑流模擬過程中，匯流累積量是基于水流方向數(shù)據(jù)計(jì)算而來的.對每一個(gè)網(wǎng)格來說，其匯流累積量的大小代表著其上游有多少個(gè)網(wǎng)格的水流方向最終匯流經(jīng)過該網(wǎng)格，匯流累積的數(shù)值越大，該區(qū)域越易形成地表徑流.

一個(gè)流域的匯流累積量計(jì)算是提取該流域河流網(wǎng)絡(luò)的基礎(chǔ).流域內(nèi)一個(gè)網(wǎng)格的匯流量反映了其匯聚水流能力的強(qiáng)弱.一個(gè)網(wǎng)格的匯流累積量越大，其匯流能力也就越強(qiáng)，該網(wǎng)格所代表的地形特征就有可能是河谷；反之，匯流累積量為零的地方則可能代表流域的分水嶺.匯流累積量計(jì)算的基本思想是，假設(shè)以規(guī)則格網(wǎng)表示的DEM上每個(gè)網(wǎng)格有一個(gè)單位水量，按照水往低處流的自然規(guī)律，根據(jù)區(qū)域地形的水流方向數(shù)據(jù)計(jì)算每點(diǎn)處所流過的水量數(shù)值，從而得到該區(qū)域匯流累積量.

通過水流方向圖能夠計(jì)算出所有網(wǎng)格的匯流累積量，輸出的匯流累積量圖上每一個(gè)網(wǎng)格單元的值就代表著注入該網(wǎng)格單元的單位水量的數(shù)量，其值越高注入水流越多.由水流方向數(shù)據(jù)到匯流累積量計(jì)算的過程如圖3所示.

ArcGIS中匯流累積量計(jì)算步驟為：打開Fill Accumulation工具，以Flow Direction工具生成的柵格數(shù)據(jù)作為輸入，即可得到匯流累積量數(shù)據(jù).其中，可選項(xiàng)“Input Weight Raster”用于施加降雨量、蒸發(fā)量等外部因素的影響，默認(rèn)賦值為1.

3.4 對匯流設(shè)置閾值

在匯流累積圖中，每一個(gè)網(wǎng)格的匯流累積量代表著能夠注入該網(wǎng)格的所有單位水量的數(shù)量, 當(dāng)網(wǎng)格的匯流累積量大于某一給定的閾值時(shí),認(rèn)為該網(wǎng)格位于水道上.匯流閾值是河系網(wǎng)絡(luò)提取的關(guān)鍵參數(shù)，也是河道認(rèn)定的門檻.利用Hydrology中定義河系的命令可以設(shè)置匯流閾值.初始分析時(shí)，可選擇自小至大的幾個(gè)閾值，以便觀察最可能的河系的大致范圍.

設(shè)置閥值的步驟為：以Flow Accumulation生成的柵格數(shù)據(jù)作為輸入，輸出單值化后的柵格數(shù)據(jù).在Map Algebra Exp ression中填入con (即Condition)語句， 語句格式為：

con ( < condition >,< true_Expresssion >,{…} ,{false_Expression} )

其中{< false_expression>}為可選項(xiàng)，當(dāng){< false_expression>}為空且語句返回的值為假時(shí)，網(wǎng)格將被賦予默認(rèn)值NULL.

閾值的選擇視對待提取的河網(wǎng)的精度、疏密的要求而定.

3.5 河流網(wǎng)絡(luò)提取

以設(shè)定的匯流閾值為標(biāo)準(zhǔn),利用流徑處理命令從匯流累積柵格圖層中提取河系柵格圖層,凡是累積量或大于等于匯流閾值的柵格即被定義為河道.

ArcGIS中河流網(wǎng)絡(luò)提取步驟為：分別以Single OutputMap Algebra生成的柵格數(shù)據(jù)和Flow Direction生成的數(shù)據(jù)作為“Input Stream Raster”以及 “ Input Flow Direction Raster”的輸入，輸出的文件即是矢量化后的河網(wǎng)文件.

4 結(jié)語

ArcGIS中的Hydrology工具包提供了強(qiáng)大的水文信息提取功能，可以快速完成河網(wǎng)提取、流域劃分、水文信息統(tǒng)計(jì)、DEM數(shù)據(jù)的裁剪及拼接等任務(wù)[6-8].本文在分析DEM原理的基礎(chǔ)上，使用Hydrology工具包進(jìn)行了流域特征中的河網(wǎng)提取的實(shí)驗(yàn)，闡述了河網(wǎng)提取的一般過程.

[1] 陸守一.地理信息系統(tǒng)[M].北京:高等教育出版社,2004:182-187.

[2] 石軍南.DEM生成及其在森林公園規(guī)劃中的應(yīng)用[J].中南林業(yè)調(diào)查規(guī)劃,2001,20(1):36-41.

[3] 李志林,朱慶.數(shù)字高程模型[M].武漢:武漢大學(xué)出版社,2000.

[4] 郭仁忠.空間分析[M].武漢:武漢測繪科技大學(xué)出版社,2000:139-145.

[5] 孫慶艷,余新曉,胡淑萍，等.基于ArcGIS環(huán)境下DEM流域特征提取及應(yīng)用[J].北京林業(yè)大學(xué)學(xué)報(bào),2008,30(2):144-147.

[6] 張超,鄭鈞,張尚弘,等.ArcGis 9.0中基于DEM的水文信息提取方法[J].水利水電技術(shù),2005,36(11):1-4.

[7] 趙健,賈忠華,羅紈.ARCGIS環(huán)境下基于DEM的流域特征提取[J].水資源與水工程學(xué)報(bào),2006,17(1):74-76.

[8] 徐亞菲,李向新,賴金富,等.基于DEM和ArcGIS的水文信息提取方法研究[J].科技情報(bào)開發(fā)與經(jīng)濟(jì),2008,18(6):135-139.