地球靜止軌道CW方程建模誤差分析

任家棟， 張大力， 曾慶雙

(哈爾濱工業(yè)大學(xué) 航天學(xué)院，黑龍江 哈爾濱 150001)

在目前的航天器相對運動框架中，動力學(xué)模型呈多樣化趨勢，具體取決于兩星的距離、運行軌道的扁率、運動遞推的運算復(fù)雜度等因素，難以用單一標(biāo)準(zhǔn)來明確某個動力學(xué)模型更加適合或更加有效[1]。相對運動描述方法分為相對軌道要素描述法和直角坐標(biāo)描述法2類[2]。本文主要針對相對導(dǎo)航制導(dǎo)設(shè)計應(yīng)用場景，如二階段濾波、非線性預(yù)測估計等濾波方法需要明確的干擾誤差特性，基于此，對直角坐標(biāo)系描述方法進(jìn)行分析。

多年來，研究人員在CW方程[3]的基礎(chǔ)上，從近圓和近距離等理想化假設(shè)條件下出發(fā)，向保留高階項、考慮軌道攝動、橢圓軌道等幾方面進(jìn)行了擴展研究。文獻(xiàn)[4-5]分別研究了圓軌道假設(shè)下二階、三階項對相對運動的影響。Lawden[6]以真近點角為參變量，給出了橢圓軌道下的相對動力學(xué)解析描述，稱為Lawden方程，已經(jīng)廣泛應(yīng)用于橢圓軌道的相對運動控制領(lǐng)域。文獻(xiàn)[7-8]采用逐次逼近的QV(Quadratic-Volterra solution)方法研究了橢圓軌道下考慮高階項的分析模型。

軌道環(huán)境攝動，如地球非球形引力攝動、大氣阻力攝動、太陽光壓攝動等均會導(dǎo)致期望軌跡產(chǎn)生偏差，是影響編隊飛行的關(guān)鍵因素?；贑W方程，文獻(xiàn)[9-10]提出了一種地球非球形引力攝動的補償模型，提升了相對運動精度。文獻(xiàn)[11]將變化幅度較小的地心距及軌道傾角平均化，推導(dǎo)了考慮J2項攝動干擾的近地軌道相對運動解析方程。文獻(xiàn)[12]采用分析法研究了三體引力對衛(wèi)星相對運動的影響，三體引力主要對對軌道面內(nèi)運動產(chǎn)生長周期作用。太陽光壓是高軌衛(wèi)星的主要攝動元素之一，對于高面質(zhì)比(high area-to-mass ratio,HAMR)目標(biāo)其影響甚至占主導(dǎo)地位。太陽光壓攝動一般考慮為非保守力，但在全軌道考慮時會呈現(xiàn)一些近似保守力的影響特性，如平均軌道半長軸不變、軌道偏心率或傾角以年或年的倍數(shù)周期進(jìn)行波動[13-14]。

目前對CW方程的各誤差要素的影響分析成果較多，但缺少地球靜止軌道小傾角、小偏心率等特性下針對相對導(dǎo)航制導(dǎo)應(yīng)用場景的分析，本文結(jié)合地球靜止軌道的特點，推導(dǎo)CW方程的建模誤差，并對其主要誤差項進(jìn)行分析，以得到誤差的規(guī)律和量化特性。

1 高軌環(huán)境分析與建模

衛(wèi)星在軌道上始終受到各種攝動力的作用，主要包括：地球非球形攝動，高層大氣的氣動力攝動，太陽、月亮的引力攝動，太陽光輻射壓力攝動等。在攝動力的作用下，衛(wèi)星軌道偏心率e、周期、升交點赤經(jīng)Ω和軌道傾角i等不斷演化，呈現(xiàn)周期或非周期特性。雖然攝動力量級通常較小，但其長期影響不容忽視，需對各攝動的特性進(jìn)行詳細(xì)數(shù)值分析。

由于大氣阻力隨軌道高度的增加而急劇減小，地球靜止軌道攝動主要考慮地球非球形引力攝動、日月引力攝動(也稱三體引力攝動)、太陽光壓攝動，其量級特性見表1[15]。

表1 地球靜止軌道各攝動力與理想地球引力比值

下面以柯氏定理為基礎(chǔ)對相對運動模型進(jìn)行推導(dǎo)。慣性系下的矢量在轉(zhuǎn)動坐標(biāo)系下求導(dǎo)為:

(1)

目標(biāo)星(target)和追蹤星(chaser)的絕對動力學(xué)方程分別為：

(2)

式中：rt、rc分別為目標(biāo)星和追蹤星在慣性坐標(biāo)系下的位置矢量；at、ac分別為目標(biāo)星和追蹤星的環(huán)境攝動加速度。

圖1 地心軌道坐標(biāo)系Fig.1 Geocentric coordinate system

相對運動通常在衛(wèi)星軌道坐標(biāo)系下研究，軌道坐標(biāo)系采用VVLH坐標(biāo)系，定義為OoXoYoZo。

軌道坐標(biāo)系OoXoYoZo不是慣性系，因此在軌道坐標(biāo)系下推導(dǎo)相對運動方程需考慮OoXoYoZo的轉(zhuǎn)動信息。為了便于研究其轉(zhuǎn)動特性，建立地心軌道系OXYZ，原點在地心，相當(dāng)于軌道坐標(biāo)系平移到地球球心。因此地心軌道系OXYZ和星體軌道系OoXoYoZo有同樣的轉(zhuǎn)動特性。

(3)

顯然，ωx=0，即軌道坐標(biāo)系不存在X軸方向的轉(zhuǎn)動分量。

根據(jù)柯氏定理求二階導(dǎo)：

(4)

式中(Fx,Fy,Fz)為軌道攝動加速度。解式(4)，并結(jié)合兩體運動軌道動力學(xué)可得瞬時軌道根數(shù)表示的瞬時軌道角速度和角加速度：

(5)

式中：a為軌道半長軸；ec為偏心率；f為緯度幅角；Fx/r是軌道攝動力所帶來的。

根據(jù)苛氏定理，在追蹤星軌道坐標(biāo)系下的兩星相對運動方程表示為：

(6)

(7)

式中：O(n2)為地球球形引力場的2階以上高階項引力差；Δatc為軌道環(huán)境攝動力差。式(7)可改寫為：

(8)

若忽略ΔF的影響，上述方程是線性方程，即著名的Hill方程或CW方程，該方程被廣泛用于交會對接動力學(xué)分析與控制問題中。為了得到CW方程誤差的量化特性，有必要對ΔF的各組成部分進(jìn)行詳細(xì)研究。

2 相對運動模型誤差推導(dǎo)

分析可知，CW方程是近似方程，其誤差體現(xiàn)在加速度項ΔF。將ΔF分解為2部分：一部分是地球球形引力場的高階項和環(huán)境攝動力差等未建模作用力ΔF1；另一部分是圓軌道假設(shè)下，由于實際偏心率的存在等帶來CW方程的未建模作用力ΔF2。

2.1 圓軌道假設(shè)下的建模誤差

首先，對地球球形引力場高階項和環(huán)境攝動力差等未建模作用力ΔF1進(jìn)行量化分析，該部分誤差可改寫為：

ΔF1=O(n2)+Δatc

(9)

式中：Δatc表示空間環(huán)境中各攝動力的綜合作用；O(n2)表示地球球形引力場兩星作用力差的二階及其高階項，若只考慮二階項，有：

(10)

其次，對于日月引力對相對運動的作用力參見式(2)。日月引力對兩星作用力的差消去共同部分，作用形式與地球球形引力場對兩星作用力的形式相同。同相對動力學(xué)中地球球形引力場的一階近似過程(7)，在太陽或月球中心慣性坐標(biāo)系下，對兩星的日月引力差取一階形式為：

(11)

式中：mt、mc表示在月球中心慣性坐標(biāo)系下目標(biāo)星和追蹤星的位置矢量；st、sc表示在太陽中心慣性坐標(biāo)系下目標(biāo)星和追蹤星位置矢量；sx、sy、sz表示衛(wèi)星在太陽軌道系下的相對位置坐標(biāo)；mx、my、mz表示衛(wèi)星在月球軌道系下的相對位置坐標(biāo)。

進(jìn)一步的，由表1可見，地球非球形引力攝動中，各攝動項形式相似，帶諧項J2項引力勢能作用力最大，較其他項高約3個數(shù)量級。在地心慣性系下，J2攝動作用力表示為：

(12)

略去r-7以上高階項，式(12)可改寫為：

(13)

(14)

傾角i對于地球靜止軌道衛(wèi)星正常情況下是小量，小于0.1°，進(jìn)行一階展開得：

y=-sinisinΩ(xt-xc)-sinicosΩ(yt-yc)

-cosi(zt-zc)≈-(zt-zc)

(15)

將式(15)代入式(14)，J2項的攝動差在追蹤星軌道坐標(biāo)系可表示為：

(16)

綜上所述，地球球形引力場的高階項和環(huán)境攝動力差等未建模作用力ΔF1可統(tǒng)一表示為：

(17)

2.2 非圓軌道的相關(guān)誤差

由圓軌道假設(shè)等帶來CW方程的未建模誤差ΔF2表達(dá)式：

(18)

其中：

(19)

(20)

ΔF2表達(dá)式中有相對速度項，為了獲得相對位置表示的統(tǒng)一表達(dá)式。將速度項用相對位置表示。

在追蹤星軌道系下，目標(biāo)星相對運動的一階閉環(huán)形式，那么垂直于軌道平面運動，即軌道法向運動可以表示為：

(21)

易見，該方向運動為正弦周期振蕩，振蕩角速度為追蹤星瞬時軌道角速度。軌道面內(nèi)是橢圓運動，可表示為：

(22)

令Δl為一軌時間內(nèi)兩星距離變化量，zc為橢圓運動Z軸的中心坐標(biāo)，表示為：

(23)

因此，X軸相對速度可表示為：

(24)

將式(5)、(21)、(24)代入ΔF2中，忽略兩階小量得軌道面法向力攝動為：

(25)

2.3 誤差綜合

統(tǒng)計各誤差的作用效果如表2，包含CW方程推導(dǎo)過程中引入的主要誤差源。地球球形引力場線性化中二次項誤差與衛(wèi)星高度的四次方成反比，與兩星的距離平方成正比；地球田諧項攝動可以忽略；帶諧項攝動中J2項誤差與衛(wèi)星高度五次方成反比，與兩星的距離成正比；日月引力的攝動作用與衛(wèi)星高度相關(guān)性弱，與兩星的距離成正比；圓軌道假設(shè)中軌道偏心率相關(guān)誤差與衛(wèi)星高度的三次方成反比，與兩星的距離成正比；軌道面法向力攝動通過影響軌道坐標(biāo)系的角速度耦合相對運動引入誤差，與衛(wèi)星高度的一次方成反比，與兩星的距離和一軌漂移量成正比。

表2 軌道攝動力帶來的誤差Table 2 The error of the orbital perturbation force

綜上所述，CW方程的誤差一方面與衛(wèi)星高度有關(guān)，會隨衛(wèi)星高度的增加急速變??；另一方面與兩星的距離有關(guān)，包含線性相關(guān)和平方相關(guān)項。

3 地球靜止軌道相對運動模型誤差特性分析

本節(jié)對地球靜止軌道相對運動模型誤差進(jìn)行特性分析，得到其誤差的量級以及與運動構(gòu)型的關(guān)系，為后面章節(jié)的相對導(dǎo)航設(shè)計提供依據(jù)。

前節(jié)分析可得：日月引力攝動、軌道面法向力攝動、J2攝動量級相同，且都與兩星的距離成正比；地球球形引力場二次項在兩星距離小于1 km時與J2等攝動同階，在兩星距離大于1 km時攝動力隨距離平方增加；偏心率相關(guān)攝動大小與軌道偏心率成正比，在衛(wèi)星軌道偏心率小于1×10-5時攝動力與J2等攝動同階。CW方程在高軌誤差的作用形式變得復(fù)雜，相關(guān)因素較多，沒有明確的主要誤差，不易補償。

但對于典型的在軌飛行軌道，通常地球靜止軌道衛(wèi)星軌道偏心率小于0.001，兩星相對距離大于1 km，因此地球球形引力場的二次項和偏心率相關(guān)攝動可以作為CW方程的主要誤差源，用以分析CW方程誤差的量級和作用形式。

利用軌道仿真軟件STK的HPOP軌道輸出，對主要攝動和總攝動的特性進(jìn)行數(shù)值分析。主要攝動項包括地球引力二次項和偏心率相關(guān)攝動，總攝動在動力學(xué)模型中按式(8)右端引出。

分析條件為：追蹤星偏心率0.001，兩星相距10 km，伴飛橢圓1 km。主要攝動項的量化特性見圖2，總攝動項的量化特性見圖3。

圖2 地球靜止軌道主要攝動累加Fig.2 The main perturbation of the geostationary orbit

圖3 地球靜止軌道CW方程各軸總攝動Fig.3 The total perturbation of each axis of the CW equation of geostationary orbit

仿真可見，上述主要攝動項在各類誤差中占主要部分，與總誤差量級相同，且振蕩和常偏特性類似，能夠反映地球靜止軌道相對運動模型CW方程的誤差特性。圖3中的波動部分是目標(biāo)進(jìn)出陰影過程中的太陽光壓攝動差。

下面對地球引力二次項和偏心率相關(guān)攝動項的量化特性進(jìn)行獨立分析。圖4是對表2中，地球球形引力二階項在兩星10 km穩(wěn)定伴飛工況下各軸作用力的量化特性。相對運動橢圓和軌道面外振蕩控制的比較小，此時地球球形引力場二階項在軌道面法向Y軸作用力很小，該項與兩星的相對距離沒有關(guān)系；地球球形引力場的二階項對X軸和Z軸的攝動作用與兩星相對距離和伴飛橢圓尺寸均相關(guān)，由于相對運動在Z軸為周期運動，X軸攝動項表現(xiàn)出周期振蕩；Z軸的攝動作用主要特點是常偏，在穩(wěn)定伴飛期間，由于Y軸振蕩小，Z項的攝動表現(xiàn)為跟兩星距離平方成正比的常偏和一個與伴飛橢圓相關(guān)的小幅振蕩的綜合。

圖4 地球球形引力場二階項攝動各軸分量Fig.4 Second order perturbation of the earth′s spherical gravitational field perturbs each axial component

圖5是對表2中偏心率相關(guān)項在兩星10 km穩(wěn)定伴飛工況下各軸作用力的仿真輸出，呈現(xiàn)周期性的振蕩特性。從相對運動方程CW方程的穩(wěn)態(tài)解析式(22)可得，各軸相對位置呈現(xiàn)周期性變化，偏心率相關(guān)項各軸攝動耦合因素多，偏心率和相對運動橢圓尺寸聯(lián)合影響振蕩的幅值，相對運動的當(dāng)前位置和當(dāng)前軌道幅角等綜合影響振蕩的相位。

圖5 偏心率相關(guān)攝動各軸分量Fig.5 Eccentricity correlative perturbation of each axis component

綜合分析表明，地球球形引力二次項和偏心率相關(guān)攝動項反映了CW方程在地球靜止軌道的主要誤差特性，可以在相對導(dǎo)航濾波設(shè)計中用以定量分析狀態(tài)方程的誤差情況，或設(shè)計干擾估計模型。地球球形引力場二次項的誤差在各軸呈現(xiàn)周期振蕩特性，并隨兩星距離變大，在Z軸產(chǎn)生常值偏差分量；偏心率相關(guān)攝動與兩星相對構(gòu)型線性相關(guān)，各軸均呈現(xiàn)為周期振蕩特性。

5 結(jié)論

1)考慮地球球形引力場高階項、地球非球形引力攝動、三體引力攝動和太陽光壓攝動，分析了地球靜止軌道環(huán)境中各攝動力的對相對運動的作用量級，對兩星相對運動方程進(jìn)行了建模分析，得到了完整的誤差要素。

2)為研究CW方程的誤差特性，對加速度項進(jìn)行深入分析，給出了各誤差源的近似表達(dá)式、作用量級和影響因素，分析表明，CW方程的誤差隨衛(wèi)星高度的增加急速變小，并且與相對距離存在線性和平方復(fù)合的相關(guān)特性。

3)對典型的相對運動型態(tài)進(jìn)行了數(shù)值仿真，結(jié)果表明地球球形引力二次項、偏心率相關(guān)項以及J2項相關(guān)誤差為CW方程的主要誤差來源，在軌道系坐標(biāo)系各軸向具體呈現(xiàn)為常值偏差和軌道周期振蕩組合的分布特性。