基于EGM2008重力场高程异常值模型在矿山测量中的应用

王 兵

（贵州地矿测绘院有限公司，贵州 贵阳 550018）

随着全球卫星定位系统GNSS（Global Navigation Satellite System，全球导航卫星系统）不断发展和民用精度的不断提高，在工程测绘和矿山测绘等领域的重要性日益明显。目前，基于GNSS的RTK（Real Time Kinematic，实时动态）实时动态差分测量不仅在传统地面测绘还是在航空摄影测量中都发挥着重要作用。而RTK实时动态差分测量所得到高程系统为以参考椭球面为基准的大地高系统，而我国在工程领域或采矿等领域采用的是以我国1985国家高程为基准的似大地水准面正常高系统，故RTK实时动态差分测量所得到大地高不能直接使用，必须通过高程异常改正为正常高后才能使用，故建立全省高程异常值模型是必要的。建立全省高程异常值模型既利于在水准点稀少的区域获得高精度的正常高，又利于在测量点大地高向正常高批量转换电算化时有比较精确的模型数据基础，提高电算化精度和计算机运行效率，对贵州范围内高程异常值变化进行应用研究和掌握贵州范围内高程异常值变化趋势和值域分布，可为工程测绘或矿山测量等领域的正常高成果精度分析提供数据基础。

1 EGM2008全球重力场概述（Earth Gravitational Model EGM2008）

EGM2008重力场模型是较EGM84、EGM96之后由美国国家地理空间情报局地球重力场研发小组(US National Geospatial-Intelligence Agency)发布的较新地球超高阶重力场模型，该模型的阶数达到了2159次，并且加球谐系数的阶扩展到2190次，在精度和分辨率上都有巨大的提升。该机构发布EGM2008重力场模型的目的在于解算更加精确的WGS84坐标系统所对应的似大地水准面[1]，如下图1所示。

图1 EGM2008似大地水准面图

2 EGM2008高程异常计算

采用EGM2008重力场模型计算参考椭球面与EGM2008似大地水准面高程异常值，计算公式为：

式中：φ为地心纬度；λ为地心经度；ρ为地心向经；GM为地心引力常数；α为参考椭球长半轴；γP为正常重力值；Cmm和Smm分别为完全规格化位系数；为规格化Legndra函数；N为EGM2008重力场模型展开的最高阶数[2]。

本次采用AllTrans EGM2008 Calculator软件进行网格计算，经查询贵州省级行政区域在WGS84坐标系下范围为东经103.60°～109.59°，北纬24.62°～29.23°，以该区域为网格区域进行采样计算（如下图2所示），采样间隔为1′×1′，计算结果选择Bi-Quadratic值，EGM文件采用Und_min1x1_egm2008_isw=82_WGS84_TideFree_SE文件1′×1′网格[3]。

图2 AllTrans EGM2008 Calculator软件进行网格计算

经过该软件计算得到1′×1′网格的EGM2008高程异常值，输出的文件格式为文本文件数据交换格式，可以用文本文件或CSV格式EXCEL文件打开查看和编辑。

3 起算面偏差值测定和偏差改正

由于EGM2008模型起算面为全球高程基准所定义的起算面，与我国通过在青岛验潮站测得的平均海水面所使用的起算面有一定的差值，即起算面偏差Δζ[4]。起算面偏差可以通过平均分布在全省的多个水准点测得，采用在已知1985高程基准高程值的水准点上架设RTK设备测得该水准点的大地高，设所测得的大地高为H大地高，水准点处EGM2008高程异常值为ξEGM2008，已知水准点1985高程基准高程值为H85高，则有：

设AllTrans EGM2008 Calculator软件进行网格计算得到1′×1′网格的EGM2008高程异常值为m×n常数矩阵

采用全省多个水准点联测求得起算面偏差改正数m×n常数矩阵

则改正后大地高至1985国家高程基准高程异常值m×n矩阵为：

4 高程异常模型建立

根据公式3虽然可得到可用于RTK实时动态差分测量所得的大地高向1985国家高程基准转换的高程异常值，但是为m×n高程异常值矩阵，采样间隔为1′×1′，矩阵在实际使用中并不利于计算机查询，且在行列号相邻处存在截断误差[5]。

为解决这一问题，笔者采用不规则三角网的数字高程模型（简称TIN）根据矩阵中记录的高程异常值集合将区域划分为若干相连的三角面网络，形成一系列互不交叉、互不重叠的连接在一起的三角形来表示高程异常表面，形成一个三维空间的高程异常线性模型（如下图3所示）。用该TIN模型线性插值法生成加密的数字栅格模型，可有减少截断误差对查询结果的影响，并且利于计算机查询和二次开发[6]。我国使用的测绘成果坐标系统目前有大地坐标系统，有高斯克吕格投影三度分带和六度分带，由于采用高斯克吕格3度分带在贵州范围存在跨带，而六度分带变形过大，故模型采用大地坐标系统经纬度投影，可实现在全省内坐标统一，使用者在查询时须将平面坐标转换为大地坐标后根据坐标点位查询相应的高程异常值。若采用电算化方式编程计算时须将坐标转换为数字栅格模型所在行列号之后进行查询，计算机返回的行列号值即为查询坐标处的高程异常值[7]。

图3 贵州范围高程异常数字栅格模型

5 实例应用

本文以我单位承担的位于贵州省范围内一矿山倾斜摄影测量为例，受业主方委托，我单位采用倾斜摄影测量方式向业主提供5cm分辨率的矿山地表OSGB格式三维模型，要求坐标系统为2000国家大地坐标系，高程系统为1985国家高程基准。矿山为已开采的砂石矿山，占地面积约20hm2。

首先，我单位工作人员在矿山范围内均匀布设像控点5个，采用购买千寻CORS位置服务账号获得2000国家大地坐标系和大地高实时动态差分信息[4]。在像控点上正确架设RTK接收机，在固定解状态下平滑10次采集获得像控点2000国家大地坐标和大地高。采用计算机编程将平面坐标转换为大地坐标，在贵州范围高程异常数字栅格模型中采用计算机编程获得每个像控点处高程异常值，并将RTK获得的大地高改正为1985国家高程基准正常高[8]，具体信息见表1所示。

表1 像控点坐标及高程表 单位：m

为验证以上通过贵州范围高程异常数字栅格模型转换正常高的准确性，我单位在省测绘资料档案馆收集到距测区约2km处的已知点YZ01，采用RTK在已知点上校正高程后实测得到5个像控点正常高，并对其进行进精度比较分析，结果如表2所示。

表2 高程精度对比表 单位：m

经过实测与采用贵州范围高程异常数字栅格模型转换得到的高程对比发现，该测区高程最大误差约为2cm，最小误差约为0.5cm，满足业主方倾斜摄影测量的精度要求。经计算转换高程中误差为±1.2cm，完全满足四等水准高程精度要求，证明贵州范围高程异常数字栅格模型在该区域内准确可靠。

6 结语

本文通过具体的流程详述了应用EGM2008重力场模型计算参考椭球面与EGM2008似大地水准面的高程异常值，对该高程异常值进行起算面偏差改正，建立贵州范围高程异常数字栅格模型。对采用三角网插值法建立贵州范围高程异常数字栅格模型的算法公式进行图解和推导，并根据我单位承担的位于贵州省范围内一矿山倾斜摄影测量为实例验证贵州范围高程异常数字栅格模型精度和可靠性，本方法适用于各类采用1985国家高程基准的工程测绘和矿山测量等图根测量和碎部测量工作中，采用RTK实时动态差分测量进行高程异常改正精度可满足于四等水准测量精度要求。该模型可直接用于相关测绘地理信息系统软件开发，可用于监测定位传感器大地高转正常高自动化程序开发等，对类似测绘以及软件开发工作具有重要的参考价值。