/*!
  \file
  \brief Bayesian Filtering for linear Gaussian models (Kalman Filter) and extensions
  \author Vaclav Smidl.

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

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

#ifndef EKFfix_H
#define EKFfix_H


#include <estim/kalman.h>
#include "fixed.h"
#include "matrix.h"
#include "matrix_vs.h"
#include "reference_Q15.h"
#include "parametry_motoru.h"
#include "mpf_double.h"
#include "fast_exp.h"
#include "ekf_mm.h"
#include "qmath.h"

using namespace bdm;

double minQ(double Q); 

void mat_to_int16(const imat &M, int16 *I);
void vec_to_int16(const ivec &v, int16 *I);
void int16_to_mat(int16 *I, imat &M, int rows, int cols);
void int16_to_vec(int16 *I, ivec &v, int len);
void UDtof(const mat &U, const vec &D, imat &Uf, ivec &Df, const vec &xref);


//! EKF for testing q44
class EKFtest: public EKF_UD{
	void bayes ( const vec &yt, const vec &cond ) {
		EKF_UD::bayes(yt,cond);
		vec D =  prior()._R()._D();
		
		if (D(3)>10) D(3) = 10;
		
		prior()._R().__D()=D;
	}
};
UIREGISTER(EKFtest);


class EKFscale: public EKFCh{
	LOG_LEVEL(EKFscale, logCh, logA, logC);
public:
	mat Tx;
	mat Ty;
	
	void from_setting ( const Setting &set ){
		EKFCh::from_setting(set);
		vec v;
		UI::get(v, set, "xmax", UI::compulsory);
		Tx = diag(1./v);
		UI::get(v, set, "ymax", UI::compulsory);
		Ty = diag(1./v);
		
		UI::get(log_level, set, "log_level", UI::optional);
	};
	
	void log_register(logger &L, const string &prefix){
		BM::log_register ( L, prefix );
		
		L.add_vector ( log_level, logCh, RV ("Ch", dimension()*dimension() ), prefix );
		L.add_vector ( log_level, logA, RV ("A", dimension()*dimension() ), prefix );
		L.add_vector ( log_level, logC, RV ("C", dimensiony()*dimension() ), prefix );
	};
	
	void log_write() const{
		BM::log_write();
	    if ( log_level[logCh] ) {
			mat Chsc=Tx*(est._R()._Ch());
			vec v(Chsc._data(), dimension()*dimension());
			if (v(0)<0)
				v= -v;
			log_level.store( logCh, round(v*32767));
		}
		if (log_level[logA]){
			mat Asc = Tx*A*inv(Tx);
			vec v(Asc._data(), dimension()*dimension());
			log_level.store( logA, round(v*32767));
		}
		if (log_level[logC]){
			mat Csc = Ty*C*inv(Tx);
			vec v(Csc._data(), dimensiony()*dimension());
			log_level.store( logC, round(v*32767));
		}
	}
	
};
UIREGISTER(EKFscale);

/*!
\brief Extended Kalman Filter with UD matrices in fixed point16 arithmetic

An approximation of the exact Bayesian filter with Gaussian noices and non-linear evolutions of their mean.
*/
class EKFfixedUD : public BM {
	public:
		LOG_LEVEL(EKFfixedUD,logU, logG, logD, logA, logP);
		
		void init_ekf(double Tv);
		void ekf(double ux, double uy, double isxd, double isyd);
				
		/* Constants - definovat jako konstanty ?? ?kde je vyhodnejsi aby v pameti byli?*/
		int16 Q[16]; /* matrix [4,4] */
		int16 R[4]; /* matrix [2,2] */
		
		int16 x_est[4]; /* estimate and prediction */
		
		int16 PSI[16]; /* matrix [4,4] */
		int16 PSIU[16]; /* matrix PIS*U, [4,4] */
		
		int16 Uf[16]; // upper triangular of covariance (inplace)
		int16 Df[4];  // diagonal covariance
		int16 Dfold[4]; // temp of D
		int16 G[16];  // temp for bierman
		
		int16 cA, cB, cC, cG, cH;  // cD, cE, cF, cI ... nepouzivane
		
		enorm<fsqmat> E;
		mat Ry;
		
	public:
		//! Default constructor
		EKFfixedUD ():BM(),E(),Ry(2,2){
			int16 i;
			for(i=0;i<16;i++){Q[i]=0;}
			for(i=0;i<4;i++){R[i]=0;}
			
			for(i=0;i<4;i++){x_est[i]=0;}
			for(i=0;i<16;i++){Uf[i]=0;}
			for(i=0;i<4;i++){Df[i]=0;}
			for(i=0;i<16;i++){G[i]=0;}
			for(i=0;i<4;i++){Dfold[i]=0;}
			
			for(i=0;i<16;i++){PSI[i]=0;}
			
			set_dim(4);
			dimy = 2;
			dimc = 2;
			E._mu()=zeros(4);
			E._R()=zeros(4,4);
			init_ekf(0.000125);
		};
		//! Here dt = [yt;ut] of appropriate dimensions
		void bayes ( const vec &yt, const vec &ut );
		//!dummy!
		const epdf& posterior() const {return E;};
		void log_register(logger &L, const string &prefix){
			BM::log_register ( L, prefix );
			
				L.add_vector ( log_level, logG, RV("G",16), prefix );
				L.add_vector ( log_level, logU, RV ("U", 16 ), prefix );
				L.add_vector ( log_level, logD, RV ("D", 4 ), prefix );
				L.add_vector ( log_level, logA, RV ("A", 16 ), prefix );
				L.add_vector ( log_level, logP, RV ("P", 16 ), prefix );
				
		};
		//void from_setting();
};

UIREGISTER(EKFfixedUD);

/*!
 * \brief Extended Kalman Filter with UD matrices in fixed point16 arithmetic
 * 
 * An approximation of the exact Bayesian filter with Gaussian noices and non-linear evolutions of their mean.
 */
class EKFfixedUD2 : public BM {
public:
	LOG_LEVEL(EKFfixedUD2,logU, logG, logD, logA, logC, logP);
	
	void init_ekf2(double Tv);
	void ekf2(double ux, double uy, double isxd, double isyd);
	
	/* Constants - definovat jako konstanty ?? ?kde je vyhodnejsi aby v pameti byli?*/
	int16 Q[4]; /* matrix [4,4] */
	int16 R[4]; /* matrix [2,2] */
	
	int16 x_est[2]; /* estimate and prediction */
	int16 y_est[2]; /* estimate and prediction */
	int16 y_old[2]; /* estimate and prediction */
	
	int16 PSI[4]; /* matrix [4,4] */
	int16 PSIU[4]; /* matrix PIS*U, [4,4] */
	int16 C[4]; /* matrix [4,4] */
	
	int16 Uf[4]; // upper triangular of covariance (inplace)
	int16 Df[2];  // diagonal covariance
	int16 Dfold[2]; // temp of D
	int16 G[4];  // temp for bierman
	
	int16 cA, cB, cC, cG, cH;  // cD, cE, cF, cI ... nepouzivane
	
	enorm<fsqmat> E;
	mat Ry;
	
public:
	//! Default constructor
	EKFfixedUD2 ():BM(),E(),Ry(2,2){
		int16 i;
		for(i=0;i<4;i++){Q[i]=0;}
		for(i=0;i<4;i++){R[i]=0;}
		
		for(i=0;i<2;i++){x_est[i]=0;}
		for(i=0;i<2;i++){y_est[i]=0;}
		for(i=0;i<2;i++){y_old[i]=0;}
		for(i=0;i<4;i++){Uf[i]=0;}
		for(i=0;i<2;i++){Df[i]=0;}
		for(i=0;i<4;i++){G[i]=0;}
		for(i=0;i<2;i++){Dfold[i]=0;}
		
		for(i=0;i<4;i++){PSI[i]=0;}
		for(i=0;i<4;i++){C[i]=0;}
		
		set_dim(2);
		dimc = 2;
		dimy = 2;
		E._mu()=zeros(2);
		E._R()=zeros(2,2);
		init_ekf2(0.000125);
	};
	//! Here dt = [yt;ut] of appropriate dimensions
	void bayes ( const vec &yt, const vec &ut );
	//!dummy!
	const epdf& posterior() const {return E;};
	void log_register(logger &L, const string &prefix){
		BM::log_register ( L, prefix );
		
		L.add_vector ( log_level, logG, RV("G2",4), prefix );
		L.add_vector ( log_level, logU, RV ("U2", 4 ), prefix );
		L.add_vector ( log_level, logD, RV ("D2", 2 ), prefix );
		L.add_vector ( log_level, logA, RV ("A2", 4 ), prefix );
		L.add_vector ( log_level, logC, RV ("C2", 4 ), prefix );
		L.add_vector ( log_level, logP, RV ("P2", 4 ), prefix );
		
	};
	//void from_setting();
};

UIREGISTER(EKFfixedUD2);

/*!
 * \brief Extended Kalman Filter with UD matrices in fixed point16 arithmetic
 * 
 * An approximation of the exact Bayesian filter with Gaussian noices and non-linear evolutions of their mean.
 */
class EKFfixedUD3 : public BM {
public:
	LOG_LEVEL(EKFfixedUD3,logU, logG, logD, logA, logC, logP);
	
	void init_ekf3(double Tv);
	void ekf3(double ux, double uy, double isxd, double isyd);
	
	/* Constants - definovat jako konstanty ?? ?kde je vyhodnejsi aby v pameti byli?*/
	int16 Q[9]; /* matrix [4,4] */
	int16 R[4]; /* matrix [2,2] */
	
	int16 x_est[3]; /* estimate and prediction */
	int16 y_est[2]; /* estimate and prediction */
	int16 y_old[2]; /* estimate and prediction */
	
	int16 PSI[9]; /* matrix [4,4] */
	int16 PSIU[9]; /* matrix PIS*U, [4,4] */
	int16 C[6]; /* matrix [4,4] */
	
	int16 Uf[9]; // upper triangular of covariance (inplace)
	int16 Df[3];  // diagonal covariance
	int16 Dfold[3]; // temp of D
	int16 G[9];  // temp for bierman
	
	int16 cA, cB, cC, cG, cF, cH;  // cD, cE, cF, cI ... nepouzivane
	
	enorm<fsqmat> E;
	mat Ry;
	
public:
	//! Default constructor
	EKFfixedUD3 ():BM(),E(),Ry(2,2){
		int16 i;
		for(i=0;i<9;i++){Q[i]=0;}
		for(i=0;i<4;i++){R[i]=0;}
		
		for(i=0;i<3;i++){x_est[i]=0;}
		for(i=0;i<2;i++){y_est[i]=0;}
		for(i=0;i<2;i++){y_old[i]=0;}
		for(i=0;i<9;i++){Uf[i]=0;}
		for(i=0;i<3;i++){Df[i]=0;}
		for(i=0;i<4;i++){G[i]=0;}
		for(i=0;i<3;i++){Dfold[i]=0;}
		
		for(i=0;i<9;i++){PSI[i]=0;}
		for(i=0;i<6;i++){C[i]=0;}
		
		set_dim(3);
		dimc = 2;
		dimy = 2;
		E._mu()=zeros(3);
		E._R()=zeros(3,3);
		init_ekf3(0.000125);
	};
	//! Here dt = [yt;ut] of appropriate dimensions
	void bayes ( const vec &yt, const vec &ut );
	//!dummy!
	const epdf& posterior() const {return E;};
	void log_register(logger &L, const string &prefix){
		BM::log_register ( L, prefix );		
	};
	//void from_setting();
};

UIREGISTER(EKFfixedUD3);

/*!
 * \brief Extended Kalman Filter with Chol matrices in fixed point16 arithmetic
 * 
 * An approximation of the exact Bayesian filter with Gaussian noices and non-linear evolutions of their mean.
 */
class EKFfixedCh : public BM {
public:
	LOG_LEVEL(EKFfixedCh,logCh, logA, logP);
	
	void init_ekf(double Tv);
	void ekf(double ux, double uy, double isxd, double isyd);
	
	/* Constants - definovat jako konstanty ?? ?kde je vyhodnejsi aby v pameti byli?*/
	int16 Q[16]; /* matrix [4,4] */
	int16 R[4]; /* matrix [2,2] */
	
	int16 x_est[4]; /* estimate and prediction */
	
	int16 PSI[16]; /* matrix [4,4] */
	int16 PSICh[16]; /* matrix PIS*U, [4,4] */
	
	int16 Chf[16]; // upper triangular of covariance (inplace)
	
	int16 cA, cB, cC, cG, cH;  // cD, cE, cF, cI ... nepouzivane
	
	enorm<chmat> E;
	mat Ry;
	
public:
	//! Default constructor
	EKFfixedCh ():BM(),E(),Ry(2,2){
		int16 i;
		for(i=0;i<16;i++){Q[i]=0;}
		for(i=0;i<4;i++){R[i]=0;}
		
		for(i=0;i<4;i++){x_est[i]=0;}
		for(i=0;i<16;i++){Chf[i]=0;}
		
		for(i=0;i<16;i++){PSI[i]=0;}
		
		set_dim(4);
		dimc = 2;
		dimy =2;
		E._mu()=zeros(4);
		E._R()=zeros(4,4);
		init_ekf(0.000125);
	};
	//! Here dt = [yt;ut] of appropriate dimensions
	void bayes ( const vec &yt, const vec &ut );
	//!dummy!
	const epdf& posterior() const {return E;};
	void log_register(logger &L, const string &prefix){
		BM::log_register ( L, prefix );
		
		L.add_vector ( log_level, logCh, RV ("Ch", 16 ), prefix );
		L.add_vector ( log_level, logA, RV ("A", 16 ), prefix );
		L.add_vector ( log_level, logP, RV ("P", 16 ), prefix );		
	};
	//void from_setting();
};

UIREGISTER(EKFfixedCh);

/*!
 * \brief Extended Kalman Filter with UD matrices in fixed point16 arithmetic
 * 
 * An approximation of the exact Bayesian filter with Gaussian noices and non-linear evolutions of their mean.
 */
class EKFfixedCh2 : public BM {
public:
	LOG_LEVEL(EKFfixedCh2,logCh, logA, logC, logP, logDet, logRem);
	
	void init_ekf2(double Tv);
	void ekf2(double ux, double uy, double isxd, double isyd);
	
	/* Constants - definovat jako konstanty ?? ?kde je vyhodnejsi aby v pameti byli?*/
	int16 Q[4]; /* matrix [4,4] */
	int16 R[4]; /* matrix [2,2] */
	
	int16 x_est[2]; /* estimate and prediction */
	int16 y_est[2]; /* estimate and prediction */
	int16 y_old[2]; /* estimate and prediction */
	
	int16 PSI[4]; /* matrix [4,4] */
	int16 PSICh[4]; /* matrix PIS*U, [4,4] */
	int16 C[4]; /* matrix [4,4] */
	
	int16 Chf[4]; // upper triangular of covariance (inplace)
	
	int16 cA, cB, cC, cG, cH;  // cD, cE, cF, cI ... nepouzivane
	
	enorm<fsqmat> E;
	mat Ry;
	
public:
	//! Default constructor
	EKFfixedCh2 ():BM(),E(),Ry(2,2){
		int16 i;
		for(i=0;i<4;i++){Q[i]=0;}
		for(i=0;i<4;i++){R[i]=0;}
		
		for(i=0;i<2;i++){x_est[i]=0;}
		for(i=0;i<2;i++){y_est[i]=0;}
		for(i=0;i<2;i++){y_old[i]=0;}
		for(i=0;i<4;i++){Chf[i]=0;}
		
		for(i=0;i<4;i++){PSI[i]=0;}
		for(i=0;i<4;i++){C[i]=0;}
		
		set_dim(2);
		dimc = 2;
		dimy = 2;
		E._mu()=zeros(2);
		E._R()=zeros(2,2);
		init_ekf2(0.000125);
	};
	//! Here dt = [yt;ut] of appropriate dimensions
	void bayes ( const vec &yt, const vec &ut );
	//!dummy!
	const epdf& posterior() const {return E;};
	void log_register(logger &L, const string &prefix){
		BM::log_register ( L, prefix );
		
		L.add_vector ( log_level, logCh, RV ("Ch2", 4 ), prefix );
		L.add_vector ( log_level, logA, RV ("A2", 4 ), prefix );
		L.add_vector ( log_level, logC, RV ("C2", 4 ), prefix );
		L.add_vector ( log_level, logP, RV ("P2", 4 ), prefix );
		L.add_vector ( log_level, logDet, RV ("Det", 1 ), prefix );
		L.add_vector ( log_level, logRem, RV ("Rem", 1 ), prefix );
		
	};
	void from_setting ( const Setting &set ){
		BM::from_setting(set);
		vec dQ,dR;
		UI::get ( dQ, set, "dQ", UI::optional );
		UI::get ( dQ, set, "dR", UI::optional );
		if (dQ.length()==2){
			Q[0]=prevod(dQ[0],15);      // 1e-3
			Q[3]=prevod(dQ[1],15);      // 1e-3
		}
		if (dR.length()==2){
			R[0]=prevod(dR[0],15);      // 1e-3
			R[3]=prevod(dR[1],15);      // 1e-3
		}
	}
};

UIREGISTER(EKFfixedCh2);


//! EKF for comparison of EKF_UD with its fixed-point16 implementation
class EKF_UDfix : public BM {
	protected:
		//! logger
		LOG_LEVEL(EKF_UDfix,logU, logG);
		//! Internal Model f(x,u)
		shared_ptr<diffbifn> pfxu;
		
		//! Observation Model h(x,u)
		shared_ptr<diffbifn> phxu;
		
		//! U part
		mat U;
		//! D part
		vec D;
				
		mat A;
		mat C;
		mat Q;
		vec R;
		
		enorm<ldmat> est;
		
		
	public:
		
		//! copy constructor duplicated 
		EKF_UDfix* _copy() const {
			return new EKF_UDfix(*this);
		}
		
		const enorm<ldmat>& posterior()const{return est;};
		
		enorm<ldmat>& prior() {
			return const_cast<enorm<ldmat>&>(posterior());
		}
		
		EKF_UDfix(){}
		
		
		EKF_UDfix(const EKF_UDfix &E0): pfxu(E0.pfxu),phxu(E0.phxu), U(E0.U), D(E0.D){}
		
		//! Set nonlinear functions for mean values and covariance matrices.
		void set_parameters ( const shared_ptr<diffbifn> &pfxu, const shared_ptr<diffbifn> &phxu, const mat Q0, const vec R0 );
		
		//! Here dt = [yt;ut] of appropriate dimensions
		void bayes ( const vec &yt, const vec &cond = empty_vec );
		
		void log_register ( bdm::logger& L, const string& prefix ){
			BM::log_register ( L, prefix );
			
			if ( log_level[logU] )
				L.add_vector ( log_level, logU, RV ( dimension()*dimension() ), prefix );
			if ( log_level[logG] )
				L.add_vector ( log_level, logG, RV ( dimension()*dimension() ), prefix );
			
		}
		/*! Create object from the following structure
		
		\code
		class = 'EKF_UD';
		OM = configuration of bdm::diffbifn;    % any offspring of diffbifn, bdm::diffbifn::from_setting
		IM = configuration of bdm::diffbifn;    % any offspring of diffbifn, bdm::diffbifn::from_setting
		dQ = [...];                             % vector containing diagonal of Q
		dR = [...];                             % vector containing diagonal of R
		--- optional fields ---
		mu0 = [...];                            % vector of statistics mu0
		dP0 = [...];                            % vector containing diagonal of P0
		-- or --
		P0 = [...];                             % full matrix P0
		--- inherited fields ---
		bdm::BM::from_setting
		\endcode
		If the optional fields are not given, they will be filled as follows:
		\code
		mu0 = [0,0,0,....];                     % empty statistics
		P0 = eye( dim );              
		\endcode
		*/
		void from_setting ( const Setting &set );
		
		void validate() {};
		// TODO dodelat void to_setting( Setting &set ) const;
		
};
UIREGISTER(EKF_UDfix);


class MPF_pmsm_red:public BM{
		double qom, qth, r;


public:
		MPF_pmsm_red(){
			dimy=2;
			dimc=2;
			qom=1e-1;
			qth=1e-6;
			r=1e-1;
		};
	void bayes ( const vec &val, const vec &cond ) {
/*		const double &isa = val(0);
		const double &isb = val(1);
		const double &usa = cond(0);
		const double &usb = cond(1);*/
		mpf_bayes((floatx)val(0),(floatx)val(1),(floatx)cond(0), (floatx)cond(1));//isa,isb,usa,usb);
	}
	
	class mp:public epdf{
		LOG_LEVEL(mp,logth);
	public:
		mp():epdf(){set_dim(3);		log_level[logth]=true;
}
		vec sample() const {return zeros(3);}
		double evallog(const vec &v) const {return 0.0;}
		vec mean() const {vec tmp(3); 
			floatx es,ec,eo;
			mpf_mean(&es, &ec, &eo); 
			tmp(0)=es;tmp(1)=ec;tmp(2)=eo;
			return tmp;
		}
		vec variance() const {return zeros(3);}
		void log_register ( bdm::logger& L, const string& prefix ) {
			epdf::log_register ( L, prefix );
			if ( log_level[logth] ) {
				int th_dim = N; // dimension - dimension of cov
				L.add_vector( log_level, logth, RV ( th_dim ), prefix );
			}
		}
		void log_write() const {
			epdf::log_write();
			if ( log_level[logth] ) {
				floatx Th[N];
				mpf_th(Th);
				vec th(N); 
				for(int i=0;i<N;i++){th(i)=Th[i];}
				log_level.store( logth, th );
			}
		}
	};
	
	mp mypdf;
	const mp& posterior() const {return mypdf;}
	
	void from_setting(const Setting &set){
	    BM::from_setting(set);
		UI::get ( log_level, set, "log_level", UI::optional );

		UI::get(qom,set,"qom",UI::optional);
		UI::get(qth,set,"qth",UI::optional);
		UI::get(r,set,"r",UI::optional);
	}
	void validate(){
		mpf_init((floatx)qom,(floatx)qth,(floatx)r);
		
	}
};
UIREGISTER(MPF_pmsm_red);

//! EKF with covariance R estimated by the VB method
class EKFvbR: public EKFfull{	
    LOG_LEVEL(EKFvbR,logR,logQ,logP);

	//! Statistics of the R estimator
	mat PsiR;
	//! Statistics of the Q estimator
	mat PsiQ;
	//! forgetting factor
	double phi;
	//! number of VB iterations
	int niter;
	//! degrees of freedom
	double nu;
	//! stabilizing element
	mat PsiR0;
	//! stabilizing element
	mat PsiQ0;
	
	void from_setting(const Setting &set){
	    EKFfull::from_setting(set);
		if (!UI::get(phi,set,"phi",UI::optional)){
			phi = 0.99;
		}
		if (!UI::get(niter,set,"niter",UI::optional)){
			niter = 3;
		}
		PsiQ = Q;
		PsiQ0 = Q;
		PsiR = R;
		PsiR0 = R;
		nu = 3;
    }	
	void log_register ( logger &L, const string &prefix ) {
		EKFfull::log_register(L,prefix);
		L.add_vector(log_level, logR, RV("{R }",vec_1<int>(4)), prefix);
		L.add_vector(log_level, logQ, RV("{Q }",vec_1<int>(4)), prefix);
		L.add_vector(log_level, logP, RV("{P }",vec_1<int>(4)), prefix);
	}
	void bayes ( const vec &val, const vec &cond ) {
		vec diffx, diffy;
		mat Psi_vbQ;
		mat Psi_vbR;
		nu = phi*nu + (1-phi)*2 + 1.0;
		
		//save initial values of posterior
		vec mu0=est._mu();
		fsqmat P0=est._R();
		vec xpred = pfxu->eval(mu0,cond);
		
		for (int i=0; i<niter; i++){
			est._mu() = mu0;
			est._R() = P0;
			
			EKFfull::bayes(val,cond);

//			diffy = val - fy._mu();
//			Psi_vbR = phi*PsiR + (1-phi)*PsiR0+ outer_product(diffy,diffy)/*+C*mat(est._R())*C.T()*/;
//			R = Psi_vbR/(nu-2);
	
			diffx = est._mu() - xpred;
			Psi_vbQ = phi*PsiQ + (1-phi)*PsiQ0+ outer_product(diffx,diffx);/*mat(est._R());  A*mat(P0)*A.T();*/
			//
			Psi_vbQ(0,1) = 0.0;
			Psi_vbQ(1,0) = 0.0;
			Q = Psi_vbQ/(nu-2);
		}
		PsiQ = Psi_vbQ;
		PsiR = Psi_vbR;
//		cout <<"==" << endl << Psi << endl << diff << endl << P0 << endl << ":" << Q;
		log_level.store(logQ, vec(Q._M()._data(),4));
		log_level.store(logR, vec(R._M()._data(),4));
		{
			mat Ch(2,2);
			Ch=chol(est._R()._M());
			log_level.store(logP, vec(Ch._data(),4));
		}
	}
};
UIREGISTER(EKFvbR);


class ekfChfix: public BM{
	ekf_data E;
public:
	ekfChfix(){
		init_ekfCh2(&E,0.000125);set_dim(2);	dimc = 2;
		dimy = 2;
	} 
	void bayes ( const vec &val, const vec &cond ) {
		int16 ux,uy;
		ux=prevod(cond[0]/Uref,15);
		uy=prevod(cond[1]/Uref,15);
		
		int16 yx,yy;
		// zadani mereni
		yx=prevod(val[0]/Iref,15);
		yy=prevod(val[1]/Iref,15);
		
		int16 detS, rem;
		ekfCh2(&E, ux,uy,yx,yy, &detS, &rem);
		
		Est._mu()=vec_2(E.x_est[0]*Uref/32768., E.x_est[1]*Uref/32768.);
		
		ll = -0.5*qlog(detS)-0.5*rem;
	}
	const epdf& posterior() const {return Est;};
	void from_setting ( const Setting &set ){
		BM::from_setting(set);
		vec dQ,dR;
		UI::get ( dQ, set, "dQ", UI::optional );
		UI::get ( dR, set, "dR", UI::optional );
		if (dQ.length()==2){
			E.Q[0]=prevod(dQ[0],15);      // 1e-3
			E.Q[3]=prevod(dQ[1],15);      // 1e-3
		}
		if (dR.length()==2){
			E.dR[0]=prevod(dR[0],15);      // 1e-3
			E.dR[1]=prevod(dR[1],15);      // 1e-3
		}
	}
	
	enorm<fsqmat> Est;
	mat Ry;
};
UIREGISTER(ekfChfix);



class ekfChfixQ: public BM{
	LOG_LEVEL(ekfChfixQ,logQ,logCh,logC,logRes)
	ekf_data E;
	int64 PSI_Q0, PSI_Q1, PSI_Q0_reg, PSI_Q1_reg;
	int Q_ni;
	int phi_Q;
public:
	ekfChfixQ(){
		init_ekfCh2(&E,0.000125);set_dim(2);	dimc = 2;
		dimy = 2;
		Q_ni = 7;
		phi_Q = ((1<<Q_ni) -1)<<(15-Q_ni);

		PSI_Q0_reg = ((int64)((1<<15)-phi_Q)*E.Q[0])<<Q_ni;
		PSI_Q1_reg = ((int64)((1<<15)-phi_Q)*E.Q[3])<<Q_ni;
		PSI_Q0 = ((int64)E.Q[0])<<Q_ni;
		PSI_Q1 = ((int64)E.Q[3])<<Q_ni;
	} 
	void bayes ( const vec &val, const vec &cond ) {
		int16 ux,uy;
		ux=prevod(cond[0]/Uref,15);
		uy=prevod(cond[1]/Uref,15);
		
		int16 yx,yy;
		// zadani mereni
		yx=prevod(val[0]/Iref,15);
		yy=prevod(val[1]/Iref,15);
		
		int16 detS, rem;
		ekfCh2(&E, ux,uy,yx,yy, &detS, &rem);
		
		Est._mu()=vec_2(E.x_est[0]*Uref/32768., E.x_est[1]*Uref/32768.);
		
		int16 xerr0 = E.x_est[0]-E.x_pred[0];
		PSI_Q0 = ((int64)(phi_Q*PSI_Q0) + (((int64)xerr0*xerr0)<<15) + PSI_Q0_reg)>>15;		
		E.Q[0] = PSI_Q0>>Q_ni;

		xerr0 = E.x_est[1]-E.x_pred[1];
		PSI_Q1 = ((int64)(phi_Q*PSI_Q1) + (((int64)xerr0*xerr0)<<15) + PSI_Q1_reg)>>15;		
		E.Q[3] = PSI_Q1>>Q_ni;
		
		cout << E.Q[0] << "," << E.Q[3] <<endl;
		ll = -0.5*qlog(detS)-0.5*rem;
	}
	const epdf& posterior() const {return Est;};
	void from_setting ( const Setting &set ){
		BM::from_setting(set);
		vec dQ,dR;
		UI::get ( dQ, set, "dQ", UI::optional );
		UI::get ( dR, set, "dR", UI::optional );
		if (dQ.length()==2){
			E.Q[0]=prevod(dQ[0],15);      // 1e-3
			E.Q[3]=prevod(dQ[1],15);      // 1e-3
		}
		if (dR.length()==2){
			E.dR[0]=prevod(dR[0],15);      // 1e-3
			E.dR[1]=prevod(dR[1],15);      // 1e-3
		}

		UI::get(Q_ni, set, "Q_ni", UI::optional);
		{ // zmena Q!!
			phi_Q = ((1<<Q_ni) -1)<<(15-Q_ni);

			PSI_Q0_reg = ((int64)((1<<15)-phi_Q)*E.Q[0])<<Q_ni;
			PSI_Q1_reg = ((int64)((1<<15)-phi_Q)*E.Q[3])<<Q_ni;
			PSI_Q0 = ((int64)E.Q[0])<<Q_ni;
			PSI_Q1 = ((int64)E.Q[3])<<Q_ni;

		}
		UI::get(log_level, set, "log_level", UI::optional);

	}
	void log_register ( logger &L, const string &prefix ) {
		BM::log_register(L,prefix);
		L.add_vector(log_level, logQ, RV("{Q }",vec_1<int>(2)), prefix);
		L.add_vector(log_level, logCh, RV("{Ch }",vec_1<int>(3)), prefix);
		L.add_vector(log_level, logC, RV("{C }",vec_1<int>(4)), prefix);
		L.add_vector(log_level, logRes, RV("{dy }",vec_1<int>(2)), prefix);
	}
	void log_write() const {
		BM::log_write();
		if ( log_level[logQ] ) {
			log_level.store( logQ, vec_2((double)E.Q[0],(double)E.Q[3] ));
		}
		if ( log_level[logCh] ) {
			log_level.store( logCh, vec_3((double)E.Chf[0],(double)E.Chf[1], (double)E.Chf[3]));
		}
		if ( log_level[logRes] ) {
			log_level.store( logRes, vec_2((double)E.difz[0],(double)E.difz[1]));
		}
		if ( log_level[logC] ) {
			vec v(4);
			for (int i=0;i<4;i++) v(i)=E.C[i];
			log_level.store( logC, v);
		}
	}
	
	enorm<fsqmat> Est;
	mat Ry;
};
UIREGISTER(ekfChfixQ);

#endif // KF_H


