離散卡爾曼濾波用于GPS動(dòng)態(tài)變形數(shù)據(jù)處理

1、離散卡爾曼濾波用于GPS動(dòng)態(tài)變形數(shù)據(jù)處理馬攀1,文鴻雁1,2(1.桂林工學(xué)院土木工程系,廣西桂林541004;2.武漢大學(xué)測(cè)繪學(xué)院,湖北武漢430079)摘要:研究了一種卡爾曼濾波算法,1998動(dòng)態(tài)變形數(shù)據(jù),采用自編軟件(用C語(yǔ)言)進(jìn)行處理2,0很相近,并,GPS動(dòng)態(tài)變形的數(shù)據(jù)精度的效果比:卡爾曼濾波;GPS;高斯-馬爾柯夫序列中圖分類(lèi)號(hào):O211164;P22814文獻(xiàn)標(biāo)識(shí)碼:A0引言目前,在變形監(jiān)測(cè)中,由于GPS、測(cè)地機(jī)器人等新技術(shù)不斷引用,積累資料多,數(shù)據(jù)量大,如何及時(shí)地有效地從大量的變形監(jiān)測(cè)信息中進(jìn)行數(shù)據(jù)挖掘,即抽取、轉(zhuǎn)換、分析和進(jìn)行建模處理,從中提取關(guān)鍵性的數(shù)據(jù),對(duì)變形進(jìn)行分析與解釋

2、,尤其是對(duì)建筑物的變形進(jìn)行預(yù)測(cè),使人們能及時(shí)地針對(duì)變形建筑物變形狀態(tài)進(jìn)行有效處理,對(duì)滑坡等災(zāi)害及時(shí)做好有關(guān)方面的防范和預(yù)防措施,從而極大地減少因建筑物變形或滑坡等地質(zhì)災(zāi)害造成的損失.所以對(duì)變形監(jiān)測(cè)信息進(jìn)行數(shù)據(jù)挖掘,特別是對(duì)大型建筑物的變形進(jìn)行分析與預(yù)測(cè)十分重要.而在GPS用于動(dòng)態(tài)變形觀測(cè)中,是僅在確定的瞬間來(lái)研究系統(tǒng)的性能,而在此瞬間所得的觀測(cè)值本身是帶有隨機(jī)噪聲的.如果能在這個(gè)系統(tǒng)中分別找出表示它們的狀態(tài)參數(shù)與其觀測(cè)值之間的函數(shù)關(guān)系,那么就可以利用卡爾曼濾波來(lái)減弱隨機(jī)噪聲的干擾,達(dá)到提高GPS變形觀測(cè)數(shù)據(jù)精度的目的.為此,本文將對(duì)卡爾曼濾波用于某大壩GPS動(dòng)態(tài)形變數(shù)據(jù)處理作一些探討.收稿日期

3、:2001-12-24;修訂日期:2002-04-03基金項(xiàng)目:廣西壯族自治區(qū)教育廳資助項(xiàng)目(1998-169)作者簡(jiǎn)介:馬攀(1977-),男,四川南部人,碩士研究生,大地測(cè)量學(xué)與測(cè)量工程專(zhuān)業(yè).第22卷第3期馬攀等:離散卡爾曼濾波用于GPS動(dòng)態(tài)變形數(shù)據(jù)處理2351離散線性系統(tǒng)的卡爾曼濾波111離散線性系統(tǒng)卡爾曼濾波的數(shù)學(xué)描述當(dāng)不考慮系統(tǒng)具有確定性輸入時(shí),可設(shè)系統(tǒng)狀態(tài)方程和量測(cè)方程分別為2),Xk=Fk/k-1Xk-1+k-1Wk-1(k=1,2),Yk=HkXk+Ek(k=1,2(1)(2)式中Xk為k時(shí)刻的n維狀態(tài)矢量,也是被估計(jì)矢量;Yk為k時(shí)刻上m維量測(cè)矢量;Fk/k-1為k-1到k時(shí)

4、刻的系統(tǒng)一步轉(zhuǎn)移矩陣(nn階);Wk-1為k-1時(shí)刻的系統(tǒng)噪聲(r維);k-1為系統(tǒng)噪聲矩陣(nr階),它表征由k-1到k時(shí)刻的各個(gè)系統(tǒng)噪聲分別影響k時(shí)刻各個(gè)狀態(tài)的程度;Hk為k時(shí)刻的量測(cè)矩陣(mn階);Ek為k時(shí)刻的m維量測(cè)噪聲.卡爾曼濾波要求Wk和Ek的白噪聲序列,有E(Wk)=0,ECov(Wk,Wj)EkEj=k,(kEj)=0,(3)(4)式中Cov(A,k,在卡;kj是Kronecker函數(shù),即kj=初始狀態(tài)的一、二階統(tǒng)計(jì)特性為E(X0)=mX0,Var(X0)=CX0,Ek都不相關(guān).(6)0(kj),1(k=j).(5)式中Var(X0)為對(duì)X0求方差的符號(hào).卡爾曼濾波要求mX0

5、和CX0為已知量,且要求X0與Wk以及為便于敘述,現(xiàn)將離散卡爾曼濾波的其它方程總結(jié)如下?tīng)顟B(tài)預(yù)測(cè)Xk/k-1=Fk/k-1Xk-1;T方差預(yù)測(cè)Dk/k-1=Fk/k-1Dk-1FTk/k-1+k/k-1Qk-1k/k-1;(7)(8)(9)(10)-1狀態(tài)估計(jì)Xk=Xk/k-1+Kk(Yk-HkXk/k-1);方差迭代Dk=(I-KkHk)Dk/k-1;T濾波增益Kk=Dk/k-1HTk(HkDk/k-1Hk+Rk);(11)(12)初始值X0=E(X0),D0=Var(X0).Yk-1求Xk的預(yù)測(cè)值的公式,稱(chēng)此式為一步預(yù)測(cè)公式;而稱(chēng)Xk/HkXk/k-1)k-1(7)和(8)式是利用Xk-1求

6、Xk的預(yù)報(bào)值Xk/k-1及其方差公式,也就是利用Y1,Y2,為一步預(yù)測(cè)值.(9)式中的(Yk-稱(chēng)為預(yù)報(bào)殘差,(9)式的物理意義是:濾波值等于預(yù)報(bào)值加一修正項(xiàng),該修正項(xiàng)由預(yù)報(bào)殘差乘增益矩陣構(gòu)成.因此,通常稱(chēng)Kk為濾波增益矩陣.112卡爾曼濾波的性質(zhì)和特點(diǎn)從111可以看到,卡爾曼濾波方程是一組遞推計(jì)算公式,其計(jì)算過(guò)程是一個(gè)不斷地預(yù)測(cè)、修正的過(guò)程.在求解時(shí)不需要貯存大量的觀測(cè)數(shù)據(jù),并且當(dāng)?shù)玫叫碌挠^測(cè)數(shù)據(jù)時(shí),可隨時(shí)算得新的濾波值,便于實(shí)時(shí)處理觀測(cè)成果.卡爾曼濾波還具有以下特點(diǎn):卡爾曼濾波遞推公式是根據(jù)廣義最小二乘原理得到的,按廣義最小二乘原理得到的結(jié)果具有無(wú)偏性和方差最小性.因此,對(duì)于高斯馬爾柯夫序列

7、模型來(lái)說(shuō),卡爾曼濾波值X(k/k)是由觀測(cè)值Y1,Y2,Yk得到Xk的最小方差無(wú)偏估計(jì),其誤差方差陣DX(k/k)就是所有估計(jì)中的最小方差陣.對(duì)于一個(gè)馬爾柯夫序列模型來(lái)說(shuō),X(k/k)則是Xk的線性最小方差無(wú)偏估計(jì)2.對(duì)于具有一致完全能觀性和一致完全能控性的線性系統(tǒng),濾波遞推充分多的步數(shù)后,濾波誤差236桂林工學(xué)院學(xué)報(bào)2002年方差陣Dk/k將不依賴(lài)于D0,濾波值也漸漸地不依賴(lài)于初始值的選取3.從卡爾曼濾波方程可以看出,增益矩陣Kk與觀測(cè)值無(wú)關(guān),因此,可以預(yù)先算出,從而減少實(shí)時(shí)處理的計(jì)算工作量.從(11)式可知:當(dāng)Rk增大時(shí),增益矩陣Kk變小,這是因?yàn)橛^測(cè)噪聲增大,則濾波的增益就應(yīng)取小一些,以

8、減弱觀測(cè)噪聲對(duì)濾波值的影響.在應(yīng)用卡爾曼濾波中,除了盡可能精確地描述動(dòng)態(tài)方程和觀測(cè)方程外,還要選取適當(dāng)?shù)腝k和Rk,這個(gè)可以通過(guò)試算獲得,也可以自適應(yīng)估計(jì)之.預(yù)報(bào)殘差(Yk-HkXk/k-1)是一個(gè)零均值的高斯馬爾柯夫序列2.卡爾曼濾波是將全部待估參數(shù)都作為正態(tài)隨機(jī)量,按照廣義最小二乘原理來(lái)求定參數(shù)的最佳估計(jì)值.當(dāng)已知參數(shù)的先驗(yàn)統(tǒng)計(jì)性質(zhì)時(shí),由于卡爾曼濾波考慮了這種性質(zhì),因此,二乘平差估值具有更高的精度.2卡爾曼濾波程序設(shè)計(jì)211從13Kk=Dk/Xk=Fk/k-1Hk(HkDk/Tk-1Hk+Rk)k-1T-1,(13a)(13b)(13c)k-1.k-1Xk-1+Kk(Yk-HkFk/KkH

9、k)Dk/k-1Xk-1),Ck=(I-Dk/k-1k-1k-1,Qk-1Tk/=Fk/k-1Ck-1FTk/+k/(13d)上式中Rk為mm階的觀測(cè)噪聲Ek的協(xié)方差陣;Xk為n維向量,第k時(shí)刻經(jīng)濾波后的估值;Kk為nm階的增益矩陣;Ck為nn階的估計(jì)誤差協(xié)方差陣;Qk為nn階的模型噪聲Wk的協(xié)方差陣.根據(jù)上述公式,可以從X0=E(X0),D0(給定)出發(fā),利用已知的矩陣Qk,Rk,Hk,Fk/k-1以).及k時(shí)刻的觀測(cè)值Yk,遞推地算出每個(gè)時(shí)刻的狀態(tài)估計(jì)Xk(k=1,2,如果線性系統(tǒng)是定常的,則有Fk/k-1=F,Hk=H,即它們都是常陣;如果模型噪聲Wk和觀測(cè)噪聲Ek都是平穩(wěn)隨機(jī)序列,則Q

10、k和Rk都是常陣.在這種情況下,常增益的離散卡爾曼濾波是漸進(jìn)穩(wěn)定的.212程序流程圖程序設(shè)計(jì)流程圖見(jiàn)圖1.3算例及數(shù)據(jù)分析311算例參數(shù)在本次大壩GPS動(dòng)態(tài)形變測(cè)量數(shù)據(jù)處理的算例中,是對(duì)我國(guó)1998年洪水期某大壩,用雙頻GPS定位儀器進(jìn)行垂直于大壩軸線方向位移觀測(cè)值的“濾波”.觀測(cè)歷元間隔為15s,每隔6h計(jì)算一次,共計(jì)65d的觀測(cè)數(shù)據(jù).)T,并取其初值為設(shè)觀測(cè)點(diǎn)的狀態(tài)向量為Xk=(X,XX0=(113,0)T.狀態(tài)轉(zhuǎn)移矩陣為6F=10tk1,其中tk為采樣間隔,在本例中取tk=6h,即第22卷第3期馬攀等:離散卡爾曼濾波用于GPS動(dòng)態(tài)變形數(shù)據(jù)處理237F=106.量測(cè)矩陣為H=(1,0).動(dòng)

11、態(tài)系統(tǒng)維數(shù)為n=2,觀測(cè)系統(tǒng)維數(shù)為m=1.初始估計(jì)誤差協(xié)方差陣取為D0=0000.模型噪聲協(xié)方差陣取為Q=1.201.取觀測(cè)向量序列長(zhǎng)度k=260.對(duì)輸出的狀態(tài)向量序列每5個(gè)打印一次.312結(jié)果分析與討論,的方差DX.由于濾波的,所以可以在濾波中求得觀測(cè)值Y的改正值V(即-的估值),以及虛擬觀測(cè)值Yx的改正值Vx(即-x的估值),2并利用V和值Vx來(lái)求得驗(yàn)后的單位權(quán)方差0,這里只考慮Dx=0和信號(hào)X的數(shù)學(xué)期望x=0(即Yx=0時(shí))的情況.2根據(jù)濾波中求驗(yàn)后單位權(quán)方差的公式20=/n,這里的=-YTPV,于是將312的計(jì)算數(shù)據(jù)代2入上式,可得算例的驗(yàn)后單位權(quán)方差0=/n=30618/260=11

12、18.又根據(jù)求中誤差估值的基本公式=/n,可得算例的觀測(cè)值中誤差1130測(cè)=mm,濾波值中誤差1107mm.濾=圖1卡爾曼濾波程序流程圖Fig11ThesequencechartofKalmanfilteringprogram由圖2不難看出,該算例成果圖形象直觀地顯示出濾波值數(shù)據(jù)曲線與原始數(shù)據(jù)曲線的變化趨勢(shì)基本一致;同時(shí),比較上面的計(jì)算數(shù)據(jù)可知,算例的驗(yàn)22后單位權(quán)方差0=1118與驗(yàn)前所取的方差0=112相近,說(shuō)明驗(yàn)前選取的方差是恰當(dāng)?shù)?并且,濾波值比觀測(cè)值的精度提高了約1717%.至此,也就是說(shuō)明該卡爾曼濾波數(shù)學(xué)模型較好地模擬了目標(biāo)系統(tǒng)的物理變化規(guī)律,表明所運(yùn)用的卡爾曼濾波程序的運(yùn)算結(jié)果比

13、較理想.4結(jié)論(1)通過(guò)對(duì)大壩GPS動(dòng)態(tài)形變測(cè)量數(shù)據(jù)的算例及其數(shù)據(jù)分析可知,運(yùn)用卡爾曼濾波方法可以明顯地改善GPS動(dòng)態(tài)變形的數(shù)據(jù)精度,并能大量地(待估狀態(tài)向量可以多維)實(shí)時(shí)處理GPS動(dòng)態(tài)變形數(shù)據(jù).由于適合在微機(jī)或小型計(jì)算機(jī)上應(yīng)用,因此其處理速度是其它方法難以比擬的.(2)在卡爾曼濾波程序用于算例的試算中,對(duì)于選取不同的模型噪聲和觀測(cè)噪聲協(xié)方差陣,獲得的各濾波值也略有差異.這說(shuō)明合理確定一個(gè)待估狀態(tài)系統(tǒng)的數(shù)學(xué)模型,選擇適當(dāng)?shù)慕?shù),更好238桂林工學(xué)院學(xué)報(bào)2002年圖2大壩于壩軸線垂直方向(Fig12Theachievementofdaminspoolofassumed)1,則卡爾曼濾波的效果更佳.(3)對(duì)于穩(wěn)定的系統(tǒng)來(lái)說(shuō),只要初始值足夠精確,就能保證以后的濾波值達(dá)到任意指定的精度.在實(shí)際應(yīng)用中僅僅只有穩(wěn)定性是不夠的,漸進(jìn)穩(wěn)定的系統(tǒng)則具有更優(yōu)的性質(zhì),即不論初始值有多大的偏差,只要濾波時(shí)間充分長(zhǎng),就能自行消除初始值偏差的影響,保證濾波值的準(zhǔn)確性.(4)由于在實(shí)踐中,卡爾曼濾波系統(tǒng)存在建模誤差及計(jì)算機(jī)有效字長(zhǎng)的限制等原因,可能導(dǎo)致濾波發(fā)散.然而由于卡爾曼濾波的廣泛應(yīng)用,克服濾波發(fā)散將是一個(gè)重要的研究方向.參考文獻(xiàn):12345ApplicationofdiscreteKalmanfilteringtoGPSdataprocessingofkinematicde