基于地面重力的衛(wèi)星重力梯度檢校方法

瞿慶亮， 常曉濤* ， 朱廣彬， 于勝文， 周苗， 劉偉

1 山東科技大學測繪科學與工程學院， 山東青島 266590 2 自然資源部國土衛(wèi)星遙感應(yīng)用中心， 北京 100048

0 引言

地球重力場衛(wèi)星可以高效探測地球重力場的精細結(jié)構(gòu)，尤其是2009年歐洲航天局發(fā)射的GOCE(gravity field and steady-state ocean circulation explorer)衛(wèi)星，作為全球首顆以重力梯度觀測模式恢復(fù)高精度高分辨率地球重力場的地球重力衛(wèi)星，能夠以100 km空間分辨率構(gòu)建1～2 cm精度的全球大地水準面模型和1 mGal精度的地球重力場模型(Bouman et al., 2011, 2016; Bouman and Fuchs, 2012; Van Der Meijde et al., 2015; Xu et al., 2017; Siemes, 2018)，為人類更深入地了解地球的內(nèi)部結(jié)構(gòu)和研究海洋環(huán)流提供了豐富的地球觀測基礎(chǔ)數(shù)據(jù).GOCE衛(wèi)星的核心載荷是星載靜電重力梯度儀，其具有高靈敏特性，可以高精度探測地球重力場的中短波信號，配合SST-hl(Satellite-to-Satellite Tracking in the high-low mode)技術(shù)獲取的地球重力場的長波信號，可實現(xiàn)高精度地球重力場的恢復(fù)(蔡林等, 2012; 朱廣彬等, 2018; Liang et al., 2020; 劉煥玲等, 2020).但梯度儀受衛(wèi)星空間環(huán)境及觀測條件影響，觀測值中存在難以避免的噪聲，這會增加重力場恢復(fù)數(shù)據(jù)處理的難度，并進一步對重力場恢復(fù)的精度和分辨率產(chǎn)生不利影響.在軌檢校是顧及梯度儀存在觀測噪聲的實際情況，通過對其在軌標定最大限度地消除觀測系統(tǒng)誤差，提高重力場恢復(fù)的精度和分辨率.因此，在衛(wèi)星總體任務(wù)確定的情況下，有效的在軌檢校顯然能降低衛(wèi)星工程系統(tǒng)對梯度儀研制的指標要求，對控制衛(wèi)星工程成本、提高觀測質(zhì)量、實現(xiàn)工程總體要求具有重要意義.

研制梯度模式地球重力場探測衛(wèi)星是我國《國家民用空間基礎(chǔ)設(shè)施中長期發(fā)展規(guī)劃(2015—2025年)》的重要內(nèi)容.突破衛(wèi)星在軌檢校相關(guān)技術(shù)，面向我國發(fā)展梯度模式重力衛(wèi)星的總體任務(wù)，研究提出針對性的衛(wèi)星梯度儀在軌檢校方法是順利完成規(guī)劃任務(wù)、確保工程技術(shù)達到國際先進水平的關(guān)鍵.在此相關(guān)研究領(lǐng)域，國際上，Arabelos和Tscherning(1998)提出利用地面重力對衛(wèi)星重力梯度觀測值進行檢校的方法之后，Wolf和Denker(2005)利用地面重力驗證了衛(wèi)星軌道高度處的衛(wèi)星重力梯度張量觀測值，并將各分量上重力梯度延拓值與觀測值進行了差異比對.此外，Eshagh(2010, 2011)研究提出利用地面重力以及大地水準面高計算衛(wèi)星重力梯度的方法，并在此基礎(chǔ)上對衛(wèi)星重力梯度張量觀測值進行檢校和驗證.prlák等(2015)在Eshagh研究的基礎(chǔ)上，提出利用衛(wèi)星測高數(shù)據(jù)計算衛(wèi)星重力梯度的方法.Yildiz等(2016)在南極地區(qū)利用航空重力數(shù)據(jù)，對GOCE重力梯度張量實測數(shù)據(jù)進行了驗證.Bouman和Fuchs(2012)利用先驗重力場模型對測量帶寬內(nèi)的重力梯度值進行了驗證.Rispens和Bouman(2011)聯(lián)合先驗重力場模型與星敏感器數(shù)據(jù)對加速度計測量值進行了檢校.國內(nèi)學者針對GOCE衛(wèi)星的在軌檢校工作進行了跟蹤研究.羅志才等(2009)在開展GOCE衛(wèi)星重力梯度測量數(shù)據(jù)預(yù)處理方法研究的基礎(chǔ)上，討論了利用重力場模型、GOCE SST數(shù)據(jù)以及地面重力數(shù)據(jù)對GOCE衛(wèi)星在軌檢校的技術(shù)方法.吳云龍(2010)利用地面重力對GOCE衛(wèi)星的在軌檢校方法進行了進一步的研究.徐新禹等(2015)對同階次的不同重力場模型、同一模型不同階次以及是否考慮漂移參數(shù)對GOCE重力梯度觀測值檢校的影響進行了系統(tǒng)分析.這些基于GOCE衛(wèi)星的在軌檢校技術(shù)方法為我國自主梯度模式重力衛(wèi)星的檢校提供了良好的技術(shù)基礎(chǔ).

總的來說，國內(nèi)外學者對利用地面實測重力或地球重力場模型實現(xiàn)梯度儀在軌檢校進行了比較廣泛和深入的研究，為GOCE衛(wèi)星的在軌檢校做出了重要貢獻，但在地面重力數(shù)據(jù)的空間尺度、精度以及梯度儀本身噪聲對在軌檢校影響方面的綜合分析研究存在明顯不足，在提高檢校效果、降低研制成本方面尚存在一定的發(fā)展空間.本文吸收國際國內(nèi)對衛(wèi)星梯度儀在軌檢校的最新技術(shù)成果，在借鑒GOCE衛(wèi)星在軌檢校方法的基礎(chǔ)上，面向中國發(fā)展的梯度模式重力衛(wèi)星的規(guī)劃任務(wù)，研究提出一種綜合考慮地面重力數(shù)據(jù)空間尺度、精度并顧及梯度儀噪聲的衛(wèi)星重力梯度在軌檢校方法，通過球諧分析與球諧綜合將地面數(shù)據(jù)延拓至衛(wèi)星軌道高度，利用旋轉(zhuǎn)矩陣建立局部指北坐標系與衛(wèi)星梯度儀坐標系的轉(zhuǎn)換關(guān)系，成功實現(xiàn)地面數(shù)據(jù)與衛(wèi)星觀測數(shù)據(jù)的比對檢校.

1 研究方法

1.1 地面重力向上延拓

地球質(zhì)體外部引力位場是一個調(diào)和場，任意一點的引力位V滿足Laplace方程.利用分離變量法，推導球近似條件下引力位的球諧級數(shù)展開式為：

(1)

正常引力位V′的球諧級數(shù)展開式為：

(2)

T=W-U=V-V′.

(3)

將式(1)、式(2)代入式(3)，得到擾動位T的球諧級數(shù)展開形式如下：

(4)

由地面重力異常Δg與擾動位的關(guān)系：

(5)

推導得到地面重力異常的球諧函數(shù)展開表達式為：

(6)

根據(jù)球諧函數(shù)的正交性，由地面重力異常計算擾動位系數(shù)公式如下：

(7)

擾動位系數(shù)加上正常位系數(shù)得到引力位系數(shù)后，利用球諧綜合將地面重力向上延拓得到衛(wèi)星高度處的重力梯度張量.局部指北坐標系(local north-oriented frame, LNOF)下，重力梯度張量Vzz分量可表示為如下形式：

重力梯度張量其他分量的表達形式可參見Petrovskaya和Vershkov(2006)，由此即可實現(xiàn)與衛(wèi)星梯度張量觀測值的比對，并基于此開展衛(wèi)星梯度儀的在軌檢校.

1.2 檢校模型構(gòu)建

衛(wèi)星重力梯度觀測值是在梯度儀坐標系(gradiometer reference frame, GRF)下給出，而基于地面重力的梯度延拓值是在LNOF下，因此需要建立LNOF與GRF之間的轉(zhuǎn)換關(guān)系，可由旋轉(zhuǎn)矩陣表示為：

VGRF=RVGRFRT，

(9)

圖1 GOCE衛(wèi)星參考系統(tǒng)：OG-XGYGZG表示梯度儀坐標系，OS-XNYNZN表示局部指北坐標系，O-XIYIZI表示慣性坐標系，O-XEYEZE表示地固坐標系Fig.1 The reference systems for GOCE satellite. OG-XGYGZG is the GRF coordinate system, OS-XNYNZN is the LNOF coordinate system, O-XIYIZI is the IRF coordinate system, O-XEYEZE is the EFRF coordinate system

重力梯度觀測值的在軌檢校是利用外部獨立數(shù)據(jù)計算的重力梯度值與衛(wèi)星重力梯度觀測值進行比較，建立函數(shù)關(guān)系式，進而求得檢校參數(shù)，即：

(10)

v=Ax-l，

(11)

利用最小二乘估計尺度因子和偏差參數(shù)，可得：

x=(ATPA)-1ATPl.

(12)

2 實驗分析

根據(jù)衛(wèi)星設(shè)計指標，對地面觀測數(shù)據(jù)精度、空間尺度以及重力梯度儀觀測噪聲等關(guān)鍵檢校要素開展實驗分析.沿軌重力梯度張量采用EGM2008重力場模型，結(jié)合GOCE衛(wèi)星精密軌道數(shù)據(jù)SST_PSO_2(precise science orbits)進行計算，數(shù)據(jù)時間從2011年2月11日至6月23日，采樣率為1 s.

2.1 地面重力向上延拓與檢校評估

為評估地面重力向上延拓精度，采用EGM2008重力場模型計算格網(wǎng)分辨率為0.25°×0.25°的全球重力異常，分別加入標準差為5 μGal、1 mGal和2 mGal的白噪聲，作為地面重力起算數(shù)據(jù)，將其向上延拓得到衛(wèi)星軌道處的重力梯度張量，與沿軌直接球諧綜合計算的EGM2008重力梯度張量進行比較，得到重力梯度各分量的延拓誤差，如表1所示.

表1 重力梯度張量各分量的延拓誤差(mE)Table 1 Standard deviations of upward continuation errors of each component of gravitational gradients (mE)

從表1可以看出，基于地面重力數(shù)據(jù)得到的重力梯度張量各分量的向上延拓誤差隨著地面觀測數(shù)據(jù)精度的降低而變大.當?shù)孛嬷亓?shù)據(jù)精度為5 μGal時，各分量延拓誤差均在10-2mE水平；當?shù)孛嬷亓?shù)據(jù)精度為1 mGal時，各分量延拓誤差在0.3～0.6 mE；當?shù)孛嬷亓?shù)據(jù)精度達到2 mGal時，各分量延拓誤差在0.4～1.2 mE.參考GOCE衛(wèi)星重力梯度張量四個高精度分量中，Vxx、Vyy分量在軌實際測量精度為10 mE/Hz1/2，Vzz、Vxz分量測量精度在20 mE/Hz1/2左右(Floberghagen et al., 2011; Rummel et al., 2011)，基于地面重力的重力梯度向上延拓誤差遠小于GOCE衛(wèi)星重力梯度觀測值的噪聲水平.圖2給出了地面重力精度分別為5 μGal、1 mGal和2 mGal時，重力梯度張量對角線分量的向上延拓誤差分布情況.

由圖2可以看出，重力梯度張量對角線分量的向上延拓誤差隨機分布在衛(wèi)星軌道上，即地面重力觀測數(shù)據(jù)的地理分布對向上延拓誤差產(chǎn)生的影響較小.當?shù)孛嬷亓τ^測精度為5 μGal時，Vxx、Vyy兩個分量的向上延拓誤差較為接近，均在±0.2 mE以內(nèi)，而Vzz分量則達到±0.4 mE；當?shù)孛嬷亓τ^測精度為1 mGal時，Vxx、Vyy兩個分量延拓誤差范圍達到±3 mE，Vzz分量增至±5 mE，相比5 μGal精度的地面重力，增幅較為明顯；當?shù)孛嬷亓τ^測精度在2 mGal時，Vxx、Vyy兩個分量延拓誤差范圍為±4 mE，Vzz分量為±7 mE.由此看出，Vzz分量向上延拓誤差大小約為Vxx、Vyy兩個分量之和，但符號相反.

圖2 重力梯度張量對角線分量的向上延拓誤差分布Fig.2 The distribution of upward continuation errors in diagonal components of gravitational gradients

利用地面重力檢校衛(wèi)星重力梯度觀測值，往往在0.005～0.1 Hz的測量帶寬內(nèi)進行(Veicherts et al., 2011).由于濾波后的重力梯度觀測值中偏差參數(shù)很小，可以忽略不計，待檢校參數(shù)一般只考慮尺度因子.因此，本文在重力梯度觀測中加入尺度因子為1.05的待檢校參數(shù)，并采用巴特沃斯帶通濾波器對地面重力延拓值與沿軌重力梯度球諧綜合值進行同樣的濾波處理.以一周為檢校周期(Bouman et al., 2011)，求定的尺度因子結(jié)果如表2所示.

表2 重力梯度觀測值尺度因子檢校中誤差Table 2 The root mean square errors of the scale factors of gravitational gradient observations

從表2可以看出，重力梯度觀測值檢校參數(shù)的確定精度與地面重力觀測精度相關(guān).其中，當?shù)孛嬷亓τ^測精度在5 μGal水平時，檢校得到重力梯度張量尺度因子的中誤差在10-5～10-6量級；地面重力觀測精度達到1 mGal時，尺度因子的中誤差在10-4～10-6量級；地面重力觀測精度為2 mGal時，尺度因子的中誤差在10-4量級.檢校后的重力梯度觀測值應(yīng)滿足Laplace方程，重力梯度張量矩陣對角線分量之和，即重力梯度張量的跡理論值為零.當重力梯度觀測值中存在誤差，重力梯度張量的跡實際值通常不為零，此時跡越小，檢校精度越高，以此對重力梯度觀測值檢校精度進行評價.利用不同精度的地面重力數(shù)據(jù)檢校后，重力梯度張量的跡如圖3所示.

由圖3可以看出，利用5 μGal精度的地面重力對重力梯度觀測值進行檢校后，重力梯度張量的跡優(yōu)于10-3mE；地面重力精度為1 mGal時，檢校后的重力梯度張量的跡仍優(yōu)于10-2mE；地面觀測精度降至2 mGal水平時，檢校后的重力梯度張量的跡偏大，在10-2～10-1mE量級.這說明，隨著地面重力觀測數(shù)據(jù)精度水平的提升，重力梯度觀測值的檢校精度會相應(yīng)升高1個數(shù)量級左右. 目前地面重力數(shù)據(jù)獲取精度在5 μGal，而航空重力數(shù)據(jù)獲取精度也能達到1 mGal左右(Forsberg and Olesen, 2010)，因此衛(wèi)星重力梯度觀測值在軌檢校時，可采用地面重力或航空重力，其數(shù)據(jù)精度能夠滿足衛(wèi)星重力梯度儀檢校要求.

圖3 地面重力數(shù)據(jù)檢校后重力梯度張量的跡Fig.3 The gravitational gradient trace calibrated by ground gravity data

2.2 重力梯度儀噪聲影響

衛(wèi)星重力梯度觀測值在測量帶寬內(nèi)呈現(xiàn)白噪聲特性，不同噪聲量級會對重力梯度儀檢校產(chǎn)生較大影響.參考GOCE衛(wèi)星重力梯度儀發(fā)射前的測量精度指標為6～7 mE，實際在軌測量精度在10～20 mE水平(Marque et al., 2010)，本文在設(shè)計仿真實驗時，在沿軌重力梯度球諧綜合值中加入1.05的尺度因子，對其進行帶通濾波處理后，加入標準差分別為6 mE、10 mE以及20 mE的白噪聲，作為待檢校的梯度觀測值.利用5 μGal和1 mGal精度的地面重力對梯度觀測值進行檢校，得到的尺度因子中誤差結(jié)果如表3所示.

表3 不同噪聲水平的重力梯度觀測值尺度因子檢校中誤差Table 3 The root mean square errors of the scale factors of gravitational gradient observations with different noise

從表3可看出，測量帶寬內(nèi)的重力梯度觀測值在不同噪聲量級影響下，利用5 μGal或者1 mGal精度的地面重力，計算得到六個分量上尺度因子的中誤差均在10-4～10-6水平.這表明，利用高精度的地面重力能夠高精度的求定不同觀測噪聲水平的重力梯度觀測值中的檢校參數(shù).以重力梯度張量的跡作為檢校精度評價指標，利用5 μGal精度的地面重力對不同噪聲水平下的梯度儀觀測值進行檢校后，重力梯度張量的跡如圖4所示.

圖4 不同噪聲水平的梯度儀觀測值檢校后重力梯度張量的跡Fig.4 The calibrated gravitational gradient trace with different noise of gradiometer observations

由圖4可以看出，檢校后重力梯度張量的跡隨著梯度儀觀測值噪聲增加而變大，不同噪聲水平下，檢校后重力梯度張量的跡分別為7.0 mE、9.0 mE和12.8 mE左右.結(jié)合表3和圖4可以得到，基于地面重力的重力梯度觀測值在軌檢?？梢杂行跤^測值中的系統(tǒng)偏差，但觀測值中的偶然誤差仍會對檢校后重力梯度張量的跡產(chǎn)生較大的影響.GOCE衛(wèi)星發(fā)射前，重力梯度張量的跡精度指標要求為11 mE(Siemes et al., 2012; Stummer et al., 2012)，考慮到未來國產(chǎn)重力梯度儀的精度要優(yōu)于該指標，因此其測量噪聲控制在10 mE以內(nèi)較為合適.

2.3 地面重力空間尺度影響

利用地面重力對衛(wèi)星重力梯度觀測值進行在軌檢校時，由于觀測成本限制，地面重力實地觀測空間尺度是有限的，可采用中心區(qū)域?qū)崪y地面重力數(shù)據(jù)、外部區(qū)域先驗重力場模型值填充的方式進行.

考慮到EGM2008重力異常精度約為10 mGal(章傳銀等, 2009; Pavlis et al., 2012)，本文研究中，在外部區(qū)域的EGM2008模型重力異常中加入標準差為10 mGal白噪聲作為填充值，在中心區(qū)域的EGM2008模型重力異常中分別加入標準差為5 μGal和1 mGal的白噪聲，作為實測地面重力，進行向上延拓后，將延拓的重力梯度值與沿軌直接球諧綜合計算的EGM2008重力梯度值進行比較.不同空間尺度下，重力梯度各分量向上延拓誤差如圖5所示.

由圖5可以看出，檢校中心區(qū)域采用5 μGal或者1 mGal精度的地面重力向上延拓時，重力梯度各分量的延拓誤差隨地面重力空間尺度的增加有逐漸變小的趨勢.相同空間尺度下，檢校中心區(qū)域采用5 μGal精度的地面重力向上延拓得到的重力梯度張量延拓誤差要小于1 mGal精度的地面重力；對于重力梯度分量而言，Vzz、Vxz及Vyz三個分量延拓誤差相較其他幾個分量略大.由于GOCE重力衛(wèi)星測量帶寬下限為0.005 Hz，相當于衛(wèi)星飛行200 s，因此在利用地面重力檢校梯度觀測值時，要求檢校中心區(qū)域至少為12°×12°.在該尺度條件下，中心區(qū)域采用5 μGal精度的地面重力，重力梯度六個分量的延拓誤差在0.2～0.4 mE；中心區(qū)域采用1 m Gal精度的地面重力，重力梯度六個分量的延拓誤差在0.2～0.9 mE.

圖5 不同空間尺度下地面重力向上延拓誤差(a) 5 μGal地面重力； (b) 1 mGal地面重力.Fig.5 The upward continuation errors for different calibration areas with the precision of ground gravity data in 5 μGal (a) and 1 mGal (b)

同樣，利用局部地面重力對重力梯度觀測值進行在軌檢校需在測量帶寬內(nèi)進行，研究中將沿軌重力梯度球諧綜合值作為真值，加入1.05的尺度因子后對其進行帶通濾波，并加入標準差為10 mE的白噪聲，作為待檢校的梯度觀測值.當?shù)孛嬷亓τ^測精度為5 μGal，空間尺度為12°×12°時，利用地面重力檢校得到衛(wèi)星重力梯度觀測值六個分量上尺度因子的統(tǒng)計結(jié)果見表4.

表4 地面觀測空間尺度為12°×12°時檢校得到的尺度因子統(tǒng)計值Table 4 The statistical values of the scale factors when the ground gravity data in a calibration area of 12°×12°

從表4可以看出，檢校得到的重力梯度觀測值尺度因子的平均值與1.05之間最大差異為0.002，尺度因子中誤差在10-3～10-4量級.這說明，利用12°×12°的地面重力觀測數(shù)據(jù)，能夠以較高精度確定出檢校參數(shù).此外，將檢校前后重力梯度張量的跡與梯度儀觀測值只含有白噪聲時的重力梯度張量的跡進行了對比，并計算了檢校前后重力梯度Vzz分量的變化，如圖6所示.

圖6 重力梯度張量的跡(a)以及Vzz分量(b)檢校前后的變化Fig.6 The differences before and after calibration of gravitational gradient trace (a) and Vzz component (b) of gravitational gradients

由圖6可看出，檢校前重力梯度張量的跡平均值在9.5 mE水平，在利用地面重力數(shù)據(jù)對梯度觀測值進行檢校后，重力梯度張量的跡平均值約為9.0 mE水平，與梯度觀測值只含有白噪聲條件下的梯度張量的跡水平一致，這說明重力梯度觀測值中白噪聲的存在是導致檢校后重力梯度張量的跡不為零的主要因素；對于重力梯度Vzz分量而言，在檢校前觀測值與真值之間差異較為明顯，其值在19.3 mE，這是由于梯度觀測值中的尺度因子和噪聲共同影響的結(jié)果，而檢校后觀測值與真值的差異在4.4 mE，相比檢校前有了明顯提升.這表明利用本文設(shè)計提出的基于地面重力的衛(wèi)星重力梯度在軌檢校方法可有效實現(xiàn)未來我國自主重力衛(wèi)星梯度儀的高精度在軌檢校.

3 結(jié)論

作為重力梯度衛(wèi)星任務(wù)的重要技術(shù)環(huán)節(jié)，重力梯度儀的在軌檢校和觀測值的在軌驗證可有效改善重力梯度衛(wèi)星數(shù)據(jù)預(yù)處理效果，提高衛(wèi)星觀測質(zhì)量.本文面向我國未來重力梯度衛(wèi)星發(fā)展規(guī)劃，研究了基于地面重力的衛(wèi)星重力梯度觀測值的在軌檢校方法，從地面觀測數(shù)據(jù)精度、空間尺度以及重力梯度儀觀測噪聲等檢校要素進行了分析研究以及實驗驗證.結(jié)果表明，在12°×12°的地面檢校區(qū)域內(nèi)，利用精度優(yōu)于1 mGal的地面重力觀測數(shù)據(jù)即可實現(xiàn)觀測噪聲低于10 mE的重力梯度儀高精度在軌檢校.在此條件下，重力梯度張量各分量延拓誤差相比梯度儀觀測噪聲低一個量級，并且能夠以10-3～10-4水平確定檢校參數(shù)，檢校后重力梯度張量的跡在8.7～9.4 mE水平，可有效實現(xiàn)我國國產(chǎn)衛(wèi)星重力梯度儀的高精度在軌檢校.

致謝感謝ESA(European Space Agency)提供的GOCE衛(wèi)星數(shù)據(jù)和ICGEM(International Centre for Global Earth Models)提供的EGM2008全球重力場模型.