高云军，刘召芹，万文辉，孔　振

(1.山东科技大学，山东 青岛 266590；2.中国科学院遥感与数字地球研究所 遥感科学国家重点实验室， 北京 100101)



DEM网格间距对重力近区地形改正精度影响及验证

高云军1,2，刘召芹2，万文辉2，孔振1,2

(1.山东科技大学，山东 青岛 266590；2.中国科学院遥感与数字地球研究所 遥感科学国家重点实验室， 北京 100101)

摘要：利用分形随机算法建立平地、丘陵和山地3种精细地形仿真场景，将DEM逐级重采样为不同网格间距，分析不同DEM网格间距对3种地形的重力近区地形改正误差影响。发现随着DEM网格间距的增加，近区地形改正误差随之增大。对于平地，使用1∶10 000的DEM，网格间距为5 m仍能够满足规范要求；对于丘陵地，使用1∶5 000的DEM，网格间距为2.5 m能够满足规范要求；对于山地，使用1∶1 000的DEM,网格间距1 m能够满足规范要求。通过消费级无人机获取丘陵地精细地形，验证地形仿真的结论，同时说明消费级无人机能够应用于重力近区地形改正。

关键词：地形仿真；DEM；网格间距；重力近区地形改正；消费级无人机

重力测点周围的近区(0～20 m)地形起伏变化，对测点重力观测值的影响最大[1]，重力近区地形改正对地形产品的精度要求较高而相应的大比例尺地形数据的缺乏，随着大比例尺重力调查的开展，重力近区地形改正成为制约重力调查精度的重要因素。测点的实际地形改正值大小无法测得，只能通过加密特征点来逼近真实的地形改正值[2]。随着仪器设备及测绘手段的发展，使得获取近区测点附近精细三维地形成为可能。D.Schiavone等[3]使用三维激光扫描仪获取测点100 m范围内网格间距为1 m的DEM(Digital Elevation Model)；邸凯昌等[4]提出一种基于全景立体视觉和摄影测量技术的快速近区重力地形改正方法，能够自动生成测点30 m范围内网格间距为1 m的DEM，这些方法均是通过测点周围的DEM进行重力近区地形改正。DEM网格间距直接关系到地形表达的精度，DEM网格间距越小，越能够逼近真实的地形，地形精度越高，重力近区地形改正的精度也越高，但过于精细的DEM会增加DEM生产的工作量，随着DEM网格间距的增大，地形改正误差也随之增大。因此，在满足重力测量规范的情况下，如何选择合适的DEM网格间距，对于快速实现重力近区地形改正具有实际的指导意义。

本文利用分形随机算法仿真模拟生成平地、丘陵和山地3种地形的DEM，网格间距为0.1 m，以此DEM为真值，通过立方卷积插值法将3种DEM逐级重采样为不同网格大小的DEM，利用“方域”形算法分别计算重力近区地改值，来评价不同网格间距对重力近区地改值的影响。为了验证仿真，利用消费级无人机在试验近区采集影像，对无人机上搭载的超宽角鱼眼相机进行了高精度标定，利用标定得到的相机参数对原始影像数据畸变校正，制作实验区高精度DEM，结合测区内重力测点位，计算不同DEM网格距时重力近区地形改正误差，与仿真结果进行验证分析。

1网格间距影响仿真分析

1.1地形仿真

采用分形随机算法仿真三维地形，通过递归算法用简单的规则来模拟复杂地形[5]。本次研究采用最广泛的中点位移法来模拟三维地形，其基本过程[6]如图1所示。

图 1　正方形细分方法

1)设置二维地形数组4个角点(点A,B,C,D)的高程值hA,hB,hC,hD，在正方形中心点(O)处生成一个随机变量rO，计算中心点的高程值

(1)

2)在正方形每条边的中心处(点E,F,G,H) 各生成一个随机量rE，rF，rG，rH，计算各点的高程值如下：

(2)

3)根据重新生成的4个正方形网格，重复1)和2)直到地形数据细分到所需的地形精细程度为止。

随机变量主要是为了增加地形起伏高度的随机性，使仿真地形更加真实。本文采用高斯分布模型，为了防止地形出现明显的周期性，每次产生的随机变量乘以一个地形变化的比例系数和分形参数，随机变量的计算公式[7]为

(3)

根据文献[8]有关地形类别划分规则，利用分形随机算法建立0.1 m网格间距的平地、丘陵和山地3类地形，生成仿真地形面积都为0.5 km×0.5 km，每类地形均匀选取64个测点，如图2所示。不同地形测点20 m范围内地形点与测点高差统计如表1所示。

表1　 地形点与测点高差统计结果　m

1.2地形改正计算和误差分析

由于地形复杂，应用基本公式计算地形改正值时，积分上下限无法给出，一般均采用近似积分的办法，将以测点为中心的四周地形分割成许多小块，计算出每一个小块地形质量对测点的重力值影响值，然后累加求和便得到该点的地形改正值[9]。考虑到仿真DEM数据为n×n的网格数据，采用方形域算法进行地形改正值的计算[10]，计算公式为

(4)

式中：σ为岩石密度，取2.67 g/cm3；G为万有引力常数取6.67×10-8cm3/(g·s2)；d为网格间距；r为节点到测点的距离；ΔH为节点与测点的高程差；ci,j为梯形积分系数。

图2　三类仿真地形及测点位分布

利用DEM进行重力近区地形改正，精度好坏无法按传统方法进行直接评价，参考DEM在中区地形改正中精度的评价方式，以网格间距为0.1 m的DEM计算地形改正值作为真值(见表2)；评价重采样后的不同DEM网格计算的地形改正值误差。

利用立方卷积插值法将0.1 m的DEM以0.1 m为增量逐级重采样为不同网格间距，计算测点在不同网格间距时的重力近区地形改正值，3种地形不同网格间距及地形改正中误差见图3。

表2　利用网格距0.1 m的DEM计算的地改值

图3　不同DEM网格间距对应的地改误差

从图3可见地形改正误差随着DEM网格间距的增加而变大，而且地形起伏越大，DEM网格对地形改正误差影响越明显，平地最小，丘陵次之，山地影响最大。

1) 平坦区由于地形起伏较小，DEM网格间距对地形改正的影响较小，当DEM网格间距为5 m时，地形改正中误差为0.005×10-5m/s2，仍能够满足文献[11]中限差0.020×10-5m/s2的精度要求；根据表3大中比例尺数字高程模型要求，对于平坦区域的重力近区地形改正，DEM比例尺为1∶1万，网格间距为5 m能够满足重力近区地形改正要求。

2) 对于丘陵地区，地形起伏变大，当DEM网格间距大于3 m后，地形改正中误差已经超出规范中的限差；对于丘陵地形，结合表3，建议在实际工作中，DEM比例尺取1∶5 000，网格间距为2.5 m。

表3　数字高程模型的格网尺寸

3) 由于山地地形起伏最大，随着DEM网格间距的增大，地形的平滑现象越明显，地形坡度信息损失越多，地形改正误差变化趋势最大，当DEM网格间距大于1.1 m后，地形改正中误差已经超过规范限差。山地地形DEM比例尺建议不小于1∶1 000，网格间距至少为1 m，方能满足大比例尺重力调查规范中的精度要求。

2无人机重力近区地形改正对比验证

无人机遥感观测技术作为一种远程观测技术，克服了地形、气候的限制,以及成本低、作业简单、效率高等独特优势，广泛应用于地形图测绘和地质灾害监测等领域[12]。消费级旋翼无人机作为一种新的数据获取平台，因其操作简单和系统成本低，极大地降低了航拍难度和成本。本次实验通过消费级无人机获取研究区域精细地形，验证不同的网格间距对重力近区地形改正的影响，同时说明消费级无人机在重力近区地形改正适用性。

2.1实验区和数据

本次实验区位于唐山市迁西县，地形以丘陵为主。实验选用轻小型消费级无人机DJI Phantom 2 Vision+(DJI2+)，搭载视场角为140°的鱼眼相机。经过实地探勘地形，设定无人机飞行高度约为120 m，东西方向规划为4条航线，共获取91张影像，覆盖面积约为0.55 km×0.54 km。无人机航拍影像位置如图4所示。

图4　无人机航拍影像位置

2.2数据处理

为了消除相机本身的系统误差和影像畸变，对相机进行了高精度标定。采用等距投影来表示鱼眼相机的成像模型[13]，利用三维相机标定场对鱼眼相机内参数进行高精度的标定，标定后反投误差为0.21个像素；利用标定得到的相机参数对原始影像数进行畸变矫正，从而得到无畸变影像数据。

利用SFM(Structure from Motion)技术[14]，解算得到实验区三维密集点云数据；结合外业实测的GPS控制点，对三维数据进行绝对定向。由于三维密集点云中不可避免地混有地物点，为了提取真实数字地面模型，采用基于坡度的滤波算法，滤去地面上的植被和房屋，根据大比例尺DEM生产规范，生成网格间距为0.1 m的DEM。

2.3DEM精度验证

(5)

误差分布如表4所示。

表4　误差统计结果　m

从表4可以得出：无人机制作的DEM内插高程与实测高程点高程中误差为0.088 m，能够完全满足1∶500数字高程模型的精度要求[16]。

2.4地形仿真对比验证

为了验证仿真结论，在测区内选取了49个测点，在测点近区范围内与测点最大高差为22.67 m，最小高差为1.03 m，平均高差为9.15 m，测区地形及点位分布如图5所示。

图5　测区地形及点位分布

为了分析不同的DEM网格距对重力近区地形改正的影响，利用立方卷积插值法将网格间距为0.1 m的DEM以0.1 m为增量逐级重采样为不同网格间距，利用方形域算法计算测点在不同网格间距时的近区地形改正值；将各测点在网格间距为0.1 mDEM时计算的地改值作为真值，经计算地形改正值最大值为0.386×10-5m/s2，最小值为0.006×10-5m/s2，平均地形改正值为0.116×10-5m/s2，计算在不同DEM网格间距时的地形改正中误差差，仿真和实测地形改正中误差与DEM网格间距关系见图6。

图6　不同DEM网格间距对应的地形改正误差

随着DEM网格间距的增大，地形改正误差逐渐增大，实测地形和仿真地形的地形改正误差的变化趋势基本一致。当DEM网格间距大于3 m后，实测地形改正误差已经超出规范中的限差[11]，实测结果和地形仿真的结论相符，对于丘陵地形DEM比例尺使用1∶5 000，网格间距为2.5 m完全能够满足规范精度要求，同时说明消费级无人机航测能够应用于重力近区地形改正中。

3结束语

随着DEM网格间距的增大，重力近区地形改正误差逐渐增大，地形起伏越大，地形改正误差对DEM网格间距的敏感度越高。对于地势起伏比较小的地形，DEM的比例尺为1：10 000，网格间距为5 m时，引起的重力近区地形改正误差在规范要求限差内；对于丘陵地为主的地形，DEM的比例尺宜采用1∶5 000，网格间距为2.5 m，对于以山地为主的地形，DEM比例尺宜采用1∶1 000，网格间距为1 m。在条件允许的情况下，应优先使用更为密集网格的高程数据计算重力近区地形改正值，从而提高重力近区地形改正精度。

消费级无人机从垂直方向成像能够从完整度上实现较优的精细地形重建，由于起飞和降落容易，能够应用于重力地形改正地形条件复杂的测区，在保证重力近区地形改正精度的前提下，极大地提高了地形改正的效率，为我国重力近区地形改正提供了一条有效途径。

参考文献：

[1]刘宽厚，耿涛，杨怀英，等.基于便携激光测距仪的重力测量近区地形改正系统[J]. 物探与化探，2012，36(3)：404-408.

[2]张国利，赵更新，匡海阳，等. 近几年地面重力调查工作方法技术的一些进展[J]. 地球物理学进展，2015，30(1)：386-390.

[3]SCHIAVONE D,CAPOLONGO D,LODDO M.High resolution dems for near-station terrain correction in gravimetry[R]．2007 EGM．In-ternational Workshop，2007．

[4]邸凯昌，吴凯，刘召芹，等.全景立体视觉的快速近区重力地形改正方法(英文)[J]. 遥感学报， 2013(4)：759-767.

[5]徐琼燕，李俊峰，吴印华.基于分形理论的三维地形生成算法研究与实现[J].计算机与现代化，2011(5)：70-73.

[6]李巍，孙义威，万文辉,等. 月面形貌仿真设计及在嫦娥三号任务中的应用[J]. 地理研究，2014，33(6)：1040-1048.

[7]俞静，石晶晶.基于分形理论的三维模拟地形仿真[J].成都大学学报(自然科学版)，2011，30(3)：249-252.

[8]杨亚斌，韩革命，梁萌.重力近区地形改正精度探讨[J].物探化探计算技术，2011，33(1)：92-96.

[9]王成锡.中国区域重力信息系统集成[D].北京：中国地质大学，2012.

[10] 刘文锦，奚家鉴，张兴雅.区域重力测量的地形改正方法[J].物探与化探，1983，7(2):77-83.

[11] 中华人民共和国地质矿产部.大比例尺重力勘查规范:DZ/T0171—1997[S].北京：中国标准出版社，1997.

[12] 高姣姣. 高精度无人机遥感地质灾害调查应用研究[D].北京：北京交通大学，2010.

[13] SCHNEIDER D,SCHWALBE E,MAAS H G. Validation of geometric models for fisheye lenses. ISPRS Journal of Photogrammetry and Remote Sensing,2009,64：259-266.

[14] 魏占玉，Arrowsmith Ramon，何宏林,等. 基于SfM方法的高密度点云数据生成及精度分析[J].地震地质，2015，37(2):636-648.

[15] 张祖勋，张剑清.数字摄影测量学[M].武汉：武汉大学出版社，2012：65-67.

[16] 国家测绘局.基础地理信息数字成果 1∶500 1∶1 000 1∶2 000 数字高程模型：CHT 9008.2-2010[S].北京：测绘出版社，2010.

[责任编辑：刘文霞]

DOI:10.19349/j.cnki.issn1006-7949.2016.08.004

收稿日期：2015-10-11

基金项目：国家973课题(2013CB733202);国家自然科学基金资助项目(41471388)

作者简介：高云军(1990-)，男，硕士研究生.

中图分类号：P631

文献标识码：A

文章编号：1006-7949(2016)08-0014-05

Influence and verification of DEM grid spacing on the accuracy ofgravity near-region terrain correction

GAO Yunjun1,2，LIU Zhaoqin2，WAN Wenhui2，KONG Zhen1,2

(1.Shandong University of Science and Technology,Qingdao 266590，China；2.State Key Laboratory of Remote Sensing Science,Institute of Remote Sensing and Digital Earth,Chinese Academy of Science,Bejing 100101,China)

Abstract：This paper first establishes a fine terrain simulation in flat region, hilly region and mountain region using random fractal algorithm. Then, the DEM is resampled in different grid spacing. Lastly, the influence of different grid spacing on three terrain of gravity near-region correction is analyzed in this way. The error of gravitynear-region correction increases as DEM grid spacing increases. In flat region, the calculated DEM data with 1∶1 000 scales and grid spacing of 5 m can meet the requirements of specifications. In hilly region, the calculated DEM data with 1∶5 000 scales and grid spacing of 2.5 m can meet the requirements of specifications. And in mountain region, the calculated DEM data with 1∶1 000 scales and grid spacing of 1m can also meet the requirements of specifications. In this research, the results of terrain simulation are verified with a fine terrain in hilly region acquiring with consumer UAV (Unmanned Aerial Vehicle). And it shows that consumer UAV can be applied to the terrain correction of gravity.

Key words：terrain simulation; DEM; grid spacing; gravity near-region terrain correction; consumer UAV