|
10 | 10 | import math |
11 | 11 | import matplotlib.pyplot as plt |
12 | 12 |
|
| 13 | +Cx = np.diag([1.0, 1.0, math.radians(30.0)])**2 |
| 14 | + |
13 | 15 | # Simulation parameter |
14 | 16 | Qsim = np.diag([0.2])**2 |
15 | 17 | Rsim = np.diag([1.0, math.radians(30.0)])**2 |
|
18 | 20 | SIM_TIME = 50.0 # simulation time [s] |
19 | 21 | MAX_RANGE = 20.0 # maximum observation range |
20 | 22 |
|
| 23 | +POSE_SIZE = 3 # [x,y,yaw] |
| 24 | +LM_SIZE = 2 # [x,y] |
| 25 | + |
21 | 26 | # Particle filter parameter |
22 | 27 | NP = 100 # Number of Particle |
23 | 28 | NTh = NP / 2.0 # Number of particle for re-sampling |
@@ -61,26 +66,43 @@ def observation(xTrue, xd, u, RFID): |
61 | 66 |
|
62 | 67 | def motion_model(x, u): |
63 | 68 |
|
64 | | - F = np.matrix([[1.0, 0, 0, 0], |
65 | | - [0, 1.0, 0, 0], |
66 | | - [0, 0, 1.0, 0], |
67 | | - [0, 0, 0, 0]]) |
| 69 | + F = np.matrix([[1.0, 0, 0], |
| 70 | + [0, 1.0, 0], |
| 71 | + [0, 0, 1.0]]) |
68 | 72 |
|
69 | 73 | B = np.matrix([[DT * math.cos(x[2, 0]), 0], |
70 | 74 | [DT * math.sin(x[2, 0]), 0], |
71 | | - [0.0, DT], |
72 | | - [1.0, 0.0]]) |
| 75 | + [0.0, DT]]) |
73 | 76 |
|
74 | 77 | x = F * x + B * u |
75 | 78 |
|
76 | 79 | return x |
77 | 80 |
|
78 | 81 |
|
| 82 | +def calc_n_LM(x): |
| 83 | + n = int((len(x) - POSE_SIZE) / LM_SIZE) |
| 84 | + return n |
| 85 | + |
| 86 | + |
| 87 | +def jacob_motion(x, u): |
| 88 | + |
| 89 | + Fx = np.hstack((np.eye(POSE_SIZE), np.zeros( |
| 90 | + (POSE_SIZE, LM_SIZE * calc_n_LM(x))))) |
| 91 | + |
| 92 | + jF = np.matrix([[0.0, 0.0, -DT * u[0] * math.sin(x[2, 0])], |
| 93 | + [0.0, 0.0, DT * u[0] * math.cos(x[2, 0])], |
| 94 | + [0.0, 0.0, 0.0]]) |
| 95 | + |
| 96 | + G = np.eye(POSE_SIZE) + Fx.T * jF * Fx |
| 97 | + |
| 98 | + return G, Fx, |
| 99 | + |
| 100 | + |
79 | 101 | def ekf_slam(xEst, PEst, u, z): |
80 | 102 | # Predict |
81 | 103 | xEst = motion_model(xEst, u) |
82 | | - # [G,Fx]=jacobF(xEst, u); |
83 | | - # PEst= G'*PEst*G + Fx'*R*Fx; |
| 104 | + G, Fx = jacob_motion(xEst, u) |
| 105 | + PEst = G.T * PEst * G + Fx.T * Cx * Fx |
84 | 106 |
|
85 | 107 | return xEst, PEst |
86 | 108 |
|
@@ -137,11 +159,11 @@ def main(): |
137 | 159 | [-5.0, 20.0]]) |
138 | 160 |
|
139 | 161 | # State Vector [x y yaw v]' |
140 | | - xEst = np.matrix(np.zeros((4, 1))) |
141 | | - xTrue = np.matrix(np.zeros((4, 1))) |
142 | | - PEst = np.eye(4) |
| 162 | + xEst = np.matrix(np.zeros((POSE_SIZE, 1))) |
| 163 | + xTrue = np.matrix(np.zeros((POSE_SIZE, 1))) |
| 164 | + PEst = np.eye(POSE_SIZE) |
143 | 165 |
|
144 | | - xDR = np.matrix(np.zeros((4, 1))) # Dead reckoning |
| 166 | + xDR = np.matrix(np.zeros((POSE_SIZE, 1))) # Dead reckoning |
145 | 167 |
|
146 | 168 | # history |
147 | 169 | hxEst = xEst |
|
0 commit comments