黃美麗，王大軼，馮 昊，田百義

(北京空間飛行器總體設(shè)計部，北京 100094)

0 引 言

對于低軌衛(wèi)星，軌道高度受大氣影響將不斷降低，并且軌道高度越低，大氣密度越大，軌道衰減速度越快[1]。因此，衛(wèi)星需要定期消耗推進劑進行軌道控制以維持軌道高度在一定的范圍內(nèi)。軌道維持的推進劑預算是衛(wèi)星軌道設(shè)計的一個重要組成部分。

衛(wèi)星軌道維持的推進劑消耗量由大氣密度及衛(wèi)星迎風面質(zhì)比決定，衛(wèi)星總體構(gòu)型以及軌道確定后，迎風面質(zhì)比也就確定了，軌道維持的推進劑消耗量取決于大氣密度[2]。

大氣密度是一個復雜的多元變量函數(shù)，因此預測難度極大[3]。幾十年來，大氣模型不斷在改進和發(fā)展，但由于高層大氣的變化極其復雜，模型仍存在約15%～30%的誤差，在空間環(huán)境擾動期間甚至可達100%或更高。因此，進一步提高大氣模型預報精度成為大氣模型研究需要解決的重要問題。目前，大氣模型預報精度的改進主要有兩種方式：一是建立新模型，增加更多更新的探測數(shù)據(jù)或者對模型的控制變量進行改進；二是在原模型基礎(chǔ)上進行結(jié)果修正。無論是哪種改進方法，其基礎(chǔ)都是要獲取大量準確的大氣密度觀測數(shù)據(jù)[4]。

衛(wèi)星的兩行根數(shù)(TLE)、我國在軌衛(wèi)星定軌數(shù)據(jù)，都是實測軌道數(shù)據(jù)，將一段時間內(nèi)的軌道數(shù)據(jù)刨除變軌因素，衛(wèi)星軌道參數(shù)的變化反映了衛(wèi)星所受攝動力的情況[5]。把其中與大氣阻力無關(guān)的攝動因素去除之后，得到了衛(wèi)星受大氣阻力攝動的軌道變化情況，再結(jié)合衛(wèi)星的外形結(jié)構(gòu)、阻力系數(shù)、姿態(tài)變化等因素，可對大氣密度進行反演，得到近實時的大氣密度。

近年來，國內(nèi)外在由TLE數(shù)據(jù)進行大氣密度反演以及在此基礎(chǔ)上的大氣密度模型修正等領(lǐng)域進行了不少研究，但多注重于反演方法的推導及個例的應用，缺乏大量數(shù)據(jù)反演獲得中長期平均大氣密度的研究，并且尚無將反演方法應用于衛(wèi)星推進劑預算的研究[6-12]。

本文針對我國典型遙感衛(wèi)星運行軌道開展中長期的平均大氣密度反演研究，并利用反演結(jié)果進行大氣密度局部修正以及推進劑預算應用研究，達到提高推進劑預算準確度、在保障衛(wèi)星安全運行基礎(chǔ)上最大限度減少推進劑預算的目的。

1 中長期大氣密度反演算法研究

1.1 低軌衛(wèi)星軌道高度衰減分析

對于低軌衛(wèi)星，大氣阻力造成軌道機械能的長期衰減。對于一般的低地球軌道，大氣阻力產(chǎn)生與衛(wèi)星速度相反的阻力，對應的阻力加速度大小為：

(1)

式中：fD為大氣阻力加速度，單位km/d2；CD為阻力系數(shù)，對于一般衛(wèi)星其取值范圍為2.0～2.3；S為衛(wèi)星有效迎風面積，單位km2；m為衛(wèi)星質(zhì)量，單位kg；ρ為衛(wèi)星飛行高度處的大氣密度，單位kg/km3；V為衛(wèi)星速度，單位km/d。

大氣阻力沿速度的反向，將造成軌道面內(nèi)參數(shù)半長軸a，偏心率e,近地點幅角ω和平近點角M的變化，對應的攝動方程如下：

(2)

式中：θ為衛(wèi)星真近點角，單位rad；E為衛(wèi)星偏近點角，單位rad；r為衛(wèi)星地心距，單位km；p=a(1-e2)為軌道半通徑，單位km；n為衛(wèi)星軌道運動平均角速度，單位rad/d；R,T為大氣阻力加速度的徑向和橫向分量。

(3)

對于橢圓的低地球軌道，其高度的衰減涉及了半長軸、遠地點高度、近地點高度和偏心率的衰減，由于衛(wèi)星在橢圓軌道上不同位置的大氣密度也是變化的，導致軌道衰減方程(2)和近似公式(3)不再適用，橢圓軌道各參數(shù)的衰減情況可參考文獻[12]。

1.2 大氣密度反演算法分析

1)積分反演算法

大氣密度反演基于實測的衛(wèi)星軌道平根數(shù)變化，采用帶有大氣攝動的軌道動力學模型(見式(2)),對衛(wèi)星所在高度中長期的等效大氣密度反演，給出大氣密度反演算法理論和流程。該方法適用于近圓軌道中長期的大氣密度反演。

衛(wèi)星軌道數(shù)據(jù)包括地心J2000坐標系下的軌道瞬根和平根數(shù)，軌道參數(shù)包含軌道半長軸、偏心率、傾角、升交點赤經(jīng)、近地點幅角和平近點角。平根數(shù)中的軌道半長軸消去了地球非球形攝動的短周期和長周期攝動項，主要受大氣阻力攝動影響。主要表現(xiàn)是半長軸在軌逐漸衰減，可根據(jù)半長軸衰減率反推衛(wèi)星當時所在高度的等效大氣密度情況?；谲壍烂鎯?nèi)參數(shù)攝動方程，大氣密度反演方案可按圖1進行。

詳細算法流程描述如下：

(2)對攝動方程(2)積分，可得到衛(wèi)星運行ΔT時間后的平軌道半長軸at。

2)解析反演算法

(4)

以495 km軌道為例，采用解析反演算法得到的大氣密度和積分反演算法相比，誤差約0.05%，兩種方法相差極小。且解析反演算法與積分反演算法相比，可極大節(jié)約計算時間，因此，本文采用解析算法進行大氣密度的反演。

上述方法用于反演衛(wèi)星中長期的大氣密度情況，需對長期實測的軌道數(shù)據(jù)進行簡單處理，獲取衛(wèi)星在某段時期內(nèi)的平均衰減率，以便對該段時間的平均大氣密度進行反演。

2 典型軌道中長期大氣密度反演分析

2.1 典型軌道選取

迄今為止，我國已成功發(fā)射了六發(fā)遙感衛(wèi)星(ZY-2(01)、ZY-2(02)、ZY-2(03)、YAOGAN5、YAOGAN12和YAOGAN21),運行于495km軌道上，并且有多發(fā)衛(wèi)星在495km附近的軌道上運行。未來，我國還將往該軌道附近陸續(xù)發(fā)射多顆遙感衛(wèi)星。該軌道上持續(xù)有衛(wèi)星運行已超過11年的一個太陽活動周期，進行大氣密度反演對后續(xù)的軌道維持推進劑預算是非常有參考意義的。

因此，本文選擇典型衛(wèi)星高度495km進行大氣密度反演，并與MSIS大氣密度模式計算的大氣密度進行對比，驗證本文的反演模型和算法。大氣密度計算主要取決于太陽活動率，太陽活動率以波長10.7的射電流量F10.7表示。MSIS00模式計算大氣密度時采用實測的F10.7數(shù)據(jù)，根據(jù)美國國家海洋和大氣管理局(NOAA)發(fā)布的實測F10.7數(shù)據(jù)，2000年至2016年，實測F10.7數(shù)據(jù)如圖2所示。

大氣密度主要依賴于太陽活動水平F10.7，F(xiàn)10.7既有約為11年的長周期變化，又有周期約為27天的短周期變化。對于衛(wèi)星推進劑預算，只需考慮長周期變化即可。

2.2 中長期大氣密度反演結(jié)果

以ZY-2(01)、ZY-2(02)、ZY-2(03)、YAOGAN5、YAOGAN12和YAOGAN21衛(wèi)星為例，衛(wèi)星平均軌道高度均為495km。

以YAOGAN12和YAOGAN21衛(wèi)星為例，2015年期間兩顆星軌道半長軸衰減情況如圖3所示。

由圖3可知，兩顆星在2015年期間進行了5次正常的軌道維持。雙星軌道高度、整星質(zhì)量和有效迎風面積較為一致，因此，雙星軌道半長軸衰減也較為一致。

大氣密度反演過程中，ΔT的取值應小于衛(wèi)星軌道維持周期。ΔT取值越大，所得大氣密度越接近ΔT時間內(nèi)的平均密度，對獲取每天的大氣密度和F10.7的關(guān)系不利；而取值越小，對抑制軌道半長軸誤差不利，所得大氣密度受軌道定軌誤差影響較大。綜上考慮，495 km軌道大氣密度反演時，取ΔT=5天。

基于實測F10.7計算所得的MSIS00大氣模式密度和本文反演算法所得的2010年至2016年反演密度情況見圖4。

本文方法所得大氣密度與MSIS大氣密度模式相比，得到如下結(jié)論：

1)MSIS00大氣密度模式與反演大氣密度的變化趨勢吻合度較高。

2)在太陽活動高年，MSIS00模型誤差較大，誤差最大可達60%。

3)在太陽活動低年，MSIS00模型誤差較小，模型精度可達10%。

3 大氣密度模型局部修正方法研究

3.1 大氣密度模型局部修正方法

本節(jié)以MSIS00大氣密度模式為例，進行大氣模型局部修正方法研究。

通過在軌數(shù)據(jù)進行大氣密度反演，獲得一段時間內(nèi)的大氣密度數(shù)據(jù)：

ρs=ρs(t)

(5)

通過大氣密度模型計算一段時間內(nèi)的大氣密度數(shù)據(jù)，計算大氣密度時采用實測的F10.7數(shù)據(jù)。

ρm=ρm(t,F10.7)

(6)

從而得到同一時間的實測大氣密度與模型密度之間的關(guān)系系數(shù)

(7)

利用衛(wèi)星長期在軌的數(shù)據(jù)，通過多項式曲線擬合，得到兩種密度之間的關(guān)系系數(shù)ρx與F10.7的關(guān)系：

ρx(F10.7)=A·f(F10.7)

(8)

式中：A為多項式系數(shù)向量。

通過上述關(guān)系系數(shù)與F10.7的擬合表達式，可以對大氣密度模型進行修正，從而再用于工程設(shè)計中，優(yōu)化衛(wèi)星推進劑預算及總體設(shè)計。

3.2 典型軌道大氣密度模型局部修正

對典型軌道高度495 km，通過式(5)～(8)可獲取大氣密度模型修正系數(shù)同F(xiàn)10.7的關(guān)系如圖5所示。

采用5次多項式擬合得出模型修正系數(shù)同F(xiàn)10.7的關(guān)系見下式：

(9)

其中多項式系數(shù)向量A取值如下：

A= [5.538039×10-11,-5.496297×10-8,

2.098230×10-5,-0.003817,0.322250,

-8.822965]

(10)

模型修正前后和反演的大氣密度情況見圖6。

由圖6可知：

1)在太陽活動高年，模型修正后的大氣密度平均值同反演的大氣密度平均值相比，最大誤差約20%，與模型修正前的50%誤差相比，模型精度有了較大提高。

2)在太陽活動低年，修正前后的大氣密度模型同反演的大氣密度均有較高的契合度。

3)通過系數(shù)修正主要提高了太陽活動高年的大氣密度模型精度。

4 大氣密度模型局部修正在推進劑預算中的應用

4.1 低軌衛(wèi)星推進劑預算

低軌衛(wèi)星推進劑預算主要包括對初軌調(diào)整、軌道機動、軌道維持以及姿態(tài)控制等消耗的推進劑進行預算并留有一定的余量。其中，對大氣阻力引起的軌道衰減進行軌道維持所需的推進劑預算是重點及難點，這主要因為其他幾項推進劑預算一般情況下均可準確預估，而軌道維持所需推進劑則與大氣密度預測密切相關(guān)，而大氣密度預測的不確定性很大，是目前整體推進劑預算的難點。

傳統(tǒng)的軌道保持推進劑預算方法主要有兩種：

1)對于壽命3年左右的衛(wèi)星，一般都按太陽活動中高年(F10.7=175)的情況，根據(jù)大氣密度模型預報的大氣密度進行推進劑預算，推進劑余量較大。

2)對于壽命5年以上的衛(wèi)星，一般按一個太陽活動周期中高、中、低年組合的形式，根據(jù)預計運行年份，以大氣密度模型預報的大氣密度進行推進劑預算，余量相對較小。

可以看出，兩種方法都受限于大氣密度模型預報精度。對于典型軌道，在具備長期(一個太陽活動周期11年以上)大氣密度反演數(shù)據(jù)的前提下，可利用反演結(jié)果對大氣密度模型進行局部修正，提高大氣模型預報精度，從而提高推進劑預算精度。

4.2 應用算例及結(jié)果分析

本文以495 km軌道高度上的5年壽命衛(wèi)星推進劑預算作為應用算例進行大氣密度模型局部修正方法應用研究。

根據(jù)前文分析，對于5年壽命衛(wèi)星，按太陽活動中高年的情況，根據(jù)大氣密度模型預報的大氣密度進行推進劑預算。取衛(wèi)星的迎風面質(zhì)比為0.008 m2/kg，分析可得大氣密度模型局部修正前后的軌道維持推進劑預算比對如表1所示[10]。

由表1可知，軌道保持推進劑消耗由修正前的148 kg減少為110 kg，比原設(shè)計結(jié)果節(jié)省38 kg推進劑，節(jié)省了近35%。

表1 大氣密度模型局部修正前后軌道維持推進劑預算比對Table 1 Propellant budget before modified model and after

5 結(jié)束語

本文基于實測軌道數(shù)據(jù)進行大氣密度反演方法、大氣模型局部修正方法及工程應用研究，研究成果具有創(chuàng)新性以及工程實用性，可在我國未來航天器的高精度軌道預報、推進劑預算、碰撞規(guī)避等任務中得以推廣和應用，可以為我國未來航天器的總體設(shè)計等提供技術(shù)支持，促進形成合理、可行、優(yōu)化的設(shè)計方案。

需要指出的是，本文采用的面質(zhì)比計算方法適用于合作目標，對于非合作目標，測軌參數(shù)里包含的B*參數(shù)是一個綜合了阻力系數(shù)、面質(zhì)比和大氣密度的綜合氣動參數(shù)，后續(xù)如何根據(jù)該參數(shù)獲取大氣密度將進行深入研究。