劉衛(wèi),劉四清,龔建村,王榮蘭

(1.中國科學(xué)院國家空間科學(xué)中心,北京100190；2.中國科學(xué)院微小衛(wèi)星創(chuàng)新研究院,上海201210；3.中國科學(xué)院大學(xué),北京100049)

1 引言

自1961年美國子午儀-4A衛(wèi)星在軌爆炸以來，可觀測空間碎片數(shù)量出現(xiàn)大規(guī)模增長[1]。1978年1月，入軌僅3個(gè)月載有核反應(yīng)堆的宇宙954衛(wèi)星再入時(shí)，在加拿大上空失控，沿約1000km的區(qū)域釋放了約30kg核輻射物[2]。1996年7月，法國櫻桃衛(wèi)星被1986年11月阿里安火箭爆炸的碎片撞擊損壞[3]；2009年在用衛(wèi)星美國銥星33和俄羅斯宇宙2251發(fā)生碰撞，產(chǎn)生可追蹤碎片超過2000塊[4，5]；1989年3月的太陽大爆發(fā)造成了大量太空編目目標(biāo)的丟失，2012年7月美國空間監(jiān)視網(wǎng) (SSN)更是創(chuàng)下目標(biāo)丟失列表紀(jì)錄[6，7]，原因也正是太陽活動峰年的熱層大氣模式誤差。近年來空間目標(biāo)頻繁發(fā)生在軌碰撞，隕落進(jìn)入地球大氣燒毀，甚至砸至地面，對人類的空間資產(chǎn)和地面安全造成威脅。熱層密度預(yù)測值的不確定性造成的空間碎片的軌道衰減，對碰撞發(fā)生時(shí)間和位置的預(yù)測有很大影響[1]。碰撞預(yù)警計(jì)算中漏警、虛警時(shí)有發(fā)生，隕落預(yù)報(bào)中誤差較大等這些問題都是熱層大氣密度預(yù)測誤差直接或間接造成的。20世紀(jì)50年代末，第一顆人造衛(wèi)星上天后觀測到大氣阻力對衛(wèi)星軌道運(yùn)動的影響?；谛l(wèi)星軌道衰變數(shù)據(jù)，Jacchia等人在20世紀(jì)60～70年代開發(fā)出第一代通用經(jīng)驗(yàn)大氣模型，如HarrisPriester、Jacchia模型[8]。經(jīng)過數(shù)十年發(fā)展，目前已有Jacchia、DTM和MSIS等系列熱層大氣密度模式，這些模式平均而言，均存在15%～20%的誤差，空間環(huán)境擾動期間誤差甚至超過100%[9，10]，無法滿足空間態(tài)勢感知等高精度的軍事需求。

Travis等人利用CHAMP和GRACE精密星歷(POE)估算大氣密度，評估大擾動期間的密度變化，發(fā)現(xiàn)使用POE推導(dǎo)密度與加速度計(jì)反演密度預(yù)測的軌道之間的誤差很小。因此Travis認(rèn)為加速度計(jì)觀測到的熱層密度高頻變化項(xiàng)，對軌道分析應(yīng)用都不重要[11]。大氣阻力可以從近地軌道衛(wèi)星的軌道信息中推斷出來，為測量質(zhì)量密度提供了一種直接的手段。CHEN等人基于半長軸的阻力攝動方程，使用CHAMP衛(wèi)星反演密度，發(fā)現(xiàn)應(yīng)用反演密度可以減少軌道預(yù)測誤差[12]。Rodney等人通過2003～2007年期間CHAMP衛(wèi)星加速度計(jì)數(shù)據(jù)，驗(yàn)證密度模型延時(shí)誤差的影響，并分析延遲對軌道預(yù)測結(jié)果的影響[13]。

在熱層大氣模型性能評測中，通常給出相對誤差標(biāo)準(zhǔn)差，如15%～20%，而忽略誤差均值問題。那么熱層大氣密度模式誤差的統(tǒng)計(jì)特性上是否為零均值，非零均值如何影響航天器軌道預(yù)報(bào)，標(biāo)準(zhǔn)差又是如何影響軌道預(yù)報(bào)的呢。典型熱層模式誤差及對軌道預(yù)報(bào)影響的相關(guān)文獻(xiàn)多見報(bào)道，但未見文獻(xiàn)給出上述問題答案。本文通過熱層模式誤差特性的統(tǒng)計(jì)和軌道預(yù)報(bào)數(shù)值試驗(yàn)，將給出上述問題的答案。期待為熱層模式的軌道預(yù)報(bào)應(yīng)用改進(jìn)、建模指引方向。

2 模式誤差特性

從熱層大氣建模的角度來看，熱層模式不應(yīng)存在系統(tǒng)性偏差，也就是相對誤差均值不為零的情況。由于空間環(huán)境變化極其復(fù)雜，累積的熱層密度探測資料在時(shí)間空間上存在不連續(xù)、分辨率低等情況，目前熱層大氣的部分變化機(jī)制仍不確定。目前的熱層大氣模式在不同空間環(huán)境下表現(xiàn)出差異化的系統(tǒng)偏差。

2.1 模式長期特性

熱層模式的長期特性是指大時(shí)間尺度 (太陽活動周)上表現(xiàn)出來的偏差、標(biāo)準(zhǔn)差等性質(zhì)。

為定量評估動態(tài)改進(jìn)前后模式性能，我們定義模式偏離實(shí)測的百分比，即相對誤差ε:

式中，ρO，ρM分別是密度實(shí)測值和模式值。

PARDINI等人通過分析四顆衛(wèi)星半長軸衰減，評估了JR-71、 MSIS-90、 MSISE00、 GOST-2004、JB2006和JB2008六種大氣模式，指出由于建模時(shí)采用高度無關(guān)的固定值阻力系數(shù)，此六種大氣模式均高估500km以下大氣密度約7%～20%[14]。平均偏差最小的模式依次是 JB08、MSISE00和GOST-2004。同時(shí)給出了第23太陽活動周峰年各模式在不同高度的平均偏差，見表1。

表1 各模式不同高度的平均偏差Tab.1 The mean bias of different models

劉衛(wèi)等人基于CHAMP衛(wèi)星2001年5月15日至2008年12月31日期間加速度計(jì)反演熱層大氣密度數(shù)據(jù)為基準(zhǔn)，對JB2008和MSISE00兩模式誤差進(jìn)行統(tǒng)計(jì)分析[15]，給出了 JB2008和MSISE00兩種模式在不同類型的空間環(huán)境下，相對誤差均值、標(biāo)準(zhǔn)差的特性，如表2所示。

表2 兩大氣模式相對誤差均值與標(biāo)準(zhǔn)差@400kmTab.2 The mean relative error and standard deviation at 400km altitude

由上述統(tǒng)計(jì)分析表明，長期而言典型模式誤差均值仍不為零，并且在不同空間環(huán)境、高度表現(xiàn)也不同。

2.2 事件期間的模式特性

空間環(huán)境事件期間的模式特性是指小時(shí)間尺度 (幾天)上，發(fā)生空間環(huán)境擾動事件時(shí)，熱層大氣模式表現(xiàn)出的誤差特性，包括誤差均值和標(biāo)準(zhǔn)差。圖1是基于CHAMP加速度計(jì)反演熱層密度，以日為區(qū)間單位統(tǒng)計(jì)給出的MSISE00模式的相對誤差均值和標(biāo)準(zhǔn)差隨太陽活動指數(shù)F10.7和地磁活動指數(shù)Ap的變化情況。圖中可見MSISE00模式相對誤差均值μ和標(biāo)準(zhǔn)差σ，隨空間環(huán)境變化，幾次磁暴均存在μ和σ陡增，μ達(dá)到100%，σ甚至可以達(dá)到300%。

圖2給出2003年10月31日磁暴事件期間JB2008和MSISE00模式相對誤差分布。圖2可見，JB2008相對MSISE00模式誤差分布更集中，數(shù)據(jù)統(tǒng)計(jì)兩模式在該磁暴事件期間的相對誤差均值分別為48%和75%。

圖1 MSISE00模式相對誤差均值和標(biāo)準(zhǔn)差隨空間環(huán)境的變化Fig.1 Variation of mean and standard deviation of MSISE00 model error with space environment

圖2 2003年10月31日磁暴事件期間JB2008和MSISE00模式相對誤差分布圖Fig.2 Relativeerror distribution of JB2008 and MSISE00 models during the 2003/10/31 event

由上述統(tǒng)計(jì)可見，無論是熱層模式的長期特性還是空間環(huán)境事件期間的模式特性，熱層大氣模式誤差均值和標(biāo)準(zhǔn)差都不為零，空間環(huán)境事件期間表現(xiàn)更加明顯。

3 軌道預(yù)報(bào)影響統(tǒng)計(jì)

3.1 軌道預(yù)報(bào)

CHAMP衛(wèi)星的快速科學(xué)軌道 (Rapid Science Orbit，RSO)數(shù)據(jù)具有較好的連續(xù)性，且軌道位置精度優(yōu)于10cm[12]。我們以CHAMP衛(wèi)星采樣間隔30s的RSO星歷數(shù)據(jù)為參考軌道，對軌道預(yù)報(bào)結(jié)果進(jìn)行比較分析，軌道預(yù)報(bào)中的主要模型和參量:

(1)99×99階地球引力場模型GGM05C；

(2)30×30階TOPEX 4.0海潮模型；

(3)CHAMP加速度計(jì)反演大氣密度；

(4)坐標(biāo)、時(shí)間及天文常數(shù)系統(tǒng)采用IERS(2010)標(biāo)準(zhǔn)；

(5)日、月引力攝動；

(6)太陽光壓，固體潮攝動。

3.2 高斯分布誤差

符合高斯分布誤差產(chǎn)生同樣是模式誤差如何影響軌道分析的關(guān)鍵。我們的高斯分布誤差產(chǎn)生是基于兩個(gè)均勻分布隨機(jī)變量的組合獲取的，具體理論如下。如果x1，x2… 是具有聯(lián)合概率分布p(x1，x2…)dx1dx2…的隨機(jī)偏離，并且如果y1，y2…每個(gè)都是所有x的函數(shù)，則y的聯(lián)合概率分布為:

考慮 (0，1)內(nèi)兩個(gè)一致偏離x1，x2與兩個(gè)變量y1，y2之間的變換:

等價(jià)地可寫成:

雅克比行列式可以如下計(jì)算[16]:

上式是y1函數(shù)與y2函數(shù)的乘積，都獨(dú)立服從正態(tài)分布。這是符合 (0，1)高斯分布的隨機(jī)變量y，那么任意(μ，σ)高斯分布的隨機(jī)變量y′可如下獲取:

圖3是實(shí)際仿真計(jì)算中產(chǎn)生的高斯分布隨機(jī)誤差，圖3中 (a)、(b)、(c)設(shè)定的均值和標(biāo)準(zhǔn)差分別是:30%、100%；30%、50%；10%、80%；實(shí)際產(chǎn)生的誤差均值和標(biāo)準(zhǔn)差分別是:29.08%、98.68%；29.54%、49.34%；10.74%、78.94%，符合要求。

式中，ρC，ρO分別是計(jì)算中使用的密度值和密度實(shí)測值。高斯分布隨機(jī)誤差由上述公式產(chǎn)生后，按照式 (8)計(jì)算熱層密度，進(jìn)行軌道預(yù)報(bào)并統(tǒng)計(jì)預(yù)報(bào)誤差情況。

4 結(jié)果分析

前節(jié)模式特性分析可見，空間環(huán)境事件期間熱層大氣模式存在較大的誤差均值和標(biāo)準(zhǔn)差。因此空間環(huán)境擾動期的航天器軌道預(yù)報(bào)是航天工程應(yīng)用的痛點(diǎn)問題。我們以3h值ap≥236為條件篩選出2001～2010年期間，具有CHAMP軌道數(shù)據(jù)的16個(gè)磁暴事件，見表3。施加符合不同均值和標(biāo)準(zhǔn)差的高斯分布熱層密度誤差，分析軌道預(yù)報(bào)誤差。發(fā)現(xiàn)16個(gè)磁暴事件中，在相同誤差參數(shù)下，軌道預(yù)報(bào)特性相似，下面選取幾個(gè)典型事例進(jìn)行展示分析。

4.1 相同均值的比較

以2003年10月29日磁暴事件為例，在軌道預(yù)報(bào)中，分別施加誤差均值μ=0%，-10%，-30%，10%，30%，每種誤差均值下又結(jié)合不同的誤差標(biāo)準(zhǔn)差σ=0%，10%，30%，50%，80%，100%進(jìn)行軌道預(yù)報(bào)誤差統(tǒng)計(jì)。圖4展示的是相對誤差均值相同，不同標(biāo)準(zhǔn)差下軌道預(yù)報(bào)誤差。圖中相同均值下，軌道預(yù)報(bào)誤差對標(biāo)準(zhǔn)差的變化不敏感，隨標(biāo)準(zhǔn)差的增大有所增加。整體上軌道預(yù)報(bào)誤差隨誤差均值絕對值增加而迅速增大，均值的正負(fù)決定預(yù)報(bào)偏差的方向。

圖3 實(shí)際產(chǎn)生的高斯隨機(jī)誤差分布圖Fig.3 Actual Gauss random error distribution

表3 篩選的空間環(huán)境事件信息Tab.3 The selected space environmental event information

圖4 不同相對誤差標(biāo)準(zhǔn)差下軌道預(yù)報(bào)誤差 (均值0%，-10%，-30%)Fig.4 Orbit prediction errors(mean values 0%，-10%，-30%，10% and 30%)under different relative error standard deviations

4.2 相同標(biāo)準(zhǔn)差比較

以2003年10月30日磁暴事件為例，在軌道預(yù)報(bào)中，分別施加誤差標(biāo)準(zhǔn)差σ=10%，30%，50%，80%，100%，每種相對誤差標(biāo)準(zhǔn)差又結(jié)合不同的均值μ=0%，-10%，-30%，10%，30%進(jìn)行軌道預(yù)報(bào)誤差統(tǒng)計(jì)。圖5、圖6展示的是相對誤差標(biāo)準(zhǔn)差相同，不同均值下軌道預(yù)報(bào)誤差。與上節(jié)結(jié)果一致，軌道預(yù)報(bào)誤差主要取決于誤差均值絕對值，對標(biāo)準(zhǔn)差的變化不敏感。誤差均值的正負(fù)決定預(yù)報(bào)偏差的方向，同時(shí)證明軌道預(yù)報(bào)結(jié)果對熱層小尺度的變化不敏感，熱層密度誤差可以通過平均進(jìn)行補(bǔ)償。

圖5 不同相對誤差均值下軌道預(yù)報(bào)誤差 (標(biāo)準(zhǔn)差10%和30%)Fig.5 Orbit prediction errors(standard deviations 10% and 30%)under different relative error mean values

圖6 不同相對誤差均值下軌道預(yù)報(bào)誤差 (標(biāo)準(zhǔn)差50%和80%)Fig.6 Orbit prediction errors(standard deviations 50% and 80%)under different relative error mean values

5 結(jié)論

本文介紹基于高斯分布的密度誤差對航天器軌道的影響分析的必要性，統(tǒng)計(jì)典型熱層模式誤差的長期特性和空間環(huán)境事件期間模式誤差特性，給出航天器軌道預(yù)報(bào)中引入高斯分布熱層密度相對誤差的方法，并進(jìn)行驗(yàn)證。篩選出2001～2010年期間，具有CHAMP軌道數(shù)據(jù)的16個(gè)磁暴事件，施加不同均值和標(biāo)準(zhǔn)差的高斯分布熱層密度誤差，分析軌道預(yù)報(bào)誤差。得到的主要結(jié)論如下:

(1)熱層大氣模式存在非零均值誤差，且該誤差隨空間環(huán)境而變化；

(2)典型熱層模式誤差均值和標(biāo)準(zhǔn)差，在空間環(huán)境擾動期最大可超過100%；

(3)軌道預(yù)報(bào)誤差主要取決于熱層模式誤差均值的絕對值，對誤差標(biāo)準(zhǔn)差不敏感，誤差均值的正負(fù)決定預(yù)報(bào)偏差的方向；

(4)軌道預(yù)報(bào)結(jié)果對熱層小尺度的變化不敏感，熱層密度誤差可以通過平均進(jìn)行補(bǔ)償。

期望本文對熱層大氣模式改進(jìn)、建模具有方向指引意義。