控制系統(tǒng)仿真課程設計--基于kalman濾波的信息融合算法設計

上傳人：奔*** IP屬地：河北更新時間：2024-03-01 格式：docx 頁數(shù)：23 大?。?84.25KB 人氣指數(shù)：12 舉報 版權申訴

控制系統(tǒng)仿真課程設計--基于kalman濾波的信息融合算法設計_第1頁

已閱讀1頁，還剩22頁未讀，繼續(xù)免費閱讀

版權說明：本文檔由用戶提供并上傳，收益歸屬內(nèi)容提供方，若內(nèi)容存在侵權，請進行舉報或認領

文檔簡介

1、　　控制系統(tǒng)仿真課程設計　?。?010級）　　控制系統(tǒng)仿真課程設計　　一、題目　　基于Kalman濾波的信息融合算法設計　　1) 學習并

2、掌握線性系統(tǒng)Kalman濾波的基本原理和基本公式；　　2) 學習并掌握一種常用的融合算法；　　3) 學習并利用Matlab軟件實現(xiàn)基本的Kalman濾波和信息融合算法的仿真。　　二、主要要求　　1) 具備基本的概率與數(shù)理統(tǒng)計知識；&l

3、t;/p>　　2) 熟悉并掌握基本的Matlab軟件編寫能力；　　3) 學習并掌握正交投影定理和矩陣求逆定理；　　4) 了解Kalman濾波的功能、來源和基本原理；　　5) 掌握Kalman濾波的推導過程和基本運行公式；　　6) 了解信息融合的

4、基本概念和方法；　　7) 掌握一種典型的多傳感器信息融合算法：分布式局部估計值加權融合。　　三、主要內(nèi)容　　一）線性系統(tǒng)的Kalman濾波　　考慮如下一類單傳感器線性動態(tài)估計系統(tǒng)<b&g

5、t;　　(1)　　(2)　　其中，是離散的時間變量；是系統(tǒng)的狀態(tài)向量，是系統(tǒng)的狀態(tài)轉(zhuǎn)移矩陣；是狀態(tài)的觀測向量，是相應的觀測矩陣；和是零均值的高斯白噪聲過程，且滿足如下條件：　　， (3)　　初始狀態(tài)

6、為一隨機向量，且滿足　　(4)　　那么，線性系統(tǒng)的Kalman濾波基本公式如下：　　計算狀態(tài)的一步預測值　　(5)　　計算一步預測誤差協(xié)方差陣</p

7、>　　(6)　　計算增益陣　　(7)　　計算狀態(tài)估計值　　(8)</

8、b>　　和估計誤差協(xié)方差陣　　(9)　　其中和為時刻的狀態(tài)估計以及相應的估計誤差協(xié)方差陣。　　那么，Kalman濾波仿真程序執(zhí)行方案如下：　　確定初

9、始狀態(tài)、初始狀態(tài)估計和相應的協(xié)方差矩陣；給定狀態(tài)轉(zhuǎn)移矩陣、過程噪聲方差、測量矩陣和測量噪聲方差（這些量均可認為是常量）　　產(chǎn)生仿真信號數(shù)據(jù)　　從開始循環(huán)（L為給定的仿真時刻長度）　　當時　　a1)

10、利用和隨機函數(shù)產(chǎn)生一個高斯白噪聲；　　a2) 根據(jù)式(1)有；　　a3) 利用和隨機函數(shù)產(chǎn)生一個高斯白噪聲；　　a4) 根據(jù)式(2)有。　　當時　　b1) 利用和隨機函數(shù)產(chǎn)生一個高斯白噪聲；

11、　　b2) 根據(jù)式(1)有；　　b3) 利用和隨機函數(shù)產(chǎn)生一個高斯白噪聲；　　b4) 根據(jù)式(2)有。　　開始Kalman濾波估計　　從開始循環(huán)（L為給定的仿真時刻長度）

12、　　當時　　a1) 根據(jù)式(5)和式(6)有　　, 　　a2) 利用式(7)-(9)計算估計和相應的估計誤差協(xié)方差矩陣。　　當時

13、b1) 根據(jù)式(5)和式(6)計算和；　　b2) 利用式(7)-(9)計算估計和相應的估計誤差協(xié)方差矩陣。　　問題：給定相應參數(shù)（也鼓勵采用其他參數(shù)），進行Kalman濾波估計算法程序的編寫，并進行繪圖和分析　　標量情形：<b&

14、gt;　　，，，，，，　　（1）請利用Matlab軟件進行Kalman濾波估計仿真程序編寫；　　%mykf.m　　%produce system　　clear;</p&g

15、t;　　A=1;　　P0=100;　　X0=10;　　C=1;　　Q=0.1;&

16、lt;/b>　　R=10;　　%real states and measure states　　for k=1:150　　W(k)=sqrt(Q)*randn(1,1);　　V(k)=s

17、qrt(R)*randn(1,1);　　if k==1　　X(k)=A*X0+W(k);　　Z(k)=C*X(k)+V(k);　　else　　X(k)=A

18、*X(k-1)+W(k);　　Z(k)=C*X(k)+V(k);　　end　　end　　%predict states and estimate states<p&g

19、t;　　for k=1:150　　if k==1　　X_yc(k)=A*X0;　　Z_yc(k)=C*X_yc(k);　　P_yc(k)=A*P0*A'+Q;　　K(k)=P_yc(k

20、)*C'/(C*P_yc(k)*C'+R);　　X_gj(k)=X_yc(k)+K(k)*(Z(k)-Z_yc(k));　　P_gj(k)=(eye(1)-K(k)*C)*P_yc(k);　　else　　X_yc(k)

21、=A*X(k-1);　　Z_yc(k)=C*X_yc(k);　　P_yc(k)=A*P_yc(k-1)*A'+Q;　　K(k)=P_yc(k)*C'/(C*P_yc(k)*C'+R);　　X_gj(k)=X_yc(k)+K(k)*(Z(k)-Z_y

22、c(k));　　P_gj(k)=(eye(1)-K(k)*C)*P_yc(k);　　end　　end　　%create figure

23、figure　　t=1:150;　　plot(t,Z(1,t),'-og')　　hold on　　plot(t,X_gj(1,t),'r')&l

24、t;/p>　　hold on　　plot(t,X(1,t),'b')　　hold off　　legend('觀測','估計','狀態(tài)')</

25、p>　　xlabel('仿真次數(shù)')　　ylabel('數(shù)值')　　figure　　plot(abs(Z-X),'-og');　　hold

26、on　　plot(abs(X_gj-X));　　hold off　　legend('預測與真實之差','估計與真實之差')　　xlabel('仿真次數(shù)')</p

27、>　　ylabel('數(shù)值')　?。?）繪出狀態(tài)預測值和狀態(tài)估計值的曲線圖；　　（3）繪出預測誤差協(xié)方差和估計誤差協(xié)方差的曲線圖；　?。?）對仿真結(jié)果進行分析。　　預測值和估計值都能夠在一定程度上反應真實值，但是估計值比觀測值更

28、接近真實值。　　狀態(tài)估計值：表明估計值是在預測值的基礎上進行優(yōu)化后得到結(jié)果，所以估計值更準確一些。　　矢量情形：　　，，，，　　，，&

29、lt;p>　?。?）請利用Matlab軟件進行Kalman濾波估計仿真程序編寫；　　%mykf.m　　%produce system　　clear;　　A=[1 1;0 1];</p&

30、gt;　　P0=[100 10;10 100];　　x0=[1;0.1];　　X_0=[10;1];　　C=[1 0;0 1];　　Q=[0.1 0;0 0.1];　　R=[10 0;0 10];</p

31、>　　%real states and measure states　　for k=1:150　　W(:,k)=sqrt(Q)*randn(2,1);　　V(:,k)=sqrt(R)*randn(2,1);　　if k==1&

32、lt;/b>　　X(:,k)=A*x0+W(:,k);　　Z(:,k)=C*X(:,k)+V(:,k);　　else 　　X(:,k)=A*X(:,k-1)+W(:,k);　　Z(:,k

33、)=C*X(:,k)+V(:,k);　　end　　end　　%predict states and estimate states　　for k=1:150

34、;　　if k==1　　X_yc(:,k)=A*X_0;　　Z_yc(:,k)=C*X_yc(:,k);　　P_yc(:,:,k)=A*P0*A'+Q;　　T_yc(k)=trace(P_yc(:,:,k));</p&g

35、t;　　K(:,:,k)=P_yc(:,:,k)*C'/(C*P_yc(:,:,k)*C'+R);　　X_gj(:,k)=X_yc(:,k)+K(:,:,k)*(Z(:,k)-Z_yc(:,k));　　P_gj(:,:,k)=(eye(2)-K(:,:,k)*C)*P_yc(:,:,k);&l

36、t;p>　　T_gj(k)=trace(P_gj(:,:,k));　　else　　X_yc(:,k)=A*X_gj(:,k-1);　　Z_yc(:,k)=C*X_yc(:,k);　　P_yc(:,:,k)=A*P_gj(:,:,

37、k-1)*A'+Q;　　T_yc(k)=trace(P_yc(:,:,k));　　K(:,:,k)=P_yc(:,:,k)*C'/(C*P_yc(:,:,k)*C'+R);　　X_gj(:,k)=X_yc(:,k)+K(:,:,k)*(Z(:,k)-Z_yc(:,k));

38、;　　P_gj(:,:,k)=(eye(2)-K(:,:,k)*C)*P_yc(:,:,k);　　T_gj(k)=trace(P_gj(:,:,k));　　end　　end

39、　%create figure　　figure 　　t=1:150;　　plot(t,X(1,t),'-or')　　hold on</p

40、>　　plot(t,X_gj(1,t),'g')　　plot(t,Z(1,t),'--')　　hold off　　legend('分量一狀態(tài)','分量一估計','分量一預測&#

41、39;)　　xlabel('仿真次數(shù)')　　ylabel('數(shù)值')　　figure 　　plot(t,X(2,t),'-or',t,X_gj(2,t),'g')<

42、/p>　　hold on　　plot(t,Z(2,t),'--')　　hold off　　legend('分量二狀態(tài)','分量二估計','分量二預測

43、9;)　　xlabel('仿真次數(shù)')　　ylabel('數(shù)值')　　figure　　plot(t,abs(Z(2,t)-X(2,t)),'-or')<p&

44、gt;　　hold on　　plot(t,abs(X_gj(2,t)-X(2,t)),'g')　　hold off　　legend('預測與真實之差','估計與真實之差')

45、;　　xlabel('仿真次數(shù)')　　ylabel('數(shù)值')　　figure　　plot(t,T_gj(t),'g',t,T_yc(t),'-or')<p&g

46、t;　　legend('估計','預測')　　xlabel('仿真次數(shù)')　　ylabel('數(shù)值')　?。?）繪出狀態(tài)預測值和狀態(tài)估計值的曲線圖（每個狀態(tài)包括兩個分量）；　　圖1.1<

47、;/b>　　圖2.1　　（3）繪出預測誤差協(xié)方差陣跡(Trace)和估計誤差協(xié)方差陣跡的曲線圖；　　圖3.1　　（4）對仿真結(jié)果進行分析。　　分量的

48、估計值比分量的觀測值更接近真實值。整個時也是估計值更準確。　　針對矢量情形，自行選取三組不同的參數(shù)進行Kalman濾波的仿真，并進行相應仿真結(jié)果的比較分析。　　改變Q變大（Q=4）　　改變R變小（R=4）　　改變H變?。℉=0.99）<p&g

49、t;　　當R的值變小時，預測值的陣跡會變得下墜更快，預測值本身的震蕩會減小，對真實值的偏離會變小。當Q的值增大時，估計值也會更加偏離真實值。當H變小時，預測值與真實值偏差變大，估計值與真實值的偏差也會變大。　　二）基于線性Kalman濾波信息融合算法　　考慮如下一類多傳感器線性動態(tài)估計系統(tǒng)　　(

50、10)　　， (11)　　其中，是離散的時間變量，為傳感器的數(shù)目；是系統(tǒng)的狀態(tài)向量，是系統(tǒng)的狀態(tài)轉(zhuǎn)移矩陣；是狀態(tài)的觀測向量，是相應的觀測矩陣；和是零均值的高斯白噪聲過程，且滿足如下條件：　　， (12)

51、　初始狀態(tài)為一隨機向量，且滿足　　(13)　　那么，對于每一個傳感器觀測均可執(zhí)行一)當中基于單個觀測的Kalman濾波估計，可得到個局部估計和相應的估計誤差協(xié)方差矩陣。從而，可利用分布式加權融合技術將上述個局部Kalman濾波估計進行融合，即：　　(14)&

52、lt;/b>　　此時，和為融合后的狀態(tài)估計和相應的融合估計誤差協(xié)方差矩陣。　　問題：給定相應參數(shù)（也鼓勵采用其他參數(shù)），進行上述分布式融合算法的仿真　　給定如下參數(shù)：　　，，，，</p

53、>　　，，　　（1）請利用Matlab軟件進行分布式融合估計算法仿真程序編寫；　　%mykf.m　　%produce system　　clear;&l

54、t;/b>　　clc;　　A=[1 1;0 1];　　P0=[100 10;10 100];　　x0=[1;0.1];　　X_0=[10;1];

55、　　C=[1 0;0 1];　　Q=[0.1 0;0 0.1];　　R=[4 0;0 4];　　%guan ce qi 1 real states and measure states　　for k=1:150　　W(:,k)=s

56、qrt(Q)*randn(2,1);　　V1(:,k)=sqrt(R)*randn(2,1);　　if k==1　　X1(:,k)=A*x0+W(:,k);　　Z1(:,k)=C*X1(:,k)+V1(:,k);

57、　　else 　　X1(:,k)=A*X1(:,k-1)+W(:,k);　　Z1(:,k)=C*X1(:,k)+V1(:,k);　　end　　end</

58、b>　　%predict states and estimate states　　for k=1:150　　if k==1　　X_yc1(:,k)=A*X_0;　　Z_yc1(:,k)=

59、C*X_yc1(:,k);　　P_yc1(:,:,k)=A*P0*A'+Q;　　T_yc1(k)=trace(P_yc1(:,:,k));　　K1(:,:,k)=P_yc1(:,:,k)*C'/(C*P_yc1(:,:,k)*C'+R);　　X_g

60、j1(:,k)=X_yc1(:,k)+K1(:,:,k)*(Z1(:,k)-Z_yc1(:,k));　　P_gj1(:,:,k)=(eye(2)-K1(:,:,k)*C)*P_yc1(:,:,k);　　T_gj1(k)=trace(P_gj1(:,:,k));　　else&l

61、t;/p>　　X_yc1(:,k)=A*X_gj1(:,k-1);　　Z_yc1(:,k)=C*X_yc1(:,k);　　P_yc1(:,:,k)=A*P_gj1(:,:,k-1)*A'+Q;　　T_yc1(k)=trace(P_yc1(:,:,k));

62、　　K1(:,:,k)=P_yc1(:,:,k)*C'/(C*P_yc1(:,:,k)*C'+R);　　X_gj1(:,k)=X_yc1(:,k)+K1(:,:,k)*(Z1(:,k)-Z_yc1(:,k));　　P_gj1(:,:,k)=(eye(2)-K1(:,:,k)*C)*P_yc1(:,:,k);</p&

63、gt;　　T_gj1(k)=trace(P_gj1(:,:,k));　　end　　end　　%guan ce qi 2 real states and measure states<p&

64、gt;　　for k=1:150　　V2(:,k)=sqrt(R)*randn(2,1);　　if k==1　　X2(:,k)=A*x0+W(:,k);　　Z2(:,k)=C*X2(:,k)+V2(:,k);&

65、lt;p>　　else 　　X2(:,k)=A*X2(:,k-1)+W(:,k);　　Z2(:,k)=C*X2(:,k)+V2(:,k);　　end　　end</b&

66、gt;　　%predict states and estimate states　　for k=1:150　　if k==1　　X_yc2(:,k)=A*X_0;　　Z_yc2(:,k)=C*

67、X_yc2(:,k);　　P_yc2(:,:,k)=A*P0*A'+Q;　　T_yc2(k)=trace(P_yc2(:,:,k));　　K2(:,:,k)=P_yc2(:,:,k)*C'/(C*P_yc2(:,:,k)*C'+R);　　X_gj2

68、(:,k)=X_yc2(:,k)+K2(:,:,k)*(Z2(:,k)-Z_yc2(:,k));　　P_gj2(:,:,k)=(eye(2)-K2(:,:,k)*C)*P_yc2(:,:,k);　　T_gj2(k)=trace(P_gj2(:,:,k));　　else<

69、/p>　　X_yc2(:,k)=A*X_gj2(:,k-1);　　Z_yc2(:,k)=C*X_yc2(:,k);　　P_yc2(:,:,k)=A*P_gj2(:,:,k-1)*A'+Q;　　T_yc2(k)=trace(P_yc2(:,:,k));&l

70、t;p>　　K2(:,:,k)=P_yc2(:,:,k)*C'/(C*P_yc2(:,:,k)*C'+R);　　X_gj2(:,k)=X_yc2(:,k)+K2(:,:,k)*(Z2(:,k)-Z_yc2(:,k));　　P_gj2(:,:,k)=(eye(2)-K2(:,:,k)*C)*P_yc2(:,:,k);

71、;　　T_gj2(k)=trace(P_gj2(:,:,k));　　end　　end　　% rong he 　　for k=1:150

72、;　　P_rh(:,:,k)=inv(P_gj2(:,:,k))+inv(P_gj1(:,:,k));　　P_rhgj(:,:,k)=inv(P_rh(:,:,k));　　T_rhgj(k)=trace(P_rhgj(:,:,k));　　X_rhgj(:,k)=P_rhgj(:,:,k)*inv(P_gj2(:,

73、:,k))*X_gj2(:,k)+P_rhgj(:,:,k)*inv(P_gj1(:,:,k))*X_gj1(:,k);　　end　　%create figure　　figure

74、　　t=1:150;　　plot(t,X2(1,t),'r',t,Z2(1,t),'g')　　hold on　　plot(t,X_gj2(1,t),'-ok',t,X_rhgj(1,t),'--

75、9;)　　hold off　　legend(' 觀測器2狀態(tài)一',',預測分量一',' 估計分量一),',' 融合分量一')　　xlabel('仿真次數(shù)')

76、;　　ylabel('數(shù)值')　　figure　　plot(t,X2(2,t),'r',t,Z2(2,t),'g')　　hold on　　plot(

77、t,X_gj2(2,t),'-ok',t,X_rhgj(2,t),'--')　　hold off　　legend(' 觀測器2狀態(tài)二',',預測分量二', '估計分量二,',' 融合分量二')<

78、;p>　　xlabel('仿真次數(shù)')　　ylabel('數(shù)值')　　figure　　plot(t,abs(X_rhgj(1,t)-X2(1,t)),'--',t,abs(X_gj2(1,t)-X2(1,t)),'

79、k')　　legend('融合與真實之差','估計與真實之差')　　xlabel('仿真次數(shù)')　　ylabel('數(shù)值')　　figure

80、　　k=1:150;　　plot(T_gj1(1,k),'g');　　hold on　　plot(T_gj2(1,k),'r');<b&g

81、t;　　hold on　　plot(T_rhgj(1,k),'-oy');　　hold off　　legend('估計一陣跡','估計二陣跡','估計融合陣跡')<

82、;p>　　xlabel('仿真次數(shù)')　　ylabel('數(shù)值')　?。?）繪出融合狀態(tài)估計和每個局部估計狀態(tài)估計值的曲線圖　?。?）繪出融合估計誤差協(xié)方差矩陣跡和每個局部估計誤差協(xié)方差陣跡曲線圖；　?。?）對上述仿真結(jié)果進行分析。&

83、lt;/p>　　融合之后的估計值與真實值的誤差比兩個分量都要小。這是因為通過多個設備對同一狀態(tài)的加權估計，得到的結(jié)果可以進一步減小與真實值的誤差。　　四實踐總結(jié)　　經(jīng)過這兩個星期的短學期，從中學到了不少以前沒有深入過的知識，對Matlab的操作更加熟悉，對kalman濾波有了一定的了解，讓

84、我獲益匪淺，拓寬了我的知識面。 　　五致謝　　非常感謝xx老師以及同學們在我上機過程為我解答疑問，對我的錯誤提出指正。　　六參考文獻　　付夢印.kalman濾波理論及其在導航系統(tǒng)中的應用.科學

眾賞文庫> 全部分類> 畢業(yè)設計

溫馨提示

1. 本站所有資源如無特殊說明，都需要本地電腦安裝OFFICE2007和PDF閱讀器。圖紙軟件為CAD,CAXA,PROE,UG,SolidWorks等.壓縮文件請下載最新的WinRAR軟件解壓。
2. 本站的文檔不包含任何第三方提供的附件圖紙等，如果需要附件，請聯(lián)系上傳者。文件的所有權益歸上傳用戶所有。
3. 本站RAR壓縮包中若帶圖紙，網(wǎng)頁內(nèi)容里面會有圖紙預覽，若沒有圖紙預覽就沒有圖紙。
4. 未經(jīng)權益所有人同意不得將文件中的內(nèi)容挪作商業(yè)或盈利用途。
5. 眾賞文庫僅提供信息存儲空間，僅對用戶上傳內(nèi)容的表現(xiàn)方式做保護處理，對用戶上傳分享的文檔內(nèi)容本身不做任何修改或編輯，并不能對任何下載內(nèi)容負責。
6. 下載文件中如有侵權或不適當內(nèi)容，請與我們聯(lián)系，我們立即糾正。
7. 本站不保證下載資源的準確性、安全性和完整性, 同時也不承擔用戶因使用這些下載資源對自己和他人造成任何形式的傷害或損失。

控制系統(tǒng)仿真課程設計--基于kalman濾波的信息融合算法設計

文檔簡介

溫馨提示

最新文檔

評論

控制系統(tǒng)仿真課程設計--基于kalman濾波的信息融合算法設計

文檔簡介

溫馨提示

最新文檔

評論

免費下載