深空探測(cè)器三程多普勒建模與算法實(shí)現(xiàn)

曹建峰，黃 勇，劉 磊，張 宇，鄭愛武，胡松杰

(1. 航天飛行動(dòng)力學(xué)技術(shù)重點(diǎn)實(shí)驗(yàn)室，北京100094；2. 北京航天飛行控制中心，北京100094；3. 中科院上海天文臺(tái)，上海 200030)

深空探測(cè)器三程多普勒建模與算法實(shí)現(xiàn)

曹建峰1,2,3，黃 勇3，劉 磊1,2，張 宇1,2，鄭愛武1,2，胡松杰1,2

(1. 航天飛行動(dòng)力學(xué)技術(shù)重點(diǎn)實(shí)驗(yàn)室，北京100094；2. 北京航天飛行控制中心，北京100094；3. 中科院上海天文臺(tái)，上海 200030)

針對(duì)高精度深空探測(cè)多普勒應(yīng)用，提出了一種間接使用中心天體、航天器的速度與加速度進(jìn)行多普勒建模的方法。避免了計(jì)算機(jī)字長(zhǎng)截?cái)嗾`差、星歷表插值誤差對(duì)建模精度的影響。使用該算法對(duì)環(huán)火星探測(cè)器進(jìn)行了連續(xù)1周的多普勒計(jì)算，通過常規(guī)的距離差分方法檢驗(yàn)了算法的正確性，對(duì)計(jì)算的1小時(shí)多普勒數(shù)據(jù)進(jìn)行10階多項(xiàng)式擬合，擬合后數(shù)據(jù)殘差小于0.002mm/s. 分析結(jié)果表明該算法計(jì)算多普勒結(jié)果正確，計(jì)算精度較常規(guī)距離差分算法在精度上得到1個(gè)量級(jí)的提升。

深空探測(cè)；多普勒；觀測(cè)建模

0 引 言

深空網(wǎng)(Deep space network, DSN)在深空探測(cè)任務(wù)中發(fā)揮著至關(guān)重要作用，測(cè)量類型主要包括多普勒測(cè)速與測(cè)距[1-2]。測(cè)距數(shù)據(jù)因?yàn)椴豢杀苊獾拇嬖谙到y(tǒng)性偏差，環(huán)繞飛行階段的航天器測(cè)軌工作更依賴于高精度的DSN多普勒測(cè)量。

美國(guó)推進(jìn)噴氣實(shí)驗(yàn)室(Jet propulsion laboratory, JPL)于1958年負(fù)責(zé)籌建了首個(gè)深空網(wǎng)，完成了早期的深空探測(cè)任務(wù)。近60年深空探測(cè)任務(wù)的不斷發(fā)展使得深空網(wǎng)測(cè)量能力與精度得到了大幅提升。早期Mariner 9、Viking 1-2的測(cè)軌數(shù)據(jù)主要為S波段雙程多普勒，60s積分雙程多普勒測(cè)量的噪聲水平約為1mm/s[3]。歐空局(European space agency, ESA)發(fā)射火星探測(cè)器——火星快車(Mars express, MEX)的測(cè)軌主要依賴于S和X波段的雙程、三程多普勒測(cè)量，其1s積分周期的S和X波段多普勒測(cè)量噪聲水平分別約為1.2mm/s和0.2mm/s[4]。日本臼田站在試驗(yàn)跟蹤中的多普勒精度優(yōu)于0.1mm/s[5]。嫦娥二號(hào)(Chang’E-2, CE-2)小行星探測(cè)試驗(yàn)后期的測(cè)軌完全依賴于DSN和中國(guó)VLBI網(wǎng)，受限于測(cè)控頻段以及星載轉(zhuǎn)發(fā)設(shè)備的能力，在小行星探測(cè)試驗(yàn)期間，CE-2在S頻段的測(cè)量精度約0.7mm/s(1 s積分周期)[6]，鄭為民在對(duì)火星快車跟蹤試驗(yàn)中獲取了0.17mm/s的多普勒測(cè)量精度[7]，這一狀況在未來的自主火星探測(cè)任務(wù)必將得到改善，測(cè)量精度有望達(dá)到ESA水準(zhǔn)。

使用多普勒測(cè)量數(shù)據(jù)對(duì)航天器開展精密定軌計(jì)算工作中，對(duì)多普勒測(cè)量進(jìn)行觀測(cè)建模是軌道計(jì)算的基礎(chǔ)，通常建模精度需高于測(cè)量數(shù)據(jù)精度1個(gè)量級(jí)以上時(shí)可以忽略建模引入的誤差。目前對(duì)于多普勒的建模主要使用距離差分建模，對(duì)于深空跟蹤數(shù)據(jù)在算法實(shí)現(xiàn)時(shí)必須采用四精度浮點(diǎn)數(shù)。哥達(dá)德分析中心指出雙精度浮點(diǎn)數(shù)引入的截?cái)嗾`差，大行星位置計(jì)算誤差都是制約多普勒建模精度的重要因素。對(duì)于積分周期較大的多普勒建模，如60積分周期，可以采用距離差分建模，文中稱為常規(guī)算法。對(duì)于較小的積分周期，必須在算法實(shí)現(xiàn)上進(jìn)行優(yōu)化處理。GEODYN軟件針對(duì)深空多普勒建模問題，在行星歷表位置計(jì)算中重新對(duì)插值引數(shù)進(jìn)行了計(jì)算[8]。歐空局在對(duì)多普勒觀測(cè)建模進(jìn)行了算法上的改進(jìn)，基于泰勒展開，但是仍然必須依賴四精度浮點(diǎn)數(shù)，其建模精度可以優(yōu)于0.01mm/s。論文對(duì)多普勒建模算法進(jìn)行了進(jìn)一步改進(jìn)，使用雙精度浮點(diǎn)數(shù)定義可以實(shí)現(xiàn)優(yōu)于0.01mm/s的建模精度。

1 三程多普勒測(cè)量模型

根據(jù)無線電信號(hào)上下行鏈路構(gòu)型的差異，地基測(cè)站對(duì)航天器的跟蹤模式可以分為單程、雙程和三程模式[9]。圖1給出了不同測(cè)量模式下的跟蹤示意圖。單程多普勒測(cè)量模式通過星載晶振產(chǎn)生下行信號(hào)，由地面臺(tái)站直接接收跟蹤信號(hào)，該模式只有下行鏈路，其測(cè)量精度受限于星載晶振的穩(wěn)定性，因此目前很少在實(shí)際深空探測(cè)任務(wù)中使用，我國(guó)“螢火一號(hào)”火星探測(cè)任務(wù)原計(jì)劃采用該測(cè)量模式[10]。雙程模式由地面臺(tái)站向探測(cè)器發(fā)射上行信號(hào)，星上轉(zhuǎn)發(fā)器接收上行信號(hào)，并產(chǎn)生相干的下行信號(hào)，再由地面跟蹤站接收。三程模式與雙程模式相類似，但是上行發(fā)射站與下行接收站不同，屬于開源測(cè)量。在觀測(cè)建模算法中，三程測(cè)量與雙程測(cè)量的原理相同，僅僅是上下行站不同，因此本文僅針對(duì)三程多普勒測(cè)量進(jìn)行建模分析，將上下行站設(shè)置為相同則得到雙程多普勒建模。

Moyer給出了三程多普勒計(jì)算公式[11]，

(1)

(2)

式中：fR為接收頻率，fS為發(fā)射頻率，M為星上轉(zhuǎn)發(fā)比，c表示光速，ΔT為接收端的積分周期，UTC為協(xié)調(diào)世界時(shí)。下標(biāo)s,e分別表示積分周期的起始與結(jié)束，下標(biāo)1表示上行站發(fā)射信號(hào)，2表示衛(wèi)星轉(zhuǎn)發(fā)信號(hào)，3表示下行站接收星上轉(zhuǎn)發(fā)信號(hào)，12表示信號(hào)上行路徑，23表示信號(hào)下行路徑(下文中的下標(biāo)意義相同)。

為便于計(jì)算，對(duì)式(1)～(2)進(jìn)行變換，三程多普勒可以表示為，

[(TDB1-UTC1)e-(TDB1-UTC1)s]-

(3)

式中：TDB表示質(zhì)心動(dòng)力學(xué)時(shí)，大括號(hào)內(nèi)第1項(xiàng)表示1個(gè)積分周期內(nèi)信號(hào)傳輸?shù)墓庑袝r(shí)差異，第2、3項(xiàng)為時(shí)間系統(tǒng)的差異。

2 三程多普勒測(cè)量模型算法實(shí)現(xiàn)

受限于計(jì)算機(jī)字長(zhǎng)、行星歷表精度，當(dāng)積分周期較小時(shí)(如積分周期為1s)，通過式(3)直接計(jì)算深空探測(cè)器多普勒精度較低，無法滿足優(yōu)于0.01mm/s的測(cè)量需求。通常做法是使用四精度浮點(diǎn)數(shù)定義來解決計(jì)算機(jī)字長(zhǎng)截?cái)鄦栴}，但仍不足以克服星歷表插值誤差，另外四精度浮點(diǎn)數(shù)較雙精度浮點(diǎn)數(shù)在數(shù)值計(jì)算效率上會(huì)大大降低。JPL在建立行星歷表DE414時(shí)發(fā)現(xiàn)使用雙精度與四精度混合定義會(huì)比使用雙精度定義的效率降低30倍[12]。

根據(jù)式(3)，將多普勒計(jì)算步驟分2部分完成：首先計(jì)算第一部分的光行時(shí)，然后計(jì)算第二部分的時(shí)間系統(tǒng)差異，并通過算法改進(jìn)提高計(jì)算精度。

2.1 光行時(shí)的計(jì)算

式(3)中大括號(hào)內(nèi)第1項(xiàng)根據(jù)其物理意義可以表示為，

(TDB3-TDB1)e-(TDB3-TDB1)s=

[(l23)e-(l23)s]+[(l12)e-(l12)s]+

[(T23)e-(T23)s]+[(T12)e-(T12)s]

(4)

式中：l表示信號(hào)傳輸?shù)闹本€距離，T表示引力延遲。積分開始時(shí)刻的直線距離可以表示為，

(5)

式中：rE表示信號(hào)接收時(shí)刻(發(fā)射時(shí)刻)地球在質(zhì)心天球參考系中的位置矢量，rsta表示該時(shí)刻測(cè)站在天球參考系位置矢量，rP表示信號(hào)轉(zhuǎn)發(fā)時(shí)刻中心天體在質(zhì)心天球參考系中的位置矢量，rsat表示衛(wèi)星相對(duì)于中心天體的質(zhì)心天球參考系位置矢量，式(5)的計(jì)算需要進(jìn)行光行時(shí)迭代。

積分結(jié)束時(shí)刻的距離為，

(6)

式中：ΔrE,ΔrP表示表示地球、中心天體從積分開始至積分結(jié)束在質(zhì)心天球參考系中的位置矢量變化，Δrsta表示測(cè)站在一個(gè)積分周期內(nèi)相對(duì)于地球的位置矢量變化，Δrsat表示航天器在1個(gè)積分周期內(nèi)相對(duì)于中心天體的位置矢量變化。

標(biāo)記，

Δr= (ΔrE+Δrsta)-(Δrsat+ΔrP)=

(7)

式中：x,y,z為r的坐標(biāo)分量，Δx,Δy,Δz為Δr的坐標(biāo)分量。則le-ls在r處可做泰勒展開，丟棄高階項(xiàng)有，

(8)

為確保計(jì)算精度，式(7)中ΔrE與ΔrP不宜直接使用行星歷表計(jì)算，而使用質(zhì)心天球參考系中地球與中心天體的速度計(jì)算，

(9)

引力延遲部分可以表示為，

(10)

2.2 時(shí)間系統(tǒng)差異的計(jì)算

式(3)中大括號(hào)內(nèi)第2、3項(xiàng)表示一個(gè)積分周期內(nèi)TDB與UTC差的變化，可以分解為,

[(TDB-UTC)e-(TDB-UTC)s]=

[(TDB-TT)e-(TDB-TT)s]+

[(TT-TAI)e-(TT-TAI)s]+

[(TAI-UTC)e-(TAI-UTC)s]

(11)

TDB與TT的轉(zhuǎn)換關(guān)系一般表示為下式，

TDB-TT=wf+wj+wt

(12)

式中：wf為Fairhead模型，wt為與地理位置相關(guān)的項(xiàng)，wj為基于行星歷表給出的改正項(xiàng)。

經(jīng)分析，行星歷表改正項(xiàng)對(duì)多普勒建模的影響小于10-7m/s，可以不予考慮。Fairhead模型項(xiàng)展開表達(dá)式包含787項(xiàng)，

(13)

ΔT表示從J2000.0起算的時(shí)間，單位為1000年，

(14)

與地理位置相關(guān)的項(xiàng)，

wt= 2.9×10-14usin(UT1+λ+L-LSa)+

1.0×10-13usin(UT1+λ-2M)+

1.33×10-13usin(UT1+λ-D)+

專名屈折變化之非本真性的意義在于殊相相對(duì)于語(yǔ)言的獨(dú)立性：雖然我們對(duì)語(yǔ)言的不同使用方式創(chuàng)造了各種不同的對(duì)象，然而在某種意義上，殊相可以被視為語(yǔ)言不得不原封不動(dòng)地接納的前語(yǔ)言“自在客體”，而作為其語(yǔ)言上的對(duì)應(yīng)物，本真專名在充當(dāng)句子之專名時(shí)在很多語(yǔ)言中都不經(jīng)歷變化。

1.33×10-13usin(UT1+λ+L-LJ)-

2.29×10-13usin(UT1+λ+2L+M)-

2.2×10-12vcos(L+M)+

5.312×10-12usin(UT1+λ-M)

-1.3677×10-11usin(UT1+λ+2L)-

1.31840×10-10vcos(L)+

3.17679×10-10usin(UT1+λ)

(15)

式中：M表示地月系質(zhì)心軌道運(yùn)動(dòng)的平近點(diǎn)角，UT1為世界時(shí)，λ表示測(cè)站地理經(jīng)度，D為日月平角距，L為太陽(yáng)平黃經(jīng)，LJ為木星的日心平黃經(jīng)，LSA為土星的日心平黃經(jīng)。u是測(cè)站距地球自轉(zhuǎn)軸的距離(單位km)，v是測(cè)站距赤道平面北向距離。

時(shí)間系統(tǒng)差異對(duì)積分周期的導(dǎo)數(shù)則可以表示為

(16)

在1950-2050時(shí)間范圍內(nèi)，考慮0.003ps/s精度需要(對(duì)應(yīng)于多普勒約0.001mm/s)，式(16)需要考慮wf的部分系數(shù)(A系數(shù)170項(xiàng)，B系數(shù)4項(xiàng)，C系數(shù)1項(xiàng))，wt的全部項(xiàng)，而wj完全不需要考慮，A,B,C系數(shù)參見SOFA文檔。

2.3TDB至UTC的迭代計(jì)算

式(11)的計(jì)算中涉及到TDB至UTC的轉(zhuǎn)換，而SOFA標(biāo)準(zhǔn)算法僅提供了TT計(jì)算TDB的算法，其說明是可以使用TDB替代TT計(jì)算TDB-TT，因而嚴(yán)格來講，TDB至UTC需要進(jìn)行迭代計(jì)算，具體迭代步驟如圖2。

圖3給出了不迭代計(jì)算與迭代1次計(jì)算的差異，可以看到兩者的差異小于0.5ps，且主要為周年變化，屬于長(zhǎng)周期項(xiàng)，對(duì)于短時(shí)間的多普勒計(jì)算，該差異可以忽略，因此在實(shí)際算中，可以直接使用TDB代替TT計(jì)算，該簡(jiǎn)化引入的誤差對(duì)計(jì)算精度的影響可以忽略。

至此，利用式(3)、 (4)、 (11)可以實(shí)現(xiàn)多普勒的理論值計(jì)算。

3 精度驗(yàn)證

使用常規(guī)的距離差分的方法已經(jīng)對(duì)觀測(cè)建模的正確性進(jìn)行了檢驗(yàn)[13-14]，本節(jié)的計(jì)算主要驗(yàn)證本文方法計(jì)算的精度。這里的常規(guī)算法是通過距離差分實(shí)現(xiàn)，其優(yōu)點(diǎn)是原理與算法簡(jiǎn)單，但由于航天器與測(cè)站距離的直接使用，計(jì)算中截?cái)嗾`差直接影響建模的精度；另外，采用JPL星歷表插值計(jì)算大行星位置時(shí)，其精度僅能達(dá)到0.05mm，該誤差并不能通過距離差分而完全消除。本文提供的算法相對(duì)復(fù)雜，避開了直接使用大行星位置，間接使用了航天器相對(duì)于中心天體的速度與加速度。

以火星快車軌道為例，基于雙精度浮點(diǎn)數(shù)，分別使用常規(guī)數(shù)值算法與本文提供的算法進(jìn)行多普勒測(cè)速的計(jì)算，計(jì)算的弧長(zhǎng)為7天，考慮積分周期為1s，為了便于比較，不考慮測(cè)站與目標(biāo)之間的遮擋關(guān)系，假設(shè)衛(wèi)星全弧段可見，然后比較兩者之間的差異(圖4)。兩者的一致性非常高，絕對(duì)差異小于0.2mm/s，且差異表現(xiàn)為噪聲。

理論上每個(gè)觀測(cè)弧段內(nèi)多普勒測(cè)速值應(yīng)該是光滑連續(xù)的，因此使用多項(xiàng)式擬合后的殘差可以在一定程度上反映多普勒測(cè)速建模的誤差，分別對(duì)兩種算法計(jì)算的結(jié)果選取1小時(shí)觀測(cè)弧段數(shù)據(jù)進(jìn)行擬合。圖5(a)給出了擬合示意圖，使用的是10階多項(xiàng)式，(B)為常規(guī)算法擬合后的殘差，最大偏差達(dá)±0.1mm/s，圖5(c)為使用本文算法計(jì)算多普勒數(shù)據(jù)擬合后的殘差，最大偏差為±0.002mm/s。比較圖5(b)與圖5(c)可知,圖中算法的差異主要由常規(guī)算法計(jì)算的誤差引起。

4 結(jié) 論

論文提出了一種適用于深空探測(cè)任務(wù)的高精度多普勒建模算法，該算法使用環(huán)繞中心天體與航天器的位置、速度、加速度進(jìn)行多普勒測(cè)速的計(jì)算，避免了直接使用測(cè)距建模引入的截?cái)嗾`差。使用常規(guī)距離差分建模算法驗(yàn)證了算法的正確性，通過多項(xiàng)式擬合的方法檢驗(yàn)了算法計(jì)算的精度。

[1] 李平，張紀(jì)生．NASA深空網(wǎng)(DSN)的現(xiàn)狀及發(fā)展趨勢(shì)[J]．飛行器測(cè)控學(xué)報(bào)，2003, 22(4): 10-17．[LiPing,ZhangJi-sheng.ThecurrentstatusandthefutureofNASAdeepspacenetwork[J].JournalofSpacecraftTT&CTechnology, 2003, 22(4): 10-17.]

[2] 于志堅(jiān)，于益農(nóng)，董光亮，等．歐空局深空網(wǎng)的現(xiàn)狀和發(fā)展趨勢(shì)[J]．飛行器測(cè)控學(xué)報(bào)，2004, 12(3)：1-6．[YuZhi-jian,YuYi-nong,DongGuang-liang,etal.StatusandtrendofdevelopmentofESAdeepspacenetwork[J].JournalofSpacecraftTT&CTechnology, 2004, 12(3): 1-6.]

[3]BornGH,ChristensenEJ,FerrariAJ,etal.ThedeterminationofthesatelliteorbitofMariner9[J].CelestialMechanics, 1974, 9: 395-414.

[4]RosenblattP,LaineyV,MaistreSL,etal.AccurateMarsexpressorbitstoimprovethedeterminationofthemassandephemerisoftheMartianmoons[J].PlanetaryandSpaceScience, 2008, 56(7): 1043-1053.

[5]KikuchiF,LiuQ,HanadaH,etal.PicosecondaccuracyVLBIofthetwosubsatellitesofSELENE(KAGUYA)usingmultifrequencyandsamebeammethods[J].RadioScience, 2009, 44(2): 413-421.

[6] 曹建峰，胡松杰，劉磊，等．“嫦娥二號(hào)”小行星探測(cè)試驗(yàn)定軌計(jì)算與分析[J]．北京航天航天大學(xué)學(xué)報(bào)，2014, 40(8): 1095-1101．[CaoJian-feng,HuSong-jie,LiuLei,etal.OrbitdeterminationandanalysisforChang'E-2asteroidexploration[J].JournalofBeijingUniversityofAeronauticsandAstronautics, 2014, 40(8): 1095-1101.]

[7] 鄭為民，馬茂莉，王文彬．深空探測(cè)器被動(dòng)式高精度多普勒測(cè)量方法與應(yīng)用[J]. 宇航學(xué)報(bào), 2013, 34(11): 1462-1467. [ZhengWei-min,MaMao-li,WangWen-bin.High-precisionpassivedopplermeasurementmethodanditsapplicationsindeepspaceexplorater[J].JournalofAstronautics, 2013, 34(11): 1462-1467.]

[8]RowlandsDD,MarshallJA,MccarthyJ,etal.GeodynIIsystemdescription.volume1 [R].HughesSTXCorp.Greenbelt,Maryland, 1993.

[9]JianN,ShangK,ZhangS,etal.Adigitalopen-loopdopplerprocessingprototypefordeep-spacenavigation[J].ScienceinChinaSeriesG:Physics,Mechanics&Astronomy, 2009, 52(12): 1849-1857.

[10] 吳季，朱光武，趙華，等．螢火一號(hào)火星探測(cè)計(jì)劃的科學(xué)目標(biāo)[J]．空間科學(xué)學(xué)報(bào)，2009, 29(5): 449-455．[WuJi,ZhuGuang-wu,ZhaoHua,etal.OverviewofscientificobjectivesofChina-RussiajointMarsexplorationprogramYH-1 [J].Chin.J.SpaceSCI, 2009, 29(5): 449-455.]

[11]MoyerTD.Formulationforobservedandcomputedvaluesofdeepspacenetworkdatatypesfornavigation[M].California:JetPropulsionLaboratory, 2000.

[12]StandishEM.JPLplanetaryephemerisDE414，IOM343R-06-002[R].California:JetPropulsionLaboratory, 2006: 1-8.

[13]CaoJF,HuangY,HuXG,etal.MarsexpresstrackingandorbitdeterminationtrialswithChineseVLBInetwork[J].ChineseScienceBulletin, 2010, 55(32): 3654-3660.

[14] 曹建峰，黃勇，胡小工，等．深空探測(cè)中多普勒的建模與應(yīng)用[J]．宇航學(xué)報(bào)，2011, 32(7): 1583-1589． [CaoJian-feng,HuangYong,HuXiao-gong,etal.Modelingandapplicationofdopplerdataindeepspaceexploration[J].JournalofAstronautics, 2011, 32(7): 1583-1589.]

通信地址：北京市5130信箱120分箱(100094)

電話：13811206835

E-mail:jfcao@foxmail.com

胡松杰(1973-)，男，博士，研究員，主要從事軌道動(dòng)力學(xué)研究。本文通信作者。

通信地址：北京市5130信箱120分箱(100094)

E-mail:husongjie@aliyun.com

(編輯：張宇平)

Modeling and Algorithm Realization of Three-Way Doppler for Deep Space Exploration

CAO Jian-feng1,2,3, HUANG Yong3, LIU Lei1,2, ZHANG Yu1,2, ZHENG Ai-wu1,2, HU Song-jie1,2

(1. Science and Technology on Aerospace Flight Dynamics Laboratory, Beijing 100094, China; 2.Beijing Aerospace Control Center, Beijing 100094, China; 3. Shanghai Astronomical Observatory, Shanghai 200030, China)

Aiming at the high precision Doppler application for deep space exploration, a new algorithm for Doppler modeling using the velocity and acceleration of the central object and the spacecraft is proposed. The algorithm avoids the truncation error induced by the direct use of ranging and the positional uncertainty of ephemerides interpolation. Continuous Doppler computation for a Martian orbiter lasting one week is performed, and the general computation algorithm using ranging difference is used to test the correctness of the new algorithm. In addition, one-hour length of the Doppler data is fitted with a polynomial of order 10, and the residual error is less than 0.002mm/s. The analysis results show that the new algorithm is correct and is approximately one order of magnitude improvement in accuracy compared with the conventional modeling algorithm.

Deep space exploration; Doppler; Observation modeling

2016-08-01;

2016-11-21

國(guó)家自然科學(xué)基金(11203003, 61573049, 11373013, 61304233, 11473056)

V412.4+1

A

1000-1328(2017)03-0304-06

10.3873/j.issn.1000-1328.2017.03.011

曹建峰(1982-)，男，博士，工程師，主要從事軌道動(dòng)力學(xué)研究。