基于ANUDEM与山脊抬升的DEM简化

高 翔，徐 柱

（西南交通大学 遥感信息工程系，四川 成都611756）

为得到粗尺度DEM，需要对精细尺度DEM进行简化［1－2］。现阶段，常用的地形简化方法有两类：基于滤波的简化和基于地形特征的简化［3］。滤波方法常使用移动窗口对DEM进行过滤运算，如表面简化［4］和稀疏采样［5］等。该类方法虽简单高效，但易造成地形的过度平滑从而丢失重要的地形特征。基于特征的简化方法首先提取地形结构特征，之后利用得到的特征信息重建粗尺度DEM，因此，基于特征的地形简化方法能够更好地保留地形特征［6－7］。

本文利用ANUDEM（Australia National University DEM）算法基于特征点以及水文要素对DEM实施地形简化［8］。该方法采用迭代有限微分内插技术和地形强化算法，自动去除伪洼地［9－11］，同时该算法可利用等高线、水文线、洼地、边界线和海岸线等要素对插值过程进行约束从而得到与现实地表高度吻合的地形表面。然而，实验表明，基于特征点与水文的ANUDEM插值方法能够充分保留地表的水文要素，却无法顾及山体的结构特征，如山脊线要素。作为地形骨架线，山脊线在很多应用中有着比水文线更为重要的作用，如山体单元划分［12－13］与冰川编目［14］。因此，在ANUDEM算法基础上，提出一种山脊抬升的DEM简化方法。该方法将特征点、水文线和山脊线同时加入到地形简化过程中，最终生成基于特征点与骨架线约束的地形简化表面。实验部分，对不同分辨率下的简化DEM进行山脊抬升并对比分析抬升前后的山脊线保留情况。实验表明，新方法在保证地形参数准确性的情况下更好地保留了初始地形的骨架信息

1 方 法

为充分保留地形骨架信息，利用特征点与骨架线进行地形简化。该方法由4部分组成：①特征点、骨架线的提取；②基于特征点与水文线的ANUDEM内插；③简化地形的山脊抬升；④山脊抬升前后山脊线保留度统计分析。本文利用Arc GIS 10．2进行DEM、等高线和地形透视图的展示；利用开源软件GAT（Geospatial Analysis Tools）完成特征点与水文线的提取；利用Arc Tool box中的Topo to Raster工具完成ANUDEM插值过程，山脊线提取与抬升、地形参数和骨架线保留度的计算在IDLE（Python GUI）平台中编程实现。

1．1 提取特征点、骨架线

1．1．1 提取特征点

常用的特征点提取方法有：基于线简化的特征点提取算法［15］、VIP算法［16］及Z最大容差法。其中，Z最大容差法是最为常用的特征点提取算法。该方法首先提取一系列点并构建TIN；接着向TIN中加入高程误差（简化后DEM与初始DEM高程差）最大的点并重构TIN；循环执行第二步直到任意点的高程误差均小于设定阈值为止［17］。其中，Z容差阈值与简化DEM分辨率的对应关系如表1所示［18］。

表1 容差值与DEM分辨率的对应关系 m

1．1．2 提取水文线与山脊线

本文采用水文模拟法提取DEM中的水文线，其一般步骤为①对初始DEM实施填洼操作②计算DEM流向；③基于流向计算汇水累积量；④基于汇水累积量计算水文线（本文选取500为累积量判别阈值，因为当阈值为500时可获得较为丰富的水系）；⑤基于流向数据，将水文转换为线状要素。之后，利用文献［14］设计的算法提取DEM中的山脊线，并利用文献［19］对地形特征线进行优化处理。

1．2 山脊抬升

完成ANUDEM插值过程后，对简化表面实施山脊抬升。若仅增加山脊区域的高程，山脊与其邻近区域将存在较大高差。为减小山脊区域高程的剧烈变化，本文设计了基于最短欧式距离的山脊抬升方案，该方法通过计算栅格至最近山脊的最短距离来确定当前栅格的抬升高度。基于抬升方法的DEM山脊修复公式为

其中：R为当前DEM分辨率，D为某栅格至其最近山脊的最小欧式距离，K为抬升率，H为抬升高程值。

图1为不同K值和不同H值下，地表抬升高程值与最小欧式距离的函数曲线，可以发现：①K值相同时，H值越大，相同距离下的抬升值越大，且随着距离的缩小，这一规律更加明显；②H值相同时，K值越大，相同距离下的抬升值越小，其变化趋势也越明显。因此，为减小初始DEM与处理后的DEM间的高程差异，H值不应过大，K值则不应过小。实验部分将通过研究地形基本参数的变化规律确定H与K的取值。

图1 抬升高程值与距离的函数关系

2 统计分析

本文选用平均坡度与表面粗糙度对山脊线抬升前后的地形复杂度进行刻画。平均坡度¯S与表面粗糙度K计算如式（2）、式（3）所示。

其中：S为坡度，A为投影面积，A′为表面面积，i为第i个栅格单元，n为DEM栅格总数。为描述山脊线抬升前后地表高程的变化规律，计算其高程均方根误差，如式（4）所示。

图2所示，式中的L为初始DEM提取出的山脊线原长（黑色线条）。为计算地形简化后山脊线保留长度，首先为初始线段创建半径r＝R／2的缓冲区，其中R为DEM分辨率；接着统计简化后的骨架线落入缓冲区内的线段长度L′i；最终基于式（5）获得保留度［18］。

图2 骨架线保留度

图2中实线为原始DEM骨架线，边长L；虚线为简化后DEM骨架线；粗线为简化后骨架线落入原始DEM骨架线缓冲区内的线段，长度为L′i。

3 实验与分析

本文选取GDEM 30 m分辨率DEM进行实验分析，DEM大小及基本参数如表2所示。

表2 DEM地形基本参数

由表2可知：实验所选DEM高差、平均坡度以及表面粗糙度都较大，这表明该区域具有一定的地表复杂度，因此适合作为实验区域进行山脊抬升前后的对比分析。最终，所选DEM提取的骨架线如图3所示。

图3 实验DEM骨架线提取结果（黑色细线为等高线，灰色粗线为水文线，黑色粗线为山脊线）

3．1 山脊抬升对地形的影响

确定式（1）中山脊抬升的基本参数K和H，研究不同取值下山脊抬升前后简化地表的均方根误差、平均坡度以及表面粗糙度的变化规律。实验数据如表3和表4所示。其中，表3为H值相同时，不同K值下的地形参数取值；表4为K值相同时，不同H值下的地形参数取值。

表3 不同K值下的地形参数取值

表4 同H值下的地形参数取值

由表中数据可知 与均方根相比 平均坡度与表面粗糙度变化极小，因此山脊抬升不会显著影响简化地形的基本形态。图4进一步描述了均方根误差随山脊抬升参数的变化规律。可以发现：①H值不变时，随着K值的增加，高程均方根误差呈幂函数递减趋势，当K值小于1时，均方根误差迅速减小，K值大于1时，均方根误差的变化趋于平缓，如图4（a）所示，因此实验部分的K值取1；②K值不变时，随着H值的增加，均方根误差呈线性增长，如图4（b）所示，因此为在一定的均方根误差内提升山脊区域DEM高程，实验部分的H值取10。

完成抬升参数的确定后，利用ANUDEM算法与山脊抬升法得到简化DEM。其中，50 m、250 m、500 m分辨率下的简化DEM如图5所示。

图4 均方根误差与抬升参数的函数关系

图5 基于ANUDEM与山脊抬升法得到的不同简化程度下的DEMs

3．2 山脊保留程度研究

3．2．1 等高线对比

图6为基于ANUDEM方法与山脊抬升法得到的50 m分辨率DEM的局部等高线，等高距为100 m。可以发现：①非山体区域，山脊抬升前后的等高线近似相同；②山体区域，山脊抬升前后的等高线具有一定的差异，与原始DEM相比，山脊抬升后的DEM更好地保持了山体特征。图6（a）丢失了高程为4 600 m的山峰，而图6（b）保留了这一地形要素，同时，图6（b）中高程值为4 700 m的等高线与原始等高线的吻合度也高于图6（a）。

3．2．2 山脊线保留度对比

为定量描述抬升前后山脊线保留程度，利用式（5）进行保留度计算，结果如表5所示。

通过表5数据，可以发现：

1）随着地形简化程度的缩小，特征点百分比逐渐増大，山脊线保留度逐渐增加。如当分辨率为500 m时，抬升后的山脊线保留度仅为43%左右，而当分辨率减小至50 m时，抬升后的山脊线保留度达到了92%。

图6 山脊抬升前后DEM等高线差异

2）当地形简化程度较高时，特征点百分比较小，简化DEM的数据量较小，对DEM实施抬升处理并不能够有效增加山脊线保留度，随着简化程度的降低，山脊抬升能够一定程度地增加山脊线保留程度。如当分辨率等于500 m时，特征点总数百分比仅为0．38%，两者的山脊线保留度均为43%左右，无明显差异。当分辨率减小至50 m时，与抬升前相比，抬升前后的山脊保留度分别为78%与92%，保留度增加了14%左右。

4 结 论

本文研究了一种基于ANUDEM与山脊抬升的DEM简化方法。该方法首先利用最大容差法提取地形特征点，接着利用水文模拟法提取地表水文线以及山脊线，之后通过ANUDEM插值法获得简化DEM，最终对简化DEM实施山脊抬升从而获得顾及骨架线的简化地表。

为确定山脊抬升方法的基本参数，计算抬升后DEM的均方根误差、平均坡度以及表面粗糙度。数据表明，H值不变时，均方根误差随K值的减小呈幂函数递减，当K值为2时，均方根误差的变化趋于稳定；K值不变时，均方根误差随H值的增加呈线性增长。与均方根误差相比，山脊抬升对平均坡度与表面粗糙度的影响较小，因此山脊抬升不会造成地形复杂度的较大变化。

为进一步研究不同简化地形下山脊抬升对山体特征的保留程度，计算山脊抬升前后的山脊线保留度。实验表明，在插值生成的简化DEM中实施山脊抬升能够一定程度的增加山体结构的保留程度，且随着地形简化程度的降低，山脊抬升法的作用也逐渐明显。

［1］ ZHANG X，DRAKE N A，WAINWRIHGT J，et al．Co mparison of Slope Esti mates fro m Low Resolution DEMs：Scaling Issues and a Fractal Method f or Their Solution［J］．Earth Surface Processes and Landf or ms，1999，24（9）：763－779．

［2］ 周小军，王光霞，薛志伟，等．基于DEM化简的等高线综合研究［J］．测绘工程，2014，23（2）：10－14．

［3］ LI Zhilin．Multi－scale Digital Terrain Modelling and Analysis［C］／／Zhou Q．，Lees B．，Tang G．Advances in Digital Terrain Analysis．Berlin：Springer－Verlag，2008：59－83．

［4］ LI Zhilin，ZHU Q．Digital Elevation model［M］．Wuhan：Wuhan University of Surveying and Mapping Press，2003．

［5］ LI Zhilin．Variation of t he Accuracy of Digital Terrain Models with Sampling Interval［J］．Photogrammetric Record，1992，14（79）：113－128．

［6］ GESCH D B．The Effects of DEM Generalization Methods on Derived Hydr ologic feat ures［C］／／Lowell K，Jaton A．Spatial Accuracy Assess ment：Land Information Uncertainty in Natural Resources．Chelsea：Ann Ar bor Press，1999：255－262．

［7］ AI T，LI J．A DEM Generalization by Minor Valley Branch Detection and Grid Filling［J］．ISPRS Journal of Photogrammetry and Remote Sensing，2010，65（2）：198－207．

［8］ HUTCHINSON M F．ANUDEM Version 4．6．2．http：／／cres．anu．edu．au／soft ware／anudem．ht ml as of July 1999．

［9］ 杨勤科，Mcvicar，T．R．，李领涛，等．ANUDEM－专业化数字高程模型插值算法及其特点［J］．干旱地区农业研究，2006，24（3）：36－41．

［10］张彩霞，杨勤科，段建军．一种高质量的数字高程模型（DEM）建立方法－ANUDEM法［J］．中国农业通报，2005，21（12）：411－415．

［11］丑述仁，姚志宏，曹佳云，等．基于Hutchinson的DEM建立及质量评价［J］．地理空间信息，2012，10（3）：127－129．

［12］肖飞，张百平，凌峰．基于DEM的地貌实体单元自动提取方法［J］．地理研究，2008，27（2）：459－467．

［13］徐静，王春，张耀民，等．规则格网DEM中平直面状特征地 形 识 别 与 提 取［J］．测 绘 科 学，2014，39（8）：163－166．

［14］郭万钦，刘时银，余蓬春，等．利用流域边界和坡向差自动提取山脊线［J］．测绘科学，2011，36（6）：210－213．

［15］DOUGLAS D H，PEUCKER T K．Algorith ms for the Reduction of the Nu mber of Points Required to Represent a Digitized Line or Its Caricature［J］．Canadian Cartographer，1973，10（2）：112－122．

［16］CHEN Z，GUEVARA J．Systematic Selection of Very Important Points（VIP）fro m Digital Terrain Model for Constructing Triangular Irregular Net works［R］．In：Proc．8thInternational Sy mposium on Computer－Assisted Cartography，Auto Carto8，Balti more，MD，1987，29 March－3 April．50－56．

［17］CHANG K T．Introduction to Geographic Infor mation Systems［M］．Singapore：Mc Graw－Hill，2008．

［18］CHEN Y，WILSON J P，ZHU Q，et al．Comparison of Drainage－constrained Methods f or DEM Generalization［J］．Co mputers and Geosciences，2012，48，41－49．

［19］冯长强，江南，张巍，等．地形特征线提取后续问题优化处理［J］．测绘工程，2014，23（3）：28－31．