阮仁桂，賈小林，馮來平，王 龍，張 奮

1. 西安測繪研究所，陜西 西安 710054； 2. 地理信息工程國家重點實驗室，陜西 西安 710054

北斗三號全球導航衛(wèi)星系統(tǒng)(BDS-3)已于2020年7月31日正式開通運行。得益于獨有的全星座星間鏈路，BDS-3的空間信號精度和定位服務性能已躋身主要GNSS的前列，其中廣播軌道的徑向精度優(yōu)于10 cm[1-2]?？梢灶A見，BDS-3將在高精度定位、時間頻率傳遞、參考框架建立與維持等領域發(fā)揮越來越大的作用[3]。

高精度的衛(wèi)星軌道是實現(xiàn)BDS-3高精度應用的前提。精準確定衛(wèi)星軌道，需要精細地考慮衛(wèi)星受到的各種攝動力[4]。BDS-3衛(wèi)星運行于中、高軌道，力模型的不確定性主要來自太陽光壓、地球輻射等非保守力攝動。為提高定軌精度，文獻[5—7]對BDS-3衛(wèi)星的軌道模型進行了探討，主要工作集中于改進BDS-3衛(wèi)星的太陽光壓模型。針對BDS-3衛(wèi)星平臺為長方體的特點，文獻[5]基于可估計參數(shù)的箱體-兩翼(box-wing,BW)模型[8-9]和頻譜分析方法構建了兩個先驗補償模型以增強五參數(shù)ECOM(extend CODE orbit model)[10]模型(ECOM1)。文獻[6]同樣基于可估參數(shù)的BW模型，針對不同廠家的BDS-3衛(wèi)星建立了分析型光壓模型。為了支持高精度應用[11]，中國衛(wèi)星導航系統(tǒng)管理辦公室發(fā)布了北斗衛(wèi)星的元數(shù)據(jù)[12]，包括幾何和物理參數(shù)。文獻[7]根據(jù)BDS-3衛(wèi)星元數(shù)據(jù)建立的BW模型聯(lián)合ECOM1可以獲得優(yōu)于九參數(shù)的擴展ECOM模型(D4B1)[13]的定軌結果。與太陽光壓相比，地球輻射和天線輻射產生的攝動力要小得多，研究表明二者會引起導航衛(wèi)星軌道徑向的系統(tǒng)性偏差[14-16]，但它們對BDS-3定軌的影響尚未被涉及。

文獻[6—7,17]研究表明，即使考慮了BW模型補償，BDS-3衛(wèi)星軌道的SLR殘差仍然有不可忽視的系統(tǒng)性偏差，這說明BDS-3衛(wèi)星的非保守力模型還有待進一步完善。本文基于BDS-3衛(wèi)星元數(shù)據(jù)和天線輻射功率信息構建地球輻射、天線輻射和BW光壓模型對典型的經驗型太陽光壓模型進行補償，考查各非保守力因素或模型對精密定軌的影響，以期優(yōu)選出適合BDS-3 MEO衛(wèi)星的非保守力建模策略。

1 BDS-3衛(wèi)星的非保守力建模

BDS-3衛(wèi)星受到的非保守力包括[4]：太陽光壓、地球輻射、天線輻射和衛(wèi)星體熱輻射。其中，衛(wèi)星體熱輻射與衛(wèi)星體表面輻射率、絕對溫度有關，難以準確建模且部分被光壓模型吸收[13]，在這里不單獨考慮。

1.1 太陽光壓攝動

太陽光壓是BDS-3衛(wèi)星承受的最大的非保守力。針對導航衛(wèi)星的太陽光壓建模，國內外學者已經積累了非常豐富的經驗[10,13,18-21]，大體上可以分為分析法和經驗法。兩種方法各有優(yōu)劣，為了獲得最佳效果，聯(lián)合采用經驗法和分析法對太陽光壓進行建模逐漸成為一種趨勢。因此，太陽光壓引起的攝動加速度aSRP可表示為

aSRP=aapr+aemp

式中，aapr表示先驗的光壓攝動加速度，通常采用分析模型進行計算，如ROCK模型[18]、BW模型[9]等，其中的模型參數(shù)需事先測定或計算得到；aemp表示經驗模型，如GSPM[20]、ECOM系列[10]、可估計參數(shù)的BW模型[8]。本文將采用BW模型和ECOM系列模型進行太陽光壓建模。

1.1.1BW模型

BW模型的核心思想是將衛(wèi)星視為一個箱體和兩翼的太陽帆板所構成，在給定各表面的幾何、物理參數(shù)及衛(wèi)星姿態(tài)信息的條件下，可以計算出太陽光壓在各個表面產生的攝動加速度，將各個面的加速度疊加就得到太陽光壓引起的總加速度。由于衛(wèi)星體和太陽帆板表面材料的熱輻射特性不同，分別用式(1)和式(2)計算衛(wèi)星體各表面元(i)和太陽帆板(sp)的太陽光壓加速度[9]

(1)

(2)

式中，a*(*=i,sp，下同)表示加速度矢量；A*表示面積；e⊙表示太陽方向矢量；en,*表示表面元(*)的法向矢量；θ*為e⊙和en,*的夾角；α*、δ*和ρ*表示表面元(*)的吸收、漫反射和鏡面反射系數(shù)；m表示衛(wèi)星質量；c為光速；Φ⊙表示太陽在衛(wèi)星所處位置的輻射通量，計算公式為

式中，Φ0≈1367 W/m2表示太陽在地球附近的長期平均輻射通量；1AU表示太陽到地球的平均距離；r⊙和r分別表示太陽和衛(wèi)星的位置矢量。

1.1.2ECOM模型

ECOM模型是歐洲定軌中心(Center of Orbit Determination at Europe,CODE)針對GPS衛(wèi)星提出的經驗型太陽光壓模型[10]，現(xiàn)已被廣泛應用于各種導航衛(wèi)星的精密定軌。其核心思想是將太陽光壓加速度在3個正交方向eD、eY和eB上進行表示

aSRP=D(u)·eD+Y(u)·eY+B(u)·eB

(3)

式中，eD為衛(wèi)星指向太陽的矢量；eY為“日-星-地”平面法向矢量；eB=eD×eY，構成右手系；D(u)、Y(u)和B(u)為對應3個方向的加速度，采用截斷傅里葉級數(shù)表示，其具體計算公式有多個不同的版本[10]，最新的擴展(或改進的)ECOM模型為

(4)

式中，[D0,Y0,B0]、[D2i,c,B2i-1,c]、[D2i,s,B2i-1,s]上進行(i=1,2,…)是擴展ECOM模型參數(shù)，需作為待估參數(shù)進行解算；Δu=u-u⊙,u⊙和u分別表示太陽在軌道面的投影和衛(wèi)星相對于軌道升交點的夾角。

式(4)中eD方向的加速度僅考慮了偶次項，eB方向僅考慮了奇次項。研究表明，對于衛(wèi)星體為正方體的GPS衛(wèi)星，僅考慮{D0,Y0,B0,Bc,Bs}五參數(shù)的簡化ECOM模型，即ECOM1模型，可獲得最優(yōu)的軌道和有關大地參數(shù)的解算結果[10]。對于非正方體的衛(wèi)星，文獻[13]推薦七參數(shù){D0,D2,c,D2,s,Y0,B0,Bc,Bs}的D2B1模型(下文用ECOM2表示)或九參數(shù){D0,D2,c,D2,s,D4,c,D4,s,Y0,B0,Bc,Bs}的D4B1模型[13,21]。

1.2 地球輻射攝動

地球輻射攝動是指地球表面受到太陽光照射后直接反射或紅外輻射對衛(wèi)星引起的輻射壓力。文獻[15]建立了一系列的地球輻射通量模型，包括分析型、數(shù)值型和數(shù)據(jù)型[15]。其中，將地球視為球形朗伯輻射體的分析模型計算簡單，與基于CERES(clouds and Earth 's radiant energy system)數(shù)據(jù)的結果相比，最大誤差不超過15%，平均誤差為1.6%。這也說明準確計算地球輻射通量非常困難。分析模型計算直接反射Φrefl和紅外輻射Φinfr通量的公式分別為

(5)

用Φrefl和Φinfr代替Φ⊙，用地球方向矢量e⊕代替e⊙，就可用式(1)和式(2)計算地球輻射引起的衛(wèi)星體各表面和太陽帆板的攝動加速度。

1.3 天線輻射攝動

天線輻射攝動是衛(wèi)星的導航信號天線持續(xù)對地發(fā)射射頻信號產生的與輻射方向相反的攝動力。天線輻射加速度aAT可由式(6)計算

(6)

式中，P為天線信號發(fā)射功率。

1.4 衛(wèi)星信息

中國衛(wèi)星導航系統(tǒng)管理辦公室提供了BDS-3衛(wèi)星表面幾何和光學參數(shù)文件[12]。各類衛(wèi)星的參數(shù)見表1，其中用CAST和SECM區(qū)分中國空間技術研究院和上海微小衛(wèi)星工程中心研制的衛(wèi)星。由于文獻[12]僅給出了各表面的吸收系數(shù)，沒有提供鏡面反射率。本文根據(jù)文獻[7]，假設各衛(wèi)星各個表面的鏡面反射率為1，僅CAST衛(wèi)星+X和-Z面的鏡面反射率為0；同時假設太陽帆板反面的光學參數(shù)與正面相同，且所有表面材料對紅外輻射的鏡面反射、漫反射和吸收系數(shù)分別為0.1、0.1和0.8。GNSS衛(wèi)星的天線輻射功率可以用大口徑天線進行測量[14]，本文進行天線輻射建模時采用IGS(international GNSS service)提供的數(shù)據(jù)，SECM衛(wèi)星和CAST衛(wèi)星的天線輻射功率分別為280、310 W。

表1 BDS-3衛(wèi)星的幾何和物理參數(shù)Tab.1 Geometric and physical parameters of BDS-3 satellites

姿態(tài)信息是建立分析型光壓模型和地球輻射模型的另一關鍵信息。BDS-3 MEO和IGSO衛(wèi)星在大部分時間里采用名義偏航(NYS)模式：即+Z軸指向地心，+X軸位于“日-星-地”平面并指向太陽側。當衛(wèi)星進入深度地影區(qū)(太陽相對于軌道面的高度角小于3°)且日、星、地三者接近共線時進入控制偏航(CYS)模式：+X軸與日-星-地平面存在一定夾角，具體的姿態(tài)模型見文獻[7]。根據(jù)上文的公式和有關數(shù)據(jù)可得：CAST和SECM衛(wèi)星的太陽光壓加速度分別約為130、70 nm/s2，地球輻射加速度分別約為1.5、0.5 nm/s2，天線輻射加速度約為1.5 nm/s2。

2 試驗和分析

為了分析不同非保守力因素對BDS-3衛(wèi)星精密定軌的影響，收集全球觀測網的數(shù)據(jù)，構建多個非保守力模型進行定軌試驗。

2.1 數(shù)據(jù)及試驗設計

收集2021年年積日001—110期間全球約75個監(jiān)測站的觀測數(shù)據(jù)，測站分布如圖1所示。采用西安測繪研究所的SPODS軟件[22-23]進行定軌解算，定軌弧長為30 h，以每天的12∶00為中心。按照IGS或iGMAS分析中心的日常數(shù)據(jù)處理模式，同時解算軌道、鐘差、測站坐標、對流層延遲等參數(shù)，具體策略見表2。

表2 數(shù)據(jù)處理算法、模型和策略Tab.2 Algorithms, models and strategies for data processing

圖1 本文所用的IGS和iGMAS測站分布Fig.1 Distribution of the selected IGS and iGMAS stations

本文首先關注經驗型太陽光壓模型ECOM1和ECOM2的優(yōu)劣，在此基礎上考查地球輻射、天線輻射和先驗光壓模型的影響。初步分析發(fā)現(xiàn)地球輻射與天線輻射主要作用于徑向方向，對定軌結果的影響相似。為此，設計了14組不同的非保守力模型配置，見表3。模型的命名采用以下規(guī)則：“B”代表BW模型;“Sn”(n=1、2)分別代表ECOM1和ECOM2；“E”代表地球輻射；“A”代表天線輻射；“*”表示“E”和“A”僅用于SECM衛(wèi)星。在不至于混淆的情況下，下文用ECOM1和ECOM2分別代表以ECOM1和ECOM2為經驗光壓模型的一類非保守力模型。

表3 非保守力模型配置表Tab.3 Configuration of the non-conservative force models

2.2 結果分析

通過考查相鄰兩天定軌結果6 h重疊段的互差進行軌道精度的內部評估；分析激光數(shù)據(jù)殘差進行外部檢核?？紤]到IGSO和GEO衛(wèi)星的區(qū)域覆蓋特性，幾何觀測條件較差，本文不對其結果進行分析。

2.2.1 軌道重疊互差分析

在試驗數(shù)據(jù)弧段內，C20、C21所在的軌道和C29、C30所在軌道分別在年積日033—040和083—094期間進入深度地影區(qū)，軌道面內的衛(wèi)星會進入CYS模式。圖2統(tǒng)計了NYS和CYS模式下不同模型得到的BDS-3 MEO衛(wèi)星軌道重疊互差在R、T和N方向的RMS。由圖2可以看出，在NYS模式下，不同力模型同一軌道方向的精度差異不大，3個方向的RMS分別在3、5和5.5 cm左右，只是BS1和BS1EA*在軌道T和N方向的精度似乎較其他模型稍差。力模型相同時，與NYS模式的結果相比，CYS模式下軌道各方向的誤差都明顯增大，增幅普遍超過1 cm；其中ECOM2在R方向的增幅超過3 cm，這說明，ECOM2對姿態(tài)模式更敏感。還可看出，考慮地球輻射、天線輻射和BW模型不能顯著降低軌道重疊互差。

圖2 不同模型得到的軌道重疊互差RMSFig.2 RMS of overlap orbit differences with different force models

2.2.2 激光殘差分析

BDS-3衛(wèi)星都安裝了激光反射棱鏡，國際衛(wèi)星激光測距組織對4顆BDS-3衛(wèi)星(C20、C21、C29和C30)的激光測距數(shù)據(jù)可以用于評估衛(wèi)星軌道徑向精度。在進行激光殘差計算時，各顆衛(wèi)星的激光反射器的質心偏差參數(shù)見表4[12]。值得注意的是，表4中的數(shù)值與文獻[5—7,17]所采用的數(shù)值不完全一致。

表4 BDS-3衛(wèi)星的激光反射器質心偏差參數(shù)Tab.4 Center of mass corrections of the laser array for BDS-3 satellites m

圖3繪出有代表性的6個力模型得到的激光殘差隨星-地-日夾角ε的變化，紅色數(shù)值為紅色擬合直線的斜率?？梢钥闯觯δＰ拖嗤瑫r，衛(wèi)星C20和C21的激光殘差具有更明顯的隨ε變化的趨勢。在經驗光壓模型基礎上，加入地球輻射和天線輻射模型不能減弱激光殘差與ε的相關性；而加入BW模型有助于降低CAST衛(wèi)星激光殘差與ε的相關性。對比圖3(a)、(d)可見，相比于ECOM1，采用ECOM2可明顯減弱激光殘差與ε的相關性。

圖3 激光殘差隨ε的變化Fig.3 Variation of SLR residuals with ε

同一廠家衛(wèi)星的激光殘差具有相似的統(tǒng)計特性，圖4統(tǒng)計了兩類衛(wèi)星激光殘差的STD。由圖4可以看出，在NYS模式下，不同模型的激光殘差STD基本都小于4.5 cm；ECOM2對應的STD更小，基本不超過3.5 cm。在CYS模式下，ECOM1對應的STD普遍略有增大，但都不超過5 cm；而ECOM2對應的STD顯著增大，增幅為2～3 cm?？梢?，在經驗光壓模型基礎上，增加BW模型似可降低激光殘差STD，而加入地球輻射和天線輻射則未顯示出顯著的有益影響。

圖4 衛(wèi)星軌道的SLR殘差STDFig.4 STD of SLR residuals of orbits with different force models

圖5統(tǒng)計了激光殘差的系統(tǒng)偏差?？梢钥闯?，在NYS模式下，僅考慮經驗光壓模型時，不同廠家衛(wèi)星的系統(tǒng)偏差符號相反；在經驗光壓模型基礎上，考慮地球輻射和天線輻射，可顯著降低SECM衛(wèi)星的系統(tǒng)偏差，但卻導致CAST衛(wèi)星的系統(tǒng)偏差顯著增大。例如，在S2EA的結果中，SECM衛(wèi)星的系統(tǒng)偏差非常接近于0，而CAST衛(wèi)星的系統(tǒng)偏差則超過7 cm。此外，增加BW模型或采用ECOM2，有助于降低激光殘差的系統(tǒng)偏差。與S1的結果相比，對于CAST和SECM衛(wèi)星，BS1EA*的系統(tǒng)偏差分別降低13%和88%，S2的系統(tǒng)偏差分別降低15%和45%，BS2EA*的系統(tǒng)偏差分別降低26%和100%。

圖5 衛(wèi)星軌道的SLR殘差的系統(tǒng)偏差Fig.5 SLR residuals biases of orbits with different force models

與NYS模式的結果相比，在CYS模式下，對于CAST衛(wèi)星，各模型得到的系統(tǒng)偏差都普遍增大；對于SECM衛(wèi)星，采用ECOM2模型對應的系統(tǒng)偏差顯著增大，ECOM1模型的結果則未見顯著惡化。

3 結 論

本文基于衛(wèi)星幾何、光學參數(shù)及天線輻射功率等信息，對BDS-3衛(wèi)星非保守力進行建模，采用全球觀測網的數(shù)據(jù)進行定軌試驗，分析在不同的經驗光壓模型(ECOM1和ECOM2)基礎上，增加BW光壓模型、地球輻射和天線輻射模型對定軌結果的影響。結果表明，經驗光壓模型是影響軌道精度的主要因素，在名義偏航模式下，ECOM2模型具有更好的表現(xiàn)，但ECOM1模型對衛(wèi)星的姿態(tài)模式更不敏感。在控制偏航模式下，ECOM2得到的軌道徑向重疊互差RMS和激光殘差STD顯著增大。在經驗光壓模型基礎上，增加地球輻射或天線輻射，可以顯著降低SECM衛(wèi)星激光殘差系統(tǒng)偏差；對CAST衛(wèi)星則適得其反，其原因還有待進一步分析。此外，增加BW模型對于提高軌道精度也是有益的。為了進一步滿足BDS-3高精度應用需求，構建更加準確的非保守力模型，需要衛(wèi)星廠家提供更詳細和準確的衛(wèi)星表面幾何和光學參數(shù)等信息。