基于滑動(dòng)窗濾波與卡爾曼濾波的信息融合

黃金虎 喬 祁

（1、中國(guó)直升機(jī)設(shè)計(jì)研究所，江西 景德鎮(zhèn) 333001 2、中航技進(jìn)出口有限責(zé)任公司，北京 102600）

1 總體設(shè)計(jì)思路

本文以無(wú)人直升機(jī)活塞發(fā)動(dòng)機(jī)的轉(zhuǎn)速信息融合為研究對(duì)象，具體實(shí)施過(guò)程如圖1 所示。首先分析發(fā)動(dòng)機(jī)的開車數(shù)據(jù)，建立對(duì)應(yīng)的離散化數(shù)據(jù)模型；其次是確定初始參數(shù)、融合過(guò)程；最后分析對(duì)比融合數(shù)據(jù)與原始數(shù)據(jù)。

圖1 設(shè)計(jì)流程圖

活塞發(fā)動(dòng)機(jī)的轉(zhuǎn)速信號(hào)主要包含三路傳感器采集，分別是頻率采集單元（脈沖信號(hào)）、發(fā)動(dòng)機(jī)采集單元TCU（數(shù)字信號(hào)）、旋翼轉(zhuǎn)速信號(hào)（脈沖信號(hào)）。本文最終目的是在飛控計(jì)算機(jī)內(nèi)具備實(shí)時(shí)融合三路傳感器為一路，同時(shí)如果有一路信號(hào)出現(xiàn)故障時(shí)，能自動(dòng)剔除該信號(hào)。

三路傳感器的之間相互獨(dú)立，額定轉(zhuǎn)速5500rpm，旋翼轉(zhuǎn)速乘傳動(dòng)比9.524，完整開車轉(zhuǎn)速信號(hào)數(shù)據(jù)如圖2 所示，離合器嚙合后，旋翼轉(zhuǎn)速與發(fā)動(dòng)機(jī)轉(zhuǎn)速重疊。在時(shí)間t=515～1672s 區(qū)間，發(fā)動(dòng)機(jī)轉(zhuǎn)速處于額定閉環(huán)控制狀態(tài)。在時(shí)間t=1672～1900s 區(qū)間，發(fā)動(dòng)機(jī)轉(zhuǎn)速處于停車過(guò)程。

圖2 各路信號(hào)數(shù)據(jù)

2 建立數(shù)學(xué)模型

發(fā)動(dòng)機(jī)轉(zhuǎn)速特性離散化的模型及觀測(cè)模型如下：

式中：x 為發(fā)動(dòng)機(jī)的轉(zhuǎn)速函數(shù)，k 為時(shí)間計(jì)數(shù)，f 為狀態(tài)轉(zhuǎn)移矩陣，B 為控制輸入矩陣，μ 為外部輸入包括油門開度、扭矩、活塞缸空氣輸入、環(huán)境溫度等，q 為發(fā)動(dòng)機(jī)本身系統(tǒng)偏差造成的高斯噪聲，包含缸體振動(dòng)、氣動(dòng)負(fù)載變化、油品、其他環(huán)境因素等。z 為各傳感器的觀測(cè)值包括發(fā)動(dòng)機(jī)轉(zhuǎn)速傳感器（脈沖信號(hào)）、TCU 采集、旋翼傳感器（脈沖信號(hào)），h 為轉(zhuǎn)速到傳感器映射矩陣，v 為各傳感器采樣時(shí)的高斯噪聲。

3 窗空間

窗函數(shù)在信號(hào)處理中是指一種除在給定區(qū)間之外取值均為0 的實(shí)函數(shù)。本文借鑒此類方法，對(duì)有限時(shí)間范圍的波形進(jìn)行分析。本文的運(yùn)算周期tick=40ms，建立一個(gè)按照先進(jìn)先出原則的樣本分析窗空間speed。其中speed 為13×3 矩陣，每列分別對(duì)應(yīng)于頻率采集單元（脈沖信號(hào)）、發(fā)動(dòng)機(jī)采集單元（數(shù)字信號(hào)）、旋翼轉(zhuǎn)速信號(hào)（脈沖信號(hào)）數(shù)值。

4 濾波參數(shù)的初始化

發(fā)動(dòng)機(jī)的輸入矩陣為μ，設(shè)μ1 為油門值、μ2 為扭矩負(fù)載值、μ3 活塞缸空氣輸入、μ4 環(huán)境溫度等，則外部輸入向量為μ=[μ1，μ2, μ3，μ4]T，發(fā)動(dòng)機(jī)轉(zhuǎn)速函數(shù)為x(μ1，μ2, μ3，μ4)，則控制陣B 為：

分析完整的開車時(shí)間歷程數(shù)據(jù)，發(fā)動(dòng)機(jī)在額定狀態(tài)下，轉(zhuǎn)速近似一條水平直線，即宏觀上發(fā)動(dòng)機(jī)轉(zhuǎn)速對(duì)任意輸入的偏導(dǎo)數(shù)為0，則控制矩陣為0。實(shí)際上，采用了PID 閉環(huán)控制的活塞發(fā)動(dòng)是一套穩(wěn)定的總距隨動(dòng)控制系統(tǒng)。在額定工況下，發(fā)動(dòng)的轉(zhuǎn)速，在任意負(fù)載下，轉(zhuǎn)速的波動(dòng)范圍為（-50，50）。引起此類偏差是由發(fā)動(dòng)機(jī)氣缸特性、氣動(dòng)負(fù)載、油品、環(huán)境造成的，屬于高斯噪聲；如果想提高模型的預(yù)測(cè)準(zhǔn)確性，實(shí)際使用中可通過(guò)黑盒法辨識(shí)出發(fā)動(dòng)機(jī)轉(zhuǎn)速模型，本文不做擴(kuò)展。實(shí)際上對(duì)于機(jī)械慣性系統(tǒng)，在很短的時(shí)間的微小輸入都可以看作是噪聲干擾。因此公式(1)可以簡(jiǎn)化為：

式中：帶下標(biāo)x_(k)為轉(zhuǎn)速在本次tick 的先驗(yàn)估計(jì)值，x(k-1)為上次tick 的轉(zhuǎn)速后驗(yàn)估計(jì)值(最優(yōu)估計(jì))，f 為狀態(tài)方程。

4.1 狀態(tài)方程f：在一個(gè)tick=40ms 下，默認(rèn)狀態(tài)下，上一tick 的采樣與本次tick 的采樣相同，轉(zhuǎn)速不變，f=1。

4.2 過(guò)程噪聲q:分析歷程數(shù)據(jù)，得到模型高斯噪聲q=1.5。

4.3 估計(jì)值協(xié)方差p：初始化p=1，需要迭代過(guò)程實(shí)時(shí)跟新。

4.4 先驗(yàn)估計(jì)協(xié)方差p_：初始化p_=1，需要迭代過(guò)程實(shí)時(shí)跟新。

4.5 窗空間列方差v：初始化v=[0 0 0]，用于存儲(chǔ)speed 每列的方差，不斷更新觀測(cè)方程的方差矩陣。

有時(shí)候我們?cè)诓牧仙险也坏矫鞔_的主題或標(biāo)題，而問(wèn)題設(shè)問(wèn)上也沒(méi)有具體對(duì)象的限制。那么我們可能遇到了“全開”式論題，我們只能根據(jù)材料信息，將零碎的信息源整合成一個(gè)對(duì)象，以此為基礎(chǔ)進(jìn)行“要義”式、“推因”式或“推果”式立論。若是零碎的材料信息可以整合成兩個(gè)對(duì)象，可以考慮在二者之間建立聯(lián)系，使之成為簡(jiǎn)約句子結(jié)構(gòu)的主語(yǔ)或賓語(yǔ)。下面請(qǐng)看例3：

4.6 測(cè)量噪聲協(xié)方差r：初始化r=[1 0 0;0 1 0;0 0 1]; 3×3 矩陣，用于求取測(cè)量殘差，需要使用v 對(duì)對(duì)角線賦值，由于各傳感器相互獨(dú)立，因此其他元素賦值0。

4.7 觀測(cè)矩陣z：初始化z=[0;0;0]; 1×3 的列向量，z1、z2、z3 分別對(duì)應(yīng)頻率采集單元（脈沖信號(hào)）、TCU（數(shù)字信號(hào)）、旋翼轉(zhuǎn)速信號(hào)（脈沖信號(hào)）最新數(shù)值。

4.8 觀測(cè)到狀態(tài)矩陣h：由公式2 得，h=[1;1;1]，是1維列向量，用于映射傳感器測(cè)試量與預(yù)測(cè)值，統(tǒng)一它們的量綱。

4.9 轉(zhuǎn)速最優(yōu)估計(jì)x：初始化x=0。

4.10 濾波增益矩陣K：初始化K=[0 0 0]，需要迭代過(guò)程實(shí)時(shí)跟新。

5 迭代運(yùn)算

融合的迭代運(yùn)算主要以卡爾曼濾波[2]的為基礎(chǔ)進(jìn)行，但系統(tǒng)噪聲及觀察噪聲做了相應(yīng)調(diào)整。

5.1 卡爾曼第一公式

式中：狀態(tài)轉(zhuǎn)移公式x_是本次的先驗(yàn)估計(jì)，f 是狀態(tài)轉(zhuǎn)移矩陣，x 是上次的最優(yōu)估計(jì)。

5.2 卡爾曼第二公式

式中：噪聲協(xié)方差矩陣的傳遞公式，預(yù)測(cè)模型帶來(lái)的噪聲q，表示不確定性在各公式的傳遞。

5.3 人工干預(yù)通道

為了避免矩陣計(jì)算機(jī)運(yùn)算過(guò)程中出現(xiàn)非法數(shù)，需要對(duì)窗空間的各列向量均方根值v 做出限定。設(shè)定上邊界v(1)Max=25000,v(2)Max=25000,v(3)Max=25000。

設(shè)定下邊界，使用拋物線函數(shù)，對(duì)稱軸由實(shí)際開車均方根值數(shù)據(jù)得到。極小值可人工設(shè)置，本文稱之為“人工干預(yù)通道”，通過(guò)調(diào)整各傳感器極小值來(lái)主動(dòng)選擇置信度，極值越小置信度越高。

5.4 余度管理與控制

當(dāng)某路傳感器出現(xiàn)背離其它傳感器信號(hào)時(shí)，采用三判二的原則，方差v 賦極大值，能使其融合權(quán)重占比降至最低。

5.5 測(cè)量殘差協(xié)方差

5.6 卡爾曼第三公式

卡爾曼系數(shù)的更新，其中K 為1×3 的列向量，每個(gè)元素對(duì)應(yīng)的數(shù)學(xué)意義分別表示頻率采集、TCU 采集單元、旋翼轉(zhuǎn)速采集在轉(zhuǎn)速融合中所占之權(quán)重。預(yù)估值在轉(zhuǎn)速信息融合中所占權(quán)重則為1-K×h。

5.7 卡爾曼第四公式

后驗(yàn)估計(jì)（最優(yōu)估計(jì)）=模型數(shù)據(jù)后驗(yàn)估計(jì)+k×（觀測(cè)殘差）。

5.8 卡爾曼第五公式

預(yù)估模型的噪聲協(xié)方差矩陣的更新。

5.9 窗空間的更新

程序等待一tick 后，更新speed 的13×3 的矩陣，以最新的轉(zhuǎn)速采樣z1k+1、z2k+1、z3k+1替換上一tick 的值，隊(duì)列其他元素依次向前，z1k-12、z2k-12、z3k-12丟棄，按照“先進(jìn)先出”的原則往前替換，如公式16、17 所示。

求取k+1 時(shí)刻speed 的每列均方根值，更新窗空間方差v，用于不斷校準(zhǔn)觀測(cè)方程的方差矩陣。程序跳轉(zhuǎn)至卡爾曼第一公式。

6 融合效果與分析

飛控計(jì)算機(jī)分別記錄下k 矩陣，轉(zhuǎn)速融合最優(yōu)估計(jì)x，畫圖分析如下。

6.1 額定階段的轉(zhuǎn)速與權(quán)重

選取時(shí)間段（1015～1050s）波動(dòng)最大部分分析，如圖3所示。

圖3 額定階段采樣值與融合值

k1、k2 、k3、k4 分別表示頻率采集、發(fā)動(dòng)機(jī)采集單元、旋翼轉(zhuǎn)速采集、預(yù)估值在轉(zhuǎn)速信息融合中所占權(quán)重。其中k1+k2+k3+k4=1。在發(fā)動(dòng)機(jī)額定狀態(tài)，頻率采集的權(quán)重約50%，TCU 采集的權(quán)重約12%，旋翼采集約18%，預(yù)估值的權(quán)重約20%。這些值的置信度可以通過(guò)調(diào)節(jié)人工干預(yù)通道進(jìn)行調(diào)節(jié)。當(dāng)某一信號(hào)波動(dòng)較大時(shí)，與之對(duì)應(yīng)的占比將相應(yīng)的自動(dòng)降低如圖4 所示。

圖4 額定階段權(quán)重對(duì)比

6.2 停車階段的轉(zhuǎn)速與權(quán)重

選擇時(shí)間t=1672～1900s 區(qū)間，分析發(fā)動(dòng)機(jī)停車的過(guò)程數(shù)據(jù)如圖5 所示，離合器率先進(jìn)行脫開的過(guò)程，旋翼轉(zhuǎn)速下降，最后發(fā)動(dòng)機(jī)停車。

圖5 停車階段采樣值與融合值

發(fā)動(dòng)機(jī)的停車過(guò)程，能模擬當(dāng)旋翼傳感器失效時(shí)，融合算法經(jīng)過(guò)若干步長(zhǎng)后，能有效剔除失效信號(hào)。信息失真度低。如圖6 所示，旋翼轉(zhuǎn)速在最終融合過(guò)程中權(quán)重被調(diào)節(jié)到0%。

圖6 停車階段權(quán)重對(duì)比

7 結(jié)論

本文基于卡爾曼濾波的預(yù)測(cè)與跟新的回歸校準(zhǔn)思想，引入了譜分析中窗函數(shù)概念，改進(jìn)了常規(guī)的“三判二”式的信號(hào)選擇方式，設(shè)計(jì)了一套既能滿足時(shí)域離散系統(tǒng)、過(guò)濾信號(hào)中高斯噪聲的，又能診斷并屏蔽失效信號(hào)的、具備人工干預(yù)、余度控制融合的方法。