-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathKalman2D.cpp
66 lines (54 loc) · 1.31 KB
/
Kalman2D.cpp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
#include "Kalman2D.h"
Kalman2D::Kalman2D()
{
//ctor
m_x = Matrix(2,1);
m_p = Matrix(2,2);
m_q = Matrix(2,2);
}
Kalman2D::~Kalman2D()
{
//dtor
}
void Kalman2D::Reset(double qx, double qv, double r, double pd, double ix){
m_q.Set(1,1,qx*qx);
m_q.Set(1,2,qx*qv);
m_q.Set(2,1,qx*qv);
m_q.Set(2,2,qv*qv);
m_r = r;
m_p.Set(1,1,pd);
m_p.Set(2,2,pd);
m_p.Set(1,2,0);
m_p.Set(2,1,0);
m_x.Set(1,1,ix);
m_x.Set(2,1,0);
}
double Kalman2D::Update(double mx,double mv, double dt){
double data_f[] = {1,0,dt,1};
double data_h[] = {1,0,0,1};
double data_ht[] = {1,0,0,1};
double data_id[] = {1,0,0,1};
double data_z[] = {mx,mv};
Matrix f = Matrix(2,2,data_f);
Matrix h = Matrix(2,2,data_h);
Matrix ht = Matrix(2,2,data_ht);
Matrix id = Matrix(2,2,data_id);
Matrix z = Matrix(2,1,data_z);
//预测系统状态,其中U={0,0}
// X = F*X + H*U
m_x = f * m_x;
//f.DisplayMtr();
//m_x.DisplayMtr();
//预测系统方差
// P = F*P*F^T + Q
m_p = f * m_p * (f.TransposeMtr()) + m_q;
//修正卡尔曼混合系数
Matrix tmp = h * m_p * ht + m_r;
Matrix sinv = tmp.InverseMtr();
Matrix k = m_p * ht * sinv;
//修正系统状态
m_x = (id - k) * m_x + k * z;
//修正系统均方差
m_p = (id - k*h) * m_p;
return m_x.Get(1,1);
}