利用Wenzel加權譜組合法拓展GOCE重力場模型

張敏利 趙德軍 陳永祥，2

1 西安測繪信息技術總站，西安市西影路36號，710054

2 大地測量與地球動力學國家重點實驗室，武漢市徐東大街430號，430077

截止到2013年底，歐空局發布了第4 代GOCE重力場模型，其采用了27 個月的GOCE衛星數據，包括兩個產品，即時域法模型GO＿CONS＿GCF＿2＿TIM＿R4（簡稱Tim4）和直接法模型GO＿CONS＿GCF＿2＿DIR＿R4（簡稱Dir4）。Tim4僅采用純粹的GOCE 衛星數據，即重力梯度數據和GPS 高低衛衛跟蹤軌道數據，而Dir4除了采用GOCE 衛星數據外，還采用了25a的LAGEOS 1／2 激 光 測 衛 數 據、8a的GRACE 低低跟蹤衛星數據和K 波段距離變率數據。采用純GOCE衛星數據恢復的Tim4模型的累積大地水準面誤差和重力異常誤差分別為3.2cm 和0.9 mGal（100km 空間分辨率），采用200階的Tim4和201～2 190階EGM08組合的模型在德國的精度達到4.5cm（通過875個GPS／水準點的外部檢核）。

GOCE任務的基本目標是以前所未有的精度恢復中高階地球重力場，但對于許多大地測量和地球物理應用來說，還需要更高階次的模型。最簡單的方法就是采用GOCE＋EGM08的譜組合法，EGM08模型含有豐富的地面重力信息，但其低頻部分僅采用不足5a的GRACE 信息，因此，GOCE模型能有效地彌補EGM08模型中低頻的不足。譜組合法簡單實用，文獻［1］通過實測GPS水準驗證簡單組合模型比EGM08模型精度至少提高12%。文獻［2］用GOCO03S 和EGM08 的簡單組合模型計算的重力大地水準面與北大西洋沿岸驗潮站確定的海平面形態一致。文獻［3］用GOCO03S和歐洲重力場模型EGG08 的加權組合計算的大地水準面，在德國精度為2.4cm，高于EGM08 精度2.6 cm 以及歐洲重力場模型EGG08精度2.8cm，而第3代時域GOCE 模型Tim3與EGG08 的組合模型精度更高，達到2.1 cm。但文獻［3］采用高斯低通濾波來確定GOCE模型和超高階重力場模型（EGM08或EGG08）的權，其濾波帶寬很難確定，需要反復試算，且沒有顧及模型自身的精度。因而，本文嘗試采用Wenzel的加權譜組合來構建超高階重力場模型，最終達到了令人滿意的效果。

1 Wenzel加權組合模型

不同組織發布的重力場模型對應的參考橢球長半軸、地球引力常量與GRS80橢球都不一致，因此，首先需要將重力場位系數“標準化”到GRS80橢球上來：

σn為重力場模型的誤差階誤差：

δCnm、δSnm為模型系數的誤差。由誤差傳播定律，組合模型系數的誤差為：

2 模型精度分析

2.1 模型頻譜估計

根據 式（4）計 算 了Dir4、Tim4、EGM08、EIGEN－6C2四個模型的系數階誤差（圖1）。從中看出，對于純GOCE 數據構造的Tim4模型在46～179階精度優于EGM08，這是因為EGM08模型中低頻采用了GRACE 數據，而GRACE 衛星主要恢復中低頻重力場，GOCE 衛星主要恢復中高頻靜態重力場。Dir4 模型精度在214 階內優于EGM08，但是在214 階之后精度低于EGM08模型，是因為EGM08高頻部分采用衛星測高和地面重力測量數據來構建，故在214階之后體現出EGM08的優勢。這也說明，GOCE 模型和EGM08 模型在頻譜上具有互補性。EIGEN－6C2是德國地學中心發布的完全階次到1 949的重力場模型，其中低頻部分采用Dir2模型作為參考場，高頻部分采用EGM08，因而其階誤差也是從214階處開始減小。

圖1 不同重力場模型的誤差階誤差Fig.1 Root of error degree variances from different geopotential models

由式（5）計算Dir4、Tim4 分別和EGM08、EIGEN－6C2的組合模型的系數誤差，再根據式（3）計算組合模型的頻譜權（圖2）。從中看出，由于Dir4中低頻精度明顯優于EGM08，所以Dir4＋EGM08組合模型的權在150階內幾乎為1，也就是說，在150 階內幾乎用Dir4 位系數替換了EGM08相應的位系數。由于Tim4 精度低于EIGEN－6C2，因而其低階部分的權幾乎為0。

圖2 Wenzel頻譜組合權Fig.2 Spectral weights derived from Wenzel’s combination schema

還可以用大地水準面的累積誤差來評估模型：

R為地球平均半徑。根據式（6）將組合模型的大地水準面的累積誤差繪制于圖3。從中看出，所有組合模型的精度都要高于EGM08模型，其中Dir4＋EGM08精度最高，在250階處為2.8cm，在低于210 階處，其精度略優于Dir4＋EIGEN－6C2組合模型，但210階之后，二者迅速分化，但這兩個組合模型精度都高于EIGEN－6C2。

圖3 不同組合模型的大地水準面累積誤差Fig.3 Cumulative geoid errors derived from different combined models

2.2 外符合精度檢驗

由重力位系數計算大地水準面的方法見文獻［5］。本文借助于美國地理空間情報局NGA 提供的harmonic＿synth＿v02.f軟件，采用美國大地測量局NGS提供的GPS水準網數據GPSBM2009（包含美國大陸18 399個點的GPS／水準數據）。GPSBM2009數據的大地基準采用北美大地基準NAD83，其定位定向與ITRF 之間存在一定的系統差，因此，首先采用NGS提供的坐標框架和歷元轉換軟件Htdp，將大地坐標從NAD83框架轉換到ITRF2008框架［6］。將實測的大地水準面和EGM08模型計算的大地水準面之差繪制于圖4，從中看出，NAVD88高程基準存在明顯的東西方向和南北方向的傾斜，東西方向傾斜達到1.5m，原因歸咎于水準網系統誤差的累積［7］。

利用曲面函數模型對大地水準面之差進行擬合，然后對擬合后的殘差進行精度統計。文獻［8］給出了6種趨勢面函數模型，本文采用其中較簡單的四參數模型：

ΔN為大地水準面之差，v為殘差，a0、a1、a2、a3為曲面函數模型系數，Bi、Li為大地坐標。根據式（7）由最小二乘法求得函數模型系數和殘差，圖4擬合后的殘差分布見圖5。

圖4 大地水準面之差Fig.4 Distribution of geoid height differences

圖5 曲面擬合后大地水準面之殘差Fig.5 Distribution of geoid height residual errors after a curve fitting

表1給出了不同組合模型經曲面擬合后的殘差統計。從中看出，所有組合模型精度都優于EGM08模型，但效果卻沒有文獻［1－2］那么顯著，即使精度最高的Dir4＋EGM08 模型，精度也只提升了8%，原因可能是構建EGM08模型時使用了美國高精度的地面重力數據，而EGM08在美國精度已經很高了，提升空間較小。Tim4＋EIGEN－6C2和EIGEN－6C2精度一致，其原因從圖2可知，Tim4＋EIGEN－6C2的組合權很小，也就是說Tim4 在組合模型中幾乎沒有貢獻。相反，Dir4＋EGM08的組合權幾乎為1，即Dir4在計算大地水準面時貢獻很大，因而其組合模型精度最高。

表1 不同模型大地水準面差異統計／cmTab.1 Statistics of the geoid height differences from differentgravity field models／cm

3 結 語

1）聯合LAGEOS、GRACE 和GOCE 數據解算的直接法模型Dir4 精度明顯優于只采用GOCE數據的時域法模型Tim4。

2）大地水準面主要是由重力場中低位系數貢獻的，因而用GOCE＋EGM08 的組合模型能提高水準面計算精度。

3）大地水準面的累積誤差和美國實測數據表明，200階的Dir4和2 190 階的EGM08 組合模型表現最優，其精度最高。

4）從文獻［1］和文獻［9］對我國西藏地區的分析來看，在缺乏地面重力信息的地區，GOCE 模型能有效地彌補重力場中高頻信號，因而拓展的GOCE超高階重力場模型在這些地區優勢更明顯。

［1］張興福.多類重力場模型的精度分析及聯合確定GPS點正常高的方法［J］.測繪學報，2013，42（1）：6－12（Zhang Xingfu，et al.The Accuracy Analysis of Multiple Earth Gravity Field Models and Normal Height Determination of GPS Points Based on Combination Method［J］.Acta Geodaetica et Cartographica Sinica，2013，42（1）：6－12）

［2］Woodworth P L.Towards Worldwide Height System Unification Using Ocean Information［J］.Journal of Geodetic Science，2012（2）：302－318

［3］Rülke A.Unification of European Height System Realizations［J］.Journal of Geodetic Science，2013（2）：343－354

［4］Wenzel H G.Zur Geoidbestimmung Durch Kombination von Schwereanomalien und Einem Kugelfunktionsmodell mit Hilfe von Integralformeln［J］.Zeitschrift für Vermessungswesen，1981，106（3）：102－111

［5］王正濤，黨亞民，晁定波.超高階地球重力位模型確定的理論與方法［M］.北京：測繪出版社，2011（Wang Zhengtao，Dang Yamin，Chao Dingbo.Theory and Methodology of Ultra－High－Degree Geopotential Model Determination［M］.Beijing：Surveying and Mapping Press，2011）

［6］Soler T.Transforming Positions and Velocities between the International Terrestrial Reference Frame of 2000 and North American Datum of 1983［J］.Journal of Surveying Engineering，2004，4：49－55

［7］Roman D R.Geodesy，Geoid，and Vertical Datums：a Perspective from the U.S.National Geodetic Survey［C］.FIG Congress 2010：Facing the Challenges－Building the Capacity，Sydney，2010

［8］Kotsakis C.Evaluation of EGM08in Greece Using GPS and Leveling Heights［R］.IAG International Symposium on Gravity Geoid and Earth Observation，Chania，Greece，2008

［9］高春春.GOCE地球重力場模型在青藏地區的評價及分析［J］.測繪學報，2014，43（1）：21－29（Gao Chunchun.Evaluation and Analysis of GOCE Earth Gravity Field Model in Qinghai－Tibet Region［J］.Acta Geodaetica et Cartographica Sinica，2014，43（1）：21－29）