李文文, 李敏,2*, 施闖,2, 趙齊樂,2

1 武漢大學(xué)衛(wèi)星導(dǎo)航定位技術(shù)研究中心， 武漢 430079 2 地球空間信息技術(shù)協(xié)同創(chuàng)新中心， 武漢 430079

?

基于GRACE星載加速度計數(shù)據(jù)的熱層密度反演

李文文1, 李敏1,2*, 施闖1,2, 趙齊樂1,2

1 武漢大學(xué)衛(wèi)星導(dǎo)航定位技術(shù)研究中心， 武漢 430079 2 地球空間信息技術(shù)協(xié)同創(chuàng)新中心， 武漢 430079

本文主要研究了利用GRACE星載加速度計數(shù)據(jù)反演熱層密度.首先聯(lián)合采用GRACE衛(wèi)星2007—2009年星載加速度計數(shù)據(jù)和星載GPS數(shù)據(jù)進(jìn)行動力學(xué)定軌并同時估計加速度計校正參數(shù)，依此對加速度計數(shù)據(jù)進(jìn)行了校正.根據(jù)Sentman稀薄空氣動力學(xué)方程計算衛(wèi)星空氣動力系數(shù)，對校正后的加速度計數(shù)據(jù)進(jìn)行處理，反演得到了該時期沿軌熱層大氣密度.為分析反演密度的精度，將本文反演得到的GRACE-A衛(wèi)星沿軌密度與Doornbos的解算結(jié)果，以及經(jīng)驗密度模型NRLMSISE00，HASDM模型進(jìn)行比較分析.統(tǒng)計結(jié)果表明，本文反演結(jié)果比Doornbos系統(tǒng)性偏大約5%～8%，二者間的標(biāo)準(zhǔn)差(STD, standard deviation)在10%以內(nèi)，具有較好的符合性.其差異主要是由于采用了不同的加速度標(biāo)校以及空氣動力系數(shù)計算策略.本文反演得到的熱層密度較HASDM模型表現(xiàn)為正的系統(tǒng)性偏差且幅度在4%以內(nèi)，而Doornbos的結(jié)果較HASDM模型約系統(tǒng)性偏小4%～7%，二者與HASDM模型的標(biāo)準(zhǔn)差均為30%左右.另外，本文反演密度與NRLMSISE00模型之間存在約30%～40%的系統(tǒng)性bias，其STD也在30%左右.

熱層密度； GRACE； 加速度計； HASDM； NRLMSISE00； 空氣動力

1 引言

星載加速度計數(shù)據(jù)包含了衛(wèi)星在軌所受非保守力信息，通過分離太陽輻射光壓(Solar Radiation Pressure, SRP)、地球反照輻射等其他非保守力后可以獲得衛(wèi)星所受的空氣動力，進(jìn)而可以用于熱層密度反演.隨著CHAMP、GRACE、GOCE等重力衛(wèi)星計劃的順利實施，星載加速度計資料跨度可達(dá)一個太陽活動周期，這為研究熱層密度提供了大量高精度、高時空分辨率的數(shù)據(jù).Bruinsma和Biancale(2003)，Bruinsma等(2004), Sutton等(2007), Doornbos等(2010), Doornbos(2012), 汪宏波和趙長印(2009)研究了采用CHAMP、GRACE星載加速度計實測數(shù)據(jù)計算熱層大氣密度.這些成果極大地促進(jìn)了熱層密度對太陽耀斑、太陽風(fēng)暴和地磁暴的響應(yīng)(Bruinsma et al., 2006, 2010; Lei et al., 2011; Liu et al., 2011; Sutton et al., 2005, 2006)、熱層大尺度行擾(Bruinsma et al., 2007a, Bruinsma and Forbes, 2008, 2009; Sutton et al., 2006)、太陽輻射變化等對熱層密度的影響(Guo et al., 2007)、熱層密度年內(nèi)周期性變化(Guo et al., 2008)等方面的研究.

利用CHAMP、GRACE星載加速度計數(shù)據(jù)反演熱層密度時首先要對加速度計數(shù)據(jù)進(jìn)行校正.通常認(rèn)為GFZ發(fā)布的CHAMP、GRACE衛(wèi)星的加速度計數(shù)據(jù)在每個方向均存在尺度誤差和偏差誤差(Case et al., 2010).Gerlach等(2003)，徐天河和楊元喜(2004) 研究了采用能量守恒法進(jìn)行加速度計校準(zhǔn)，該方法通過能量守恒方程并采用現(xiàn)有重力場模型估計標(biāo)校參數(shù)，但是嚴(yán)重依賴于先驗重力場模型精度.Helleputte等(2008, 2009)在精密定軌(Precise Orbit Determination, POD)過程中以加速度計數(shù)據(jù)代替非保守力，通過估計加速度計標(biāo)校參數(shù)實現(xiàn)對CHAMP和GRACE衛(wèi)星加速度計數(shù)據(jù)的校正.其認(rèn)為如果加速度計校正參數(shù)精確確定，則聯(lián)合高精度的保守力模型以及校正后的加速度計數(shù)據(jù)進(jìn)行軌道積分后得到的軌道應(yīng)接近于“真實軌道”，因此加速度計校正的好壞可以通過軌道精度來體現(xiàn).由于采用該校正方法時尺度參數(shù)與偏差參數(shù)強(qiáng)反相關(guān)，van Helleputte等(2009)對尺度和偏差因子引入了強(qiáng)約束，另外為了補(bǔ)償力學(xué)模型誤差還估計了經(jīng)驗力參數(shù).結(jié)果表明該方法定軌結(jié)果與JPL精密科學(xué)軌道(Precision Scientific Orbit, PSO)相比三維RMS(Root Mean Square)為3.5 cm，SLR檢核RMS為2.3 cm，且將尺度因子固定為常數(shù)時可估計獲得更為平滑的偏差序列.熊永清等(2011a)提出了利用CHAMP衛(wèi)星精密軌道數(shù)據(jù)作為偽觀測數(shù)據(jù)，通過軌道擬合的方式估計衛(wèi)星初始狀態(tài)向量和加速度計參數(shù).該方法能夠避免處理GPS數(shù)據(jù)的復(fù)雜性.為消弱尺度因子和線性偏差的相關(guān)性，熊永清等(2011b)在以上方法基礎(chǔ)上聯(lián)合多個單天法方程進(jìn)行二次校標(biāo)，其多天只解算一組共同的尺度因子，而線性偏差每個單天解算一組.

衛(wèi)星空氣動力系數(shù)CA也是限制反演精度的重要因素之一.CA通常與衛(wèi)星的幾何形狀、表面材質(zhì)物理性質(zhì)(對大氣粒子的吸附、反射或其他交互作用性質(zhì))以及大氣狀態(tài)，如：大氣溫度，大氣中各種成分的濃度等因素有關(guān)，其可通過Cook模型(Cook 1965)、Sentman稀薄空氣動力學(xué)方程(Sentman, 1961)等模型進(jìn)行計算.Cook模型中忽略了粒子隨機(jī)熱運(yùn)動對衛(wèi)星的作用.Koppenwallner(2008)，Sutton (2009)的分析表明，對于CHAMP、GRACE這類衛(wèi)星，其與氣流方向平行的側(cè)面面板較為細(xì)長，大氣粒子的熱運(yùn)動會導(dǎo)致CA顯著增大.因此這種情況下采用Cook模型計算的空氣動力系數(shù)偏低，會導(dǎo)致反演得到的熱層密度系統(tǒng)性偏大.Sentman方程(Sentman, 1961)考慮了入射粒子的不規(guī)則熱運(yùn)動，并假設(shè)反射粒子全漫反射分布，被認(rèn)為是目前適合用于計算CHAMP、GRACE衛(wèi)星CA系數(shù)的模型(Sutton 2009).由于低軌衛(wèi)星在軌飛行中，衛(wèi)星表面會吸附一層氧原子，從而改變了衛(wèi)星表面材料對大氣粒子吸附、反射等性質(zhì)，Moe等(2004, 2005)在Sentman方程基礎(chǔ)上引入能量調(diào)節(jié)系數(shù)α.能量調(diào)節(jié)系數(shù)決定了粒子是否保持其平均熱能(α=0)或者獲得衛(wèi)星的表面溫度(α=1)，對CA有顯著影響.考慮到大氣中各種成分對CA的貢獻(xiàn)是非線性的，Doornbos (2012)在其研究中分別計算各種大氣成分對CA的影響并最終取其濃度加權(quán)平均值.

Bruinsma等(2003,2004)最早采用CHAMP衛(wèi)星加速度計數(shù)據(jù)研究熱層密度反演，其將校正后的加速度計數(shù)據(jù)扣除SRP等其他非保守力后得到“大氣阻力觀測值”，以該阻力觀測值和由經(jīng)驗?zāi)Ｐ兔芏扔嬎愕淖枇Φ谋戎底鳛槌叨纫蜃樱瑢δＰ兔芏冗M(jìn)行改正從而得到了“觀測密度”.在此基礎(chǔ)上，Bruinsma等(2004)分析了由太陽、地磁活動以及潮汐引起的熱層密度變化.Sutton等(2007)則系統(tǒng)研究了利用CHAMP、GRACE加速度數(shù)據(jù)反演熱層密度的方法，利用校正后的星固系X軸和Y軸加速度計數(shù)據(jù)解算熱層密度和大氣風(fēng)速，并分析了加速度計校準(zhǔn)誤差、大氣阻力/升力系數(shù)誤差、先驗風(fēng)速模型誤差、SRP誤差等對反演結(jié)果精度的影響，其結(jié)果表明反演密度誤差一般小于15%，但是在地磁活動劇烈期間會由于大氣風(fēng)速誤差增加而增至30%.Doornbos(2012)在以上研究基礎(chǔ)上，利用Angara軟件對衛(wèi)星表面力進(jìn)行精化建模，并分析了多種密度反演方式.國內(nèi)汪宏波和趙長印(2008)處理了CHAMP衛(wèi)星2001年至2002年星載加速度計數(shù)據(jù)，并反演得到了沿軌熱層密度.其利用反演密度對經(jīng)驗大氣模型DTM94和MSIS90的精度進(jìn)行了檢核，結(jié)果表明在太陽活動峰年DTM94模型誤差為30%～35%，MSIS90模型為25%～30%(汪宏波和趙長印, 2008).

上述研究中，Bruinsma等(2003,2004)，Sutton等(2007), 汪宏波和趙長印(2008)采用單天估計一組偏差和尺度因子的策略進(jìn)行加速度計校正.然而，根據(jù)van Helleputte等(2008,2009)的結(jié)果，利用該策略校正的加速度計數(shù)據(jù)進(jìn)行軌道積分的精度僅在10 cm水平，表明其標(biāo)校精度相對較低.Doornbos(2012)根據(jù)van Helleputte等(2009)的研究結(jié)果進(jìn)行加速度計校正，其校正時為了補(bǔ)償保守力模型誤差每600 s估計一組切向、法向和徑向經(jīng)驗力參數(shù).然而，頻繁調(diào)節(jié)的經(jīng)驗力參數(shù)同時還會吸收部分非保守力，從而導(dǎo)致校正參數(shù)可能出現(xiàn)偏差.Bruinsma等(2003,2004)，汪宏波和趙長印(2008)在其研究中均采用Cook模型計算CHAMP衛(wèi)星空氣動力系數(shù)，這可能導(dǎo)致其反演的密度偏大.

本文在以上基礎(chǔ)上研究了利用GRACE星載加速度計數(shù)據(jù)進(jìn)行熱層密度反演.在反演中，采用了與Van Helleputte等 (2009)相同的方式對GRACE衛(wèi)星加速度計進(jìn)行標(biāo)校，但是本文在不引入經(jīng)驗力參數(shù)的前提下，通過調(diào)整定軌中的動力學(xué)策略，達(dá)到了與前者一致的軌道精度水平，同時也保證了標(biāo)校參數(shù)的完整性.通過與已有結(jié)果和模型的比較，對本文反演密度的精度進(jìn)行了評價.

2 數(shù)據(jù)源與模型

本文選擇GRACE衛(wèi)星2007—2009年L1b級的數(shù)據(jù)進(jìn)行處理，其包括星載GPS數(shù)據(jù)(GPS)，加速度計數(shù)據(jù)(ACC)，姿態(tài)四元數(shù)數(shù)據(jù)(SCA)，衛(wèi)星質(zhì)量數(shù)據(jù)(MAS)，衛(wèi)星機(jī)動數(shù)據(jù)(THR)等.數(shù)據(jù)的具體說明可參見文獻(xiàn)(Case et al.,2010).在進(jìn)行加速度標(biāo)校時，本文采用IGS提供的最終GPS精密星歷和30 s采樣率鐘差產(chǎn)品作為數(shù)據(jù)處理的時空基準(zhǔn)，其軌道、鐘差精度分別在2.5 cm(1D RMS)和75 ps(RMS)量級(http:∥igs.org/products).反演過程中，衛(wèi)星姿態(tài)通過恒星照相儀SCA數(shù)據(jù)計算，衛(wèi)星在軌質(zhì)量通過MAS數(shù)據(jù)獲取.為維持KBR儀器相互對準(zhǔn)，GRACE衛(wèi)星在軌運(yùn)行時經(jīng)常需要進(jìn)行機(jī)動調(diào)整，由此引入的機(jī)動力需在反演時扣除，否則會引起密度反演誤差.本文在處理時則依據(jù)THR數(shù)據(jù)剔除了機(jī)動時期的數(shù)據(jù).

為了對反演密度的精度進(jìn)行評價，本文將處理結(jié)果與Doornbos、Sutton博士反演的結(jié)果進(jìn)行了比較.Doornbos、Sutton博士的數(shù)據(jù)成果均可以在線公開獲取，其具體的處理方法與策略可以參考Sutton等(2007)與Doornbos(2012).同時還對比了經(jīng)驗密度模型，包括：NRLMSISE00模型(Picone et al., 2002)和HASDM模型(Bowman and Storz, 2003).其中，NRLMSISE00模型由美國海軍實驗室在MSISE系列經(jīng)驗密度模型的基礎(chǔ)上進(jìn)行改進(jìn)得到，是目前低軌空間物體軌道確定和預(yù)報，以及熱層-電離層耦合等空間天氣研究中常用的經(jīng)驗?zāi)Ｐ?其誤差一般在10%～30%且隨高度增加而增大，在400～500 km高度處約在25%水平(Macros et al., 2006).HASDM模型則由美國空軍研究中心開發(fā)，其使用75～80個低軌空間物體的雷達(dá)跟蹤數(shù)據(jù)在確定軌道的同時估計大氣密度的改正值，并以球諧函數(shù)的形式每3 h估計一次Jacchia71模型中的拐點溫度和夜間散逸層最低溫度(Bowman and Storz, 2003).HASDM模型在200～800 km高度范圍的精度可達(dá)6%～8%(Storz et al., 2005).

3 加速度計校正

3.1 聯(lián)合GPS與加速度計數(shù)據(jù)的校正方法

GRACE衛(wèi)星上搭載的SuperSTAR星載加速度計能以極高的精度測量衛(wèi)星所受到的非保守力.SuperSTAR的X、Y軸(指向星固系下的Z、X軸)是高敏感軸，而Z軸(指向星固系Y軸)為低敏感軸，三軸的測量精度都認(rèn)為可以達(dá)到10-9m·s-2(Frommknecht et al., 2003).然而由于星上溫度變化，儀器老化以及儀器開關(guān)切換等原因，會導(dǎo)致加速度計測量的加速度數(shù)據(jù)包含了偏差和尺度變化.因此在使用加速度計數(shù)據(jù)時，首先要對其進(jìn)行校準(zhǔn).由于大氣阻力的作用，切向的加速度計校正參數(shù)比較容易準(zhǔn)確確定.法向與徑向的尺度與偏差參數(shù)與衛(wèi)星的初始狀態(tài)向量強(qiáng)相關(guān)，且法向是約束最弱的方向(van Helleputte et al., 2008).

GRACE衛(wèi)星搭載的GPS接收機(jī)能夠提供大量連續(xù)的高精度觀測數(shù)據(jù)，本文聯(lián)合星載GPS和加速度計數(shù)據(jù)進(jìn)行動力學(xué)定軌，并同時估計加速度計校正參數(shù).定軌中動力學(xué)模型策略主要參照GFZ最新發(fā)布的Level 2數(shù)據(jù)處理標(biāo)準(zhǔn)0005版本(Dahle et al., 2012)并進(jìn)行了一定調(diào)整.保守力方面選擇精度盡量高的模型，其中地球非球形引力攝動采用EIGEN-6C模型計算，且靜態(tài)重力場截至130階，時變部分則截至50階，海潮模型采用FES2004模型并截至40階，固體潮、極潮、相對論等其他保守力則均采用IERS Conventions(2003)的推薦模型(McCarthy and Petit, 2004).非保守力模型則以加速度計數(shù)據(jù)替代.在使用加速度計數(shù)據(jù)時，先由GFZ發(fā)布的初始校正參數(shù)(Bettadpur, 2004)進(jìn)行預(yù)校正，然后在定軌過程中同時估計X、Y、Z三方向的偏差參數(shù)，尺度參數(shù)則不予估計.這樣處理實際上是將尺度因子固定為常數(shù)，而僅估計偏差因子.其原因是一方面GRACE衛(wèi)星星上溫度控制較為穩(wěn)定，其尺度參數(shù)可以認(rèn)為是常數(shù)而不隨時間變化(Bruinsma et al., 2007)；另一方面，由于偏差參數(shù)和尺度參數(shù)強(qiáng)相關(guān)，估計的偏差參數(shù)會補(bǔ)償尺度參數(shù)先驗值的誤差.Bezděk(2010)，熊永清等(2011b)的研究中也表明了這一現(xiàn)象.偏差參數(shù)一般采用分段常數(shù)模型估計.本文測試了多組不同分段間隔情況下定軌的精度，發(fā)現(xiàn)估計間隔為60 min時軌道精度最佳，因此采用每60 min估計一組偏差參數(shù)的策略.這也與GFZ的數(shù)據(jù)處理標(biāo)準(zhǔn)(Dahle et al., 2012)中的策略一致.對于GPS觀測值，采用無電離層組合偽距(PC)和無電離層組合相位(LC)作為基本觀測模型，觀測誤差主要考慮GRACE和GPS衛(wèi)星天線相位偏心(phase center offset, PCO)以及GPS衛(wèi)星天線相位中心變化(phase center variation, PCV)、相對論效應(yīng)、相位纏繞誤差等.模糊度參數(shù)未進(jìn)行固定，而是直接采用浮點解.具體的數(shù)據(jù)處理策略見表1.3.2 校正結(jié)果評價

3.2.1 GRACE衛(wèi)星軌道精度

如前所述，利用校正后的加速度計數(shù)據(jù)結(jié)合保守力模型進(jìn)行軌道積分，所得軌道的精度可以反映加速度計校正的精度.依據(jù)表1中的數(shù)據(jù)處理策略，對GRACE衛(wèi)星2007—2009年數(shù)據(jù)進(jìn)行精密定軌處理，并將解算軌道與JPL PSO比較，統(tǒng)計每個定軌弧段的軌道差值RMS，其時序圖見圖1.由圖中可見，以JPL PSO為參考軌道，本文解算的GRACE-A/B衛(wèi)星軌道在切向(Along-track) 與法向(Cross-track) 的RMS均主要在1.0～2.0 cm范圍變化，而徑向(Radial) RMS序列變化相對較為穩(wěn)定，對于GRACE-A、B星均在1.5 cm左右變化.表2為三方向軌道差值的平均RMS，其中徑向平均RMS較法向與切向均略差，這可能是由于JPL PSO對GRACE衛(wèi)星的天線PCV進(jìn)行了改正，而本文并未考慮該影響.GRACE-A，B衛(wèi)星三年軌道差值的平均3D RMS分別為2.6 cm和2.3 cm.

SLR技術(shù)常作為外部檢核的方式用于軌道精度評價.采用NPT格式SLR數(shù)據(jù)對本文解算的GRACE衛(wèi)星軌道進(jìn)行檢核，統(tǒng)計每個弧段SLR觀測值殘差(Observed Minus Computed, OMC)的平均值與標(biāo)準(zhǔn)差.統(tǒng)計時剔除高度角小于10°以及OMC絕對值大于50 cm的觀測值，最終所有剔除的觀測值約占總觀測數(shù)量的2.5%，所有弧段的SLR OMC平均值與標(biāo)準(zhǔn)差時間序列見圖2.可見OMC平均值主要在±3 cm以內(nèi)變化，未見明顯系統(tǒng)偏差.標(biāo)準(zhǔn)差(Standard Deviation, STD)則均主要在5 cm以內(nèi).表3為所有弧段SLR OMC平均值，STD和RMS的均值.同時對這三年JPL PSO軌道進(jìn)行了SLR檢核和精度統(tǒng)計.對比可見，SLR檢核本文軌道的殘差平均值均小于2 mm，JPL PSO則存在約1 cm的系統(tǒng)偏差.二者SLR檢核RMS均在2.5 cm水平.

表1 利用GRACE星載GPS和加速度計數(shù)據(jù)進(jìn)行精密定軌 以及加速度計校正的數(shù)據(jù)處理策略Table 1 Data processing strategy for POD and accelerometer data calibration using GPS observations and accelerometer data

表2 GRACE-A/B衛(wèi)星與JPL PSO軌道 比較的平均RMS統(tǒng)計Table 2 RMS statistics of orbit differences of the GRACE-A/B satellite orbits w.r.t. JPL PSO orbits

表3 本文解算的GRACE-A/B衛(wèi)星軌道 以及JPL PSO軌道的SLR檢核精度統(tǒng)計Table 3 Mean, STD, RMS statistics of the calculated GRACE-A/B orbits from this paper and the JPL PSO

3.2.2 加速度計校正參數(shù)

圖3為2007—2009年間GRACE-A和GRACE-B衛(wèi)星加速度計偏差參數(shù)在星固系X,Y,Z方向的偏差時間序列.圖中個別弧段估計得到的偏差出現(xiàn)明顯跳變，其與衛(wèi)星控制器切換以及衛(wèi)星溫度突變等原因相關(guān)(Bruinsma et al., 2007; Doornbos, 2012).對于GRACE雙星而言三個方向的偏差參數(shù)變化都非常平緩，且存在隨時間的緩慢飄移，其主要由加速度計的老化導(dǎo)致，這也與Bruinsma等(2007)的結(jié)果一致.根據(jù)估計的偏差時序，對于GRACE-A/B衛(wèi)星其在X方向偏差參數(shù)的年均變化分別僅為1×10-8m·s-2和5×10-9m·s-2.而在Y方向的偏差參數(shù)年均變化相比X方向明顯較大，分別為9×10-8m·s-2和1.84×10-7m·s-2.

圖1 2007—2009年GRACE-A/B衛(wèi)星與JPL PSO軌道比較的RMS序列Fig.1 RMS series of orbit differences with respect to the JPL PSO during 2007—2009

圖2 2007—2009年GRACE-A/B衛(wèi)星SLR殘差平均值(a)與標(biāo)準(zhǔn)差(b)序列Fig.2 The mean, STD series of the SLR OMC of GRACE-A/B satellites during 2007—2009

圖3 GRACE-A/B衛(wèi)星2007—2009年加速度計偏差序列(A) GRACE-A衛(wèi)星； (B) GRACE-B衛(wèi)星.Fig.3 Time series of estimated bias parameters of GRACE-A/B satellites(A) GRACE-A satellite; (B) GRACE-B satellite.

4 基于加速度計數(shù)據(jù)的熱層密度反演方法與策略

本節(jié)介紹了利用X軸加速度計數(shù)據(jù)提取熱層密度的基本方法，并分析了Cook模型和Sentman模型用于計算空氣動力系數(shù)的差異，以及其他非保守力對反演的影響.基于這些分析，本節(jié)確定了密度反演的具體數(shù)據(jù)處理策略.

4.1 利用X軸加速度計數(shù)據(jù)反演熱層密度

4.1.1 反演方法

(1)

式中,ρ為大氣密度，vr和vr分別為衛(wèi)星與大氣相對速度及其單位矢量，m為衛(wèi)星質(zhì)量，CA為整星總的空氣動力系數(shù)，其在vr方向的分量一般稱為阻力系數(shù)CD，垂直于vr方向的分量則稱為升力系數(shù)CL，Aref這里是參考面積可以為任一給定正值，uL方向為升力方向單位矢量，CD,i，CL,i系數(shù)分別為為面板i的大氣阻力和升力系數(shù).在計算時通常將每塊面板CD，CL系數(shù)除以Aref以使其無量綱.

GRACE衛(wèi)星姿態(tài)相對穩(wěn)定，其三軸旋轉(zhuǎn)歐拉角較小一般在1°以內(nèi).由于大氣阻力的作用，加速度計X軸分量信號強(qiáng)度大且校正精度較高(星固系X軸指向衛(wèi)星沿軌方向)，非常利于提取大氣密度.利用校正后的加速度計數(shù)據(jù)扣除SRP等非保守力并投影至X軸方向，得到大氣動力在X軸的分量aobs,X：

(2)

(3)

其中，acal為校正后的加速度計數(shù)據(jù)(慣性系下)，asrp為太陽輻射壓力，aalb為地球反照輻射壓力，aIR為地球紅外輻射壓力，eX為慣性系轉(zhuǎn)星固系旋轉(zhuǎn)矩陣中X方向的單位向量.根據(jù)大氣動力、熱層密度和CA之間的函數(shù)關(guān)系(式1)反解熱層密度：

(4)

由式(2)和(4)可知，利用校正后的加速度計數(shù)據(jù)反演密度時，需要對除空氣動力外的其他非保守力進(jìn)行準(zhǔn)確建模和扣除，并準(zhǔn)確計算衛(wèi)星的空氣動力系數(shù).

4.1.2 空氣動力系數(shù)模型

CD,CL系數(shù)反映了衛(wèi)星外表面與大氣粒子之間的交互作用的強(qiáng)弱，并與衛(wèi)星外表面材質(zhì)物理特性以及大氣粒子的成分等有關(guān).精確確定CD，CL系數(shù)是目前反演熱層密度的難點.由于CHAMP、GRACE這類衛(wèi)星側(cè)面面板較長，在Cook模型中被忽略的分子熱運(yùn)動會對其空氣動力系數(shù)有較大影響.Sutton(2009) 研究表明對于CHAMP衛(wèi)星，相比采用Cook模型，采用Sentman方程后反演的熱層密度與經(jīng)驗密度之間的標(biāo)準(zhǔn)差降低了約36%，偏差也由約56.4%～61.3%降低至21.2%～23.3%.

根據(jù)Bettadpur(2012)給出的GRACE衛(wèi)星星體幾何模型，圖4顯示了GRACE-A衛(wèi)星面板與星固系示意圖.其中1, 2號面板為衛(wèi)星前后面，3, 4, 5和6號面板是衛(wèi)星的側(cè)面，7, 8號分別為衛(wèi)星底面和頂面，9號面板為衛(wèi)星圓柱狀底部的平面投影.對于9號面板因其法向方向不唯一(可沿衛(wèi)星X-Y平面任意方向)且面積較小，故而未參與整星空氣動力系數(shù)的計算.圖5為分別利用Sentman方程和Cook模型計算的GRACE-A衛(wèi)星各個面板的CD系數(shù).由圖4中可見，由于1，8號面板外法向與相對速度方向夾角θ>90°，因此在Cook模型中其CD系數(shù)始終為0.2號面板外法向正對大氣氣流方向(θ接近0°)，其對衛(wèi)星的空氣動力系數(shù)貢獻(xiàn)最大，Cook模型計算的結(jié)果與Sentman模型相差約為0.05.相差最大的面板為7號面板，Sentman模型計算的CD系數(shù)約高出Cook模型0.4～0.6，這是由于7號面板是GRACE衛(wèi)星面積最大的面板，兩種模型之間的差異通過面板面積被放大.對所有面板的空氣動力系數(shù)求和，采用Cook模型和Sentman方程計算的整星空氣動力系數(shù)分別為2.37和3.50，其差異可達(dá)32%(以Sentman模型結(jié)果為參考).

圖4 GRACE-A衛(wèi)星面板與星固系示意圖(Doornbos, 2012)Fig.4 Spacecraft body-fixed directions and panels of GRACE-A satellite

4.1.3 GRACE衛(wèi)星非保守力分析

本文采用Box-wing模型計算SRP.地球反照輻射壓力采用Bhanderi and Bak(2005)提出的模型計算，所需的地表反照輻射系數(shù)采用CERES網(wǎng)站提供的1°×1°格網(wǎng)分辨率的月平均地表反射系數(shù)，而地球長波輻射則參考Rim(1992)給出的經(jīng)驗?zāi)Ｐ?由于熱層密度與太陽輻射強(qiáng)度相關(guān)，圖6分別給出了在太陽活動高年(2003年1月3日)和低年時期(2009年1月2日)GRACE-A衛(wèi)星所受的太陽輻射壓力(SRP)，地球輻射壓力(ALB+IR)以及空氣動力(AERO)的比較.在熱層反演中，最重要的是大氣阻力在X軸方向的信號強(qiáng)度.圖中可見，在太陽活動高年，由于熱層密度較大，其所受大氣阻力遠(yuǎn)大于太陽輻射壓力，在X軸方向可以達(dá)到(0.8～3.2)×10-7m·s-2，比太陽輻射壓力高一個量級.而在太陽活動低年，大氣阻力在X方向稍大于太陽輻射壓力，二者均在10-8m·s-2量級.并且圖中可見在大氣阻力X分量極小處，SRP的X分量與其大小幾乎相等.因此，衛(wèi)星所受SRP需要精確計算.而對于地球輻射壓力，其量級通常較小在(1-2)×10-8m·s-2水平且主要作用于衛(wèi)星徑向(星固系Z軸方向)，因此在實際反演計算中可以忽略.

4.2 策略分析

根據(jù)4.1節(jié)可知，加速度計校正，地球輻射壓力，空氣動力系數(shù)，衛(wèi)星質(zhì)量，大氣風(fēng)速等均能影響反演結(jié)果.大氣阻力主要作用于衛(wèi)星沿軌方向(即星固系X軸)，從而使得加速度計數(shù)據(jù)X軸分量信號強(qiáng)度大得多也更容易校正，這為反演熱層密度提供了有利條件.第3節(jié)中軌道精度的分析結(jié)果也驗證了本文校正的結(jié)果具有較高精度.而地球輻射壓力主要作用于衛(wèi)星徑向，在X軸的分量要比大氣阻力低2～3量級，因此在反演中可以忽略.GRACE衛(wèi)星質(zhì)量年均變化1 kg左右，以GFZ公布的480 kg為參考質(zhì)量，1 kg質(zhì)量誤差引起的密度偏差約為0.2%.在反演時也可采用GFZ發(fā)布的精確的衛(wèi)星質(zhì)量數(shù)據(jù)(MAS數(shù)據(jù)).大氣風(fēng)速誤差(相對于旋轉(zhuǎn)大氣)對密度反演的影響比較復(fù)雜，其導(dǎo)致衛(wèi)星與大氣的相對速度在量級與方向上均會產(chǎn)生誤差，從而引起衛(wèi)星空氣動力系數(shù)以及密度計算誤差.Bruinsma等(2004) 通過仿真分析表明大氣風(fēng)速與反演誤差幾乎呈線性關(guān)系，緯度分量在赤道地區(qū)每100 m·s-1誤差約引起4%的誤差，在兩極地區(qū)影響可以忽略，而經(jīng)線分量相反，在兩極地區(qū)每100 m·s-1誤差約引起5%的誤差，在赤道地區(qū)則較小為2%.大氣風(fēng)速一般在幾百米每秒量級，但是在空間環(huán)境擾動等時期，大氣風(fēng)速可以增大至千米每秒量級，因此在反演時應(yīng)盡量減小其影響.

圖5 GRACE-A衛(wèi)星各面板大氣阻力系數(shù)(α=0.93, 2003年1月3日)Fig.5 Drag coefficients of each panel of GRACE-A satellite based on Cook and Sentman models on 2003 Jan-3rd (fixed α=0.93)

圖6 不同太陽活動強(qiáng)度時期GRACE衛(wèi)星空氣動力與太陽輻射壓力(A)、地球輻射壓力(B)比較Fig.6 Comparison between aerodynamic, SRP and Earth radiation during different solar activity periods

根據(jù)以上分析，本文按照表4中的策略反演熱層密度.其中，加速度計數(shù)據(jù)采用3.1節(jié)估計得到的每60 min一組的標(biāo)校參數(shù)進(jìn)行校正，SRP采用Box-wing模型計算，Albedo與IR均不予考慮.空氣動力系數(shù)采用Sentman模型計算，因其對大氣溫度與衛(wèi)星面板溫度變化并不敏感，本文直接采用NRLMSISE00模型計算大氣溫度，衛(wèi)星表面溫度則采用Doornbos(2012)的推薦值273 K.調(diào)節(jié)系數(shù)是影響空氣動力系數(shù)的重要參量，本文研究中采用Sutton博士的推薦值0.93.大氣風(fēng)速采用HWM93模型計算.此外，為避免衛(wèi)星姿態(tài)調(diào)整機(jī)動時期引入的機(jī)動力的影響，本文根據(jù)THR產(chǎn)品(包含衛(wèi)星機(jī)動記錄數(shù)據(jù))剔除了衛(wèi)星點火機(jī)動時期的數(shù)據(jù).

表4 熱層密度反演數(shù)據(jù)處理策略Table 4 Data processing strategy of thermosphere density derivation

5 熱層密度反演結(jié)果與統(tǒng)計分析

采用第4節(jié)提出的策略處理得到了GRACE-A和GRACE-B衛(wèi)星2007—2009年間的沿軌熱層密度.為了對反演結(jié)果的精度進(jìn)行評價，本節(jié)將反演得到的GRACE-A衛(wèi)星2007—2009年沿軌密度與Doornbos(2012)等的反演結(jié)果，以及NRLMSISE00，HASDM等經(jīng)驗?zāi)Ｐ瓦M(jìn)行比較分析.

5.1 與Doornbos反演結(jié)果比較

圖7為Doornbos, Sutton博士反演的2007年1月28日GRACE-A衛(wèi)星沿軌密度與本文反演的結(jié)果進(jìn)行對比，其中ρobs為本文反演結(jié)果，Density_Direct為Doornbos采用單軸法反演結(jié)果，Density_Sutton為Sutton博士反演得到的結(jié)果.Density_Direct的反演策略與表4類似，其差異主要為：(1) 加速度計校正參數(shù)采用van Helleputte等(2009)中的結(jié)果；(2) 計算空氣動力學(xué)系數(shù)時分別考慮各種氣體粒子單獨(dú)的貢獻(xiàn).Sutton博士的反演策略也與本文類似，但是其將反演結(jié)果根據(jù)衛(wèi)星軌道緯度進(jìn)行了合并(Sutton et al., 2007)，因此其數(shù)據(jù)點時刻并不規(guī)則.由圖中可見三者反演密度的變化趨勢以及幅度均基本一致.

將本文反演結(jié)果ρobs與Density_Direct進(jìn)行比較，對二者密度比值按年份進(jìn)行對數(shù)正態(tài)分布擬合并統(tǒng)計bias和標(biāo)準(zhǔn)差；為了便于對照分析，將Density_Direct與HASDM模型也進(jìn)行了比較，結(jié)果分別見圖8和表5.其中年均比較數(shù)據(jù)點數(shù)約為2400000，統(tǒng)計密度比值的采樣間隔設(shè)置為0.002.由圖8可見兩組密度比值擬合的相關(guān)系數(shù)均達(dá)到了0.95以上，且Density_Direct與HASDM密度比值的擬合相關(guān)系數(shù)稍大一些.Density_Direct與ρobs，Density_Direct與HASDM的密度比值均主要分布在0.5～2.0區(qū)間內(nèi)，但是Density_Direct與ρobs的比值分布更為集中，表明二者結(jié)果更為接近.由表5可見，Density_Direct較HASDM模型、ρobs均偏小約4%～8%.Density_Direct與HASDM模型密度比值的STD在30%左右，與ρobs之間STD則均小于10%，表明本文反演結(jié)果與Doornbos的結(jié)果之間符合性較好.這主要是由于二者采用了相同的數(shù)據(jù)進(jìn)行反演，且反演中采用的模型類似.

ρobs與Density_Direct之間的差異主要由計算空氣動力系數(shù)與加速度計標(biāo)校的策略不一致所造成.Doornbos在計算CA時分別考慮了每種氣體粒子的貢獻(xiàn)，而本文在計算時則是直接采用大氣粒子的平均分子質(zhì)量.圖9為2007年1月1—15日半個月內(nèi)二者CA比值的頻數(shù)統(tǒng)計，其中統(tǒng)計的空氣動力系數(shù)數(shù)量為129600，采樣間隔為0.005.對該頻數(shù)分布進(jìn)行對數(shù)正態(tài)擬合的平均值與標(biāo)準(zhǔn)差分別為0.977和1.016，擬合的相關(guān)系數(shù)為0.915.值得注意的是計算空氣動力系數(shù)時所需要輸入的數(shù)據(jù)或模型參數(shù)均可以很容易保持一致，因此空氣動力系數(shù)差異的變化不大，其造成的反演密度差異也應(yīng)在此水平.Doornbos(2012)通過仿真分析了加速度計校正誤差對反演密度的影響，表明當(dāng)X軸加速度計數(shù)據(jù)存在10 nm·s-2的偏差時，其造成的密度反演偏差和STD分別為-2.7%和1.5%.本文研究中將尺度參數(shù)固定為常數(shù)并每60 min估計一組偏差參數(shù)，而van Helleputte等(2009)則單天估計一組尺度與偏差參數(shù)且在切向、徑向與法向均引入了經(jīng)驗力參數(shù)，二者標(biāo)校參數(shù)的差異對密度反演造成的影響值得進(jìn)一步研究.

圖7 Doornbos，Sutton等研究與本文反演結(jié)果對比Fig.7 Densities comparison of this paper with respect to results from Doornbos and Sutton

圖8 2007(A), 2008(B)和2009(C)年Density_Direct與HASDM模型以及ρ obs的密度比值的頻數(shù)統(tǒng)計以及對數(shù)正態(tài)分布擬合Fig.8 Frequency number count and log-normal fit of ratios of Density_Direct over ρobs and HASDM model densities during year 2007 (A), 2008 (B) and 2009 (C), respectively

圖9 Doornbos與本文計算的空氣動力系數(shù)差異的頻數(shù)統(tǒng)計和對數(shù)正態(tài)分布擬合Fig.9 Frequency count and log-normal fit of the differences of the aerodynamic coefficients between Doornbos and this paper

圖10 GRACE-A衛(wèi)星2007年1月28日沿軌反演密度，NRLMSISE00和HASDM模型密度以及密度比值(a) 沿軌ρobs，NRLMSISE00和HASDM密度；(b) ρobs和NRLMSISE00分別與HASDM密度比值；(c) 星下點緯度; (d) 軌道高度； (e) 星下點地方時.Fig.10 The accelerometer-derived, NRLMSISE00 and HASDM model densities of GRACE-A satellite on January 16, 2007, as well as their density ratios(a) Accelerometer-derived, NRLMSISE00 and HASDM model densities; (b) Ratio series of accelerometer-derived and NRLMSISE00 densities over HASDM model densities; (c) Orbit latitude; (d) Orbit altitude; (e) Local time of ground tracks.

圖11 2007(A), 2008(B)和2009(C)年ρ obs和 NRLMSISE00模型與HASDM模型密度比值的頻數(shù)統(tǒng)計以及對數(shù)正態(tài)分布擬合Fig.11 Frequency number count and log-normal fit of ratios of ρobs and NRLMSISE00 over HASDM model densities during year 2007 (A), 2008 (B) and 2009 (C) respectively

RatioofDensity_DirectoverHASDMRatioofDensity_Directoverρobs200720082009200720082009Bias-0.070-0.043-0.068-0.057-0.072-0.074STD1.2811.3281.3311.0741.0841.089

5.2 與NRLMSISE00和HASDM經(jīng)驗?zāi)Ｐ捅容^

圖10為2007年1月28日GRACE-A衛(wèi)星沿軌的反演密度ρobs，NRLMSISE00和HASDM模型密度以及密度比值.該天為地磁平靜日，其全天Ap指數(shù)均小于10，平均Ap指數(shù)為5，F(xiàn)10.7指數(shù)和81天平均F10.7指數(shù)分別為79.0和84.5 sfu，表明太陽活動強(qiáng)度較低.為了便于比較，圖中還給出了衛(wèi)星軌道緯度，高度和當(dāng)?shù)貢r.圖中可見軌道高度約在475 km左右.總體上ρobs與經(jīng)驗?zāi)Ｐ兔芏鹊淖兓厔莼疽恢?，但與HASDM更為接近，二者在密度峰值處大小非常接近.而NRLMSISE00模型預(yù)報密度較ρobs與HASDM均系統(tǒng)性偏大.ρobs與HASDM比值主要分布在0.5～1.5范圍內(nèi)，且在兩極地區(qū)時出現(xiàn)較大幅度的跳變，其原因可能是在極地地區(qū)大氣風(fēng)速很容易增大至千米每秒量級，采用HWM93模型計算的大氣風(fēng)速可能會存在較大的誤差，從而導(dǎo)致反演結(jié)果的誤差增大.圖 10中同時可見ρobs和NRLMSISE00模型密度在午夜時期均出現(xiàn)一個甚至多個峰值，這反映的是午夜密度極大(Midnight Density Maximum, MDM)現(xiàn)象.HASDM模型則未能反映出該變化.因此可以發(fā)現(xiàn)在午夜期間，ρobs、NRLMSISE00模型與HASDM模型的密度比值序列均存在峰值現(xiàn)象.

將本文反演的密度ρobs與NRLMSISE00模型預(yù)報密度以及HASDM模型比較并進(jìn)行統(tǒng)計.圖11為2007—2009年密度比值的頻數(shù)分布統(tǒng)計，以及對數(shù)正態(tài)分布擬合.圖中可見，ρobs和NRLMSISE00模型以及HASDM模型之間擬合的相關(guān)系數(shù)均達(dá)到了0.99以上.表6為密度比值的對數(shù)正態(tài)分布擬合的bias與STD.可見，ρobs與HASDM模型之間bias較小，均在4%以內(nèi)，而與NRLMSISE00模型之間存在約30%～40%的年平均bias，其原因目前認(rèn)為主要是由于2007—2009年太陽活動異常寧靜從而導(dǎo)致NRLMSISE00模型產(chǎn)生了尺度性bias(Emmert et al., 2010; Solomon et al., 2010; Liu et al.,2011).

表6 2007, 2008和2009年ρobs和NRLMSISE00與HASDM模型密度比值的對數(shù)正態(tài)分布擬合統(tǒng)計Table 6 Statistics of log-normal fit of ratios of ρobsand NRLMSISE00 over HASDM model densities during year 2007, 2008 and 2009 respectively

STD統(tǒng)計結(jié)果來看，ρobs與NRLMSISE00以及HASDM的STD均在30%左右，但是其與HASDM模型之間的STD均要稍小于NRLMSISE00模型.值得注意的是，ρobs較HASDM模型系統(tǒng)性偏大，而Doornbos反演結(jié)果系統(tǒng)性偏小約4%～7%，二者與HASDM模型之間的STD均在30%左右.

以上的統(tǒng)計分析表明本文采用加速度計反演的熱層密度與HASDM模型之間bias較小，而NRLMSISE00模型與本文反演密度以及HASDM均存在較大的bias.一般而言，不同密度模型或數(shù)據(jù)之間的bias主要來自于空氣動力系數(shù)誤差，以及對熱層密度變化物理機(jī)制的錯誤模型化等，因此其常常與高度、緯度、地方時、太陽輻射和地磁活動強(qiáng)度等因素有關(guān)(Doornbos et al., 2008; 陳旭杏等, 2013; Shi et al., 2015).本文利用加速度計數(shù)據(jù)反演密度時，分段加速度計校正參數(shù)的誤差雖然會導(dǎo)致密度反演出現(xiàn)bias，但是從長時間來看加速度計標(biāo)校誤差是隨機(jī)的(Sutton et al., 2007)，而衛(wèi)星空氣動力系數(shù)誤差則會導(dǎo)致顯著的bias.一般而言，能量調(diào)節(jié)系數(shù)對空氣動力系數(shù)CA的影響最為明顯，調(diào)節(jié)系數(shù)越大CA越小.對于GRACE衛(wèi)星而言，其調(diào)節(jié)系數(shù)一般認(rèn)為在0.8以上，但是精確確定其調(diào)節(jié)系數(shù)目前仍然是一項挑戰(zhàn)(Sutton, 2009).類似的，對于NRLMSISE00，HASDM等經(jīng)驗?zāi)Ｐ?，其使用的衛(wèi)星阻力觀測數(shù)據(jù)等同樣存在由彈道系數(shù)誤差所導(dǎo)致的bias.

另一方面，ρobs與目前公認(rèn)精度最高的模型HASDM之間的標(biāo)準(zhǔn)差可達(dá)30%.其很大部分可能是由于HASDM采用的低階球諧模型僅能在相對較大時空尺度上對熱層密度的平均誤差進(jìn)行校正，而ρobs能夠反映熱層密度更高時空分辨率上的變化.Lechtenberg 等(2013) 將CHAMP和GRACE加速度計反演得到的密度與HASDM模型進(jìn)行比較時發(fā)現(xiàn)二者總體上非常符合，然而前者能夠反應(yīng)熱層密度變化的諸多高頻變化現(xiàn)象，如：午夜密度極大，地磁極尖區(qū)密度增加，大氣行擾等，而HASDM模型無法探測到以上變化.本文圖10中也顯示出了類似的現(xiàn)象.其次，本文所研究的數(shù)據(jù)時段為2007—2009年，屬于太陽活動低年時期，而Bowman和Storz(2003) 對HASDM模型精度進(jìn)行評價時采用的是2001—2002年的數(shù)據(jù)，屬于太陽活動高年.在太陽活動低年時期由于太陽輻射較小導(dǎo)致熱層密度為整個太陽活動周期內(nèi)最低，因此大氣阻力信號強(qiáng)度較弱.在考慮其他非保守力如SRP等誤差不變的情況下，大氣阻力減小從而會使得反演熱層密度的相對誤差增大.而另一方面，HASDM模型采用低軌空間物體阻力觀測值來校正Jacchia71模型，其在太陽活動低年也同樣存在因阻力信號減弱所致的阻力觀測值相對誤差增大的問題.這也是導(dǎo)致本文統(tǒng)計結(jié)果中ρobs與HASDM模型之間STD較大的原因之一.因此，采用整個太陽活動周期內(nèi)的反演密度數(shù)據(jù)將有利于更準(zhǔn)確地評價反演精度，特別是在不同太陽輻射強(qiáng)度條件下的反演密度精度.另外，對于反演得到的結(jié)果需要進(jìn)行嚴(yán)格的質(zhì)量控制和粗差剔除.現(xiàn)階段本文僅僅根據(jù)THR產(chǎn)品對衛(wèi)星點火時期的數(shù)據(jù)進(jìn)行了刪除，然而對于極區(qū)大氣風(fēng)速誤差等導(dǎo)致的反演粗差仍然參與了比較，這也可能會對統(tǒng)計結(jié)果產(chǎn)生干擾.

6 總結(jié)

本文介紹了利用加速度計數(shù)據(jù)反演熱層密度的模型與方法，并對GRACE衛(wèi)星2007—2009年的數(shù)據(jù)進(jìn)行了處理與分析.首先利用星載GPS數(shù)據(jù)進(jìn)行動力學(xué)精密定軌，并以加速度計數(shù)據(jù)代替非保守力，估計得到了該時段內(nèi)GRACE衛(wèi)星的加速度計標(biāo)校參數(shù).在此基礎(chǔ)上利用校正后的加速度計數(shù)據(jù)扣除太陽輻射光壓，并根據(jù)Sentman方程計算衛(wèi)星空氣動力系數(shù)，反演得到了2007—2009年沿GRACE衛(wèi)星軌道的熱層大氣密度.本文研究的主要結(jié)論如下：

(1) 對不同太陽活動強(qiáng)度時期的衛(wèi)星所受空氣動力、太陽輻射壓力，地球反照輻射等非保守力進(jìn)行比較，結(jié)果表明在太陽活動高年衛(wèi)星阻力信號較強(qiáng)，其在X軸分量比SRP高出一個量級；而太陽活動低年時期由于衛(wèi)星阻力信號變?nèi)酰赡軐?dǎo)致反演密度誤差增加.

(2) 利用GRACE衛(wèi)星2007—2009年的GPS與加速度計數(shù)據(jù)進(jìn)行動力學(xué)定軌處理，通過軌道精度分析對加速度計校正精度進(jìn)行了評價.結(jié)果表明，采用本文處理策略，利用GPS和加速度計數(shù)據(jù)定軌結(jié)果與JPL PSO在同一精度量級，GRACE-A/B星與JPL PSO軌道差值3D RMS平均值分別為2.6 cm和2.3 cm，SLR軌道檢核的RMS平均值分別為2.4 cm和2.6 cm.通過該方法估計得到的三個方向的偏差參數(shù)變化都非常平緩，X方向偏差參數(shù)的年均變化在10 nm·s-2量級.

(3) 分析了利用GRACE衛(wèi)星加速度計進(jìn)行密度反演的數(shù)據(jù)處理策略，對反演結(jié)果進(jìn)行分析與評估.對GRACE衛(wèi)星2007—2009年的加速度計數(shù)據(jù)進(jìn)行處理得到了該時期內(nèi)沿GRACE衛(wèi)星軌道的熱層大氣密度.將本文反演的結(jié)果與Doornbos解算的結(jié)果，經(jīng)驗密度模型NRLMSISE00和HASDM進(jìn)行對比分析.結(jié)果表明本文反演結(jié)果與Doornbos結(jié)果之間存小于8%的系統(tǒng)性bias，二者STD統(tǒng)計值在10%以內(nèi)，具有較好的符合性.二者之間的差異主要是由不同的加速度標(biāo)校和空氣動力系數(shù)計算策略造成.本文反演結(jié)果較HASDM偏大，bias在4%以內(nèi)，而Doornbos結(jié)果較HASDM模型值偏小，bias在7%以內(nèi)，二者與HASDM模型比較的STD均在30%左右.本文結(jié)果與NRLMSISE00模型之間存在約30%～40%的系統(tǒng)bias，這主要是由于NRLMSISE00模型在2007—2009年太陽活動異常寧靜期間存在顯著模型偏差所致.

致謝 感謝ISDC提供GRACE衛(wèi)星相關(guān)數(shù)據(jù)，以及IGS提供GPS精密軌道和鐘差產(chǎn)品.感謝Doornbos和Sutton博士提供基于GRACE衛(wèi)星加速度計反演的密度數(shù)據(jù).感謝各基金項目對論文研究的支持.

Bettadpur S V. 2004. Recommendation for a-priori bias & scale parameters for level-1B ACC data (Release 00). GRACE TN-04-02.

Bettadpur S V. 2012. GRACE Product Specification Document, GRACE 327-720. The University of Texas at Austin.

Bezděk A. 2010. Calibration of accelerometers aboard GRACE satellites by comparison with POD-based nongravitational accelerations.JournalofGeodynamics, 50(5): 410-423.

Bhanderi D D V, Bak T. 2005. Modeling Earth Albedo for Satellites in Earth Orbit.∥ AIAA Guidance, Navigation, and Control Conference and Exhibit. San Francisco, California: AIAA.

Bowman B R, Storz M F. 2003. High Accuracy Satellite Drag Model (HASDM) review.AdvancesintheAstronauticalSciences, 116: 1943-1952. Bruinsma S, Biancale R. 2003. Total densities derived from accelerometer data.JournalofSpacecraftandRockets, 40(2): 230-236. Bruinsma S, Tamagnan D, Biancale R. 2004. Atmospheric densities derived from CHAMP/STAR accelerometer observations.PlanetaryandSpaceScience, 52(4): 297-312. Bruinsma S, Forbes J M, Nerem R S, et al. 2006. Thermosphere density response to the 20-21 November 2003 solar and geomagnetic storm from CHAMP and GRACE accelerometer data.JournalofGeophysicalResearch:SpacePhysics, 111(A6): A06303.

Bruinsma S, Biancale R, Perosanz F. 2007. Calibration parameters of the CHAMP and GRACE accelerometers. Poster presented at the GRACE Science Team Meeting.

Bruinsma S L, Forbes J M. 2007. Global observation of traveling atmospheric disturbances (TADs) in the thermosphere.GeophysicalResearchLetters, 34(14): L14103. Bruinsma S L, Forbes J M. 2008. Medium- to large-scale density variability as observed by CHAMP.SpaceWeather, 6(8): S08002. Bruinsma S L, Forbes J M. 2009. Properties of traveling atmospheric disturbances (TADs) inferred from CHAMP accelerometer observations.AdvancesinSpaceResearch, 43(3): 369-376.

Bruinsma S L, Forbes J M. 2010. Anomalous behavior of the thermosphere during solar minimum observed by CHAMP and GRACE.JournalofGeophysicalResearch:SpacePhysics, 115(A11): A11323.

Case K, Kruizinga G L H, Wu S C. 2010. GRACE Level 1B Data Product User Handbook, JPL Publication D-22027.

Chen X X, Hu X, Xiao C Y, et al. 2013. Correction method of the low earth orbital neutral density prediction based on the satellites data and NRLMSISE-00 model.ChineseJ.Geophys. (in Chinese), 56(10): 3246-3254, doi: 10.6038/cjg20131003.

Cook G E. 1965. Satellite drag coefficients.PlanetaryandSpaceScience, 13(10): 929-946.

Dahle C, Flechtner F, Gruber C, et al. 2012. GFZ GRACE Level-2 Processing Standards Document for Level-2 Product Release 0005. Deutsches GeoForschungsZentrum GFZ.

Doornbos E, Klinkrad H, Visser P. 2008. Use of two-line element data for thermosphere neutral density model calibration.AdvancesinSpaceResearch, 41(7): 1115-1122.

Doornbos E, van Den Ijssel J, Luehr H, et al. 2010. Neutral density and crosswind determination from arbitrarily oriented multiaxis accelerometers on satellites.JournalofSpacecraftandRockets, 47(4): 580-589.

Doornbos E. 2012. Thermospheric Density and Wind Determination from Satellite Dynamics. Berlin Heidelberg: Springer.

Emmert J, Lean J, Picone J. 2010. Record-low thermospheric density during the 2008 solar minimum.GeophysicalResearchLetters, 37(12): 245-269.Frommknecht B, Oberndorfer H, Flechtner F, et al. 2003. Integrated sensor analysis for GRACE-development and validation.AdvancesinGeosciences, 1: 57-63. Gerlach C, Sneeuw N, Visser P, et al. 2003. CHAMP gravity field recovery with the energy balance approach: first results.∥ First CHAMP Mission Results for Gravity, Magnetic and Atmospheric Studies. Berlin Heidelberg: Springer, 134-139.

Guo J P, Wan W X, Forbes J M, et al. 2007. Effects of solar variability on thermosphere density from CHAMP accelerometer data.JournalofGeophysicalResearch:SpacePhysics, 112(A10): A10308.

Guo J P, Wan W X, Forbes J M, et al. 2008. Interannual and latitudinal variability of the thermosphere density annual harmonics.JournalofGeophysicalResearch:SpacePhysics, 113(A8): A08301.

Koppenwallner G. 2008. Comment on special section: New perspectives on the satellite drag environments of earth, mars, and venus.JournalofSpacecraftandRockets, 45(6): 1324-1327. Lechtenberg T F, McLaughlin C A, Locke T, et al. 2013. Thermospheric density variations: observability using precision satellite orbits and effects on orbit propagation.SpaceWeather, 11(1): 34-35.

Lei J H, Thayer J P, Lu G, et al. 2011. Rapid recovery of thermosphere density during the October 2003 geomagnetic storms.JournalofGeophysicalResearch:SpacePhysics, 116(A3): A03306.

Liu H X, Doornbos E, Yamamoto M, et al. 2011. Strong thermospheric cooling during the 2009 major stratosphere warming.GeophysicalResearchLetters, 38(12): L12102. Marcos F A, Bowman B R, Sheehan R E. 2006. Accuracy of Earth′s thermospheric neutral density models.∥ AIAA/AAS Astrodynamics Specialist Conference and Exhibit, Guidance, Navigation, and Control and Co-located Conferences. Keystone, CO.: AIAA.

McCarthy D D, Petit G. 2004. IERS conventions 2003. IERS Technical Note No. 32. Frankfurt am Main: Verlag des Bundesamts für Kartographie und Geod?sie.

Moe K, Moe M M, Rice C J. 2004. Simultaneous analysis of multi-instrument satellite measurements of atmospheric density.JournalofSpacecraftandRockets, 41(5): 849-853.

Moe K, Moe M M. 2005. Gas-surface interactions and satellite drag coefficients.PlanetaryandSpaceScience, 53(8): 793-801.

Montenbruck O, Gill E. 2000. Satellite Orbits: Models, Methods and Applications. Heidelberg: Springer-Verlag.

Picone J M, Hedin A E, Drob D P, et al. 2002. NRLMSISE-00 empirical model of the atmosphere: Statistical comparisons and scientific issues.JournalofGeophysicalResearch:SpacePhysics, 107(A12): SIA 15-1-SIA 15-16.

Rim H J. 1992. TOPEX Orbit Determination Using GPS Tracking System [Ph. D. thesis]. Austin: The University of Texas.

Sentman L H. 1961. Free Molecule Flow Theory and its Application to the Determination of Aerodynamic Forces. DTIC Document.

Shi C, Li W W, Li M, et al. 2015. Calibrating the scale of the NRLMSISE00 model during solar maximum using the two line elements dataset.AdvancesinSpaceResearch, 56(1): 1-9. Solomon S C, Woods T N, Didkovsky L V, et al. 2010. Anomalously low

solar extreme-ultraviolet irradiance and thermospheric density during solar minimum.GeophysicalResearchLetters, 37(16): L16103. Storz M F Bowman B R, Branson M J I, et al. 2005. High accuracy satellite drag model (HASDM).AdvancesinSpaceResearch, 36(12): 2497-2505.

Sutton E K, Forbes J M, Nerem R S. 2005. Global thermospheric neutral density and wind response to the severe 2003 geomagnetic storms from CHAMP accelerometer data.JournalofGeophysicalResearch:SpacePhysics, 110(A9): A09S40.

Sutton E K, Forbes J M, Nerem R S, et al. 2006. Neutral density response to the solar flares of October and November, 2003.GeophysicalResearchLetters, 33(22): L22101.

Sutton E K, Nerem R S, Forbes J M. 2007. Density and winds in the thermosphere deduced from accelerometer data.JournalofSpacecraftandRockets, 44(6): 1210-1219.

Sutton E K. 2009. Normalized force coefficients for satellites with elongated shapes.JournalofSpacecraftandRockets, 46(1): 112-116.

Van Helleputte T, Visser P. 2008. GPS based orbit determination using accelerometer data.AerospaceScienceandTechnology, 12(6): 478-484. Van Helleputte T, Doornbos E, Visser P. 2009. CHAMP and GRACE accelerometer calibration by GPS-based orbit determination.AdvancesinSpaceResearch, 43(12): 1890-1896. Wang H B, Zhao C Y. 2008. Use CHAMP/STAR accelerometer data to evaluate atmospheric density models during solar maximum year.ActaAstronomicaSinica(in Chinese), 49(2): 168-178. Xiong Y Q, Wang H B, Zhao C Y, et al. 2011a. A rapid and effective calibration method for CHAMP accelerometer.ProgressinAstronomy(in Chinese), 29(2): 228-237. Xiong Y Q, Wang H B, Xu X L, et al. 2011b. The calibration research for GRACE accelerometer.ScientiaSinicaPhysica,Mechanica&Astronomica(in Chinese), 41(11): 1319-1327.

Xu T H, Yang Y X. 2004. Calibration for CHAMP accelerometry data based on known Earth gravity field model.ActaGeodaeticaetCartographicSinica(in Chinese), 33(3): 200-204.

附中文參考文獻(xiàn)

陳旭杏, 胡雄, 肖存英等. 2013. 基于衛(wèi)星數(shù)據(jù)和NRLMSISE-00模型的低軌道大氣密度預(yù)報修正方法. 地球物理學(xué)報, 56(10): 3246-3254, doi: 10.6038/cjg20131003.

汪宏波, 趙長印. 2008. 用CHAMP加速儀數(shù)據(jù)校驗太陽活動峰年的大氣模型精度. 天文學(xué)報, 49(2): 168-178.

熊永清, 汪宏波, 趙長印等. 2011a. CHAMP加速儀資料的快速校標(biāo)研究. 天文學(xué)進(jìn)展, 29(2): 228-237.

熊永清, 汪宏波, 許曉麗等. 2011b. GRACE加速儀資料定標(biāo)研究. 中國科學(xué): 物理學(xué)、力學(xué)、天文學(xué), 41(11): 1319-1327.

徐天河, 楊元喜. 2004. 利用現(xiàn)有重力場模型求定CHAMP衛(wèi)星加速度計修正參數(shù). 測繪學(xué)報, 33(3): 200-204.

(本文編輯 汪海英)

Thermosphere mass density derivation using on-board accelerometer observations from GRACE satellites

LI Wen-Wen1, LI Min1,2*, SHI Chuang1,2, ZHAO Qi-Le1,2

1ResearchCenterofGNSS,WuhanUniversity,Wuhan430079，China2CollaborativeInnovationCenterforGeospatialTechnology,Wuhan430079,China

This paper mainly discusses the thermosphere mass density derivation using the GRACE on-board accelerometer data. The precision of the accelerometer calibration as well as the satellite aerodynamics model have significant impact on the derivation. In this paper the accelerometer calibration is realized through precise orbit determination (POD) process by using the on-board GPS and accelerometer data, while the aerodynamic coefficients of the GRACE twin satellites are calculated using the Sentman equations. The GRACE accelerometer data during year 2007—2009 are processed and the densities along the satellites′ orbit during this period are obtained. To assess the precision of the accelerometer-derived densities, densities from the GRACE-A satellite during 2007—2009 are compared to the results from Doornbos (2012), as well as the empirical NRLMSISE00 and HASDM models. The results indicate that the densities from this paper are systematically about 8% larger than those from Doornbos (2012), while the STD between these two independent dataset is under 10%. Their differences are mainly attributed to the different strategies adopted for accelerometer calibration and aerodynamic coefficients calculation. The accelerometer-derived densities in this paper are systematically larger than the HASDM model with a bias about 4%, while the densities from Doornbos are 4%～7% lower. The STDs of the derived densities and the results of Doornbos with respect to HASDM model are both around 30%. The bias and STD between the accelerometer-derived densities and the empirical NRLMSISE00 model are 30%～40% and around 30%, respectively.

Thermosphere mass density; GRACE; Accelerometer; HASDM; NRLMSISE00; Aerodynamic

10.6038/cjg20160903.

國家自然科學(xué)基金項目(41274049，41375041，41274049，41375041，41231174)，國家高技術(shù)研究發(fā)展計劃(2014AA123101)，地球空間環(huán)境與大地測量教育部重點實驗室開放研究基金(15-02-05)資助.

李文文，男，1990年生，在讀博士，主要研究方向包括星載GPS低軌衛(wèi)星精密定軌，利用低軌空間物體數(shù)據(jù)反演高層大氣密度等.

E-mail:cheeselee@whu.edu.cn

10.6038/cjg20160903

P228

2015-08-05，2016-07-19收修定稿

李文文, 李敏, 施闖等. 2016. 基于GRACE星載加速度計數(shù)據(jù)的熱層密度反演. 地球物理學(xué)報，59(9):3159-3174,

Li W W, Li M, Shi C, et al. 2016. Thermosphere mass density derivation using on-board accelerometer observations from GRACE satellites.ChineseJ.Geophys. (in Chinese),59(9):3159-3174,doi:10.6038/cjg20160903.

*通訊作者 李敏，男，1983年生，2011年畢業(yè)于武漢大學(xué)，獲博士學(xué)位，副教授，目前研究方向包括利用低軌空間物體數(shù)據(jù)反演高層大氣密度，低軌衛(wèi)星精密軌道確定，GNSS精密定位等.E-mail：limin@whu.edu.cn