辛寧 邱樂德 張立華 劉乃金

(1 中國空間技術(shù)研究院通信衛(wèi)星事業(yè)部, 北京 100094) (2 航天東方紅衛(wèi)星有限公司, 北京 100094)

?

一種重力測量衛(wèi)星質(zhì)心在軌標(biāo)定改進(jìn)算法

辛寧1邱樂德1張立華2劉乃金1

(1 中國空間技術(shù)研究院通信衛(wèi)星事業(yè)部, 北京 100094) (2 航天東方紅衛(wèi)星有限公司, 北京 100094)

針對重力測量衛(wèi)星質(zhì)心在軌標(biāo)定算法中存在的陀螺精度不高且易失效的問題，提出了一種應(yīng)用星敏感器和靜電加速度計(jì)的質(zhì)心標(biāo)定改進(jìn)算法。在該算法中，將一個明顯大于干擾力矩的周期性磁力矩，作用于衛(wèi)星上；將星敏感器觀測量代入預(yù)測濾波器中，估計(jì)衛(wèi)星的角速度及角加速度；利用靜電加速度計(jì)信息設(shè)計(jì)擴(kuò)展卡爾曼濾波，從而實(shí)現(xiàn)衛(wèi)星質(zhì)心的標(biāo)定。數(shù)學(xué)仿真結(jié)果表明：此算法能夠?qū)πl(wèi)星的角速度及質(zhì)心位置進(jìn)行實(shí)時估計(jì)，衛(wèi)星質(zhì)心三軸最佳標(biāo)定精度均優(yōu)于0.05 mm，可實(shí)現(xiàn)陀螺失效情況下衛(wèi)星質(zhì)心較為精確的標(biāo)定。

重力測量衛(wèi)星；質(zhì)心在軌標(biāo)定；預(yù)測濾波；擴(kuò)展卡爾曼濾波；星敏感器

1 引言

為得到高精度的地球重力場數(shù)據(jù)，重力測量衛(wèi)星要利用靜電加速度計(jì)測量非重力因素對衛(wèi)星軌道的影響[1]。靜電加速度計(jì)中檢驗(yàn)質(zhì)量塊的質(zhì)心與衛(wèi)星質(zhì)心的位置偏差，將直接影響非重力因素的測量精度。產(chǎn)生上述偏差的主要因素包括：地面衛(wèi)星質(zhì)心校準(zhǔn)存在的誤差；衛(wèi)星質(zhì)量的損失，如衛(wèi)星在軌飛行中燃料的消耗，衛(wèi)星在空間真空電離環(huán)境中長期運(yùn)行時部分物質(zhì)的蒸發(fā)等；其他因素的影響，如衛(wèi)星發(fā)射段振動[2]。因此，開展衛(wèi)星質(zhì)心的在軌標(biāo)定算法研究，為質(zhì)心在軌修正任務(wù)提供依據(jù)，具有重要的意義[3]。

目前，重力測量衛(wèi)星質(zhì)心的在軌標(biāo)定算法，必須利用陀螺來獲取衛(wèi)星的角速度信息，且存在一定的局限性。陀螺的誤差漂移特性會使測量精度降低，而且陀螺有可能失效，從而導(dǎo)致質(zhì)心在軌標(biāo)定算法無法應(yīng)用。文獻(xiàn)[4]中詳細(xì)介紹了“重力校正和氣候?qū)嶒?yàn)”(GRACE)衛(wèi)星的質(zhì)心在軌標(biāo)定算法，該算法除了需要陀螺的角速度數(shù)據(jù)外，還要精確已知控制輸入以及衛(wèi)星姿態(tài)動力學(xué)特性，這些參數(shù)在質(zhì)心標(biāo)定期間往往存在不確定性。文獻(xiàn)[5]中提出將陀螺數(shù)據(jù)代入最小模型誤差估計(jì)算法中估計(jì)衛(wèi)星角加速度，利用擴(kuò)展卡爾曼濾波估計(jì)衛(wèi)星質(zhì)心。最小模型誤差估計(jì)算法是一種離線估計(jì)算法，因此無法對質(zhì)心位置進(jìn)行實(shí)時估計(jì)。為了提高重力測量衛(wèi)星質(zhì)心在軌標(biāo)定算法的可靠性，避免陀螺的不利影響，本文在目前已有研究的基礎(chǔ)上，提出了一種應(yīng)用高精度星敏感器和靜電加速度計(jì)的質(zhì)心標(biāo)定改進(jìn)算法，在保證質(zhì)心標(biāo)定精度的同時，可為重力測量衛(wèi)星設(shè)計(jì)提供一種新的解決方案。

本文提出的質(zhì)心標(biāo)定改進(jìn)算法，首先利用高精度星敏感器，結(jié)合非線性預(yù)測濾波估計(jì)衛(wèi)星的角速度及角加速度；然后通過擴(kuò)展卡爾曼濾波估計(jì)實(shí)現(xiàn)衛(wèi)星質(zhì)心的在軌標(biāo)定。該算法中的濾波器不需要衛(wèi)星的控制輸入以及動力學(xué)特性信息，能夠?qū)崟r估計(jì)衛(wèi)星質(zhì)心、角速度及角加速度，具有快速的收斂性能和良好的穩(wěn)定精度，可實(shí)現(xiàn)無陀螺條件下的衛(wèi)星質(zhì)心在軌標(biāo)定。

2 衛(wèi)星質(zhì)心在軌標(biāo)定改進(jìn)算法

重力測量衛(wèi)星質(zhì)心在軌標(biāo)定改進(jìn)算法的基本原理如圖1所示。根據(jù)地球磁場模型及衛(wèi)星的軌道位置，給衛(wèi)星施加一定頻率的磁力矩，使衛(wèi)星產(chǎn)生相應(yīng)的轉(zhuǎn)動；將靜電加速度計(jì)、星敏感器作為檢測傳感器，對星敏感器觀測矢量進(jìn)行預(yù)測濾波處理，估計(jì)出角速度及角加速度；對靜電加速度計(jì)觀測矢量進(jìn)行擴(kuò)展卡爾曼濾波處理，分離出具有線性特征的非重力加速度分量，由此來確定衛(wèi)星質(zhì)心的位置。

圖1 衛(wèi)星質(zhì)心在軌標(biāo)定流程Fig.1 Flow of on-orbit calibration of satellite center of mass

2.1 系統(tǒng)模型

重力測量衛(wèi)星本體坐標(biāo)系與靜電加速度計(jì)坐標(biāo)系的關(guān)系如圖2所示，衛(wèi)星本體坐標(biāo)系的原點(diǎn)OS1為衛(wèi)星質(zhì)心，XS1軸垂直于衛(wèi)星軌道平面，YS1軸表示衛(wèi)星飛行方向，ZS1軸表示衛(wèi)星對地指向。靜電加速度計(jì)坐標(biāo)系的原點(diǎn)為檢驗(yàn)質(zhì)量塊質(zhì)心位置OA1，理論上與OS1重合，YA1軸和ZA1軸分別平行于衛(wèi)星本體坐標(biāo)系的ZS1軸和XS1軸。

圖2 衛(wèi)星本體坐標(biāo)系與靜電加速度計(jì)坐標(biāo)系Fig.2 Satellite body-fixed coordinate system and electrostatic accelerometer coordinate system

靜電加速度計(jì)在衛(wèi)星本體坐標(biāo)系下輸出的非保守力可表示為[6]

Fobs=MAB(s1Fout)+Ab+An

(1)

式中：MAB為靜電加速度計(jì)坐標(biāo)系到衛(wèi)星本體坐標(biāo)系的轉(zhuǎn)換矩陣；s1為靜電加速度計(jì)的標(biāo)度因子；Fout為靜電加速度計(jì)在靜電加速度計(jì)坐標(biāo)系下的實(shí)際輸出，見式(2)；Ab為零偏誤差；An為噪聲誤差。

(2)

由文獻(xiàn)[6]中可知，由重力梯度和非保守力因素引起的加速度分量，在較短的標(biāo)定時間t內(nèi)基本保持線性關(guān)系，因此質(zhì)心在軌標(biāo)定的測量方程為

(3)

(4)

顯然，在構(gòu)建上述測量方程時，必須要利用衛(wèi)星的角速度和角加速度信息。由于陀螺精度不高，難以準(zhǔn)確獲取角速度及角加速度信息，這成為重力測量衛(wèi)星質(zhì)心在軌標(biāo)定的一個難點(diǎn)。為此，本文將星敏感器測量的衛(wèi)星姿態(tài)四元數(shù)代入到預(yù)測濾波中，以精確地獲取衛(wèi)星角速度及角加速度信息；再利用擴(kuò)展卡爾曼濾波對質(zhì)心位置進(jìn)行標(biāo)定。

2.2 角速度預(yù)測濾波估計(jì)

(5)

(2)求取正弦基模參數(shù)S(i)。

(6)

(3)利用正弦基模參數(shù)S(i)和時間序列i的載體姿態(tài)四元數(shù)q(i+1)，可以求得i+1時的轉(zhuǎn)換矩陣Tq，ω(i+1)。

(7)

式中：q0，q1，q2，q3為姿態(tài)四元數(shù)q的4個分量。

(8)

式中：Nu為由姿態(tài)四元數(shù)q的標(biāo)部產(chǎn)生的標(biāo)部值。

由式(7)和式(8)可得角速度估計(jì)值的狀態(tài)方程為

(9)

式中：d(i)為模型誤差；I為單位向量。

使用四元數(shù)表示的衛(wèi)星姿態(tài)運(yùn)動微分方程為

(10)

將式(10)表示為離散形式，根據(jù)四元數(shù)定義可得角速度估計(jì)值的觀測方程為

q(i)=

q(i-1)+v(i)

(11)

式中：v(i)為零均值高斯白噪聲。

預(yù)測濾波的物理概念[7]是：設(shè)已經(jīng)獲得了i時的狀態(tài)估計(jì)值，采樣間隔時間為Δi，在接到i+Δi時的測量信息q(i+1)后，利用模型最小誤差逼近可預(yù)測[i,i+Δi]范圍內(nèi)的模型誤差d1(i)，然后將d1(i)的預(yù)測值代入系統(tǒng)狀態(tài)方程和觀測輸出方程中，將狀態(tài)估計(jì)值傳播到i+1時，從而得到較高的估計(jì)精度。

為了估計(jì)出模型誤差，首先將q(i+1)采用泰勒級數(shù)法展開，直到模型誤差d1(i)出現(xiàn)，有

q(i)+y(i)+Λ(Δi)S(i)d1(i)

(12)

(13)

根據(jù)預(yù)測濾波理論，由式(12)可得

d1(i)=-{[Λ(Δi)S(i)]TR-1[Λ(Δi)S(i)]+W}-1·

[Λ(Δi)S(i)]TR-1[y(i)-q(i+1)+q(i)]

(14)

式中：W為3×3階加權(quán)矩陣；R為量測噪聲協(xié)方差矩陣。

上述過程的實(shí)質(zhì)是先預(yù)測下一時刻的模型誤差校正量，然后再對狀態(tài)估計(jì)值進(jìn)行誤差修正，從而有效地提高了角速度的估計(jì)精度。

2.3 角加速度預(yù)測濾波估計(jì)

衛(wèi)星角加速度的表達(dá)式為

(15)

式中：Jp為衛(wèi)星轉(zhuǎn)動慣量；M為衛(wèi)星力矩。

將式(15)等號右邊的狀態(tài)信息全部定義為模型誤差d2(i)，由此獲得的角加速度的系統(tǒng)狀態(tài)方程和觀測方程為

(16)

2.4 質(zhì)心的擴(kuò)展卡爾曼濾波估計(jì)

在衛(wèi)星角速度及角加速度估計(jì)值已知的條件下，質(zhì)心估計(jì)參數(shù)X=[lαβ]T，進(jìn)行估計(jì)的系統(tǒng)狀態(tài)方程及觀測方程為

(17)

待估狀態(tài)變量一步預(yù)測估計(jì)均值為

(18)

式中：k為時間序列標(biāo)志符；狀態(tài)轉(zhuǎn)移矩陣Φ(tk)=I9×9。

一步估計(jì)均方誤差矩陣為

P(tk+1/k)=Φ(tk)P(tk)ΦT(tk)

(19)

反饋增益矩陣為

K(tk+1)=P(tk+1/k)HT(tk+1)·

[H(tk+1)P(tk+1/k)H(tk+1)+Q(tk+1)]-1

(20)

估計(jì)均方誤差矩陣為

P(tk+1)=[I-K(tk+1)H(tk+1)]P(tk+1/k)

(21)

狀態(tài)估計(jì)值為

X(tk+1)=X(tk+1/k)+K(tk+1)·

(22)

3 仿真實(shí)例

3.1 仿真條件

重力測量衛(wèi)星質(zhì)心在軌標(biāo)定改進(jìn)算法的數(shù)學(xué)仿真參數(shù)，如表1所示。

為了使衛(wèi)星質(zhì)心在軌標(biāo)定改進(jìn)算法能夠得到合適的精度，分別在俯仰、偏航和滾動方向施加周期性的磁偶極矩Mm=30sin(0.2πt)A·m2，為了檢驗(yàn)所設(shè)計(jì)的磁偶極矩信號的持續(xù)時間與質(zhì)心標(biāo)定精度的關(guān)系，將信號的時間間隔分別設(shè)定為100s，200s，300s，400s。為了分析得到數(shù)據(jù)采集的軌道相位與質(zhì)心標(biāo)定精度的關(guān)系，分別選在赤道附近及南極附近進(jìn)行衛(wèi)星質(zhì)心在軌標(biāo)定。

表1 仿真參數(shù)

3.2 仿真流程

衛(wèi)星質(zhì)心在軌標(biāo)定改進(jìn)算法仿真流程見圖3。具體步驟為：①當(dāng)衛(wèi)星運(yùn)行到地球赤道上空時，將磁力矩器作為姿態(tài)控制元件，分別在俯仰、偏航、滾動方向產(chǎn)生較大的周期性磁偶極矩，獲得的磁力矩作用于衛(wèi)星上，使衛(wèi)星進(jìn)行姿態(tài)機(jī)動，姿態(tài)機(jī)動完成后恢復(fù)正常工作姿態(tài)；②將姿態(tài)機(jī)動過程中星敏感器輸出的姿態(tài)四元數(shù)進(jìn)行基于預(yù)測濾波的四元數(shù)-角速度算法解算，獲得衛(wèi)星的角速度；③根據(jù)解算的角速度，采用預(yù)測濾波解算出衛(wèi)星的角加速度；④利用靜電加速度計(jì)數(shù)據(jù)與衛(wèi)星質(zhì)心位置的線性關(guān)系，帶入衛(wèi)星的角速度、角加速度，采用擴(kuò)展卡爾曼濾波進(jìn)行解算，即可獲得衛(wèi)星質(zhì)心在衛(wèi)星本體坐標(biāo)系下的位置。

圖3 仿真流程Fig.3 Simulation flow

3.3 數(shù)值仿真結(jié)果

在赤道上空進(jìn)行姿態(tài)機(jī)動，在衛(wèi)星滾動(X軸)、俯仰(Y軸)和偏航(Z軸)方向上的磁偶極矩如圖4所示。磁偶極矩是幅度為30 A·m2、頻率為0.1 Hz的正弦波，機(jī)動總時間為300 s，間隔為100 s。圖5為衛(wèi)星的磁力矩，三軸的磁力矩最大幅值分別為0.4 mN·m，0.75 mN·m，0.75 mN·m，形式仍為正弦波，頻率仍為0.1 Hz。采用基于預(yù)測濾波的四元數(shù)-角速度算法，得到的衛(wèi)星角速度估計(jì)值如圖6所示，可以看出，基于預(yù)測濾波的四元數(shù)-角速度算法具有很好的魯棒性。三軸角速度估計(jì)誤差的平均值低于3×10-8rad/s，標(biāo)準(zhǔn)差低于5×10-7rad/s，具有很好的估計(jì)精度，作為參考角速度完全符合要求。衛(wèi)星角加速度估計(jì)值如圖7所示。改進(jìn)算法將全部狀態(tài)信息定義為模型誤差，僅利用量測信息便能精確地估計(jì)出衛(wèi)星質(zhì)心標(biāo)定期間的角加速度。靜電加速度計(jì)的輸出結(jié)果如圖8所示，X軸方向的線性加速度近似效果相對較差，主要原因是X軸方向?yàn)樾l(wèi)星的飛行方向，大氣阻力在此方向上的量級遠(yuǎn)遠(yuǎn)大于其他兩個方向。

圖4 赤道上空時衛(wèi)星的磁偶極矩Fig.4 Satellite magnetic dipole moment over equator

圖5 赤道上空時衛(wèi)星的磁力矩Fig.5 Satellite magnetic torque over equator

圖6 赤道上空時衛(wèi)星的角速度估計(jì)值Fig.6 Estimation of angular rate over equator

圖8 赤道上空時靜電加速度計(jì)的輸出Fig.8 Electrostatic accelerometer output over equator

衛(wèi)星質(zhì)心在軌標(biāo)定精度如表2所示，根據(jù)美國GRACE衛(wèi)星的指標(biāo)，質(zhì)心三軸標(biāo)定精度須在0.05 mm以內(nèi)[8]。

由表2可以看出，X軸、Y軸和Z軸的最佳標(biāo)定精度分別為0.022 mm，0.005 mm，0.008 mm。因此，本文提出的質(zhì)心標(biāo)定改進(jìn)算法滿足標(biāo)定精度需求。在相同條件下，在赤道上空進(jìn)行姿態(tài)機(jī)動，質(zhì)心標(biāo)定結(jié)果要略優(yōu)于在南極上空機(jī)動的結(jié)果。其主要原因是：赤道上空激勵起的三軸角加速度較為均勻，這樣利于質(zhì)心標(biāo)定精度的提高；而南極上空激勵的三軸角加速度幅值差異較大。通過不同的激勵時長看，300～900 s的總機(jī)動時間是可以接受的，如果激勵時間再長，質(zhì)心測量精度反而降低。這是因?yàn)椋簶?biāo)定中假設(shè)的非重力加速度項(xiàng)在質(zhì)心標(biāo)定期間保持線性關(guān)系，而隨著標(biāo)定時間的延長，這種假設(shè)的誤差越來越大，反而不利于質(zhì)心標(biāo)定精度的提高。

表2 衛(wèi)星質(zhì)心在軌標(biāo)定精度

4 結(jié)束語

重力測量衛(wèi)星質(zhì)心在軌標(biāo)定的核心問題是，在缺乏直接敏感元件的情況下，如何高精度地確定衛(wèi)星角速度及角加速度信息。針對這一問題，本文提出了一種預(yù)測濾波與擴(kuò)展卡爾曼濾波相結(jié)合的質(zhì)心在軌標(biāo)定改進(jìn)算法。該算法通過設(shè)計(jì)估計(jì)角速度及角加速度的預(yù)測濾波與質(zhì)心的擴(kuò)展卡爾曼濾波，實(shí)現(xiàn)衛(wèi)星質(zhì)心的標(biāo)定，模型結(jié)構(gòu)比較簡單，且利用非線性預(yù)測濾波建模非常方便。仿真結(jié)果表明：在不使用陀螺或者陀螺失效的情況下，衛(wèi)星質(zhì)心在軌標(biāo)定仍具有較高的精度，因而改進(jìn)算法為重力測量衛(wèi)星質(zhì)心在軌標(biāo)定提供了一種新的解決思路。

References)

[1]周澤兵，白彥崢，祝竺，等. 衛(wèi)星重力測量中加速度計(jì)在軌參數(shù)校準(zhǔn)方法研究[J].中國空間科學(xué)技術(shù)，2009，29(6)：74-80

Zhou Zebing，Bai Yanzheng，Zhu Zhu，et al. In-orbit cali-bration methods of accelerometer parameters on satellite-borne gravimetry [J]. Chinese Space Science and Technology，2009，29(6)：74-80 (in Chinese)

[2]佘世剛，王鍇，周毅，等. 高精度星間微波測距技術(shù)[J].宇航學(xué)報(bào)，2006，27(3)：403-406

She Shigang， Wang Kai，Zhou Yi，et al. The technology of high accuracy inter-satellite microwave ranging [J]. Journal of Astronautics， 2006，27(3)：403-406 (in Chinese)

[3]陳光鋒，唐富榮，薛大同．重力衛(wèi)星在軌質(zhì)心修正原理[J]．宇航學(xué)報(bào)，2005，26(5):567-570

Chen Guangfeng，Tang Furong，Xue Datong. The trim principle of center of mass of gravity satellite during orbit flight [J]. Journal of Astronautics，2005，26(5)：567-570 (in Chinese)

[4]Wang F, Tapley B D. Precise determination of accelero-meter proof mass [C]//Proceedings of AAS/AIAA Space Flight Mechanics Meeting. Washington D.C.：AIAA，2002：739-752

[5]王本利,廖鶴,韓毅. 基于MME/EKF算法的衛(wèi)星質(zhì)心在軌標(biāo)定[J].宇航學(xué)報(bào)，2010，31(9)： 2150-2156

Wang Benli, Liao He, Han Yi. On-orbit calibration of satellite center of mass based on MME/EKF algorithm [J]. Journal of Astronautics，2010，31(9)：2150-2156 (in Chinese)

[6]Tangyin S, Williams T. Mass property estimation using coasting maneuvers [J]. Journal of Guidance, Control and Dynamic, 1997，20(4)：625-632

[7]Crassidis J L. Predictive filtering for attitude estimation without rate sensors [J]. Journal of Guidance， Control and Dynamic， 1997，20(3):522-527

[8]Kim J R. Simulation study of a low-low satellite-to-satellite tracking mission [D].Austin: University of Texas at Austin, 2000

[9]Furun Wang. Study on center of mass calibration and K-band ranging system calibration of the GRACE mission [D]. Austin: University of Texas at Austin, 2003

(編輯：夏光)

An Improved Algorithm for On-orbit Calibration of Gravity Measurement Satellite Center of Mass

XIN Ning1QIU Lede1ZHANG Lihua2LIU Naijin1

(1 Institute of Telecommunication Satellite， China Academy of Space Technology， Beijing 100094， China) (2 DFH Satellite Co. Ltd.， Beijing 100094， China)

In view of low precision and vulnerability of gyros, an improved algorithm for on-orbit calibration of gravity measurement satellite center of mass based on star sensor and electrostatic accelerometer is proposed．A periodic magnetic torque which is larger than other external torques is applied on satellite．The observation of star sensor is used to estimate angular rate and acceleration based on predictive filter．The extended Kalman filter is used to estimate the center of mass of satellite by using electrostatic accelerometer data．The numerical simulation results indicate that the angular acceleration and center of mass can be real-time estimated synchronously， and the calibration accuracy of three axes is better than 0.05mm, which can demonstrate the effectiveness of the proposed algorithm even without gyros．

gravity measurement satellite； on-orbit calibration of center of mass； predictive filter； extended Kalman filter； star sensor

2014-11-14；

2015-01-28

國家自然科學(xué)基金(91438205)

辛寧，男，博士，研究方向?yàn)樾l(wèi)星系統(tǒng)總體設(shè)計(jì)。Email：xinning7@sina.com。

P228

A

10.3969/j.issn.1673-8748.2015.04.008