root/tests/testKF.cpp @ 24

Revision 22, 1.5 kB (checked in by smidl, 16 years ago)

upravy Kalmana

Line 
1#include <itpp/itbase.h>
2#include <estim/libKF.h>
3
4using namespace itpp;
5
6//These lines are needed for use of cout and endl
7using std::cout;
8using std::endl;
9
10int main() {
11
12
13        // Klaman filter
14        mat A, B,C,D,R,Q,P0;
15        vec mu0;
16        mat Mu0;;
17        // input from Matlab
18        it_file fin( "testKF.it" );
19
20        mat Dt, Xt,Xt2,XtE;
21        int Ndat;
22
23        bool xxx= fin.seek( "d" );
24        fin >>Dt;
25        fin.seek( "A" ); 
26        fin >> A;
27        fin.seek( "B" ); 
28        fin >> B;
29        fin.seek( "C" ); 
30        fin >> C;
31        fin.seek( "D" ); 
32        fin >> D;
33        fin.seek( "R" ); 
34        fin >> R;
35        fin.seek( "Q" ); fin >> Q;
36        fin.seek( "P0" ); fin >> P0;
37        fin.seek( "mu0" ); fin >> Mu0; 
38        mu0=Mu0.get_col(0);
39       
40        Ndat = Dt.cols();
41        Xt=zeros( 2,Ndat );
42        Xt2=zeros( 2,Ndat );
43        XtE=zeros( 2,Ndat );
44
45        Kalman<fsqmat> KF( A,B,C,D,fsqmat(R),fsqmat(Q),fsqmat(P0),mu0 );
46//      cout << KF;
47        KalmanFull KF2( A,B,C,D,R,Q,P0,mu0 );
48        RV rx("1","{x}","2","0","0");
49        RV ru("2","{u}","1","0","0");
50        RV ry("3","{y}","1","0","0");
51        bilinfn fxu(rx,ru,A,B);
52        bilinfn hxu(rx,ru,C,D);
53        EKF<fsqmat> KFE(&fxu,&hxu,Q,R,mu0,P0);
54
55        Xt.set_col( 0,KF.mu );
56        Xt2.set_col( 0,KF.mu );
57        XtE.set_col( 0,KF.mu );
58        for ( int t=1;t<Ndat;t++ ) {
59                KF.bayes( Dt.get_col( t ));
60                KF2.bayes( Dt.get_col( t ));
61                KFE.bayes( Dt.get_col( t ));
62                Xt.set_col(t,KF.mu);
63                Xt2.set_col(t,KF2.mu);
64                XtE.set_col(t,KFE.mu);
65//              Kmu = KF.mu;
66//              cout <<  "t:" <<t<< "  " << dt<<"  "<<Kmu <<endl;
67        }
68
69        it_file fou( "testKF_res.it" );
70        fou << Name("xth") << Xt;
71        fou << Name("xth2") << Xt2;
72        fou << Name("xthE") << XtE;
73        //Exit program:
74        return 0;
75
76}
Note: See TracBrowser for help on using the browser.