GPS高程拟合方法3.1等值线图示法等值线图示法是最直接的求算高程异常的方法。这种方法的核心思想就是内插的思想,绘制高程异常的等值线图,然后采用内插法来确定未知点的高程异常值。具体操作十分的简单,在测区内制定分布均匀的GPS点,用水准测量的方法来测定这些点的水准高,根据公式ζ=H-Hr求出这些点的高程异常,选择适当的比例尺按照已知点的平面坐标展会在图纸内,对已知点标注出高程异常值,再确定等高距,绘制出高程异常值的等值线图。之后就可以内插出待测点的高程异常值,进而求出待测点的正常高。这种方法只适用地形相对平坦的地方,在此种测区内采用这种方法拟合的高程精度可达到厘米级。测区的地形相对复杂内插出的高程异常值就不准确,而且这种内插法的精度往往取决于两个方面,分别是测区内GPS点的分布密度和已知点大地高的精确度。首先GPS点的分布比较密集,那么内插精度就相对较高,如果比较稀疏这时候就要借助于此测区的重力测量资料,提高内插精度。且还要注意GPS点间高程异常的非线性变化。另外就是水准点的精度,联测时尽量选取高精度的正常高,尽可能使得出的高程异常值准确,进而才能内插出待测点高精度的高程异常值。这种方法虽然简单易操作,但是有其弱点,就是精度不高,只有当对拟合精度要求不高的时候才使用此种方法(注:等值线法不需构造数学模型)。3.2狭长带状区域线性拟合解析内插法作为拟合高程最常用的方法,主要思想是把似大地水准面用数学曲面近似拟合,建立所在测区内最为接近似大地水准面的数学模型,以此来计算测区内任意点的高程异常值,从而计算出正常高。这种方法计算出的高程异常值的精度是由所采用的数学模型和似大地水准面的拟合程度所决定的。解析内插法在选择数学模型时,首先要考虑的就是GPS点的分布情况。GPS点的分布情况可分为带状分布和面状分布。若GPS点是呈线状布设,而且是以沿线似大地水准面为一条连续且光滑的曲线,这时就可以采用相对于狭长带状区域的解析内插法来内插出待定点的高程异常值,从而求出待定点的正常高。这种线状分布的内插原理是:测区内已知水准点,用GPS测出其GPS高程,计算出已知水准点高程异常值,根据已知点的平面坐标和计算得出的高程异常值,构造出一个插值函数,这个函数是用来拟合GPS分布线上的似大地水准面的。用这个函数内插出位置点的高程异常值。下面是两种用来拟合线状分布的GPS高程的内插法。3.2.1多项式曲线拟合法多项式曲线拟合是线状分布拟合的主要方法。多项式拟合顾名思义其插值函数是一个m次的代数多项式,若高程控制点的高程异常为ζ,坐标为𝑥𝑙(或𝑦𝑖或𝑑𝑖或拟合坐标或𝑥𝑖-𝑥0或𝑦𝑖-𝑦0)的函数关系为下式:𝜁𝑙=𝑎0+𝑎1𝑥𝑙1+𝑎2𝑥𝑙2+𝑎3𝑥𝑙3+……𝑎𝑚𝑥𝑙𝑚3-1各高程控制点的已知高程异常与其拟合值之差为下式所示:𝑟𝑖=𝜁𝑙(x)-𝜁𝑖(i=0,1,2…n)3-2上式我们称之为离差。(3-1)中𝑥𝑙是拟合点到参考点(𝑥0,𝑦0)的直线距离,𝑥0,𝑦0为设定的常数值。在一般情况下都认为,𝑥0,𝑦0就是测区内已知点坐标的均值。多项式曲线拟合使用起来非常方便,但是它有自身的局限性,即是使用这种方法的时候,所测路线不能太长,要限制控制点到测点的距离不能太远,通常把距离控制在300米以内。这个要求是因为使用多项式曲线方法拟合似大地水准面,如果它拟合的范围太大,点位的高程异常变化就越复杂,削高补低的方法不能满足我们所要求的精度。随着多项式阶数的增大,也会使拟合出的曲线振荡的更厉害,从而造成拟合的误差增大。这些造成了多项式曲线拟合的缺陷,但是在路线较短的情况下,这种方法有足够的精度来拟合GPS点的正常高程。在式(3-1)中用m次多项式拟合似大地水准面,这个m的值如何取定,一般情况下如果测区不是很长,地形相对平坦,那么我们通常取m取为3。也就是说多项式为三次多项式。若测区比较长或者是测区地形比较复杂就要依情况而定,增加多项式的次数,提高拟合精度。依上述分析m的取值主要和测区长度以及测区的复杂程度有关。3.2.2三次样条曲线拟合法三次样条曲线拟合法针对测线长,已知点多的测区GPS高程拟合问题。由上述可以知道,当测线比较长已知点较多的时候,就需要构造高次的拟合多项式,当m值比较高的情况下,会出现不稳定的现象,对求解高程异常值会有比较大的影响,并且最小二乘法在求多项式系数中也会增大削高补低的误差,因此为了避免测线长、已知点多这种情况下所出现的问题,通常采用分段拟合的方法,采用三次样条函数拟合数学模型。这种方法很好的解决了因测线长而引起的问题。三次样条曲线的实质就是一个拼接而成的连续函数,在把测线分为多段的情况下,每段设为三次多项式函数,然后将这些多项式函数组成三次样条函数。为了计算准确,应用中要求这种构造的曲线不仅在连接点处函数要连续。而且还要求这个函数的一级导数还有二级导数全部要是连续的,才能保证在分段之后构造的三次样条函数后期运算中能够计算出准确的高程异常值。设过n个已知点,𝜁𝑖和𝑥𝑖(或𝑦𝑖或拟合坐标)在区间[𝑥𝑖,𝑥𝑖+1](i=1,2,…,n-1)上有三次样条函数关系:ζ(x)=ζ(𝑥𝑖)+(x-𝑥𝑖)ζ(𝑥𝑖,𝑥𝑖+1)+(x-𝑥𝑖)(𝑥−𝑥𝑖+1)ζ(x,𝑥𝑖,𝑥𝑖+1)3-3式中,x为待定点坐标,ζ(𝑥𝑖,𝑥𝑖+1)为一阶差商,ζ(𝑥𝑖,𝑥𝑖+1)=(𝜁𝑖+1−𝜁𝑖)/(𝑥𝑖+1-𝑥𝑖);ζ(x,𝑥𝑖,𝑥𝑖+1)为二阶差商,ζ(x,𝑥𝑖,𝑥𝑖+1)=16[𝜁,,(𝑥𝑖)+𝜁,,(x)+𝜁,,(𝑥𝑖+1),而𝜁,,(x)(i=1,2,…,n-1),满足系数矩阵为对称三角阵的线性方程组(𝑥𝑖-𝑥𝑖+1)𝜁,,(𝑥𝑖−1)+2(𝑥𝑖+1,𝑥𝑖−1)𝜁,,(𝑥𝑖)+(𝑥𝑖+1,𝑥𝑖)𝜁,,(𝑥𝑖+1)=6[𝜁,,(x,𝑥𝑖+1)-ζ(x,𝑥i)]3-4ζ(𝑥0)=𝜁,,(𝑥𝑛)=0用追赶法解上面方程组,可求出𝜁,,(𝑥𝑖)和ζ(𝑥𝑖,𝑥𝑖+1),而𝜁,,(x)=𝜁,,(𝑥𝑖)+(x-𝑥𝑖)𝜁,,(𝑥i,𝑥𝑖+1)3-5这种做法有诸多好处,其中优点有三点:其一计算简便,其二保留了多项式的优点,其三克服了多项式的缺点。多项式的缺点是单个多项式会有不灵活不稳定的现象。由于三次样条曲线的种种优点,往往在实际中当遇到测线长已知点多的情况下采用此方法拟合高程。3.3曲面拟合法曲面拟合法是用于GPS点的分布在一定区域的时候,且可以选择数学曲面拟合该区域的似大地水准面,构造适当的数学模型,计算该区域内的高程异常值,然后求出正常高。这种拟合法的主体思想和曲线拟合法异曲同工的。具体思想是:已知测区的若干已知水准点,并且用GPS测定这些点的高程,利用公式求得这些点的高程异常,有了已知点的高程异常,已知点的平面坐标是已知的,所以利用其平面坐标(x,y)和高程异常值ζ构造出来的数学模型拟合最为接近于该测区的似大地水准面,然后内插出未知点的高程异常值ζ,进而求出正常高。3.3.1多项式曲面拟合法测站点的大地高H与正常高h之间有如下关系:h=H-ζ3-6多项式函数拟合法的基本思想是在小区域GPS网内,将似大地水准面看成曲面(或平面),将高程异常表示为平面坐标(x,y)的函数,通过网中起算点(既进行了GPS测量又进行了几何水准联测的点)已知的高程异常确定测区的似大地水准面形状,求出其余各点的高程异常,然后根据式(3-6)求出其他点的正常高,其数学模型为:ζ=f(x,y)+ε3-7式中f(x,y)是拟合的似大地水准面;ε是拟合误差,f(x,y)=𝑎0+𝑎1x+𝑎2y+𝑎3𝑥2+𝑎4xy+𝑎5𝑦2+…3-8x=B-0y=L-00=1𝑛∑0=1𝑛∑其中:n为GPS网中点的数量,(B,L)为已知点的大地坐标,𝑎0,𝑎1,𝑎2,𝑎3,𝑎4,𝑎5……为拟合待定参数;x,y为各GPS点的平面坐标的近似值,一般取起算点的平面坐标减去网中全部点平面坐标的均值。〔1〕二次曲面拟合取(3-8)式中的一、二次项后将大地水准面拟合为:f(x,y)=𝑎0+𝑎1x+𝑎2y+𝑎3𝑥2+𝑎4xy+𝑎5𝑦23-9即得二次曲面拟合模型:ζ=[𝑎0𝑎1𝑎2𝑎3𝑎4𝑎5][1𝑥𝑦𝑥2𝑥𝑦𝑦2]T+ε3-10每一个起算点可以组成一个上式,若共存在m个这样的起算点,则可列出m个方程:𝜁1=𝑎0+𝑎1𝑥1+𝑎2𝑦1+𝑎3𝑥12+𝑎4𝑥1𝑦1+𝑎5𝑦11+1𝜁2=𝑎0+𝑎1𝑥2+𝑎2𝑦2+𝑎3𝑥22+𝑎4𝑥2𝑦2+𝑎5𝑦21+2……𝜁𝑚=𝑎0+𝑎1𝑥𝑚+𝑎2𝑦𝑚+𝑎3𝑥𝑚2+𝑎4𝑥𝑚𝑦𝑚+𝑎5𝑦𝑚1+𝑚3-11从而组成误差方程:V=-Bx+L3-12上式中,B=[1𝑥1𝑦1𝑥12𝑦12𝑥1𝑦11𝑥2𝑦2𝑥22𝑦22𝑥2𝑦21𝑥𝑚𝑦𝑚𝑥𝑚2𝑦𝑚2𝑥𝑚𝑦𝑚],x=[𝑎0𝑎1𝑎2𝑎3𝑎4𝑎5],L=[𝜁1𝜁2𝜁𝑚],解得x=(𝑃)−1PL3-13解算出𝑎𝑖即可求出网中其余点的高程异常,并利用式(3-6)求出各未知点的正常高h。〔2〕多项式平面拟合在小范围或平原地区,可以认为大地水准面趋近于平面。此时,可选用公式(3-8)的前三项,将大地水准面拟合为:f(x,y)=a0+a1x+a2y3-14拟合模型为:ζ=[𝑎0𝑎1𝑎2][1𝑥𝑦]T+ε3-15其中,ai(i=0,1,2)为未知参数,此时要求公共点至少3个。〔3〕多项式相关平面拟合:也叫做四参数曲面拟合,若选用公式(3-8)的前三项和第五项进行拟合,则拟合曲面的表达式变为:f(x,y)=a0+a1x+a2y+a3xy3-16拟合模型为:ζ=[𝑎0𝑎1𝑎2𝑎5][1𝑥𝑦𝑥𝑦]T+ε3-17其中,ai(i=0,1,2,3)为未知参数,此时需要公共点至少4个。3.3.2移动曲面拟合法移动曲面拟合法是一种局部逼近法,其基本思想是以每一个内差点为中心,利用内差点周围数据点的值,建立一个拟合曲面,使其到各个数据点的距离之加权平方和为极小,而这个曲面在内插点上的值就是所求的内插值。设P为内插的点,下面对P构造相应的曲面。本文取如下的二次多项式曲面为例:f(x,y)=𝑎0+𝑎1x+𝑎2y+𝑎3𝑥2+𝑎4xy+𝑎5𝑦23-18设选取数据点的坐标为(𝑥𝑖,𝑦𝑖),i=1,2,…,n;n≥6且设内插点P的坐标为(𝑥𝑝,𝑦𝑝)。将(𝑥𝑖,𝑦𝑖)改化到以P为原点的局部坐标系中,即:𝑥𝑖,=𝑥𝑖-𝑥𝑝𝑦𝑖,=𝑦𝑖-𝑦𝑝3-19形成新的坐标(𝑥𝑖,,𝑦𝑖,),为移动坐标。任一点数据(𝑥𝑖,𝑦𝑖)假设距离d的递减函数为:ω(d)=ω√(𝑥𝑖−𝑥𝑝)2+(𝑦𝑖−𝑦𝑝)23-20将ω(d)作为权函数,对每个数据点赋予权𝜔𝑖,这里𝜔𝑖不是代表数据点的观测精度,而是反映该点与内插点的相关程度的大小。因此,权𝜔𝑖确定的原则应该与该数据点和内插点的距离𝑑𝑖有关,𝑑𝑖越小,它对内插点的影响越大,则权应越大。相反,𝑑𝑖越小,它对内插点的影响越小,则权应越小。最后,由最小二乘法解如下带权的极小值问题:Min丨∑𝜔𝑖𝑣𝑖2𝑛𝑖=1丨=R3-21为了给出二次多项式曲面,式(3-18)的系数,那么这时需要选取P点周围的数据点。当点数不够多时,则应扩大R的取值。现在这里由n个数据点的值,可得到如下的方程式:𝑣𝑖=𝑎0+𝑎1𝑥𝑖,+𝑎2𝑦𝑖,+𝑎3𝑥,𝑖2+𝑎4𝑥𝑖,𝑦𝑖,+𝑎5𝑦,𝑖2-𝑓𝑖(i=1,2,…,n)3-