张 亮，杨勤科，2，兰 敏

（1.中国科学院 水利部 水土保持研究所，陕西 杨凌712100；2.西北大学 城市与环境学院，西安710069；3.中煤西安设计工程有限责任公司，西安710054）

数字高程模型（Digital Elevation Models，DEMs）是一组表示景观内基于某一基准面的高程点在空间上的分布状况的数字阵列［1］。DEM包括规则格网、不规则三角网和等高线三种模型，规则格网DEM数据结构简单的优点使其在数字地形分析中的应用最为广泛［2－3］。基于地形测量成果数据（等高线、高程点和河流等）通过插值建立规则格网DEM，是目前最为成熟和经济实用的DEM建立方法［4－6］。在基于等高线插值建立格网DEM过程中，只有设置合理的分辨率才能保证地形信息不损失或少损失，又不引入过多的数据冗余［7－9］。曲率是表示曲线弯曲程度的参数，它可以直接反映等高线的复杂程度，是等高线信息量的重要指标［10］。因而计算等高线曲率，可以通过衡量等高线建立DEM过程中信息的转换率，为准确设置DEM分辨率提供参考。为满足实际应用的需求，前人研究者们提出了多种地形曲率指标及其算法［11］。然而这些算法均基于格网DEM，通过局部窗口的曲面拟合函数，来模拟地形表面的方式实现；由于等高线节点数据分布不规则等原因，如何基于等高线直接计算地形曲率则缺乏必要的研究。另外，格网DEM本身是等高线的派生数据，基于其计算的地形曲率则是二次派生数据，精度比直接基于等高线计算的结果低，所以研究用等高线计算地形曲率有其必要性。本文将讨论矢量格式等高线曲率的计算方法，以期为DEM分辨率的确定提供科学依据，对数字地形分析中地形曲率的研究也具有理论意义。

1 研究方法

1.1 等高线曲率算法

等高线曲率，是指通过地表上某一点的等高面（水平面）与地表的交线（即等高线）在该点的曲率，曲率是衡量曲线弯曲程度的指标；等高线曲率表达了地表物质运动的汇合和发散模式。

为了量化表示等高线的弯曲程度，本文用等高线上相邻三个节点拟合曲率圆，来计算等高线曲率。不在一条直线上的三点可以唯一确定一个圆，中间一点处的曲率用半径的倒数表示。如三点共线（即等高线的局部为直线），则直接赋曲率值为零。

如图1所示，设r为三点（Pi－1，Pi，Pi＋1）所组成三角形的外接圆半径，那么1／r即可认为是点Pi处的曲率。根据海伦公式、正弦定理、三角形面积公式可以推导出，圆拟合法计算等高线曲率的公式为：

式中：Cc——节点 Pi处的等高线曲率；a，b，c——线段的Pi－1Pi，PiPi＋1，Pi＋1Pi－1长度，等高线曲率单位为弧度每米（rad／m）。由于该方法以单条等高线上的节点为基础来计算曲率，因而暂不能表示地形的凹凸变化。本文主要讨论等高线曲率，为简便起见，下文中未指明的曲率一律特指等高线曲率。

图1 等高线曲率计算示意图

1.2 等高线曲率的计算

本研究基于Microsoft.NET平台用ArcEngine组件库开发定制程序，用于计算等高线曲率。图2为等高线曲率计算的数据流程图，表1为数据点的存储结构，主要步骤包括：（1）对等高线数据（CNT）做拓扑和连接等方面的错误修复，拓扑错误包括等高线叠合、交叉、自相交等，修复后的等高线图层标记为CNT＿conn；（2）使用ArcGIS的Feature Vertices To Points工具将CNT＿conn转换成点图层，并用Add XY coordinates工具提取等高线节点坐标，结果标记为CNT＿points；（3）编写程序，以CNT＿points中的“ORIG＿ID”字段（记录了所在等高线的索引号）为依据读取同一条等高线上相邻三节点的坐标数据，用公式（1）计算等高线曲率并存储于CNT＿points属性表的curvature字段中。

图2 等高线曲率计算数据流程

表1 等高线节点数据结构

1.3 试验数据

本文选取黄土丘陵区、黄土塬区和东北漫岗丘陵区三个典型试验样区。黄土丘陵样区（图3a，样区1）位于陕北安塞县县南沟小流域，平均坡度29.9°，梁地形比较简单平缓，沟坡则较复杂，海拔1 070.3～1 338.5m，面积2.41km2。黄土塬样区（图3b，样区2）位于陕西长武县，平均坡度19.2°，塬面地形比较简单，沟坡地形比较复杂，海拔范围977.9～1 197.0m，面积2.44km2。东北漫岗丘陵区（图3c，样区3）位于黑龙江省拜泉县罗家屯沟小流域，平均坡度3.9°，整体呈波状缓坡丘陵，少部分较陡峭，海拔范围254.8～305.8m，面积2.40km2。研究选用1∶1万等高线作为试验数据，其中黄土丘陵区等高距5m（CNT1），黄土塬区等高距10m（CNT2），东北漫岗丘陵区等高距2.5m（CNT3）。

图3 试验样区等高线

1.4 地形类型划分及相关地形指标提取

地形类型划分：用手工的方式分别勾绘黄土高原两个试验样区的沟沿线和塬边线，将黄土丘陵区地形简单区分为梁和沟壑；黄土塬区分为现代沟坡、古代沟坡及塬面三种；根据坡度特征将漫岗丘陵样区划分为较缓坡和较陡坡丘陵（图4）。这种划分并非严格的地貌类型划分，仅用于分析比较不同地形中等高线曲率的统计分析。

图4 地形分区图

相关地形指标提取：用样区1等高线数据（5m等高距）建立2.5m分辨率格网DEM（DEM1），样区2等高线（10m等高距）建立5m分辨率DEM（DEM2），样区3等高线（2.5m 等高距）建立1m分辨率DEM（DEM3）。然后分别提取其坡度（SLOPE1，SLOPE2，SLOPE3）、坡 长 （SLPLEN1，SLPLEN2，SLPLEN3）、剖 面 曲 率 （PROCURV1，PROCURV2，PROCURV3）和 平 面 曲 率 （PLANCURV1，PLANCURV2，PLANCURV3），用于比较分析等高线曲率与相关地形指标之间的空间格局分布关系。

2 结果与分析

2.1 空间格局分析

等高线曲率的计算结果（图5）可见：（1）等高线比较弯曲的地方曲率值大，等高线比较平直的地方曲率值小；（2）在山谷线和山脊线附近等高线曲率计算值比较大，而一般斜坡上曲率值小。可见本文建议的圆拟合曲率算法计算的等高线曲率，结果是正确的。但仔细观察可见，在比较平值的等高线上，曲率计算值略有波动，这可能是由于等高线上的线段节点的局地误差所致，对此有待进一步研究。

图5 试验样区等高线曲率图（局部）

按照前述地形类型划分（图4）对等高线曲率进行地形分区统计（表2），结果表明：（1）黄土丘陵区沟壑地形等高线曲率要普遍高于梁地形；（2）黄土塬区现代沟坡等高线曲率比古代沟坡和塬面均大；（3）漫岗丘陵区较陡坡地形也是比较缓坡地形等高线曲率大。而这种大小变化，均由等高线弯曲程度不同所致，这进一步说明本文提出的算法，计算结果符合实际情况。

2.2 统计分布分析

黄土丘陵区、黄土塬区和漫岗丘陵区等高线曲率相关统计数据（图6）表明：（1）等高线曲率值集中在0.2以下（92%～98%），而且曲率集中于小值的趋势依黄土丘陵区、黄土塬区、漫岗丘陵区的次序增强；（2）等高线曲率的平均值和范围方面，黄土丘陵区均大于黄土塬区和漫岗丘陵区。这是由于黄土丘陵区等高线弯曲程度普遍较大，而黄土塬区在塬面和古代沟坡上，特别是漫岗丘陵区等高线相对比较平直。

2.3 与常见地形指标的比较

提取常见地形指标，包括坡度、坡长、剖面曲率和平面曲率如图7—9所示（等高线曲率图空白为背景色），其与等高线曲率的对比可以看出：（1）等高线曲率与平面曲率的空间格局最为相似，特征线条沿地形特征线方向，都是在山脊线和山谷线处曲率大；但数值上并不完全相同。（2）在丘陵区的沟壑、塬区的现代沟和漫岗丘陵区的较陡坡区，等高线曲率、剖面曲率和坡度等数值均较大，但从空间上并不具有十分精确的对应关系。（3）剖面曲率反映的是地形在垂直方向上二阶导数，与坡度的变率相关，而坡度主要由等高线之间的距离决定，因此剖面曲率图和坡度图在空间格局和纹理特征上都最为相似。（4）等高线曲率与坡长的空间分布关系似乎不太大。

表2 试验样区等高线曲率分区统计 rad／m

图6 等高线曲率频率曲线

图7 黄土丘陵区各地形因子

图8 黄土塬区各地形因子

图9 漫岗丘陵区各地形因子

3 结论与讨论

3.1 结 论

（1）本文设计的圆拟合算法，从曲率的数学含义入手，用曲线拟合的方法，可实现对等高线曲率的计算。结果表明在等高线比较弯曲的地方，等高线曲率比较大；从地形来看，沟壑地带（沟沿线到坡脚线之间）等地形比较复杂的地方等高线曲率比较大。等高线曲率计算值的空间分布与地形特征具有良好的对应关系，结果合理可信。

（2）不同复杂程度的样区，曲率计算结果的空间和统计分布有所不同，地形越复杂等高线曲率值越大；黄土丘陵区的沟壑处等高线曲率最大，黄土塬区的现代沟坡次之，古代沟坡区较小，塬面和缓坡丘陵区则最小。等高线曲率可直接、客观地反映等高线的复杂程度，在一定程度上也可以反映地形的复杂程度。

3.2 讨 论

本文初步实现了对等高线曲率的计算，但是还有一些问题有待深入研究，包括：（1）已经实现的等高线曲率计算方法是基于单条等高线来计算曲率，没有使用等高线拓扑关系信息，因而尚不能区分曲率的正负，也就不能区分地形的的凸凹变化。（2）等高线曲率与栅格数据计算的坡度、坡长、平面曲率和剖面曲率等地形指标之间，从概念、空间格局、统计分布等方面，具体的相同处和不同处，有什么量化关系。（3）比较平滑的等高线上，本文计算的曲率值存在的局部变异，有待通过对算法容错能力的提高，或者通过对等高线的修复处理来解决。

［1］ Moore I D，Grayson R B，Ladson A R.Digital terrain modelling：a review of hydrological，geomorphological，and biological applications［J］.Hydrological Processes，1991，5（1）：3－30.

［2］ 周启鸣，刘学军.数字地形分析［M］.北京：科学出版社，2006：57－63.

［3］ 程峥，姚志宏，张婧，等.南方红壤农作区数字流域建立和应用［J］.水土保持研究，2012，19（3）：12－15.

［4］ Hutchinson M F.A new procedure for gridding elevation and stream line data with automatic removal of spurious pits［J］.Journal of Hydrology，1989，106（3）：211－232.

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

［6］ Wilson J P.Digital terrain modeling［J］.Geomorphology，2012，137（1）：107－121.

［7］ 沈晶玉，史明昌.DEM网格尺寸对地形因子影响研究：以北京市延庆县八达岭小流域为例［J］.水土保持研究，2006，13（5）：66－69.

［8］ 杨勤科，张彩霞，李领涛，等.基于信息含量分析法确定DEM分辨率的方法研究［J］.长江科学院院报，2006，23（5）：21－23.

［9］ 王峰，王春梅.地形因子与DEM分辨率关系的初步研究：以蒙阴县为例［J］.水土保持研究，2009，16（4）：225－229.

［10］ 陈其明.在自动制图中对曲线信息量的数学描述与压缩［J］.测绘学报，1981，10（3）：239－254.

［11］ Evans I S.General geomorphometry，derivatives of altitude，and descriptive statistics［M］∥Chorley R J.Spatial Analysis in Geomorphology.London：Methuen＆Co.，1972.