Révision 546
Ajouté par Khalid YAZID il y a presque 3 ans
branch/yazid/sp4a3/sp4a3_kalman.c | ||
---|---|---|
debug=0; ///Mettre à 1 pour afficher les matrices.
|
||
///Ajouter votre code ci-dessous///
|
||
// Kalman
|
||
double Mobs[2][1] ={{xobs},{yobs}};
|
||
double Mobs[2][1] ={{xobs},{yobs}}; //{{xobs},{yobs}}
|
||
// X = F*X
|
||
Mul_Mat_Mat(4,4,F, 4,1,X, X1);
|
||
|
||
Plot_Mat(X1," X(k+1|k) = ");
|
||
|
||
//P = F*P*F'+Q;
|
||
Mul_Mat_Mat(4,4,F, 4,4, P,P1);
|
||
Mul_Mat_Mat(4,4,P1, 4,4, FT,P2);
|
||
Add_Mat_Mat(4,4,P2,4,4,Q, P3);
|
||
Mul_Mat_Mat(4,4,F, 4,4, P,P1); //P = F*P
|
||
Mul_Mat_Mat(4,4,P1, 4,4, FT,P2); //P = F*P*F'
|
||
Add_Mat_Mat(4,4,P2,4,4,Q, P3); //P = F*P*F'+Q;
|
||
|
||
|
||
Plot_Mat(P3,"P(k+1|k) = F.P(k|k).FT + Q = ");
|
||
|
||
// K = P*H' / ( H*P*H' + R);
|
||
Mul_Mat_Mat(2,4,H, 4,4, P3,K1);
|
||
Mul_Mat_Mat(2,4,K1, 4,2, HT,K2);
|
||
Mul_Mat_Mat(2,4,H, 4,4, P3,K1); // H*P
|
||
Mul_Mat_Mat(2,4,K1, 4,2, HT,K2); // H*P*H'
|
||
Add_Mat_Mat(2,2,K2,2,2,R, K3);
|
||
Inverse_Mat_22(2,2,K3,K4);
|
||
Mul_Mat_Mat(4,4,P3, 4,2, HT,K5);
|
||
Mul_Mat_Mat(4,2,K5, 2,2, K4,K6);
|
||
Mul_Mat_Mat(4,4,P3, 4,2, HT,K5); //P*H'
|
||
Mul_Mat_Mat(4,2,K5, 2,2, K4,K6); // K = P*H' / ( H*P*H' + R);
|
||
|
||
Plot_Mat(K6,"K = ");
|
||
|
||
//X = X + K*([xobs(i);yobs(i)]-H*X);
|
||
|
||
Mul_Mat_Mat(2,4,H, 4,1, X1,X2);
|
||
Sub_Mat_Mat(2,1,Mobs,2,1,X2,X3) ;
|
||
Mul_Mat_Mat(2,4,H, 4,1, X1,X2); //H*X
|
||
Sub_Mat_Mat(2,1,Mobs,2,1,X2,X3) ; //[xobs(i);yobs(i)]-H*X
|
||
Plot_Mat(X3,"DELTA = Obs - H.X(k+1|k)");
|
||
Mul_Mat_Mat(4,2,K6, 2,1,X3,X4);
|
||
Add_Mat_Mat(4,1,X1,4,1,X4,X5);
|
||
Mul_Mat_Mat(4,2,K6, 2,1,X3,X4); //K*([xobs(i);yobs(i)]-H*X);
|
||
Add_Mat_Mat(4,1,X1,4,1,X4,X5); //X = X + K*([xobs(i);yobs(i)]-H*X);
|
||
|
||
Plot_Mat(X5," X(k+1|k+1) = X(k+1|k) + K.Delta = ");
|
||
|
||
// P = P - K*H*P;
|
||
Mul_Mat_Mat(4,2,K6, 2,4, H,P4);
|
||
Mul_Mat_Mat(4,4,P4, 4,4, P3,P5);
|
||
Sub_Mat_Mat(4,4,P3,4,4,P5, P6);
|
||
Mul_Mat_Mat(4,2,K6, 2,4, H,P4); //K*H
|
||
Mul_Mat_Mat(4,4,P4, 4,4, P3,P5); //K*H*P
|
||
Sub_Mat_Mat(4,4,P3,4,4,P5, P6); // P = P - K*H*P;
|
||
|
||
Plot_Mat(P6," P(k+1|k+1) = P(k+1|k) - K.H.P(k+1|k) = ");
|
||
|
Formats disponibles : Unified diff
Indentation et commentaire programme finale