肖 杰,張 錦
(1.中國科學院測量與地球物理研究所,湖北武漢430077;2.中國科學院動力大地測量學重點實驗室,湖北武漢430077;3.中國科學院大學,北京100049;4.太原理工大學測繪科學與技術系,山西太原030024)
建立礦區似大地水準面模型的主要目的是要將GPS觀測所獲取的大地高轉化為我國工程常用的正常高,從而替代勞動強度大、效率低的水準測量。將某一點的大地高轉化為正常高,其關鍵是要知道該點的高程異常值。
高程異常的確定方法有重力法和幾何法。重力法需要精度和分布較好的重力數據及地形數據,但由于我國重力數據較為稀缺及重力施測的難度和精度問題,限制了此方法的使用。幾何法則是通過模型擬合高程異常曲面,然后用內插的方法獲得待定點的高程異常值。由于所選模型的不同,產生了不同的擬合方法,基本可分為函數模型和統計模型。函數模型有多項式擬合[1-2]、多面函數法[3]及多種模型組合擬合[4]等;統計模型有加權平均法[5]、Kriging 方法[6]等,其實質都是對局部區域內似大地水準面的逼近。本文利用某礦區GPS/水準數據,通過對多種擬合模型擬合精度的比較分析,建立適合礦區的似大地水準面擬合模型,方便礦區測繪生產。
加權平均模型的基本原理是:設有n個控制點,則高程異常ζ的計算公式為

式中,pi是已知點的權函數。
Shepard于1964年提出了一種局部逼近模型,選定R>0,并規定權函數為[7-8]


此方法是一種改進的加權平均法,需要合理選擇R,使恰當數量的點(xi,yi)落入圓域中。權函數的光滑性和衰減性越好,曲面的擬合效果就越好。
多項式擬合模型在水準面擬合中是應用最多也是最普遍的,三次多項式擬合模型為[1]

式中,α0,α1,…,α9為多項式擬合系數;ζ為高程異常值,即GPS所測大地高與水準高程之差,單位為m;xm、ym為擬合區的中心平面坐標,單位為m。
多面函數法[9-10]是Hardy于1977年提出的,該方法認為任何一個圓滑的數學表面總可用一系列有規則的數學表面的總和以任意精度逼近。
多面函數方程的一般形式為

式中,ζ可看做水準面擬合的高程異常值;ai為待定系數;q(x,y,xi,yi)為核函數,通常為 x、y 的二次函數,其表達式為

式中,(x,y)為待求點坐標;(xi,yi)為中心點坐標,即核函數結點坐標;δ2為光滑因子,用來對核函數進行調整;β一般可選某個非零實數,如 0.5、1、-0.5等,選0.5 時為正雙曲面,選 -0.5 時為倒雙曲面。
多面函數法中的核函數和光滑因子δ2對最終的擬合效果有非常重要的影響,文獻[10]對光滑因子的選取作了詳細介紹,認為相同的光滑因子正雙曲面擬合精度高于倒雙曲面擬合精度。除此之外,核函數結點的選擇也是非常重要的,一般的原則是選擇分布均勻:區域內具有代表性的特征點作為結點,不同的結點坐標和結點密度將直接影響最終的擬合效果。
當確定了以上各個參數,即可列出誤差方程,用最小二乘法進行待定系數求解ai,其矩陣表達式用A表示,矩陣A的求解表達式可表示為

當計算出系數矩陣A后,即可通過式(3)求取各個點位的高程異常值。
圖1為某礦區控制點點位分布圖,共191個控制點,均進行了不少于2 h的GPS靜態觀測和三等或四等水準觀測。整個礦區控制網覆蓋范圍約1000 km2,礦區地形起伏較大,最大高差約517m。
為保證有足夠精度的大地高,采用CGCS2000坐標系下約束平差后的大地高(用4個分布均勻的國家高等級控制點的CGCS2000坐標進行約束平差)代替WGS-84坐標系下的大地高。根據控制點的大地高和水準高程可以計算出各點的高程異常值。下面分別用上述4種擬合模型計算高程異常值。其中多面函數法經反復試算,選用62個分布均勻的控制點作為結點坐標。加權平均模型、Shepard插值模型、三次多項式模型亦選用此62個控制點作為已知點參與計算,其余點作為外部檢核點。圖2為參與解算的62個已知點點位分布圖。

圖1 礦區控制點點位分布圖

圖2 62個控制點點位分布圖
圖3為三次多項式模型內符合精度差值圖,圖4為加權平均模型、Shepard插值模型和三次多項式模型外符合精度差值圖。

圖3 三次多項式內符合精度差值
圖5為多面函數法內符合精度插值圖。多面函數法的擬合計算選某礦區的182個GPS/水準點參與擬合參數計算。圖2為多面函數結點點位分布圖,其中,光滑因子δ取0.01,β取0.5。經計算后,內符合檢核點高程異常值計算結果與實測高程異常值(GPS控制點的大地高與水準觀測的正常高之差,以下同)計算結果之差如圖5所示,內符合檢核點中差值的絕對值大于4 cm的點數有35個,差值最大的數值高達0.3m。

圖4 加權平均模型、Shepard插值模型和三次多項式模型外符合精度差值

圖5 多面函數法高程異常值計算結果與實測值之差
從以上4種擬合方法的計算結果來看,加權平均模型插值效果最差,Shepard模型較加權平均模型插值效果有極大提高。多面函數法是GPS/水準擬合中應用較多的方法,其結點坐標,δ、β值的選擇對計算結果影響較大,需反復試算。4種方法的均方根誤差(RMSE)見表1。

表1 4種方法擬合均方根誤差 m
從表1中可以看出,三次多項式擬合方法的均方根誤差最小。
通過上述實例數據可以看出,在類似于某礦區這種大范圍的、地形較復雜的地區單純用一種方法來擬合整個礦區的似大地水準面,其擬合精度不夠理想。根據似大地水準面的特性,可以把高程異常曲面分為中長波項和短波項,用不同的擬合方法分別進行擬合,然后組合在一起,這樣可以發揮各種方法的優勢,彌補其劣勢,達到提高擬合精度的目的。
通常曲面的中長波項(即趨勢項)是比較平滑的曲面,用多項式擬合即可很好地去除,剩余的部分即為短波項影響,一般由于地形的高低起伏等因素引起,用Shepard插值法或多面函數法即可很好地對剩余高程異常進行擬合。圖6為用三次多項式去除趨勢項,用Shepard插值模型擬合的結果與實測高程異常之差值比較,均方根誤差±0.018m。

圖6 組合法擬合結果與實測值之差
文中首先分別用加權平均模型、Shepard插值模型、多項式擬合模型及多面函數法4種方法對某礦區控制點的GPS/水準數據進行似大地水準面擬合,對4種方法進行對比分析和精度評定,在這種單獨擬合的情況下,三次多項式擬合精度較其他3種方法較高。根據局部似大地水準面的物理和幾何特性,用組合擬合法進一步改善了擬合精度,擬合均方根誤差不超過±2 cm,表明此模型可以在地形起伏較大的礦區加以應用,結合GPS觀測成果可以取代四等及以下幾何水準測量。
[1]中華人民共和國國家質量監督檢驗檢疫總局,中國國家標準化管理委員會.GB/T 23709—2009區域似大地水準面精化基本技術規定[S].北京:中國標準出版社,2009.
[2]劉長建,柴洪洲,吳洪舉,等.GPS水準多項式擬合自動優選算法[J].測繪科學技術學報,2009(1):49-51.
[3]陶本藻,姚宜斌,趙美超.論多面函數推估與協方差推估[J].測繪通報,2002(9):4-6.
[4]鐘波,羅志才.GPS水準綜合模型的應用研究[J].測繪通報,2007(6):5-7.
[5]張玲彬,武雁剛,張坤軍.GPS水準加權綜合擬合模型的應用研究[J].黃河水利職業技術學院學報,2008(3):37-39.
[6]李明,高星偉,文漢江,等.Kriging方法在GPS水準擬合中的應用[J].測繪科學,2009(1):106-107.
[7]黃友謙.曲線曲面的數值表示和逼近[M].上海:上海科學技術出版社,1984.
[8]史利民,王仁宏.幾種基于散亂數據擬合的局部插值方法[J].數學研究與評論,2006(2):283-291.
[9]劉經南,施闖,姚宜斌,等.多面函數擬合法及其在建立中國地殼平面運動速度場模型中的應用研究[J].武漢大學學報:信息科學版,2001,26(6):500-503.
[10]馬洪濱,董仲宇.多面函數GPS水準高程擬合中光滑因子求定方法[J].東北大學學報:自然科學版,2008,29(8):1176-1178.