基于MATLAB似大地水準面精化的區域應用

封華梅

(陜西省一八五煤田地質有限公司，陜西 榆林 719000)

0 引言

似大地水準面指的是從地面點沿正常重力線量取正常高所得端點構成的封閉曲面，與大地水準面的差值為正常高與正高之差[1]。大地高與正高或正常高之間的差距，分別稱為大地水準面差距和高程異常[2]。求解大地水準面差距或高程異常的主要方法有幾何法[3]、重力法[4]及組合法等[5]。利用GNSS/水準擬合法是幾何確定似大地水準面的一種方法，其中多項式擬合是建立擬合模型的一種。其原理就是在區域范圍內的GNSS網中，存在若干個聯測了水準高程的GNSS點(稱GNSS水準重合點)，根據GNSS水準重合點上已知的大地高和正常高確定GNSS水準重合點上的高程異常，然后再由GNSS水準重合點的平面坐標與高程異常的函數關系構造某一曲面來逼近似大地水準面，進而推算出測區中未進行水準聯測的GNSS點的高程異常，從而得到測區內所需GNSS點的正常高[6]。

1 函數模型理論

函數模型擬合的最大優點是對于趨勢性變化的擬合效果較好，但需要合適的函數模型形式，并需要確定合適的參數個數。若函數模型的針對性強，則擬合與推估的效果好。函數模型的形式一經確定，計算過程中一般不再變化。基于選定的函數模型及其相應的擬合原則，很容易求得模型中的待估參數，如此即完成了函數模型的構造與擬合。基于擬合的函數模型，即可求解未知點的信號。然而，任何待定的模型均難以“準確”地表征逼近場的一切特征。于是難免有部分“剩余誤差”；另一方面，任何函數模型逼近均基于足夠的已知點信息，這類已知信息的任何誤差都將對函數模型的待估參數的估計帶來影響。這種影響又必將對后續未知點值的推估帶來誤差[7]。

1.1 多項式擬合模型

多項式函數的數學模型是

ξ(x，y)=F(x，y，β)

(1)

式中:ξ(x，y)—點(x，y)處的高程異常值；F(x，y，β)—所選擇的函數形式。

將(1)式表示成

ξ(x，y)=a0+a1x+a2y+a3x2+a4y2+a5xy

(2)

式中:a0～a5—待定系數。

通常采用較多的模型是二次多項式模型，以此推導平差過程：

利用公共點上GNSS測定的大地高和水準測量的正常高計算出該點上的高程異常，若測區內有n(n>6)個這樣的公共點，則可算出高程異常ξi(i=1，2，…，n)，根據模型可以列出n個誤差方程

(3)

寫成矩陣形式為

V=BX-L

(4)

其中：

V=[v1v2…vn]T

X=[a0a1a2a3a4a5]T

L=[ξ1ξ2…ξn]T

根據最小二乘法原理VTPV=min可得

X=(BTB)-1BTL

(5)

就可解算出擬合方程的待定系數，然后根據模型解算出未知點上的高程異常值。

1.2 GNSS/水準擬合的精度分析

目前GNSS/水準擬合的精度主要通過內符合精度和外符合精度來評定[8]。

內符合精度：根據參與擬合的已知點的實測高程異常ξi與擬合高程異常ξ＇i的差值vi進行計算，若參與擬合的已知點數為n，則內符合精度σ內為

(6)

外符合精度:根據檢核點的實測高程異常ξi與擬合高程異常ξ＇i的差值vi進行計算，若檢核點的數為m，則外符合精度σ外為

(7)

2 基于MATLAB程序的實例分析計算

為了分析GNSS/水準擬合法在實際生產中應用的可靠性，在此將利用小保當煤礦工業廣場和風井廣場的實測控制點進行驗證分析，采用以下評定精度的標準：最小擬合誤差Vmin、最大擬合誤差Vmax、內符合精度σ內、外符合精度σ外。

2.1 測區概況介紹

該區域位于小保當煤礦工業廣場和風井廣場之間，形狀為長方形，東西長約3.6 km，南北寬約1.4 km，面積約5 km2。本區域西高東低，局部均相對平坦。在該區域均勻分布GNSS/水準重合點21個，在相對獨立坐標系的成果見表1。

2.2 參與計算的已知點選取

根據選點原則，在測區均勻選擇部分GNSS/水準重合點進行參與計算，通過分析選擇XG02、XG03、XG04、XG05、XG06、XG07、XG08、XG10、XG11、XG12、XGJ8、YF04、YF06、EF03等14個點參與計算。

表1 GNSS/水準重合點成果表

2.3 基于MATLAB程序計算求取擬合系數

通過Matlab程序對數據進行計算，部分代碼顯示如下

format long g；

m=length(P(：，1))；

B=zeros(m，6)；

for i=1：m；

B(i，：)=[1，P(i，1)，P(i，2)，P(i，1)×P(i，1)，P(i，2)×P(i，2)，P(i，1)×P(i，2)]；

end

ss=P(：，3)；

cs=pinv(B′×B)×B′×ss。

將數據代入計算出二次多項式擬合模型系數為

a0=-0.000 000 233 687 411 332 634

a1= 0.000 104 563 005 212 421

a2=-0.000 261 039 683 129 869

a3= 0.000 000 004 476 835 742 444 78

a4=0.000 000 054 012 643 907 004 7

a5=-0.000 000 029 420 139 308 138 1

2.4 通過計算的模型系數進行精度評定

通過Matlab程序對數據計算精度，部分代碼顯示如下

m=length(P(：，1))；

fori=1：m；

v(i，：)=[cs(1)+cs(2)×P(i，1)+cs(3)×P(i，2)+cs(4)×P(i，1)×P(i，1)+cs(5)×P(i，2)×P(i，2)+cs(6)×P(i，1)×P(i，2)]；

end

v

利用Matlab程序進行計算，并分別給出計算結果。參與計算的點殘差值見表2，未參與計算的點殘差值見表3。

表2 參與計算點殘差值表

表3 未參與計算點殘差值表

根據上述(6)式計算，內符合精度為σ內=0.030 m；上述(1～10)式計算，外符合精度為σ外=0.058 m。

3 結論

(1)各種GNSS/水準擬合法的精度都達到了厘米級的精度，選用合適的模型可以代替三、四等水準測量。

(2)在平緩地區GNSS/水準擬合法的精度大于起伏地區，所以地形對GNSS/水準擬合法精度的影響較大。

(3)選擇不同擬合點數對精度的影響不是很大，通過比較并不是點數選擇的越多其擬合精度就越好。