CN104613932A - 利用垂线偏差与重力异常确定似大地水准面模型的方法 - Google Patents

Abstract

本发明涉及利用垂线偏差与重力异常确定似大地水准面模型的方法，有效提高似大地水准面的建模精度，解决大地测量中高程基准面的确定问题，采用高程异常的微分关系，由格网垂线偏差与重力异常计算格网相对高程异常差；通过GPS/水准点施加约束，得到格网节点的高程异常，利用每个节点的垂线偏差、重力异常、数字高程数据，可以计算出相邻节点间的相对高程异常差；计算相邻格网之间的高程异常差，采用参数平差模型进行平差计算，确定似大地水准面模型；选择3个百万图幅4°×6°区域，分辨率为的垂线偏差、重力异常以及地形格网数据，分别计算建立3个区域的格网似大地水准面模型，分析精度，本发明方法新颖独特，效率高，精度高，是大地测量上的一大创新。

Description

利用垂线偏差与重力异常确定似大地水准面模型的方法

技术领域

本发明涉及大地测量领域，特别是一种利用垂线偏差与重力异常确定似大地水准面模型的方法。

背景技术

似大地水准面是从地面点量算正常高所构成的一个几何曲面，可以看作是正常高系统中地面点高程的起算面。传统的高程基准是以水准测量建立并以水准标石所体现，但在当今，这种传统方式正在改变。通过重力场数据建立的格网大地水准面所体现的数字高程基准已成为发展趋势。它为应用卫星定位技术进行三维测量提供了必要基础，也避免了水准标石的维护困难。建立高分辨率高精度似大地水准面模型已成为当前大地测量主要任务之一。

我国是世界上重力场最复杂的地区，但陆地重力测量目前大约只有80余万点，主要是地质部门测量的，用于大地测量的重力点不足20万点，而且点的分布非常不均匀，东部地区较为密集，西部地区和山区比较稀疏，约40％的5′×5′格网无实测重力点。

目前通用的似大地水准面建模方法是按照Molodensky理论联合地球重力场模型、重力、地形等数据直接计算一点的高程异常再用GPS/水准拟合系统差，由于似大地水准面不能有效地吸收地形等高频信息并且难以利用GPS/水准数据对似大地水准面作真正地实际控制，因此很难再大幅度提高似大地水准面精度。似大地水准面的确定方法主要还是依赖密集的GPS/水准，这不但成本高，而且由于采用的重力数据格网分辨率低，导致似大地水准面模型的分辨率也低。

由于垂线偏差对高频地形信息比较敏感，因此首先可以依据Molodensky理论联合地球重力场模型、重力、地形等数据建立高精度的垂线偏差数值模型，然后利用垂线偏差与重力异常，基于天文重力水准的方法计算得到相对高程异常差，再对其施加GPS/水准控制点予以约束，最后采用平差的方法获得高程异常，从而建立高精度高分辨率的似大地水准面模型，但至今未见有利用垂线偏差与重力异常确定似大地水准面模型方法的公开报导。

发明内容

针对上述情况，为克服现有技术之缺陷，本发明之目的就是提供一种利用垂线偏差与重力异常确定似大地水准面模型的方法，可有效提高似大地水准面的建模精度，用于解决大地测量中高程基准面的确定问题。

本发明解决的技术方案是，由以下步骤实现：

一、格网相对高程异常差计算：采用高程异常的微分关系，由格网垂线偏差与重力异常计算格网相对高程异常差；

二、GPS/水准控制及格网平差：通过GPS/水准点施加约束，得到格网节点的高程异常，利用每个节点的垂线偏差、重力异常、数字高程数据，可以计算出相邻节点间的相对高程异常差；

三、平差模型：计算相邻格网之间的高程异常差，采用参数平差模型进行平差计算，确定似大地水准面模型；

四、分析区域格网似大地水准面模型的精度：选择3个百万图幅4°×6°区域，分辨率为1′×1′的垂线偏差、重力异常以及地形格网数据，分别计算建立3个区域的格网似大地水准面模型，分析精度。

本发明方法新颖独特，利用垂线偏差与重力异常确定相对高程异常差，然后施加GPS/水准点进行控制，采用平差计算得到高程异常，从而确定似大地水准面模型，有效用于大地测量中高程基准面的确定，效率高，精度高，是大地测量上的一大创新。

附图说明

图1为本发明利用垂线偏差与重力异常确定似大地水准面模型框示图。

图2为本发明GPS/水准控制网图。

图3为本发明的平差计算格网角点的高程异常图。

图4为本发明利用垂线偏差计算相对高程异常差图。

图5为本发明区域一GPS/水准点点位分布情况图。

图6为本发明计算高程异常与实测高程异常之间的差值图(单位：米)。

图7为本发明区域二GPS/水准点点位分布情况图。

图8为本发明计算高程异常与实测高程异常之间的差值图(单位：米)。

图9为本发明区域三GPS/水准点点位分布情况图。

图10为本发明计算高程异常与实测高程异常之间的差值图(单位：米)。

图11为本发明分区域检核点计算误差曲线图。

图12为本发明分区域误差曲线示意图。

具体实施方式

以下结合附图和具体情况对本发明的具体实施方式作详细说明。

由图1所示，本发明在具体实施中，包括设置GPS/水准点、高程异常控制点，格网高程数据、格网垂直偏差、格网重力异常的计算，获得相对高程异常差，与高程异常控制点经平差计算获得高程异常，从而构建似大地水准面模型，具体由以下步骤实现：

一、格网相对高程异常差计算：按照高程异常的微分关系由格网垂线偏差与重力异常计算格网相对高程异常差，其计算公式为：

$\mathrm{d\ζ}=-(\mathrm{\ξ}\cos A+\mathrm{\η}\sin A)\mathrm{ds}-\frac{\mathrm{\Δg}}{\stackrel{\‾}{\mathrm{\γ}}}\mathrm{dh}$               式(1)

式中，ξ为垂线偏差子午分量，η为垂线偏差卯酉分量，Δg为重力异常，dh为两点间高差，为重力异常改正项，为平均正常重力值；

二、GPS/水准控制及格网平差：采用垂线偏差和重力异常计算只能获得格网相对高程异常差，要想获得高程异常，须通过GPS/水准点施加约束，方法是，设m行、n列的分辨率为1′×1′的区域，格网中的A、B、C、D点表示GPS/水准控制点，每50公里一个GPS/水准控制点，实测的GPS/水准控制点往往不在格网的节点处，通过插值拟合的方法归算到控制点所在格网的四个角点上，由这四个角点作为控制点(如图2所示)；

计算平差中控制点的高程异常，方法是：

1、平差计算GPS/水准控制点所在网格四个角点的高程异常；

(a)计算控制点到该网格角点的高程异常差，计算公式为：

$\mathrm{\Δ\ζ}=-\mathrm{\υ\Δs}-\frac{\mathrm{\Δg}}{\stackrel{\‾}{\mathrm{\γ}}}\mathrm{dh}$            式(2)

其中υ＝ξ_中cosα+η_中sinα，υ为格网节点与GPS/水准点的垂线偏差的平均值，ξ_中为格网节点与GPS/水准点子午方向垂线偏差的平均值，η_中为格网节点与GPS/水准点卯酉方向垂线偏差的平均值，α为该方向的方位角；

(b)以GPS/水准控制点的高程异常作为已知点，平差得到该格网四个角点的高程异常

以GPS/水准高程异常ζ₀为已知点，按水准网平差的方法组建误差方程组，求解法方程，得到网格四个角点的高程异常ζ₁,ζ₂,ζ₃,ζ₄(如图3所示)；

2、同理可以求得该区域中其它GPS/水准点所在格网四个角点的高程异常，将这些平差求得的角点作为控制点，对整个网进行约束，平差求得其它节点的高程异常值，进而得到该区域的似大地水准面模型；

利用每个节点的垂线偏差、重力异常、数字高程数据，计算出相邻节点间的相对高程异常差，相对高程异常差相当于水准网中的相邻点间的高差，而控制点的高程异常则可认为是水准网中的已知点高程，利用垂线偏差计算相对高程异常差(如图4)，子午方向和卯酉方向相邻节点的高程异常差分别为：

Δζ_子午＝ζ_上-ζ_下＝-1/2(ξ_上+ξ_下)Δs＝-ξ_中Δs              式(3)

Δζ_卯酉＝ζ_右-ζ_左＝-1/2(η_右+η_左)Δs＝-η_中Δs              式(4)

考虑重力异常改正项，各节点之间的高程异常差为：

                          式(5)

                          式(6)

式中ξ_中，η_中为相邻节点子垂线偏差午分量和卯酉分量的平均值，Δg_中为相邻两节点的重力异常平均值，Δh为相邻两节点的高程差，为平均正常重力值i，j分别为节点行、列序号，一个m×n个格网区域共包含(m-1)×n个子午方向高程异常差和m×(n-1)个卯酉方向高程异常差；

三、平差模型

利用式(5)与式(6)计算相邻格网之间的高程异常差，以GPS/水准点提供已知点的高程异常，作为起算点进行控制，采用平差计算得到格网中每个节点的高程异常，为便于编程实现，采用参数平差模型进行平差计算；

设为格网节点高程异常的平差改正数，为格网节点高程异常的近似值，对于m行×n列格网沿纬线方向，误差方程为：

$V_{i,j}^h=\mathrm{\δ}{\hat{X}}_{i,j+1}-\mathrm{\δ}{\hat{X}}_{i,j}+X_{i,j+1}^0-X_{i,j}^0-L_{i,j}^h$                       式(7)

式中，为卯酉方向观测量(垂线偏差等数据计算得到的高程异常差)的残差，为卯酉方向观测量(由垂线偏差等数据计算的卯酉方向相邻两结点的高程异常差)；

同样，对于相邻两条纬线之间的每条经线方向，列出下列误差方程：

$V_{i,j}^l=\mathrm{\δ}{\hat{X}}_{i+1,j}-\mathrm{\δ}{\hat{X}}_{i,j}+X_{i+1,j}^0-X_{i,j}^0-L_{i,j}^l$                        式(8)

式中，为子午方向观测量(垂线偏差等数据计算得到的高程异常差)的残差，为子午方向观测量(由垂线偏差等数据计算的子午方向相邻两节点的高程异常差)；

组成误差方程：

$\left\{\begin{array}{c}{V}_{\mathrm{1,1}}^h=\mathrm{\δ}{\hat{X}}_{\mathrm{1,2}}-\mathrm{\δ}{\hat{X}}_{\mathrm{1,1}}+X_{\mathrm{1,1}}^0-L_{\mathrm{1,1}}^h\\ .\\ .\\ .\\ V_{1,j}^h=\mathrm{\δ}{\hat{X}}_{1,j+1}-\mathrm{\δ}{\hat{X}}_{1,j}+X_{1,j+1}^0-X_{1,j}^0-L_{1,j}^h\\ .\\ .\\ .\\ V_{1,n-1}^h=\mathrm{\δ}{\hat{X}}_{1,n}-\mathrm{\δ}{\hat{X}}_{1,n-1}+X_{1,n}^0-X_{1,n-1}^0-L_{1,n-1}^h\\ V_{2,1}^h=\mathrm{\δ}{\hat{X}}_{\mathrm{2,2}}-\mathrm{\δ}{\hat{X}}_{\mathrm{2,1}}+X_{\mathrm{2,2}}^0-X_{\mathrm{2,1}}^0-L_{\mathrm{2,1}}^h\\ .\\ .\\ .\\ V_{i,j}^h=\mathrm{\δ}{\hat{X}}_{i,j+1}-\mathrm{\δ}{\hat{X}}_{i,j}+X_{i,j+1}^0-X_{i,j}^0-L_{i,j}^h\\ .\\ .\\ .\\ V_{m,n1}^h=\mathrm{\δ}{\hat{X}}_{m,n}-\mathrm{\δ}{\hat{X}}_{m,n-1}+X_{m,n}^0-X_{m,n-1}^0-L_{m,n-1}^h\\ V_{\mathrm{1,1}}^l=\mathrm{\δ}{\hat{X}}_{\mathrm{2,1}}-\mathrm{\δ}{\hat{X}}_{\mathrm{1,1}}+X_{\mathrm{2,1}}^0-X_{\mathrm{1,1}}^0-L_{\mathrm{1,1}}^l\\ V_{\mathrm{1,2}}^l=\mathrm{\δ}{\hat{X}}_{\mathrm{2,2}}-\mathrm{\δ}{\hat{X}}_{\mathrm{1,2}}+X_{\mathrm{2,2}}^0-X_{\mathrm{1,2}}^0-L_{\mathrm{1,2}}^l\\ .\\ .\\ .\\ V_{i,j}^l=\mathrm{\δ}{\hat{X}}_{i+1,j}-\mathrm{\δ}{\hat{X}}_{i,j}+X_{i+1,j}^0-X_{i,j}^0-L_{i,j}^l\\ .\\ .\\ .\\ V_{m-1,n}^l=\mathrm{\δ}{\hat{X}}_{m,n}-\mathrm{\δ}{\hat{X}}_{m-1,n}+X_{m,n}^0-X_{m-1,n}^0-L_{m-1,n}^l\end{array}\right.$        式(9)

即

$V=\mathrm{A\δ}\hat{X}+l$                          式(10)

其中i＝1,2，...m，j＝1,2，...n，其法方程为

$\mathrm{N\δ}\hat{X}+U=0$                         式(11)

以上为自由网平差模型，缺少基准，因此需要施加一定数量的GPS/水准点进行控制，才能计算出每个节点的高程异常，从而确定似大地水准面模型，为不破坏系数阵的规律性以及易于编程实现，将控制点作为“未知点”，但是赋予其大的权，即对法方程系数阵N中与控制点相对应的元素加上一个大的数10e30，从而在平差计算中控制点的值不会产生变化，这样既顾及系数阵的规律性又利用自由网平差的模型，方便易行；

同时为提高计算速度与精度，采用下式计算格网节点高程异常的平差值

$\mathrm{\δ}\hat{X}=-N/U$                       式(12)

${\hat{X}}_{i,j}=\mathrm{\δ}\hat{X}+X_{i,j}^0$                      式(13)

其中，N＝A^TPA，U为法方程自由项向量，U＝A^TPl，P为观测量的权阵，为格网结点高程异常的平差值；

四、分析区域格网似大地水准面模型的精度

选择3个百万图幅4°×6°区域，分辨率为1′×1′的垂线偏差、重力异常以及地形格网数据，分别计算建立3个区域的格网似大地水准面模型，第一区域内共布测212个实测GPS/水准控制点，其中选择37个在区域内均匀分布的点作为计算控制点，实测GPS/水准控制点和计算控制点采用混合穿插分布，平差计算与水准面拟合采用相同的控制点，剩余175个点作为计算检核点，首先利用构建好的格网似大地水准面模型，内插得到剩余175个GPS/水准点的高程异常值，然后与其实测的高程异常值求差，对差值进行统计计算，分析评定第一区域格网似大地水准面模型的精度；

第二区域内布测292个实测GPS/水准控制点，其中选择42个在区域内均匀分布的点作为计算控制点，实测GPS/水准控制点和计算控制点采用混合穿插分布，平差计算与水准面拟合采用相同的控制点，剩余250个点作为计算检核点，首先利用构建好的格网似大地水准面模型，内插得到剩余250个GPS/水准点的高程异常值，然后与其实测的高程异常值求差，对差值进行统计计算，分析评定第二区域格网似大地水准面模型的精度

第三区域内共布测了256个实测GPS/水准控制点，其中选择45个在区域内均匀分布的点作为计算控制点，实测GPS/水准控制点和计算控制点采用混合穿插分布，平差计算与水准面拟合采用相同的控制点，剩余211个点作为计算检核点，首先利用构建好的格网似大地水准面模型，内插得到剩余211个GPS/水准点的高程异常值，然后与其实测的高程异常值求差，对差值进行统计计算，分析评定第三区域格网似大地水准面模型的精度。

本发明在经实际应用，取得了非常好的效果，精度高，精度在10厘米左右，为了证明本发明的有效性，选取1个百万图幅(4°×6°)区域，分辨率为1′×1′的垂线偏差、重力异常以及地形格网数据，按本发明方法计算高程异常差，进行平差统计，结果如表1、表2所示。

表1相对高程异常差的平差残差值及中误差

同一区域，加入GPS/水准点进行控制

表2相对高程异常差的平差残差值及中误差

由表1、表2可以看出，加入GPS/水准控制点后，相对高程异常差的残差的绝对值较自由网平差要大，参数平差值的改正数最大绝对值为分米级，这是因为利用垂线偏差与重力异常等数据计算获得的似大地水准面与由GPS/水准计算得到的似大地水准面之间存在系统误差，但是由中误差可以看出，其内符合精度都比较高。

同时，选择3个百万图幅(4°×6°)区域，分辨率为1′×1′的垂线偏差、重力异常以及地形格网数据，采用本发明方法，分别计算建立3个区域的格网似大地水准面模型，实验采用Matlab 2012b平台，计算结果分别如下：

区域一内共布测了212个实测GPS/水准控制点，其中选择37个在区域内均匀分布的点作为计算控制点(平差计算与水准面拟合采用相同的控制点)，剩余175个点作为计算检核点，点位分布情况如图5(星星为控制点，圆点为检核点)：

为了评价依据本发明所构建的格网似大地水准面模型的精度，首先利用构建好的格网似大地水准面模型，内插得到剩余175个GPS/水准点的高程异常值，然后与其实测的高程异常值求差，差值结果如图6所示，最后对差值进行统计计算，精度评定结果列于表3。

表3基于垂线偏差与高程异常建立的似大地水准面模型的精度(单位：厘米)

区域二内共布测292个实测GPS/水准控制点，其中选择42个在区域内均匀分布的点作为计算控制点(平差计算与水准面拟合采用相同的控制点)，剩余250个点作为计算检核点，点位分布情况如图7(星星为控制点，圆点为检核点)，具体实施方案及精度评定方法同区域一，计算高程异常与实测高程异常之间的差值及差值的统计结果见图8与表4。

表4基于垂线偏差与高程异常建立的似大地水准面模型的精度(单位：厘米)

区域三内共布测了256个实测GPS/水准控制点，其中选择45个在区域内均匀分布的点作为计算控制点(平差计算与水准面拟合采用相同的控制点)，剩余211个点作为计算检核点，点位分布情况如图9(星星为控制点，圆点为检核点)，具体实施方案及精度评定方法同区域一，计算高程异常与实测高程异常之间的差值及差值的统计结果见图10与表5。

表5基于垂线偏差与高程异常建立的似大地水准面模型的精度(单位：厘米)

以上结果清楚的表明，本发明所采用的利用垂线偏差和高程异常计算相对高程异常差，进而求解高程异常确定似大地水准面的方法，与实测点进行比较其精度在10厘米左右。

为了进一步说明确定似大地水准面的精度，采用同样的数据源以及精度评价方法，利用两种不同技术方案建立了我国的似大地水准面模型。一种是依据目前通用的基于Molodensky理论直接计算高程异常，然后利用GPS/水准点拟合修正，建立的似大地水准面模型精度如图11及表6所示；另一种是先依据Molodensky理论计算垂线偏差，然后按照本发明方法利用垂线偏差与重力异常计算相对高程异常差，通过GPS/水准予以控制，最后通过平差计算建立我国的似大地水准面模型，建立的似大地水准面模型的精度如图12及表7所示。

表6分区域误差统计特征(cm)

表7分区域误差统计特征(cm)

由图11、图12以及表6、表7可以看出，本发明模型的精度优于传统模型精度，具有实际的应用价值。

由上述可以看出，本发明是提出基于垂线偏差与重力异常确定我国似大地水准面的方法。利用垂线偏差对高频地形信息敏感的特点，以高分辨率地形数据弥补重力场数据的不足，以确定相对高程异常差代替以往的直接确定高程异常的方法，类似于GPS从直接定位到差分定位的转变，而且这种格网相对高程异常构成了严密的几何条件，可以通过平差消除矛盾提高精度，比以往采用数值拟合的简单方法提高了严密度。根据估计，不需要专门测量GPS/水准控制点，只是在现有国家三级(C级)GPS/水准网的控制下，相对控制点的高程异常精度可以优于±3cm，有效用于大地测量领域高程基准面的建立，精度高，具有很强的实际应用价值，是大地测量上的一大创新。

Claims (2)

1.一种利用垂线偏差与重力异常确定似大地水准面模型的方法，其特征在于，由以下步骤实现：

一、格网相对高程异常差计算：采用高程异常的微分关系，由格网垂线偏差与重力异常计算格网相对高程异常差；

二、GPS/水准控制及格网平差：通过GPS/水准点施加约束，得到格网节点的高程异常，利用每个节点的垂线偏差、重力异常、数字高程数据，可以计算出相邻节点间的相对高程异常差；

三、平差模型：计算相邻格网之间的高程异常差，采用参数平差模型进行平差计算，确定似大地水准面模型；

四、分析区域格网似大地水准面模型的精度：选择3个百万图幅4°×6°区域，分辨率为1′×1′的垂线偏差、重力异常以及地形格网数据，分别计算建立3个区域的格网似大地水准面模型，分析精度。

2.根据权利要求1所述的利用垂线偏差与重力异常确定似大地水准面模型的方法，其特征在于，由以下步骤实现：

一、格网相对高程异常差计算：按照高程异常的微分关系由格网垂线偏差与重力异常计算格网相对高程异常差，其计算公式为：

$\mathrm{d\ζ}=-(\mathrm{\ξ}\cos A+\mathrm{\η}\sin A)\mathrm{ds}-\frac{\mathrm{\Δg}}{\stackrel{\‾}{\mathrm{\γ}}}\mathrm{dh}$   式(1)

式中，ξ为垂线偏差子午分量，η为垂线偏差卯酉分量，Δg为重力异常，dh为两点间高差，为重力异常改正项，为平均正常重力值，可以取9.797644656m/s²；

二、GPS/水准控制及格网平差：采用垂线偏差和重力异常计算只能获得格网相对高程异常差，要想获得高程异常，必须通过GPS/水准点施加约束，方法是，设m行、n列的分辨率为1′×1′的区域，格网中的A、B、C、D点表示GPS/水准控制点，每50公里一个GPS/水准控制点，实测的GPS/水准控制点往往不在格网的节点处，通过插值拟合的方法归算到控制点所在格网的四个角点上，由这四个角点作为控制点；

计算平差中控制点的高程异常，方法是：

(1)、平差计算GPS/水准控制点所在格网四个角点的高程异常；

(a)计算控制点到该格网角点的高程异常差，计算公式为：

$\mathrm{\Δ\ζ}=-\mathrm{\υ\Δs}-\frac{\mathrm{\Δg}}{\stackrel{\‾}{\mathrm{\γ}}}\mathrm{dh}$   式(2)

其中υ＝ξ_中cosα+η_中sinα，υ为格网节点与GPS/水准点的垂线偏差的平均值，ξ_中为格网节点与GPS/水准点子午方向垂线偏差的平均值，η_中为格网节点与GPS/水准点卯酉方向垂线偏差的平均值，α为该方向的方位角；

(b)以GPS/水准控制点的高程异常作为已知点，平差得到该格网四个角点的高程异常

以GPS/水准高程异常ζ₀为已知点，按水准网平差的方法组建误差方程组，求解法方程，得到网格四个角点的高程异常ζ₁,ζ₂,ζ₃,ζ₄；

(2)、同理可以求得该区域中其它GPS/水准点所在格网四个角点的高程异常，将这些平差求得的角点作为控制点，对整个网进行约束，平差求得其它节点的高程异常值，进而得到该区域的似大地水准面模型；

利用每个节点的垂线偏差、重力异常、数字高程数据，计算出相邻节点间的相对高程异常差，相对高程异常差相当于水准网中的相邻点间的高差，而控制点的高程异常则可认为是水准网中的已知点高程，利用垂线偏差计算相对高程异常差，子午方向和卯酉方向相邻节点的高程异常差分别为：

Δζ_子午＝ζ_上-ζ_下＝-1/2(ξ_上+ξ_下)Δs＝-ξ_中Δs  式(3)

Δζ_卯酉＝ζ_右-ζ_左＝-1/2(η_右+η_左)Δs＝-η_中Δs  式(4)

考虑重力异常改正项，各节点之间的高程异常差为：

  式(5)

  式(6)

式中ξ_中，η_中为相邻节点子垂线偏差午分量和卯酉分量的平均值，Δg_中为相邻两节点的重力异常平均值，Δh为相邻两节点的高程差，为平均正常重力值i，j分别为节点行、列序号，一个m×n个格网区域共包含(m-1)×n个子午方向高程异常差和m×(n-1)个卯酉方向高程异常差；

三、平差模型

利用式(5)与式(6)计算相邻格网之间的高程异常差，以GPS/水准点提供已知点的高程异常，作为起算点进行控制，采用平差计算得到格网中每个节点的高程异常，为便于编程实现，采用参数平差模型进行平差计算；

设为格网节点高程异常的平差改正数，为格网结点高程异常的近似值，对于m行×n列格网沿纬线方向，误差方程为：

$V_{i,j}^h=\mathrm{\δ}{\hat{X}}_{i,j+1}-\mathrm{\δ}{\hat{X}}_{i,j}+X_{i,j+1}^0-X_{i,j}^0-L_{i,j}^h$   式(7)

式中，为卯酉方向观测量的残差，为卯酉方向观测量，由垂线偏差等数据计算的卯酉方向相邻两结点的高程异常差；

同样，对于相邻两条纬线之间的每条经线方向，列出下列误差方程：

$V_{i,j}^l=\mathrm{\δ}{\hat{X}}_{i+1,j}-\mathrm{\δ}{\hat{X}}_{i,j}+X_{i+1,j}^0-X_{i,j}^0-L_{i,j}^l$   式(8)

式中，为子午方向观测量的残差，为子午方向观测量，由垂线偏差等数据计算的子午方向相邻两节点的高程异常差；

组成误差方程：

$\begin{array}{c}{V}_{\mathrm{1,1}}^h=\mathrm{\δ}{\hat{X}}_{\mathrm{1,2}}-\mathrm{\δ}{\hat{X}}_{\mathrm{1,1}}+X_{\mathrm{1,2}}^0-X_{\mathrm{1,1}}^0-L_{\mathrm{1,1}}^h\\ .\\ .\\ .\\ V_{1,j}^h=\mathrm{\δ}{\hat{X}}_{1,j+1}-\mathrm{\δ}{\hat{X}}_{1,j}+X_{1,j+1}^0-X_{1,j}^0-L_{1,j}^h\\ .\\ .\\ .\\ V_{1,n-1}^h=\mathrm{\δ}{\hat{X}}_{1,n}-\mathrm{\δ}{\hat{X}}_{1,n-1}+X_{1,n}^0-X_{1,n-1}^0-L_{1,n-1}^h\\ V_{\mathrm{2,1}}^h=\mathrm{\δ}{\hat{X}}_{\mathrm{2,2}}-\mathrm{\δ}{\hat{X}}_{\mathrm{2,1}}+X_{\mathrm{2,2}}^0-X_{\mathrm{2,1}}^0-L_{\mathrm{2,1}}^h\\ .\\ .\\ .\\ V_{i,j}^h=\mathrm{\δ}{\hat{X}}_{i,j+1}-\mathrm{\δ}{\hat{X}}_{i,j}+X_{i,j+1}^0-X_{i,j}^0-L_{i,j}^h\\ .\\ .\\ .\\ V_{m,n-1}^h=\mathrm{\δ}{\hat{X}}_{m,n}-\mathrm{\δ}{\hat{X}}_{m,n-1}+X_{m,n}^0-X_{m,n-1}^0-L_{m,n-1}^h\\ V_{\mathrm{1,1}}^l=\mathrm{\δ}{\hat{X}}_{\mathrm{2,1}}-\mathrm{\δ}{\hat{X}}_{\mathrm{1,1}}+X_{\mathrm{2,1}}^0-X_{\mathrm{1,1}}^0-L_{\mathrm{1,1}}^l\\ V_{\mathrm{1,2}}^l=\mathrm{\δ}{\hat{X}}_{\mathrm{2,2}}-\mathrm{\δ}{\hat{X}}_{\mathrm{1,2}}+X_{\mathrm{2,2}}^0-X_{\mathrm{1,2}}^0-L_{\mathrm{1,2}}^l\\ .\\ .\\ .\\ V_{i,j}^l=\mathrm{\δ}{\hat{X}}_{i+1,j}-\mathrm{\δ}{\hat{X}}_{i,j}+X_{i+1,j}^0-X_{i,j}^0-L_{i,j}^l\\ .\\ .\\ .\\ V_{m-1,n}^l=\mathrm{\δ}{\hat{X}}_{m,n}-\mathrm{\δ}{\hat{X}}_{m-1,n}+X_{m,n}^0-X_{m-1,n}^0-L_{m-1,n}^l\end{array}$   式(9)

即

$V=\mathrm{A\δ}\hat{X}+l$   式(10)

其中i＝1,2，…m，j＝1,2，…n，其法方程为

$\mathrm{N\δ}\hat{X}+U=0$   式(11)

以上为自由网平差模型，缺少基准，因此需要施加一定数量的GPS/水准点进行控制，才能计算出每个节点的高程异常，从而确定似大地水准面模型，为不破坏系数阵的规律性以及易于编程实现，将控制点作为“未知点”，但是赋予其大的权，即对法方程系数阵N中与控制点相对应的元素加上一个大的数10e30，从而在平差计算中控制点的值不会产生变化，这样既顾及系数阵的规律性又利用自由网平差的模型，方便易行；

同时为提高计算速度与精度，采用下式计算格网节点高程异常的平差值

$\mathrm{\δ}\hat{X}=-N/U$   式(12)

${\hat{X}}_{i,j}=\mathrm{\δ}\hat{X}+X_{i,j}^0$   式(13)

其中，N＝A^TPA，U为法方程自由项向量，U＝A^TPl，P为观测量的权阵，为格网结点高程异常的平差值；

四、分析区域格网似大地水准面模型的精度：

选择3个百万图幅4°×6°区域，分辨率为1′×1′的垂线偏差、重力异常以及地形格网数据，分别计算建立3个区域的格网似大地水准面模型，第一区域内共布测212个实测GPS/水准控制点，其中选择37个在区域内均匀分布的点作为计算控制点，实测GPS/水准控制点和计算控制点采用混合穿插分布，平差计算与水准面拟合采用相同的控制点，剩余175个点作为计算检核点，首先利用构建好的格网似大地水准面模型，内插得到剩余175个GPS/水准点的高程异常值，然后与其实测的高程异常值求差，对差值进行统计计算，分析评定第一区域格网似大地水准面模型的精度；

第二区域内布测292个实测GPS/水准控制点，其中选择42个在区域内均匀分布的点作为计算控制点，实测GPS/水准控制点和计算控制点采用混合穿插分布，平差计算与水准面拟合采用相同的控制点，剩余250个点作为计算检核点，首先利用构建好的格网似大地水准面模型，内插得到剩余250个GPS/水准点的高程异常值，然后与其实测的高程异常值求差，对差值进行统计计算，分析评定第二区域格网似大地水准面模型的精度；

第三区域内共布测了256个实测GPS/水准控制点，其中选择45个在区域内均匀分布的点作为计算控制点，实测GPS/水准控制点和计算控制点采用混合穿插分布，平差计算与水准面拟合采用相同的控制点，剩余211个点作为计算检核点，首先利用构建好的格网似大地水准面模型，内插得到剩余211个GPS/水准点的高程异常值，然后与其实测的高程异常值求差，对差值进行统计计算，分析评定第三区域格网似大地水准面模型的精度。