顧及EGM2008和殘差地形模型的GPS高程轉(zhuǎn)換方法研究

谷延超，范東明

（西南交通大學(xué) 測量工程系，四川 成都 610031）

為了克服水準(zhǔn)測量在山區(qū)和丘陵地區(qū)實(shí)施困難的缺點(diǎn)，GPS大地高轉(zhuǎn)換為正常高以代替水準(zhǔn)測量則是山區(qū)水準(zhǔn)測量的有效方法之一［1－4］。GPS高程轉(zhuǎn)換即將以參考橢球面為基準(zhǔn)面的GPS大地高通過某種特定的函數(shù)轉(zhuǎn)換為以似大地水準(zhǔn)面為基本面的正常高，即精確確定高程異常。確定高程異常大致可分為數(shù)學(xué)擬合法和重力法兩大類［2－4］。數(shù)學(xué)擬合法，其本質(zhì)為利用似大地水準(zhǔn)面的連續(xù)光滑性，使用一定的數(shù)學(xué)曲面來擬合高程異常，但是數(shù)學(xué)擬合并不能夠真實(shí)地反映似大地水準(zhǔn)面的變化。由于在重山區(qū)似大地水準(zhǔn)面變化不平穩(wěn)，數(shù)學(xué)擬合法使用效果并不理想［2］。重力法是指利用局部或者全球重力數(shù)據(jù)求解高程異常，將GPS大地高轉(zhuǎn)換為正常高。重力法充分利用了地球物理場，但由于重力數(shù)據(jù)缺乏和重力場模型精度對于工程應(yīng)用而言較低，所以工程建設(shè)一般都不采用重力法直接進(jìn)行GPS高程轉(zhuǎn)換。

由于地球重力場模型能以一定的分辨率精確表達(dá)地球重力場中長波信號，但小于該分辨率的重力場短波信號則遺失，使得利用地球重力場求解高山區(qū)高程異常存在信號遺失［5］。而地形數(shù)據(jù)則能夠較細(xì)致地表達(dá)地球重力場短波信號，但其對地球重力場中長波信息表達(dá)欠佳。故可結(jié)合地球重力場和地形數(shù)據(jù)精確求解高程異常。

本文首先利用SRTM數(shù)據(jù)（Shuttle Radar To－pography Mission）求得殘差地形模型（RTM）數(shù)據(jù)，獲取地球重力場短波信號，然后利用重力場模型求解地球重力場中長波信號，恢復(fù)地球重力場全波段信號，最后利用合理的數(shù)學(xué)擬合法消除系統(tǒng)誤差。實(shí)驗(yàn)證明該方法能較好地恢復(fù)地球重力場全波段信號，在地形起伏較大的山區(qū)能提高GPS高程轉(zhuǎn)換精度，并且不受數(shù)據(jù)限制，能應(yīng)用于各類工程，體現(xiàn)該方法的有效性、實(shí)用性和優(yōu)越性。

1 EGM2008高程異常

EGM2008模型是由美國國家地理空間情報局地球重力場研發(fā)小組釋放的全球超高階地球重力場模型，EGM2008模型達(dá)到2 159階次（球諧系數(shù)可擴(kuò)展至2 190，階次為2 159），因此，它能以約5′的空間分辨率計(jì)算擾動位T和重力相關(guān)量如重力異常Δg和高程異常ζ等［5－6］。

地球重力場球諧函數(shù)模型例如EGM2008由一系列完全正?；那蛑C函數(shù)系數(shù)組成，根據(jù)球諧函數(shù)連續(xù)疊加可計(jì)算擾動位T，利用Bruns公式可求得地面上任意一點(diǎn)的高程異常［3］

式中：（r，θ，λ）為以地球質(zhì)心為坐標(biāo)原點(diǎn)、Z軸與地球自轉(zhuǎn)軸重合的坐標(biāo)系中的球坐標(biāo)（其中r為向徑，θ為余緯，λ為經(jīng)度），γ為該點(diǎn)的正常重力值，Pnm（cosθ）為完全正?；兆尩戮喗Y(jié)函數(shù)，Cnm，Snm為完全正?；瘮_動位系數(shù)，m為實(shí)際地球引力位減去正常引力位的系數(shù)，a＝6.378 136 46×106m為參考橢球長半軸，nmax為系列地球重力場模型展開的球諧函數(shù)最高階數(shù)，fM＝3.986 004 415×1014m3s－2為地心引力常數(shù)。

2 RTM高程異常

由于球諧函數(shù)展開截止至2 159階，EGM2008不能夠體現(xiàn)地球重力場中波長小于5′的高頻信息。這將在求解擾動位時產(chǎn)生信號遺失，同時高程異常求解時也將產(chǎn)生相似的信號遺失誤差［5］。

SRTM是美國航天飛機(jī)飛行中利用雷達(dá)測圖技術(shù)得到的數(shù)字地表模型（DSM），可獲取分辨率為3″（約90m）的數(shù)字高程數(shù)據(jù)HSRTM［5］。SRTM高程數(shù)據(jù)的精度已經(jīng)由多位研究者加以驗(yàn)證，通過與各國的國家高程坐標(biāo)數(shù)據(jù)庫進(jìn)行對比，精度在平坦地區(qū)可以達(dá)到4～6m，在山地差距將達(dá)到11～14m［5，7］。

地球重力場模型結(jié)合數(shù)字地形模型能夠更有效地利用地形和重力場的短波信號的強(qiáng)相關(guān)性，故殘差地形模型（RTM）數(shù)據(jù)可以用來重建地球重力場短波信號。

2.1 SRTM構(gòu)建RTM數(shù)據(jù)

RTM數(shù)據(jù)由DTM減去一個平滑的參照網(wǎng)格。如果參考網(wǎng)格選擇與球諧函數(shù)的最大階次保持一致，RTM數(shù)據(jù)將增加球諧函數(shù)最大階次以上信號。該方法可使球諧函數(shù)高程異常的波段信號增加，同時減小信號遺失誤差。本文使用SRTM通過高通濾波求解平滑參考面，利用SRTM數(shù)據(jù)和平滑參考面構(gòu)建RTM數(shù)據(jù)［5］。高通濾波求解平滑參考面。

式中：H為DTM數(shù)據(jù)矩陣，（k，l）為平滑參考面格網(wǎng)中心對應(yīng)的DTM矩陣的坐標(biāo)。參數(shù)n可通過濾波窗口的大小ω和DTM分辨率Δx計(jì)算。

為了確保建立平滑參考面和EGM2008的球諧函數(shù)展開的最高階次相同，窗口寬度ω應(yīng)與球諧函數(shù)展開的空間分辨率保持一致，為2 159階次對應(yīng)5′的窗口。

根據(jù)式（2）、式（3）求得平滑參考面，將SRTM數(shù)據(jù)與平滑面相減可得殘余地形模型

2.2 RTM高程異常計(jì)算

根據(jù)位理論和Bruns公式可得［1，3］

式中：Tr為殘余地面模型引起的擾動位，r為積分主體到地面點(diǎn)的距離，G為萬有引力常數(shù)（G＝6.673×10－11m3s－2g－1），ρ為 質(zhì) 量 密 度（ρ＝2.670×103kgm－3），hr為參考面高程，h為流動點(diǎn)高程。利用RTM數(shù)據(jù)恢復(fù)地球重力場短波信號時，用hr＝0和h＝zRTM。

由于顧及EGM2008和RTM數(shù)據(jù)求解的高程異常與實(shí)際工程基準(zhǔn)面存在高程基準(zhǔn)不一致等系統(tǒng)偏差，因此，需要在上述基礎(chǔ)上利用多項(xiàng)式擬合系統(tǒng)誤差。

式中：a0，a1，a2為二元一次曲面擬合參數(shù)；（x，y）為擬合點(diǎn)的重心化坐標(biāo)。當(dāng)GPS/水準(zhǔn)點(diǎn)個數(shù)多于3個時，可通過最小二乘方法求解擬合參數(shù)。

3 實(shí)例分析

本文以某隧道洞外GPS控制網(wǎng)為例進(jìn)行方法驗(yàn)證，隧道全長約7km，在隧道進(jìn)出口以及斜井處布設(shè)15個GPS點(diǎn)，并按照三等水準(zhǔn)測量的要求對各點(diǎn)進(jìn)行聯(lián)測。圖1為GPS控制點(diǎn)分布圖，該隧道區(qū)域內(nèi)，隧道穿越的山嶺最高處可達(dá)1 030m，最低處約為540m，高差可達(dá)400m。從圖1可以看出該區(qū)黃土沖溝發(fā)育，該隧道山區(qū)呈“V”形。控制點(diǎn)GPS7847和GPS7848位于谷底周圍，高程較低，其余各點(diǎn)的高程大致相等。

圖1 控制點(diǎn)分布示意圖

首先根據(jù)式（1）求解各GPS控制點(diǎn)的EGM2008高程異常，然后按照式（2）～（4）構(gòu)建RTM數(shù)據(jù)，根據(jù)式（5）計(jì)算RTM高程異常。下面以控制點(diǎn)“1×4”為例介紹RTM高程異常。

圖2為控制點(diǎn)“1×4”為中心的SRTM高程數(shù)據(jù)，由圖可知：控制點(diǎn)位于山腳，控制點(diǎn)周圍高程變化比較明顯，高差約為200m；在半徑200km的范圍內(nèi)西部變化較為平穩(wěn)，而東部變化劇烈，高差可達(dá)1 500m，為明顯的山區(qū)地貌。

圖2 SRTM高程數(shù)據(jù)

圖3為由SRTM數(shù)據(jù)構(gòu)建的控制點(diǎn)“1×4”RTM數(shù)據(jù)，經(jīng)過去除參考面高程的RTM高程變化平緩，但地形急劇變化區(qū)域RTM變化較大，圖中可發(fā)現(xiàn)僅在河谷和重山區(qū)高程變化較大，高差可達(dá)到1 200m，表明RTM高程代表地球重力場小于EGM2008的最高分辨率（約為5′）的短波信號在山區(qū)不能夠忽視，需加以考慮。

圖3 RTM高程數(shù)據(jù)

圖4為控制點(diǎn)“1×4”的RTM高程異常中心局部圖，可得出控制點(diǎn)附近的地形對高程異常影響較大。RTM高程異常并不隨著距離的增大而消失，故確定合理的積分半徑也將影響RTM高程異常［2，5］。

圖4 RTM高程異常局部圖

圖5為控制點(diǎn)“1×4”的RTM高程異常數(shù)據(jù)與積分半徑間的關(guān)系?？砂l(fā)現(xiàn)：積分半徑達(dá)到160km時，RTM高程異常趨于穩(wěn)定。試驗(yàn)表明其他控制點(diǎn)在積分半徑達(dá)到160km時RTM高程異常也趨于穩(wěn)定，故本試驗(yàn)選取積分半徑R＝200km。

圖5 積分半徑與RTM高程異常

對所有控制點(diǎn)數(shù)據(jù)進(jìn)行EGM2008高程異常和RTM高程異常求解后，做以下分析：首先在大地高的基礎(chǔ)上進(jìn)行EGM2008高程異常改正，稱為方法1；然后繼續(xù)進(jìn)行RTM高程異常改正，稱為方法2；最后采用二元一次曲面擬合消除系統(tǒng)誤差，本試驗(yàn)固定GPS201、進(jìn)2、1×4和2×4四個控制點(diǎn)，計(jì)算擬合參數(shù)，進(jìn)而求解其他點(diǎn)的擬合高程異常改正數(shù)，稱為方法3。將3種方法求解的正常高與實(shí)際正常高比較，見表1，并計(jì)算各種方法的外符合精度

表1 3種方法的比較

式中：V表示檢查點(diǎn)求解的正常高與實(shí)測值的差值，n為GPS控制點(diǎn)的數(shù)目，t為固定點(diǎn)個數(shù)。

通過表1分析可得：

1）EGM2008在該山區(qū)的精度較高，與該區(qū)域內(nèi)的真實(shí)高程異常的中誤差為±4.2cm，表明EGM2008地球重力場模型在我國適用。

2）加入RTM高程異常，在一定程度上恢復(fù)了EGM2008的遺失信號，EGM2008高程異常和RTM高程異常組合與該區(qū)域內(nèi)的真實(shí)高程異常的均方根誤差為±2.6cm，較僅進(jìn)行EGM2008高程異常改正提高約40%，表明RTM高程異常能夠有效地提高GPS高程轉(zhuǎn)換的精度。

3）利用二元一次曲面擬合存在的系統(tǒng)誤差，可以有效地提高GPS高程轉(zhuǎn)換的精度，實(shí)現(xiàn)基準(zhǔn)轉(zhuǎn)換，使得求解結(jié)果與我國高程系統(tǒng)一致。

4 結(jié)束語

本文分析了EGM2008地球重力場的信號遺失問題，利用RTM數(shù)據(jù)恢復(fù)遺失信號；重點(diǎn)通過SRTM高程數(shù)據(jù)建立RTM高程數(shù)據(jù)，計(jì)算RTM高程異?；謴?fù)高程異常短波信號；利用二元一次曲面擬合系統(tǒng)誤差。通過實(shí)驗(yàn)對比：顧及EGM2008和殘差地形模型的高程異常較僅使用EGM2008求解的高程異常的精度有較大提高，利用二元一次曲面擬合系統(tǒng)誤差可使GPS高程轉(zhuǎn)換達(dá)到厘米級精度。同時利用EGM2008和殘差地形模型也為局部或區(qū)域性的似大地水準(zhǔn)面精化提供參考。該方法基于公共數(shù)據(jù)（EGM2008和SRTM），數(shù)據(jù)資源公共開放，使得其應(yīng)用于各類工程實(shí)踐成為可能。

［1］徐紹銓，李振洪，吳云孫.GPS高程擬合系統(tǒng)的研究［J］.武漢測繪科技大學(xué)學(xué)報，1999，24（4）：337－339.

［2］張同剛，岑敏儀，馮義從，等.地形起伏對GPS工程控制網(wǎng)高程異常的影響［J］.鐵道學(xué)報，2005，27（2）：79－80.

［3］游為，范東明，付淑娟，等.GPS高程轉(zhuǎn)換的新方法研究［J］.工程勘察，2009，37（3）：60－62.

［4］熊永良，路伯祥.高山區(qū)GPS網(wǎng)正常高求解方法研究［J］.西南交通大學(xué)學(xué)報，1997，32（2）：154－155.

［5］Hirt C.Prediction of Vertical Deflections From High－degree Spherical Harmonic Synthesis and Residual Terrain Model Data［J］.Journal of Geodesy，2010，84（1）：179－190.

［6］http://earth－info.nima.mil/GandG/wgs84/gravitymod/egm2008/index.html.

［7］陳俊勇，對SRTM3和GTOPO30地形數(shù)據(jù)質(zhì)量的評估［J］.武漢大學(xué)學(xué)報：信息科學(xué)版，2005，30（11）：941－942.