GIS环境下基于GDEM的丹江口库区流域自动提取方法的研究

李 辉，黄进良，王立辉，赵国松

（1.中国科学院 测量与地球物理研究所，武汉 430077；2.中国科学院 研究生院，北京 100049）

丹江口水库位于鄂西北与豫西南的交界处，是我国正在实施的南水北调工程中线的水源地.丹江口库区主要包括湖北省十堰市所辖的张湾区、茅箭区、丹江口市、郧县和郧西县以及河南省南阳市所辖的淅川县和西峡县，共7个县（市、区）.库区生态环境质量的好坏对于南水北调中线工程的安全运营和充分发挥效益具有决定性的意义.

以某个流域为研究区域，提取该流域的水文特征，得到该流域的次一级流域即为小流域.以小流域为单元，范围比较小，同时又是一个相对完整的流域，研究流域内的生态环境治理和建设，既节省研究的人力和物力，又保持了其科学性，同时又可以推广到整个库区的流域生态环境的研究中.

流域特征提取的方法包括确定流向的单向流法和多向流法；提取河网的识别谷点法和基于流向提取河网的方法等［1］.GIS技术将一些方法应用到了软件中.“3S”技术的引进，为水文科学注入了新的血液［2］.GIS在数据管理、空间分析及可视性方面的功能［3］，已经使目前水文模型研究的重点从流域统计模型转到GIS与成熟水文模型相结合的分布式水文模型的研究方面.美国环境系统研究所为ArcGIS开发了一个水文模块，用于河网提取和分析.国内外运用此模块已经对提取水系研究地表水文特征等进行了大量工作［4－7］.本文使用30m分辨率的DEM，通过模拟水流方向，汇流分析，河网生成，出水口的确定子集边界的划分等对库区的水文信息进行了提取.

1 DEM数据来源

数字高程模型（DEM）是目前用于流域地形分析的主要数据.基于DEM提取流域的数字特征有多种比较成熟的算法［8］.DEM数据中包含了丰富的地形、地貌、水文信息，它能够反映各种分辨率的地形特征，可以提取大量的地表形态信息，如流域网格单元的坡向、坡度以及单元格之间的关系等［9－10］.

DEM获取有3种方法.第一，航空摄影测量采用的比较广泛.便于快捷的获取大面积的数据，有利于保持数据的现势性.第二，利用GPS和全站仪等测量仪器在实地采集以获取数据.多用于有限范围的大比例尺DEM建模.工作量较大.第三，地形图上数据采集.此种方法的特点是信息丰富，价格低，效率高.

本文使用的DEM数据来源于日本METI和美国NASA合作，使用ASTER数据生成的全球数字高程模型产品GDEM.ASTER GDEM基本的单元按1°×1°分片，空间分辨率1弧度秒（约30 m），垂直精度20m，水平精度30m.

2 ArcGIS的流域提取方法

基于ArcGIS的流域自动提取的命令在水文分析工具里.根据自动提取的原理，要按照一定的程序进行流域提取.图1简要描述了流域提取的基本流程.

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

2.1 数据的预处理

流域自动提取的有关命令集中在Hydrology菜单中，使用方便.由于水流在地势较低的区域的不连续性以及在平地区域的不确定性和随机性等原因，刚获得的DEM数据，一般会有洼地（Sink）和尖峰（Pink）.洼地是指某个单元格的高程低于相邻8个单元格的高程的情况；尖峰则是某个单元格的高程高于和它相邻的8个单元格的高程的情况.而这两种情况会影响到流域提取的质量，洼地单元格会导致水流可以流入却不能流出，尖峰单元格会导致水流可以流出却不能流入，因此需要对DEM首先进行填洼和削峰的处理.使用Hydroloy工具包中的Fill命令.计算机扫描每个单元格，如果是洼地，那么该洼地单元格的高程会被赋值为相邻8个单元格中高程最低的那个高程；遇到尖峰也会用临近的单元格的高程来代替.用Fill之后的数据与获得的未经处理的数据进行相减，就可以得到填洼和削峰的位置和量［11］.图2为填洼和削峰之后的DEM图.

另外，使用DEM进行流域提取还必须进行平坦区域的处理.平坦区域是指DEM本身存在的平坦区域和经过洼地填平后形成的平坦区域.Hydrology的D8算法，在平坦区域内无法产生，这种算法是通过将平地两端边缘的水流汇流栅格点连接起来生成与实际不符的平直河道.而平坦区域水流的流动方向有随机性，而自然水系通常是弯曲而成不规则的形状，于是有些栅格单元的水流方向就是多流向的，D8算法是单流向的，因此就无法处理平坦区域.处理这种情况通常采用高程增量迭加的算法来设定流向.实际中就是将平坦区域的单元格的高程进行微调（增高），生成相对合理的汇流水系.

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

2.2 流域自动提取

2.2.1 确定栅格单元水流方向 首先对DEM进行流向分析.水流方向是指栅格单元内水流流出该单元格的方向.水流总是从地势高处向地势低处流动，中间经过下渗，截留，蒸发，散发等过程，最后通过流域出水口排出流域.要提取流域，首先得确定水流在每个栅格单元里的流动方向.采用FlowDirection函数来确定水流方向.此命令将每个栅格单元赋予1个1～255的整数值来表示水流方向.单元格周围有8个相邻的单元格.ArcGIS是通过比较中心单元格与周围8个相邻单元格的高程大小，然后将高程下降最多的单元格方向作为该中心单元格的流向的.丹江口库区流域水流方向见图3.

2.2.2 汇流能力分析 汇流分析是利用FlowAccumulation函数来计算流入该单元格的上游单元格的总数目来得到汇流能力的栅格图.算法从每个单元格出发扫描水流方向矩阵，也就是FlowDirection命令的运行结果.沿水流方向一直追踪到DEM的边界.可以理解为假设在流域内每个单元格降水一个单位，然后按照水流方向向下移动，每经过一个单元格就使其积累流量增加一个单位，由此每个单元格都可以计算出其所累积的上游的流量值.流量值即为积累上游的单元格数，乘以单元格面积就可以得到每个单元格的汇流面积.汇流分析可以确定河流的网络，然后在河流网络的基础上确定出流域的边界.

图3 丹江口库区流域水流方向Fig.3 Flow direction of the drainage in Danjiangkou Reservoir

2.2.3 水流网络提取 对流向的积累栅格设置阈值来提取河网，采用Map Algebra里的Single Output Map Algebra命令，以FlowAccumulation生成的结果作为输入，输出单值化的栅格数据.在Map Algebea expression中填入Condition语句，语句 格 式 为 Con（〈congdition〉，〈true＿expression〉，〈false＿expression〉），其中，〈false＿expression〉是可选项.它为空时，如果Condition语句返回值为假，单元格将被赋值为NULL.

例如，Con（［flow＿accumulation］＞Const，1，0）表示将栅格数据flow＿accumulation中汇流值大于某个常数const的单元格赋值为1，其余的单元格则赋值为0.汇流阈值是河网提取的关键因子［12］.可以尝试几个不同的阈值来确定最终的阈值.以阈值为标准，利用流经处理命令，在流量累积数据中提取河网图层，凡是累积流量大于所设置的阈值的即被定义为河道.提取的河流网络见图4，其中图例1代表提取出的河道.

图4 丹江口库区河网提取Fig.4 Extraction of river network in Danjiangkou Reservoir

2.2.4 子流域提取 首先，上一步骤得到了河网图层，但是河道的标识都为1.所以先采用Stream Link命令对河网进行编号.输入已生成的河网数据和FlowDirection生成的水流方向数据，输出得到河网河道记有编号的栅格数据.

然后，以有河道编号的栅格数据和FlowDirection生成的水流方向数据为输入数据，用 Watershed命令，输出划分了子流域的栅格数据图层.丹江口库区子流域提取见图5.

图5 丹江口库区子流域提取Fig.5 Extraction of watersheds in Danjiangkou Reservoir

3 结果与分析

按照前面描述的方法，将DEM数据经过填洼，削峰，平地抬升，水流方向的确定，汇流能力分析，流域网络的提取，子流域提取等步骤，最后输出了流域提取的栅格数据.从提取的流域结果来看，通过生成的子流域图层与DEM叠加显示，在山地丘陵地区，流域界线与分水岭非常吻合.但是平坦区域由于等高线稀疏，算法并不考虑随机因素对河流的影响，提取的水系平直，与实际弯曲的水系有一定的误差.但总体上本文的研究方法是可行的.

需要指出，本文的研究区域是湖北省十堰市及河南省南阳市所辖的七个完整的行政县（市、区），而行政区划的边界并不是按照小流域的边界划分的，因此在研究区域的周边，有很小的部分区域没有包括在流域提取的范围内.

4 结论与讨论

ArcGIS的Hyrology工具提供了强大的水文信息提取功能，操作界面比较直观，可以根据需要创建工具包应用于其它研究，灵活方便.

本文使用30m分辨率的DEM数据对丹江口库区进行了流域提取.结果表明提取的河网同实际情况相比误差较小，基本符合实际的水网分布.同时由DEM直接提取河网，能够大大减少人力和物力的消耗.

小流域是一个相对完整的生态过程单元，以小流域为研究单元，可以从中观察认识库区的自然分异规律，进而推广应用于库区流域管理、水资源高效利用及生态建设与环境保护.

同时，本文所介绍的方法是基于地表径流模型，仅考虑地形因素对水系的影响.实际中，人为的堤坝、道路等都对水系产生影响，且DEM的精度、大小、范围、所描述的流域及其本身的特征都直接影响提取效果［13］.在河网提取过程中需要设置一个汇流阈值，根据这个值的确定可以提取出不同级别的子流域.用何种方法确定汇流阈值来提取河道，还有待进一步的研究.

［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 的分布式水文模式构建方法［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栅格系统的集水农业地表产流模拟分析［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］李 翀，杨大文.基于栅格数字高程模型DEM的河网提取及实现［J］.中国水利水电科学研究院学报，2004，2（3）：208－214.

［11］孙庆艳，余新晓，胡淑萍，等.基于ArcGIS环境下的DEM流域特征提取及应用［J］.北京林业大学学报，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.