1
+ # please do not change this file
2
+ import numpy as np
3
+ import matplotlib .pyplot as plt
4
+
5
+
6
+ def clamp (n , minn , maxn ):
7
+ return max (min (maxn , n ), minn )
8
+
9
+
10
+ def wrap2pi (a ):
11
+ return (a + np .pi ) % (2 * np .pi ) - np .pi
12
+
13
+
14
+ def addGaussianNoise (a ,u ,sigma ):
15
+ # add Gaussian noise to a scalar a, with mean u and covariance sigma
16
+ a += np .random .normal (u , sigma , 1 )
17
+ return a [0 ]
18
+
19
+
20
+ class vehicle :
21
+ # parameters
22
+ lr = 1.7
23
+ lf = 1.1
24
+ Ca = 15000.0
25
+ Iz = 3344.0
26
+ f = 0.01
27
+ m = 2000.0
28
+ g = 10
29
+
30
+ def __init__ (self , state ):
31
+ # initialize the vehicle
32
+ self .state = state
33
+ self .observation = vehicle .observation (
34
+ state .xd , state .yd , state .phid , state .delta , state .X , state .Y , state .phi )
35
+
36
+ class state :
37
+ deltaMax = np .pi / 6
38
+ deltaMin = - np .pi / 6
39
+ xdMax = 100.0
40
+ xdMin = 0.0
41
+ ydMax = 10.0
42
+ ydMin = - 10.0
43
+
44
+ def __init__ (self , xd = 0.0 , yd = 0.0 , phid = 0.0 , delta = .0 , X = .0 , Y = .0 , phi = .0 ):
45
+ self .xd = xd
46
+ self .yd = yd
47
+ self .phid = phid
48
+ self .delta = clamp (delta , self .deltaMin , self .deltaMax )
49
+ self .X = X
50
+ self .Y = Y
51
+ self .phi = wrap2pi (phi )
52
+ def showState (self ):
53
+ print ('xd\t ' , self .xd , 'yd\t ' , self .yd , 'phid\t ' , self .phid , 'delta\t ' , self .delta ,
54
+ 'X\t ' , self .X , 'Y\t ' , self .Y , 'phi\t ' , self .phi )
55
+
56
+ class observation (state ):
57
+ def __init__ (self , xd = 0.0 , yd = 0.0 , phid = 0.0 , delta = .0 , X = .0 , Y = .0 , phi = .0 ):
58
+ vehicle .state .__init__ (self , xd , yd , phid , delta , X , Y , phi )
59
+ self .addNoise ()
60
+
61
+ def copyState (self , state ):
62
+ self .xd = state .xd
63
+ self .yd = state .yd
64
+ self .phid = state .phid
65
+ self .delta = state .delta
66
+ self .X = state .X
67
+ self .Y = state .Y
68
+ self .phi = state .phi
69
+
70
+ def addNoise (self ):
71
+ self .xd = addGaussianNoise (self .xd , 0 , 0.5 )
72
+ self .yd = addGaussianNoise (self .yd , 0 , 0.5 )
73
+ self .phid = addGaussianNoise (self .phid , 0 , 0.05 )
74
+ self .delta = clamp (addGaussianNoise (self .delta , 0 , 0.05 ), self .deltaMin , self .deltaMax )
75
+ self .X = addGaussianNoise (self .X , 0 , 1 )
76
+ self .Y = addGaussianNoise (self .Y , 0 , 1 )
77
+ self .phi = wrap2pi (addGaussianNoise (self .phi , 0 , 0.5 ))
78
+
79
+ class command :
80
+ # F: N
81
+ # deltad: rad/s
82
+ deltadMax = np .pi / 6.0
83
+ deltadMin = - np .pi / 6.0
84
+ FMax = 10000.0
85
+ Fmin = - 10000.0
86
+
87
+ def __init__ (self , F_ = 0.0 , deltad_ = 0.0 ):
88
+ self .F = clamp (F_ , self .Fmin , self .FMax )
89
+ self .deltad = clamp (deltad_ , self .deltadMin , self .deltadMax )
90
+
91
+ def showCommand (self ):
92
+ print ('F:\t ' , self .F , 'deltad:\t ' , self .deltad )
93
+
94
+ def update (self , command ):
95
+ # time step
96
+ dt = 0.05
97
+ # update state
98
+ Ff = np .sign (self .state .xd )* self .f * self .m * self .g
99
+ Ftotal = command .F - Ff if np .abs (command .F )>= np .abs (Ff ) else 0
100
+ # print(Ftotal)
101
+ ax = 1 / self .m * Ftotal
102
+ if np .abs (self .state .xd ) <= 0.5 :
103
+ Fyf = 0.0
104
+ Fyr = 0.0
105
+ else :
106
+ Fyf = 2.0 * self .Ca * (self .state .delta - (self .state .yd + self .lf * self .state .phid )/ self .state .xd )
107
+ Fyr = 2.0 * self .Ca * (- (self .state .yd - self .lr * self .state .phid )/ self .state .xd )
108
+
109
+ xdd = self .state .phid * self .state .yd + ax
110
+ ydd = - self .state .phid * self .state .xd + 1.0 / self .m * (Fyf * np .cos (self .state .delta ) - Fyr )
111
+ phidd = 1.0 / self .Iz * (self .lf * Fyf - self .lr * Fyr )
112
+ Xd = self .state .xd * np .cos (self .state .phi ) - self .state .yd * np .sin (self .state .phi )
113
+ Yd = self .state .xd * np .sin (self .state .phi ) + self .state .yd * np .cos (self .state .phi )
114
+ self .state .xd += xdd * dt
115
+ self .state .yd += ydd * dt
116
+ self .state .phid += phidd * dt
117
+
118
+ self .state .delta += command .deltad * dt
119
+ # self.state.delta = command.deltad
120
+
121
+
122
+ self .state .X += Xd * dt
123
+ self .state .Y += Yd * dt
124
+ self .state .phi += self .state .phid * dt
125
+ self .applyConstrain ()
126
+ # update observation
127
+ self .observation .copyState (self .state )
128
+ self .observation .addNoise ()
129
+ # self.state.showState()
130
+
131
+
132
+ def applyConstrain (self ):
133
+ # phi should be between +- pi
134
+ self .state .phi = wrap2pi (self .state .phi )
135
+ # state constraint
136
+ self .state .delta = clamp (self .state .delta ,self .state .deltaMin ,self .state .deltaMax )
137
+ self .state .xd = clamp (self .state .xd , self .state .xdMin , self .state .xdMax )
138
+ self .state .yd = clamp (self .state .yd , self .state .ydMin , self .state .ydMax )
139
+
140
+ def showState (self ):
141
+ self .state .showState ()
142
+
143
+
144
+ if __name__ == "__main__" :
145
+ a = vehicle (vehicle .state (Y = 0.0 ,xd = 1 ))
146
+ n = 1000
147
+
148
+ X = []
149
+ Y = []
150
+ delta = []
151
+ xd = []
152
+ yd = []
153
+ phi = []
154
+ phid = []
155
+ for i in range (n ):
156
+ # if i% 1 ==0:
157
+ X .append (a .state .X )
158
+ Y .append (a .state .Y )
159
+ delta .append (a .state .delta )
160
+ xd .append (a .state .xd )
161
+ yd .append (a .state .yd )
162
+ phid .append (a .state .phid )
163
+ phi .append (a .state .phi )
164
+ if a .state .xd > 3 :
165
+ c = a .command (deltad_ = np .sin (i / 10 ), F_ = - 10000 )
166
+ else :
167
+ c = a .command (deltad_ = np .sin (i / 10 ), F_ = 10000.0 )
168
+ a .update (command = c )
169
+
170
+ plt .subplot (321 )
171
+ plt .title ('position' )
172
+ plt .plot (X ,Y ,'r' )
173
+
174
+ plt .subplot (322 )
175
+ plt .title ('delta' )
176
+ plt .plot (delta , 'r' )
177
+
178
+ plt .subplot (323 )
179
+ plt .title ('xd' )
180
+ plt .plot (xd , 'r' )
181
+
182
+ plt .subplot (324 )
183
+ plt .title ('yd' )
184
+ plt .plot (yd , 'r' )
185
+
186
+ plt .subplot (325 )
187
+ plt .title ('phi' )
188
+ plt .plot (phi , 'r' )
189
+
190
+ plt .subplot (326 )
191
+ plt .title ('phid' )
192
+ plt .plot (phid , 'r' )
193
+
194
+ plt .show ()
0 commit comments