利用Lunar Prospector Doppler數(shù)據(jù)解算并分析月球重力場(chǎng)

黃昆學(xué)，常曉濤，孫　玉

(1. 武漢大學(xué)測(cè)繪學(xué)院，湖北 武漢 430079； 2. 國(guó)家測(cè)繪地理信息局衛(wèi)星測(cè)繪應(yīng)用中心，北京 100830；

3. 代爾夫特理工大學(xué)，荷蘭 代爾夫特)

Using Lunar Prospector Doppler Data to Recover and Analyze

Lunar Gravity Field

HUANG Kunxue，CHANG Xiaotao，SUN Yu

?

利用Lunar Prospector Doppler數(shù)據(jù)解算并分析月球重力場(chǎng)

黃昆學(xué)1，常曉濤2，孫玉3

(1. 武漢大學(xué)測(cè)繪學(xué)院，湖北 武漢 430079； 2. 國(guó)家測(cè)繪地理信息局衛(wèi)星測(cè)繪應(yīng)用中心，北京 100830；

3. 代爾夫特理工大學(xué)，荷蘭 代爾夫特)

Using Lunar Prospector Doppler Data to Recover and Analyze

Lunar Gravity Field

HUANG Kunxue，CHANG Xiaotao，SUN Yu

摘要：Doppler技術(shù)可高精度測(cè)定衛(wèi)星之間的距離變化率，可以和其他技術(shù)聯(lián)合用于月球衛(wèi)星定軌。本文利用Lunar Prospector月球項(xiàng)目4個(gè)月的Doppler跟蹤數(shù)據(jù)，解算了75階次月球重力場(chǎng)模型。首先介紹了數(shù)據(jù)預(yù)處理方法、力模型的使用及參數(shù)設(shè)置等關(guān)鍵技術(shù)，然后把解算的重力場(chǎng)模型與同階次的其他月球重力場(chǎng)模型進(jìn)行比較，可以看出低階位系數(shù)較LP75G月球重力場(chǎng)模型稍差，但高階位系數(shù)精度有所提高，可能是僅利用了低軌單星數(shù)據(jù)的原因。

關(guān)鍵詞：Lunar Prospector；月球重力場(chǎng)；Doppler；數(shù)據(jù)解算

Lunar Prospector(LP)是美國(guó)NASA的月球探索任務(wù)，衛(wèi)星于1998年1月6日發(fā)射升空，經(jīng)過一系列調(diào)整之后進(jìn)入月球極地近圓軌道。LP繞月軌道被設(shè)計(jì)為過兩極的圓軌道,軌道高度設(shè)計(jì)為100 km,運(yùn)行速度約為5517 km/h,運(yùn)行周期為118 min。LP的設(shè)計(jì)初衷是用于探測(cè)、研究月球全球構(gòu)成等，重點(diǎn)是探測(cè)南極水冰和分析月球的重力場(chǎng)及磁場(chǎng)。為完成不同的科學(xué)任務(wù),LP任務(wù)分成3個(gè)階段。正常任務(wù)階段從衛(wèi)星進(jìn)入月球軌道至1999年1月，歷時(shí)1年。其軌道為極地圓軌道，平均高度為100 km[1]。之后，LP進(jìn)入短暫的過渡階段，平均軌道高度為40 km。過渡階段的主要任務(wù)是驗(yàn)證軌道預(yù)測(cè)的能力，為其后的擴(kuò)展任務(wù)階段做準(zhǔn)備。LP擴(kuò)展任務(wù)階段開始于1999年1月29日，直到7月31日撞擊月球表面為止，平均軌道高度約為30 km。采用較低的環(huán)月軌道是為了得到高精度、高分辨率的月球重力場(chǎng)模型。

美國(guó)JPL已經(jīng)利用LP數(shù)據(jù)研制了一系列高精度、高階次月球重力場(chǎng)模型，包括LP75D、LP75G、LP100J、LP100K、LP150和LP165等[2-3]。此外，也有學(xué)者利用GEODYN和SOLVE軟件處理了公布在PDS上的LP Doppler跟蹤數(shù)據(jù)，從而實(shí)踐了月球重力場(chǎng)解算技術(shù)流程。文獻(xiàn)[4]利用正常任務(wù)階段數(shù)據(jù)模制了75階月球重力場(chǎng)模型。也有學(xué)者利用擴(kuò)展任務(wù)階段數(shù)據(jù)，以LP100J為初始重力場(chǎng)模型解算了100階月球重力場(chǎng)模型。但上述文獻(xiàn)主要介紹了基本解算原理和方法，均未涉及LP數(shù)據(jù)預(yù)處理等具體信息。本文詳盡介紹了數(shù)據(jù)預(yù)處理方法、力模型的使用及參數(shù)設(shè)置等關(guān)鍵技術(shù)，然后對(duì)解算的重力場(chǎng)模型進(jìn)行了比較分析。

一、數(shù)據(jù)處理方法和流程

本文所采用的數(shù)據(jù)來自美國(guó)PDS的LP軌道數(shù)據(jù)文件(ODF)。該數(shù)據(jù)文件包含美國(guó)深空網(wǎng)(DSN)對(duì)LP的Doppler跟蹤觀測(cè)數(shù)據(jù)，包括距離變率(測(cè)速)和距離(測(cè)距)兩種數(shù)據(jù)。距離變率數(shù)據(jù)采樣間隔為10 s，精度較高(約為0.3 mm/s)；距離數(shù)據(jù)采樣間隔為33 s，精度平均約為0.5 m。由于最終要使用GEODYN Ⅱ軟件處理觀測(cè)數(shù)據(jù)，因此需要將ODF中的數(shù)據(jù)記錄轉(zhuǎn)換為GEODYN Ⅱ可識(shí)別的格式。首先根據(jù)數(shù)據(jù)格式說明文件TRK-2-18將ODF由二進(jìn)制文件轉(zhuǎn)為ASCII文件。ODF中的距離變率數(shù)據(jù)以Hz為單位,距離數(shù)據(jù)則以range unit(RU)為單位，為了使其符合GEODYN Ⅱ?qū)τ跀?shù)據(jù)格式的特殊要求，需要按下式將二者單位分別轉(zhuǎn)換為m/s和m

Doprate=0.5×C/(T×Fref)×Dopfreq

(1)

Doprange=0.5×C/(FRU×Fref)×Doprange

(2)

式中，C表示光速；T為轉(zhuǎn)換器上下行信號(hào)比，S波段應(yīng)取為240/221；Fref為基準(zhǔn)頻率，可由數(shù)據(jù)文件中讀取，其值接近于2093 MHz；FRU取值為0.5。由于ODF實(shí)際上為Ramp格式數(shù)據(jù)，其基準(zhǔn)頻率通常是隨時(shí)間線性變化的。但在實(shí)際計(jì)算過程中，注意到絕大部分情況下基準(zhǔn)頻率的變率實(shí)際上為零，因此可以將其作為普通的Doppler數(shù)據(jù)來處理，即基準(zhǔn)頻率為常量。

本文研究所采用的衛(wèi)星定軌及大地參數(shù)解算軟件為GEODYN Ⅱ和SOLVE。GEODYN Ⅱ軟件主要由3部分構(gòu)成，即TDF、IIS和IIE。TDF是格式轉(zhuǎn)換器，能夠?qū)⒛撤N數(shù)據(jù)記錄轉(zhuǎn)換為GEODYN II可處理的格式。但由ODF轉(zhuǎn)換得到的ASCII數(shù)據(jù)記錄格式仍然不能被TDF識(shí)別，因此需要利用程序再次進(jìn)行轉(zhuǎn)換。經(jīng)過TDF轉(zhuǎn)化后，觀測(cè)數(shù)據(jù)文件的格式變?yōu)镚EODYN II binary格式，該格式可以被GEODYN II所識(shí)別。根據(jù)情況進(jìn)行弧段劃分(Arc1，Arc2，…，Arcn)。由于月球背面存在重力信息空白，軌道弧段不宜過長(zhǎng)，一般取為2 d。還要注意同一弧段中不應(yīng)包括對(duì)衛(wèi)星的人工操作，即衛(wèi)星機(jī)動(dòng)的時(shí)間段應(yīng)排除在連續(xù)弧段之外。上述數(shù)據(jù)預(yù)處理過程如圖1所示。SOLVE嚴(yán)格意義上應(yīng)該屬于GEODYN Ⅱ的一部分，其主要功能是進(jìn)行矩陣求逆運(yùn)算。圖1中的ReadODF與Trans2TDF為數(shù)據(jù)格式轉(zhuǎn)換程序。

圖1　LP 數(shù)據(jù)預(yù)處理流程

二、月球重力場(chǎng)模型解算

本文解算使用1998年1月11日至5月1日正常任務(wù)階段約4個(gè)月的Doppler跟蹤雙程、三程數(shù)據(jù)，所使用的LP數(shù)據(jù)稍多于LP75G(3個(gè)月)。解算流程如圖2所示，圖2中Ematrix表示法方程矩陣。所采用的力學(xué)模型、星歷及參數(shù)設(shè)置等見表1。用于解算GLGM2的數(shù)據(jù)主要來自軌道高度約為200km的Clementine任務(wù)。以該模型作為初始重力場(chǎng)模型對(duì)LP進(jìn)行軌道積分時(shí)，顯然精度不夠，即使積分弧段為2d仍然會(huì)導(dǎo)致部分弧段不收斂。對(duì)于不收斂的弧段，先縮短其積分長(zhǎng)度直至收斂[5-6]。利用所有積分弧段解算重力場(chǎng)模型，將更新后的重力場(chǎng)模型作為先驗(yàn)重力場(chǎng)模型替換GLGM2[7-12]。更新后的重力場(chǎng)模型包含LP探測(cè)的重力場(chǎng)信息，可使所有2d積分弧段均收斂。

圖2　利用GEODYN II和SOLVE解算月球重力場(chǎng)位系數(shù)的具體流程及所需文件

ftn01行星星歷文件ftp:∥navigator.jpl.nasa.govJPLftn02集成了地磁、太陽(yáng)通量等文件ftp:∥ftp.ngdc.noaa.govNOAAftn12重力場(chǎng)模型GLGM2(70階)PDSLemoineetal.1997ftn05控制文件太陽(yáng)光壓模型:簡(jiǎn)單柱狀模型多體攝動(dòng)模型:太陽(yáng)及其他主要模型精密初始軌道:利用spicetoolkit軟件提取并可同時(shí)進(jìn)行坐標(biāo)框架轉(zhuǎn)換 PDSNAIFDSN測(cè)站坐標(biāo):DSN16、24、27、54、61、66、34、42、46Folkner1997測(cè)站速度模型:ITRF93框架Folkner1997測(cè)站坐標(biāo)改正:固體潮、海潮、極潮ftp:∥geodesy.gsfc.nasa.govLOVE:h2=0.609;k2=0.0852

續(xù)表1

首先，根據(jù)所選初始力模型(月球重力場(chǎng)模型、三體攝動(dòng)模型、太陽(yáng)光壓模型等)對(duì)LP軌道進(jìn)行數(shù)值積分；然后便可得到線性觀測(cè)方程，該觀測(cè)方程建立了觀測(cè)數(shù)據(jù)(距離、距離變率數(shù)據(jù))與重力場(chǎng)位系數(shù)、初始狀態(tài)向量等未知參數(shù)之間的函數(shù)關(guān)系，并可形成法方程矩陣(Ematrix)；最終利用SOLVE對(duì)該法方程矩陣求逆，便可得到待求參數(shù)[13-14]。

三、模型分析

LPD75的模型階方差曲線與LP75G比較接近(如圖3所示)，但在低階次(7～20階)出現(xiàn)了明顯誤差，這在其誤差階方差曲線中也體現(xiàn)了出來。主要原因是由于GLGM2及LP75G都融合了多個(gè)歷史探月任務(wù)的數(shù)據(jù)，而本文解算的LPD75僅利用了LP的Doppler跟蹤數(shù)據(jù)，所恢復(fù)的月球重力場(chǎng)模型低階位系數(shù)精度較差。圖4顯示了LPD75和LP75G模型與月球地形模型STM359_grid-03的相關(guān)性。LPD75與地形的相關(guān)性在前50階次與LP75G相比總體上略差，50階次以上則優(yōu)于LP75G。顯然，軌道較低(正常任務(wù)階段月100 km)的LP任務(wù)對(duì)月球重力場(chǎng)短波信息較為敏感，所增加的一個(gè)月的LP數(shù)據(jù)，有助于提高模型高階位系數(shù)的精度，但對(duì)50階以下的重力場(chǎng)位系數(shù)沒有改進(jìn)。

圖3　模型階方差及誤差階方差

由LPD75解算的月球自由空氣重力異常 (如圖5(a)所示)、月球大地水準(zhǔn)面起伏(如圖6(a)所示)與LP75G(如圖5(b)、圖6(b)所示)可以看出，兩個(gè)模型的差別主要集中在月球背面和極區(qū)。LPD75的高頻噪聲明顯少于LP75G。

圖4　重力場(chǎng)模型與地形模型相關(guān)性

圖5　月表自由空氣重力異常及其差值

圖6　月球大地水準(zhǔn)面起伏及其差值

四、結(jié)束語(yǔ)

本文基于LP Doppler數(shù)據(jù)，利用GEODYN Ⅱ和SOVLE軟件解算了月球重力場(chǎng)模型，得到了合理可靠的結(jié)果，驗(yàn)證了數(shù)據(jù)處理、軟件使用的正確性。文中重點(diǎn)說明了其他學(xué)者忽略的LP數(shù)據(jù)預(yù)處理方法和過程，詳細(xì)闡述了月球重力場(chǎng)解算的流程，給出了必要文件及參數(shù)設(shè)置信息。除LP數(shù)據(jù)預(yù)處理過程外，文中所述GEODYN Ⅱ與SOLVE的解算過程均為通用過程，對(duì)于其他任務(wù)和類型的數(shù)據(jù)，只需要進(jìn)對(duì)預(yù)處理過程和控制文件(ftn05)作相應(yīng)改動(dòng)即可，可對(duì)今后利用GEODYN Ⅱ和SOLVE軟件進(jìn)行月球重力場(chǎng)模型解算提供有借鑒。

參考文獻(xiàn)：

[1]BINDER A B. Lunar Prospector: Overview [J]. Science, 1998, 281(5382): 1475-1476.

[2]KONOPLIV A S, BANERDT W B, SJOGREN W L. Venus Gravity: 180th Degree and Order Model [J]. Icarus, 1999, 139(1): 3-18.

[3]KONOPLIV A S, ASMAR S W, CARRANZA E, et al. Recent Gravity Models as a Result of the Lunar Prospector Mission [J]. Icarus, 2001, 150(1): 1-18.

[4]GOOSSENS S, MATSUMOTO K. Lunar Satellite Orbit Determination Analysis and Quality Assessment from Lunar Prospector Tracking Data and Selene Simulations [J]. Advances in Space Research, 2007, 40(1): 43-50.

[5]豐海，李建成，李大煒，等.月球Airy均衡狀態(tài)與月殼厚度估計(jì)[J]. 測(cè)繪學(xué)報(bào),2012 ,41 (4): 543-548.

[6]PING J S, HUANG Q, YAN J G, et al. Lunar Topographic Model Cltm-s01 from Chang’e-1 Laser ALtimeter [J]. Science in China Series G: Physics Mechanics and Astronomy, 2009, 52(7): 1105-1114.

[7]羅志才,李瓊,鐘波. 利用GRACE時(shí)變重力場(chǎng)反演黑河流域水儲(chǔ)量變化[J]. 測(cè)繪學(xué)報(bào)，2012,41 (5): 676-681.

[8]GUO J Y, SUN Y, CHANG X T,et al. Lunar Deflections of the Vertical and Their Distribution [J]. Natural Science, 2011, 3(5): 339-343.

[9]ARAKI H, TAZAWA S, NODA H, et al. Lunar Global Shape and Polar Topography Derived from Kaguya-LALT Laser Altimetry [J]. Science, 2009, 323(5916): 897-900.

[10]李森，安智明. 克里金模型在區(qū)域重力數(shù)據(jù)插值中的應(yīng)用[J]. 測(cè)繪通報(bào), 2013 (10): 63-66.

[11]MAZARICO E, LEMOINE F, HAN S C et al. GLGM-3: A Degree-150 Lunar Gravity Model from the Historical Tracking Data of NASA Moon Orbiters [J]. Journal of Geophysical Research, 2010(115): E05001.

[12]FOLKNER W M, WILLIAMS J G, BOGGS D H. The Planetary and Lunar Ephemeris de 421 [J]. IPN Progress Report, 2008:172-178.

[13]曹金國(guó)，王建斌，戴山嶺. 重力場(chǎng)擾動(dòng)位泰勒展開式低階項(xiàng)應(yīng)用研究[J].測(cè)繪通報(bào)，2014(1):54-55.

[14]GOOSSENS S, MATSUMOTO K. Lunar Satellite Orbit Determination Analysis and Quality Assessment from Lunar Prospector Tracking Data and Selene Simulations [J]. Advances in Space Research, 2007, 40(1): 43-50.

引文格式： 黃昆學(xué)，常曉濤，孫玉. 利用Lunar Prospector Doppler數(shù)據(jù)解算并分析月球重力場(chǎng)[J].測(cè)繪通報(bào)，2016(1)：19-22.DOI:10.13474/j.cnki.11-2246.2016.0005.

作者簡(jiǎn)介：黃昆學(xué)(1985—)，男，博士生，主要研究方向?yàn)樵虑蛑亓?chǎng)解算。E-mail：shulang2003@163.com

基金項(xiàng)目：國(guó)家自然科學(xué)基金 (41204007); 國(guó)家高技術(shù)研究發(fā)展計(jì)劃(2013AA122502); 國(guó)家重點(diǎn)基礎(chǔ)研究發(fā)展計(jì)劃(2013CB733302); 中國(guó)地質(zhì)調(diào)查局地質(zhì)調(diào)查工作項(xiàng)目(12120113019100)

收稿日期：2014-11-25； 修回日期： 2015-10-16

中圖分類號(hào)：P228

文獻(xiàn)標(biāo)識(shí)碼：B

文章編號(hào)：0494-0911(2016)01-0019-04