基于梯田DEM水流路径模拟算法的比较

冯园,赵牡丹*,冯瑶,王翊人,张鹏

1.西北大学城市与环境学院,陕西西安710100

2.西安航天天绘数据技术有限公司,陕西西安710100

基于梯田DEM水流路径模拟算法的比较

冯园1,赵牡丹1*,冯瑶2,王翊人1,张鹏1

1.西北大学城市与环境学院,陕西西安710100

2.西安航天天绘数据技术有限公司,陕西西安710100

大面积梯田的修筑，极大地改变了地表微形态，使地表坡度、径流路径和汇流方式发生改变，进而影响到土壤侵蚀强度。基于5 m分辨率的DEM数据，利用真实田坎信息构建梯田DEM，采用AGREE算法对梯田DEM进行预处理，通过比较分析不同集水面积阈值提取的梯田DEM沟道信息，初步探索合适的填洼算法和流向算法应用于梯田DEM的水流路径模拟。结果表明：利用M&G算法填洼处理，MFD-md算法确定流向，能够更好的模拟梯田区域的水流路径。

梯田DEM;AGREE算法;路径模拟;集水面积阈值

DEM是进行数字地形分析的基础数据，从中可以提取坡度、坡长、地表粗糙度等因子[1,2]。近年来，DEM在水土保持、生态学应用、水文分析等诸多领域得到了较广泛的应用[3-6]。但目前，我国虽然已经完成了从1:1万到1:100万比例尺地形数据库的建设，然而各类DEM均只反映连续光滑的自然地表，不能表达诸如田坎、梯田这类普遍存在的自然或者人工突变地形信息[7]，而在黄土高原水土流失区，经多年的治理和建设，大面积梯田的修筑使得地表坡度减缓、坡长缩短，起到了很好蓄水保土的作用[8]。同时，地表微形态的改变也使得水流的汇流路径和汇流方式都发生改变，因而对于嵌入梯田信息的DEM（本研究称为梯田DEM）进行水流路径模拟研究是十分必要的。

将基于高精度遥感影像，利用真实田坎信息构建梯田DEM，通过AGREE算法对梯田DEM进行预处理，探索分析不同填洼算法和流向算法对梯田DEM水流路径模拟的差异性。进而完善梯田DEM研究的相关理论以及方法，也为基于过程的土壤侵蚀定量研究奠定基础。

1 研究区概况及主要数据源

1.1 研究区域概况

遵循典型性、科学性原则，选取陕西省黄龙县县城西北部的梯田区域作为实验样区。黄龙县位于陕西省延安市东南部，处于黄土高原丘陵沟壑区，属于陕西秦岭以北大陆性半湿润气候区。地处东经109°38′49″～110°16′49″，北纬35°24′09″～36°02′11″之间，南北长约64 km，县城西北部梯田分布广泛。

1.2 主要数据源

研究所用数据包括：（1）国家基础比例尺1:1万DLG数据，其中包括等高线、高程点数据层以及水系数据层；（2）1:1万DEM数据；（3）2.5 m SPOT高精度遥感影像数据；（4）GPS RTK野外实测的田面高程点数据。研究区域由三块完整梯田组成的小流域。

2 研究方法

2.1 梯田DEM的建立

为了构建出的梯田更加逼近真实情况，选用李慧[9]提出的基于真实田坎信息构建梯田DEM的方法。梯田具有独特的剖面形态特征以及适合构建数学模型的几何特征[10]，因而构建梯田DEM之前，须对梯田的数学模型及相关参数进行简单介绍（如图1）。

图1 梯田断面示意图Fig.1 The section diagram of terrace

设W为梯田坡面宽，L为梯田田面水平投影即梯田田面宽，θ为梯田对应的原始地表坡度，δ为内侧田坎坡度，l为梯田田坎水平投影宽，H为梯田两相邻田面高程之差，h为田坎高程差即田坎高。由此可以得到关系式：

由（1）式子变换可得到：

由此可知当θ=0，梯田的类型为水平梯田，由水平梯田的三角函数关系可知：

研究样区是一个完整的小流域（本研究的流域面积为1.5 km2），有典型的水平梯田且现有的1:1万DEM数据可以覆盖此流域。经筛选后最终确定的样区内有3块完整的梯田，分别命名为梯田1，梯田2，梯田3。整体的构建思路可归纳为：首先，在高分辨率遥感影像上提取梯田的台沿线，根据梯田构建的数学模型，求出台沿线偏移量L，其次再根据两组赋有高程值的特征线（即台沿线和台沿偏移线）作为内插属性，构建TIN，再将TIN转换成栅格DEM，最终得到梯田DEM（如图2）。

图2 梯田DEM光照渲染图(a:梯田1；b:梯田2；c：梯田3)Fig.2 Illumination rendering of terraces DEM(a:terrace 1;b:terrace 2;c:terrace 3)

2.2 梯田DEM预处理

2.2.1 DEM数据前期处理对于整个流域而言，由于构建好的梯田DEM均为1 m分辨率，而非梯田部分为5 m分辨率，所以需要对非梯田部分进行重采样处理，然后对其进行镶嵌处理。为了更好的显示镶嵌效果，使用山体阴影效果显示整个研究流域（如图3）。

图3 嵌入梯田后的研究区效果Fig.3 Embedded terraces in the study area

2.2.2 AGREE算法修复梯田DEM流域基于真实田坎信息构建的梯田DEM，尤其是梯田类型为水平梯田时，由于方法自身局限性，构建出的梯田田面近似于水平，无法反映出降雨时水流的汇流方式和汇流路径，所以需要对构建好的梯田DEM进行一定的预处理，保证水平梯田田面水流能够汇流并形成连贯的流路。AGREE算法[11]可以根据输入的流路信息修复DEM数据，使原始DEM数据的高程值在输入河道位置得以修改，为后续提取精确的流路信息奠定基础。基于上述思想，本研究首先根据野外实地测量的梯田田面水流沟道点作为基础数据，然后将RTK点转换成向量流路线，进而得到完整的流路；最终将得到的流路作为修复梯田DEM的向量资料，参数设置为：Stream buffer为1，Smooth drop/raise为1，Sharp drop/raise为100。

2.3 梯田DEM水流路径模拟

2.3.1 洼地填充DEM中存在洼地、伪洼地和平原区，有些是真实地形实际存在的，但大多数是由于数据采集和处理时造成的。洼地的存在往往会导致水流无法到达集水区边界，使得最后得到的水流路径不连续。且经过AGREE算法预处理后的DEM，很有可能产生新的洼地，因而须进行洼地填充，得到校正后的梯田DEM。目前，国内外学者已提出很多种填洼算法，本研究采用的是J&D[12]算法和M&G[13]算法。J&D算法的基本思想是先对洼地进行填平处理，然后逐步垫高平地，最终使流域内所有格网中的水都能从出口流出；M&G算法对于洼地，在预处理中就将其填为平坦区，对于平坦区，利用一定的算法对平坦区域内的高程进行特殊的修正，使其能尽可能贴合真实的地貌地形。

2.3.2 水流流向确定水流方向的确定是提取流域水系的关键，它决定了地表径流的方向以及栅格单元间流量的分配。流向算法根据其模型中每个像元点的流向数量可以分为两类：单流向算法（Single Flow Direction algorithm，SFD）和多流向算法（Multiple Flow Direction algorithm，MFD）。本研究采用了最为经典的单流向算法—D8算法[14]以及秦承志提出的一种新的多流向算法—MFD-md算法[15]。D8算法是假设每一栅格单元只有8种可能的流向并分别比较中心格网与其邻近八个格网之间高程值，高程值差距最大者则被认定为水流流向；MFD-md是算法假设水流分配随着下坡坡度变化，同时保留了单流向算法和多流向算法的优势，更为合理地模拟了地表微地形的变化对水流分配的影响。

3 结果与分析

为了反映水流路径的模拟效果，并能分析和比较不同算法模拟水流汇流，形成流路的差别，本研究将通过横向比对即对在相同集水面积阈值下不同算法提取的流路信息的丰富性，以及纵向比对即通过不同算法在各自合理集水面积阈值下所提取的流路信息的真实性，以此比较不同算法在实际应用中的效果。

3.1 横向比较分析

图4、图5分别为D8算法和MFD-md算法在阈值为500、1000及阈值为2000时，梯田DEM提取的水流路径。从图4和图5均可以看出，随着集水面积阈值的增大，提取的流路也变得稀疏，在梯田田面上表象尤为明显。图4与图5（即D8算法提取结果与MFD-md算法提取结果）相比较，MFD-md算法生成的流路更为曲折，与真实流路更为契合，并且在相同阈值条件下，D8算法和MFD-md算法提取的流路密度是不同的，MFD-md算法会生成更多短小河道，也就说明MFD-md算法在相同的集水阈值条件下能更丰富表达梯田田面的流路信息。

图4 D8算法不同阈值得到的水流路径Fig.4 The different flow paths from different threshold values of D8 algorithm

图5 MFD-md算法不同阈值得到的水流路径Fig.5 The different flow paths from different threshold values of MFD-md algorithm

3.2 纵向比较分析

以GPS RTK实测得到的梯田田面沟道点为基础，经处理后得到的本研究的“真”值水系线，其总河长为25.6173 km，流域面积为1.5 km2，所以其河网密度=25.6173/1.5=170.38（1/（10-1km））。对于D8算法提取的水流路径，集水面积阈值=400时，其河网密度为170.22（1/10-1km），对于MFD-md算法提取的水流路径，集水面积阈值=1800时，其河网密度为170.46（1/10-1km）。两者与水系线计算出河网密度=170.38（1/（10-1km））最为接近，所以可以认为利用D8算法，集水面积阈值选取400时，提取的水流路径较为合理；利用MFD-md算法，集水面积阈值选取1800时，提取的水流路径较为合理（表1）。

表1 不同算法提取河网密度的统计结果Table1 Statistical results of drainage network density extracted by different algorithms

将集水面积阈值为400时，D8算法提取的水流路径与“真”值水流路径叠合显示，D8算法提取结果用灰色线条置于顶层表示，“真”值水流路径置于底层用黑色线条表示，得到图6。将集水面积阈值为1800时，MFD-md算法提取的水流路径与“真”值水流路径叠合显示，MFD-md算法提取结果用灰色线条置于顶层表示，“真”值水流路径置于底层用黑色线条表示，得到图7。

图6 阈值=400时，D8提取流路与“真”实流路对比Fig.6 The comparison between D8 and"true"flow path when threshold is 400

图7 阈值=1800时，MFD-md提取流路与“真”实流路对比Fig.7 The comparison between MFD-md and"true"flow path when threshold is 1800

由图6可以看出，D8算法的提取结果与真实河道整体上基本相吻合，保留了大部分真实的田面流路信息。但在河道汇集繁杂的地方，与真实河道出入较大；并且在梯田田面上细小流路的模拟效果不是非常理想，会出现“短缺”现象，尤其是在水流汇集的起始点处“损失”较为严重。由图7可以看出，MFD-md算法的提取结果，整体上看，基本的流路与真实河道相吻合，几乎保留了真实河道流路的全部信息。在河道汇集繁杂的地方，与真实河道也更为契合；在梯田田面上细小流路模拟效果比D8算法模拟的结果要好，虽然在水流汇集的起始点出现“短缺”现象，但“误差”在可以接受的范围内。

4 结论与讨论

基于真实田坎信息构建梯田DEM时，对田坎信息的提取与高程赋值的方法比原有方法更合理，因而能够更精确的表达梯田DEM信息；采用AGREE算法，通过Buffer、Smooth、Sharp三个参数对水平梯田DEM田面上的高程值进行预处理，为处理梯田DEM田面因“平整”水流无法汇集和流出的难题，提供了一种可行的新方法；通过不同集水面积阈值条件下提取梯田DEM沟道信息，分析比较得出：经过M&G算法填洼处理后，再采用MFD-md算法计算流向提取的流路结果要优于J&D算法填洼处理后的D8算法，结果更接近“真实”流路。

此外，还有诸多地方待以改进，首先集水面积阈值的准确性完全取决于水系线的精度，实验结果对数据精度的依赖性太强；其次对于梯田DEM的洼地与平坦区处理，及流向确定算法，研究中只探讨了J&D算法填洼后，运用D8算法确定流向模拟结果，和M&G算法填洼后，运用MFD-md算法确定流向的模拟结果，并未做交叉分析比较。在今后的研究中应多作对比与改进，以便探寻更好的模拟水流路径的方法，完善梯田DEM水流路径模拟研究的理论体系，为后续梯田DEM在土壤侵蚀方面的研究提供更好的理论与技术支撑。

[1]李志林,朱庆.数字高程模型[M].武汉:武汉测绘科技大学出版社,2000:2-15

[2]周启鸣,刘学军.数字地形分析[J].测绘与空间地理信息,2007,30(3):146-146

[3]周启鸣,刘学军.数字地形分析[M].北京:科技出版社,2006:4-6

[4]杨勤科,Mcvicar TR,Vanniel TG,等.用ANUDEM建立水文地貌关系正确DEM的方法研究[J].测绘科学,2006,31(6):155-157

[5]郭伟玲,杨勤科,汪翠英,等.适用于地形复杂地区水土流失评价的高分辨率DEM建立方法[J].干旱地区农业研究,2008,26(3):246-252

[6]宋晓猛,张建云,占车生,等.基于DEM的数字流域特征提取研究进展[J].地理科学进展,2013,32(1):31-40

[7]杨昕,汤国安,刘学军,等.数字地形分析的理论,方法与应用[J].地理学报,2009,64(9):1058-1070

[8]杨封科.半干旱黄土丘陵区梯田集水增产效应研究[J].水土保持学报,2006,20(5):130-132

[9]李慧.梯田DEM构建方法研究[D].西安:西北大学,2014

[10]张宝国,吴倩,陈桂英,等.关于低山丘陵区水平梯田断面优化设计的探讨[J].水土保持应用技术,2007(2):31-33

[11]彭培,林爱文.基于AGREE算法的河流水系提取[J].水电能源科学,2015(4):27-29

[12]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

[13]Garbrecht J,Martz LW.The Assignment of Drainage Direction Over Flat Surfaces in Raster Digital Elevation Models[J].Journal of Hydrology,1997,193(1):204-213

[14]O'callaghan JF,Mark DM.The Extraction of Drainage Networks From Digital Elevation Data[J].Computer Vision, Graphics,and Image Processing,1984,28(3):323-344

[15]秦承志,李宝林,朱阿兴,等.水流分配策略随下坡坡度变化的多流向算法[J].水科学进展,2006,17(4):450-456

The Comparison of the Simulation Algorithms for Water Flow Path on the Terrace DEM

FENG Yuan1,ZHAO Mu-dan1*,FENG Yao2,WANG Yi-ren1,ZHANG Peng1
1.College of Urban and Environmental Science/Northwest University,Xi’an710100,China
2.Xi’an Aerospace Technology Co.,Ltd.,Xi’an710100,China

Large areas of terraces greatly changed the surface micro-morphology so as to affect the terraced surface water flow path and soil erosion.Based on 5 m resolution DEM data,the paper adopted the actual ridge-based method to build terrace DEM and used AGREE algorithm to pre-process terrace DEM,then compared the information of terrace DEM flow path under different catchment area threshold to explore proper fill sinks algorithm and flow direction algorithm applied to the terraces DEM flow path simulation.Results showed that using the M&G algorithm to fill sinks and MFD-md algorithm to determine flow could better simulate the terraced fields flow path.

Terrace DEM;AGREE algorithm;flow simulation;catchment area threshold

P208

:A

:1000-2324(2017)01-0041-05

2016-04-04

:2016-05-30

国家自然科学基金项目:梯田对坡度坡长因子的扰动特征研究(41271284)

冯园(1991-),女,硕士研究生,主要从事GIS空间分析.E-mail:15299193518@163.com

＊通讯作者:Author for correspondence.E-mail:zmudan@nwu.edu.cn