肖永強,王宏力,馮 磊,由四海,何貽洋,許 強

(1.火箭軍工程大學導彈工程學院,陜西 西安 710025;2.青州高新技術(shù)研究所測試控制系,山東 青州 262500)

0 引 言

X射線脈沖星導航是一種新興的自主導航方式,其利用脈沖星輻射的信號對航天器進行定位、授時等服務。其自主性強、可靠性高,有望擺脫航天器對全球定位系統(tǒng)(global positioning system,GPS)等衛(wèi)星導航的依賴,在民用和軍事領域有著巨大的發(fā)展前景[1-10]。但是目前要將其真正應用到實際工程中仍存在很大的差距,其中一個很重要的原因就是脈沖星方位誤差的影響[11-16]。

在脈沖星導航過程中,0.001″的脈沖星方位誤差會造成幾百米的定位誤差[17]。因此,為了提高脈沖星方位誤差估計精度,國內(nèi)學者采用了基于信標衛(wèi)星的估計[18]、魯棒濾波估計[19-21]、組合導航[22-24]等方法對脈沖星方位誤差進行估計。此外,許強在利用信標衛(wèi)星進行估計時考慮了實際存在的衛(wèi)星位置誤差的影響[25]。但是上述研究都沒有考慮航天器的時鐘鐘差影響。雖然孫守明等研究了脈沖星導航以及脈沖星與慣性、多普勒等組合導航中的鐘差修正[26-28]問題,但是在利用信標衛(wèi)星進行脈沖星方位誤差估計時,依然會受到時鐘鐘差的影響。然而從目前公開的文獻來看,利用信標衛(wèi)星進行脈沖星方位誤差估計中的鐘差修正問題尚未有學者研究。當利用信標衛(wèi)星進行脈沖星方位誤差估計時,由于航天器長時間運行,時鐘會發(fā)生漂移,進而會造成系統(tǒng)偏差,從而使脈沖星方位誤差估計精度降低。

為了解決脈沖星方位誤差估計時的鐘差影響,本文通過分析得出時鐘鐘差對脈沖星方位誤差估計的影響為緩變的過程,可看作系統(tǒng)偏差,常用的處理方式是利用增廣狀態(tài)法。但是若利用增廣狀態(tài)法處理,當系統(tǒng)狀態(tài)與系統(tǒng)偏差維數(shù)相近時,會使濾波計算量增加,且易出現(xiàn)數(shù)值不穩(wěn)定的問題[29]。因此,為了解決上述問題,本文提出了一種考慮鐘差修正的兩步卡爾曼濾波(two-stage Kalman filter,TSKF)算法。

1 脈沖星方位誤差估計原理

圖1 脈沖星方位誤差估計原理Fig.1 Principle of pulsar position error estimation

轉(zhuǎn)換過程[1]為

(1)

(2)

式中,tsat為脈沖到達航天器的真實時間;n為真實的脈沖星單位方向矢量。設脈沖星的赤經(jīng)為α,赤緯為β,則滿足:

(3)

設(Δα,Δβ)為脈沖星方位誤差,則真實的脈沖星方位與帶誤差的脈沖星方位滿足

(4)

將式(4)代入式(3),進行泰勒展開并忽略二階及以上小項可得

(5)

式中,帶誤差的脈沖星方向矢量為

(6)

記脈沖星單位方向矢量誤差Δn為

(7)

則可以將式(5)表示為

n-Δn=n′

(8)

(9)

取狀態(tài)變量為X=[ΔαΔβ]T,則傳統(tǒng)脈沖星方位誤差估計算法的狀態(tài)方程和觀測方程為

(10)

(11)

式中,狀態(tài)轉(zhuǎn)移矩陣為

觀測矩陣為

式中，Wk和ηk為系統(tǒng)噪聲和量測噪聲。

2 時鐘鐘差影響分析

由于航天器時鐘頻率和相位的漂移,脈沖到達航天器的真實時間和航天器時鐘測得的脈沖到達航天器的時間之間存在時鐘鐘差,設航天器的時鐘鐘差為δt,則滿足

(12)

時鐘鐘差可通過估計相對于標準時間的鐘差、鐘差漂移率和鐘差漂移率的變化率獲得,因此可將衛(wèi)星時鐘鐘差模型[30]表示為

(13)

式中,x1、x2和x3分別為鐘差、鐘差漂移率和鐘差漂移率變化率;τ為時間間隔。取狀態(tài)變量為x=[x1，x2，x3]T,則鐘差離散過程模型可表示為

(14)

式中,狀態(tài)轉(zhuǎn)移矩陣為

(15)

式中,q1、q2和q3為噪聲的功率譜密度。

仿真取時鐘誤差漂移率為3.637 979×10-12,時鐘漂移率的變化率為6.66×10-18/s;根據(jù)銣原子鐘模型[26],取時鐘的噪聲譜密度分別為q1=1.11×10-22s,q2=2.22×10-32/s和q3=6.66×10-45/s3。給定時鐘初始時刻的鐘差為0,取時間間隔為1 s,則可得到鐘差隨時間的變化如圖2所示。

圖2 鐘差隨時間的變化Fig.2 Change of clock error with time

由圖2可得,雖然鐘差漂移較小,但隨著航天器的長時間運行,鐘差接近5×10-6s,該值與光速相乘理論上會引起1 500 s左右的定位誤差,將會對量測方程式(11)造成較大的影響,進而會對脈沖星方位誤差估計造成嚴重的影響。

使用傳統(tǒng)脈沖星方位誤差估計算法進行仿真,分別驗證其在有無鐘差影響下的估計效果。以脈沖星B0531+21作為觀測脈沖星,其參數(shù)如表1所示。

表1 脈沖星B0531+21參數(shù)Table 1 Parameters of pulsar B0531+21

其中,P為脈沖周期,W為脈沖寬度,Fx為脈沖輻射光子流量,pf為一個脈沖周期內(nèi)脈沖輻射流量與平均輻射流量之比。

脈沖星的觀測噪聲方差[2]為

(16)

式中,A為探測器有效面積,本仿真中設為1 m2;Bx=0.005 ph·cm-2·s-1,為宇宙背景噪聲;d為脈沖寬度W與脈沖周期P之比;tobs為觀測時間,設為1 000 s;則可計算[25]得到σ=(77.69 m)2。脈沖星方位誤差設為(2 mas,2 mas),初始狀態(tài)設為0。

使用同一顆衛(wèi)星,分別在有無鐘差影響的情況下進行仿真,衛(wèi)星軌道參數(shù)如表2所示,仿真結(jié)果如圖3所示。

表2 衛(wèi)星軌道參數(shù)Table 2 Parameters of orbit

圖3 有無鐘差的脈沖星方位誤差估計結(jié)果Fig.3 Estimation results of pulsar position error without correction of clock error

對比分析圖3可得,當無鐘差影響時,傳統(tǒng)脈沖星方位誤差估計算法能較為精確地估計方位誤差。但當存在鐘差時,估計結(jié)果誤差較大,尤其是赤緯誤差,接近80 mas;而且赤經(jīng)和赤緯誤差估計結(jié)果都沒有收斂到一個定值。可見,時鐘鐘差會對脈沖星方位誤差估計精度產(chǎn)生較大的影響,因此非常有必要對航天器時鐘鐘差進行修正。

3 修正鐘差的脈沖星方位誤差估計算法

當考慮鐘差時,將式(1)和式(2)相減可得新的觀測模型為

(17)

由第1.2節(jié)分析可知,鐘差變化是一個緩慢過程,常用的處理方法是增廣狀態(tài)法,即將系統(tǒng)狀態(tài)和鐘差組合為新的狀態(tài)變量進行濾波解算。但是這樣就會使系統(tǒng)狀態(tài)由2維變?yōu)?維,不僅會增加運算負擔,還容易造成數(shù)值不穩(wěn)定的問題[29]。因此,為了解決上述問題,本文采用TSKF算法。

TSKF算法是由Friedland提出用于處理系統(tǒng)常值偏差問題[31],Ignagni 在這基礎上將其應用到處理緩變偏差上[32]。該算法解耦狀態(tài)估計與系統(tǒng)偏差估計,可有效減小濾波過程中的計算量,提高數(shù)值計算穩(wěn)定性。王奕迪等在前期研究中也將其應用到了脈沖星導航中[33-34]。

根據(jù)TSKF算法原理取第一步濾波的狀態(tài)量為X=[Δα,Δβ]T,第二步濾波狀態(tài)量b為鐘差δt,可將方位誤差估計的TSKF算法方程寫為

一步濾波方程:

Xk/k-1=AkXk-1

(18)

(19)

(20)

Xk=Xk/k-1+Kk(Zk-HkXk/k-1)

(21)

Pk=(I-KkHk)Pk/k-1

(22)

二步濾波方程:

Uk=Ak-1Vk-1+Bk-1

(23)

Sk=HkUk+Ck

(24)

(25)

Vk=Uk-KkSk

(26)

(27)

(28)

最終估計結(jié)果為

(29)

式中,Qk為系統(tǒng)噪聲Wk的方差;Pk為狀態(tài)量的協(xié)方差;Rk為觀測噪聲ηk的方差;Bk為δt在式(11)中的驅(qū)動方程;Ck為δt在觀測方程(12)中的驅(qū)動方程。由上述分析可得,Bk=0,Ck=c·Φ。

仿真條件同第2節(jié),估計結(jié)果如圖4所示。

圖4 TSKF算法估計結(jié)果Fig.4 Estimation results of TSKF algorithm

為了進一步驗證所提算法的有效性,選取不同的鐘差對比傳統(tǒng)算法與本文所提算法的估計結(jié)果,如表3所示。

表3 不同鐘差條件下TSKF算法估計偏差Table 3 TSKF estimation bias under different clock errors

由圖4和表3可得,當不考慮鐘差的傳統(tǒng)算法估計偏差較大時,本文所提出的考慮鐘差的TSKF算法估計精度較高,赤經(jīng)和赤緯誤差估計精度均能保持在0.1 mas以內(nèi),且都能收斂到2 mas左右?？梢姳疚乃崴惴苡行Ц綦x時鐘鐘差的影響,使估計精度保持在無鐘差影響下的水平。

為比較TSKF與增廣狀態(tài)濾波算法(augmented state Kalman filter algorithm,ASKF)的精度,仿真條件同第2節(jié),得到仿真結(jié)果如圖5所示。由圖5可得,在10～50天內(nèi),TSKF赤經(jīng)估計結(jié)果振幅略高于ASKF,二者都能很快收斂到2 mas左右,最終估計結(jié)果精度相當,均在0.01 mas以內(nèi)。0～50天,TSKF赤緯估計結(jié)果偏差大于ASKF,但處于可控范圍內(nèi),且TSKF收斂速度略快于ASKF,最終估計結(jié)果TSKF精度略低于ASKF,但偏差在0.06 mas以內(nèi),仍在可控范圍內(nèi)。其他條件不變,對比兩種算法在不同鐘差條件下的估計結(jié)果,如表4所示。表4也驗證了TSKF算法估計結(jié)果接近ASKF,使方位誤差估計精度保持在較高的水平。為比較兩者的計算量,統(tǒng)計上述仿真實驗中兩種算法Matlab程序中的浮點計算次數(shù),得到TSKF算法133 480 860次,ASKF算法316 957 544次,可得TSKF運算量遠小于ASKF,運算量減小了57.89%,有效地減小了計算負擔。

圖5 兩種算法估計結(jié)果對比Fig.5 Estimated results comparison of two algorithms

表4 不同鐘差條件下兩種算法估計偏差Table 4 Estimation bias of two algorithms under different

由上述分析可知,本文提出的TSKF算法估計精度較高,雖略低于ASKF算法估計精度,但處于可控范圍內(nèi),而TSKF算法計算量顯著低于ASKF算法,具有更高的運算效率。

為驗證TSKF算法的魯棒性,分析當仿真用的噪聲方差與數(shù)據(jù)真實方差不一致時的估計結(jié)果。仿真時將真實值的觀測噪聲增大3倍,其他條件不變,即仿真中觀測噪聲方差R=(0.077 69 km)2,而真實值的觀測噪聲為0.233 07 km,估計結(jié)果如圖6所示。

分析圖6可知,當仿真用的噪聲方差與數(shù)據(jù)真實方差不一致時,在0～50天赤經(jīng)估計誤差略大于圖4估計結(jié)果,但估計結(jié)果也能在第80天左右收斂到2 mas左右。在0～100天赤緯估計結(jié)果偏差大于圖4赤緯估計結(jié)果,但在第100天后也能收斂到2 mas左右?？梢?TSKF算法在仿真用的噪聲方差與數(shù)據(jù)真實方差不一致時,仍能使估計結(jié)果收斂到無鐘差影響的狀態(tài)下。

圖6 噪聲不一致時的估計結(jié)果Fig.6 Estimated results with inconsistent noise

由于時鐘鐘差主要由初始鐘差、鐘差漂移率和鐘差漂移率變化率確定,因此為進一步驗證本文算法的魯棒性,在原鐘差模型的基礎上添加周期項:

(30)

式中,A0和B0為幅值;Ts為衛(wèi)星軌道周期。仿真時,鐘差真實值采用上述帶周期項的模型,幅值A0和B0都取為1×10-6s,其他條件不變,得到估計結(jié)果如圖7所示。

圖7 增加周期項的估計結(jié)果Fig.7 Estimated results of increasing period item

其他條件不變,選取不同幅值A0和B0的大小進行仿真實驗,每次仿真實驗進行20次,并選取最后6天的估計結(jié)果取平均值,估計結(jié)果如表5所示。

表5 不同幅值條件下的估計偏差Table 5 Estimation bias under different amplitude conditions

分析上述圖表可得,當鐘差中加入周期變化值而與模型值不一致時,雖然仿真前期偏差較大,但隨著時間的推移,估計結(jié)果仍能收斂到2 mas左右。雖然與模型一致時會在2 mas附近上下抖動,但這種抖動幅度是微小的,可以通過取平均值的方法來提高精度,也驗證了TSKF算法具有較強的魯棒性。

4 結(jié) 論

(1) 利用航天器(信標衛(wèi)星)進行脈沖星方位誤差估計時,航天器的時鐘鐘差會造成系統(tǒng)偏差,對估計結(jié)果產(chǎn)生影響,且這種偏差和影響是不容忽略的。

(2) 考慮鐘差修正的TSKF算法能有效隔離時鐘鐘差對脈沖星方位誤差估計的影響,在保證濾波解算穩(wěn)定性的同時,能使估計精度保持在無鐘差影響下的水平。

(3) 與ASKF算法相比,TSKF算法具有相近的估計精度和更高的運算效率,且TSKF算法具有較強的魯棒性。