劉瑛,張德明

(太原市勘察測繪研究院,山西太原 030002)

基于EGM2008和多面函數(shù)的太原市似大地水準(zhǔn)面精化

劉瑛?,張德明

(太原市勘察測繪研究院,山西太原 030002)

采用EGM2008重力場模型和趨勢面擬合與基于正交最小二乘思想篩選節(jié)點(diǎn)的多面函數(shù)法,運(yùn)用移去恢復(fù)原理建立了太原市區(qū)的似大地水準(zhǔn)面模型,擬合得到的模型內(nèi)符合與外符合中誤差均小于2 cm,完全可以滿足太原市區(qū)內(nèi)常規(guī)測量的要求。

EGM2008;趨勢面擬合;多面函數(shù);似大地水準(zhǔn)面

1　引 言

建立一個(gè)區(qū)域的似大地水準(zhǔn)面模型,比較標(biāo)準(zhǔn)的方法是移去恢復(fù)法,即將測區(qū)范圍內(nèi)離散的點(diǎn)進(jìn)行GPS水準(zhǔn)聯(lián)測后,把求出采樣點(diǎn)的高程異常分解成中長波分量、短波分量和殘余部分[1]。此時(shí)用模型大地水準(zhǔn)面來表示中長波部分,用足夠分辨率的數(shù)字地形模型計(jì)算短波部分,以上兩部分可認(rèn)為是趨勢項(xiàng),在計(jì)算時(shí)先移去此部分。再對殘余部分用數(shù)學(xué)或統(tǒng)計(jì)的方法,如多項(xiàng)式、多面函數(shù)、最小二乘配置等方式進(jìn)行擬合。最后再將擬合結(jié)果和趨勢向相疊加,恢復(fù)最初的高程異常進(jìn)而擬合出該區(qū)域的似大地水準(zhǔn)面模型。

通常趨勢項(xiàng)部分由重力場模型計(jì)算求出,在精化似大地水準(zhǔn)面方面,許多學(xué)者對各個(gè)重力模型進(jìn)行了對比,均顯示EGM2008模型精度要優(yōu)于其他模型,并肯定了其在中國大陸地區(qū)的適用性。因此在進(jìn)行似大地水準(zhǔn)面精化時(shí),EGM2008是計(jì)算模型大地水準(zhǔn)面的首選地球重力場模型[2]。

在對殘余部分進(jìn)行擬合時(shí),采用多面函數(shù)是比較普遍的做法。為了突出反映高程異常的變化,在進(jìn)行多面函數(shù)擬合之前,可先進(jìn)行趨勢面分析,再對剩余部分進(jìn)行曲面擬合;而在多面函數(shù)擬合中,節(jié)點(diǎn)的選擇對最終的結(jié)果影響較大,因此本文采用正交最小二乘的思想來篩選節(jié)點(diǎn)[3]。

2　基于EGM2008重力場模型和趨勢面正交多面函數(shù)的移去恢復(fù)法基本思想

高程異常ζ可分為長波分量、短波分量和殘余部分,即:

長波分量ζGM可由重力場模型求出[4],短波分量ζT通過數(shù)字高程模型計(jì)算[4],移去的過程即求出殘余部分的高程異常:

如果沒有數(shù)字高程模型,也可以將短波部歸入殘余高程異常:

由于重力場模型、GPS大地高、水準(zhǔn)高均存在測量誤差、重力模型高頻分量以及不同基準(zhǔn)面等的影響,模型大地水準(zhǔn)面和實(shí)測水準(zhǔn)面之間還會存在系統(tǒng)性偏差[5,6]。而從擬合的角度來說,殘余的高程異常仍可認(rèn)為由趨勢項(xiàng)和隨機(jī)信號兩部分組成,因此可先對趨勢面進(jìn)行分析,再對剩余部分進(jìn)行擬合操作。將高程異常的趨勢面表示為與坐標(biāo)的函數(shù):

式中f(x,y)常為多項(xiàng)式函數(shù),可由最小二乘法求出系數(shù)。判斷最終擬合的結(jié)果是否顯著,可進(jìn)行F檢驗(yàn)[7]。如果存在異常的觀測值,則需采用抗差估計(jì)來控制其對參數(shù)求解的影響[3,8]。

在經(jīng)過趨勢面分析之后,采用多面函數(shù)進(jìn)行最終的擬合,多面函數(shù)可將某一點(diǎn)的高程異常ζ(x,y)表示為:

式中,αi為待求參數(shù),n為參數(shù)個(gè)數(shù),q為核函數(shù)。

擬合的關(guān)鍵就是核函數(shù)的確定,這其中又以節(jié)點(diǎn)的確定最為關(guān)鍵,可從觀測點(diǎn)中選擇均勻分布的點(diǎn)作為結(jié)點(diǎn),但此種做法帶有一定的隨機(jī)性。為了能用最少的參數(shù)逼近到滿足要求的精度,通常需要對節(jié)點(diǎn)進(jìn)行篩選,可通過顯著性t檢驗(yàn)[9]、逐步回歸[10]以及正交最小二乘的方式[3,11]進(jìn)行選擇。

采用正交最小二乘的方式來篩選,若將所有觀測點(diǎn)作為節(jié)點(diǎn),并表示成列向量的形式,誤差方程的系數(shù)矩陣可表示為:

式中,βi=[q(x1,y1,x0i,y0i) q(x2,y2,x0i,y0i)… q(xn,yn,x0i,y0i)]T,將βi相互正交化,則觀測向量Z和βi的夾角θ為:

θ值越大,則對應(yīng)的βi,即參數(shù)αi對應(yīng)的節(jié)點(diǎn)對擬合的貢獻(xiàn)越大。因此,可以通過判斷正交化的系數(shù)矩陣與觀測值夾角大小的方式來篩選結(jié)點(diǎn),從而以較少的參數(shù)達(dá)到較高的擬合度,同時(shí)又能有足夠的多余觀測來保證推估的可靠性。

3　太原市區(qū)似大地水準(zhǔn)面精化案例分析

3．1概況

為了滿足太原市城市建設(shè)和規(guī)劃發(fā)展及現(xiàn)代化測繪生產(chǎn)的需要,本文利用了太原市區(qū)的控制測量成果和EGM2008重力場模型對市區(qū)的高程異常進(jìn)行了擬合?？刂泣c(diǎn)均為E級GPS點(diǎn),共194個(gè),較為均勻地分布在了太原市建成區(qū)及清徐縣的平坦地區(qū)。所有點(diǎn)位均聯(lián)測了二等水準(zhǔn),獲得了在WGS84坐標(biāo)系下的大地高以及在太原高程基準(zhǔn)大沽高程系下的正常高。求得了這些點(diǎn)的高程異常值采用MATLAB雙調(diào)和樣條插值內(nèi)插出的高程異常等值線如圖1。從194個(gè)點(diǎn)中,選擇146個(gè)點(diǎn)進(jìn)行擬合計(jì)算,其余48個(gè)點(diǎn)作為檢查點(diǎn)用于外部質(zhì)量評價(jià),點(diǎn)位分布如圖2所示。

圖1　高程異常等值線

圖2　擬合點(diǎn)與檢查點(diǎn)分布

3．2模型高程異常的移去

本文利用1′×1′的EGM2008模型文件計(jì)算出了太原市區(qū)的格網(wǎng)模型高程異常,如圖3所示。

圖3　EGM2008下的太原市區(qū)高程異常

同時(shí)根據(jù)控制點(diǎn)的經(jīng)緯度計(jì)算了各個(gè)點(diǎn)位的模型高程異常,并通過GPS水準(zhǔn)獲取的高程異常求出了殘余高程異常。從圖1和圖3可以看出,通過控制點(diǎn)的高程異常繪制的等值線和模型高程異常值均呈現(xiàn)由東北向西南減小的趨勢的各控制點(diǎn)殘余高程異常如圖4所示。3．3 趨勢面及多面函數(shù)擬合

圖4　殘余高程異常值

扣除模型高程異常后的殘余高程異常ζ仍有系統(tǒng)性偏差,本文采用四參數(shù)模型[5]進(jìn)行趨勢面擬合:

ζRES=a0+a1cosφcosλ+a2cosφsinλ+a3sinφ+v(8)

擬合過后再對殘余的高程異常進(jìn)行多面函數(shù)擬合。首先基于正交最小二乘篩選節(jié)點(diǎn),將所有的146個(gè)高程異常點(diǎn)作為節(jié)點(diǎn)插值的情況求出系數(shù)矩陣并進(jìn)行Gram-Schmidt正交化,再根據(jù)式(7)求出各個(gè)參數(shù)的對應(yīng)系數(shù)矩陣向量與高程異常向量的夾角,從中篩選出39個(gè)點(diǎn)作為節(jié)點(diǎn)。再通過簡單試算確定曲面類型和平滑因子,即可確定核函數(shù),并求出每個(gè)節(jié)點(diǎn)對應(yīng)的參數(shù)。測區(qū)范圍內(nèi)扣除模型高程異常的殘余部分經(jīng)過趨勢面和多面函數(shù)擬合之后的結(jié)果如圖5所示。

圖5　殘余高程異常擬合結(jié)果

3．4結(jié)果及精度評定

將各控制點(diǎn)的高程異常移去模型高程異常并擬合出殘余高程異常的參數(shù)后,太原市區(qū)的任意一點(diǎn)高程異常值即可通過簡單的恢復(fù)過程求出,擬合的最終結(jié)果如圖6所示。

圖6　最終擬合結(jié)果

參與計(jì)算的146個(gè)控制點(diǎn)擬合殘差以及48個(gè)檢查點(diǎn)的不符值如圖7、圖8所示。

圖7　擬合殘差

圖8　檢查點(diǎn)殘差

最終的統(tǒng)計(jì)結(jié)果如表1所示。

擬合結(jié)果統(tǒng)計(jì)　表1

從圖7、圖8、表1看出,采用EGM2008重力場模型以及趨勢面和多面函數(shù)進(jìn)行擬合的結(jié)果,除一個(gè)點(diǎn)的殘差為3．51 cm外,其余擬合殘差均小于3 cm,內(nèi)符合精度為1．62 cm;外符合精度的最大不符值為3．14 cm,其余結(jié)果均小于2．5 cm,外符合精度為1．13 cm,總體上有著較好的擬合精度。

4　結(jié) 論

本文利用移去恢復(fù)原理,將高程異常分成模型高程異常,殘余高程異常的趨勢項(xiàng)部分和殘余高程異常的信號部分,分別利用EGM2008重力場模型、趨勢面擬合和正交多面函數(shù)法建立了太原市區(qū)的似大地水準(zhǔn)面模型。最終結(jié)果的內(nèi)符合精度為1．6 cm,外符合精度為1．1 cm,該擬合模型可以滿足太原市區(qū)常規(guī)的測繪生產(chǎn)需求。

[1] 寧津生,羅志才,楊沾吉等．深圳市1Km高分辨率厘米級高精度大地水準(zhǔn)面的確定[J]．測繪學(xué)報(bào),2003,32(2): 102～107．

[2] 劉東順．基于EGM2008的似大地水準(zhǔn)面精化方法研究[D]．沈陽:東北大學(xué),2010．

[3] 張菊清,劉平芝．抗差趨勢面與正交多面函數(shù)結(jié)合擬合DEM數(shù)據(jù)[J]．測繪學(xué)報(bào),2008(4):526～530．

[4] GB/T 23709-2009．區(qū)域似大地水準(zhǔn)面精化基本技術(shù)規(guī)定[S]．

[5] 肖杰,張錦,鄧增兵等．礦區(qū)似大地水準(zhǔn)面精化方法研究[J]．測繪通報(bào),2015(2):14～18．

[6] 榮敏,周巍,陳春旺．重力場模型EGM2008和EGM96在中國地區(qū)的比較與評價(jià)[J]．大地測量與地球動力學(xué), 2009,29(6):123～125．

[7] 劉大杰,陶本藻．實(shí)用測量數(shù)據(jù)處理方法[M]．北京:測繪出版社,2000．

[8] 趙麗華．區(qū)域地殼運(yùn)動模型實(shí)現(xiàn)的理論與方法研究[D]．西安:長安大學(xué),2011．

[9] 陶本藻,余井泉．變形模型運(yùn)動參數(shù)的篩選和運(yùn)算[J]．測繪學(xué)報(bào),1988(4):241～246．

[10] 趙承坤,黃立人．速率面擬合法中核函數(shù)中心點(diǎn)的選擇[J]．地殼形變與地震,1991(2):48～54．

[11] S Chen,CN Cowan,PM Grant．Orthogonal Least Squares Learning Algorithms for Radial Basis Function Networks[J]．IEEE Trans Neural Networks,1991,2(2):302～309．

[12] 黃立人,陶本藻,趙承坤．多面函數(shù)擬合在地殼垂直運(yùn)動研究中的應(yīng)用[J]．測繪學(xué)報(bào),1993(1):25～32．

Taiyuan Quasi-Geoid Refinement Based on EGM2008 and Multi-Surface Function

Liu Ying,Zhang Deming

(Taiyuan Institute of Surveying and Mapping,Taiyuan 030002,China)

A quasi-geoid model in Taiyuan urban areas established by using remove-restore method based on EGM2008 gravity field model,trend surface fitting and multi-surface function which choosing nodes by using the orthogonal least squares．The result shows that root mean square error of internal and external quality evaluation are both less than 2 cm,completely meet the conventional surveying production demand in the urban area of Taiyuan．

EGM2008;trend surface fitting;multi-surface function;quasi-geoid model

1672-8262(2016)01-20-04

P223

B

?2015—10—09

劉瑛(1979—),男,工程師,主要從事工程測量和測量項(xiàng)目管理工作。