❶ 擴展卡爾曼濾波(EKF)演算法詳細推導及模擬(Matlab)
姓名:王柯禕
學號:20021110373T
轉自 :https://blog.csdn.net/gangdanerya/article/details/105105611
【嵌牛導讀】介紹擴展卡爾曼濾波(EKF)演算法的詳細推導,局限性和MATLAB模擬。
【嵌牛鼻子】擴展卡爾曼濾波(EKF)
【嵌牛正文】
擴展卡爾曼濾波演算法 是解決非線性狀態估計問題最為直接的一種處理方法,盡管EKF不是最精確的」最優「濾波器,但在過去的幾十年成功地應用到許多非線性系統中。所以在學習非線性濾波問題時應該先從EKF開始。
EKF演算法是將非線性函數進行泰勒展開,然後省略高階項,保留展開項的一階項,以此來實現非線性函數線性化,最後通過卡爾曼濾波演算法近似計算系統的狀態估計值和方差估計值。
一、EKF演算法詳細推導
【注】EKF推導參考的是黃蔚的博士論文「CKF及魯棒濾波在飛行器姿態估計中的應用研究」,論文中EKF,UKF和CKF等演算法講解的都很詳細,值得一看。
我們把KF與EKF演算法拿出來對比可以發現:
二、EKF演算法局限性:
該演算法線性化會引入階段誤差從而導致濾波精度下降,同時當初始狀態誤差較大或系統模型非線性程度較高時,濾波精度會受到嚴重影響甚至發散。
需要計算雅克比矩陣,復雜,計算量大,影響系統的實時性,還會導致EKF演算法的數值穩定性差。
當系統存在模型失配,量測干擾,量測丟失,量測延遲或狀態突變等復雜情況時,EKF演算法魯棒性差。
三、Matlab模擬:
clear all;clc; close all;
tf = 50;
Q = 10;w=sqrt(Q)*randn(1,tf);
R = 1;v=sqrt(R)*randn(1,tf);
P =eye(1);
x=zeros(1,tf);
Xnew=zeros(1,tf);
x(1,1)=0.1;
Xnew(1,1)=x(1,1);
z=zeros(1,tf);
z(1)=x(1,1)^2/20+v(1);
zjian=zeros(1,tf);
zjian(1,1)=z(1);
for k = 2 : tf
%%%%%%%%%%%%%%%模擬系統%%%%%%%%%%%%%%%
x(:,k) = 0.5 * x(:,k-1) + (2.5 * x(:,k-1) / (1 + x(:,k-1).^2)) + 8 * cos(1.2*(k-1)) + w(k-1);
z(k) = x(:,k).^2 / 20 + v(k);
%%%%%%%%%%%%%%%EKF開始%%%%%%%%%%%%%%%
Xpre = 0.5*Xnew(:,k-1)+ 2.5*Xnew(:,k-1)/(1+Xnew(:,k-1).^2) + 8 * cos(1.2*(k-1));
zjian =Xpre.^2/20;
F = 0.5 + 2.5 * (1-Xnew.^2)/((1+Xnew.^2).^2);
H = Xpre/10;
PP=F*P*F'+Q;
Kk=PP*H'*inv(H*PP*H'+R);
Xnew(k)=Xpre+Kk*(z(k)-zjian);
P=PP-Kk*H*PP;
end
t = 2 : tf;
figure; plot(t,x(1,t),'b',t,Xnew(1,t),'r*'); legend('真實值','EKF估計值');
模擬結果:
❷ 卡爾曼濾波器的演算法
在這一部分,我們就來描述源於Dr Kalman 的卡爾曼濾波器。下面的描述,會涉及一些基本的概念知識,包括概率(Probability),隨機變數(Random Variable),高斯或正態分配(Gaussian Distribution)還有State-space Model等等。但對於卡爾曼濾波器的詳細證明,這里不能一一描述。首先,我們先要引入一個離散控制過程的系統。該系統可用一個線性隨機微分方程(Linear Stochastic Difference equation)來描述:X(k)=A X(k-1)+B U(k)+W(k)再加上系統的測量值:Z(k)=H X(k)+V(k)上兩式子中,X(k)是k時刻的系統狀態,U(k)是k時刻對系統的控制量。A和B是系統參數,對於多模型系統,他們為矩陣。Z(k)是k時刻的測量值,H是測量系統的參數,對於多測量系統,H為矩陣。W(k)和V(k)分別表示過程和測量的雜訊。他們被假設成高斯白雜訊(White Gaussian Noise),他們的covariance 分別是Q,R(這里我們假設他們不隨系統狀態變化而變化)。對於滿足上面的條件(線性隨機微分系統,過程和測量都是高斯白雜訊),卡爾曼濾波器是最優的信息處理器。下面我們來用他們結合他們的covariances 來估算系統的最優化輸出(類似上一節那個溫度的例子)。首先我們要利用系統的過程模型,來預測下一狀態的系統。假設現在的系統狀態是k,根據系統的模型,可以基於系統的上一狀態而預測出現在狀態:X(k|k-1)=A X(k-1|k-1)+B U(k) ……….. (1)式(1)中,X(k|k-1)是利用上一狀態預測的結果,X(k-1|k-1)是上一狀態最優的結果,U(k)為現在狀態的控制量,如果沒有控制量,它可以為0。到現在為止,我們的系統結果已經更新了,可是,對應於X(k|k-1)的covariance還沒更新。我們用P表示covariance:P(k|k-1)=A P(k-1|k-1) A』+Q ……… (2)式(2)中,P(k|k-1)是X(k|k-1)對應的covariance,P(k-1|k-1)是X(k-1|k-1)對應的covariance,A』表示A的轉置矩陣,Q是系統過程的covariance。式子1,2就是卡爾曼濾波器5個公式當中的前兩個,也就是對系統的預測。現在我們有了現在狀態的預測結果,然後我們再收集現在狀態的測量值。結合預測值和測量值,我們可以得到現在狀態(k)的最優化估算值X(k|k):X(k|k)= X(k|k-1)+Kg(k) (Z(k)-H X(k|k-1)) ……… (3)其中Kg為卡爾曼增益(Kalman Gain):Kg(k)= P(k|k-1) H』 / (H P(k|k-1) H』 + R) ……… (4)到現在為止,我們已經得到了k狀態下最優的估算值X(k|k)。但是為了要令卡爾曼濾波器不斷的運行下去直到系統過程結束,我們還要更新k狀態下X(k|k)的covariance:P(k|k)=(I-Kg(k) H)P(k|k-1) ……… (5)其中I 為1的矩陣,對於單模型單測量,I=1。當系統進入k+1狀態時,P(k|k)就是式子(2)的P(k-1|k-1)。這樣,演算法就可以自回歸的運算下去。卡爾曼濾波器的原理基本描述了,式子1,2,3,4和5就是他的5 個基本公式。根據這5個公式,可以很容易的實現計算機的程序。