00001 
00013 #include <stat/loggers.h>
00014 #include <estim/libKF.h>
00015 #include "simulator.h"
00016 #include "pmsm.h"
00017 
00019 class pmsmDS : public DS {
00020 
00021 protected:
00023         int L_x, L_ou, L_oy, L_iu, L_optu;
00025         vec profileWw;
00027         vec profileMz;
00029         double dt_prof;
00031         int Dt;
00033         bool opt_modu;
00035 public:
00037         pmsmDS ()   {Dt=125; Drv=RV ( "{o_ua o_ub o_ia o_ib t_ua t_ub }" );}
00038         void set_parameters ( double Rs0, double Ls0, double Fmag0, double Bf0, double p0, double kp0, double J0, double Uc0, double DT0, double dt0 ) {
00039                 pmsmsim_set_parameters ( Rs0, Ls0, Fmag0, Bf0, p0, kp0, J0, Uc0, DT0, dt0 );
00040         }
00042         void set_options ( string &opt ) {
00043                 opt_modu = ( opt.find ( "modelu" ) ==string::npos );
00044         }
00045         void getdata ( vec &dt ) {dt=vec ( KalmanObs,6 );}
00046         void write ( vec &ut ) {}
00047 
00048         void step() {
00049                 static int ind=0;
00050                 static double dW; 
00051                 static double Ww; 
00052                 static double Mz; 
00053                 if ( t>=dt_prof*ind ) {
00054                         ind++;
00055                         
00056                         if ( ind<profileWw.length() ) {
00057                                 
00058                                 if ( profileWw.length() ==1 ) {
00059                                         Ww=profileWw ( 0 ); dW=0.0;
00060                                 }
00061                                 else {
00062                                         dW = profileWw ( ind )-profileWw ( ind-1 );
00063                                         dW *=125e-6/dt_prof;
00064                                 }
00065                         }
00066                         else {
00067                                 dW = 0;
00068                         }
00069                         
00070                         if ( ind<profileMz.length() ) {
00071                                 
00072                                 Mz = profileMz(ind);
00073                         }
00074                         else {
00075                                 Mz = 0;
00076                         }
00077                 }
00078                 Ww += dW;
00079                 
00080                 for ( int i=0;i<Dt;i++ ) {      pmsmsim_step ( Ww , Mz);}
00081         };
00082 
00083         void log_add ( logger &L ) {
00084                 L_x = L.add ( rx, "x" );
00085                 L_oy = L.add ( ry, "o" );
00086                 L_ou = L.add ( ru, "o" );
00087                 L_iu = L.add ( ru, "t" );
00088                 
00089                 if ( opt_modu ) {
00090                         L_optu = L.add ( ru, "model" );
00091                 }
00092         }
00093 
00094         void logit ( logger &L ) {
00095                 L.logit ( L_x, vec ( x,4 )      );
00096                 L.logit ( L_oy, vec_2 ( KalmanObs[2],KalmanObs[3] ) );
00097                 L.logit ( L_ou, vec_2 ( KalmanObs[0],KalmanObs[1] ) );
00098                 L.logit ( L_iu, vec_2 ( KalmanObs[4],KalmanObs[5] ) );
00099                 if ( opt_modu ) {
00100                         double sq3=sqrt ( 3.0 );
00101                         double ua,ub;
00102                         double i1=x[0];
00103                         double i2=0.5* ( -i1+sq3*x[1] );
00104                         double i3=0.5* ( -i1-sq3*x[1] );
00105                         double u1=KalmanObs[0];
00106                         double u2=0.5* ( -u1+sq3*KalmanObs[1] );
00107                         double u3=0.5* ( -u1-sq3*KalmanObs[1] );
00108 
00109                         double du1=0.7* ( double ( i1>0.1 ) - double ( i1<-0.1 ) ) +0.05*i1;
00110                         double du2=0.7* ( double ( i2>0.1 ) - double ( i2<-0.1 ) ) +0.05*i2;
00111                         double du3=0.7* ( double ( i3>0.1 ) - double ( i3<-0.1 ) ) +0.05*i3;
00112                         ua = ( 2.0* ( u1-du1 )- ( u2-du2 )- ( u3-du3 ) ) /3.0;
00113                         ub = ( ( u2-du2 )- ( u3-du3 ) ) /sq3;
00114                         L.logit ( L_optu , vec_2 ( ua,ub ) );
00115                 }
00116                 
00117         }
00118 
00119         void set_profile ( double dt, const vec &Ww, const vec &Mz ) {dt_prof=dt; profileWw=Ww; profileMz=Mz;}
00120 };
00121 
00123 class pmsmCRB : public EKFCh{
00124         protected:
00125                 vec interr;
00126                 vec old_true;
00127                 vec secder;
00128                 int L_CRB;
00129                 int L_err;
00130                 int L_sec;
00131         public:
00133         pmsmCRB():EKFCh(){old_true=zeros(6);}
00134 
00135         void bayes(const vec &dt){
00136                 
00137                 vec true_state=vec(x,4); 
00138                 est.set_mu(true_state);
00139                 
00140                 
00141                 old_true(4)=KalmanObs[4];
00142                 old_true(5)=KalmanObs[5];
00143                 interr = (true_state - pfxu->eval(old_true));
00144                 
00145                 
00146                 IMpmsm2o* pf = dynamic_cast<IMpmsm2o*>(pfxu);
00147                 if (pf) {secder=pf->eval2o(vec_2(KalmanObs[4],KalmanObs[5]));}
00148                                 
00149                 EKFCh::bayes(dt);
00150                 old_true.set_subvector(0,true_state);
00151         }
00152         
00153         void log_add(logger &L, const string &name="" ){
00154                 L_CRB=L.add(rx,"crb");
00155                 L_err=L.add(rx,"err");
00156                 L_sec=L.add(rx,"d2");
00157         }
00158         void logit(logger &L){
00159                 L.logit(L_err, interr);
00160                 L.logit(L_CRB,diag(est._R().to_mat()));
00161                 L.logit(L_sec,secder);
00162         }
00163 };