基于UKF的FADS／INS融合大氣數(shù)據(jù)估計

蔣保睿，劉 鵬，肖地波

（成都信息工程大學控制工程學院，成都 610225）

1 引 言

對于高速飛行器而言，大氣數(shù)據(jù)包括空速（V）、攻角（α）、側(cè)滑角（β）、總壓（q）、靜壓（p）等參數(shù)，其測量精度和可靠性影響飛行器控制系統(tǒng)的正常工作及性能發(fā)揮，關(guān)系到飛行安全和質(zhì)量。目前普遍使用探針傳感技術(shù)測量大氣數(shù)據(jù)，凸出的空速管和風向標不僅在高速飛行時測量數(shù)據(jù)不夠準確，還會影響飛行器的機動和隱身性能。 為了解決上述問題，人們開始尋找更為可靠的測量方法。嵌入式大氣數(shù)據(jù)傳感系統(tǒng)（Flush Air Data Sensing System，F(xiàn)ADS）采用嵌入在飛行器前端不同位置上的壓力傳感器陣列測量飛行器表面的壓力分布，根據(jù)測得的壓力值，按照一定的算法解算得到大氣數(shù)據(jù)。美國Dryden 飛行研究中心在20 世紀60年代開始對FADS 進行研究，在理論上證明了該系統(tǒng)可滿足現(xiàn)代飛行器對飛行速度、機動性等方面的要求。

慣性導航系統(tǒng)（Inertial Navigation System，INS）短時間的導航精確度高于FADS，但其突出缺點是受到漂移誤差影響，導航的精度隨時間的增長而降低。 為解決這些問題，將大氣數(shù)據(jù)系統(tǒng)、慣導系統(tǒng)組合使用是較為常見的一種方式。 上世紀末，F(xiàn)ADS 被首次應(yīng)用在飛行器上，標志著系統(tǒng)的可行性；近年來，陸辰等人提出了一種基于FADS 的擴展卡爾曼濾波（Extended Kalman Filtering，EKF）實時估計大氣參數(shù)的方法，驗證了大氣數(shù)據(jù)系統(tǒng)的測量范圍和可靠性；Karlgaard 等人在研究火星探測器著陸階段的空氣數(shù)據(jù)估計時，使用全狀態(tài)（速度和姿態(tài)）慣性量的輔助，與原有的大氣數(shù)據(jù)處理算法相比，精度得到了提升。

為了提升算法在復雜飛行狀態(tài)時對參數(shù)的估計能力、進一步降低誤差，本文提出一種基于無跡卡爾曼濾波（Unscented Kalman Filtering，UKF）的大氣數(shù)據(jù)和全狀態(tài)慣導數(shù)據(jù)估算真空速、攻角和側(cè)滑角的方法。 該方法利用INS 測量慣性數(shù)據(jù)，同時FADS 測量大氣數(shù)據(jù)，最后通過UKF 濾波，解算出精度更高的攻角、側(cè)滑角、真空速等大氣參數(shù)。 另外，本文FADS 解算部分提出了一個快速求解動靜壓的方法，在減少迭代步驟、提高計算速度上取得了良好的效果。

2 系統(tǒng)設(shè)計

基于UKF 的FADS／INS 融合大氣數(shù)據(jù)估計系統(tǒng)包括慣性導航系統(tǒng)、用于測量氣壓值的FADS 系統(tǒng)和信息融合估計算法。 本章介紹慣導系統(tǒng)和FADS 的數(shù)據(jù)解算、UKF 濾波算法以及根據(jù)濾波輸出值求解大氣數(shù)據(jù)的方法。

2．1 慣導系統(tǒng)數(shù)據(jù)解算

其中，旋轉(zhuǎn)矩陣的值完全由當前時刻的姿態(tài)角決定，求解公式為

2．2 FADS 模型解算

FADS 通過飛行器表面特定區(qū)域的壓力分布反推得到飛行參數(shù)，解算順序一般為：動靜壓、攻角側(cè)滑角、馬赫數(shù)和氣壓高度。

已知測壓孔的位置參數(shù)和每個孔的壓力值，可列寫如下方程求解動靜壓：

式中：σ——第i 個孔對應(yīng)的氣流方向，p——第i 個測壓孔在當前時刻所測壓力，f——一個單調(diào)增函數(shù)，其自變量為動靜壓之比，因變量為修正系數(shù)，用于修正不同馬赫數(shù)下測壓孔數(shù)據(jù)，可由實驗或牛頓理論確定。 當孔的數(shù)量為5 時，式中n ＝5，即用五個方程組成的方程組求解動壓q與靜壓P兩個未知數(shù)。 如果直接利用計算機求數(shù)值解，不僅耗時長、計算繁瑣，而且精度不高。 利用偽逆，公式（3）可化簡為

經(jīng)迭代即可求解q與P，經(jīng)測試，在結(jié)果精度相同的情況下，改進算法的運行時長比原算法減少約70%。

分別選擇豎直對稱軸和水平對稱軸的三個測壓孔數(shù)據(jù)，得到局部攻角與局部側(cè)滑角，具體表達式為

其中，A、B 用于求解攻角，利用垂直的測壓孔數(shù)據(jù)，與豎直分布的測壓孔數(shù)據(jù)有關(guān)；同理A′、B′、C′用于求解側(cè)滑角，與水平分布的測壓孔數(shù)據(jù)有關(guān)。

2．3 UKF 方程設(shè)計

在適用于非線性系統(tǒng)的改進卡爾曼濾波算法中，常用方法有EKF 和UKF。 相較于EKF，UKF 采用無跡變換（Unscented Transformation，UT），在狀態(tài)預測時采用非線性傳播方式，避免了EKF 在使用線性逼近非線性變換時的截斷誤差，濾波結(jié)果更加準確，且避免了繁瑣的雅克比求解工作。

選取飛機的真空速、加速度、姿態(tài)角、姿態(tài)角角速度作為狀態(tài)量，即

因為大氣紊流的時間短、速度及其改變量小，

所以可以認為地速與真空速的加速度一致，即

經(jīng)整理可得狀態(tài)方程

式中：u ＝（A，A，A，ω，ω，ω）——輸入量；w——系統(tǒng)噪聲向量。

由公式（4） ～公式（8）可得真空速V、攻角α、側(cè)滑角β，基于此三個大氣數(shù)據(jù)，按公式（6）求得地面坐標系下的空速分量，即觀測方程為

2．4 濾波計算

其中，a 的選擇控制了采樣點的分布，通常取一個不大于1 的正數(shù)，此處取0.2；b 是一個非負的權(quán)系數(shù)，對于高斯分布，可以從理論上證明b ＝2 是最優(yōu)的；n 為維數(shù)，根據(jù)系統(tǒng)狀態(tài)量X 的形式，此處取n ＝12；λ ＝n（a－1）。

大氣數(shù)據(jù)中的真空速V、攻角α、側(cè)滑角β 均可由濾波的輸出值經(jīng)公式（6）計算得到。 馬赫數(shù)是真空速與當?shù)芈曀俚谋戎?，?/p>

3 仿真實驗

本文設(shè)計FADS／INS 的UKF 算法目的在于提高飛行器大氣數(shù)據(jù)的估計精度，開展仿真實驗驗證算法性能，將濾波前測量結(jié)果的誤差、EKF 算法誤差和本文提出的UKF 濾波的誤差進行對比。 模型由數(shù)據(jù)生成、單傳感器信息預處理、EKF／UKF 傳感器融合、飛行參數(shù)解算等部分組成。 整體框圖如圖1 所示，在對比仿真中，除濾波算法，模型的參數(shù)設(shè)置均相同。

圖1 仿真模型整體框圖Fig．1 Overall block diagram of simulation model

飛行仿真條件如下：飛行高度起始值為3 km，上升段上升至20 km，巡航300 s 后下降至3 km；飛行馬赫數(shù)上升段段從0.1 提升至7，下降段再由7下降為0.1；初始位置記為地面坐標系的（0，0，0）點，地面坐標系的X 軸朝向正東，Y 軸指向正北，由右手坐標系確定Z 軸垂直向下。 初始航向角0°，初始時刻飛行器坐標系與地面坐標系朝向相同，飛行器坐標系0 點位于地面坐標系的（0，0， －1000 m）處。 飛行總時長600 s，采樣間隔1 Hz。

對慣性元件的誤差建模，使用一階高斯——馬爾科夫過程，并且設(shè)定相關(guān)時間系數(shù)和標準差兩個參數(shù)。 相關(guān)時間系數(shù)越小，一階高斯—馬爾科夫過程越趨近于白噪聲表現(xiàn)（變化迅速）。 反之，系數(shù)越大，越趨近于隨機游走表現(xiàn)（變化緩慢）。 標準差在一定程度上由測量設(shè)備的精度決定，仿真測試時設(shè)置加速度計的相關(guān)時間系數(shù)為0. 9，標準差為0.02 m／s，陀螺儀對角加速度的測量誤差包括刻度因數(shù)誤差、零偏、隨機噪聲與耦合誤差，刻度因數(shù)誤差和零偏誤差可在裝置內(nèi)部進行一定精度的補償，耦合誤差因為數(shù)值很小通常忽略不計。 各孔壓力數(shù)值由馬赫數(shù)、攻角、側(cè)滑角等參數(shù)的實際值求解并添加噪聲而得，噪聲服從一階高斯——馬爾科夫過程，相關(guān)時間系數(shù)0.5，標準差約為20 Pa。

3．1 真空速各軸分量估計結(jié)果

運行程序，EKF 與UKF 直接估計飛行器的三軸速度分量結(jié)果如圖2 所示。

圖2 濾波前后三軸速度分量對比Fig．2 Comparison of three?axis velocity components before and after filtering

兩種濾波方法對速度均有很好的量測性和收斂性，隨著迭代次數(shù)的增加，UKF 的實時性一直可以得到保證，但是EKF 的輸出結(jié)果抖動較大，是在線性化過程中泰勒展開導致強非線性部分數(shù)據(jù)丟失的體現(xiàn)。

對比濾波前（僅INS）、EKF 和UKF 輸出速度的大小，在600 s 的采樣時間內(nèi)，濾波前各軸速度平均誤差之和為10.18 m／s，EKF 各軸速度平均誤差之和為0. 98 m／s，UKF 各軸速度平均誤差之和為0.38 m／s，不同方法濾波后各軸速度分量的誤差隨時間變化曲線如圖3 所示。

圖3 濾波后速度誤差分量對比Fig．3 Comparison of velocity error components after filtering

因為飛行速度基數(shù)大，所以分析X 軸的速度——時間曲線，發(fā)現(xiàn)EKF 和UKF 的估計誤差均較大，但兩者相比，UKF 的誤差比EKF 的小。 在Y和Z 軸分量上，UKF 的濾波效果明顯好于EKF。

3．2 攻角、側(cè)滑角估計誤差

濾波后輸出的各軸速度值由公式（6）計算可得攻角α 和側(cè)滑角β 的估計值。 圖4 和圖5 給出了采用EKF 與UKF 得到的攻角和側(cè)滑角隨時間變化曲線，在整個仿真階段，由UKF 估計的攻角和側(cè)滑角誤差小，且振蕩幅度低。

圖4 不同方案估計攻角的值與誤差Fig．4 Estimated angle of attack values and errors of different schemes

圖5 不同方案估計側(cè)滑角的值與誤差Fig．5 Values and errors of side slip angle estimated by different schemes

α 和β 誤差統(tǒng)計特性見表1。

表1 不同方案估計α、β 的誤差結(jié)果對比（單次）Tab．1 Comparison of error of different schemes for estimating α and β（single time）

設(shè)定每次飛行時間500 s，運行100 次飛行仿真程序，不同方案估計α、β 的誤差結(jié)果對比見表2。

表2 不同方案估計α、β 的誤差結(jié)果對比（100 次）Tab．2 Comparison of error of different schemes for estimating α and β（100 times）

從表1 和表2 數(shù)據(jù)可以得到結(jié)論：INS 直接估計的各項指標都是最差的，從數(shù)值上判斷其估計效果最差，UKF 估計的效果最好，EKF 介于兩者之間。

4 結(jié)束語

本文針對超音速飛行器在飛行時估計飛行數(shù)據(jù)時遇到的精度差、可靠性低等問題，提出了一種基于多傳感器融合的大氣數(shù)據(jù)估計算法。 經(jīng)數(shù)據(jù)預處理、信息融合、飛行數(shù)據(jù)計算等步驟，輸出馬赫數(shù)、攻角和側(cè)滑角。 仿真結(jié)果表明，采用UKF 的FADS／INS 融合估計方法與單傳感器估計和EKF 估計方法相比，具有較高的精度和良好的穩(wěn)定性，可為精確估計飛行器的大氣參數(shù)提供支持，避免因參數(shù)估計錯誤而引起的事故。

但是本UKF 估計方法也具有一定的缺陷：沒有分析非球形或楔形飛行器機頭FADS 的模型結(jié)算方法，針對特殊形狀的機頭應(yīng)單獨進行空氣動力學建模；其次，系統(tǒng)對風速建模時選擇的是定常風，若大氣紊流和突風較大，可能會影響算法的精確度，若能對模型進行改進，將更加有效地適用于各種飛行狀態(tài)下的大氣數(shù)據(jù)估計。