卡爾曼濾波粗差修正方法應(yīng)用

齊公玉,邱衛(wèi)寧,花向紅

(1.武漢大學(xué) 測(cè)繪學(xué)院,湖北 武漢 430079;2.武漢大學(xué) 災(zāi)害監(jiān)測(cè)與防治研究中心,湖北 武漢 430079)

卡爾曼濾波粗差修正方法應(yīng)用

齊公玉1,2,邱衛(wèi)寧1,2,花向紅1,2

(1.武漢大學(xué) 測(cè)繪學(xué)院,湖北 武漢 430079;2.武漢大學(xué) 災(zāi)害監(jiān)測(cè)與防治研究中心,湖北 武漢 430079)

基于標(biāo)準(zhǔn)卡爾曼濾波,提出對(duì)含有粗差的觀測(cè)值進(jìn)行整體修正和局部修正。最后,通過工程實(shí)例驗(yàn)證,兩種修正方法均能夠有效抵抗粗差的影響,并可獲得比標(biāo)準(zhǔn)卡爾曼濾波更準(zhǔn)確的濾波結(jié)果。

卡爾曼濾波;粗差修正;應(yīng)用

卡爾曼濾波理論作為一種重要的最優(yōu)估計(jì)理論被廣泛用于各種動(dòng)態(tài)數(shù)據(jù)處理中。其采用線性狀態(tài)方程描述系統(tǒng),結(jié)合觀測(cè)方程,采用遞推形式計(jì)算系統(tǒng)的最優(yōu)估值[1]。它無需存儲(chǔ)大量觀測(cè)數(shù)據(jù),便可進(jìn)行參數(shù)估計(jì)。但是,卡爾曼濾波容錯(cuò)能力很差,觀測(cè)值中有很小的粗差,會(huì)使濾波發(fā)散,使濾波估值發(fā)生大的偏移。所以,如何有效地探測(cè)和修正動(dòng)態(tài)測(cè)量數(shù)據(jù)中粗差觀測(cè)值是非常值得研究的問題。

目前,粗差的定位常用的方法有以均值漂移模型為基礎(chǔ)的粗差探測(cè)法和以污染誤差模型理論為基礎(chǔ)的抗差估計(jì)[2-3]等方法,但是,這些方法要求觀測(cè)網(wǎng)必須具有一定的幾何圖形強(qiáng)度和多余觀測(cè)數(shù)。在實(shí)際工程中,如水準(zhǔn)檢測(cè)網(wǎng)中,一般多余觀測(cè)數(shù)很少,如果將含有粗差的觀測(cè)值剔除將使一些狀態(tài)參數(shù)無法得到估計(jì)。這樣就需要對(duì)含有粗差觀測(cè)值進(jìn)行修正,而不能簡(jiǎn)單的舍棄。本文將粗差修正分為兩類—整體修正和局部修正,并根據(jù)實(shí)際的工程數(shù)據(jù),分析兩種方法均可以獲得比標(biāo)準(zhǔn)卡爾曼濾波更加準(zhǔn)確的濾波結(jié)果,因此,這種方法的研究具有實(shí)際的工程應(yīng)用價(jià)值。

1 標(biāo)準(zhǔn)卡爾曼濾波

設(shè)一線性系統(tǒng)的狀態(tài)方程和觀測(cè)方程為

式中:Xk為狀態(tài)向量,Lk為觀測(cè)向量,Φk,k-1為狀態(tài)轉(zhuǎn)移矩陣,Uk-1為控制向量,一般不考慮,Γk,k-1,Bk為系數(shù)矩陣,Ωk-1為系統(tǒng)動(dòng)態(tài)噪聲向量,Δk為觀測(cè)噪聲向量,其隨機(jī)模型為

2 粗差探測(cè)方法

令

由統(tǒng)計(jì)性質(zhì)可知,Vk是均值為0的高斯隨機(jī)量,根據(jù)式(4),可以得到其協(xié)方差陣為

式中:Dx(k/k-1)為預(yù)測(cè)誤差協(xié)方差陣,DΔ(k)為觀測(cè)誤差協(xié)方差陣。

假定濾波模型是正確的,需要判斷觀測(cè)值是否存在粗差,根據(jù)文獻(xiàn)[4]中,構(gòu)造檢驗(yàn)統(tǒng)計(jì)量 Tk=VkTDV(k)Vk～χ2(nk,0),對(duì)于給定的顯著水平α,當(dāng)Tk＜χ2a時(shí),認(rèn)為濾波正常,反之認(rèn)為濾波不正常,即觀測(cè)值存在粗差。

3 粗差修正

在判斷觀測(cè)值存在粗差后,粗差的修正可以分為2種方法,一種為整體修正,一種為局部修正。整體修正為對(duì)k期的觀測(cè)值統(tǒng)一進(jìn)行迭代運(yùn)算;而局部修正是對(duì)探測(cè)到的個(gè)別粗差觀測(cè)值進(jìn)行迭代修正。

3.1 整體修正法

設(shè)k期的觀測(cè)值Lk,探測(cè)到存在粗差觀測(cè),假設(shè) k期以前濾波均正常,并且已經(jīng)求得^X(k-1/k-1),則加入粗差修正過程的第k期卡爾曼濾波過程為

1)根據(jù)式(1)的狀態(tài)方程,可知 k期的預(yù)測(cè)值X^(k/k-1),并計(jì)算出其協(xié)方差陣;

5)迭代結(jié)束后,選擇修正后的L′k計(jì)算參數(shù)的濾波估值。

3.2 局部修正法

在文獻(xiàn)[5]中提出根據(jù)式(5)判別粗差觀測(cè)值的判別式為

式中:i,i為對(duì)角線上第i個(gè)元素;Vi(k)為V(k)第 i個(gè)分量。C為常數(shù),可根據(jù)實(shí)際物理背景確定。根據(jù)式(6),探查到含有粗差觀測(cè)值的具體位置。

設(shè)k期的觀測(cè)值Lk,探查到第 j個(gè)觀測(cè)值Lkj存在粗差,假設(shè) k期以前濾波均正常,并且已經(jīng)求得^X(k-1/k-1),則加入粗差修正過程的第 k期卡爾曼濾波局部修正過程如下:

1)根據(jù)式(1)的狀態(tài)方程,得到 k期的預(yù)測(cè)值X^(k/k-1),并計(jì)算出其協(xié)方差陣;

2)按照標(biāo)準(zhǔn)卡爾曼濾波進(jìn)行計(jì)算,得到濾波結(jié)果^X(k/k);

Lk的修正值,將第 j個(gè)觀測(cè)值換為L′k中第j個(gè)值;

4)對(duì)步驟2)、3)進(jìn)行迭代運(yùn)算,并計(jì)算每次迭代 后 的 單 位 權(quán) 方 差 因 子 σ^2

5)迭代結(jié)束后,選擇最后修正后的 L′k計(jì)算參數(shù)的濾波估值。

4 實(shí)驗(yàn)分析

為驗(yàn)證本方法的計(jì)算精度,文中采用某一地表沉降監(jiān)測(cè)項(xiàng)目中的對(duì)動(dòng)態(tài)水準(zhǔn)監(jiān)測(cè)網(wǎng)數(shù)據(jù)進(jìn)行處理。監(jiān)測(cè)網(wǎng)中共有12個(gè)未知地表點(diǎn)D1～D12和1個(gè)已知點(diǎn)。數(shù)據(jù)采樣時(shí)間為3～4 d,共監(jiān)測(cè)20期。本文對(duì)3期數(shù)據(jù)進(jìn)行處理,根據(jù)前2期確定的初值,對(duì)第3期數(shù)據(jù)進(jìn)行濾波。選取未知點(diǎn)的高程和高程的變化率為狀態(tài)參數(shù),其中Δt=4。

濾波參數(shù)初值的選取:X0選取第2期高程平差結(jié)果。X˙0選取為第1期和第2期平差結(jié)果之差與數(shù)據(jù)采樣間隔的比值。X(0/0)=X0;DΩ=E;D(0/0)=E。狀態(tài)方程和觀測(cè)方程為

數(shù)據(jù)處理過程:首先,根據(jù)僅有偶然誤差時(shí)的觀測(cè)值進(jìn)行粗平差得到結(jié)果作為真值。其次,在觀測(cè)值 H2,H9中分別加入5 cm的觀測(cè)誤差,用標(biāo)準(zhǔn)卡爾曼濾波及整體修正和局部修正法分別進(jìn)行求解,并分別求出與真值的差值。最后,對(duì)數(shù)據(jù)結(jié)果進(jìn)行比較,得出結(jié)論。

假設(shè)預(yù)測(cè)值準(zhǔn)確,通過MA TLAB程序解算得到表1中的數(shù)據(jù)。

表1 數(shù)據(jù)計(jì)算結(jié)果比較

通過表1中數(shù)據(jù)比較,可以得到如下結(jié)論:

1)整體和局部修正方法都可以有效地抵抗粗差的影響,獲得比較接近平差值的濾波結(jié)果。

2)整體修正法比局部修正法收斂較快。

3)局部修正結(jié)果比整體修正總體上結(jié)果要準(zhǔn)確。但是局部修正需要提前進(jìn)行粗差定位。經(jīng)檢驗(yàn),如果粗差定位時(shí)存在未被探測(cè)到的粗差觀測(cè),則濾波結(jié)果會(huì)變壞。但是如果定位粗差觀測(cè)值多于實(shí)際的粗差觀測(cè)值,濾波結(jié)果仍然較好。

5 結(jié)束語

本文提出卡爾曼濾波在進(jìn)行有粗差觀測(cè)值的濾波時(shí)發(fā)散的問題,基于此問題在標(biāo)準(zhǔn)卡爾曼濾波的基礎(chǔ)上,并在假定濾波模型正確后,將卡爾曼濾波粗差修正分為整體修正和局部修正兩方法,并給出各自的計(jì)算步驟。最后,將兩種方法應(yīng)用于工程實(shí)際數(shù)據(jù)的解算中,通過與標(biāo)準(zhǔn)卡爾曼濾波比較計(jì)算可以得出,兩種修正方法在含有粗差時(shí)的濾波中均能夠有效抵抗粗差觀測(cè)值的影響,獲得比較接近真值的濾波結(jié)果。但從結(jié)果分析可以發(fā)現(xiàn),兩種方法又有各自的優(yōu)缺點(diǎn)。計(jì)算結(jié)果表明,卡爾曼濾波粗差修正方法具有實(shí)際的應(yīng)用價(jià)值和研究?jī)r(jià)值。但如何快速準(zhǔn)確的粗差定位仍需要進(jìn)一步的研究。

[1]王秋平,陳 娟,王顯利,等,一種新的自適應(yīng)非線性卡爾曼濾波算法[J].光電工程,2008,35(7):17-21.

[2]朱建軍,曾卓喬.污染誤差模型下的測(cè)量數(shù)據(jù)處理理論[J].測(cè)繪學(xué)報(bào),1999(3):215-220.

[3]宋力杰,楊元喜.均值漂移模型粗差探測(cè)法與LEGE法的比較[J].測(cè)繪學(xué)報(bào),1999(4):295-300.

[4]許國輝,張新長.卡爾曼濾波模型粗差的探測(cè)及其在施工變形中的測(cè)量[J].中山大學(xué)學(xué)報(bào):自然科學(xué)版,2003,42(3):89-91.

[5]祝轉(zhuǎn)民,秋宏興,李濟(jì)生,等.動(dòng)態(tài)測(cè)量數(shù)據(jù)野值的辨識(shí)與剔除[J].系統(tǒng)工程與電子技術(shù),2004,26(2):147-149.

[6]劉大杰 ,于正林,陶本藻.形變測(cè)量動(dòng)態(tài)數(shù)據(jù)的處理方法[M].北京:測(cè)繪出版社,1992:185-193.

[7]HAMPEL.Robust Statistics App roach Based on Influence[M].John Wiley Press,1986.

Themethod of grosserror modification in Kalman Filtering and itsapplication

Q IGong-yu1,2,Q IU Wei-ning1,2,HUA Xiang-hong1,2
(1.School of Geodesy and Geomatics,Wuhan University,Wuhan 430079,China;2.Hazard Monitoring and Prevention Research Center,Wuhan University,W uhan 430079,China)

M ethods of integralmodification and part modification are raised to modify gross error based on the standard Kalman Filtering.In the end,the resultsof an engineering examp le indicate that themethods can both effectively control the effects of gross error and obtain the more accurate results.

Kalman Filtering;gross error modification;app lication

P207

A

1006-7949(2010)01-0050-03

2009-04-22

國家863計(jì)劃資助項(xiàng)目(2008AA 12Z308);國土資源大調(diào)查項(xiàng)目(1212010914015);武漢市科技計(jì)劃資助項(xiàng)目(200760323113)

齊公玉(1986-),女,碩士研究生.

[責(zé)任編輯李銘娜]