星間測量自主導航改進狀態相關黎卡提方程濾波

，

中國空間技術研究院 錢學森空間技術實驗室，北京 100094

分布式衛星系統是由物理上互不相連、共同協作完成同一航天任務的多衛星構成的空間系統[1]，由于其具有對地面依賴性小、自主性強、靈活可變、綜合成本低等特點，成為衛星技術發展的新趨勢。分布式衛星系統實現自主協同需具備的核心和基礎能力是自主導航[2]，其自主導航問題與普通衛星存在較大區別：既需要確定每個衛星的精確空間位置和速度，又需要確定衛星之間的相對位置速度關系。在分布式衛星系統應用任務中，衛星間的相對信息相對于絕對信息更容易獲得，近年來，利用星間相對觀測的自主導航技術備受關注。這些相對觀測量包括視線角、星間相對距離和相對位置矢量等[3-9]。能否基于相對測量序列求解各自衛星的絕對位置和速度，在理論上取決于系統的能觀性能，在應用上涉及高性能定軌濾波算法。

文獻[10-11]研究了基于星間測距的雙星自主導航方法，指出沒有先驗信息，僅依靠相對距離測量的導航系統是虧秩的，并給出了嚴格的數學證明。文獻[12]針對雙星編隊首次提出了僅利用地心慣性坐標系中星間相對位置矢量序列作為觀測量的自主絕對導航方法，并采用局部線性化的方法分析了其能觀性，文章指出此方法在半長軸相等的共面圓軌道情況下不適用，但未給出嚴格的理論證明。文獻[13]對文獻[12]提出的導航方法進行了物理角度的闡述，把雙星系統視作一個巨型加速度計，并對定軌精度進行了分析，并說明該方法除幾個0°傾角的特殊軌道外，在大部分情況下是能觀的。文獻[14]首次對基于星間位置矢量測量的雙星自主導航系統使用非線性系統局部能觀性理論進行分析，通過數學證明指出在兩衛星瞬時地心距不相等的情況下可以定軌。文獻[2]在上述研究基礎上提出了非線性系統局部k階能觀性的概念，并證明了理想二體軌道模型下該自主導航系統局部4～6階能觀的充要條件及其與非線性系統能觀性秩條件成立的關系。文獻[15]考慮能觀性強弱對濾波算法的影響，給出了一個基于相對位置矢量自主導航系統能觀性秩條件的能觀度指標，并據此對無跡濾波(Unscented Kalman Filter，UKF)算法進行改進，提高了導航定位的精度。但由于在算法中，指標只能基于位置速度的估計值而非其真值得到，因此不能真實反映系統能觀度的大小。

衛星動力學模型是個典型的非線性系統，一般采用非線性濾波算法對軌道進行估計。常用的方法有擴展卡爾曼濾波[16-19](Extended Kalman Filter，EKF)和無跡濾波[20-22]等。狀態相關黎卡提方程(State-Dependent Riccati Equation，SDRE)方法[23]通過構造線性結構解決非線性系統的最優控制問題，隨之應用于非線性系統的估計[24]，構成的濾波器被稱為狀態相關黎卡提方程濾波(State-Dependent Riccati Equation Filter，SDREF)。在無隨機游走的確定性系統中，SDREF可以全部獲得系統的非線性信息。這使得可以將非線性系統作為線性系統來處理，不存在線性化誤差，且無需計算Jacobian矩陣。

在本文的研究中，SDREF被用作利用地心慣性坐標系中星間相對位置矢量序列作為觀測量的自主導航估計器，并與EKF算法在濾波步驟、計算負擔和濾波精度等方面進行了比較。本文推導并定義了一種新的在線反映能觀性程度的特征量，提出了基于該特征量對濾波增益陣進行在線調整的改進SDREF算法。數學仿真表明,該改進SDREF算法在保持SDREF計算量小的優點下，導航估計精度明顯優于采用EKF和SDREF算法的結果。

1 SDREF原理

SDREF是一種參數化的濾波方法，它將非線性系統改寫為線性化結構，系統矩陣和觀測矩陣的參數由當前狀態決定，稱為狀態相關系數(State-Dependent Coefficients，SDC)。SDREF與EKF唯一的不同之處在于EKF是在估計狀態附近展開進行線性化近似，而SDREF是將非線性系統做參數化處理，這一轉換過程稱為偽線性化。

考慮如下的連續非線性系統:

(1)

式中：X為n維狀態向量；Y為m維觀測向量；f(·)，h(·)分別為n維和m維非線性函數；w和v分別為n維和m維高斯白噪聲向量，各分量間互不相關，且滿足

(2)

由文獻[25]可知，如果f(·),h(·)可微且滿足f(0)=0，h(0)=0，則非線性系統(1)可以轉化為如下具有SDC的線性形式：

(3)

其狀態方程和觀測方程的系數矩陣以狀態向量為參數。由于沒有進行線性化近似，式(1)轉化為式(3)的過程中不會損失任何信息，可以運用線性系統標準Kalman濾波方法對式(3)進行處理。

首先對式(3)進行離散化，得到

(4)

式中：Φ(Xk)為狀態轉移矩陣；H(Xk)為觀測矩陣；wk和vk分別為離散化的狀態噪聲和觀測噪聲，當采樣間隔Δt較小時滿足

狀態轉移矩陣Φ(Xk)如下：

(5)

在計算時一般取二階近似：

(6)

步驟1：選定濾波初值

步驟2：狀態一步預測

(7)

步驟3：計算一步預測均方誤差

(8)

步驟4：計算濾波增益

(9)

步驟5：計算估計均方誤差

Pk+1=(I-Kk+1Hk+1)Pk+1/k

(10)

步驟6：狀態更新

(11)

從上述過程中可以看出，與EKF算法相比，雖然SDREF沒有線性化損失信息量，但是f(X)和h(X)在改寫時拆解出了一個右乘狀態向量，其對狀態方程和觀測方程的描述有所欠缺。因此SDREF和EKF兩種算法對模型的近似程度高低難以給出明確的定論，后文會通過理論分析和數值仿真對SDREF和EKF算法在相對位置矢量測量自主導航系統中的表現進行探討。

2 系統模型

考慮基于星間相對位置矢量測量的雙星自主導航系統，將兩顆衛星記為衛星A和衛星B，測量關系如圖1所示。

O是地球中心，rA和rB分別是地心到兩衛星的矢量，rAB是衛星A到衛星B的位置矢量。已知衛星A搭載了星敏感器和星間相對測量設備(包括無線電或激光測距儀、光學或雷達測向儀)，可以測量自身在慣性空間中的姿態、衛星A與衛星B間的距離以及衛星B在衛星A本體坐標系中的方位，通過坐標轉換可以得到兩衛星在慣性空間中的相對位置矢量。

本文選取地心赤道慣性坐標系為參考坐標系，坐標原點為地球質心，基準面為赤道面，X軸指向春分點，Z軸指向北極，Y軸與X軸和Z軸構成右手直角坐標系。含J2攝動項的衛星軌道動力學模型為[26]：

(12)

(13)

選取地心赤道慣性坐標系下兩衛星的位置速度向量為系統狀態變量，記

X為12維狀態變量，則帶J2攝動項的系統狀態方程為：

(14)

其中系統噪聲w包括未考慮的高階地球非球形引力、太陽光壓、大氣阻力等攝動力。

觀測方程為：

Y=h(X)+v=rA-rB+v=HX+v

(15)

3 SDREF與EKF計算量對比

在基于相對位置矢量觀測的雙星導航系統中，由于觀測矩陣是線性定常的，SDREF與EKF唯一的區別在狀態轉移矩陣Φk的計算上，前者由當前時刻的參數化矩陣求得，而后者由Jacobian矩陣求得。

對于SDREF，狀態方程(14)可以改寫為：

(16)

記I3×3為I3，則離散后的狀態轉移矩陣為

(17)

對于EKF，狀態方程非線性函數的Jacobian矩陣

(18)

式中：

狀態轉移矩陣

(19)

分析式(17)中SDREF狀態轉移矩陣的形式，可以得到以下定理。

定理1 在以式(14)為系統方程、式(15)為觀測方程的導航系統中，若SDREF滿足以下條件：

1)初始時刻估計向量協方差矩陣P0為對角4×4分塊陣；

2)觀測矩陣H為對角1×4分塊陣；

3)狀態誤差和觀測誤差為各分量互不相關的白噪聲，即Qk與Rk為對角陣。

則后續每一時刻的狀態誤差協方差矩陣Pk都是對角4×4分塊陣。

證明：用數學歸納法進行證明。

1)當k=1時，顯然成立。

一步預測均方誤差陣Pg+1/g為：

(20)

Pg+1/g是對角4×4的分塊陣。

(21)

則濾波增益矩陣Kg+1為：

(22)

Kg+1是對角4×1分塊陣。

將式(20)和式(22)帶入式(10)中，由于Kg+1，Hg+1和Rg+1都是對角分塊陣，所以Pg+1是對角4×4分塊陣。

由式(1)(2)，如果SDREF滿足定理1的條件1)～3)，則每一時刻的誤差協方差矩陣Pk都是對角4×4分塊陣。

在本文分析的自主導航系統中，定理1的條件是比較容易滿足的，此時Pk，Φk，Kk均為對角分塊陣，12維矩陣運算可簡化為4維矩陣運算，計算量約是EKF的12.5%。

4 基于能觀度的改進濾波算法

能觀度是系統能觀性強弱的度量，表征了系統狀態可由輸出反映的程度。非線性系統能觀度的大小會隨著狀態變化而波動，從而對濾波算法產生影響。在能觀度較大的時刻，系統由觀測反映狀態的能力較強，濾波估計誤差較小；反之，能觀度較小的時刻，由觀測反映狀態的能力較弱，濾波估計誤差較大。

目前常用的能觀度評價準則有基于能觀性矩陣的評價準則和基于濾波誤差協方差陣P的評價準則。基于能觀性矩陣的評價準則需要計算矩陣的奇異值，求解比較繁瑣，且要基于狀態真值進行計算。由于在實際濾波中真值是未知的，因此這種評價標準并不能準確反映系統能觀度的大小。基于濾波誤差協方差陣P的評價準則不需要求解能觀性矩陣，計算量比較小，但是在非線性濾波算法中，P陣不能反映真實的估計誤差協方差大小，所以此準則也不能真實反映系統能觀度的情況。

本文從能觀度的物理意義出發，基于濾波算法量測更新過程給出了非線性系統能觀度的另外一種表征方式。

(23)

矩陣的Frobenius范數定義為：若A=(aij)∈Rm×n，規定

(24)

已知Frobenius范數具有如下性質：

式中：λi(ATA)為矩陣ATA的第i個特征值。

依據以上分析，在濾波的時間更新過程中對一步預測估計誤差協方差陣的計算公式進行修改，引入一個縮放參數ρ，將式(8)改寫為：

(25)

由此，依據當前時刻系統能觀度的大小對誤差協方差陣進行了修正，在適當擴大P的基礎上兼顧到能觀度的影響，增強了濾波的穩定性，提高了濾波精度。由式(7)、式(9)～式(11)、式(25)構成的算法本文稱之為改進的狀態相關黎卡提方程濾波算法(Modified State-Dependent Riccati Equation Filter, MSDREF)。

5 仿真與分析

考慮分布式衛星系統在衛星編隊飛行、自主空間在軌服務和非合作目標測量等應用場景中，星間相對距離較近，因此選取兩個軌道相近的衛星。兩衛星的初始軌道參數如表1所示。

表1 兩衛星初始軌道參數

衛星A搭載了星敏感器和星間相對測量設備(包括無線電或激光測距儀、光學或雷達測向儀)，可以測量自身在慣性空間中的姿態、衛星A與衛星B間的距離以及衛星B在衛星A本體坐標系中的方位。利用上述測量量，相對位置矢量在慣性空間中的坐標可由下式解得：

(26)

仿真時間為40 000 s，星敏感器觀測數據更新率為0.5 Hz，星間測距數據采樣率為0.2 Hz，設定濾波周期為10 s，在同樣的初始條件下分別用EKF、SDREF和MSDREF算法對基于星間相對位置矢量測量的自主導航系統進行仿真，MSDREF的參數γ取0.3。

定義位置誤差和速度誤差分別為:

將本算例中自主導航系統在3種算法下的估計結果進行對比，如表2所示。可以看出EKF的計算時間約是SDREF和MSDREF的7～8倍，這與前述對3種算法計算量的分析吻合。在穩態誤差方面，SDREF由于軌道模型非線性程度不高且不存在線性化誤差，其估計精度與EKF基本相當；而MSDREF根據系統能觀度的變化對濾波協方差陣進行在線調整，改善了EKF和SDREF對能觀度敏感的問題，位置和速度估計的穩態誤差是EKF的21%。

表2 EKF、SDREF和MSDREF估計結果對比

6 結束語

本文針對分布式衛星系統高精度自主導航的任務需求，選取典型的基于星間相對位置觀測的雙星模型，采用SDREF算法對自主導航系統進行估計，通過數學證明了在滿足一定條件時，SDREF的計算量顯著小于EKF的結論。文章提出了一種改進的SDREF算法，根據系統能觀度強弱對標準SDREF算法進行在線修正，并選取一組近距離雙星編隊任務進行數學仿真，對EKF、SDREF和MSDREF的估計結果進行了對比。在相同條件下，MSDREF的估計精度明顯優于SDREF和EKF，且保持了SDREF計算量小的優點。本文研究結果可以為未來分布式衛星工程實施中高精度自主導航算法的設計提供重要參考。