/*!
  \file
  \brief TR 2525 file for testing Wishart random walk on PMSM simulator
  \author Vaclav Smidl.

  -----------------------------------
  BDM++ - C++ library for Bayesian Decision Making under Uncertainty

  Using IT++ for numerical operations
  -----------------------------------
*/



#include <estim/particles.h>
#include <estim/ekf_template.h>
#include <stat/functions.h>

#include "user_info.h"
#include "../pmsmDS.h"

using namespace bdm;

int main ( int argc, char* argv[] ) {
	const char *fname;
	if ( argc>1 ) {fname = argv[1];	}
	else { fname = "pmsm_wishart.cfg"; }
	UIFile F ( fname );

	double h = 1e-6;
	int Nsimstep = 125;

	// Kalman filter
	double k;
	double l;
	F.lookupValue ( "k",k);
	F.lookupValue ( "l",l);

	int Ndat;
	int Npart;
	F.lookupValue ( "ndat", Ndat );
	F.lookupValue ( "Npart", Npart );

	vec Qdiag;
	vec Rdiag;
	UI::get( Qdiag, F, "dQ" ); //( "1e-6 1e-6 0.001 0.0001" ); //zdenek: 0.01 0.01 0.0001 0.0001
	UI::get( Rdiag, F, "dR" );// ( "1e-8 1e-8" ); //var(diff(xth)) = "0.034 0.034"
		
	pmsmDS* DS = UI::build<pmsmDS>(F,"system");

	// internal model
	IMpmsm fxu;
	//                  Rs    Ls        dt       Fmag(Ypm)    kp   p    J     Bf(Mz)
	fxu.set_parameters ( 0.28, 0.003465, Nsimstep*h, 0.1989, 1.5 ,4.0, 0.04, 0.0 );
	// observation model
	OMpmsm hxu;

	vec mu0= "0.0 0.0 0.0 0.0";
	chmat Q ( Qdiag );
	chmat R ( Rdiag );
	EKFCh KFE ;
	KFE.set_parameters ( &fxu,&hxu,Q,R );
	KFE.set_est ( mu0, chmat ( 0.01*eye( 4 ) ) );
	KFE.set_rv ( rx );

	RV rQ ( "{Q }","16" );
	EKFCh_chQ KFEp ;
	KFEp.set_parameters ( &fxu,&hxu,Q,R );
	KFEp.set_est ( mu0, chmat (0.01* eye( 4 ) ) );

	rwiWishartCh* evolQw = new rwiWishartCh; 
	evolQw->set_parameters(4, k, sqrt(Qdiag),l);
	MPF<EKFCh_chQ> M;
	M.set_parameters ( evolQw,evolQw,Npart );
	// initialize
	chmat Ch0(diag(Qdiag));
	evolQw->condition ( vec(Ch0._Ch()._data(),16) ); //Zdenek default
	M.set_statistics ( evolQw->_e() , &KFEp );
	//

	M.set_rv ( concat ( rQ,rx ) );

	dirfilelog *L = UI::build<dirfilelog> ( F, "logger" );// ( "exp/mpf_test",100 );
	
	KFE.log_level[logbounds] = true;
	KFE.log_add ( *L,"KF" );
	M.log_level[logbounds] = true;
	M.log_add ( *L,"M" );
	DS->log_add(*L);
	L->init();

	// SET SIMULATOR
	pmsmsim_set_parameters ( 0.28,0.003465,0.1989,0.0,4,1.5,0.04, 200., 3e-6, h );
	vec dt ( 2 );
	vec ut ( 2 );
	vec xt ( 4 );
	vec xtm=zeros ( 4 );
//	double Ww=0.0;
	
	vec obs(6);
	for ( int tK=1;tK<Ndat;tK++ ) {
		//Number of steps of a simulator for one step of Kalman
		DS->step();
		DS->getdata(obs);
		
		dt=obs.left(2);
		ut=obs.right(2);
		//estimator
		KFE.bayes ( concat ( dt,ut ) );
		M.bayes ( concat ( dt,ut ) );

		DS->logit (*L);

		KFE.logit ( *L );
		M.logit ( *L );
		L->step();
	}
	L->finalize();
	//Exit program:

	delete L;
	return 0;
}