顧及多方向觀測值權(quán)比反演地球重力場的動力積分法

羅志才， 周浩， 鐘波,2*， 李瓊

1 武漢大學測繪學院， 武漢 430079 2 武漢大學地球空間環(huán)境與大地測量教育部重點實驗室， 武漢 430079 3 武漢大學測繪遙感信息工程國家重點實驗室， 武漢 430079

?

顧及多方向觀測值權(quán)比反演地球重力場的動力積分法

羅志才1,2,3， 周浩1， 鐘波1,2*， 李瓊1

1 武漢大學測繪學院， 武漢 430079 2 武漢大學地球空間環(huán)境與大地測量教育部重點實驗室， 武漢 430079 3 武漢大學測繪遙感信息工程國家重點實驗室， 武漢 430079

考慮到不同坐標系下各個方向觀測值對反演地球重力場的頻譜貢獻不同，建立了顧及多方向觀測值權(quán)比的動力積分法，并利用該方法反演了高精度的GOCE HL-SST衛(wèi)星重力場模型.首先，分析了不同坐標系下各個方向觀測值與地球重力場信息的響應關(guān)系，其中慣性系(IRF)下X、Z方向的觀測值分別對扇諧系數(shù)、帶諧系數(shù)最為敏感，Z方向的解算精度在全頻段均略高于X、Y方向；地固系(EFRF)下各個方向的獨立解算精度均與能量守恒法的解算精度相當；局部指北坐標系(LNOF)下X、Z和Y三個方向的解算精度依次遞減，且Y方向在47階附近有明顯“駝峰”現(xiàn)象.其次，比較了不同坐標系下顧及三個方向觀測值權(quán)比的加權(quán)解算模型，其中加權(quán)聯(lián)合解算模型精度在20至70階次均明顯優(yōu)于等權(quán)解算模型，在帶諧項和共振階次精度提升明顯，且LNOF下的加權(quán)聯(lián)合解算精度要優(yōu)于IRF和EFRF.最后，比較了GOCE和CHAMP衛(wèi)星的模型解算精度，采用本文計算方法，僅利用2個月GOCE軌道觀測值解算的模型精度優(yōu)于包含更長觀測時段信息的AIUB-CHAMP01S和EIGEN-CHAMP03S模型，且略優(yōu)于ASU-GOCE-2months模型.

動力積分法； 加權(quán)； GOCE； 地球重力場

1 引言

21世紀以來，隨著CHAMP、GRACE和GOCE衛(wèi)星計劃的成功實施，人類對地球重力場的認知達到了前所未有的高度.其中，各個重力衛(wèi)星均搭載了高頻GPS接收機，不僅為各個衛(wèi)星計劃的有效載荷提供了高精度的位置信息，又可應用于改善地球重力場中長波信息(Bock et al.，2011；J?ggi et al.，2011a；Bezděk et al.，2014).

基于衛(wèi)星軌道獲取地球重力場模型一直是衛(wèi)星重力學的研究熱點，特別是GOCE衛(wèi)星成功發(fā)射之后，國內(nèi)外不僅致力于梯度數(shù)據(jù)確定地球重力場的研究(Pail et al., 2011；鄭偉等，2011；萬曉云等，2012)，也掀起了基于GOCE高低衛(wèi)-衛(wèi)跟蹤數(shù)據(jù)(HL-SST,High-Low Satellite-to-Satellite Tracking)反演地球重力場模型研究的新熱潮，研究內(nèi)容主要包括兩個方面：一方面是基于不同方法的研究，另一方面是基于不同衛(wèi)星任務的比較.Baur等(2012)詳細討論了采用GOCE幾何軌道恢復地球重力場模型的加速度法，分析了GOCE衛(wèi)星獲取地球重力場長波信息的能力；J?ggi等(2011b)基于天體力學法反演GOCE HL-SST地球重力場模型，其結(jié)果表明幾何軌道的協(xié)方差信息對恢復地球重力場模型的影響較??；Baur等(2014)基于GOCE幾何軌道數(shù)據(jù)分析了天體力學法、短弧法、平均加速度法、瞬時加速度法和能量守恒法恢復HL-SST地球重力場信息的能力，結(jié)果表明除能量守恒法之外的其他方法解算精度均相當；J?ggi 等(2011a)、Bezděk(2014)等分別基于天體力學法、加速度法比較了CHAMP、GRACE和GOCE的反演精度，結(jié)果表明GOCE衛(wèi)星的解算精度優(yōu)于CHAMP和GRACE.特別地，Visser等(2014)和J?ggi等(2015)探討了GOCE衛(wèi)星軌道數(shù)據(jù)恢復時變重力場模型的可行性，也是近年來探尋GRACE衛(wèi)星任務和GRACE Follow-On衛(wèi)星計劃之間空白時期時變信號獲取方法的重要研究內(nèi)容.此外，鐘波等(2013)、游為等(2012)、Su和Fan(2013)、Huang和Fan(2012)、黃強等(2013)也分別利用加速度法、短弧積分法和能量守恒法等方法反演了多組GOCE HL-SST模型.

動力積分法雖然計算較為復雜，但其理論嚴密，具有可同時獲取高精度地球重力場模型和高精度動力學軌道等優(yōu)點(鄒賢才，2007).由于GPS衛(wèi)星觀測的固有特性，各個方向觀測值的精度各不相同，但在現(xiàn)有研究中，尚未討論動力積分法中不同坐標系下各個方向觀測值與地球重力場信號的響應關(guān)系和貢獻權(quán)比等關(guān)鍵問題.同時，Bruinsma等(2013)基于動力積分法，聯(lián)合激光測衛(wèi)、GRACE和GOCE等觀測數(shù)據(jù)獲取了GOCE衛(wèi)星直接解模型，但未對GOCE HL-SST模式下的動力積分法進行詳細討論.除此之外，國內(nèi)外尚無關(guān)于動力積分法獲取GOCE HL-SST地球重力場模型的詳細研究.鑒于此，本文將建立顧及多方向觀測值權(quán)比的動力積分法，并以GOCE HL-SST為例，詳細討論不同坐標系下各個方向觀測值與地球重力場信息的響應關(guān)系，并引入多方向觀測值的聯(lián)合解算公式，比較加權(quán)和等權(quán)的聯(lián)合解算模型精度；最后與CHAMP衛(wèi)星任務恢復的地球重力場模型以及國際上采用GOCE HL-SST數(shù)據(jù)恢復的最新重力場模型進行比較.

2 動力積分法反演地球重力場模型的基本原理

動力積分法是將衛(wèi)星星歷擾動作為地球重力場信息的泛函，通過建立軌道攝動與地球重力場位系數(shù)之間的關(guān)系，精密獲取地球重力場模型的方法.動力積分軌道的精度不僅受到初始狀態(tài)的影響，還會受到各種先驗力模型參數(shù)不確定性的影響，因此在地球重力場模型解算過程中通常一并求解.

2.1 基于HL-SST技術(shù)的動力積分反演方法

根據(jù)動力積分法的基本原理可知，積分軌道誤差主要源于初始狀態(tài)誤差ΔX0和先驗力模型參數(shù)誤差Δp，因此將軌道攝動ΔX(t)表示為二者的線性組合(周旭華，2005；肖云，2006；鄒賢才，2007；張興福，2007；游為，2011)：

(1)

利用公式(1)建立的觀測方程，可通過最小二乘法直接獲得位系數(shù)修正值、初始狀態(tài)誤差和力模型誤差等參數(shù).考慮到力模型誤差的累積影響，動力積分法通常分弧段進行：根據(jù)公式(2)將分弧段形成的法方程進行約化、疊加、求解(公式(2)和公式(3)中N均為法方程陣分塊，W為伴隨矩陣分塊)；利用求得的全局變量Xe修正地球重力場位系數(shù)，然后根據(jù)公式(3)更新初始狀態(tài)向量和加速度計校準參數(shù)等局部變量Xl，代入下一次循環(huán)中，直至獲得各類參數(shù)的最佳估值.

(2)

(3)

基于上述方法，即可根據(jù)HL-SST技術(shù)獲取初始狀態(tài)、加速度校準因子以及位系數(shù)等參數(shù)的最佳估值，實現(xiàn)加速度計的精密校準和動力學軌道的精密確定，并更新地球重力場模型.

2.2 不同坐標系下的動力積分法

CHAMP、GRACE和GOCE等重力衛(wèi)星提供的是地固系EFRF(Earth-Fixed Reference Frame)下的衛(wèi)星幾何軌道，而現(xiàn)有大部分學者均是在慣性系IRF(Inertial Reference Frame)下實現(xiàn)的動力積分法.為了分析不同坐標系下各個觀測值分量與地球重力場信息的響應關(guān)系，本文分別實現(xiàn)了地固系EFRF、慣性系IRF和局部指北坐標系LNOF(Local North-Oriented Frame)下的動力積分法，計算中需要用到的轉(zhuǎn)換關(guān)系如公式(4)所示：

(4)

2.3 多方向觀測值的聯(lián)合反演

HL-SST技術(shù)給出了三個方向的觀測值，在不同坐標系下均可依次表示為X、Y、Z三個正交方向.由于天線相位中心變化以及GPS觀測的固有特點，低軌衛(wèi)星的幾何軌道精度在各個方向不一致(Bock et al.，2011；J?ggi et al.，2009)，因此需要單獨考慮各個方向觀測值對反演結(jié)果的影響.首先將公式(2)簡化為：

(5)

其中，Nc為聯(lián)合解的法方程陣，是關(guān)于各個獨立解算模型法方程陣Ni的加權(quán)平均；i表示參與聯(lián)合解算的觀測值方向，對應于不同坐標系下的X、Y、Z三個正交方向；Ri為分辨率矩陣(Jackson，1972).

特別地，由公式(6)可知，聯(lián)合解算模型Xc可表示為關(guān)于獨立計算模型Xi的加權(quán)平均，各個分量權(quán)值為Ri.Ri為對角線占優(yōu)矩陣，可以利用Ri的對角線元素分析各個正交方向觀測值對各個聯(lián)合解算位系數(shù)的相對貢獻.考慮到：

R1+R2+…+Ri=I,

(7)

也可以直接使用Ri的對角線元素之和分析各個分量整體上對聯(lián)合解算的平均貢獻.

3 數(shù)值解算與分析

本文基于動力積分法編寫了HL-SST模式下的地球重力場反演軟件，所有解算過程均采用聯(lián)合OpenMP和MPI編寫的并行函數(shù)庫(周浩等，2011, 2015)，軌道積分和變分方程解算均采用Gauss-Jackson數(shù)值積分器(羅志才等，2013)，攝動力模型包括三體攝動力、固體潮、海潮、大氣潮、極潮和相對論效應，非保守力由加速度計觀測數(shù)據(jù)提供，具體描述見表 1.

表1 各項攝動力模型及其主要參數(shù)Table 1 Summary of perturbation models and their key elements

(8)

為分析三個方向觀測值的貢獻權(quán)比以及對聯(lián)合解的影響，本文將以GOCEHL-SST數(shù)據(jù)為例進行研究.其中，基于GOCEHL-SST技術(shù)恢復地球重力場模型需要Level2的幾何軌道SST_PKI_2、簡化動力學軌道SST_PRD_2以及Level1b的加速度計數(shù)據(jù)EGG_NOM_1b(Frommknechtetal.，2011；Gruberetal.，2010).精密軌道數(shù)據(jù)由瑞士伯爾尼大學天文研究所(AIUB) 提供，精度為2cm；三個加速度計數(shù)據(jù)可以由梯度儀3個坐標軸上對稱分布的6個加速度計給出.為減小計算量，在不影響解算精度的前提下，將剔除粗差的幾何法軌道和加速度計觀測數(shù)據(jù)降采樣至5s.需要說明的是，雖然GOCE衛(wèi)星引入無阻力控制技術(shù)，可以抵消大部分的大氣阻力等非保守力的作用，但仍有部分殘余，因此在表 1中依然引入了加速度計的校準因子.下文的數(shù)值解算與分析均采用上述預處理后的觀測數(shù)據(jù)完成.

3.1 不同坐標系下各方向觀測值的解算精度

以2009年12月的觀測數(shù)據(jù)為例，采用不同方向的觀測值解算了最大截斷階次為100的地球重力場模型，解算結(jié)果如圖1和圖2所示.其中，圖1反映了解算模型相對于“真實重力場模型”的位系數(shù)誤差譜，圖1a、1b和1c分別表示僅采用慣性系下X、Y、Z方向觀測值計算的位系數(shù)誤差譜，圖1d是三方向等權(quán)解算的位系數(shù)誤差譜；圖2至圖4給出了在不同坐標系下采用各個觀測分量及其聯(lián)合解的階誤差RMS.

在人衛(wèi)軌道研究中，慣性坐標系IRF的XY平面為地球赤道面，X軸指向J2000.0歷元的春分點，Z軸指向地球平均旋轉(zhuǎn)軸.在IRF下，衛(wèi)星可看作是沿著一個橢圓面飛行的，而地球相對于這一軌道面旋轉(zhuǎn)，旋轉(zhuǎn)速度為地球自轉(zhuǎn)平均角速度.因此，若軌道傾角為90°，X、Y方向的軌道攝動完全反映了東西向的重力場信息，而Z方向的軌道攝動完全反映了南北向的重力場信息.已有研究結(jié)果表明，GRACE飛行模式對帶諧系數(shù)部分較為敏感，而SWARM飛行模式對扇諧系數(shù)部分較為敏感(Wang，2011；Zhou et al.，2014).根據(jù)上述討論，若軌道傾角為90°，理論上IRF下X、Y方向的觀測值對扇諧系數(shù)敏感，而Z方向的觀測值對帶諧系數(shù)敏感.GOCE衛(wèi)星軌道傾角為96.7°，相當于繞X軸作6.7°旋轉(zhuǎn).因此，X方向的觀測值主要反映了東西向的重力場信息，對扇諧系數(shù)最為敏感，對帶諧系數(shù)最不敏感，導致扇諧系數(shù)解算精度最高，帶諧系數(shù)解算精度最低(圖1a)，反映在圖2所示的階誤差RMS，表現(xiàn)為其低階次誤差最大；Y方向同時反映了東西向和南北向的重力場信息，25階之前的帶諧系數(shù)與Z方向的精度相當，60階之前的扇諧系數(shù)與X方向的精度相當(圖1b)；Z方向也同時反映了東西向和南北向的重力場信息，但以南北向為主，導致帶諧系數(shù)精度最高，扇諧部分精度較低，特別是75階之后誤差增加較為明顯(圖1c).特別地，Y方向50階之后的田諧系數(shù)誤差增加明顯，Z方向在50至75階的近帶諧項誤差增加明顯，具體原因有待進一步研究.綜合上述討論，反映在圖2所示的階誤差RMS上，其結(jié)果是：在73階之前，Y方向的解算精度優(yōu)于X方向的解算精度；73階之后，Y方向的解算精度不及X方向的解算精度；Z方向的解算精度在全頻段均較高.

圖1 IRF下各種解的位系數(shù)誤差譜(a) X方向； (b) Y方向； (c)Z方向； (d)等權(quán)聯(lián)合解.Fig.1 Coefficient error spectrum in IRF(a)—(c) Solutions solely derived from X, Y, Z component respectively; (d) Solution derived from the combination of three components with equal weight.

圖2 IRF下不同方向的獨立解和聯(lián)合解模型的精度Fig.2 Accuracy of individual and combined solutions in IRF

圖3 EFRF下不同方向的獨立解和聯(lián)合解模型的精度Fig.3 Accuracy of individual and combined solutions in EFRF

圖4 LNOF下不同方向的獨立解和聯(lián)合解模型的精度Fig.4 Accuracy of individual and combined solutions in LNOF

3.2 不同坐標系下顧及權(quán)比的聯(lián)合解精度

基于3.1節(jié)的討論可知，不同坐標系下各個方向觀測值與地球重力場信息的頻譜響應各不相同；同時，考慮到GPS觀測的固有特性，本文計算了顧及各個方向觀測權(quán)比的聯(lián)合解.首先，基于公式(6)得到了各個方向觀測值對應的分辨率矩陣Ri，三者之和與單位陣I之差均達到了10-12量級，可認為符合公式(7)的檢核條件；其次，將各個分辨率矩陣的對角線元素相加，可以得到X、Y、Z三個方向觀測值對聯(lián)合解的貢獻權(quán)比；最后，基于該權(quán)比計算了IRF、EFRF和LNOF下的加權(quán)聯(lián)合解，其精度分別見圖2的IRF(optimal)、圖3的EFRF(optimal)和圖4的LNOF(optimal).顯然地，采用等權(quán)處理方式的聯(lián)合解精度均優(yōu)于采用單方向觀測值的解算精度，而加權(quán)聯(lián)合解精度均優(yōu)于等權(quán)聯(lián)合解.

為了更加清晰地比較兩種不同處理方式下的解算精度，圖5給出了LNOF下等權(quán)和加權(quán)聯(lián)合解相對于“真實模型”的位系數(shù)誤差譜，對比圖5a和圖5b可知，顧及各個方向觀測值的權(quán)比后，其聯(lián)合解算模型的帶諧項系數(shù)明顯優(yōu)于等權(quán)模式；由圖5c中二者的差值可知，顧及各個方向觀測值的權(quán)比也能夠提高共振階次附近的部分位系數(shù)精度.此外，圖6給出了相對于極區(qū)“真實模型”的重力異常.由于GOCE衛(wèi)星軌道傾角為96.7°，其星下點軌跡在兩極均存在6.7°的極空白，所以利用解算模型計算的該區(qū)域重力異常的誤差較大(圖6中黑色圓內(nèi))，即使采用加權(quán)處理可提高該部分的精度，為了進一步提高極區(qū)重力異常精度，需要引入其他類型觀測數(shù)據(jù)聯(lián)合求解；而在極空白影響區(qū)域外，采用加權(quán)處理方式，也能夠提高部分區(qū)域的重力異常精度.

圖7給出了在不同坐標系下，采用等權(quán)處理方式得到的聯(lián)合解算模型之間的精度差異.由于旋轉(zhuǎn)矩陣的不變特性，三種坐標系下的解算精度差值均在10-16量級以下，可認為在等權(quán)情況下(各個方向觀測值的貢獻權(quán)比均為33.3%)，三種坐標系下的解是等價的.

將公式(6)中各個分辨率矩陣Ri的對角線元素相加，可以得到X、Y和Z三個方向觀測值對聯(lián)合解的貢獻權(quán)比.由于各個坐標系下不同方向觀測值與地球重力場的頻譜響應不同，基于公式(6)計算的三個方向觀測值對聯(lián)合解的貢獻比例各不相同，具體比值見表 2.其中，不論在何種坐標系下，X方向的權(quán)最大，均超過了40%；Y和Z方向的權(quán)比均小于33.3%，即在基于公式(6)的加權(quán)解算中，Y和Z方向的觀測值均經(jīng)過了降權(quán)處理.

表2 不同坐標系下各個分量對聯(lián)合解的貢獻權(quán)比Table 2 Contribution of each component for the combined solutions in different frame

為了便于比較不同坐標系下顧及各個方向觀測值權(quán)比的聯(lián)合解精度，圖8給出了三種坐標系下加權(quán)聯(lián)合解的階誤差RMS.不同于等權(quán)聯(lián)合解，加權(quán)聯(lián)合解在不同坐標系下的精度各不相同：在前10階和70階次之后，三種坐標系下的聯(lián)合解與等權(quán)解的精度相當；在10至70階次，EFRF與IRF下的加權(quán)聯(lián)合精度幾乎一致，只有在極少數(shù)階次略優(yōu)；而在10至70階次，LNOF下的加權(quán)聯(lián)合解精度略優(yōu)于EFRF和IRF下的聯(lián)合解精度.

由圖2可知，在IRF下僅使用X方向觀測值得到的解算模型的階誤差RMS最大，根據(jù)誤差傳播定律，理論上應該減少該方向的權(quán)值；而由表2可知，IRF系下的X方向權(quán)比為42.0%，即基于公式(6)的聯(lián)合處理方式對X方向觀測值做了升權(quán)處理.為檢驗本文加權(quán)的合理性，在解算中將表 2中IRF系下X方向和Z方向的權(quán)比互換，即對X方向做降權(quán)處理，對Z方向做升權(quán)處理，得到的解算模型精度如圖8中IRF(XZinverse)所示.采用該聯(lián)合解算方式得到的精度不但低于加權(quán)聯(lián)合解，在高階次還不及等權(quán)聯(lián)合解.根據(jù)上述討論，基于公式(6)的聯(lián)合解，不僅考慮了各個方向獨立解的精度，還考慮了各個方向觀測值之間的相互關(guān)系.因此，在后續(xù)研究中建議均采用本文的加權(quán)聯(lián)合解算方式；而在GOCE HL-SST反演中，則建議在LNOF下實現(xiàn)加權(quán)聯(lián)合解算.

3.3 顧及權(quán)比的GOCE解精度評定

CHAMP、GRACE和GOCE衛(wèi)星雖然任務各有不同，但均搭載了高精度GPS接收機，可同時實現(xiàn)衛(wèi)星精密定軌和地球重力場低頻信息的獲取.其中，GRACE和CHAMP衛(wèi)星軌道傾角分別為89.0°和87.3°，可以不考慮極空白的影響，二者均搭載了BlackJack接收機，采用了Chockering天線.考慮到上述相似性，且CHAMP衛(wèi)星軌道高度(450 km)低于GRACE衛(wèi)星(500 km)，因此本文僅比較了CHAMP與GOCE衛(wèi)星任務的HL-SST解算模型精度.

圖5 LNOF下各種解的位系數(shù)誤差譜(a)等權(quán)；(b)加權(quán)；(c)二者之差.Fig.5 Coefficient error spectrum in LNOF(a)—(b) Solutions derived from the combination of three components with equal weight and different weight respectively; (c) Differences of these two solutions.

圖6 等權(quán)解和加權(quán)解計算的極區(qū)重力異常差異Fig.6 Difference of polar gravity anomaly computed from equal and different weight

圖7 不同坐標系下等權(quán)解算結(jié)果的差異Fig.7 Comparison of the equal weight solutions derived from different frame

圖8 不同坐標系下加權(quán)聯(lián)合解算結(jié)果的差異Fig.8 Comparison of the weighted solutions derived from different frame

AIUB-CHAMP01S和EIGEN-CHAMP03S模型分別采用1年和2.8年的CHAMP觀測數(shù)據(jù)反演，其中AIUB系列采用天體力學法計算，EIGEN系列采用動力積分法計算，且EIGEN-CHAMP03S從70階次開始做正則化處理，各個模型的階誤差RMS如圖9所示.由于HL-SST模式恢復地球重力場模型的能力隨著高度增加而衰減，GOCE衛(wèi)星的軌道高度(260 km)遠低于CHAMP衛(wèi)星(450 km)，因此，在顧及三個方向觀測值權(quán)比的前提下，僅采用2個月觀測值解算的WHU-GOCE-2months模型精度優(yōu)于AIUB-CHAMP01S和EIGEN-CHAMP03S.

ASU-GOCE-2months模型是Bezděk等(2014)基于加速度法計算的重力場模型，計算中采用了2009年11月至12月共計兩個月的GOCE HL-SST數(shù)據(jù)，截斷階次為75.對比采用同時期觀測數(shù)據(jù)計算得到的WHU-GOCE-2months模型可知，二者在前60階的階誤差RMS趨于一致；采用本文顧及多方向觀測權(quán)比的動力積分法，60階之后的WHU-GOCE-2months模型優(yōu)于加速度法解算的ASU-GOCE-2months模型.此外，對比分別采用2個月和3個月觀測數(shù)據(jù)得到的WHU-GOCE-2months和WHU-GOCE-3months可知，在觀測精度一致的前提下，增加多余觀測可有效提升解算精度.根據(jù)上述討論，采用共計4年的GOCE HL-SST觀測數(shù)據(jù)，解算精度有望實現(xiàn)全頻段提升.

圖9 WHU-GOCE系列模型與其他模型的精度比較Fig.9 Comparison between WHU-GOCE and other typical solutions

4 結(jié)論

本文在基于HL-SST技術(shù)恢復地球重力場基本原理的基礎(chǔ)上，給出了不同坐標系下的動力積分法，并引入了多方向觀測值的聯(lián)合解算公式，通過計算多組GOCE HL-SST衛(wèi)星重力場模型，得到如下結(jié)論：

(1)不同坐標系下各個方向觀測值與地球重力場信號的響應關(guān)系各不相同，IRF下X方向主要反映了東西向的重力場信息，對扇諧系數(shù)最為敏感；Y方向同時反映了東西向和南北向的重力場信息；Z方向主要反映了南北向的重力場信息，對帶諧系數(shù)最為敏感；Z方向的階誤差RMS在全頻段優(yōu)于X、Y方向；EFRF下各個方向的解算精度差別較小，且與能量守恒法的解算精度相當；LNOF系下X方向的階誤差RMS最小，Z方向次之，Y方向最差，且在共振階次47階出現(xiàn)最大的“駝峰”現(xiàn)象.

(2)顧及三個方向觀測值權(quán)比的聯(lián)合解算方法能夠提高解算精度，且在帶諧系數(shù)和共振階次部分表現(xiàn)最為明顯；在等權(quán)聯(lián)合解算中，不同坐標系下的解算精度相同，而由于顧及了三方向觀測值的相關(guān)性， LNOF系下的加權(quán)聯(lián)合解精度最高.

(3)采用2個月 GOCE衛(wèi)星HL-SST數(shù)據(jù)，基于本文方法計算的WHU-GOCE-2months重力場模型精度與基于加速度法計算的ASU-GOCE-2month模型精度大體相當，且在60階次之后更優(yōu)；WHU-GOCE-2months模型精度優(yōu)于分別采用1年和2.8年CHAMP衛(wèi)星觀測值解算的AIUB-CHAMP01S和EIGEN-CHAMP03S模型.

總體而言，采用本文的顧及多方向觀測值權(quán)比的動力積分法解算精度較高，建議在后續(xù)的研究中采用本文的處理方法.目前，對不同精度觀測值的加權(quán)處理方式較多，且方法之間的數(shù)學原理各不相同，同時局限于科研人員對衛(wèi)星重力數(shù)據(jù)的數(shù)學物理特性的認知，針對衛(wèi)星重力數(shù)據(jù)的加權(quán)處理方式仍需要進一步展開.

致謝 感謝ESA提供解算所需的GOCE衛(wèi)星觀測數(shù)據(jù)，感謝Bezděk博士對重力場解算過程中相關(guān)問題的幫助，感謝兩位匿名審稿專家提出的寶貴修改意見.

Baur O, Reubelt T, Weigelt M, et al. 2012. GOCE orbit analysis: Long-wavelength gravity field determination using the acceleration approach.AdvancesinSpaceResearch, 50(3): 385-396.

Baur O, Bock H, H?ck E, et al. 2014. Comparison of GOCE-GPS gravity fields derived by different approaches.JournalofGeodesy, 88(10): 959-973.

Bock H, J?ggi A, Meyer U, et al. 2011. GPS-derived orbits for the GOCE satellite.JournalofGeodesy, 85(11): 807-818.

Bruinsma S L, F?rste C, Abrikosov O, et al. 2013. The new ESA satellite-only gravity field model via the direct approach.GeophysicalResearchLetters, 40(14): 3607-3612.

Frommknecht B, Lamarre D, Meloni M, et al. 2011. GOCE level 1b data processing.JournalofGeodesy, 85(11): 759-775.

Gruber T, Rummel R, Abrikosov O, et al. 2010. GOCE level 2 product data handbook. GO-MA-HPF-GS-0110. Huang Q, Fan D M, You W. 2013. Recovery of earth′s gravitational field model based on GOCE satellite orbits.GeomaticsandInformationScienceofWuhanUniversity(in Chinese), 38(8): 907-910.

Huang Q, Fan D M. 2012. The average acceleration approach applied to gravity coefficients recovery based on GOCE orbits.GeodesyandGeodynamics, 3(4): 18-22.

Jackson D D. 1972. Interpretation of inaccurate, insufficient and inconsistent data.GeophysicalJournalInternational, 28(2): 97-109.

J?ggi A, Dach R, Montenbruck O, et al. 2009. Phase center modeling for LEO GPS receiver antennas and its impact on precise orbit determination.JournalofGeodesy, 83(12): 1145-1162.

J?ggi A, Bock H, Prange L, et al. 2011a. GPS-only gravity field recovery with GOCE, CHAMP, and GRACE.AdvancesinSpaceResearch, 47(6): 1020-1028.

J?ggi A, Prange L, Hugentobler U. 2011b. Impact of covariance information of kinematic positions on orbit reconstruction and gravity field recovery.AdvancesinSpaceResearch, 47(9): 1472-1479.

J?ggi A, Bock H, Meyer U, et al. 2015. GOCE: assessment of GPS-only gravity field determination.JournalofGeodesy, 89(1): 33-48. Jekeli C. 1999. The determination of gravitational potential differences from satellite-to-satellite tracking.CelestialMechanicsandDynamicalAstronomy, 75(2): 85-101.

Luo Z C, Zhou H, Zhong B, et al. 2013. Analysis and validation of Gauss-Jackson integral algorithm.GeomaticsandInformationScienceofWuhanUniversity(in Chinese), 38(11): 1364-1368. McCarthy D D, Petit G. 2004. IERS conventions (2003). International Earth Rotation and Reference Systems Service, GERMANY.Pail R, Bruinsma S, Migliaccio F, et al. 2011. First GOCE gravity field models derived by three different approaches.JournalofGeodesy, 85(11): 819-843.

Su Y, Fan D M. 2013. Gravity field recovery from GOCE orbits using the energy conservation approach.GeodesyandGeodynamics, 4(2): 40-46.Visser P N A M, van der Wal W, Schrama E J O, et al. 2014. Assessment of observing time-variable gravity from GOCE GPS and accelerometer observations.JournalofGeodesy, 88(11): 1029-1046.

Wan X Y, Yu J H, Zeng Y Y, et al. 2012. Frequency analysis and filtering processing of gravity gradients data from GOCE.ChineseJ.Geophys. (in Chinese), 55(9): 2909-2916, doi: 10.6038/j.issn.0001-5733.2012.09.010.

Wang X X. 2011. Time-variable gravity field determination from satellite constellations [Ph. D. thesis]. Technische Universit?t München.

Wang Z T. 2005. Theory and methodology of earth gravity field recovery by satellite-to-satellite tracking data [Ph. D. thesis] (in Chinese). Wuhan: Wuhan University.

Xiao Y. 2006. Research on the earth gravity field recovery from satellite-to-satellite tracking data [Ph. D. thesis] (in Chinese). Zhengzhou: PLA Information Engineering University.Yi W Y, Rummel R, Gruber T. 2013. Gravity field contribution analysis of GOCE gravitational gradient components.StudiaGeophysicaetGeodaetica, 57(2): 174-202.

You W. 2011. Theory and methodology of earth′s gravitational field model recovery by LEO data [Ph. D. thesis] (in Chinese). Chengdu: Southwest Jiaotong University.

You W, Fan D M, He Q B. 2012. Recovering earth′s gravitational field model using GOCE satellite orbits.GeomaticsandInformationScienceofWuhanUniversity(in Chinese), 37(3): 294-297.

Zhang X F. 2007. The earth′s field model recovery on the basis of satellite-to-satellite tracking missions [Ph. D. thesis] (in Chinese). Shanghai: Tongji University.Zheng W, Shao C G, Luo J, et al. 2006. Numerical simulation of Earth′s gravitational field recovery from SST based on the energy conservation principle.ChineseJ.Geophys. (in Chinese), 49(3): 712-717.Zheng W, Hsu H T, Zhong M, et al. 2011. Accurate and rapid determination of GOCE Earth′s gravitational field using time-space domain method associated with Kaula regularization.ChineseJ.Geophys. (in Chinese), 54(1): 14-21, doi: 10.3969/j.issn.0001-5733.2011.01.003.

Zheng W. 2015. Theory and Approach of Satellite Gravity Inversion on the Basis of Energy Balance Principle (in Chinese). Beijing: Science Press.

Zhong B. 2010. Study on the determination of the earth′s gravity field from satellite gravimetry mission GOCE [Ph. D. thesis] (in Chinese). Wuhan: Wuhan University.

Zhong B, Luo Z C, Zhou H. 2013. Gravity field recovery using satellite average accelerations derived from high-low satellite-to-satellite tracking.GeomaticsandInformationScienceofWuhanUniversity(in Chinese), 38(8): 902-906.

Zhou H, Zhong B, Luo Z C, et al. 2011. Application of parallel algorithms based on OpenMP to satellite gravity field recovery.JournalofGeodesyandGeodynamics(in Chinese), 31(5): 123-127.

Zhou H, Luo Z C, Zhong B, et al. 2015. MPI parallel algorithm in satellite gravity field model inversion on the basis of least square method.ActaGeodaeticaetCartographicaSinica(in Chinese), 44(8): 833-839.

Zhou H, Luo Z C, Zhong B, et al. 2014. A simulation study of Earth gravity field model inversion from SWARM mission.EGUGeneralAssemblyConferenceAbstracts, 16: 5888.

Zhou X H. 2005. Study on satellite gravity and its application [Ph. D. thesis] (in Chinese). Wuhan: Institute of Geodesy and Geophysics, Chinese Academy of Sciences.

Zou X C. 2007. Theory of satellite orbit and earth gravity field determination [Ph. D. thesis] (in Chinese). Wuhan: Wuhan University.

附中文參考文獻

黃強, 范東明, 游為. 2013. 利用GOCE衛(wèi)星軌道數(shù)據(jù)反演地球重

力場模型. 武漢大學學報·信息科學版, 38(8): 907-910.

羅志才, 周浩, 鐘波等. 2013. Gauss-Jackson積分器算法分析與驗證. 武漢大學學報·信息科學版, 38(11): 1364-1368.

萬曉云, 于錦海, 曾艷艷. 2012. GOCE引力梯度的頻譜分析及濾波. 地球物理學報, 55(9): 2909-2916, doi: 10.6038/j.issn.0001-5733.2012.09.010.

王正濤. 2005. 衛(wèi)星跟蹤衛(wèi)星測量確定地球重力場的理論與方法[博士論文]. 武漢: 武漢大學.

肖云. 2006. 基于衛(wèi)星跟蹤衛(wèi)星數(shù)據(jù)恢復地球重力場的研究[博士論文]. 鄭州: 中國人民解放軍信息工程大學.

游為. 2011. 應用低軌衛(wèi)星數(shù)據(jù)反演地球重力場模型的理論和方法[博士論文]. 成都: 西南交通大學.

游為, 范東明, 賀全兵. 2012. 利用GOCE衛(wèi)星軌道反演地球重力場模型. 武漢大學學報·信息科學版, 37(3): 294-297.

張興福. 2007. 應用低軌衛(wèi)星跟蹤數(shù)據(jù)反演地球重力場模型[博士論文]. 上海: 同濟大學.

鄭偉, 邵成剛, 羅俊等. 2006. 基于衛(wèi)-衛(wèi)跟蹤觀測技術(shù)利用能量守恒法恢復地球重力場的數(shù)值模擬研究. 地球物理學報, 49(3): 712-717.

鄭偉, 許厚澤, 鐘敏等. 2011. 基于時空域混合法利用Kaula正則化精確和快速解算GOCE地球重力場. 地球物理學報, 54(1): 14-21, doi: 10.3969/j.issn.0001-5733.2011.01.003.

鄭偉. 2015. 基于能量守恒原理的衛(wèi)星重力反演理論與方法. 北京: 科學出版社.

鐘波. 2010. 基于GOCE衛(wèi)星重力測量技術(shù)確定地球重力場的研究[博士論文]. 武漢: 武漢大學.

鐘波, 羅志才, 周浩. 2013. SST-HL模式下利用衛(wèi)星均值加速度反演地球重力場. 武漢大學學報·信息科學版, 38(8): 902-906.

周浩, 鐘波, 羅志才等. 2011. OpenMP并行算法在衛(wèi)星重力場模型反演中的應用. 大地測量與地球動力學, 31(5): 123-127.

周浩, 羅志才, 鐘波等. 2015. 利用最小二乘直接法反演衛(wèi)星重力場模型的MPI并行算法. 測繪學報, 44(8): 833-839.

周旭華. 2005. 衛(wèi)星重力及其應用研究[博士論文]. 武漢: 中國科學院測量與地球物理研究所.

鄒賢才. 2007. 衛(wèi)星軌道理論與地球重力場模型的確定[博士論文]. 武漢: 武漢大學.

(本文編輯 何燕)

Dynamic integral approach based on weighted multi-direction observations for inversion of the earth′s gravity field

LUO Zhi-Cai1,2,3， ZHOU Hao1， ZHONG Bo1,2*， LI Qiong1

1SchoolofGeodesyandGeomatics,WuhanUniversity,Wuhan430079,China2KeyLaboratoryofGeospaceEnvironmentandGeodesy,MinistryofEducation,WuhanUniversity,Wuhan430079,China3StateKeyLaboratoryofInformationEngineeringinSurveying,MappingandRemoteSensing,WuhanUniversity,Wuhan430079,China

The precise GPS high-low satellite-to-satellite tracking (HL-SST) data is the key supplement to recover the long-wavelength information of the earth′s gravity field for the gravity field and steady-state ocean circulation explorer (GOCE) mission. The dynamic integral approach is one of the efficient techniques to determine the spherical harmonic coefficients (SHCs) from GOCE HL-SST observations. However, the traditional dynamic integral approach is based on the inertial reference frame (IRF), and the discrepancy of spectral contribution of multi-direction observations are not fully considered in the SHCs determination.To analyze the contribution of observations in different directions, we build a new dynamic integral approach which takes the weighted ratios of multi-direction observations into consideration. Meanwhile, a dynamic integral approach at different reference frames is also established, which is used to analyze the response relationship between the observations at different frames. The high precision GOCE HL-SST data are used to evaluate the new dynamic integral approach. In the SHCs determination, the conservative perturbing forces are modeled by precise background models, and the non-conservative perturbing forces are observed by onboard instruments.Using the new dynamic integral approach, the response relationship between multi-direction observations in different frames and earth gravity information is analyzed firstly. In the IRF, the observations fromXandZdirections are most sensitive to sectorial and zonal coefficients, respectively, and the accuracy of the model recovered fromZdirection is higher than those recovered fromXandYdirections in the whole frequency band. In the earth-fixed reference frame (EFRF), all of the solutions recovered from individual components are of similar accuracy over all harmonic degrees, and approximately equal to the solution recovered from the energy balance approach (EBA). In the local north-oriented frame (LNOF), the accuracy of solutions decreases fromX,ZtoYdirections, and there is a hump peaking around degrees 47 in the solution recovered solely fromYdirection. Secondly, in terms of the weighted solutions which take the different contribution of each component into consideration, their accuracies between 20 to 70 degrees are higher than the equal weighted ones, and the coefficients in the zonal and resonant area are improved obviously. At the north pole and south pole, the gravity anomalies derived from weighted solutions also present better performance than equal weighted solutions. Although the equal weighted solutions have a good consistency in different frames due to the rotation invariant features, the weighted solutions, which adequately take the relationships between all components, have separate accuracy curves, and they are decreased from IRF, EFRF to LNOF. Thirdly, we compare the gravity field solutions based on 2 months of GOCE HL-SST with those obtained from the challenge mini-satellite payload (CHAMP) mission. Because of considering the contribution of observations in different directions, the GOCE HL-SST solution determined with our dynamic integral approach is slightly better than ASU-GOCE-2months model. Although only using 2 months of GOCE orbits, our solutions are much better than AIUB-CHAMP01S and EIGEN-CHAMP03S models, which are determined from 1 year and 2.8 year of CHAMP data, respectively.From the comparison with equal weighted dynamic integral approach, the numerical results indicate the new dynamic integral approach, which fully considers the contribution of multi-directions, is more suitable for the SHCs determination than the traditional approach. The numerical results also demonstrate that the new dynamic integral approach performance in LNOF is superior to that in IRF and EFRF. Therefore, it is more preferable to adopt the weighted dynamic integral approach in LNOF for SHCs determinations with HL-SST data.

Dynamic integral approach; Weighting; GOCE; Earth gravity field model

羅志才， 周浩， 鐘波等. 2015. 顧及多方向觀測值權(quán)比反演地球重力場的動力積分法.地球物理學報，58(9):3061-3071,

10.6038/cjg20150904.

Luo Z C, Zhou H, Zhong B, et al. 2015. Dynamic integral approach based on weighted multi-direction observations for inversion of the earth′s gravity field.ChineseJ.Geophys. (in Chinese),58(9):3061-3071,doi:10.6038/cjg20150904.

10.6038/cjg20150904

P223

2015-02-09，2015-06-12收修定稿

國家自然科學基金項目(41131067，41174020，41374023，41474019)；地理信息工程國家重點實驗室開放基金項目(SKLGIE2013-M-1-3)；地球空間環(huán)境與大地測量教育部重點實驗室開放基金項目(13-02-05)；大地測量與地球動力學國家重點實驗室開放基金項目(SKLGED2015-1-3-E)資助.

羅志才，男，1966年生，工學博士，教授，博士生導師，現(xiàn)主要從事物理大地測量學和衛(wèi)星重力學研究.E-mail:zhcluo@sgg.whu.edu.cn*通訊作者 鐘波，男，1980年生，工學博士，講師，主要從事衛(wèi)星重力數(shù)據(jù)處理及應用研究.E-mail:bzhong@sgg.whu.edu.cn