基于實時中值估計的分布式穩(wěn)健濾波融合方法

黃 河，陳慶良，張 旭

(中國洛陽電子裝備試驗中心，河南 洛陽 471000)

0 引言

在電子信息裝備多源測量數(shù)據(jù)融合處理工作中，通常會遇到復(fù)雜電磁環(huán)境造成的連續(xù)粗差。這時，基于白噪聲假設(shè)的傳統(tǒng)方法往往表現(xiàn)不好，需要研究數(shù)據(jù)融合的穩(wěn)健方法。

關(guān)聯(lián)、配準(zhǔn)、檢擇、濾波、加權(quán)是數(shù)據(jù)融合中相互聯(lián)系、相互影響的幾個過程，相關(guān)的處理方法直接關(guān)系到最終的數(shù)據(jù)質(zhì)量，本文主要關(guān)注其中的穩(wěn)健濾波、穩(wěn)健加權(quán)問題。對此，很多文獻(xiàn)進(jìn)行了討論，提出了許多不同方法。魏喜慶等[1]將高斯過程回歸引入到容積卡爾曼濾波之中，文獻(xiàn)[2—4]也介紹了一些穩(wěn)健濾波方法，但卻均未對樣本崩潰點進(jìn)行分析。加權(quán)融合方面，文獻(xiàn)[5—8]將權(quán)值最優(yōu)分配原則與卡爾曼濾波結(jié)合應(yīng)用，但未考慮穩(wěn)健性問題；文獻(xiàn) [9]采用M估計實現(xiàn)參數(shù)穩(wěn)健估計和權(quán)值最優(yōu)分配原則下的加權(quán)融合，但只適用于靜態(tài)數(shù)據(jù)融合，未討論運動目標(biāo)測量數(shù)據(jù)的穩(wěn)健濾波融合方法。

本文針對電子信息裝備運動目標(biāo)多源測量數(shù)據(jù)加權(quán)融合中的復(fù)雜電磁環(huán)境影響，在以往工作[10-12]的基礎(chǔ)上，研究提出了基于實時中值估計的分布式穩(wěn)健濾波融合方法。

1 權(quán)值最優(yōu)分配原則應(yīng)用中的問題

1.1 融合權(quán)值最優(yōu)分配原則

令加權(quán)因子分別為W1，W2，…，WS，則

(1)

可知總均方誤差為：

(2)

可見，總均方誤差是關(guān)于各加權(quán)因子的多元二次函數(shù)，因此必然存在最小值。根據(jù)多元函數(shù)求極值理論，可求出總均方誤差最小時所對應(yīng)的加權(quán)因子為：

(3)

則有

(4)

其對應(yīng)的最小均方誤差為：

(5)

1.2 問題提出

電子信息裝備多源測量數(shù)據(jù)融合處理中，在應(yīng)用上述精度構(gòu)建融合權(quán)重的權(quán)值最優(yōu)分配原則時，由于動態(tài)測量中設(shè)備精度的不平穩(wěn)性，以及沒有絕對意義上的真值，需要對設(shè)備精度進(jìn)行動態(tài)的近似估計，從而進(jìn)行自適應(yīng)的加權(quán)融合。這種情況下，由于常用的標(biāo)準(zhǔn)差形式的精度估計方法不具備穩(wěn)健性，以及復(fù)雜電磁環(huán)境下粗差的影響，需要研究穩(wěn)健的濾波方法和穩(wěn)健的精度估計方法。

本文構(gòu)建了一種高樣本崩潰點的穩(wěn)健濾波方法，并用穩(wěn)健濾波殘差刻度估計作為動態(tài)精度構(gòu)建融合權(quán)重，從而實現(xiàn)穩(wěn)健的自適應(yīng)加權(quán)融合，這里稱作穩(wěn)健刻度法。

2 基于實時中值估計的穩(wěn)健濾波方法

2.1 中值估計的穩(wěn)健性

中值(中位數(shù))估計是一種基本的穩(wěn)健估計方法。對采樣序列yi，其中值為：

(6)

式(6)中，y(j)表示對數(shù)據(jù){y1,y2,…,ym}按從大到小排序后的第j個數(shù)值。

中值估計是按極小化極大準(zhǔn)則的一種最優(yōu)估計，其影響函數(shù)有界，樣本崩潰點接近50%，因此中值估計有良好的穩(wěn)健性，可以提供工程應(yīng)用中十分重要的高樣本崩潰點這一性質(zhì)。因此，本文研究提出以下基于實時中值估計的穩(wěn)健遞推濾波方法。

2.2 基于雙向鏈表排序的實時中值估計算法

由于中值估計需要對采樣序列進(jìn)行排序，作者從節(jié)省內(nèi)存、減少運算的角度出發(fā)，在考慮了數(shù)組、單鏈表、雙向鏈表、隊列等數(shù)據(jù)結(jié)構(gòu)后，選定帶順序號索引的雙向鏈表作為算法實現(xiàn)的基本數(shù)據(jù)結(jié)構(gòu)，排序算法選用插入排序，因為它更適合向有序表中添加元素。

應(yīng)用中應(yīng)先對雙向鏈表進(jìn)行相應(yīng)聲明，創(chuàng)建、初始化頭節(jié)點，并根據(jù)實際需要設(shè)定雙向鏈表節(jié)點數(shù)目。節(jié)點刪除時，應(yīng)依據(jù)先入先出的原則，按照插入順序索引尋找應(yīng)刪除的節(jié)點；節(jié)點插入時，按照x值大小進(jìn)行排序，并對相應(yīng)節(jié)點(新節(jié)點及大小順序發(fā)生變化的原有節(jié)點)的大小順序索引進(jìn)行賦值。排序完成后，即可得到中值。實時中值估計主要算法流程如圖1所示。

2.3 穩(wěn)健濾波

2.3.1預(yù)測

(7)

(8)

圖1 排序雙向鏈表節(jié)點插入與刪除算法流程圖Fig.1 The main flow for inserting and deleting node of the ordered doubly-linked list

2.3.2濾波

(9)

(10)

(11)

至此，對勻速直線運動，我們構(gòu)造了基于一階差分及其誤差刻度中值估計的實時穩(wěn)健濾波方法，它實際上共包含三重中值估計，因而是一種高樣本崩潰點的穩(wěn)健遞推濾波方法。對于勻加速直線運動，可以應(yīng)用類似的方法用二階差分構(gòu)造相應(yīng)的穩(wěn)健濾波。

實際應(yīng)用中還應(yīng)考慮目標(biāo)機(jī)動對算法的影響，據(jù)參考文獻(xiàn)[13],目標(biāo)緩慢轉(zhuǎn)彎往往引起長達(dá)1 min的相關(guān)加速，目標(biāo)作躲避機(jī)動時的相關(guān)加速時間為10～30 s，而由大氣湍流或飛機(jī)自身隨機(jī)因素引起的相關(guān)加速度時間只有1～2 s。綜合起來考慮，為避免將機(jī)動當(dāng)作干擾處理，算法的設(shè)計樣本容量在時間尺度上應(yīng)不超過2 s。

3 融合方法及算法流程

3.1 分布式融合方法

(12)

3.2 算法流程

在工程應(yīng)用中，穩(wěn)健濾波融合方法中的實時中值估計和系統(tǒng)誤差配準(zhǔn)[11]均可應(yīng)用上述基于雙向鏈表排序的實時中值估計算法(見2.2節(jié))。圖2為基于實時中值估計的穩(wěn)健濾波算法流程。

圖2 基于實時中值估計的穩(wěn)健遞推濾波算法流程Fig.2 The robust filtering algorithm based on real time median estimation

圖3為基于穩(wěn)健濾波殘差刻度估計的分布式自適應(yīng)加權(quán)融合算法流程。

圖3 基于穩(wěn)健濾波殘差刻度的加權(quán)融合算法流程Fig.3 The distributed weighted fusion algorithm based on robust filtering residual scale estimation

4 測試結(jié)果分析

在A、B兩雷達(dá)參加的某單目標(biāo)動態(tài)測量中，目標(biāo)大部分航段作勻速直線飛行，部分航段(圖4—圖7中5 500～7 000航段)作弱機(jī)動飛行，A雷達(dá)部分航段含連續(xù)干擾值，在將兩雷達(dá)測量數(shù)據(jù)轉(zhuǎn)至同一直角坐標(biāo)系并配準(zhǔn)后，用本文方法對兩雷達(dá)數(shù)據(jù)進(jìn)行穩(wěn)健濾波和加權(quán)融合測試。

將本文方法對A雷達(dá)的濾波結(jié)果與卡爾曼濾波進(jìn)行比較(分別采用勻速模型和當(dāng)前統(tǒng)計模型卡爾曼濾波進(jìn)行對比，勻速模型卡爾曼濾波中，過程噪聲Q分別取1，9，20倍的Q0，Q0為一確定過程噪聲)，并以更高精度的GNSS測量數(shù)據(jù)作為相對真值，得到濾波融合測試結(jié)果。

結(jié)果顯示(表1，圖4—圖6)，與卡爾曼濾波相比，本文方法濾波結(jié)果受粗差的影響較小，具有較高的樣本崩潰點。在將濾波樣本容量設(shè)為39時，本文穩(wěn)健濾波方法能夠濾除長達(dá)15的連續(xù)干擾值。

在對A雷達(dá)含干擾段數(shù)據(jù)的濾波測試中，勻速模型卡爾曼濾波發(fā)生了視在發(fā)散，其中圖4(b)，(c)，(d)的濾波發(fā)散是由于量測噪聲偏離假定模型(粗差)引起的，圖5(b)，(c)，(d)的濾波發(fā)散是由于弱機(jī)動航段系統(tǒng)狀態(tài)偏離假定模型引起的，兩種情形下濾波發(fā)散程度均受所設(shè)定過程噪聲的影響。當(dāng)前統(tǒng)計模型卡爾曼濾波對機(jī)動段處理結(jié)果較好(圖4(e),圖5(e)中5 500～7 000航段)，但當(dāng)遇到連續(xù)粗差時，處理結(jié)果較差(容易將連續(xù)干擾值判為機(jī)動)。

表1 A雷達(dá)含干擾航段濾波測試有關(guān)數(shù)據(jù)Tab.1 Relevant data in the filtering test of the disturbed measurement data of radar A

可見，由于卡爾曼濾波需要確知狀態(tài)方程，并受過程噪聲、量測噪聲的影響[14]，在復(fù)雜電磁環(huán)境下的電子信息裝備測量數(shù)據(jù)處理工作中，當(dāng)卡爾曼濾波假定模型與測量數(shù)據(jù)真實模型不相符時，濾波效果受到較大影響。

在弱機(jī)動條件下，本文濾波方法也具有一定的機(jī)動濾波能力，在目標(biāo)作弱機(jī)動前后仍是有效的(圖6(a)，(b)，(c)，(d)，圖4(f), 圖5(f))；同時加權(quán)融合測試也達(dá)到了較好的穩(wěn)健效果，融合結(jié)果受粗差影響較小(圖7)。

圖4 A雷達(dá)Y方向濾波前后的一次差Fig.4 Error of Radar A in Y direction before and after filtering

圖5 A雷達(dá)Z方向濾波前后的一次差Fig.5.Error of Radar A in Z direction before and after filtering

圖6 兩雷達(dá)某測量序列濾波融合前后的一階差分Fig.6 The first difference of a certain measurement sequence before and afterdata processing

5 結(jié)論

本文提出了基于實時中值估計的分布式穩(wěn)健濾波融合方法，該方法以穩(wěn)健統(tǒng)計理論為基礎(chǔ)，針對電子信息裝備多源測量數(shù)據(jù)加權(quán)融合中的復(fù)雜電磁環(huán)境干擾和傳統(tǒng)方法的不穩(wěn)健性，給出基于實時中值估計和勻速模型的穩(wěn)健濾波方法，并將穩(wěn)健濾波殘差刻度應(yīng)用于權(quán)值最優(yōu)分配原則下的融合權(quán)重動態(tài)計算，最終得到穩(wěn)健加權(quán)融合結(jié)果。

典型運動目標(biāo)實測數(shù)據(jù)測試結(jié)果表明，本文方法具有較高的樣本崩潰點，能有效處理含一定長度連續(xù)粗差的測量數(shù)據(jù)，且能用于弱機(jī)動情形。與卡爾曼濾波采用最小均方誤差為估計準(zhǔn)則相比，本文濾波方法所包含的的中值估計屬于穩(wěn)健估計，這也是本文濾波方法相較卡爾曼濾波更具穩(wěn)健性的根本原因。由于本文方法是在標(biāo)量(或解耦)意義下針對勻速直線運動目標(biāo)位置數(shù)據(jù)得到的，所以在以下條件下是通用的：測量數(shù)據(jù)為等周期時間序列；狀態(tài)方程近似于勻速模型，或當(dāng)系統(tǒng)狀態(tài)模型與勻速模型發(fā)生偏離時，偏離應(yīng)較小或較慢，以確保過程噪聲的影響小于量測噪聲的影響。在上述條件下，本文方法有以下特點：樣本崩潰點較高，能夠濾除干擾等因素造成的連續(xù)野值；無需精確的狀態(tài)方程、過程噪聲、量測噪聲等先驗統(tǒng)計信息；濾波過程實現(xiàn)了對量測噪聲的自適應(yīng)。

當(dāng)然，本文對穩(wěn)健濾波融合方法的研究還不夠充分，給出了一種勻速模型和標(biāo)量(或解耦)意義下的穩(wěn)健濾波融合方法，作者將在以后的工作中結(jié)合其他穩(wěn)健方法開展進(jìn)一步研究。