張艷兵,田林亞,畢繼鑫,王永敏

(河海大學(xué) 地球科學(xué)與工程學(xué)院,江蘇 南京 211100)

基于丹麥法計(jì)算權(quán)因子的GPS高程擬合方法研究

張艷兵,田林亞,畢繼鑫,王永敏

(河海大學(xué) 地球科學(xué)與工程學(xué)院,江蘇 南京 211100)

針對(duì)GPS高程擬合這一問題,介紹了等值線法、解析擬合法等傳統(tǒng)GPS高程擬合方法,并選用丹麥法計(jì)算權(quán)因子,構(gòu)造等價(jià)權(quán)矩陣,結(jié)合總體最小二乘法抵抗粗差點(diǎn)對(duì)GPS高程擬合結(jié)果的影響。結(jié)果表明:等值線法與其它兩種方法對(duì)不含粗差的等精度觀測數(shù)據(jù)擬合的效果相當(dāng);基于丹麥法計(jì)算權(quán)因子的抗差總體最小二乘法對(duì)不等精度觀測數(shù)據(jù)的偶然誤差處理效果更佳,且通過權(quán)因子可以有效處理粗差點(diǎn)對(duì)高程擬合的影響,擬合效果更佳。

GPS高程擬合;等值線法;解析擬合法;丹麥法;抗差總體最小二乘法

0 引 言

傳統(tǒng)高程測量是指利用水準(zhǔn)儀進(jìn)行測量,獲取所需高程,即正常高[1]。隨著現(xiàn)代化進(jìn)程的推進(jìn),社會(huì)建設(shè)的規(guī)模越來越大,對(duì)高程系統(tǒng)快速獲取要求更高,使得傳統(tǒng)測量手段在大范圍、地形起伏較大、通視性不良的區(qū)域難以實(shí)施。而GPS的高精度定位和快速獲取使其得到了廣泛關(guān)注,結(jié)合水準(zhǔn)測量體系,形成了GPS高程測量理論[2],降低了水準(zhǔn)測量在地形起伏較大區(qū)域的作業(yè)難度。

GPS直接獲取的高程是大地高,而在具體工程應(yīng)用中所使用的高程系統(tǒng)是正常高,兩者的參考面并不是同一基準(zhǔn)面[3],所以GPS高程要在獲取的大地高基礎(chǔ)上,減去兩個(gè)參考面的差異值。由于地球密度分布不均勻,使得似大地水準(zhǔn)面成為一個(gè)連續(xù)變化的不規(guī)則面,進(jìn)而導(dǎo)致高程異常值的不規(guī)則變化,難于準(zhǔn)確獲取[4-5]。如圖1所示,由于θ很小,所以控制點(diǎn)的高程異常值ξi可以根據(jù)ξi=Hi-hi直接獲取,然后采用某種方法建立控制點(diǎn)高程異常值與待求點(diǎn)高程異常值的相互關(guān)系,最后根據(jù)這種相互關(guān)系確定待求點(diǎn)的高程異常值估值。因此,如何精確測定控制點(diǎn)與待求點(diǎn)的相互關(guān)系是提高GPS高程測量精度的一個(gè)關(guān)鍵問題。

圖1 觀測點(diǎn)大地高與正常高的示意圖

1 等值線法

由于高程異常值隨著平面坐標(biāo)不斷變化,故可以將其轉(zhuǎn)化為三維模型,根據(jù)繪制等高線原理,繪制該區(qū)域高程異常值的等值線畫圖,然后根據(jù)待求點(diǎn)平面坐標(biāo)在等值線畫圖中拾取高程異常值估值。

具體過程為首先在測區(qū)選取適當(dāng)?shù)腉PS控制網(wǎng),但需要注意:1) GPS控制點(diǎn)需均勻分布;2) GPS控制點(diǎn)要布設(shè)在測區(qū)最高點(diǎn)、最低點(diǎn)、地性線上,使其可以充分反映測區(qū)范圍內(nèi)高程異常面的起伏變化;3) 在地勢起伏較大的區(qū)域,需加密GPS控制點(diǎn)的數(shù)量。其次選取其中m個(gè)控制點(diǎn)作為已知點(diǎn),并對(duì)其進(jìn)行水準(zhǔn)聯(lián)測,獲取已知點(diǎn)的高程異常值,其它控制點(diǎn)作為檢核點(diǎn),然后利用matlab畫出該區(qū)域高程異常值的等值線畫圖,并在圖上標(biāo)出檢核點(diǎn)的位置。最后利用matlab中點(diǎn)拾取工具,拾取檢核點(diǎn)的高程異常值信息,將其作為檢核點(diǎn)的高程異常值估值。

2 解析擬合法

解析擬合法,就是將已知點(diǎn)的高程異常值與檢核點(diǎn)的高程異常值的相關(guān)關(guān)系用數(shù)學(xué)知識(shí)精確表達(dá)為函數(shù)模型,然后根據(jù)函數(shù)模型計(jì)算檢核點(diǎn)的高程異常值估值。常見的解析擬合法有多項(xiàng)式擬合法和抗差總體最小二乘法。

2.1多項(xiàng)式擬合法

多項(xiàng)式擬合法的基本原理是首先利用GPS測量數(shù)據(jù)與水準(zhǔn)聯(lián)測高程計(jì)算出已知點(diǎn)的高程異常值ξi,然后建立高程異常值ξi與平面坐標(biāo)x,y的函數(shù)模型f(x,y),求得函數(shù)模型參數(shù),最后根據(jù)函數(shù)模型與待求點(diǎn)的平面坐標(biāo),計(jì)算待求點(diǎn)的高程異常值估值。由于二次多項(xiàng)式擬合模型[6-9]精度較高,且所需的已知點(diǎn)數(shù)量較少,在實(shí)際應(yīng)用中得到了廣泛使用,故本文將采用二次多項(xiàng)式模型擬合GPS高程異常值。

設(shè)已知點(diǎn)的高程異常ξ與其平面坐標(biāo)x,y之間關(guān)系為

ξ=f(x,y)+ε，

(1)

式中:ξ為高程異常值;f(x,y)為ξ的近似表達(dá)式;ε為誤差值。

二次多項(xiàng)式模型可表示為

f(x,y)=a0+a1x+a2y+a3xy+

a4x2+a5y2，

(2)

寫成矩陣形式:

ξ=BX+ε，

(3)

改寫為誤差方程:

V=BX-l.

(4)

要使擬合曲面最佳,只需滿足∑vivi=min,即:

VTPV=min .

(5)

在聯(lián)測的n個(gè)控制點(diǎn)中,選取m(m

2.2抗差總體最小二乘法

系數(shù)矩陣A和觀測向量ξ都是根據(jù)觀測數(shù)據(jù)獲得,當(dāng)觀測數(shù)據(jù)誤差較大甚至含有粗差時(shí),A和ξ也會(huì)受到很大影響,使平差結(jié)果出現(xiàn)較大偏差,擬合效果不佳。而抗差總體最小二乘法針對(duì)這一特性,引入權(quán)因子,優(yōu)化了平差模型,降低了觀測數(shù)據(jù)中粗差對(duì)擬合效果的影響。為便于比較抗差總體最小二乘法與多項(xiàng)式擬合法的區(qū)別,繼續(xù)采用二次多項(xiàng)式構(gòu)建總體最小二乘模型。

采用總體最小二乘模型[10],其誤差方差為

ξ+eξ=(A+EA)X,

(6)

式中:A∈Rn×m,ξ∈Rn,X∈Rm,EA為系數(shù)陣A的誤差,eξ為觀測向量ξ的觀測誤差。

若令eA=vec(EA),QA和Qξ分別為eA和eξ的協(xié)因數(shù)陣,則有:

(7)

要求式(6)的最佳估計(jì)值,需要滿足:

(8)

由于觀測值中含有粗差,所以要對(duì)總體最小二乘模型進(jìn)行改正,本文通過丹麥法計(jì)算權(quán)因子[11],抵制觀測值中粗差的影響。

(9)

計(jì)算步驟:

1) 計(jì)算初始值;

v0=0,x0=(ATPξ0A)-1ATPξ0ξ.

(10)

2) 迭代過程,計(jì)算中間值;

Qξi=QξI-1ωi,QAi=QAi-1ωi,

(11)

(12)

(13)

3) 判斷迭代終止條件為

‖Xi+1-Xi‖<ε.

(14)

根據(jù)式(14)判斷迭代過程是否結(jié)束,若不成立,重復(fù)迭代;若成立,迭代結(jié)束,計(jì)算單位權(quán)中誤差,其估值為

(15)

式中:r為自由度;λ=(Qξ+(XTQ0X)QA)-1(ξ-AX)。

3 算例分析

為了驗(yàn)證抗差總體最小二乘法在GPS高程擬合中的有效性,選取云南某地采集數(shù)據(jù)與某城市控制網(wǎng)數(shù)據(jù)進(jìn)行GPS高程異常值擬合計(jì)算。如圖2所示,數(shù)據(jù)1:云南某地區(qū)共布設(shè)了27個(gè)GPS控制點(diǎn),GPS網(wǎng)為C級(jí)網(wǎng),控制點(diǎn)全部為三等水準(zhǔn)點(diǎn),剔除3個(gè)帶有明顯粗差的控制點(diǎn),選取其中的2、6、12、14號(hào)點(diǎn)作為檢核點(diǎn),其它點(diǎn)作為已知點(diǎn)。數(shù)據(jù)2:某中型城市的城市控制網(wǎng),共已知20個(gè)GPS-E級(jí)控制點(diǎn),為研究GPS擬合原理,并對(duì)以上控制點(diǎn)進(jìn)行了3等水準(zhǔn)測量。選取其中的2、15、20號(hào)點(diǎn)作為檢核點(diǎn),其它點(diǎn)作為已知點(diǎn)。

圖2 兩組實(shí)驗(yàn)數(shù)據(jù)點(diǎn)位空間分布圖

3.1等值線法

由兩組實(shí)驗(yàn)數(shù)據(jù)分別生成間距為10 cm和20 cm的等值線畫圖,如圖3所示,然后根據(jù)檢核點(diǎn)平面坐標(biāo)在圖3上拾取其高程異常值估值,結(jié)果如表1、表2所示。

圖3 兩組實(shí)驗(yàn)數(shù)據(jù)生成的等值線畫圖

點(diǎn)號(hào)X坐標(biāo)/mY坐標(biāo)/m高程異常值/m高程異常值估值/m偏差/m22763781.594573220.44032.61632.5970.01962768522.200593641.08432.03232.064-0.032102766456.717613160.71432.04832.135-0.103122776204.448553188.40332.47132.2430.228

表2 數(shù)據(jù)2中檢核點(diǎn)高程異常值計(jì)算結(jié)果

從表1與表2中可以看出,由等值線法內(nèi)插出的檢核點(diǎn)高程異常值與真值(這里把聯(lián)測的水準(zhǔn)高程與大地高高差作為檢核點(diǎn)高程異常值的真值)之差最小可以達(dá)到cm級(jí)精度,最大可以達(dá)到厘米級(jí),所以用等值線法來內(nèi)插高程異常值的結(jié)果并不穩(wěn)定,且要求內(nèi)插點(diǎn)盡可能分布在中心位置,同時(shí),該方法也缺少精度評(píng)價(jià)指標(biāo),不易判斷生成的等值線正確與否,也不能判斷內(nèi)插點(diǎn)的高程異常值正確與否。

3.2解析擬合法

由于等值線法缺少精度衡量指標(biāo),很難判斷檢核點(diǎn)內(nèi)插出的高程異常值是否正確,因此,引入多項(xiàng)式擬合法與抗差總體最小二乘法來對(duì)兩組數(shù)據(jù)進(jìn)行擬合,并在數(shù)據(jù)1中添加粗差點(diǎn)(2 762 850.290 8,566 913.499 9,32.324 8)與(2 777 857.3 74 6,594 527.091 8,32.276 4),分析多項(xiàng)式擬合法與抗差總體最小二乘法的擬合效果。

利用多項(xiàng)式擬合法與抗差總體最小二乘法分別對(duì)數(shù)據(jù)1與數(shù)據(jù)2中的已知點(diǎn)進(jìn)行擬合處理,得到各點(diǎn)高程異常值ξi的殘差值vi,如圖4所示。其中,圖4(a)是利用多項(xiàng)式擬合法對(duì)數(shù)據(jù)1的處理結(jié)果,從中可以看出粗差點(diǎn)14,26不能得到很好的處理,對(duì)其它已知點(diǎn)的擬合結(jié)果造成了一定程度上的影響;圖4(b)是利用抗差總體最小二乘法對(duì)數(shù)據(jù)1的處理結(jié)果,從中可以粗差點(diǎn)對(duì)其它已知點(diǎn)的擬合結(jié)果影響較小;圖4(c)、4(d)示出了利用兩種方法對(duì)數(shù)據(jù)2的處理結(jié)果,可以看出兩種方法擬合效果相差不大。

圖4 兩組實(shí)驗(yàn)數(shù)據(jù)高程異常值擬合殘差 (a)多項(xiàng)式擬合法擬合殘差(數(shù)據(jù)1);(b) RWTLS法擬合殘差(數(shù)據(jù)1);(c)多項(xiàng)式擬合法擬合殘差(數(shù)據(jù)2);(d) RWTLS法擬合殘差(數(shù)據(jù)2)

表3 兩組實(shí)驗(yàn)數(shù)據(jù)中檢核點(diǎn)在不同擬合方法中高程異常值差值計(jì)算結(jié)果

根據(jù)已知點(diǎn)的殘差值與檢核點(diǎn)的高程異常值差值進(jìn)行精度分析,其中,內(nèi)符合精度計(jì)算公式為

(16)

式中,v為已知點(diǎn)殘差值向量,r為自由度。

外符合精度計(jì)算公式為

(17)

式中： Δ為檢核點(diǎn)的高程異常值真值與擬合值之差；n為檢核點(diǎn)的個(gè)數(shù)。

表4示出了兩組實(shí)驗(yàn)數(shù)據(jù)在不同擬合方法下的精度。

表4 兩組實(shí)驗(yàn)數(shù)據(jù)在不同擬合方法下的計(jì)算精度

4 結(jié)束語

本文采用等值線法、多項(xiàng)式擬合法、抗差總體最小二乘法分別對(duì)數(shù)據(jù)1與數(shù)據(jù)2進(jìn)行GPS高程擬合。當(dāng)觀測數(shù)據(jù)不含粗差,且觀測數(shù)據(jù)質(zhì)量較高時(shí),可以看出多項(xiàng)式擬合法與總體最小二乘擬合法都能取得很好的擬合效果,同時(shí),通過檢核點(diǎn)在三種方法下擬合結(jié)果的比較,可以得到等值線法內(nèi)插出的未知點(diǎn)高程異常值精度也很高,且該方法計(jì)算方便,便于提取區(qū)域中的任意點(diǎn)的高程異常值,此時(shí)等值線法更加適用;當(dāng)觀測數(shù)據(jù)中含有少量粗差時(shí),前兩種方法的擬合精度有所下降,且越靠近粗差點(diǎn),擬合效果越差,而抗差總體最小二乘法擬合效果依然可以保持穩(wěn)定,綜上,抗差總體最小二乘法適用范圍更廣,擬合效果也比較穩(wěn)定。

擬合計(jì)算中,本文只采用丹麥法來計(jì)算權(quán)因子,且取c=σ0,并未就c的取值對(duì)抗差效果進(jìn)行分析;也未引入其他模型與丹麥法計(jì)算權(quán)因子進(jìn)行比較,分析其對(duì)抗差總體最小二乘法擬合效果的影響,有待進(jìn)一步的研究。

[1] 孔祥元. 大地測量學(xué)基礎(chǔ)[M]. 武漢:武漢大學(xué)出版社, 2010.

[2] 陳為民,張旭東,符華年,等. GPS高程測量代替等級(jí)水準(zhǔn)測量的應(yīng)用研究[J]. 武漢大學(xué)學(xué)報(bào)(信息科學(xué)版),2013(7):828-831.

[3] 王珍,程壘. GPS高程擬合方法的研究和應(yīng)用[J]. 測繪通報(bào), 2015(s1):102-105.

[4] 田曉,鄭洪艷,許明元,等. 一種改進(jìn)的適用于不同地形的GPS高程擬合模型[J]. 測繪通報(bào), 2017(1):35-38,64.

[5] 吳迪軍,熊偉. 跨海橋梁GPS高程擬合方法[J]. 測繪科學(xué),2017(6):1-6.

[6] 張潘,余代俊,張玉剛,等. GPS高程擬合方法研究及精度對(duì)比試驗(yàn)[J]. 測繪通報(bào), 2015(9):54-56.

[7] 李宏博,史先琦,萬奇靈,等. GPS高程擬合在丘陵地區(qū)地形測量中的應(yīng)用研究[J]. 測繪地理信息, 2015, 40(4):44-47.

[8] 劉斌,郭際明,史俊波,等. 利用EGM2008模型與地形改正進(jìn)行GPS高程擬合[J]. 武漢大學(xué)學(xué)報(bào)(信息科學(xué)版), 2016, 41(4):554-558.

[9] 王昶,王旭. GPS高程擬合分區(qū)方法研究[J]. 大地測量與地球動(dòng)力學(xué), 2016, 36(1):26-29.

[10] 劉亞彬,鄭南山,張旭,等. GPS高程擬合的加權(quán)總體最小二乘抗差估計(jì)[J]. 大地測量與地球動(dòng)力學(xué), 2016, 36(1):30-34.

[11] 林國標(biāo),劉立龍,蔡成輝,等. 一種基于丹麥法的改進(jìn)型雙步M估計(jì)[J]. 大地測量與地球動(dòng)力學(xué),2015(2):235-238，247.

ResearchonGPSHeightFittingMethodBasedonDanishMethod

ZHANGYanbing,TIANLinya,BIJixin,WANGYongmin

(CollegeofEarthScienceandEngineering,HohaiUniversity,Nanjing211100，hina)

In this paper, the traditional GPS height fitting methods, such as contour method and analytic fitting method, are introduced for GPS elevation fitting. The Danish method is used to calculate the weight factor for the first time to construct the equivalent weight matrix, and combined with the overall least squares method to resist the influence of the rough point on the GPS elevation fitting results. The results show that the effect of the contour method is equivalent to the other two methods for fitting the equal precision observation data without gross error. Based on the weight factor calculated by Danish method, The deal result of the robust overall least squares method is better in dealing with the accidental error of the unequal accuracy observation data,and the effect of the rough points on the elevation fitting can be effectively dealt by using the weighting factor, and the fitting effect is better.

GPS elevation fitting; contour method; analytic fitting method; Danish method; robust overall least squares method

10.13442/j.gnss.1008-9268.2017.04.017

P228.4

A

1008-9268(2017)04-0096-06

2017-04-06

聯(lián)系人: 張艷兵E-mail:1821708477@qq.com

張艷兵(1992-),男,碩士研究生，研究方向?yàn)闇y量數(shù)據(jù)處理理論與方法。

田林亞(1963-),男,教授,博士,研究方向?yàn)榫芄こ虦y量理論與技術(shù)。

畢繼鑫(1994-),男,碩士研究生,研究方向?yàn)闇y量數(shù)據(jù)處理理論與方法。