CN115406401A - 减小矿区测量高程异常差值的方法 - Google Patents

Abstract

本申请公开了一种减小矿区测量高程异常差值的方法，包括以下步骤：步骤S1：从观测重力异常Δg中去除模型重力异常Δg_GM和地形重力效应δN_T得到残差重力异常，将残差重力异常代入Stokes公式，按下式计算残差大地水准面δN_r：步骤S2：在缺乏重力资料的情况下，根据高程异常的波谱结构特点分解大地水准面高N；步骤S3：在未知点上，由EGM2008地球重力场模型计算出未知点的模型高程异常，在所得结果上加上未知点剩余高程异常ξ_C，得到未知点最终的高程异常值，求得未知点的正常高。通过实践，基于“移去—恢复”法的多面函数拟合方法在大红山铜矿矿区的拟合效果较好，拟合后正常高精度可达到±2cm的要求。

Description

减小矿区测量高程异常差值的方法

技术领域

本申请涉及测绘工程技术领域，特别是一种减小矿区测量高程异常差值的方法。

背景技术

GPS技术作为现代空间大地测量的代表技术，具有精度高、成本低、快速、灵活等特点，GPS测量中需要结合高精度、高分辨率的(似)大地水准面模型，才能准确测定正高或正常高。

现有似大地水准面模型相对地面、参考椭球面相对位置关系如图1所示。目前，GPS定位技术可以精确，快速的测量得到点位三维坐标。然而，GPS所测得的高程是测站相对于WGS-84椭球面的大地高，而我国所采用的高程系统是相对于似大地水准面的正常高系统。

正常高是地面点沿着铅垂线方向到似大地水准面的距离h，而大地高是地面点沿法线方向到参考椭球面的距离H₀，H₀与h之差为高程异常。

若已知地面点的高程异常值ξ，则大地高H_大地与正常高H_正常的关系为：

H_大地＝H_正常+ξ

如何把大地高转化为正常高，成为在具体地理环境下，GPS高程应用的一个重要问题，现有解决高程转换问题的方法主要有：

一、综合利用GPS测量资料和高程异常资料确定点的高程；

二、综合利用GPS测量资料和水准测量资料确定似大地水准面的高程。

第一种方法需要足够的重力测量资料，但由于目前的重力测量资料严重不足、受制于实施重力测量难度和测量精度，更多情况下采用第二种方法：即结合GPS测量资料与水准测量资料，选择合适数学模型，通过拟合方法来获得所需要的正常高。

高程拟合方法的基本思想是：在GPS网中联测一些水准点，然后利用这些点上的正高和大地高，求出这些点位的高程异常值，再根据这些点上的高程异常值与坐标的关系，用最小二乘的方法计算出测区的似大地水准面，利用拟合出的似大地水准面，内插出其他GPS点的高程异常，从而求出各个未知点的正常高。

第二种方法在平坦或者地形不复杂的测区已经得到有效应用，并获得了令人满意的结果。

而在如大红山铜矿矿区这种地形复杂或者相邻区域地形变化较大的区域，以国家C级GPS控制点为起算基准点，在矿区中央、边缘合理选点和布网，均匀埋设11个控制点组成平面控制网。以国家一等水准点为起算基准点，沿矿区公路布设水准路线，利用电子水准仪进行三等水准测量，将高精度的水准高程引入矿区，到参与水准面模型计算的GPS控制点上。在大高差、几何水准路线难于到达的GPS点采用电磁波测距高程导线进行测量。

布设HZ01、HL(横里)、HP(红坡)、HZ12、3521、HZ04、G1、DHS5、35182、3508、HZ02-1等由11个组成的GPS C级控制网，最高点(HL)高程：1230m；最低点(35182)高程：695m；平均高程：970m。最短边：620m；最长边：4659m；已平均边长：2135m。并利用三等水准测量的方式，通过平差计算，解算出各控制点在CGCS2000坐标系下的高斯平面坐标(中央子午线东经102°)、1985年国家程基准框架下的正常高。以所选取的这11个点位构成矿区大地水准面模型建设的基础。其所提供的CGCS2000平面坐标、CGCS2000大地坐标与大地高、正常高等作为矿区进行大地水准精化的已知数据，并由这些数据通过模型构建，建立大红铜矿大地水准面基准。成果数据见表1所示。

表1 GPS控制网与高程控制测量成果表

为了对所建立的矿区大地水准面模型进行实地检验与检测，针对矿区高差大的特点，有代表性选择和布设了23个点位，作为检测点。主要利用GPS静态定位后处理技术(部份点位利用GPS RTK技术)测定了各点的GCS2000平面坐标、CGCS2000大地坐标与大地高。利用三等水准测量方法，测量并计算出这23个检测点的正常高，同时计算各点的高程异常值。检测点的成果见表2所示。

表2检测点成果表

由表2高程异常值项可以看出，各检测点的高程异常值均大于35m，说明按现有技术采用现有基准点获取的矿区大地水准面模型与复杂地理环境的实际高程差别较大，所得模型无法直接使用。如果不对现有方法进行改进，直接采用上述方法所得模型结果，会导致各区域基准平面偏离正常值，致使在GPS大地高到正常高的转换中精度严重损失，未能以相应的精度求解点位的高程。特别在高差较大的山区，大地水准面高程异常影响尤其明显，导致高程信息的精度较低。影响矿区生产的正常运行。

发明内容

本申请提供了一种减小矿区测量高程异常差值的方法，用于解决在高差起伏变化大、地形切割剧烈的矿区，直接采用国家一等水准点为起算基准点获取的GPS控制点，检测结果高程异常值高达35m以上，导致在地形复杂区域无法使用GPS技术进行准确定位的技术问题。

本申请提供了一种减小矿区测量高程异常差值的方法，包括以下步骤：

步骤S1：从观测重力异常Δg中去除模型重力异常Δg_GM和地形重力效应δN_T得到残差重力异常，将残差重力异常代入Stokes公式，按下式计算残差大地水准面δN_r：

其中，π为常数，G为观测重力异常，h_p为计算点的地形高，γ为平均正常重力，l为计算点至流动点的距离，ρ为地壳密度，；

进行大地水准面的恢复，分别计算N_GM和δN_T，其中N_GM由位系数计算得到；

步骤S₂：在缺乏重力资料的情况下，根据高程异常的波谱结构特点，按下式将大地水准面高N分解为：

ξ＝ξ_GM+ξ_ΔG+ξ_T 式(2)

其中，ξ为地球重力模型异常，ξ_GM为长波部分模型高程异常，ξ_ΔG为中波部分残差高程异常，ξ_T为短波部分地形高正；

在没有数字高程模型数据的情况下，把ξ_ΔG和ξ_T部分合在一起用数学模型逼近的方法表征；

将GPS点的高程异常按下式分为两部分求解，即：

ξ＝ξ_GM+ξ_C 式(3)

其中，ξ_C为剩余高程异常。

步骤S3：在未知点上，由EGM2008地球重力场模型计算出未知点的模型高程异常，在所得结果上加上未知点剩余高程异常ξ_C，得到未知点最终的高程异常值，求得未知点的正常高。

优选地，步骤S1包括以下步骤：

重力异常观测值也可以分解为相应的三个部分，即由位模型计算得到的模型重力异常，地形重力效应和残差重力异常，如果不考虑地形的影响，那么以上分解也可以分成两部分，即由位模型确定的部分和残余部分。这两种分解可以分别表示为：

N＝N_GM+δN_T+δN_r

Δg＝Δg_GM+δΔg_T+δΔg_r

步骤S11：移去过程，即将模型重力异常Δg_GM和地形重力效应δN_T从观测重力异常Δg中去除得到残差重力异常:

δΔg_r＝Δg-Δg_GM-δg_T

式中，Δg_GM利用位系数由下面的公式计算:

上式中，ψ，λ为计算点的纬度和经度，GM是地心引力常数，R为地球平均半径，

为完全规格化位系数，

为完全规格化缔合Legendre函数，N为位模型球谐展开的最大阶数；

δΔg_r为地形质量对计算点的引力，可以仅考虑相对于布格片的局部地形起伏的影响，即局部地形改正δΔg_TC，可以考虑相对于大地水准面的地形起伏的影响，其计算公式为：

步骤S12：将残差重力异常代入Stokes公式，计算残差大地水准面δN_r。

步骤S13：进行大地水准面的恢复，分别计算N_GM和δN_T，其中N_GM仍由位系数计算得到:

式中，ρ为计算点的地心距离，γ为平均正常重力。δN_T的计算可采用以下公式计算：

这里δN_T为局部地形影响。式中ρ为地壳密度，hp为计算点的地形高，h为流动点的地形高，l为计算点至流动点的距离。

优选地，步骤S2包括以下步骤：

在缺乏重力资料的情况下，根据高程异常的波谱结构特点，大地水准面高N可以分解为：

ξ＝ξ_GM+ξ_ΔG+ξ_T

ξ_GM——表示长波部分，称为模型高程异常，可以通过EGM2008重力场模型计算得到：

ξ_ΔG——表示中波部分，称为残差高程异常，可以通过求解剩余重力异常的边值得到：

ξ_T——表示短波部分，通过求解地形高正得到：

类似的，在没有数字高程模型数据情况下，把ξ_ΔG和ξ_T部分合在一起用数学模型逼近的方法表征。将GPS点的高程异常分为两部分求解，即：

ξ＝ξ_GM+ξ_C

式中，ξ_GM——重力场求得的高程异常；

ξ_C——剩余高程异常。

在无重力、无数字地形资料的情况下，通过若干已知的GPS/水准点，利用“移去-恢复”法求得其它未知点的高程异常，最终得出未知点的正常高。

优选地，步骤S3包括以下步骤：

结合各已知点的CGCS2000坐标值，将剩余高程异常值作为已知值，分别代入到加权平均模型、多项式曲线模型、平面模型、二次曲面模型、三次曲面模型、移动曲面模型和多面函数模型中，建立起以CGCS2000坐标值为参数，求解各未知点的剩余高程异常值函数模型。各未知点基于EGM2008模型，计算出各点的模型高程异常，最后得到各未知点的高程异常值。计算公式如下：

式中，ξ(x_i,y_i)——未知点的高程异常值；

ξ_C(x_i,y_i)——未知点基于函数模型的剩余高程异常插值；

ξ_GM(x_i,y_i)——未知点基于EGM2008的模型高程异常；

F_fit(x_i,y_i)——剩余高程异常值插值函数；

F_EGM2008(x_i,y_i)——EGM2008模型高程异常插值函数。

为了对各模型的效果进行测试，并利用上述方法计算出EMG2008模型高程异常，求出各检测点的剩余高程异常。

本申请能产生的有益效果包括：

1)本申请所提供的减小矿区测量高程异常差值的方法，结合大红山铜矿矿区高差大、地形变化激烈的实际情况，计算和比较了已知公共点和检核点数据，依据高程异常的可叠加性，采用“移去—恢复”法技术，首先移去有重力场模型计算的高程异常的中长波部分，求出的短波部分，然后选用数学模型拟合剩余的残余高程异常。通过实践，基于“移去—恢复”法的多面函数拟合方法在大红山铜矿矿区的拟合效果较好，拟合后正常高精度可达到±2cm的要求。

2)本申请所提供的减小矿区测量高程异常差值的方法，大红山铜矿矿区GPS测高基于建立的矿区大地水准面模型，以满足工程现场数据采集和内业批量数据处理的双重需求。实现GPS技术间接测定该地区任意站点的正常高，全面改善各类工程建设中高程测量的技术条件，实现GPS技术的三维定位功能。基于“移去—恢复”法的多面函数拟合方法在大红山铜矿矿区的拟合效果较好，拟合后正常高精度可达到±2cm的要求。

3)本申请所提供的减小矿区测量高程异常差值的方法，采用EGM2008重力大地水准面模型与高精度的几何大地水准面进行拟合，移去-恢复法是利用数学模型在高程转换前，首先移去用重力场模型计算的高程异常中的中、长波部分；然后对剩余的高程异常进行拟合和内插；再在内插点上利用重力场模型把移去的部分恢复，最终得到该点的高程异常。

4)本申请所提供的减小矿区测量高程异常差值的方法，采用EGM2008模型的"移去-恢复"法在大高差矿区大地水准面精化应用中的测高精度和适用性，实现了在EGM2008模型支持下的矿区地球重力场模型高程异常值计算，基于矿区高程异常的几何和物理特性，提出了建立厘米级矿区似大地水准面的模型与方法，代替常规的山区几何水准测量技术。

附图说明

图1为本申请提供的基于地球重力场模型的“移去-恢复”法示意图；

图2为本申请提供的顾及地球重力场模型的“移去-恢复”法示意图；

图3为本申请提供的EGM2008模型的全球误差分布图；

图4为基于EGM2008模型“移去恢复”法检核数据残差图；

具体实施方法

为使本发明实施方式的目的、技术方案和优点更加清楚，下面将结合本发明实施方式中的附图，对本发明实施方式中的技术方案进行清楚、完整地描述，显然，所描述的实施方式是本发明一部分实施方式，而不是全部的实施方式。通常在此处附图中描述和示出的本发明实施方式的组件可以以各种不同的配置来布置和设计。

因此，以下对在附图中提供的本发明的实施方式的详细描述并非旨在限制要求保护的本发明的范围，而是仅仅表示本发明的选定实施方式。基于本发明中的实施方式，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施方式，都属于本发明保护的范围。

本申请中未详述的且并不用于解决本申请技术问题的技术手段，均按本领域公知常识进行设置，且多种公知常识设置方式均可实现。

参见图1，下面结合附图对本发明做进一步的说明：

步骤S1：重力异常观测值分解为：由位模型计算得到的模型重力异常，地形重力效应和残差重力异常，如果不考虑地形的影响，那么以上分解也可以分成两部分，即由位模型计算得到的模型重力异常部分和残余部分，这两种分解方式可以分别表示为：

Δg＝Δg_GM+δΔg_T+δΔg_r

其中，Δg为观测重力异常，Δg_GM为模型重力异常，δΔg_T为局部地形改正值，δΔg_r为地形质量对计算点的引力；

N＝N_GM+δN_T+δN_r

其中，N为重力异常观测值，N_GM为模型重力异常，δN_T为地形重力效应，δN_r为残差重力异常；

步骤S11：移去过程，按下式将模型重力异常Δg_GM和地形重力效应δN_T从观测重力异常Δg中去除得到残差重力异常δΔg_r:

δΔg_r＝Δg-Δg_GM-δg_T

式中，模型重力异常Δg_GM采用位系数法按下式计算得到:

上式中，ψ，λ为计算点的纬度和经度，GM是地心引力常数，R为地球平均半径，

为完全规格化位系数，

为完全规格化缔合Legendre函数，N为位模型球谐展开的最大阶数；a为拟合系数；

δΔg_r为地形质量对计算点的引力，可以仅考虑相对于布格片的局部地形起伏的影响，即局部地形校正δΔg_TC，可以考虑相对于大地水准面的地形起伏的影响，计算公式为：

其中，δN_total为大地水准面地形起伏影响值，G_ρ为计算点的地心距离，γ为平均正常重力，h_p为计算点的地形高，

步骤S12：将残差重力异常代入Stokes公式，按下式计算残差大地水准面δN_r：

其中，R为地球平均半径，γ为平均正常重力，Δg_r为地形质量对计算点的引力，S(ψ)为计算点的纬度；

步骤S13：进行大地水准面的恢复，分别计算N_GM和δN_T，其中N_GM由位系数计算得到：

式中，ρ为计算点的地心距离，γ为平均正常重力；

δN_T按下式计算：

其中，δN_T为局部地形影响，ρ为地壳密度，h_p为计算点的地形高，h为流动点的地形高，l为计算点至流动点的距离，G为地心距离。

步骤S2：在缺乏重力资料的情况下，根据高程异常的波谱结构特点，将大地水准面高N分解为：

ξ＝ξ_GM+ξ_ΔG+ξ_T

其中，ξ_ΔG为中波部分，称为残差高程异常，ξ_T为短波部分地形高；ξ_GM为长波部分，称为模型高程异常，采用EGM2008重力场模型按下式计算ξ_GM：

其中，ψ，λ为计算点的纬度和经度，GM是地心引力常数，R为地球平均半径，

为完全规格化位系数，

为完全规格化缔合Legendre函数，N为位模型球谐展开的最大阶数；a为拟合系数；ρ为计算点的地心距离，γ为平均正常重力，

ξ_ΔG表示中波部分，称为残差高程异常，通过求解剩余重力异常的边值得到：

其中，R为地球平均半径，γ为平均正常重力，

为计算点的纬度、经度及地心距离值，

为高程异常；

为点的纬度；

ξ_T——表示短波部分，通过求解地形高正得到：

其中，ρ为计算点的地心距离，γ为平均正常重力，h_p为计算点的地形高，l为计算点至流动点的距离，h为流动点的地形高，G为点的地心距离；

类似的，在没有数字高程模型数据情况下，将ξ_ΔG和ξ_T部分合并后用数学模型逼近的方法表征，将GPS点的高程异常分为两部分求解，即：

ξ＝ξ_GM+ξ_C

式中，ξ_GM——重力场求得的高程异常；

ξ_C——剩余高程异常。

步骤S3：将高程异常移去EGM2008模型大地水准面高的剩余值认为是剩余高程异常值ξ_C。基于EGM2008模型的高程异常值，并以GPS大地高和三等水准测量得到的正常高求差，得到各已知点的真实高程异常值，剩余高程异常值由模型高程异常值与真实高程异常值求差得到，并认为是各已知点的真实剩余高程异常值。结合各已知点的CGCS2000坐标值，将剩余高程异常值作为已知值，分别代入到加权平均模型、多项式曲线模型、平面模型、二次曲面模型、三次曲面模型、移动曲面模型和多面函数模型中，建立起以CGCS2000坐标值为参数，求解各未知点的剩余高程异常值函数模型。各未知点基于EGM2008模型，计算出各点的模型高程异常，最后得到各未知点的高程异常值。计算公式如下：

式中，ξ(x_i,y_i)——未知点的高程异常值；

ξ_C(x_i,y_i)——未知点基于函数模型的剩余高程异常插值；

ξ_GM(x_i,y_i)——未知点基于EGM2008的模型高程异常；

F_fit(x_i,y_i)——剩余高程异常值插值函数；

F_EGM2008(x_i,y_i)——EGM2008模型高程异常插值函数。

为了对各模型的效果进行测试，并利用上述方法计算出EMG2008模型的高程异常，求出各检测点的剩余高程异常。

实施例

在大红山铜矿矿区，通过静态GPS观测，三等几何水准测量，得到了11个已知的GPS/水准点，按照步骤S1及步骤S2过程进行点位高程计算。

其中，δN_T为局部地形影响，ρ为地壳密度，h_p为计算点的地形高，h为流动点的地形高，l为计算点至流动点的距离，G为点的地心距离。

ξ＝ξ_GM+ξ_C

式中，ξ_GM——重力场求得的高程异常；

ξ_C——剩余高程异常。

则可求得这11个点的高程异常ξ_k＝h_k-H_k(k＝1,2…11)。在这些点上用地球重力场模型计算出高程异常ξ_GM，最后得到剩余高程异常ξ_C。计算结果如表3所示。

表3 11个已知点的高程异常计算表

以11个已知点的残差高程异常ξ_C作为已知数据，用常规拟合方法计算出拟合模型的拟合系数，再内插出未知点的剩余高程异常ξ_C。拟合函数如前所述的函数拟合模型，分别为二次曲面拟合，三次曲面拟合，多面函数(正双曲面)拟合、多面函数(倒双曲面)拟合，多面函数(锥面函数)拟合，多面函数(三次曲面)拟合。

按照步骤S3过程，

式中，ξ(x_i,y_i)——未知点的高程异常值；

ξ_C(x_i,y_i)——未知点基于函数模型的剩余高程异常插值；

ξ_GM(x_i,y_i)——未知点基于EGM2008的模型高程异常；

F_fit(x_i,y_i)——剩余高程异常值插值函数；

F_EGM2008(x_i,y_i)——EGM2008模型高程异常插值函数。

在未知点上，由EGM2008地球重力场模型计算出未知点的模型高程异常，再加上未知点剩余高程异常ξ_C，得到未知点最终的高程异常值，从而求得未知点的正常高。

为了对各模型的效果进行测试，在矿区布设了23个检测点，并利用上述方法计算出EMG2008模型高程异常，求出各检测点的剩余高程异常。计算结果如表4所示。

表4 23个检测点的高程异常计算表

选取均匀分布于大红山铜矿矿区的23个点进行实际测试。通过直接测量，得到检核点的大地高及正常高，求差得到检核点的高程异常ξ_i值。再运用上述基于EGM2008模型“移去-恢复”法计算得到检核点的高程异常拟合值ξ′_i，对其求差得到各个函数的残差值。基于EGM2008模型“移去恢复”法外部检核数据残差统计如下表5所示。EGM2008重力场模型结合GPS/水准测量数据，选择用移去—恢复法确定大地水准面，通常是在高程异常变化复杂，且GPS/水准点的数量足够多的情况下。采用重力场模型与GPS/水准相结合的方法确定局部似大地水准面能真实地反映大地水准面的精度，取得了良好的效果。基于EGM2008模型“移去恢复”法外部检核数据残差图如下图4所示。

由图4可知，通过本申请提供方法处理后，所得地理位置结果，逐渐趋近……，说明该方法能有效校正各测试点的误差，尤其适用于大红山矿区，高程差起伏过大的区域。结合大红山铜矿矿区高差大、地形变化激烈的实际情况，计算和比较了已知公共点和检核点数据，拟合后正常高精度可达到±2cm的要求。

表5基于EGM2008模型“移去恢复”法检核数据残差统计表

由表3～5可见，基于EGM2008模型“移去—恢复”法中，二次曲面拟合，正双曲面为核函数的多面函数拟合、锥面函数为核函数的多面函数拟合，三次曲面为核函数的多面函数拟合方法在大红山铜矿矿区的拟合效果较好，拟合后正常高精度可达到±2cm的要求。

尽管参照前述实施例对本发明进行了详细的说明，对于本领域的技术人员来说，其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换,凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (4)

1.一种减小矿区测量高程异常差值的方法，其特征在于，包括以下步骤：

步骤S1：从观测重力异常Δg中去除模型重力异常Δg_GM和地形重力效应δN_T得到残差重力异常，将残差重力异常代入Stokes公式，按下式计算残差大地水准面δN_r：

其中，π为常数，G为观测重力异常，h_p为计算点的地形高，γ为平均正常重力，l为计算点至流动点的距离，ρ为地壳密度，；

进行大地水准面的恢复，分别计算N_GM和δN_T，其中N_GM由位系数计算得到；

步骤S₂：在缺乏重力资料的情况下，根据高程异常的波谱结构特点，按下式将大地水准面高N分解为：

ξ＝ξ_GM+ξ_ΔG+ξ_T 式(2)

其中，ξ为地球重力模型异常，ξ_GM为长波部分模型高程异常，ξ_ΔG为中波部分残差高程异常，ξ_T为短波部分地形高正；

在没有数字高程模型数据的情况下，把ξ_ΔG和ξ_T部分合在一起用数学模型逼近的方法表征；

将GPS点的高程异常按下式分为两部分求解，即：

ξ＝ξ_GM+ξ_C 式(3)

其中，ξ_C为剩余高程异常；

步骤S3：在未知点上，由EGM2008地球重力场模型计算出未知点的模型高程异常，在所得结果上加上未知点剩余高程异常ξ_C，得到未知点最终的高程异常值，求得未知点的正常高。

2.根据权利要求1所述的减小矿区测量高程异常差值的方法，其特征在于，步骤S1包括以下步骤：

重力异常观测值也可以分解为相应的三个部分，即由位模型计算得到的模型重力异常，地形重力效应和残差重力异常，如果不考虑地形的影响，那么以上分解也可以分成两部分，即由位模型确定的部分和残余部分；

这两种分解可以分别表示为：

N＝N_GM+δN_T+δN_r

Δg＝Δg_GM+δΔg_T+δΔg_r

步骤S11：移去过程，即将模型重力异常Δg_GM和地形重力效应δN_T从观测重力异常Δg中去除得到残差重力异常:

δΔg_r＝Δg-Δg_GM-δg_T

式中，Δg_GM利用位系数由下面的公式计算:

上式中，ψ，λ为计算点的纬度和经度，GM是地心引力常数，R为地球平均半径，

为完全规格化位系数，

为完全规格化缔合Legendre函数，N为位模型球谐展开的最大阶数；

δΔg_r为地形质量对计算点的引力，可以仅考虑相对于布格片的局部地形起伏的影响，即局部地形改正δΔg_TC，可以考虑相对于大地水准面的地形起伏的影响，其计算公式为：

步骤S12：将残差重力异常代入Stokes公式，计算残差大地水准面δN_r；

步骤S13：进行大地水准面的恢复，分别计算N_GM和δN_T，其中N_GM仍由位系数计算得到:

式中，ρ为计算点的地心距离，γ为平均正常重力；δN_T的计算可采用以下公式计算：

这里δN_T为局部地形影响；式中ρ为地壳密度，hp为计算点的地形高，h为流动点的地形高，l为计算点至流动点的距离。

3.根据权利要求1所述的减小矿区测量高程异常差值的方法，其特征在于，步骤S2包括以下步骤：

在缺乏重力资料的情况下，根据高程异常的波谱结构特点，大地水准面高N可以分解为：

ξ＝ξ_GM+ξ_ΔG+ξ_T

ξ_GM——表示长波部分，称为模型高程异常，可以通过EGM2008重力场模型计算得到：

ξ_ΔG——表示中波部分，称为残差高程异常，可以通过求解剩余重力异常的边值得到：

ξ_T——表示短波部分，通过求解地形高正得到：

类似的，在没有数字高程模型数据情况下，把ξ_ΔG和ξ_T部分合在一起用数学模型逼近的方法表征；

将GPS点的高程异常分为两部分求解，即：

ξ＝ξ_GM+ξ_C

式中，ξ_GM——重力场求得的高程异常；

ξ_C——剩余高程异常；

在无重力、无数字地形资料的情况下，通过若干已知的GPS水准点，利用“移去-恢复”法求得其它未知点的高程异常，最终得出未知点的正常高。

4.根据权利要求1所述的减小矿区测量高程异常差值的方法，其特征在于，步骤S3包括以下步骤：

结合各已知点的CGCS2000坐标值，将剩余高程异常值作为已知值，分别代入到加权平均模型、多项式曲线模型、平面模型、二次曲面模型、三次曲面模型、移动曲面模型和多面函数模型中，建立起以CGCS2000坐标值为参数，求解各未知点的剩余高程异常值函数模型；

各未知点基于EGM2008模型，计算出各点的模型高程异常，最后得到各未知点的高程异常值；计算公式如下：

式中，ξ(x_i,y_i)——未知点的高程异常值；

ξ_C(x_i,y_i)——未知点基于函数模型的剩余高程异常插值；

ξ_GM(x_i,y_i)——未知点基于EGM2008的模型高程异常；

F_fit(x_i,y_i)——剩余高程异常值插值函数；

F_EGM2008(x_i,y_i)——EGM2008模型高程异常插值函数；

为了对各模型的效果进行测试，并利用上述方法计算出EMG2008模型高程异常，求出各检测点的剩余高程异常。

Images