授之以渔： 卡尔曼滤波器 ...大泻蜜 ...

594961689

顶！！！

7998

这个，，，看起来挺费劲的

zouyf12

卡尔曼原来是需要建模的啊，  这个很麻烦啊

hanahimiEX

多谢啦

hanahimiEX

谢谢啦

我是传奇

matlab是个好东西，建议兄弟们好好学习下

晨景的风

好吧，看不懂

aaaalook

mark一下，哈哈，很高深啊。
http://bbs.21ic.com/icview-292853-1-1.html

aaaalook

有木有？
第一个例子中的代码

1. Ts = 1;                    采样时间
   
2. s = tf('s');
   
3. sysc = 100/(60*s +1);      真实系统的传递函数
   
4. sysd = c2d(sysc, Ts);
   
5. 
   
6. t = 0:TsTs*300);
   
7. u = ones(1, length(t));     系统输入
   
8. ys = lsim(sysd, u, t, 0);   系统输出
   
9. ys = ys + 10*rand(size(ys));  测量量 = 系统输出 + 噪声
   
10. 
    
11. A = [1, 0; 0, 1];
    
12. B = [1; 0];
    
13. C = [1, 0];
    
14. N = length(A);          系统维数
    
15. Len = length(u);
    
16. 
    
17. yout = zeros(Len, 1);   滤波器输出
    
18. Xh = zeros(N, 1);       状态变量
    
19. P = eye(N);              
    
20. Q = 2 * eye(N);         系统噪声
    
21. R = 50;                 测量噪声
    
22. for i = 1 : Len
    
23.     Xh = A*Xh + B*u(i);
    
24.     P = A*P*A' + Q;
    
25.     K = P*C'* inv(C*P*C' + R);
    
26.     Xh = Xh + K*(ys(i) - C*Xh);
    
27.     P = P - K*C*P;
    
28.     yout(i) = C*Xh;
    
29. end
    
30. Ts = 1;                    采样时间
    
31. s = tf('s');
    
32. sysc = 100/(60*s +1);      真实系统的传递函数
    
33. sysd = c2d(sysc, Ts);
    
34. 
    
35. t = 0:Ts:(Ts*300);
    
36. u = ones(1, length(t));     系统输入
    
37. ys = lsim(sysd, u, t, 0);   系统输出
    
38. ys = ys + 10*rand(size(ys));  测量量 = 系统输出 + 噪声
    
39. 
    
40. A = [1, 0; 0, 1];
    
41. B = [1; 0];
    
42. C = [1, 0];
    
43. N = length(A);          系统维数
    
44. Len = length(u);
    
45. 
    
46. yout = zeros(Len, 1);   滤波器输出
    
47. Xh = zeros(N, 1);       状态变量
    
48. P = eye(N);              
    
49. Q = 2 * eye(N);         系统噪声
    
50. R = 50;                 测量噪声
    
51. for i = 1 : Len
    
52.     Xh = A*Xh + B*u(i);
    
53.     P = A*P*A' + Q;
    
54.     K = P*C'* inv(C*P*C' + R);
    
55.     Xh = Xh + K*(ys(i) - C*Xh);
    
56.     P = P - K*C*P;
    
57.     yout(i) = C*Xh;
    
58. end
    
59. 
    
60. plot(t, ys, t, yout);

复制代码

aaaalook

使用最小二乘法例子的代码：

1. 
   
2. 
   
3.     Ts = 1;               采样时间
   
4.     s = tf('s');
   
5.     sysc = 100/(60*s +1);   真实系统的传递函数
   
6.     sysd = c2d(sysc, Ts);
   
7. 
   
8.     t = 0:Ts:(Ts*300);
   
9.     u = ones(1, length(t));   系统输入
   
10.     ys = lsim(sysd, u, t, 0); 系统输出
    
11.     ys = ys + 10*rand(size(ys)); 测量量 = 系统输出 + 噪声
    
12. 
    
13.     [numd, dend] = LeastSquare(u, ys, 3);   最小二乘法获取预测系统模型
    
14.     [A, B, C, D] = tf2ss(numd, dend);       变换到状态空间形式
    
15.     Len = length(u);
    
16.     N = length(A);         系统维数
    
17. 
    
18.     yout = zeros(Len, 1);  滤波器输出
    
19.     Xh = zeros(N, 1);      状态变量
    
20.     P = eye(N);
    
21.     Q = 0.02 * eye(N);     系统噪声
    
22.     R = 50;                测量噪声
    
23.     for i = 1 : Len
    
24.         Xh = A*Xh + B*u(i);
    
25.         P = A*P*A' + Q;
    
26. 
    
27.         K = P*C'* inv(C*P*C' + R);
    
28.         Xh = Xh + K*(ys(i) - C*Xh);
    
29.         P = P - K*C*P;
    
30.         yout(i) = C*Xh;
    
31.     end
    
32. 
    
33.     plot(t, ys, t, yout);

复制代码