数字高程模型(digital elevation model，DEM)是用一组有序数值阵列形式表示地面高程的一种实体地面模型，描述的是地面高程信息，在测绘、水文、气象、地貌、地质、土壤、工程建设、通讯、军事等国民经济和国防建设以及人文和自然科学领域有着广泛应用。

为了更精细地反映地形特征及满足各行业对高精度的DEM的迫切要求，根据《广东省第一次全国地理国情普查实施方案》要求，广东省需要完成行政区域内1:10 000标准分幅2 m网格的数字高程模型精细化生产任务^[1-4]。本文以广东省国土资源技术中心承担的省级任务区DEM精细化为生产对象，重点阐述了基于数字线划地图(digital line graphic, DLG)数据的整个生产流程，各个步骤的具体方法和关键技术以及元数据的批量制作技巧。

1 DEM生产总体实施方案 1.1 数据源分析

由于生产的基础数据源的不同，DEM的生产需要分别采用相应的作业方法。

1) 采用1:10 000基础地理数据信息数据库的整合升级后并完成验收的DEM。对于坡度小于45°范围面积占图幅80%以上的图幅，直接利用已有5 m格网DEM对其重采样进行生产。

2) 广东省基础地理信息数据库中拥有覆盖任务区的1:10 000的DLG，并已经按省级要求完成1:10 000基础地理信息数据库整合升级工作，第一次全国地理国情普查DEM精细化生产采用整合升级后并完成验收的DLG数据。对于坡度小于45°范围面积占图幅80%及以下的图幅，主要利用1:10 000 DLG为资料源进行生产，然后对多种作业方法成果进行接边处理。

1.2 基于DLG的DEM精细化 1.2.1 作业流程

DEM精细化总体技术流程如图 1所示。

GeoTIN是武大吉奥信息技术有限公司推出的专门处理数字高程模型的软件，包含数据差错、三角网构建、特征点线编辑、数据裁切、数据接边等功能^[5]。在1:10 000 DLG数据中，利用等高线、高程点及水系等要素以及已有的地形特征数据，按照2 m×2 m的格网间距，基于充足的特征点采用不规则三角网(triangulated irregular network, TIN)内插方法生成DEM数据。GeoTIN生产DEM总体流程如图 2所示。

其具体步骤如下。

1) 精细化DEM生产前，要从数据源中提取等高线、高程点和水体等要素，并进行质量检查，主要检查高程点和等高线高程赋值正确性。

2) 选取双线河以及面积大于400 m²的静止水域，如湖泊、水库和池塘等，并根据相邻等高线或高程点估读其水涯线高程并赋值^[6]。

3) 对无高程点的山头、凹地或垭口等处，参考已有等高线高程，内插一个高程特征点，对狭长而缓慢的沟谷或山脊内插一条高程特征线，以提高DEM数据对真实地形的仿真度，优化DEM数据的相对高程精度。

4) 套合差超限处理。

5) 不同任务区、不同数据源、省界之间接边方案设计。

1.2.2 套合差检查

套合差检查是DEM精细化最为主要的检查工作，常见的套合差超限有局部高程值异常，与DLG套合超限而强行接边使高程异常^{[6, 7]}。

1) 通过反生等高线检查是否存在高程异常或其他不合理现象。DLG数据中部分河渠关系处理不当，会造成DEM编辑不合理。此时，应返回检查原始数据，造成矛盾的高程点(如拦水坝、主要堤、一般堤、河岸边的高程点)、落水的高程点应归属至人工高程点。河、湖、海岛上应增加人工特征点以保持地形地貌。双线河的高程应注意在有等高线穿过时应断开，分别赋上不同的高程值。要特别注意较大双线河的上下游高程。

2) 与原始等高线叠合显示，检查同名等高线的偏离情况，将偏离量超出1/2等高距的区域识别出来，对与该区域关联的DEM格网点的高程值进行合理编辑，使得所关联的DEM格网点与最近的原始等高线或高程点的高程差值的绝对值小于1/2等高距，以满足技术要求。

1.3 接边

DEM接边涉及不同任务区、不同数据源及省界之间接边。项目采用北京四维公司的PixelGrid-DEMEditDSM/DEM/DOM交互式编辑处理系统进行接边处理。1:10 000分幅DEM拼接示意图如图 3所示；PixelGrid中显示数据高程模型，并套入接边线如图 4所示；PixelGrid中DEM反生等高线如图 5所示；在融合区域内使用PixelGrid进行DEM的接边编辑如图 6所示。

为了确保DEM的生产精度，将周围8幅相邻DEM进行拼接，载入接边线，以等高线、水系面、高程点和检查点为基本数据构建TIN进而内插生成DEM并反生等高线接边时，沿着接边线设置缓冲重叠区域，在该区域内进行DEM的平滑接边编辑，PixelGrid中等高线平滑方法有均值、中值滤波、三角网内插、距离加权内插、曲面内插和置平等^[8]。这种接边方法能够只在局部小范围内对两边已经生产完成的DEM进行修改，最大程度地保证原有各任务区内部的DEM完整。接边后再重新裁切分幅，修改元数据的相关属性。

1.4 DEM裁切

按照国务院第一次全国地理国情普查领导小组办公室的要求，精细化DEM的最终成果是以DEM格网中心点坐标为定位点。为了便于DEM拼接以及带间接边处理，DEM的范围根据格网点计算公式所计算的结果再向四周增加3行3列高程格网点。格网点中心计算公式如下：

$\left\{ \begin{array}{l} {X_{\rm{起}}} = {\mathit{X}_{{\rm{max}}}} = \left[ {{\mathop{\rm int}} \left[ {\max ({\mathit{X}_{\rm{1}}}, {\mathit{X}_{\rm{2}}}, {\mathit{X}_{\rm{3}}}, {\mathit{X}_{\rm{4}}})/\Delta \mathit{d}} \right] + 1} \right]\Delta d - \Delta d/2\\ {Y_{\rm{起}}} = {Y_{\min }} = {\mathop{\rm int}} \left[ {\min ({\mathit{Y}_1}, {\mathit{Y}_{\rm{2}}}, {\mathit{Y}_{\rm{3}}}, {\mathit{Y}_{\rm{4}}})/\Delta \mathit{d}} \right] \times \Delta \mathit{d + }\Delta \mathit{d/}{\rm{2}}\\ {\mathit{X}_{\rm{止}}} = {\mathit{X}_{{\rm{min}}}} = {\mathop{\rm int}} \left[ {\min ({\mathit{X}_{\rm{1}}}, {\mathit{X}_{\rm{2}}}\mathit{, }{\mathit{X}_{\rm{3}}}, {\mathit{X}_{\rm{4}}})/\Delta \mathit{d}} \right] \times \Delta \mathit{d + }\Delta \mathit{d/}{\rm{2}}\\ {Y_{\rm{止}}} = {\mathit{Y}_{{\rm{max}}}}{\rm{ = }}\left[ {{\mathop{\rm int}} \left[ {\max ({\mathit{Y}_{\rm{1}}}, {\mathit{Y}_{\rm{2}}}, {\mathit{Y}_{\rm{3}}}\mathit{, }{\mathit{Y}_{\rm{4}}})/\Delta \mathit{d}} \right] + 1} \right] \times \Delta \mathit{d - }\Delta d/{\rm{2}} \end{array} \right. $

(1)

式中，(X₁, Y₁), …, (X₄, Y₄)为4个图廊坐标(+X指北，+Y指东)，单位为m；X_max, Y_min为DEM起始格网中心点坐标，单位为m；X_min，Y_max为DEM终止格网中心点坐标，单位为m；Δd为格网间距。

起止像元角点坐标计算公式(外扩3行3列后)为：

$\left\{ \begin{array}{l} {X_{\rm{起}}} = {\mathit{X}_{{\rm{max}}}} = {\mathop{\rm int}} \left[ {\max ({\mathit{X}_{\rm{1}}}, {\mathit{X}_{\rm{2}}}, {\mathit{X}_{\rm{3}}}, {\mathit{X}_{\rm{4}}})/2} \right]\times 2 + 7\\ {Y_{\rm{起}}} = {Y_{\min }} = {\mathop{\rm int}} \left[ {\min ({\mathit{Y}_1}, {\mathit{Y}_{\rm{2}}}, {\mathit{Y}_{\rm{3}}}, {\mathit{Y}_{\rm{4}}})/2} \right] \times 2 - 5\\ {\mathit{X}_{\rm{止}}} = {\mathit{X}_{{\rm{min}}}} = {\mathop{\rm int}} \left[ {\min ({\mathit{X}_{\rm{1}}}, {\mathit{X}_{\rm{2}}}\mathit{, }{\mathit{X}_{\rm{3}}}, {\mathit{X}_{\rm{4}}})/2} \right] \times 2 - 5\\ {Y_{\rm{止}}} = {\mathit{Y}_{{\rm{max}}}}{\rm{ = }}{\mathop{\rm int}} \left[ {\max ({\mathit{Y}_{\rm{1}}}, {\mathit{Y}_{\rm{2}}}, {\mathit{Y}_{\rm{3}}}\mathit{, }{\mathit{Y}_{\rm{4}}})/2} \right] \times 2 + 7 \end{array} \right. $

(2)

目前国内常用的处理软件通常只能以格网角点作为起点进行DEM数据的制作，因此采用栅格形式存储的DEM，起止像元角点坐标范围应在式(2)的基础上再外扩1 m(DEM格网间距2 m)。

1.4.1 单幅裁切

利用GeoTIN软件可以对单幅DEM进行标准分幅裁切，根据公式计算图廓点坐标(*.TIC文件)裁剪其外接矩形范围内的DEM文件作为结果存储。*.TIC文件格式如下(可以使用记事本编辑，西南角开始，顺序记录4个角点坐标)：

① 620 239 2 720 413

② 620 239 2 725 097

③ 626 623 2 725 097

④ 626 623 2 720 413

基于已有5 m DEM重采样2 m DEM时，格网角点坐标的定位可以利用ArcGIS栅格投影工具中的配准点选项输入西北角起点格网点角点坐标来确定。

1.4.2 批量裁切

针对多幅DEM融合后的区域，本文设计了一种批量裁切工具使用ArcGIS自带的Python脚本程序实现DEM 1:10 000分幅批量快速裁切。

1.5 DEM精细化结果

对生成的DEM进行分析，可以了解到广东省(省级任务区内)地势基本情况。广东省覆盖全地形，全省以山地为主，粤西湛江、茂名地区及粤东潮汕、梅州等地覆盖80%丘陵地；粤北韶光、河源地区高山地占30%~40%。广东省地形分类统计如表 1所示。

1) 高程中误差。采集高精度的高程检测点与DEM内差点高程比较，计算高程较差^{[9, 10]}，经过分析剔除粗差后，统计计算出DEM内插点高程中误差在设计书限差范围内, 如表 2所示。

2) 相邻图幅高程接边误差。按照设计要求，此次DEM精细化必须无缝接边，相邻图幅同名格网点高程值相同，接边精度100%，接边误差0%。

3) 套合差。将DEM成果反生成的等高线，与DLG原始等高线比较，两者套合差小于1/2等高距(仅对基于DLG的作业图幅)。

2 元数据批量制作

DEM元数据由数据基本情况、数据源情况、生产过程情况、数据质量情况和数据分发信息5部分组成，共97项。

在本项目中，开发了一款基于C#语言的批处理程序，程序主要有以下4个功能：

1) 根据IMG影像文件批量提取文件数据量、影像行列数，图幅号。结果格式如下：

文件名                 称行列数    数据量

F49G041044    2 320×3 231    8.2

…

2) 根据XLS文件批量提取指定单元格的数值，例如图幅号、图幅名称。结果格式如下：

F50G001001小正镇

…

3) 根据XLS文件由关键字段(图幅号)批量赋值各类数据至元数据指定单元格(图 7)。

4) 批量修改元数据指定单元格数值(图 8)。

各功能模块的实现调用自定义的*.TXT配置文件，配置文件遵守C#内部定义算法法则。在填写元数据之前按照标准图幅号一一对应数据项内容，输入程序即可实现自动化填写元数据批量制作，程序可以批量精确完成元数据制作中重复工作，提高工作效率，减少人工作业带来的错误风险，保证了DEM数据和元数据的一致性。

3 结束语

本文通过结合广东省地理国情的具体实践，详细论述了地理国情普查DEM精细化的工作内容和具体流程，重点阐述了基于DLG数据的整个生产流程各个步骤的具体内容，主要包括DLG数据预处理、DEM数据制作、构建、接边及元数据批量制作。在精细化DEM生产过程中，借助GeoTIN、PixelGrid辅助生产，提高DEM数据对真实地形的仿真度并保证了图幅边缘DEM的精度；元数据的制作过程中，研制了一款批处理工具实现自动化填写元数据批量制作。经过反复质量检查，广东省省级任务区内的4 247幅DEM精细化成果各项技术指标均符合国家标准。