GIS環(huán)境下基于GDEM的丹江口庫區(qū)流域自動提取方法的研究

李 輝，黃進良，王立輝，趙國松

（1.中國科學院 測量與地球物理研究所，武漢 430077；2.中國科學院 研究生院，北京 100049）

丹江口水庫位于鄂西北與豫西南的交界處，是我國正在實施的南水北調(diào)工程中線的水源地.丹江口庫區(qū)主要包括湖北省十堰市所轄的張灣區(qū)、茅箭區(qū)、丹江口市、鄖縣和鄖西縣以及河南省南陽市所轄的淅川縣和西峽縣，共7個縣（市、區(qū)）.庫區(qū)生態(tài)環(huán)境質(zhì)量的好壞對于南水北調(diào)中線工程的安全運營和充分發(fā)揮效益具有決定性的意義.

以某個流域為研究區(qū)域，提取該流域的水文特征，得到該流域的次一級流域即為小流域.以小流域為單元，范圍比較小，同時又是一個相對完整的流域，研究流域內(nèi)的生態(tài)環(huán)境治理和建設(shè)，既節(jié)省研究的人力和物力，又保持了其科學性，同時又可以推廣到整個庫區(qū)的流域生態(tài)環(huán)境的研究中.

流域特征提取的方法包括確定流向的單向流法和多向流法；提取河網(wǎng)的識別谷點法和基于流向提取河網(wǎng)的方法等［1］.GIS技術(shù)將一些方法應(yīng)用到了軟件中.“3S”技術(shù)的引進，為水文科學注入了新的血液［2］.GIS在數(shù)據(jù)管理、空間分析及可視性方面的功能［3］，已經(jīng)使目前水文模型研究的重點從流域統(tǒng)計模型轉(zhuǎn)到GIS與成熟水文模型相結(jié)合的分布式水文模型的研究方面.美國環(huán)境系統(tǒng)研究所為ArcGIS開發(fā)了一個水文模塊，用于河網(wǎng)提取和分析.國內(nèi)外運用此模塊已經(jīng)對提取水系研究地表水文特征等進行了大量工作［4－7］.本文使用30m分辨率的DEM，通過模擬水流方向，匯流分析，河網(wǎng)生成，出水口的確定子集邊界的劃分等對庫區(qū)的水文信息進行了提取.

1 DEM數(shù)據(jù)來源

數(shù)字高程模型（DEM）是目前用于流域地形分析的主要數(shù)據(jù).基于DEM提取流域的數(shù)字特征有多種比較成熟的算法［8］.DEM數(shù)據(jù)中包含了豐富的地形、地貌、水文信息，它能夠反映各種分辨率的地形特征，可以提取大量的地表形態(tài)信息，如流域網(wǎng)格單元的坡向、坡度以及單元格之間的關(guān)系等［9－10］.

DEM獲取有3種方法.第一，航空攝影測量采用的比較廣泛.便于快捷的獲取大面積的數(shù)據(jù)，有利于保持數(shù)據(jù)的現(xiàn)勢性.第二，利用GPS和全站儀等測量儀器在實地采集以獲取數(shù)據(jù).多用于有限范圍的大比例尺DEM建模.工作量較大.第三，地形圖上數(shù)據(jù)采集.此種方法的特點是信息豐富，價格低，效率高.

本文使用的DEM數(shù)據(jù)來源于日本METI和美國NASA合作，使用ASTER數(shù)據(jù)生成的全球數(shù)字高程模型產(chǎn)品GDEM.ASTER GDEM基本的單元按1°×1°分片，空間分辨率1弧度秒（約30 m），垂直精度20m，水平精度30m.

2 ArcGIS的流域提取方法

基于ArcGIS的流域自動提取的命令在水文分析工具里.根據(jù)自動提取的原理，要按照一定的程序進行流域提取.圖1簡要描述了流域提取的基本流程.

圖1 基于DEM的流域自動提取流程圖Fig.1 Process of automated extraction of drainages based on DEM

2.1 數(shù)據(jù)的預(yù)處理

流域自動提取的有關(guān)命令集中在Hydrology菜單中，使用方便.由于水流在地勢較低的區(qū)域的不連續(xù)性以及在平地區(qū)域的不確定性和隨機性等原因，剛獲得的DEM數(shù)據(jù)，一般會有洼地（Sink）和尖峰（Pink）.洼地是指某個單元格的高程低于相鄰8個單元格的高程的情況；尖峰則是某個單元格的高程高于和它相鄰的8個單元格的高程的情況.而這兩種情況會影響到流域提取的質(zhì)量，洼地單元格會導致水流可以流入?yún)s不能流出，尖峰單元格會導致水流可以流出卻不能流入，因此需要對DEM首先進行填洼和削峰的處理.使用Hydroloy工具包中的Fill命令.計算機掃描每個單元格，如果是洼地，那么該洼地單元格的高程會被賦值為相鄰8個單元格中高程最低的那個高程；遇到尖峰也會用臨近的單元格的高程來代替.用Fill之后的數(shù)據(jù)與獲得的未經(jīng)處理的數(shù)據(jù)進行相減，就可以得到填洼和削峰的位置和量［11］.圖2為填洼和削峰之后的DEM圖.

另外，使用DEM進行流域提取還必須進行平坦區(qū)域的處理.平坦區(qū)域是指DEM本身存在的平坦區(qū)域和經(jīng)過洼地填平后形成的平坦區(qū)域.Hydrology的D8算法，在平坦區(qū)域內(nèi)無法產(chǎn)生，這種算法是通過將平地兩端邊緣的水流匯流柵格點連接起來生成與實際不符的平直河道.而平坦區(qū)域水流的流動方向有隨機性，而自然水系通常是彎曲而成不規(guī)則的形狀，于是有些柵格單元的水流方向就是多流向的，D8算法是單流向的，因此就無法處理平坦區(qū)域.處理這種情況通常采用高程增量迭加的算法來設(shè)定流向.實際中就是將平坦區(qū)域的單元格的高程進行微調(diào)（增高），生成相對合理的匯流水系.

圖2 填洼和削峰后的DEMFig.2 DEM without sinks and peaks

2.2 流域自動提取

2.2.1 確定柵格單元水流方向 首先對DEM進行流向分析.水流方向是指柵格單元內(nèi)水流流出該單元格的方向.水流總是從地勢高處向地勢低處流動，中間經(jīng)過下滲，截留，蒸發(fā)，散發(fā)等過程，最后通過流域出水口排出流域.要提取流域，首先得確定水流在每個柵格單元里的流動方向.采用FlowDirection函數(shù)來確定水流方向.此命令將每個柵格單元賦予1個1～255的整數(shù)值來表示水流方向.單元格周圍有8個相鄰的單元格.ArcGIS是通過比較中心單元格與周圍8個相鄰單元格的高程大小，然后將高程下降最多的單元格方向作為該中心單元格的流向的.丹江口庫區(qū)流域水流方向見圖3.

2.2.2 匯流能力分析 匯流分析是利用FlowAccumulation函數(shù)來計算流入該單元格的上游單元格的總數(shù)目來得到匯流能力的柵格圖.算法從每個單元格出發(fā)掃描水流方向矩陣，也就是FlowDirection命令的運行結(jié)果.沿水流方向一直追蹤到DEM的邊界.可以理解為假設(shè)在流域內(nèi)每個單元格降水一個單位，然后按照水流方向向下移動，每經(jīng)過一個單元格就使其積累流量增加一個單位，由此每個單元格都可以計算出其所累積的上游的流量值.流量值即為積累上游的單元格數(shù)，乘以單元格面積就可以得到每個單元格的匯流面積.匯流分析可以確定河流的網(wǎng)絡(luò)，然后在河流網(wǎng)絡(luò)的基礎(chǔ)上確定出流域的邊界.

圖3 丹江口庫區(qū)流域水流方向Fig.3 Flow direction of the drainage in Danjiangkou Reservoir

2.2.3 水流網(wǎng)絡(luò)提取 對流向的積累柵格設(shè)置閾值來提取河網(wǎng)，采用Map Algebra里的Single Output Map Algebra命令，以FlowAccumulation生成的結(jié)果作為輸入，輸出單值化的柵格數(shù)據(jù).在Map Algebea expression中填入Condition語句，語句 格 式 為 Con（〈congdition〉，〈true＿expression〉，〈false＿expression〉），其中，〈false＿expression〉是可選項.它為空時，如果Condition語句返回值為假，單元格將被賦值為NULL.

例如，Con（［flow＿accumulation］＞Const，1，0）表示將柵格數(shù)據(jù)flow＿accumulation中匯流值大于某個常數(shù)const的單元格賦值為1，其余的單元格則賦值為0.匯流閾值是河網(wǎng)提取的關(guān)鍵因子［12］.可以嘗試幾個不同的閾值來確定最終的閾值.以閾值為標準，利用流經(jīng)處理命令，在流量累積數(shù)據(jù)中提取河網(wǎng)圖層，凡是累積流量大于所設(shè)置的閾值的即被定義為河道.提取的河流網(wǎng)絡(luò)見圖4，其中圖例1代表提取出的河道.

圖4 丹江口庫區(qū)河網(wǎng)提取Fig.4 Extraction of river network in Danjiangkou Reservoir

2.2.4 子流域提取 首先，上一步驟得到了河網(wǎng)圖層，但是河道的標識都為1.所以先采用Stream Link命令對河網(wǎng)進行編號.輸入已生成的河網(wǎng)數(shù)據(jù)和FlowDirection生成的水流方向數(shù)據(jù)，輸出得到河網(wǎng)河道記有編號的柵格數(shù)據(jù).

然后，以有河道編號的柵格數(shù)據(jù)和FlowDirection生成的水流方向數(shù)據(jù)為輸入數(shù)據(jù)，用 Watershed命令，輸出劃分了子流域的柵格數(shù)據(jù)圖層.丹江口庫區(qū)子流域提取見圖5.

圖5 丹江口庫區(qū)子流域提取Fig.5 Extraction of watersheds in Danjiangkou Reservoir

3 結(jié)果與分析

按照前面描述的方法，將DEM數(shù)據(jù)經(jīng)過填洼，削峰，平地抬升，水流方向的確定，匯流能力分析，流域網(wǎng)絡(luò)的提取，子流域提取等步驟，最后輸出了流域提取的柵格數(shù)據(jù).從提取的流域結(jié)果來看，通過生成的子流域圖層與DEM疊加顯示，在山地丘陵地區(qū)，流域界線與分水嶺非常吻合.但是平坦區(qū)域由于等高線稀疏，算法并不考慮隨機因素對河流的影響，提取的水系平直，與實際彎曲的水系有一定的誤差.但總體上本文的研究方法是可行的.

需要指出，本文的研究區(qū)域是湖北省十堰市及河南省南陽市所轄的七個完整的行政縣（市、區(qū)），而行政區(qū)劃的邊界并不是按照小流域的邊界劃分的，因此在研究區(qū)域的周邊，有很小的部分區(qū)域沒有包括在流域提取的范圍內(nèi).

4 結(jié)論與討論

ArcGIS的Hyrology工具提供了強大的水文信息提取功能，操作界面比較直觀，可以根據(jù)需要創(chuàng)建工具包應(yīng)用于其它研究，靈活方便.

本文使用30m分辨率的DEM數(shù)據(jù)對丹江口庫區(qū)進行了流域提取.結(jié)果表明提取的河網(wǎng)同實際情況相比誤差較小，基本符合實際的水網(wǎng)分布.同時由DEM直接提取河網(wǎng)，能夠大大減少人力和物力的消耗.

小流域是一個相對完整的生態(tài)過程單元，以小流域為研究單元，可以從中觀察認識庫區(qū)的自然分異規(guī)律，進而推廣應(yīng)用于庫區(qū)流域管理、水資源高效利用及生態(tài)建設(shè)與環(huán)境保護.

同時，本文所介紹的方法是基于地表徑流模型，僅考慮地形因素對水系的影響.實際中，人為的堤壩、道路等都對水系產(chǎn)生影響，且DEM的精度、大小、范圍、所描述的流域及其本身的特征都直接影響提取效果［13］.在河網(wǎng)提取過程中需要設(shè)置一個匯流閾值，根據(jù)這個值的確定可以提取出不同級別的子流域.用何種方法確定匯流閾值來提取河道，還有待進一步的研究.

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

［2］Kwan T L.Generating design hydrographs by DEM assisted geomorphic runoff simulation：a case study［J］.Journal of the American Water Resources Association，1998，34（2）：375－384.

［3］王中根，劉昌明，左七亭，等.基于DEM 的分布式水文模式構(gòu)建方法［J］.地理科學進展，2002，21（5）：432－439.

［4］Jenson S K.Application of hydrologic information automatically extracted from digital elevation models［J］.Hydrological process，1991（5）：31－44.

［5］Tarboton D G，Bras R L，Rodriguez I I.On the extraction of channel network from digital elevation data［J］.Hydrological Processes，1991（5）：81－100.

［6］梁天剛，張勝雷，戴若蘭，等.基于GIS柵格系統(tǒng)的集水農(nóng)業(yè)地表產(chǎn)流模擬分析［J］.水利學報，1998（7）：26－29

［7］陳永良，劉大有，虞強源.從DEM中自動提取自然水系［J］.中國圖象圖形學報，2002，7（1）：91－96.

［8］TRIBE A.Automated recognition of valley lines and drainage networks from grid digital elevation models：A review and a new method［J］.Journal of Hydrology，1992，139 （1／4）：263－293.

［9］Band L E.Topographical partition of watersheds with digital elevation models［J］.Water Resources Research，1986（2）：15－24.

［10］李 翀，楊大文.基于柵格數(shù)字高程模型DEM的河網(wǎng)提取及實現(xiàn)［J］.中國水利水電科學研究院學報，2004，2（3）：208－214.

［11］孫慶艷，余新曉，胡淑萍，等.基于ArcGIS環(huán)境下的DEM流域特征提取及應(yīng)用［J］.北京林業(yè)大學學報，2008，30（2）：144－147.

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

［13］Carbrecht Martzlw.Digital Elevation Model Issues in Water Resource Modeling［C］.In Proceedings of 19th ESRI International User Conference.San Diego Califomia，1999：13－22.