基于最大熵卡爾曼濾波的機載SINS星敏感器輔助空中標定技術

王 瑞，張志亮，鄭百東，戴洪德

(海軍航空大學 航空基礎學院，山東 煙臺 264001)

0 引言

由于慣性器件性能參數及系統特性隨著時間變化，慣性導航系統的精度會隨時間降低，在一定時間之后，需要對慣性導航系統的誤差進行標定，以便進行補償從而提高導航精度[1]。慣性導航系統標定方法根據標定時是直接采用慣性器件輸出參數還是導航參數可以分為分立式標定[2]和系統級標定[3]。分立式標定需要將慣性導航系統從飛機等載體上卸載，并通過轉臺精準模擬飛機運動，對慣性導航系統誤差進行標定補償。分立式標定的精度較高，但是需要繁瑣的裝卸操縱、特殊的轉臺設備，且時間效率差[4]，影響正常使用。系統級標定則可以直接在載體上在線完成，通過對比慣性導航系統輸出的導航參數和其他系統較為準確的導航參數，得出差值作為觀測量，通過相應的濾波算法估計慣性導航系統的誤差，對慣性導航系統進行補償與標定。系統級標定具有方便、快捷、效率高的特點，已成為慣性導航系統標定方法的主要研究熱點。

在系統級標定中，根據借助的外部導航參數不同，慣性導航系統能夠標定的參數的效果是不同的。利用里程計可采用位置誤差進行標定[5-6]，利用GPS信號可采用速度誤差進行標定[7]，利用星敏感器可采用姿態誤差進行標定[8]。由于位置是由速度積分而來，且位置信息的可觀測性較差，利用速度誤差和姿態誤差對慣性導航系統進行標定的效果優于位置誤差的標定效果。單獨采用一項外部設備標定慣性導航系統有可能造成其他誤差不可觀測，待標度誤差的估計收斂速度也較慢[9]。在濾波算法中，卡爾曼濾波(KF，Kalman filter)對誤差的估計效果較好，但非高斯噪聲的存在限制了卡爾曼濾波對誤差估計精度的提升。

本文針對機載捷聯慣性導航系統(SINS, strapdown inertial navigation system)空中標定問題，結合星敏感器和GPS差分信號，對比得到姿態偏差和速度偏差，以其作為觀測量，并用最大熵卡爾曼濾波(MCKF, maximum correntropy Kalman filter)方法解決觀測值中存在非高斯噪聲的問題，以提升慣性導航系統誤差標定的速度和精度。

1 星敏感器輔助SINS自標定技術

1.1 總體方案

星敏感器是天文定姿設備，通過拍攝天空恒星圖像，與存儲的星圖對比，確定飛機姿態，其精度較高，且誤差不累積。能夠輸出高精度的姿態信息，且不隨時間漂移，自主性強。將星敏感器輸出的姿態信息與SINS的姿態信息對比，得出姿態偏差。為提高系統整體可觀測性，同時引入GPS差分信號，與SINS輸出的速度信號對比，得出速度偏差。以姿態偏差和速度偏差作為系統級標定模型的觀測值，經由卡爾曼濾波器估計出系統真實誤差項，對SINS進行標定，整體方案如圖1所示。為增加姿態可觀測度，設定飛機進行機動飛行。

圖1 星敏感器輔助慣性導航系統自標定方框圖

1.2 星敏感器工作原理

星敏感器主要由3部分構成，分別為成像器件、圖像傳感器和數據解算電路，涉及到光學、機械、電子等領域。按照其工作原理進行劃分，又可分為成像系統和圖像處理系統。星敏感器工作流程如圖2所示，首先將光學鏡頭朝向天空，對鏡頭范圍內的恒星進行成像處理，成像于圖像敏感器芯片上，將光學信息轉換為電信號，經支持電路將電信號放大為電壓模擬信號，隨后由A/D轉換器轉換為數字信號，此時，星敏感器采集的光學成像成功轉換為可以存儲的數字信號，且存在一一對應的關系。然后，數據處理系統掃描上面存儲的數字星圖，找出其中的恒星，并計算出觀測到恒星的亮度信息和在像平面上的位置信息。最后，經星圖識別確定星圖內恒星與導航星庫中恒星的對應關系，通過姿態解算方法確定星敏感器的三軸姿態，進而為載體提供高精度姿態信息，實現其導航的目的。

圖2 星敏感器工作流程框圖

其中，星敏感器成像原理如圖3所示，O-xsyszs代表星敏感器坐標系，O-xvyvzv代表CCD面陣坐標系，光軸OOs與zv軸、zs軸均重合，f代表光學透鏡的焦距。p點為恒星于CCD面陣上的點，其坐標為(x,y)，此時，星光矢量在O-xsyszs中的方向向量為：

(1)

式中，Vs代表觀測噪聲。

圖3 星敏感器成像原理圖

假設星敏感器鏡頭內能夠觀測到n顆恒星，根據上式可得星光矢量在O-xsyszs中的方向向量分別為(xs1,ys1,zs1)(xs2,ys2,zs2)(xs3,ys3,zs3)…(xsn,ysn,zsn)，搜索導航星庫后可知，這n顆恒星在慣性空間中方位分別為(xi1,yi1,zi1)(xi2,yi2,zi2)(xi3,yi3,zi3)…(xin,yin,zin)，則星敏感器姿態轉移矩陣滿足下式：

(2)

也可寫為:

S=GA

(3)

式(2)和式(3)相互對應，其中A為載體的姿態轉移矩陣。

當n=3時，

A=G-1S

(4)

當n>3時，根據最小二乘法原理得：

A=(GTG)-1(GTS)

(5)

由上述分析可知，在慣性空間內，恒星相對于地球的方位始終保持不變，星敏感器正是以此來計算載體的高精度姿態信息。星敏感器的測量精度會受到地球的章動、自身誤差等影響，但總體影響較小，遠小于1″，也不會隨時間不斷積累。

綜上，星敏感器可以為本文自標定系統提供姿態信息，從而提高機載SINS自標定系統的精度。

1.3 陀螺與加速度計誤差模型

機載SINS自標定為現場在線標定，沒有拆卸操作，不影響安裝誤差。在代表機體坐標系b和慣性坐標系i下，陀螺儀誤差模型[10]可以用式(6)表示:

(6)

同樣在機體坐標系下，加速度計誤差模型如式(7)所示:

δfb=δKafb+▽a

(7)

1.4 SINS誤差方程

取n為導航坐標系，機載接連慣性導航系統的姿態誤差方程如式(8)所示，速度誤差方程如式(9)所示，位置誤差方程如式(10)～(12)所示:

(8)

(9)

(10)

(11)

(12)

1.5 狀態方程

根據SINS誤差方程，可以構建系統自標定誤差模型，但由于不進行慣導及慣導內部器件的重新拆裝，即可認為安裝誤差不發生改變。所以主要對陀螺漂移和加速度計零偏，及相應的標度因數誤差進行標定。捷聯慣性導航系統在標定時，需對姿態角誤差φE、φN、φU，速度誤差δVE、δVN，位置誤差δL、δλ、δh，陀螺儀漂移、陀螺儀標度因數誤差、加速度計零偏、加速度計標度因數誤差等進行估計與補償，將以上誤差項作為自標定模型的狀態變量，構建系統狀態方程如式(13)所示:

(13)

式中，狀態變量為X(t)。

X(t)=[φEφNφUδVEδVNδVUδLδλδhεbx

εbyεbz▽ax▽ay▽azδKaxδKayδKazδKgxδKgy

δKgz]T，系統隨機噪聲W(t)=[wgxwgywgzwaxwaywaz]T，由陀螺儀隨機噪聲wgxwgywgz和加速度計隨機噪聲waxwaywaz兩部分組成。F(t)和G(t)分別為系統矩陣和噪聲輸入矩陣:

(14)

式中，FN為9×9維慣導系統誤差矩陣，012×21為元素為0的12×21階矩陣。FS和FM分別如下所示：

其中:元素矩陣F1和F2是陀螺儀和加速度計輸出在導航坐標系的投影。

式中，C11、C12、C13等是載體系到導航系轉移矩陣的元素。

1.6 量測方程

采用GPS輔助SINS空中自標定方法僅以速度誤差作為量測量，標定天向陀螺時漂移估計較慢，效果較差。星敏感器能夠求出高精度姿態信息，但對加速度計誤差不可觀測。所以，在此處將星敏感器與SINS的姿態偏差和差分GPS與SINS的速度偏差，作為自標定模型的觀測值量，構成組合量測模型，其量測方程如式(15)所示:

Z(t)=HX(t)+V(t)

(15)

式中，觀測量Zav(t)=[φEφNφUδVEδVN]T，Hav=[I5×505×16]5×21，V(t)為量測噪聲，其中I5×5為5×5階單位陣。可以看出，量測方程采用的是“姿態+速度”匹配量測模型，與位置匹配兩側模型相比，可以提升自標定的速度。

2 基于最大熵的卡爾曼濾波

2.1 最大熵準則

相關熵是一種兩個隨機變量間的廣義相似度量，對于聯合分布函數為FXY(x,y)的兩個隨機變量X和Y是實域向量，其相關熵定義如下[12]：

(16)

式中，E[κ(X,Y)]表示κ(X,Y)的期望，κ(.,.)表示Mercer核函數，此處取高斯核函數，即:

(17)

式中，e=x-y且σ>0表示核函數的帶寬[13]。

一般情況下，由于已知數據有限且FXY(x,y)未知，可由樣本平均估計值計算相關熵：

(18)

由估算值可知，V(X,Y)為有界正數，當X=Y時，V(X,Y)取得最大值。對高斯κ(.,.)泰勒級數展開[14]：

(19)

由式(19)可以看出，變量X和Y的相關熵是對X-Y的偶次冪在取均值后，再加權求和。相關熵中有高階矩信息，若噪聲為非高斯分布，最大相關熵的魯棒性可顯著改善。

2.2 最大熵卡爾曼濾波

基于最小均方差準則的經典卡爾曼濾波器只能精確到誤差的二階統計特性，但是如采用最大相關熵準則作為卡爾曼濾波器的濾波準則，則可以應用誤差的高階統計特性，從而改善非高斯噪聲條件下的濾波器魯棒性。

SINS標定模型離散化后，式(18)和式(20)聯立改寫為：

(15)

(21)

式中，I為單位陣，

(22)

Dk=WkXk+Ek

(23)

由于Ek為白噪聲，取損失函數JL為：

(24)

式中，di,k是Dk的第i行元素，L=n+m是Dk的維數。

由最大相關熵準則可知[16]，Xk的最優估計為使損失函數JL達到最大時的值：

(25)

式中，ei,k=di,k-wi,kXk。

損失函數的最大值可通過求偏導數為0的Xk：

(26)

化簡得：

(27)

代入ei,k得：

(28)

以矩陣形式表示：

(29)

其中:

(30)

Cx,k= diag(Gσe1,k,…,Gσen,k)

(31)

Cy,k= diag} (Gσen + 1,k,…,Gσen + m,k)

(32)

將式(29)化簡后得[14]：

(33)

式中，

(34)

由此整理可得最大熵卡爾曼濾波器的算法流程：

4)更新增益矩陣;

5)狀態更新;

6)協方差陣更新;

由此，最大相關熵卡爾曼濾波通過高階矩信息，提高了非高斯噪聲條件下的魯棒性，可以解決機載捷聯慣導系統在空中標定時難以解決的非高斯噪聲問題。

3 仿真實驗與分析

3.1 仿真條件

圖4 SINS解算原理圖

SINS輸入的加速度計和陀螺儀數據與飛機運動模型相一致，并疊加各誤差項，為自標定算法提供數據支撐。

3.2 仿真方法

不同于傳統標定方法，通過采用轉臺轉動激勵出SINS相關誤差參數，實現慣導誤差參數的標定，本文基于常規起飛過程的機載SINS自標定方法，僅利用飛機常規起飛過程中的機動，激勵出慣性導航系統內部慣性器件的誤差，依據系統級標定原理，以機場跑道附近的高精度差分GPS和慣性導航輸出速度差值作為觀測量，應用Kalman濾波對慣導系統的狀態量進行估計，實現機載SINS的誤差標定。

通過飛行軌跡的設計，能夠使得系統更好地激勵出系統狀態參數，從而使機載SINS元器件誤差能被更好的標定出來。

由于在飛行過程中GPS存在精度要求無法長時間保證的問題，則充分利用起飛階段飛機的機動以及機場跑道附件能夠布設高精度GPS差分站的優勢，將空中自標定的軌跡設置在機場附近。結合常規起飛階段的仿真分析，設計恰當的飛行軌跡，且保證充分利用機場附近高精度的差分GPS。根據某型飛機的起飛機動特性以及機場的數據，在盡可能不增加飛行員額外操作的原則下設計起飛標定動作共5個階段，將起飛階段和盤旋階段相結合，共588 s。

飛行軌跡如圖5所示，假設飛機位于5 km的高度，初始速度為200 m/s，常規的機動動作主要分為5個階段：

1)勻速階段：飛機保持初始狀態飛行。

2)加速階段：飛機以2.5 m/s2加速度將飛機速度從200 m/s提升到275 m/s。

3)爬升階段：飛機在機場上方5 km高度處，以275 m/s速度進行爬升，將飛機以2/s的俯仰角拉升至30°；隨后保持角度爬升一段時間；最后將飛機改平，爬升階段完成。

4)滾轉階段：飛機在機場上方5 km高度處，以275 m/s速度進行滾轉，首先以5/s左滾轉到45，隨后保持該狀態飛行一段時間，最后以5/s右滾轉將飛機調整到水平，保持原狀態繼續飛行。

5)轉彎階段：飛機在機場上方5 km高度處，以275 m/s速度進行轉彎，首先以9/s右轉彎到90，隨后保持該狀態飛行一段時間，最后以9/s左轉彎將飛機調整到原方向，保持原狀態繼續飛行。

通過飛機機動激勵SINS輸出參數信息，利用最大相關熵卡爾曼濾波器對各項誤差參數進行標定，驗證設計的方法有效性。

3.3 仿真參數

機載捷聯慣導系統的仿真參數如表1所示。

表1 仿真參數設定

星敏感器的姿態角誤差φE=10″，φN=10″，φU=20′，GPS位置誤差1 m，速度誤差0.2 m/s。

基于以上誤差參數，設定協方差初始陣為：

P=diag([10″；10″；20′；0.2 m/s；0.2 m/s；0.2 m/s；1 m；1 m；1 m；0.03°/h；0.03°/h；0.03°/h；100 μg；100 μg；100 μg；400×10-6；400×10-6；400×10-6；400×10-6；400×10-6；400×10-6]×10)2

噪聲方差初始矩陣為：

由于采用了星敏感器和GPS信號作為空中標定的輔助信號，所以在觀測量中引入姿態和航向三軸角度信號和水平兩個速度信號。在非高斯噪聲條件下，設定核函數的帶寬σ=0.8，量測噪聲為:

仿真時，飛機航跡由爬升、轉彎、直線飛行等段組成，如圖5所示。飛機的初始狀態為：高度5 km，速度200 m/s，保持勻速直線飛行，以黑色原點代表初始位置。之后，進入分別進行爬升、轉彎、滾轉等機動動作，仿真時間共588 s。分別采用最大相關熵卡爾曼濾波和經典卡爾曼濾波依據GPS信號和星敏感器信號，對機載捷聯慣導系統進行標定。仿真結果如圖6和圖7所示，在圖中，虛線表示采用經典卡爾曼濾波后慣性器件的誤差估計曲線，實線表示采用最大熵卡爾曼濾波后慣性器件的誤差估計曲線，直實線代表各誤差項的設定值。

圖5 飛機仿真航跡

圖6 陀螺儀漂移ε與加速度計零偏▽估計曲線

圖7 加速度計標度因數誤差δKa與陀螺儀標度因數δKg誤差估計曲線

3.4 仿真結果與分析

由圖6可以看出，采用了最大熵卡爾曼濾波的星敏感器輔助SINS自標定，在三軸陀螺漂移的估計中均優于經典卡爾曼濾波，MCKF的收斂速度和估計精度均高于經典卡爾曼濾波。尤其是在z軸陀螺儀漂移的估計上，MCKF的精度顯著高于經典卡爾曼濾波，并且克服了無星敏感器輔助時，z軸陀螺儀漂移估計發散的問題[17]。但在加速度計零偏的估計上，MCKF的優勢較小，尤其是在z軸上的估計，難以克服慣導系統高度通道發散的問題。但是，對于SINS垂直通道，通常方法是借助于大氣數據參數來修正高度，因此，在后續研究中，應借助大氣數據系統的高度參數，提高對垂直通道的觀測度和實現標定效果的提升。同時對比文獻[17]可以看出，星敏感器對加速度計零偏估計的影響較小，這是由于星敏感器主要用于測量姿態數據，其對加速度的檢測影響較小。

由圖7可以看出，MCKF算法對加速度計標度因數誤差的估計精度比較高，但是收斂速度與KF相當；對陀螺儀標度因數誤差的估計精度和收斂速度優于KF。同時，和不采用星敏感器的自標定技術相比[17]，估計精度略有提升。標度因數誤差反應慣性測量器件的一致性，反應的是器件本身的測量特性，非高斯噪聲對標定因數誤差影響較小。因此，可以看出，不管是否采用最大熵原則，對標定因數的標定的影響較小。

不同濾波技術下，陀螺儀漂移、加速度計零偏以及相應的標度因數誤差的標定結果如表2所示。

表2 非高斯條件下不同濾波技術的誤差標定結果

由表2可知，在非高斯噪聲條件下，MCKF的標定結果明顯優于經典卡爾曼濾波。但是對于MCKF條件下，星敏感器的采用對標度因數的誤差標定提升不夠明顯。

4 結束語

對于機載SINS的空中自標定，為提升陀螺儀漂移和加速度計零偏的標定精度，結合前期的SINS誤差模型的可觀性分析，分別采用了星敏感器和GPS信號輔助機載SINS的空中自標定。由于空中飛行時，飛機機動、機體振動等因素引起的附加噪聲是非高斯噪聲，在此環境下，為消除非高斯噪聲對自標定精度的影響，采用了最大相關熵卡爾曼濾波方法。在最大相關熵原則中，引入了高階矩信息，可以提升濾波器的魯棒性，從而顯著提升了機載SINS誤差的自標定精度和速度。為驗證所提自標定方法的有效性，以SINS模型和飛機運動模型為基礎，通過設計飛機機動航跡，激勵SINS各項誤差，采用最大相關熵卡爾曼濾波算法來估計誤差項。同時和未采用星敏感器輔助的SINS空中標定得到的誤差項進行了對比。通過對比分析可以看出，在星敏感器輔助下，最大相關熵卡爾曼濾波算法對陀螺儀漂移和加速度零偏的標定精度有了顯著提升，尤其是天向陀螺儀的漂移，標定速度和精度都得到了提升。但是高度通道存在發散情況，因此需借助其他方法提升垂直加速度計零偏的標定精度。同時，本文在建立誤差模型時，僅對慣性元器件的常值誤差、刻度因數誤差和安裝誤差進行分析討論，未對慣性元器件溫度、高階項誤差等進行深入分析，在后續工作中可建立更全面的誤差模型，進一步提高誤差項的標定精度。