开源GRASS软件等高线生成DEM的算法改进

谢中凯，陈振杰，李飞雪

（1.江苏省地理信息技术重点实验室，江苏 南京 210023）

传统的由人工解译地形获取的等高线数据是建立DEM的重要数据源[1]。矢量等高线生成规则格网DEM的方法主要有3种[2,3]：①等高线构建TIN法，矢量等高线生成不规则三角网TIN，然后由TIN内插生成DEM；②等高线离散化法，矢量等高线离散化为数据点，然后采用一般插值方法（如IDW、Kriging插值法等）内插生成DEM[4]；③等高线内插法，矢量等高线按照预定方向直接内插生成DEM，或首先进行栅格化，然后由栅格化后的等高线内插生成DEM。

不同的GIS软件在利用等高线生成DEM时采用了不同的方法，ArcGIS采用等高线构建TIN法生成DEM，该方法是利用等高线生成DEM的常用方法，但如果在TIN生成时没有人工加入地面特征点，生成的TIN会有平三角形，这样插值生成的DEM会产生“平山顶”现象和阶梯状地形。地理资源分析支持系统GRASS是开源GIS软件[5]，采用等高线内插法生成DEM，分为矢量栅格化和栅格插值2个步骤，其中栅格插值采用线性内插方法，通过邻近的等高线进行比例内插。本文分析了GRASS软件中等高线生成DEM算法中存在的问题，并对算法进行了改进。

1 GRASS中DEM生成算法存在的问题

1）在矢量栅格化过程中，如果多条等高线通过一个栅格单元（图1中为2条），该栅格单元的取值取决于最后进行栅格化的等高线的高程值（图中红色单元为待插值栅格单元，黑色单元为等高线栅格化后的栅格单元（称为等高线栅格单元））。

图1 2条等高线过同一栅格单元

2）在栅格插值过程中，没有采用最近邻的2条等高线进行内插，而是通过四方向Flood Fill算法寻找近邻等高线。四方向Flood Fill算法的特点是可以避免在搜索时穿越等高线，具体思路为：对红色栅格单元的邻域按照“下、右、上、左”的顺序进行搜索判断（如图2），对搜索到的非等高线栅格单元重复此操作，直到找到2个不同高程值的等高线栅格单元。当算法终止时（如图3），搜索结果为具有不同高程值的蓝色等高线栅格单元，但显然紫色等高线栅格单元较上侧的蓝色栅格单元距离红色待插值栅格单元更近，因此在插值时应当采用紫色等高线栅格单元。

图2 四方向Flood Fill算法图

图3 四方向Flood Fill搜索结果图

3）内插生成的DEM存在“平山顶”现象。对于图4中的普通情况，A处所在的待插值栅格单元可通过Flood Fill算法找到2条不同高程值的等高线，采用式（1）进行线性内插。但是对于图5中的山顶情况，由于Flood Fill算法搜索时不能穿越等高线，搜索结束时只能找到1条等高线，采用式（2）确定A点的高程值，那么被等高线Z1围成的山顶将具有相同的高程值，生成的DEM存在“平山顶”现象。

图4 普通情况

图5 山顶情况

2 算法改进及实现

针对GRASS算法存在的3个问题，分别提出如下解决方案：

1）在矢量栅格化过程中，统计经过每个栅格单元的等高线的条数，并记录其高程值，然后取高程平均值作为该栅格单元的取值。

2）在栅格插值过程中，采用链表记录搜索到的等高线，每条等高线的信息作为一个链表结点，按距离升序排列，内插时取最近的2条等高线。在搜索过程中，当找到第2条等高线后，开始用前2条等高线的最大值作为搜索范围的限制，减少栅格单元搜索的个数，加快搜索速度。

3）在栅格插值过程中，对于等高线围成的山顶，采取二轮搜索策略，第二轮搜索时允许跨越一次等高线，使得山顶采用最近邻的2条等高线进行趋势内插[6]。在进行趋势内插时，待插值栅格单元采用式（3）计算，并且应设定一个较小的正值，采用式（4）计算，避免山顶的高度过高。

3 算法实现及评价

等高线生成DEM的改进算法利用Visual Studio 2010平台和GDAL库，以C++语言实现，其中矢量栅格化采用Bresenham算法，栅格插值采用改进后的GRASS线性内插算法。算法测试选取面积约4 km2的山地区，原始DEM分辨率为5 m×5 m。首先由原始DEM获取10 m等高线，然后用不同的算法由等高线重新插值生成DEM。算法评价通过3组DEM之间的比较，3组DEM分别为：ArcGIS软件生成的DEM（未加入人工特征点）、GRASS软件生成的DEM和本文改进算法生成的DEM。

在研究区内随机生成500个随机采样点，进行DEM的精度评估[7]。图6中，横轴为原始DEM的采样值，纵轴为不同算法生成的DEM的采样值，可以看出3种算法得到的拟合直线斜率均接近1，拟合度R2值达到0.99以上。但是本文提出的算法拟合直线的斜率要高于其他2种，而截距要小于其他2种，说明改进后的算法取得了较好的效果，提高了插值精度。

对插值效果的评价，采用等高线回放法[8]。如图7所示，红色线为原始等高线，绿色线为不同算法生成的DEM回放得到的等高线。图7a、b中回放得到的等高线与实际等高线相比，均在山顶处出现了失真现象，图7a中还在等高线曲率较大的地方出现了失真现象，而图7c中回放得到的等高线与实际等高线相比，误差较小。图7a中出现失真是因为在山顶等特殊地形处生成了平三角形，导致生成的DEM具有“平山顶”现象和阶梯状地形。图7b中失真是因为山顶处插值时直接采用最近等高线赋值，导致生成的DEM具有“平山顶”现象。

图6 不同算法生成的DEM精度评价

图7 等高线回放比较

在特殊地形处的插值效果方面，选取山顶处DEM进行比较。图8中上图为不同算法生成的DEM，下图为相应的山体阴影模型。图8a、b中可以看到明显的“平山顶”现象，而图8c中这种现象得到了消除，与原始DEM相接近。

图8 山顶处DEM生成比较

4 结 语

对开源GRASS软件中等高线生成DEM的算法进行改进后能提高DEM生成的精度，解决了山体陡峭地区等高线生成DEM存在的“平山顶”问题。由于现实中地形是多变的，因此等高线的分布也是各异的，对于多条等高线经过同一栅格单元，下一步的研究可以考虑统计栅格单元内不同高程值的等高线的长度，然后以长度为权重计算单元取值。

[1]宋敦江, 岳天祥, 杜正平.由等高线建立DEM的YUEHASM方法研究[J].地球信息科学学报, 2009,11(3): 325-332

[2]张凯选, 潘梦清, 方辉, 等.利用等高线生成 DEM 方法的研究[J].测绘工程, 2007,16(3): 15-18

[3]杨晓云, 唐咸远, 梁鑫.基于等高线生成DEM的内插算法及其精度分析[J].测绘工程, 2006,15(2):37-39

[4]徐潇, 谭衢霖, 王浩宇, 等.复杂地貌地形图等高线内插DEM算法的精度分析[J].遥感信息, 2013,28(6):111-115

[5]吴超.地理资源分析支持系统GRASS的研究与应用[D].武汉：华中科技大学，2005

[6]吕建峰, 刘定生，焦伟利, 等.DEM生成算法并行化研究[J].中国图像图形学报,2002, 7(5):506-512

[7]寇程, 柯长青.地形平坦地区DEM生成算法的比较研究[J].测绘与空间地理信息,2013, 36(7): 33-37

[8]唐新明, 林宗坚, 吴岚.基于等高线和高程点建立DEM的精度评价方法探讨[J].遥感信息,1999,14(3):7-10