不同插值方法下重力场模型的似大地水准面构建精度分析

Accuracy Analysis of Quasi Geoid Construction for Gravity Field Models under Different Interpolation Methods

摘要: 因高程系统和GPS系统的坐标基准不同，导致GPS高程无法直接应用到水准成果。常规的方法是通过数学模型进行高程异常拟合，但纯数学模型难以满足三等和四等水准的精度要求。使用EGM2008重力场模型的三种插值方法计算中长波，使用数据模型拟合短波，能够有效提高高程异常拟合精度。经过计算，三种插值方法的双二次插值法精度最高，4个检核点的绝对误差达到3.3 cm以内。

Abstract: Due to the different coordinate benchmarks between the elevation system and the GPS system, GPS elevation cannot be directly applied to leveling results. The conventional method is to fit elevation anomalies through mathematical models, but pure mathematical models are difficult to meet the accuracy requirements of the third and the fourth order levels. Using three interpolation methods of the EGM2008 gravity field model to calculate medium to long waves and fitting short waves with a data model can effectively improve the accuracy of altitude anomaly fitting. After calculation, the biquadratic interpolation method of the three interpolation methods has the highest accuracy, with an absolute error of within 3.3 cm for the four checkpoints.

1. 引言

随着GPS测量的广泛应用，高程测量脱离了需要步步设点的常规模型，但我国高程系统使用正常高，而GPS系统则使用正高，使得大地水准面和似大地水准面之间存在高程异常插值 [1] - [8] 。数学模型是高程异常拟合的常规方法，包括二次多项式、双三次多项式、移动三角面多项式、多面函数、BP神经网络、Shepard插值等，部分数学模型能得到四等水准的精度要求 [9] - [15] 。高程异常是由中长波和短波组成的，其中短波与当地地表地形密切关联，而中长波则可以通过重力场模型计算 [16] [17] [18] [19] [20] 。EGM2008是常用的重力场模型，其存在三种主要的插值，通过计算对比选择适用于试验区的精度最高的插值方法，结合多项式模型，进行基于移去–恢复法的似大地水准面构建，使得GPS测量数据能够直接应用到现有水准成果中。

2. 用于高程异常拟合的多项式模型

用于高程异常拟合的多项式拟合模型如式(1)：

$\xi =f\left(B,L\right)+\epsilon$ (1)

式中， $\xi$ 为高程异常， $f\left(B,L\right)=a_0+a_{1}B+a_{2}L+a_3{B}^2+a_4{L}^2+a_{5}BL+\cdots \cdots \cdots$ ， $\epsilon$ 为拟合残差，a_i为待拟合参数，根据不同的多项式次数确定i的数值，B，L为大地经纬度。

将式(1)改写为误差方程式如式(2)：

$V=BX-L$ (2)

式中， $V={\left[\begin{array}{cccc}{V}_0& V_1& \cdots & V_5\end{array}\right]}^T$ ， $X={\left[\begin{array}{cccc}{a}_0& a_1& L& a_5\end{array}\right]}^T$ ， $L={\left[\begin{array}{cccc}{\xi }_0& {\xi }_1& \cdots & {\xi }_5\end{array}\right]}^T$ ，为拟合系数矩阵， $B=\left[\begin{array}{cccccc}1& B_{11}& L_{12}& B_{13}^2& L_{14}^2& B_{15}{L}_{15}\\ ⋮& ⋮& ⋮& ⋮& ⋮& ⋮\\ 1& B_{51}& L_{52}& B_{53}^2& L_{54}^2& B_{55}{L}_{55}\end{array}\right]$ 。

二次多项式拟合中，至少需要六组已知点的高程异常及大地经纬度，收集已知数据后，将用于拟合模型的控制点数值带入，使用式(3)计算待估拟合参数X：

$X={(B^{T}B)}^{-1}{B}^{T}L$ (3)

计算X后，就是求得了a_i，回代式(1)，获得拟合模型。

3. 结合重力场模型的移去–恢复法

重力场模型中，主要有斯托克斯理论和莫洛金斯基理论，两种重力场模型都是通过建模进行扰动位的计算，斯托克斯理论将重力观测值归算至大地水准面，假设归算面外没有质量，解算归算面外部的扰动位，确定重力场。扰动位T如式(4)：

$T=\frac{R}{4\text{π}}{\displaystyle {\iint }_{\sigma }\Delta gS\left(\phi \right)\text{d}\sigma }$ (4)

大地水准高N见式(5)：

$N=\frac{T}{\gamma }=\frac{R}{4\text{π}\gamma }{\displaystyle {\iint }_{\sigma }\Delta gS\left(\phi \right)\text{d}\sigma }$ (5)

$S\left(\phi \right)=\frac{1}{S}-6S-4+10S^2-3\left(1-2S^2\right)\ln \left(S+S^2\right)$ (6)

$S=\sin \left(\frac{\phi }{2}\right)$ (7)

$\phi$ 为球面距离， $\gamma$ 为地球平均正常重力值， $\sigma$ 为单位球面， $S\left(\phi \right)$ 为斯托克斯函数。

莫洛金斯基理论将数据归算到参考椭球面而非大地水准面上，T如式(7)：

$T=T_0+T_1=\frac{R}{4\text{π}}{\displaystyle {\iint }_{\sigma }S\left(\phi \right)}\left(\Delta g+G_1\right)\text{d}\sigma$ (8)

$T_0=\frac{1}{4\text{π}R}{\displaystyle \int \Delta gS\left(\phi \right)\text{d}\sigma }$ ， $T_1=\frac{1}{4\text{π}R}{{\displaystyle \int }}^{\text{}}{G}_{1}S\left(\phi \right)\text{d}\sigma$ (9)

式中， $\sigma$ 为单位球面， $S\left(\phi \right)$ 为斯托克斯函数，G₁为地形改正。

Figure 1. Flow chart of removal restoration method combined with gravity field model

图1. 结合重力场模型的移去–恢复法流程图

计算重力场模型后，进行高程异常拟合，由于缺少高分辨率的DEM数据，将短波和残差项一起进行多项式拟合，移去–恢复法的步骤为，首先将控制点坐标分别拟合点和检核点，根据EGM2008重力场模型进行控制点的中长波项高程异常计算，使用拟合点的已知高程异常减去其中长波项高程异常，得到拟合点的待拟合部分，将待拟合部分作为式(1)中的 $\xi$ ，带入对应的控制点坐标，使用多项式模型计算待定系数X，然后回带式(1)形成高程异常拟合模型，使用高程异常拟合模型计算检核点的高程异常，加上检核点的中长波项高程异常计算值，得到经过拟合的检核点的高程异常值，与真值进行对比，验证精度，具体过程如图1所示。

4. 三种插值法的实例分析

实例使用的数据为22组重合点，选取4个检核点进行模型精度的验证，得到4个检核点的已知值以及三种插值方法计算的拟合值折线图如图2。

Figure 2. Line chart of fitted values calculated by three interpolation methods

图2. 三种插值方法计算的拟合值折线图

如图2，可知三种插值方法存在精度区别，其中与已知值最近的是双二次插值法，最远的是双线性插值法，计算三种插值方法的中误差值，双二次为3.3 cm，最临近为3.4 cm，双线性为3.6 cm。实例中进行了GPS联网的高程测量，实际水准路线达到了10公里，按照四等水准精度为每公里5 mm的要求，10公里限差为5 cm，即基于重力场模型构建的似大地水准面精度优于四等水准测量。

5. 结论

通过基于重力场模型的移去–恢复法，发现EGM2008模型的三种插值方法精度区别不是很大，但是以双二次为最优，在实际应用中，双二次插值法的每公里误差达到0.33 cm，远低于四等水准的限差值0.5 cm，结合多项式拟合能够满足四等水准测量的要求，在控制点要求精度不是很高的前提下，可以使用GPS测量配合重力场模型进行高程异常拟合。移去–恢复法能极大地缩小待拟合值的范围，使得大部分高程异常不需要进行拟合，进对短波部分进行拟合，提高了似大地水准面的拟合精度。下一步研究中，考虑使用高精度的DEM影像计算一部分短波异常，则能进一步缩小待拟合值的范围，获得更好的似大地水准面。