基于斜率約束插值的燃油流量特征提取方法

孫文柱，曲建嶺，周勝明，邸亞洲，高 峰（海軍航空工程學院 青島校區(qū)，山東 青島266041）

0 引 言

準確的飛機狀態(tài)監(jiān)控可以發(fā)現(xiàn)飛機的狀態(tài)變化趨勢，從而有效保障飛行安全?，F(xiàn)有飛機廣泛安裝有飛參記錄系統(tǒng) （flight data recorder，F(xiàn)DR），利用FDR 記錄的飛參數(shù)據(jù)監(jiān)控飛機狀態(tài)能夠?qū)崿F(xiàn)以較低的成本長時間較完整地監(jiān)測每架飛機，是一種可行性較高的飛機狀態(tài)監(jiān)控手段。

燃油流量是描述飛機發(fā)動機運行狀態(tài)的主要參數(shù)之一，通過監(jiān)測同等狀態(tài)下的燃油流量的變化趨勢可以監(jiān)控飛機發(fā)動機性能退化程度。然而目前大多數(shù)機型雖然安裝有燃油流量計，卻不實時記錄燃油流量值。因此利用飛參數(shù)據(jù)中記錄的剩余燃油量獲取燃油流量特征是一種有效的手段。目前已有許多學者對基于飛參數(shù)據(jù)中發(fā)動機狀態(tài)參數(shù)和環(huán)境參數(shù)的燃油流量提取方法進行了研究，例如：文獻 ［1］采用多元線性回歸方法對波音737－700飛機巡航段的燃油流量進行了回歸建模，適用于發(fā)動機巡航穩(wěn)定階段燃油流量的計算，但該方法無法完成整個飛行架次的燃油流量提取；文獻 ［2］以發(fā)動機運行參數(shù)為輸入，采用BP 神經(jīng)網(wǎng)絡(luò)燃油流量模型辨識的方法預測飛行全過程的油耗，但使用該方法獲得的燃油流量與發(fā)動機運行參數(shù)存在相關(guān)性，不適用于飛機性能的綜合監(jiān)控。本文提出一種飛參數(shù)據(jù)中的剩余燃油量中提取燃油流量的方法，該方法通過對飛參數(shù)據(jù)中的離散剩余燃油量進行數(shù)據(jù)清洗及插值得到連續(xù)的剩余燃油量曲線，而后求取剩余燃油量曲線對時間的變化率提取出燃油流量。

飛參數(shù)據(jù)中可用于飛機狀態(tài)監(jiān)控的數(shù)據(jù)可以定義為3種類型：累積量、瞬時狀態(tài)量和瞬時變化量。瞬時狀態(tài)量指在不受外力情況下保持不變的量，如速度、角速度等；瞬時狀態(tài)量的一次積分為累積量，如高度、距離、角度等；瞬時狀態(tài)量的一階導數(shù)為瞬時變化量，如各向加速度。累積量一般不可直接用于監(jiān)控飛機狀態(tài)；瞬時狀態(tài)量常用于多參數(shù)組合監(jiān)控，如航空發(fā)動機轉(zhuǎn)子固定轉(zhuǎn)速下的排氣溫度和滑油壓力；瞬時變化量可以用于參數(shù)的極值監(jiān)控，如飛機接地時刻過載監(jiān)控。在基于飛參數(shù)據(jù)的飛機狀態(tài)監(jiān)控過程中，為了充分挖掘飛參數(shù)據(jù)參數(shù)中的有用信息，發(fā)揮每個飛參數(shù)據(jù)參數(shù)的作用，常需要進行3類參數(shù)間的轉(zhuǎn)化。將瞬時狀態(tài)量轉(zhuǎn)化為累積量或?qū)⑺矔r變化量轉(zhuǎn)化為瞬時狀態(tài)量的問題叫做一次積分量特征提取問題。將累積量轉(zhuǎn)化為瞬時狀態(tài)量或?qū)⑺矔r狀態(tài)量轉(zhuǎn)化為瞬時變化量的問題叫做一階導量特征提取問題。本文致力于解決一階導量特征提取問題，并對一種典型的該類問題——從剩余燃油量飛參數(shù)據(jù)中提取燃油流量進行了詳細描述。

1 剩余燃油量的數(shù)據(jù)清洗

本文以飛參記錄系統(tǒng)記錄的剩余燃油量為研究對象，此剩余燃油量參數(shù)的采樣頻率為2Hz，量程0～12000 升，每個采樣值以9個二進制位存儲，采樣分辨率約為25升。某架次完整飛行過程的剩余燃油量記錄值如圖1所示，圖1中可以看出，由于采樣分辨率較大，剩余燃油量數(shù)據(jù)基本呈階梯狀單調(diào)遞減的趨勢。但在個別位置會出現(xiàn)偏離真實值較大的野值點 （6500s左右）和二值間反復變化的波動點（775s左右）。野值點的產(chǎn)生原因一般為線路間接觸不良或較大的電磁干擾，而波動點是由剩余燃油量的原理性誤差產(chǎn)生的。該型飛機剩余燃油量G 的計算公式為

式中：VI——燃油 初 始 值；VC——燃油耗 量；ρ——燃 油密度。

圖1 飛參數(shù)據(jù)中的剩余燃油量數(shù)據(jù)

VI需要在每飛行架次前人工裝定，因此對于一個架次來說，VI是恒定的。VC由葉輪式耗量傳感器輸出脈沖累積產(chǎn)生，因此VC一定是單調(diào)增加的。ρ為大氣靜溫T 的函數(shù)，即ρ＝Θ（T），而T 由飛行過程中的全壓和靜壓解算產(chǎn)生，因此ρ是實時變化的。實時變化的ρ會導致計算所得的剩余燃油量數(shù)值偶然增大，經(jīng)采樣后就形成了剩余燃油量中的波動點。

真實剩余燃油量一定是隨飛行時間的增加而減少的，即剩余燃油量數(shù)據(jù)應(yīng)具有單調(diào)遞減性。但野值點和波動點破壞了剩余燃油量的單調(diào)遞減性，因此需要對剩余燃油量數(shù)據(jù)進行清洗以去除野值點和波動點，增加剩余燃油量數(shù)據(jù)的精度［3］。對于野值點的處理采用文獻 ［4］所述參數(shù)變化率失真法進行野值點判別，判別后利用前值替補后值法對野值點進行修正。波動點的去除采用中值濾波算法，具體內(nèi)容參見文獻 ［5］。由中值濾波的原理可知，要保證中值濾波后的數(shù)據(jù)單調(diào)，只需要保證窗口滑動過程中進入窗口的數(shù)據(jù)值小于等于流出窗口的數(shù)據(jù)值。剩余燃油量在除少數(shù)波動點外，趨勢上是單調(diào)遞減的，只要使用窗口寬度足夠大，就可以保證窗口滑動過程中進入窗口的點小于等于流出窗口的點，從而完全去除剩余燃油量數(shù)據(jù)中的非單調(diào)點。為檢驗不同窗口寬度下中值濾波去除剩余燃油量波動點的效果，對現(xiàn)有1280架次飛參數(shù)據(jù)中的剩余燃油量數(shù)據(jù)進行統(tǒng)計分析，剩余燃油量數(shù)據(jù)中波動點占所有采樣點的比例隨中值濾波窗口寬度的變化規(guī)律如圖2所示。從圖2中可以看出，采樣點中波動點所占的比例隨中值濾波窗口寬度的增加而減少，當窗口寬度增大到107 時，波動點所占的比例降為0。因此，在去除剩余燃油量波動點過程中，只需要將中值濾波的窗口寬度設(shè)置為大于等于107的奇數(shù)就可以保證完全去除剩余燃油量中的波動點。

圖2 剩余燃油量數(shù)據(jù)中波動點占所有采樣點的比例隨窗口寬度的變化規(guī)律

剩余燃油量數(shù)據(jù)中帶有野值點和波動點的數(shù)據(jù)進行清洗的結(jié)果如圖3所示，從圖3中可以看到，604s處附近有野值點存在，經(jīng)過去除野值點運算后野值點消除，618～632s處、646～650s處和673～693s處均存在不同程度的波動點，使用窗口寬度為107的中值濾波算法進行濾波后波動點消除，經(jīng)過數(shù)據(jù)清洗的剩余燃油量數(shù)據(jù)呈嚴格階梯遞減狀。

圖3 剩余燃油量的數(shù)據(jù)清洗

2 剩余燃油量的插值

記數(shù)據(jù)清洗后剩余燃油量中數(shù)值序列為

因此G 的形式為階梯狀遞減，而真實的剩余燃油量應(yīng)為連續(xù)遞減，因此可以在每個階梯上選擇一個臨界點以代替此階梯上的所有值。

定義 如果樣本點滿足

則稱該樣本點為臨界點 （見圖3）。

圖3標明了臨界點的位置，臨界點位于每一個階梯的末尾，其包含了剩余燃油量遞減過程中的所有信息量，因此可以舍棄臨界點以外的所有采樣點，只對臨界點插值即可得到連續(xù)的剩余燃油量估計曲線。

插值方法的選擇會較大程度上影響特征提取結(jié)果的精度。經(jīng)總結(jié)，對臨界點插值的方法需滿足如下要求：

（1）插值所得的曲線需可導；

（2）插值曲線中的高頻分量應(yīng)盡量低以減少特征提取結(jié)果的震蕩幅度；

（3）插值方法需保持臨界點的單調(diào)性，且插值曲線斜率的絕對值應(yīng)大于某常數(shù)，其物理意義為燃油流量始終大于飛機發(fā)動機最小穩(wěn)定工作狀態(tài)時的供油量。

綜合以上3點要求，在臨界點的插值過程中應(yīng)選擇一種低階多項式的Hermite插值方法，且該插值方法應(yīng)具有約束插值斜率的功能?，F(xiàn)有對曲線斜率可調(diào)插值方法的研究主要集中于保單調(diào)插值方法［6－13］，能夠進一步限制斜率取值范圍的方法國內(nèi)外未見發(fā)表。本文在研究現(xiàn)有保單調(diào)插值方法的基礎(chǔ)上對二次樣條插值方法進行改進，提出一種適用于剩余燃油量臨界點插值的斜率約束二次樣條Hermite插值方法。

首先來討論一下兩點間的二次樣條Hermite插值問題，這是二次樣條Hermite插值問題的基礎(chǔ)。

問題1：給定兩點P1＝（t1，z1）和P2＝（t2，z2），（t1＜t2）。假設(shè)s1和s2分別為點P1和P2處的插值斜率。尋找一個函數(shù)f ∈C1［t1，t2］，使得

以下定理表明，在特定條件下，問題1 可以使用二次多項式樣條解決。

定理1 問題1可以使用二次多項式樣條解決的充要條件為

式中：Δ＝（z2－z1）／（t2－t1）。

證明：結(jié)論是顯然的。

當式 （5）成立時，兩結(jié)點間的插值方式如圖4 （a）所示，區(qū)間［t1，t2］上的插值函數(shù)為

其中，A ＝z1，B ＝s1，C ＝（s2－s1）／2（t2－t1）。

當式 （5）不成立時，兩結(jié)點間的插值方式如圖4 （b）所示，需要對區(qū)間增加一個細分點ξ，珋s為細分點處的斜率

其中，a＝s1－Δ，b＝s2－Δ。

區(qū)間［t1，ξ］上的插值函數(shù)為

其中，A ＝z1，B ＝s1，C ＝－s1）／2（ξ－t1）。區(qū)間［ξ，t2］上的插值函數(shù)為

其中，A ＝z1＋（ξ－t1）（珋s＋s1）／2，B ＝珋s，C ＝（s2－）／

2（t2－ξ）。

在點列的插值問題中，每相鄰的兩點可以看做一個兩點間的插值問題，并用定理1解決。

具有斜率約束的插值可以描述為以下問題：

問題2：給定點列

滿足

式中：K——斜率約束系數(shù)。

尋找一組插值函數(shù)f ∈C1［t1，tn］使得

以下定理可以解決問題2。

圖4 兩點間的二次樣條Hermite插值

定理2 當插值點處的斜率滿足時，問題2可以在定理1框架內(nèi)解決。式 （14）中，si為點Pi處的斜率。

證明：不失一般性，我們假設(shè)Pn單調(diào)增，令

（1）當a＋b＝0時，式 （5）成立，區(qū)間［ti，ti＋1］內(nèi)無細分點，因此

由式 （14）可得si≥K，所以f'（t）≥K；

（2）當a＋b≠0且a·b≥0時，將式 （14）代入式（8）可得

因此，f'（t）≥min（si，，si＋1）≥K；

基于以上討論，給出一種斜率約束的二次樣條Hermite插值算法如下：算法1：

步驟1 讀入點列；

步驟2 由式 （14）計算si；

步驟3 由式 （7）～式 （10）計算pi和pi＋1間的插值函數(shù)fi；

將算法1應(yīng)用于臨界點構(gòu)成的點列GC上可以得到連續(xù)的剩余燃油量曲線，從而燃油流量估計值為

3 實例分析

在某一架次飛機地面試車過程中對燃油流量進行了精確測量，選取該架次實測燃油流量可以對本文提出的方法進行驗證。已知該型飛機燃油流量最低值約為0.245L／s，因此斜率約束二次樣條Hermite插值算法中取斜率約束系數(shù)K＝－0.245，數(shù)據(jù)清洗過程中中值濾波窗口寬度為107。剩余燃油量數(shù)據(jù)清洗及臨界點插值過程如圖5所示，從圖5中可以看出：數(shù)據(jù)清洗過程去除了剩余燃油量中的野值點和波動點，插值后的剩余燃油量曲線連續(xù)平滑并單調(diào)遞減。燃油流量估計值與燃油流量實測值的比較如圖6所示，實線為燃油流量實測值，虛線為燃油流量估計值，從圖6中可以看出：發(fā)動機自423s發(fā)動機開車后，燃油流量估計值取值均大于0.245L／s，且燃油流量在0.245L／s附近時燃油流量估計值與實測值擬合程度很高，說明斜率約束的二次樣條Hermite插值有效保證了插值斜率范圍和燃油流量估計值范圍。在880s到1000s區(qū)間內(nèi)燃油流量值較高，此時燃油流量估計值與實測值擬合程度稍差。由此可以看出，本方法通過約束剩余燃油量的最小斜率，增加了燃油流量較低時的預測精度。

燃油流量的估計精度可以由相對預測誤差 （relative prediction error，RPE）來衡量，相對預測誤差的計算其公式為

通過計算得到：本實驗中燃油流量特征提取方法的RPE為4.13%，其中燃油流量值低于1.0L／s部分RPE 為1.53%，高于1.0L／s部分RPE 為8.61%。當燃油流量值較高時，燃油流量估計值與燃油流量實測值相差較大，這是由于剩余燃油量在測量時即存在誤差 （原理誤差、安裝誤差、傳輸誤差、干擾誤差等），這部分誤差在此方法中不能被有效濾除。當燃油流量值較低時，燃油流量估計值與燃油流量實測值相差較小。這說明斜率約束的燃油流量插值方法有效抑制了燃油流量值較低部分的誤差，提高了燃油流量估計的精度。因此，利用該方法得到的燃油流量值較適合于的狀態(tài)監(jiān)控。

4 結(jié)束語

本文研究了從參數(shù)中提取其一階導量的特征提取問題，并在從剩余燃油量中提取燃油流量過程中進行了實現(xiàn)，形成了一種燃油流量特征提取方法。在剩余燃油量的數(shù)據(jù)清洗過程中，引入了中值濾波理論，實驗結(jié)果表明當窗口寬度達到一定值時，中值濾波可以完全去除剩余燃油量中的波動點。在插值環(huán)節(jié)中，本文設(shè)計了斜率約束二次樣條Hermite插值方法有效保證了所求燃油流量的取值范圍。實例分析結(jié)果表明，該估計方法得出的燃油流量值與實測燃油流量值在發(fā)動機低轉(zhuǎn)速條件下具有較高的重合程度，估計精度較高，適用于所有飛參記錄系統(tǒng)中記錄了剩余燃油量的機型，具有較高的應(yīng)用價值。同時，本文提出的方法對其他一階導量特征提取問題具有借鑒作用。

［1］GENG Hong，JIE Jun.Aviation fuel flow regression model of aircraft cruise based on QAR data ［J］.Aeroengine，2008，34（4）：46－50 （in Chinese）.［耿宏，揭俊.基于QAR 數(shù)據(jù)的飛機巡航段燃油流量回歸模型 ［J］.航空發(fā)動機，2008，34 （4）：46－50.］

［2］LIU Jing.The aircraft fuel estimation model based on flight data analysis［D］.Nanjing：Nanjing University of Aeronautics and Astronautics，2010 （in Chinese）. ［劉婧.基于飛行數(shù)據(jù)分析的飛機燃油估計模型 ［D］.南京：南京航空航天大學，2010.］

［3］LI Yingying，TAN Guangyu，QU Jianling.Prehandling research of flight parameter strange points based on matlab digital filtering technique ［J］.Electric Machines and Control，2009，13 （1）：135－138 （in Chinese）. ［李映穎，譚光宇，曲建嶺.飛行記錄參數(shù)奇異點的處理 ［J］.電機與控制學報，2009，13（1）：135－138.］

［4］WU Jiangang，CHEN Zhiwei，LI Shulin，et al.Study on computer processing of recorded flight data［J］.Computer Simulation，2007，24 （2）：18－21 （in Chinese）.［吳建剛，陳志偉，李曙林，等.飛參記錄數(shù)據(jù)計算機處理的有關(guān)問題研究 ［J］.計算機仿真，2007，24 （2）：18－21.］

［5］Gonzalez R C，Woods R E.Digital image processing ［M］.3rd ed.New Jersey：Pearson，2010.

［6］Hussain M Z，Sarfraz M，Hussain M.Scientific data visualization with shape preserving C1rational cubic interpolation ［J］.Eur J Pure Appl Math，2010，3 （2）：194－212.

［7］Volkov Yu S，Bogdanov V V，Miroshnichenko V L，et al.Shape－preserving interpolation by cubic splines［J］.Mathematical Notes，2010，88 （5）：798－805.

［8］Sarfraz Muhammad H M Z.Hussain maria shape－preserving curve interpolation ［J］.International Journal of Computer Mathematics，2012，89 （1）：35－53.

［9］ZHU Xiaolin，GE Chuanfeng.C～3convex T－B interpolating curve and an algorithm for shape preserving interpolation ［J］.Journal of Engineering Graphics，2009，30 （6）：76－82 （in Chinese）.［朱曉臨，葛傳豐.連續(xù)的保凸T－B插值曲線及保形插值算法 ［J］.工程圖學學報，2009，30 （6）：76－82.］

［10］Feng Renzhong L F.A shape－preserving quasi－interpolation operator satisfying quadratic polynomial reproduction property to scattered data ［J］.Journal of Computational and Applied Mathematics，2009，225 （2）：594－601.

［11］Yong J P，Min C.A survey on shape preserving interpolation［C］／／International Conference on Digital Object Identifier，2012：2521－2523.

［12］Chen Z，Su K.Research on a class of shape－preserving quasiinterpolation functions ［C］／／2nd International Conference on Digital Object Identifier，2010：1442－1445.

［13］Meng T.A shape preserving rational cubic interpolation spline［C］／／2nd International Conference on Digital Object Identifier，2010：1－4.