擴展卡爾曼濾波確定飛機失速速度的方法

徐 駿 馳

(民航上海航空器適航審定中心 飛行性能室, 上海 200335)

0 引言

失速速度作為飛機的基本性能數據是確定起飛速度、進近速度等關鍵參數的基礎。目前,大部分飛機型號合格審定申請人所采用的方法是:用拖錐提供基準靜壓源來修正機載靜壓傳感器的位置誤差以獲得校正空速。采用該方法可以在大部分速度范圍內提供較精確的靜壓參數,但是在低速段由于出現氣流分離,飛機可能會出現抖振、飛機高度不斷降低等不穩定狀態,導致由拖錐直接測得的數據不理想。為此,部分申請人采用了結合DGPS高度和靜溫對拖錐測得的氣壓高度進行修正,進而得到相對精確的失速速度。該方法對數據的精確性有一定提升,但需要在試驗點前后保持一段時間的穩定飛行狀態,而且不能存在逆溫層。

本文提出了一種不依賴于拖錐數據,采用擴展卡爾曼濾波算法估算靜壓值從而得到失速速度的方法。該方法經過飛行試驗驗證,并與拖錐法的結果進行了對比。試驗表明,在對試驗數據進行濾波處理后,可以成功地獲取較為精確的結果。

1 擴展卡爾曼濾波

卡爾曼濾波的狀態估計算法能較為準確地對目標的運動位置和速度做出預測,但簡單卡爾曼濾波器僅適用于線性的隨機系統,不能直接應用到非線性問題中。通過在當前估計值附近將非線性函數進行一階或二階泰勒展開,舍棄高階分量,從而將非線性模型簡化為線性函數,就可以繼續使用卡爾曼濾波器,這種方法即擴展卡爾曼濾波(Extended Kalman Filter,EKF)[1]。

可以將民用大型運輸類飛機和測量系統視為非線性系統,由此可以通過非線性隨機微分方程建立相應的飛機動態模型和測量系統模型。

(1)

式中,f∈Rn,h∈Rm為非線性函數;x∈Rn為狀態矢量;u∈Rn為輸入/控制矢量;y∈Rm為輸出矢量;z∈Rm為離散點上的樣本點,具有一致的采樣頻率;矢量w和v為隨機的過程噪聲及測量噪聲。盡管認為過程噪聲是隨機變量,但通過對模型進行近似處理可以確定一些動態參數。認為該噪聲的平均值為0,是高斯白噪聲。假定這些噪聲相互獨立,受初始狀態x(0)的影響,有:

(2)

式中,R為測量噪聲;Q為系統隨機噪聲。

對方程組(1)求偏導數,得到雅克比矩陣如下:

(3)

可以將擴展卡爾曼濾波的算法分為兩步。第一步為預測:

(4)

式中,Φ=eAΔt。在上面的方程中輸入測量值完成預測值計算，然后導入濾波器。

第二步為更新:

(5)

如果選用EKF來估算靜壓的話,需要建立狀態和測量模型來構成濾波器。測量值可以由機載傳感器直接獲得。

2 相關數學模型

2.1 動力學模型

本文僅對縱向機動進行分析,將與橫航向相關的狀態作為輸入。縱向動力學模型為:

(6)

式中,u,v和w為體軸系上飛機的速度;p,q和r為飛機轉動矢量ω在體軸系上的投影;θ和φ為俯仰角和滾轉角;ax,ay和az為體軸系上的三向加速度;Ts為外界靜溫;TISA為標準大氣溫度;h為幾何高度;hp為氣壓高度;wpi為隨機輸入噪聲;bi為待估計的偏差。

對標準大氣,其溫度的計算方程為:

TISA=Tsl(1-2.25574212×10-5hp)

(7)

式中,Tsl為海平面標準大氣溫度;hp的單位為m。該公式在11000 m以下有效。

可以用下式將氣壓高度轉換為靜壓:

ps=101325(1-2.25574212×10-5hp)5.255863

(8)

式中,ps的單位為Pa。

可以從Kiel管獲得精確的總壓用于計算校正空速(vc,m/s),公式如下:

(9)

式(6)和式(7)構成了動力學模型的公式。將狀態矢量x和輸入矢量u定義為:

(10)

其他的輸入方程為:

(11)

2.2 測量模型

測量模型是關于測量傳感器的模型,需要考慮傳感器的類型和校正公式。對于失速速度的測量一般需要考慮以下四個方面的因素:(1)由DGPS測量的速度;(2)由慣導測量的俯仰角;(3)由DGPS測量的幾何高度;(4)ADC計算的校正氣壓高度。值得注意的是,由GPS測得的高度和速度應考慮偽距和多普勒漂移的問題。

測量模型的方程為:

(12)

以上方程中直接的測量參數有3個,wmi為噪聲。可以寫成矩陣的形式:

Z=[vm,θm,hm,hpm]T

(13)

在確定了動態模型和測量模型后,需要以適當的形式進行EKF濾波。

2.3 增益矩陣

EKF方程允許對模型參數和增益進行同步估算。如假定方程組(6)中的參數為常數,則將其轉換為狀態方程如下:

(14)

由此,得到增益狀態方程為:

xaug=[u,w,θ,h,hp,bax,baz,bq]T

(15)

3 飛行試驗及結果分析

為驗證該計算方法的可行性,在一架支線飛機上進行了飛行試驗。試驗機安裝了拖錐和Kiel管[2]。在進入失速前飛機保持穩定平飛20 s,然后拉桿減速直至進入自然失速。在失速改出后爬升回原配平高度保持穩定平飛20 s。共試飛了兩個架次,一個架次用于驗證卡爾曼濾波系數,命名為試驗點A;另一個架次采用DGPS綜合法,即基于DGPS數據綜合計算靜壓和空速,命名為試驗點B。試飛的幾何高度變化曲線如圖1和圖2所示(橫坐標S為采樣點)。

圖1 試驗點A幾何高度變化Fig.1 Geometric altitude variation for test point A

圖2 試驗點B幾何高度變化Fig.2 Geometric altitude variation for test point B

盡管數據采集系統可以提供更高的采樣頻率,為與DGPS一致采用了每秒鐘12.5個樣本的采樣率。

EKF的相關參數有測量噪聲矩陣R和系統隨機噪聲矩陣Q。測量噪聲R通過之前的傳感器變化曲線進行計算,Q基于殘差和理論分析并與參考方法的結果進行對比。

通過DGPS綜合法得到的最小校準空速(CAS)和通過EKF濾波計算得到的結果在表1中進行了對比。兩種校準空速的差異很小,絕對值小于0.2 kn或DGPS綜合法數值的0.2%。

表1 失速速度結果對比Table 1 Comparison of stall speed results

通過對靜壓數據進行研究可以進一步分析該方法的可靠性。圖3和圖4為試驗點A和B對應的EKF和ADC測量值相對DGPS綜合法的差異。

圖3 試驗點A壓力高度對比Fig.3 Comparison of pressure altitude for test point A

圖4 試驗點B壓力高度對比Fig.4 Comparison of pressure altitude for test point B

圖中,卡爾曼濾波法顯示出了明顯的回歸趨勢。初始時EKF依據靜壓測量值進行估算,當濾波器有了足夠的數據來提高其預測能力后,其估算值與DGPS綜合法的差值減小到了5 m以內。如圖所示,經歷了2000次采樣或160 s后,差值保持小于3 m。在失速階段ADC靜壓誤差明顯較大,接近10～15 m;而同一階段的EKF估測誤差僅1 m左右。

衡量EKF的一個主要指標是各狀態參數誤差協方差的時間歷程曲線。經過卡爾曼濾波后這些誤差協方差值均呈指數級遞減,變化曲線如圖5所示。

濾波器的更新特性意味著濾波器估算傳感器測量值的能力,其受限于理論上的方差值“HPHT+R”[3],該方差可以視為高斯白噪聲。就本試驗而言,濾波器的性能已經足夠優化,可以得到較理想的結果。

本文提出的方法可以應用到任何飛機上,但是需要謹慎選擇過程和測量噪聲協方差矩陣方程,并對濾波器進行適當調整。除使用拖錐外,還有其他有效地進行靜壓源位置誤差校準的方法,如塔側飛行法或使用DGPS。如果與本文中的方法結合,可以在整個空速校準中避免使用拖錐。

4 結束語

傳統的飛機靜壓源誤差校準試驗方法在確定飛機失速速度時不夠精確,需要對其進行改進。本文提出了一種基于擴展卡爾曼濾波通過計算確定靜壓,進而得到空速的方法。通過在一架支線飛機上進行實際飛行試驗并將該方法的結果與DGPS綜合法的結果進行對比,發現該方法的誤差在0.2 kn或0.2%以內。說明該方法可以作為一種有效的進行失速速度確定的輔助或替代方法。為得到更高精度的結果可以考慮對更新方程中白噪聲作進一步修正。同時,當系統存在較大模型誤差和計算誤差時,協方差矩陣有可能會產生畸變,容易導致結果發散。如果出現發散,可以考慮采用記憶衰減卡爾曼濾波以抑制濾波發散,提高數據的可靠性。

參考文獻：

[1] 張建帆,丁勝,陳軍,等.拓展卡爾曼濾波算法在目標跟蹤中的改進及應用[J].現代計算機,2012,(32):11-14.

[2] Celso Braga de Mendonca.Aircraft stall speed determination with stochastic filtering approach[C]//39th Annual International Symposium.Texas USA:Society of Flight Test Engineers,2008.

[3] Arthur Gelb.Applied optimal estimation[M].Massachusetts: Massachusetts Institute of Technology Press,1974.