root/tests/pmsm_unkQ.cpp @ 37

Revision 37, 2.4 kB (checked in by smidl, 16 years ago)

Matrix in Cholesky decomposition, Square-root Kalman and many bug fixes

Line 
1/*
2  \file
3  \brief Models for synchronous electric drive using IT++ and BDM
4  \author Vaclav Smidl.
5
6  -----------------------------------
7  BDM++ - C++ library for Bayesian Decision Making under Uncertainty
8
9  Using IT++ for numerical operations
10  -----------------------------------
11*/
12
13#include <itpp/itbase.h>
14#include <estim/libKF.h>
15#include <estim/libPF.h>
16#include <stat/libFN.h>
17
18#include "pmsm.h"
19
20using namespace itpp;
21/*
22// PMSM with Q on Ia and Ib given externally
23class EKF_unQ : public EKF<chmat> , public BMcond {
24public:
25        EKF_unQ( rx,ry,ru,rQ):EKF<chmat>(rx,ry,ru),BMcond(rQ){};
26        void condition(const vec &Q0){};
27};*/
28
29
30int main() {
31        // Kalman filter
32        int Ndat = 1000;
33
34//      cout << KF;
35        // internal model
36        IMpmsm fxu;
37        //                  Rs    Ls        dt       Fmag(Ypm) kp   p    J     Bf(Mz)
38        fxu.set_parameters ( 0.28, 0.003465, 20*1e-6, 0.1989,   1.5 ,4.0, 0.04, 0.0 );
39        // observation model
40        OMpmsm hxu;
41
42        vec mu0= "100 100 100 1";
43        vec Qdiag ( "0.1 0.1 0.01 0.01" );
44        vec Rdiag ( "0.02 0.02" );
45        vec vQ = "0.01:0.01:100";
46        chmat Q ( Qdiag );
47        chmat R ( Rdiag );
48        EKFCh KFE ( rx,ry,ru );
49        KFE.set_est ( mu0, chmat ( 1000*ones ( 4 ) ) );
50        KFE.set_parameters ( &fxu,&hxu,Q,R );
51
52        mat ll(100,Ndat);
53
54        EKFCh* kfArray[100];
55
56        for ( int i=0;i<100;i++ ) {
57                vec Qid ( Qdiag );
58                Qid ( 0 ) = vQ ( i ); Qid ( 1 ) = vQ ( i );
59                kfArray[i]= new EKFCh ( rx,ry,ru );
60                kfArray[i]->set_est ( mu0, chmat ( 1000*ones ( 4 ) ) );
61                kfArray[i]->set_parameters ( &fxu,&hxu,chmat ( Qid ),R );
62        }
63
64        epdf& KFEep = KFE._epdf();
65
66        //simulator values
67        vec dt ( 2 );
68        vec wt ( 2 );
69        vec ut ( 2 );
70        vec xt=mu0;
71        vec et ( 4 );
72
73        mat Xt=zeros ( 4,Ndat );
74        mat XtE=zeros ( 4,Ndat );
75        Xt.set_col ( 0,KFEep.mean() );
76
77        for ( int t=1;t<Ndat;t++ ) {
78                //simulator
79                UniRNG.sample_vector ( 2,wt );
80
81                if ( rem ( t,500 ) <200 ) ut = rem ( t,500 ) *ones ( 2 );
82                else
83                        ut=zeros ( 2 );
84
85                NorRNG.sample_vector ( 4,et );
86                NorRNG.sample_vector ( 2,wt );
87                xt = fxu.eval ( xt,ut ) + Q.sqrt_mult ( et );
88                dt = hxu.eval ( xt,ut ) + R.sqrt_mult ( wt );
89
90                //estimator
91                KFE.bayes ( concat ( dt,ut ) );
92                for ( int i=0;i<100;i++ ) {kfArray[i]->bayes( concat ( dt,ut ) );ll(i,t)=ll(i,t-1) + kfArray[i]->_ll();
93                }
94               
95                Xt.set_col ( t,xt );
96                XtE.set_col ( t,KFEep.mean() );
97        }
98
99        it_file fou ( "pmsm.it" );
100
101        fou << Name ( "xth" ) << Xt;
102        fou << Name ( "xthE" ) << XtE;
103        fou << Name ( "ll" ) << ll;
104        //Exit program:
105        return 0;
106
107}
Note: See TracBrowser for help on using the browser.