，，，

1. 長(zhǎng)沙理工大學(xué) 計(jì)算機(jī)與通信工程學(xué)院 綜合交通運(yùn)輸大數(shù)據(jù)智能處理湖南省重點(diǎn)實(shí)驗(yàn)室，長(zhǎng)沙 410004 2. 國(guó)防科學(xué)技術(shù)大學(xué) 理學(xué)院，長(zhǎng)沙 410073 3. 中國(guó)人民解放軍91550部隊(duì)，大連116023

隨著航天技術(shù)的發(fā)展及深空探測(cè)等任務(wù)的興起，傳統(tǒng)的依靠地面站的航天器導(dǎo)航方式不能滿足當(dāng)前任務(wù)的需求[1]。航天器自主導(dǎo)航，是指僅利用航天器自身搭載的敏感器信息進(jìn)行導(dǎo)航，是航天器自主運(yùn)行和管理的核心技術(shù)之一，是提高其生存能力和降低運(yùn)營(yíng)成本的關(guān)鍵技術(shù)，受到了國(guó)內(nèi)外的廣泛關(guān)注。

導(dǎo)航信息的獲取是航天器實(shí)現(xiàn)自主導(dǎo)航的前提，目前主要有慣性導(dǎo)航、衛(wèi)星導(dǎo)航、地磁導(dǎo)航、傳統(tǒng)天文導(dǎo)航等。慣性導(dǎo)航技術(shù)具有完全自主、可靠性強(qiáng)、適用范圍廣等優(yōu)點(diǎn)，但是不能單獨(dú)在長(zhǎng)期飛行的航天器上應(yīng)用[2]。衛(wèi)星導(dǎo)航目前應(yīng)用廣泛，具有成本低、精度高等優(yōu)點(diǎn)，但由于依賴人工信標(biāo)，易被摧毀和受干擾，不能夠完全自主[3]。地磁導(dǎo)航具有全天時(shí)、全天候、全地域等優(yōu)點(diǎn)，但受制于地磁模型精度及磁強(qiáng)計(jì)測(cè)量精度，目前難以達(dá)到航天器自主導(dǎo)航的精度要求[4-6]。傳統(tǒng)天文導(dǎo)航誤差不隨時(shí)間積累，自主性強(qiáng)，但易受氣象條件、敏感器性能等的影響[7-8]?？紤]到太陽(yáng)、地球、月球作為自然天體，具有運(yùn)動(dòng)規(guī)律清楚、光學(xué)特征突出等特點(diǎn)，容易識(shí)別，導(dǎo)航資源不受限制[9-10]。目前，太陽(yáng)敏感器、地球敏感器及月球敏感器技術(shù)已比較成熟，利用敏感器可以得到航天器與日地月之間的夾角信息，可用于航天器的自主導(dǎo)航。

濾波算法是自主導(dǎo)航處理測(cè)量數(shù)據(jù)求解軌道信息的主要方法，目前主要采用擴(kuò)展卡爾曼濾波(EKF)和無(wú)跡卡爾曼濾波(UKF)算法作為實(shí)時(shí)濾波算法。由于航天器自主導(dǎo)航系統(tǒng)中狀態(tài)模型和天文導(dǎo)航模型均為非線性的，EKF采用一階近似的泰勒逼近方法，損失部分精度，影響導(dǎo)航精度的提高；UKF雖然適用于非線性濾波，但當(dāng)估計(jì)維數(shù)大時(shí)，精度也會(huì)受到損失。CKF[11]是近年來(lái)提出的一種新的非線性狀態(tài)估計(jì)方法，適合于高維狀態(tài)下的實(shí)時(shí)濾波系統(tǒng)[12]。

本文研究了基于日地月天文信息的航天器全弧段自主導(dǎo)航方法。為克服導(dǎo)航系統(tǒng)的非線性特性及日月天文信息的不完整性等缺點(diǎn)，將自主導(dǎo)航分為日月可見弧段和星蝕階段兩部分?？紤]到航天器全弧段全過程的自主導(dǎo)航，其系統(tǒng)模型非線性程度高，待估參數(shù)較多，而容積卡爾曼濾波(CKF)的估計(jì)精度在高維狀態(tài)下優(yōu)于UKF算法，因此，本文采用了CKF方法進(jìn)行導(dǎo)航計(jì)算。在日月可見弧段，利用太陽(yáng)、地球、月球敏感器和測(cè)高儀提供的天文信息及航天器上存儲(chǔ)的日月星歷表，確定航天器與日地月之間的夾角及在地心赤道慣性坐標(biāo)系下的位置矢量，并以之作為觀測(cè)量，采用CKF方法進(jìn)行航天器自主定軌；在星蝕階段，利用日月可見弧段信息，用軌道預(yù)報(bào)的方式進(jìn)行導(dǎo)航，從而實(shí)現(xiàn)航天器全弧段高精度自主導(dǎo)航。

1 導(dǎo)航系統(tǒng)數(shù)學(xué)模型

1.1 航天器軌道動(dòng)力學(xué)模型

考慮地球形狀、大氣阻力、日月引力、太陽(yáng)光壓攝動(dòng)和潮汐攝動(dòng)，在地心赤道慣性坐標(biāo)系下，建立航天器的軌道動(dòng)力學(xué)模型：

(1)

式中：X=(r,v)T為航天器的狀態(tài)矢量；w=(wr,wv)T為系統(tǒng)過程噪聲:f為航天器在軌攝動(dòng)，可以表示為

(2)

式中：μe為重力場(chǎng)參數(shù)；U為地球吸引力攝動(dòng)位函數(shù)；Fa為大氣阻力攝動(dòng)力加速度；Fs為太陽(yáng)對(duì)航天器的加速度；Fm為月球?qū)教炱鞯募铀俣?；Fsr為太陽(yáng)光壓加速度；Ft為潮汐攝動(dòng)加速度[13]。令

(3)

則航天器導(dǎo)航系統(tǒng)的狀態(tài)模型可簡(jiǎn)化為:

(4)

1.2 觀測(cè)模型

如圖1所示，通過太陽(yáng)、地球、月球敏感器及測(cè)高儀,可得到某一時(shí)刻在星體坐標(biāo)系中太陽(yáng)、月球的方向單位矢量VSc0、VMc0以及地心至航天器的方向單位矢量VEc0和航天器到地球表面的距離Lc[12]。

圖1 航天器在軌位置矢量Fig.1 On-orbit position vector for spacecraft

圖2 月球、地球與航天器的幾何關(guān)系Fig.2 Geometric relationship between moon, earth and spacecraft

利用航天器在軌存儲(chǔ)的日月星歷表，可得地心赤道慣性坐標(biāo)系中月地矢量VEM和日地矢量VES。由月地距離的有限性，可知航天器-月球矢量和地球-月球矢量不平行。假設(shè)地球?yàn)榍蝮w，Re為地球半徑，則如圖2所示，由VEc=(Lc+Re)VEc0，|VEMc|=|VEM|，VMc0·VEc0=cosu1，VMc=VEc+VEMc,VMc=|VMc|VMc0，就可得到VEMc。

考慮到太陽(yáng)距離地球和航天器非常遠(yuǎn)，可將航天器-太陽(yáng)矢量和地球-太陽(yáng)矢量近似看作平行。設(shè)A為星體坐標(biāo)系到地心慣性坐標(biāo)系的轉(zhuǎn)換矩陣，則有

VEM=AVEMc,VEC=AVEC,VEM×VEC≠0

(5)

式中：VEC為地心慣性坐標(biāo)系下地航天器到地心的矢量。利用觀測(cè)量在地心慣性坐標(biāo)系中建立新的正交坐標(biāo)系[12]，設(shè)各坐標(biāo)軸的單位矢量分別為：

(6)

同樣，在星體坐標(biāo)系中也建立一個(gè)正交坐標(biāo)系，設(shè)各坐標(biāo)軸的單位矢量分別為:

(7)

令

(8)

則由式(3)～(5)有

U=AW

(9)

從而可得

A=UW-1=UWT

(10)

故在地心慣性坐標(biāo)系下太陽(yáng)、月球及地球的方向單位矢量VS0、VM0、VE0分別為：

VS0=AVSc0,VM0=AVMc0,VE0=AVEc0

(11)

則地球-航天器-月球形成的夾角為：

u1=arccos(VE0·VM0)

(12)

地球-航天器-太陽(yáng)形成的夾角為：

u2=arccos(VE0·VS0)

(13)

地心距為：

u3=Lc+Re

(14)

將航天器在地心赤道慣性坐標(biāo)系下的夾角信息及地心距作為觀測(cè)量，得觀測(cè)模型為：

Z=u+ε

(15)

式中：u=(u1,u2,u3)T；ε為觀測(cè)噪聲,且與導(dǎo)航系統(tǒng)過程噪聲w獨(dú)立。

令

u=h(X)

(16)

則導(dǎo)航系統(tǒng)的觀測(cè)模型可簡(jiǎn)化為：

Z=h(X)+ε

(17)

2 航天器軌道的自主確定

如圖1所示，根據(jù)航天器與日地月的空間幾何關(guān)系，可將航天器的自主軌道確定分為兩部分，日月可見弧段的自主定軌和星蝕弧段的自主定軌。

2.1 日月可見弧段的CKF濾波

對(duì)于日月可見弧段，可利用航天器在軌的太陽(yáng)、地球、月球敏感器及測(cè)高儀直接獲得導(dǎo)航觀測(cè)數(shù)據(jù)。對(duì)于由式(4)和式(17)所構(gòu)成的非線性導(dǎo)航系統(tǒng)：

(18)

考慮到航天器日月可見弧段自主導(dǎo)航系統(tǒng)模型的非線性程度高，待估參數(shù)較多，而CKF濾波的估計(jì)精度在高維狀態(tài)下優(yōu)于UKF算法。為降低非線性因素對(duì)導(dǎo)航濾波性能的影響，本文采用CKF濾波算法進(jìn)行航天器自主導(dǎo)航計(jì)算。

CKF濾波采用三階容積法則，用數(shù)值積分來(lái)近似高斯加權(quán)積分，其核心是求解積分形式為“非線性函數(shù)*高斯密度函數(shù)”的積分[14]，利用一組等權(quán)值容積點(diǎn)加權(quán)求和來(lái)代替加權(quán)高斯問題，也就是在計(jì)算隨機(jī)變量經(jīng)非線性變換后的均值和協(xié)方差時(shí)采用容積數(shù)值積分原則。針對(duì)式(18)的自主導(dǎo)航系統(tǒng)，CKF濾波算法具體流程如下。

(1)CKF濾波初始化

(19)

(20)

(2)時(shí)間更新

1)對(duì)于k-1時(shí)刻的狀態(tài)濾波誤差陣，將其因式分解，

(21)

2)估計(jì)容積點(diǎn):

(22)

其中：

ξi=

m=2nX，nX為非線性導(dǎo)航系統(tǒng)方程中的狀態(tài)向量的維數(shù)。

3)容積點(diǎn)傳播:

(23)

4)估計(jì)狀態(tài)預(yù)測(cè)值:

(24)

5)估計(jì)誤差協(xié)方差預(yù)測(cè)值:

(25)

式中：Qk-1為第k-1時(shí)刻的系統(tǒng)過程噪聲的協(xié)方差陣。

(3)觀測(cè)更新

1)因式分解：

(26)

2)估計(jì)容積點(diǎn)：

(27)

3)容積點(diǎn)傳播：

Zi,k|k-1=h(Xi,k|k-1)

(28)

4)估計(jì)觀測(cè)預(yù)測(cè)值：

(29)

5)估計(jì)新息協(xié)方差矩陣：

(30)

式中：Rk為第k時(shí)刻的測(cè)量噪聲的協(xié)方差陣。

6)估計(jì)互協(xié)方差矩陣:

[Xi,k|k-1-Zi,k|k-1]T

(31)

7)估計(jì)卡爾曼增益:

(32)

8)估計(jì)狀態(tài)更新:

(33)

9)估計(jì)相應(yīng)的誤差協(xié)方差:

(34)

2.2 星蝕時(shí)段的軌道預(yù)報(bào)

星蝕包括日蝕、月蝕和朔月。該觀測(cè)弧段內(nèi)，由于地球遮擋，太陽(yáng)敏感器和月球敏感器不能獲取當(dāng)前的測(cè)量數(shù)據(jù)，因此需要利用日月可見弧段得到的航天器軌道參數(shù)，通過軌道預(yù)報(bào)方法對(duì)星蝕時(shí)段的航天器進(jìn)行導(dǎo)航。

本文直接選擇軌道動(dòng)力學(xué)模型作為預(yù)報(bào)模型，由于各攝動(dòng)力的影響程度不同，除引力場(chǎng)攝動(dòng)帶諧項(xiàng)J2為10-3量級(jí)外，引力場(chǎng)高階攝動(dòng)帶諧項(xiàng)及其他攝動(dòng)力最大為10-5量級(jí)[13]，因此本文只考慮帶諧項(xiàng)J2攝動(dòng)，即預(yù)報(bào)模型為

(35)

其中:

J21P21(sinφ)·cos(λ-λ21)+J22P22(sinφ)·

cos2(λ-λ22)]}

(36)

式中符號(hào)表示見文獻(xiàn)[13]。將日月可見弧段得到的軌道參數(shù)作為軌道預(yù)報(bào)初值代入軌道預(yù)報(bào)模型，采用龍格庫(kù)塔數(shù)值積分方法，即可進(jìn)行星蝕時(shí)段的軌道預(yù)報(bào)，步驟如下：

X(t0)=X0

(37)

式中：X0為日月可見弧段得到的軌道參數(shù)。則采用四階龍格庫(kù)塔方法的軌道預(yù)報(bào)公式為：

(38)

式中：Xn和Xn+1分別為前一步點(diǎn)和當(dāng)前步點(diǎn)的值；h為積分步長(zhǎng)。

3 仿真試驗(yàn)

3.1 仿真參數(shù)設(shè)置

對(duì)由式(18)給出的航天器自主導(dǎo)航系統(tǒng)進(jìn)行全弧段全過程導(dǎo)航仿真試驗(yàn)，分別利用UKF算法[15]和本文采用的CKF算法對(duì)日月可見弧段的航天器定軌。仿真參數(shù)設(shè)置如下。

1)坐標(biāo)系：J2000.0地心赤道慣性坐標(biāo)系。

2)軌道參數(shù)：

①仿真初始?xì)v元2017年1月1日0時(shí)0分0秒；

②衛(wèi)星高度400 km，軌道傾角98°。

3)航天器動(dòng)力學(xué)模型：

①航天器理論軌道仿真模型采用的受力模型是地球引力非球形引力模型取JGM-3(70×70項(xiàng))，大氣阻力、太陽(yáng)引力、月球引力、太陽(yáng)光壓攝動(dòng)及潮汐攝動(dòng)；

②航天器導(dǎo)航系統(tǒng)狀態(tài)模型采用的受力模型是地球引力，地球引力非球形引力模型取JGM-3(20×20項(xiàng))，大氣阻力(大氣阻力因子Cd待估)。

4)導(dǎo)航敏感器的測(cè)量精度設(shè)置：太陽(yáng)敏感器5′(3σ)，月球敏感器5′(3σ)，地球敏感器5′(3σ)，測(cè)高儀1 km。

5)濾波初始誤差：位置50 km；速度5 km/s。

日月可見弧段，數(shù)據(jù)更新頻率為0.5 Hz，仿真時(shí)長(zhǎng)7 h。

3.2 仿真結(jié)果分析

對(duì)于如上的仿真模型，分別對(duì)日月可觀測(cè)弧段和星蝕時(shí)段進(jìn)行航天器自主導(dǎo)航，在日月可觀測(cè)弧段，分別采用CKF濾波和UKF濾波進(jìn)行導(dǎo)航解算，而在星蝕時(shí)段，則采用數(shù)值積分進(jìn)行軌道預(yù)報(bào)。圖3和圖4分別給出了單次仿真下的基于CKF濾波的航天器自主導(dǎo)航位置和速度方向上的濾波偏差，即位置誤差和速度誤差。為了進(jìn)一步得到航天器位置和速度各分量的導(dǎo)航精度，利用蒙特卡洛仿真，表1和表2是對(duì)50次仿真得到的位置誤差和速度誤差的統(tǒng)計(jì)結(jié)果。其中表1統(tǒng)計(jì)的是基于CKF濾波的航天器各時(shí)段的導(dǎo)航精度，表2統(tǒng)計(jì)的是日月可觀測(cè)弧段內(nèi)CKF濾波和UKF濾波的自主導(dǎo)航精度比對(duì)。

圖3 日月可見弧段位置估計(jì)誤差Fig.3 Position estimation error in the observable period

圖4 日月可見弧段速度估計(jì)誤差Fig.4 Position estimation error in the eclipse period

從仿真結(jié)果，可以看出，利用日地月天文信息進(jìn)行航天器自主導(dǎo)航，采用CKF濾波算法能很快收斂，CKF導(dǎo)航位置精度穩(wěn)定在2 km以內(nèi)。當(dāng)出現(xiàn)星蝕時(shí)，敏感器不能捕獲天體目標(biāo)信息，導(dǎo)航算法自動(dòng)切換到預(yù)報(bào)程序，從圖3和圖4可以看出，星蝕階段依賴預(yù)報(bào)算法可以實(shí)現(xiàn)目標(biāo)的正常導(dǎo)航，導(dǎo)航精度在短時(shí)間星蝕階段不受影響。此外，從表2中可以看出，利用CKF濾波算法進(jìn)行日月可觀測(cè)弧段內(nèi)航天器自主導(dǎo)航，濾波精度要優(yōu)于UKF濾波器。

表1 各時(shí)段基于CKF濾波的航天器自主導(dǎo)航精度

表2 日月可觀測(cè)弧段不同濾波算法的導(dǎo)航性能分析

4 結(jié)束語(yǔ)

本文提出了以日地月天文信息為觀測(cè)數(shù)據(jù)的航天器自主導(dǎo)航算法。對(duì)于日月可見弧段和星蝕時(shí)段分別采用導(dǎo)航濾波器和軌道數(shù)值預(yù)報(bào)的方式進(jìn)行導(dǎo)航，基于CKF濾波的導(dǎo)航精度在2 km以內(nèi)，且性能優(yōu)于UKF導(dǎo)航濾波器。仿真試驗(yàn)結(jié)果表明這種自主導(dǎo)航方法是可行的，而且得到的定位精度能夠滿足航天器自主導(dǎo)航的要求。但是，導(dǎo)航敏感器測(cè)量精度是影響自主導(dǎo)航精度的主要因素之一，而本文在運(yùn)用導(dǎo)航濾波時(shí)僅考慮導(dǎo)航敏感器的隨機(jī)誤差，忽略了導(dǎo)航系統(tǒng)的模型誤差或不確定性誤差，如敏感器測(cè)量系統(tǒng)的常值偏差、動(dòng)力學(xué)模型誤差等，因此后續(xù)將進(jìn)一步研究存在模型誤差或不確定性誤差下的航天器自主導(dǎo)航方法。此外，對(duì)于星蝕時(shí)段，一方面可以直接用軌道動(dòng)力學(xué)進(jìn)行數(shù)值預(yù)報(bào)，還可結(jié)合慣性導(dǎo)航器件進(jìn)行軌道估計(jì)，下一步將比較分析慣性導(dǎo)航器件不同誤差特性下的星蝕時(shí)段內(nèi)的軌道估計(jì)精度。

References)

[1] 房建成, 寧曉琳. 深空探測(cè)器自主天文導(dǎo)航方法[M]. 西安: 西安工業(yè)大學(xué)出版社, 2010: 124-138.

FANG J C, NING X L. Autonomous celestial navigation for deep space explorer [M]. Xi′an: Xi′an Industry University Press, 2010: 124-138(in Chinese).

[2] 周瓊峰, 唐國(guó)華, 陳雪芬. MEMS陀螺旋轉(zhuǎn)導(dǎo)航誤差的分析與仿真[J]. 中國(guó)空間科學(xué)技術(shù), 2016, 36(4): 74-80.

ZHOU Q F, TANG G H, CHEN X F. Analysis and simulation for rotatinginertial navigation error of MEMS[J]. Chinese Space Science and Technology, 2016, 36(4): 74-80(in Chinese).

[3] WANG Y Q, LUO Y, WANG P. A novel joint navigation state error discriminator based on iterative maximum likelihood estimation[J]. Science China (Information Sciences), 2015, 58(5): 125-129.

[4] 王鵬, 張迎春. 基于磁強(qiáng)計(jì)/太陽(yáng)敏感器的自主導(dǎo)航方法[J]. 系統(tǒng)工程與電子技術(shù), 2013, 35(1): 132-137.

WANG P, ZHANG Y C. Research on autonomous navigation method based on magnetometer and sum sensor[J]. Systems Engineering and Electronics, 2013, 35(1): 132-137(in Chinese).

[5] DEUTSCHMAN J, ITZHACK I Y. Evaluation of attitude and orbit estimation using actual earth magnetic field data[J]. Journal of Guidance, Control and Dynamics, 2001, 24(3): 616-626.

[6] COLLINS J T, CONGER R E. MANS: Autonomous navigation and orbit control for communication satellites[J]. AIAA, CP, 1994:94-112.

[7] WERNER E. Galileo: impact on spacecraft navigation system[J]. Journal of Global Positioning Systems, 2003, 2(2): 135-138.

[8] 楊博. 航天器星敏感器自主定位方法及精度分析[J]. 宇航學(xué)報(bào), 2002, 23(3): 81-84.

YANG B. The method of autonomous position determination using star sensor and analysis of precision for spacecraft[J]. Journal of Astronautics,2002, 23(3): 81-84(in Chinese).

[9] 黃翔宇, 崔平遠(yuǎn), 崔祜濤. 深空自主導(dǎo)航系統(tǒng)的可觀性分析[J]. 宇航學(xué)報(bào), 2006, 27(3): 332-337.

HUANG X Y, CUI P Y, CUI H T. Observability analysis of deep space autonomous navigation system[J]. Journal of Astronautics, 2006, 27(3): 332-337(in Chinese).

[10] 劉萍, 王大軼, 黃翔宇. 環(huán)月探測(cè)器自主天文導(dǎo)航系統(tǒng)的可觀度分析[J]. 中國(guó)空間科學(xué)技術(shù), 2007, 27(6): 12-18.

LIU P, WANG D Y, HUANG X Y. Analysis of observability degree in autonomous celestial navigation system of lunar orbiter[J]. Chinese Space Science and Technology,2007, 27(6):12-18(in Chinese).

[11] ARASARATNAM I, HAYKIN S. Cubature Kalman filters[J]. IEEE Transactions on Automatic Control, 2009, 54(6): 1254-1269.

[12] 朱蓉. CKF濾波算法及其在航天器自主導(dǎo)航中的應(yīng)用[D]. 長(zhǎng)沙: 國(guó)防科學(xué)技術(shù)大學(xué), 2015.

ZHU R. CKF and its application on spacecraft autonomous navigation[D]. Changsha: National University of Defense Technology, 2015(in Chinese).

[13] 李濟(jì)生，王家松.航天器軌道確定[M]. 北京: 國(guó)防工業(yè)出版社, 2003:106-188.

[14] ARASARATNAM I, HAYKIN S, HURD T R. Cubature Kalman filtering for continuous-discrete systems: theory and simulations[J]. IEEE Transactions on Signal Processing, 2010, 58(10): 4977-4993.

[15] JULIER S, UHLMANN J. Unscented filtering and nonlinear estimation[J]. IEEE Transactions on Aerospace and Electronic Systems, 2004, 92(3): 401-422.