2323

1. 北京郵電大學(xué)，北京 100876 2. 北京航天飛行控制中心，北京 100094 3. 航天飛行動力學(xué)技術(shù)重點實驗室，北京 100094

衛(wèi)星姿態(tài)描述的是衛(wèi)星本體坐標(biāo)系在參考系的指向[1-2]。在軌航天器均有特定任務(wù)，在飛行中對飛行姿態(tài)及精度均有明確要求。已有文獻分析顯示，衛(wèi)星軌道姿態(tài)變化對星載合成孔徑雷達(SAR)成像效果有較大影響。在一定的成像條件下，衛(wèi)星三軸姿態(tài)變化需滿足相應(yīng)的精度指標(biāo)，才有可能達到既定的SAR成像精度[3]。通信衛(wèi)星則需要基于姿態(tài)信息將星載天線始終指向地面接收天線。隨著中國航天和通信事業(yè)的發(fā)展，點波束天線以其靈活性好、增益高、安全保密性強等特點，獲得了越來越廣泛的應(yīng)用。但是，點波束地面覆蓋區(qū)域是由點波束天線指向確定的。如果衛(wèi)星姿態(tài)俯仰角發(fā)生1°的偏差，點波束中心經(jīng)度偏移可達十幾度，點波束將完全偏離原覆蓋區(qū)，從而直接導(dǎo)致任務(wù)失敗[4]。因此，航天器姿態(tài)是在軌航天器正常工作的基礎(chǔ)性信息。

目前，航天器的姿態(tài)確定由姿態(tài)敏感器硬件裝置和和姿態(tài)確定軟件算法結(jié)合實現(xiàn)[5-7]。常用器件包括激光陀螺儀、星敏感器及太陽敏感器等。航天器姿態(tài)敏感器在空間中的指向方位隨星體的運行而不斷改變, 因此常會受到太陽光的干擾,導(dǎo)致姿態(tài)敏感器在受干擾期間無法正常工作[8]。并且無論激光陀螺儀，或是星敏感器、太陽敏感器，均屬于星載設(shè)備，一旦發(fā)生故障則難以維修，甚至造成整星失效。開展地基航天器姿態(tài)測量技術(shù)研究，則可以避免星載敏感器失靈或者損壞后將無法繼續(xù)進行航天器姿態(tài)測量的風(fēng)險。

本文利用地基同波束干涉測量技術(shù)[9-11]，獲取航天器不同位置下行天線到地面測站的高精度距離差測量結(jié)果，進而反演解算航天器姿態(tài)信息。同波束干涉測量的基本原理是利用天線的主波束同時接收角距很小的兩個或多個探測信號，得到兩個或多個探測信號的相關(guān)相位并進行差分解算。由于探測角距相近，可以消除傳播路徑以及觀測裝置的絕大部分影響，得到的差分時延數(shù)據(jù)比傳統(tǒng)的VLBI時延測量精度更高。從測量原理上，同波束VLBI對兩個探測信號的相對位置有較強的約束。同波束干涉測量技術(shù)已經(jīng)在深空導(dǎo)航相對測量中得到了充分應(yīng)用[12-14]。

基于地基同波束干涉測量，首先建立“星載發(fā)射天線-地面接收天線”干涉距離差與航天器姿態(tài)的函數(shù)關(guān)系，給出基于干涉距離的航天器姿態(tài)方程及解算方法，最后通過數(shù)值仿真評估基于地基同波束干涉測量的航天器姿態(tài)解算精度。該方法不但適用于大型在軌航天器姿態(tài)測量，亦適用于不允許或不能安裝傳感器、采用地基觀測設(shè)備開展測量的飛行器姿態(tài)測量。

1 基于同波束干涉測量的航天器姿態(tài)測量

1.1 同波束干涉測量原理

假設(shè)星載天線為n=1,2,…,N，地面天線為m=1,2,…,M。由于星載天線位置變化將帶來星載天線與地面天線間距離，即信號傳播幾何延遲的變化，因此定義星載天線n與地面天線m間幾何延遲為τnm；由于同一衛(wèi)星發(fā)射的無線電信號傳播至同一測站中，傳播介質(zhì)基本相同，因此定義所有星載天線與地面天線m間電離層延遲、大氣延遲、設(shè)備延遲分別為τionom、τtropm、τeqm；針對同一地面天線m，定義其接收鐘差為τclockm。

假設(shè)星載天線n的發(fā)射信號為Xn(t)，則經(jīng)過傳播至地面天線m的幾何延遲τnm、電離層延遲τionom、大氣延遲τtropm、設(shè)備延遲τeqm、測站鐘差τclockm等，地面天線m接收到的信號為：

Xn(t-τnm-τionom-τtropm-τeqm-τclockm)

星載天線n=1,2至地面天線m=1,2的時延差可以利用干涉測量處理，得到星載天線1到地面天線1，2之間的干涉時延差為：

τ1=(τ12+τiono2+τtrop2+τeq2+τclock2)-

(τ11+τiono1+τtrop1+τeq1+τclock1)(1)

星載天線2到地面天線1，2之間的干涉時延差為：

τ2=(τ22+τiono2+τtrop2+τeq2+τclock2)-

(τ21+τiono1+τtrop1+τeq1+τclock1)(2)

其中利用單星載發(fā)射天線下行信號開展干涉測量的原理如圖 1所示，利用雙星載發(fā)射天線下行信號開展同波束干涉測量的原理如圖 2所示。因此星載天線n=1,2至地面天線m=1,2的時延差的差分，即同波束干涉時延為：

Δτ=(τ12-τ11)-(τ22-τ21)(3)

圖1 單星載發(fā)射天線的干涉測量原理Fig.1 Interferometry with single transmitting antenna

圖2 雙星載發(fā)射天線的同波束干涉測量原理Fig.2 Same-bean interferometry with two transmitting antennas

1.2 基于同波束干涉測量的航天器姿態(tài)測量數(shù)學(xué)模型

不失一般性，基于3個星載天線和3個地面測站開展同波束干涉測量，建立基于干涉時延的航天器姿態(tài)解算數(shù)學(xué)模型。記地固坐標(biāo)系為直角坐標(biāo)系O-XYZ。在該坐標(biāo)系中，已知地面測站坐標(biāo)分別為A、B、C，航天器坐標(biāo)為O′。以O(shè)′為原點建立航天器軌道直角坐標(biāo)系O′-X′Y′Z′，其中X′軸指向軌道運行方向，Z′軸指向地心，Y′軸與X、Z′軸成右手直角坐標(biāo)。以O(shè)′為原點建立航天器本體坐標(biāo)系，其中星載天線位置分別為R′、S′和T′。航天器歐拉角姿態(tài)的偏航角、俯仰角、翻滾角均為0時，軌道直角坐標(biāo)系與本體直角坐標(biāo)系重合。記航天器歐拉角姿態(tài)的偏航角、俯仰角、翻滾角分別為x、y、z。航天器經(jīng)過偏航角、俯仰角、翻滾角旋轉(zhuǎn)后，星載天線在地固坐標(biāo)系中的坐標(biāo)為R、S、T。利用同波束干涉測量，得到R、S到達A、B的距離差a，R、T到達A、B的距離差b，R、S到達A、C的距離差c，由此解算R′、S′、T′旋轉(zhuǎn)至R、S、T的偏航角x、俯仰角y、翻滾角z。

根據(jù)該原理，建立數(shù)學(xué)模型為：

(4)

式中：

R=Ro×Rx×Ry×Rz×R′+O′

S=Ro×Rx×Ry×Rz×S′+O′

T=Ro×Rx×Ry×Rz×T′+O′

其中：Ro為軌道直角坐標(biāo)系至地固直角坐標(biāo)系的旋轉(zhuǎn)矩陣，由航天器運行坐標(biāo)與速度唯一確定；Rx、Ry、Rz分別為由偏航角、俯仰角、翻滾角帶來的本體直角坐標(biāo)系至軌道直角坐標(biāo)系的旋轉(zhuǎn)矩陣，

1.3 基于干涉時延的航天器姿態(tài)解算

根據(jù)上述分析，利用同波束干涉測量結(jié)果解算航天器姿態(tài)信息的等效為式(4)模型的解算。對于非線性方程(4)，本節(jié)給出牛頓迭代法解算方法[15-16]。

首先利用泰勒級數(shù)對式(4)在旋轉(zhuǎn)角θ=[x,y,z]T處進行展開。為表示方便，令f(θ)=f(x,y,z)=[f1(x,y,z),f2(x,y,z),f3(x,y,z)]T，其中:

f1(x,y,z)=(|R-A|-|R-B|)-

(|S-A|-|S-B|) (5)

f2(x,y,z)=(|R-A|-|R-B|)-

(|T-A|-|T-B|) (6)

f3(x,y,z)=(|R-A|-|R-C|)-

(|S-A|-|S-C|) (7)

給定旋轉(zhuǎn)角初始估計θ0=[x0,y0,z0]T，將f(θ)在初始值θ0處進行一階泰勒級數(shù)展開可得

f(θ)=f(θ0)+J(θ-θ0)+O(‖θ-θ0‖)(8)

f(θ0)+J(θ-θ0)為f(θ)的一階近似，為雅克比矩陣J，O(‖θ-θ0‖)表示高階項，其中

忽略高階項，可以根據(jù)距離差方程組得一階近似方程

f(θ)=f(θ0)+J(θ-θ0)=d(9)

經(jīng)整理可得

J(θ-θ0)=d-f(θ0)(10)

對方程(10)兩邊同時左乘矩陣J的轉(zhuǎn)置可得

JTJ(θ-θ0)=JT[d-f(θ0)](11)

由方程(11)可得

(12)

其中旋轉(zhuǎn)角初始估計θ0=[x0,y0,z0]T及f(θ0)、距離差矢量d=[a,b,c]T已知。問題的關(guān)鍵是如何解算出雅克比矩陣J，詳細(xì)推導(dǎo)過程如下：

由定義可知距離項|R-A|，|R-B|，|R-C|，|S-A|，|S-B|，|S-C|，|T-A|，|T-B|是關(guān)于旋轉(zhuǎn)角θ=x,y,zT的復(fù)合函數(shù)，需要利用鏈?zhǔn)角髮?dǎo)法則對這些項進行求導(dǎo)。令ρA,B=(A-B)/|A-B|，Ψ=RxRyRz，可得

其中:

由此可得

(13)

其中:

RAB=(ρR,A-ρR,B)T[ΨxR′,ΨyR′,ΨzR′]

SAB=(ρS,A-ρS,B)T[ΨxS′,ΨyS′,ΨzS′]

TAB=(ρT,A-ρT,B)T[ΨxT′,ΨyT′,ΨzT′]

RAC=(ρR,A-ρR,C)T[ΨxR′,ΨyR′,ΨzR′]

SAC=(ρS,A-ρS,C)T[ΨxS′,ΨyS′,ΨzS′]

至此，解算得到雅克比矩陣J。

1.4 方程可解性分析

理論上，方程可解的充分條件是矩陣JTJ可逆(滿秩)，而JTJ滿秩的條件是雅克比矩陣J必須滿秩。受地面測站和航天器下行天線相對位置影響，雅克比矩陣J存在不滿秩的可能，將導(dǎo)致解算過程中出現(xiàn)迭代發(fā)散。同時，如果雅克比矩陣J條件數(shù)較大時，亦存在解算過程中精度因子過大的情況。

通過分析可知，雅克比矩陣J滿秩須滿足航天器上R，S，T中任意兩點與測站A，B，C中任意兩點不能共線等條件；而雅克比矩陣J條件數(shù)理論分析則較為復(fù)雜，工程實現(xiàn)中可以通過設(shè)定精度因子門限選擇應(yīng)用。

1.5 解算精度因子

在姿態(tài)方程可解條件下，由式(12)可知航天器姿態(tài)解算誤差協(xié)方差矩陣為(JTJ)-1，令

則航天器姿態(tài)解算精度(即精度因子)可以定義為

1.6 解算精度與最小觀測俯仰角關(guān)系理論分析

以星載天線R，S到達地面天線A，B的距離差d為例，航天器姿態(tài)變化引起的距離差d的變化值為Δd。在觀測俯仰角較大時，距離差d很小，變化值Δd相對于距離差d的比值Δd/d較大，因此在解算過程中能夠以較高的精度得到航天器姿態(tài)的變化，即所獲得航天器姿態(tài)解算誤差較小；相反，在觀測俯仰角較小時，由于距離差d很大，所獲得航天器姿態(tài)誤差較大。實際測量中，在不同的觀測俯仰角由姿態(tài)變化引起的距離差變化值Δd存在差異，但是這種差異引起的姿態(tài)變化誤差相對于距離差d引起的誤差可以忽略不計。

2 地基航天器姿態(tài)測量數(shù)值仿真

首先模擬在軌航天器的軌道仿真，開展同波束干涉測量仿真實驗；然后由基于干涉時延的航天器姿態(tài)測量數(shù)學(xué)模型和航天器姿態(tài)解算方法，迭代解算航天器姿態(tài)信息。

不失一般性，假設(shè)航天器軌道高度為400 km，軌道傾角為40°，偏心率為0°；航天器理論姿態(tài)(歐拉角)為偏航角5°、俯仰角-3°、翻滾角2°。星載天線在航天器本體直角坐標(biāo)系中的坐標(biāo)為R′=[50,0,0]，S′=[0,50,0]，T′=[0,-50,0]。

地面主站坐標(biāo)為：東經(jīng)110°、北緯25°、高程0 m，基線長度1000 km。根據(jù)航天器軌道和主站坐標(biāo)位置信息，得到航天器運行星下點及主站位置示意如圖3所示。

圖3 航天器運行星下點及主站位置示意Fig.3 Diagram of star points in spacecraft operation and main station coordinates

利用航天器軌道信息，得到星載天線至地面天線距離差模擬值，模擬值由理論值疊加測量噪聲得到，其中同波束干涉頻率設(shè)定為S波段，相位測量噪聲取方差為1/50周的高斯白噪聲。

航天器姿態(tài)信息解算時，設(shè)定航天器偏航角、俯仰角、翻滾角初值均為0°，進行航天器姿態(tài)迭代解算。迭代次數(shù)超過10或迭代改進量小于0.01°時迭代強制結(jié)束。由于雅克比矩陣不滿秩情況的存在，在解算過程中會出現(xiàn)迭代發(fā)散、不斷進行迭代和精度因子過大的情況，需要對其進行判斷。

2.1 迭代收斂判斷

設(shè)定最小觀測俯仰角為5°、測量時間間隔為1 s，則在2天時間內(nèi)，地面測站對航天器的觀測時長約2 850 s。該觀測時間中，航天器姿態(tài)解算誤差與迭代次數(shù)仿真結(jié)果如圖 4所示。由于雅克比矩陣存在不滿秩或條件數(shù)較大的情況，因此解算過程中可能迭代發(fā)散(由于本仿真迭代次數(shù)超過10將強制結(jié)束，因此迭代次數(shù)超過10時可以理解為迭代發(fā)散)或迭代次數(shù)較多。根據(jù)圖 4可以知道，剔除迭代次數(shù)超過10(迭代發(fā)散)的結(jié)果，亦存在部分較大的解算誤差。這是因為雅可比矩陣條件數(shù)較大導(dǎo)致。因此需進一步考慮精度因子約束，進一步提高解算結(jié)果的有效性。

圖4 迭代次數(shù)與解算誤差仿真結(jié)果Fig.4 Simulation results of calculation error and iteration numbers

2.2 精度因子與航天器姿態(tài)解算誤差關(guān)系

剔除迭代次數(shù)超過10的解算結(jié)果，設(shè)定不同精度因子限制，分別得到航天器姿態(tài)解算誤差與精度因子的關(guān)系。譬如，設(shè)置精度因子小于0.4(如圖 5所示)，航天器姿態(tài)部分解算結(jié)果仍存在較大誤差，其原因是雅可比矩陣條件數(shù)較大導(dǎo)致解算結(jié)果偏差較大。

經(jīng)試驗測試，當(dāng)精度因子小于0.19時，解算誤差小于0.04°(如圖 6所示)。而且此時解算誤差與精度因子基本成正比關(guān)系。因此可將0.19作為精度因子與解算誤差關(guān)系的一個臨界經(jīng)驗值。

圖5 精度因子σ與解算誤差仿真結(jié)果(σ<0.4)Fig.5 Simulation results of calculation error and precision factor (σ<0.4)

圖6 精度因子σ與解算誤差仿真結(jié)果(σ<0.19)Fig.6 Simulation results of calculation error and precision factor (σ<0.19)

2.3 最小觀測俯仰角與航天器姿態(tài)解算誤差關(guān)系

選擇迭代次數(shù)不超過10且精度因子小于0.19的解算結(jié)果，得到航天器姿態(tài)解算誤差與最小觀測俯仰角(即3個測站觀測航天器俯仰角的最小值)的關(guān)系如圖 7所示。

可以看出，航天器姿態(tài)解算誤差和最小觀測俯仰角基本呈反比關(guān)系，即3個測站觀測航天器的最小俯仰角越大，航天器姿態(tài)(歐拉角)誤差越小。

2.4 測量誤差精度

經(jīng)仿真試驗測試，當(dāng)干涉頻率為S頻段(約2.2 GHz)、干涉相位提取精度為1/100周時，航天器姿態(tài)解算精度可以達到0.001°量級；干涉相位提取精度為1/50周時，航天器姿態(tài)解算精度可以達到0.01°量級。

3 結(jié)束語

針對航天器姿態(tài)的地基測量需求，本文建立了基于同波束干涉時延的航天器姿態(tài)測量數(shù)學(xué)模型。針對該數(shù)學(xué)模型，推導(dǎo)了航天器姿態(tài)解算方法，給出了解算精度因子，初步分析了該模型的可解性。通過模擬在軌航天器軌道運行過程，針對航天器姿態(tài)測量進行了數(shù)值仿真試驗和誤差分析。結(jié)果顯示，利用3個地基測站針對常見空間站軌道運行的航天器3個下行天線進行同波束干涉測量，輔助精度因子約束，可以有效解算該航天器姿態(tài)信息。而且3個測站觀測航天器的最小俯仰角越大，航天器姿態(tài)解算誤差越小。作為地基測量，本文方法可以作為現(xiàn)有星基測量方法的備份。

需要說明的是，基于同波束干涉時延的航天器姿態(tài)測量數(shù)學(xué)模型中由于雅克比矩陣不滿秩導(dǎo)致解算迭代發(fā)散的理論分析還有待進一步開展。此外，本文在航天器姿態(tài)解算中采用了精度因子約束的經(jīng)驗值，亦需要后續(xù)理論分析論證。