Index: /applications/pmsm/Barcelona/QandR.m
===================================================================
--- /applications/pmsm/Barcelona/QandR.m (revision 342)
+++ /applications/pmsm/Barcelona/QandR.m (revision 349)
@@ -23,5 +23,5 @@
     1e-7];
 
-Q= StdQ.^2
+Qu = (dt/Ls)^2*(2*3)^2/12;
 
 %R
Index: /applications/pmsm/Barcelona/real.cfg
===================================================================
--- /applications/pmsm/Barcelona/real.cfg (revision 346)
+++ /applications/pmsm/Barcelona/real.cfg (revision 349)
@@ -2,10 +2,10 @@
 
 experiment={
-	ndat = 32000;
+	ndat = 8000;
 };
 
 system= {
 	type = "ITppFileDS";
-	filename = "tec0013ALL.it";
+	filename = "tec0009ALL.it";
 	varname="D";
 	rv = {type="rv"; names=["o_ia","o_ib","o_ua","o_ub"];};  // observed u
@@ -23,4 +23,5 @@
 	 dR=[0.025, 0.025];
 	 dP0=[1., 1., 1., 10., 1e-10];
+	 mu0=[1., 1., 1., 2.4, 1e-10];
 	 //connect
 	 drv = {type="rv"; names=["o_ia","o_ib","o_ua","o_ub"];};  // observed u
@@ -33,5 +34,5 @@
 logger= {
 	type="itpplog";
-	maxlen = 32000;
+	maxlen = 8000;
 	filename = "real.it";
 };
Index: /applications/pmsm/Barcelona/Mz_mx.m
===================================================================
--- /applications/pmsm/Barcelona/Mz_mx.m (revision 346)
+++ /applications/pmsm/Barcelona/Mz_mx.m (revision 349)
@@ -10,2 +10,12 @@
 hold on
 plot(M.crb(:,4),'r');
+
+figure(2)
+hold off
+plot(M.crbz(:,3));
+hold on
+plot(M.crb(:,3),'r');
+
+ind=100:size(M.crb,1)-100;
+[min(M.crb(ind,3)), max(M.crb(ind,3))]
+[min(M.crbz(ind,3)), max(M.crbz(ind,3))]
Index: /applications/pmsm/Barcelona/pwm_mx.cfg
===================================================================
--- /applications/pmsm/Barcelona/pwm_mx.cfg (revision 349)
+++ /applications/pmsm/Barcelona/pwm_mx.cfg (revision 349)
@@ -0,0 +1,35 @@
+// settings for pmsm_estim.cpp
+
+experiment={
+	ndat = 60000;
+};
+
+system= {
+	type = "pmsmDS";
+	params = {type = "external";filename = "../cfg/zcu.cfg";path = "pmsm107";};
+	tstep = 1.0; // steps for profile in [s]
+	profileW = [1, 5, 10, 100, 3, 3, 0, 0];
+	//profileM = [0, 0, 10,  10, 10, 0, 10, 0];
+	options = "modelu";
+};
+
+estimator = (
+	{type="pmsmCRB";
+	IM = {type="pmsmIM";
+	  params = {type = "external";filename = "../cfg/zcu.cfg";path = "pmsm107";};
+	  variant = "2o";// 2nd order Taylor
+	};
+	OM = {type="pmsmOM";};
+	dQ=[0.017, 0.017, 0.00000004, 1e-12];
+	dR=[0.0000000001, 0.000000001];
+	//connect
+	drv = {type="rv"; names=["o_ia","o_ib","o_ua","o_ub"];};  // observed u
+	//drv = {type="rv"; names=["o_ia","o_ib","t_ua","t_ub"];}; // true u
+	rv = {type="rv"; names=["EKF"]; sizes=[4];};
+	}
+);
+	
+logger= {
+	type="mexlog";
+	maxlen = 59999;
+};
Index: /applications/pmsm/Barcelona/real_mx.cfg
===================================================================
--- /applications/pmsm/Barcelona/real_mx.cfg (revision 346)
+++ /applications/pmsm/Barcelona/real_mx.cfg (revision 349)
@@ -2,10 +2,10 @@
 
 experiment={
-	ndat = 32000;
+	ndat = 8000;//32000;
 };
 
 system= {
 	type = "ITppFileDS";
-	filename = "tec0013ALL.it";
+	filename = "tec0009ALL.it";
 	varname="D";
 	rv = {type="rv"; names=["o_ia","o_ib","o_ua","o_ub"];};  // observed u
@@ -19,7 +19,7 @@
 	 };
 	 OM = {type="pmsmOM";};
-//	 dQ=[0.017, 0.017, 0.00000004, 1e-14, 1e-30];
-	 dQ=[0.017, 0.017, 0.00000007, 1e-14, 1e-3];
-	 dR=[0.025, 0.025];
+	 dQ=[6e-3, 6e-3, 5e-7, 3e-13, 1e-2];
+	 dR=[0.0025, 0.0025];
+	 mu0=[0., 0., 0.,3.4, 1e-10];
 	 dP0=[1., 1., 1., 10., 1e-10];
 	 //connect
@@ -31,10 +31,10 @@
 	 IM = {type="pmsmIM";
 		params = {type = "external";filename = "../cfg/zcu.cfg";path = "pmsm107";};
-		variant = "Stat";
+		variant = "2o";
 	 };
 	 OM = {type="pmsmOM";};
-//	 dQ=[0.017, 0.017, 0.00000004, 1e-14, 1e-30];
-	 dQ=[0.017, 0.017, 1.0, 1e-8];
-	 dR=[0.025, 0.025];
+	 dQ=[6e-3, 6e-3, 0.3, 5e-9];
+	 dR=[0.0025, 0.0025];
+	 mu0=[0., 0., 0.,3.4];
 	 dP0=[1., 1., 1., 10.];
 	 //connect
@@ -47,4 +47,4 @@
 logger= {
 	type="mexlog";
-	maxlen = 31999;
+	maxlen = 7999;
 };
Index: /applications/pmsm/Barcelona/Mz_mx.cfg
===================================================================
--- /applications/pmsm/Barcelona/Mz_mx.cfg (revision 346)
+++ /applications/pmsm/Barcelona/Mz_mx.cfg (revision 349)
@@ -9,7 +9,7 @@
 	params = {type = "external";filename = "../cfg/zcu.cfg";path = "pmsm107";};
 	tstep = 1.0; // steps for profile in [s]
-	profileW = [0, 1, 5, 10, 10, 3, 0, 0, -3, -3, -10, -10, -4, 0, 0];
-	profileM = [0, 10, 0, 0, 10, 0, 0, 0, 0, 0, 10, 0, 0];
-	options = "modelu";
+	profileW = [0, 3, 10, 10, 3, 0, 0, -3, -3, -3, -10, -10, -4, 0, 0];
+	profileM = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0];
+//	options = "modelu";
 };
 
@@ -26,5 +26,5 @@
 	cols = 4;
 	};*/
-	dQ=[0.0017, 0.0017, 0.00000004, 1e-14, 1e-1];
+	dQ=[0.017, 0.017, 0.00000004, 1e-14, 1e-1];
 	 dR=[0.025, 0.025];
 	 dP0=[0., 0., 0., 0., 0.];
Index: /applications/pmsm/Barcelona/real_pmsm.m
===================================================================
--- /applications/pmsm/Barcelona/real_pmsm.m (revision 346)
+++ /applications/pmsm/Barcelona/real_pmsm.m (revision 349)
@@ -1,5 +1,5 @@
 M=pmsm_estim_mex('real_mx.cfg');
 
-ids = [1:size(M.EKFMz_mean,1)];
+ids = [1:size(M.mean_EKFMz,1)];
 kftime=[ids-ids(1)]*0.000125;
 
@@ -7,20 +7,29 @@
 subplot(2,1,1);
 hold off
-plot(kftime,M.EKFMz_mean(ids,3)/2/pi);
+plot(kftime,M.mean_EKFMz(ids,3)/2/pi,'-.');
 hold on
-plot(kftime,M.EKF_mean(ids,3)/2/pi,'r');
+plot(kftime,M.mean_EKF(ids,3)/2/pi,'r');
+
+title('Measured and estimated \omega_{me}');
+xlabel('time [s]')
 
 subplot(2,1,2);
-plot(kftime, M.EKFMz_mean(ids,5));
+plot(kftime, M.mean_EKFMz(ids,5));
+title('Estimated load T_L');
+
 
 % ----- 
 
-load ~/work/vs/Pohony/DATA_MERENI/DATA_MCU/data1.dat -ascii
+load ~/work/vs/Pohony/DATA_MERENI/DATA_MCU/data9.dat -ascii
+
+data=data9;
 
 subplot(2,1,1);
 setime=[0:799]*0.004666; % fitted on currents
-plot(setime,data13([0:799]+711,6),'--k','LineWidth',3);
-plot(setime,data13([0:799]+711,11),'--g','LineWidth',3);
+plot(setime,data([0:799]+711,6),'--k','LineWidth',3);
+%plot(setime,data([0:799]+711,11),'--g','LineWidth',3);
 %plot(setime,2*pi*data1([0:799]+620,6),'--k','LineWidth',3);
 
+legend('EKF T_L estimate','EKF estimate','Sensor readings')
 
+set(gca,'XLim',[0,kftime(end)])
Index: /applications/pmsm/Barcelona/pwm.m
===================================================================
--- /applications/pmsm/Barcelona/pwm.m (revision 349)
+++ /applications/pmsm/Barcelona/pwm.m (revision 349)
@@ -0,0 +1,84 @@
+function pwm(compute)
+if nargin<1, 
+    compute=1;
+end
+
+if compute
+    M=pmsm_estim_mex('pwm_mx.cfg');
+    save pwm M
+else
+    load pwm
+end
+
+time =0.000125*[0:size(M.xia,1)-1];
+
+%%%%%%%%%%%%%% U
+figure(1);
+detail=14000:20000;
+% ua
+subplot(2,2,1);
+plot(time(detail),M.oua(detail)-M.tua(detail));
+title('Difference between simulated and non-compensated observed u_{s,\alpha} (detail)')
+
+subplot(2,2,2);
+hist(M.oua-M.tua,30);
+title('Histogram of difference on the full data record of 8s')
+
+%ub
+subplot(2,2,3);
+plot(time(detail),M.modelub(detail)-M.tub(detail));
+title('Difference between simulated and compensated observed u_{s,\beta} (detail)')
+
+subplot(2,2,4);
+hist(M.modelub-M.tub,30);
+title('Histogram of difference on the full data record of 8s')
+
+dt=0.000125;
+Ls=0.003465; 
+
+bnds = [max(abs(M.modelua-M.tua)), max(abs(M.modelub-M.tub)) 0 0];
+Qpwm = (dt/Ls)^2*(2*bnds).^2/12
+
+%%%%%%%%%%%%%%% Integration
+
+figure(2)
+subplot(4,2,1); plot(time,M.erria);xlabel('time [s]');
+title('Integration error on state variables')
+ylabel('Error on i_{s,\alpha}')
+subplot(4,2,2); hist(M.erria,30);
+title('Histogram of errors on state variables')
+
+subplot(4,2,3); plot(time,M.errib);xlabel('time [s]');
+ylabel('Error on i_{s,\beta}')
+subplot(4,2,4); hist(M.errib,30);
+%title('Histogram of error on i_{s,\beta}')
+
+subplot(4,2,5); plot(time,M.errom);xlabel('time [s]');
+ylabel('Error on \omega_{me}')
+subplot(4,2,6); hist(M.errom,30);
+%title('Histogram of error on \omega_{me}')
+
+subplot(4,2,7); plot(time,M.errth); xlabel('time [s]');
+ylabel('Error on \vartheta_{me}')
+subplot(4,2,8); hist(M.errth,30);
+%title('Histogram of error on \vartheta_{me}')
+
+Qint=[var(M.erria), var(M.errib), var(M.errom), var(M.errth)]
+
+bnds=[max(abs(M.erria)), max(abs(M.errib)), max(abs(M.errom)), max(abs(M.errth))]
+Qint = [(2*bnds).^2/12]
+
+%%%%%%%%%%%%%%%% Mz prop
+Fmag= 0.1989; 
+p = 4.0;
+J = 0.04;   
+
+Tl=10;
+dom = p/J*Tl*dt
+Qrat=(2*[Fmag/Ls*dt,Fmag/Ls*dt, 1, dt]).^2/12
+
+Qtot=Qpwm+Qint
+Qtot+Qrat
+
+%%%%%%%%%%%%%%%% R
+R=(2*0.085)^2/12*[1 1]
Index: /applications/pmsm/pmsmDS.h
===================================================================
--- /applications/pmsm/pmsmDS.h (revision 342)
+++ /applications/pmsm/pmsmDS.h (revision 349)
@@ -13,5 +13,5 @@
 #include <stat/loggers.h>
 #include <estim/libKF.h>
-#include "simulator.h"
+#include "simulator_zdenek/simulator.h"
 #include "pmsm.h"
 
@@ -149,5 +149,6 @@
 		vec true_state=vec(x,4); // read from pmsm
 		E.set_mu(true_state);
-		
+		mu=true_state;
+			
 		//integration error
 		old_true(4)=KalmanObs[4];
@@ -199,5 +200,5 @@
 	
 		void log_add(logger &L, const string &name="" ){
-			L_CRB=L.add(concat(rx,RV("Mz",1,0)),"crb");
+			L_CRB=L.add(concat(rx,RV("Mz",1,0)),"crbz");
 		}
 		void logit(logger &L){
Index: /applications/pmsm/pmsm.h
===================================================================
--- /applications/pmsm/pmsm.h (revision 332)
+++ /applications/pmsm/pmsm.h (revision 349)
@@ -30,4 +30,7 @@
 
 	void modelpwm(const vec &x0, const vec u0, double &ua, double &ub){
+/*		ua=u0[0];
+		ub=u0[1];
+		return;*/
 		double sq3=sqrt ( 3.0 );
 		double i1=x0(0);
Index: /applications/pmsm/bm_pmsm.sh
===================================================================
--- /applications/pmsm/bm_pmsm.sh (revision 349)
+++ /applications/pmsm/bm_pmsm.sh (revision 349)
@@ -0,0 +1,2 @@
+#!/bin/sh
+mex -l itpp -l bdm -l blas -l lapack -lpmsmsim -I/home/smidl/bin/matlab2008a/extern/include -I/home/smidl/work/git/mixpp/bdm -I/home/smidl/work/git/mixpp/system/linux -L/home/smidl/work/git/mixpp/bdm -L/home/smidl/work/git/mixpp/system/linux/itpp/.libs -L /usr/lib64 -L simulator_zdenek $1 
Index: /applications/pmsm/pmsm_estim_mex.cpp
===================================================================
--- /applications/pmsm/pmsm_estim_mex.cpp (revision 349)
+++ /applications/pmsm/pmsm_estim_mex.cpp (revision 349)
@@ -0,0 +1,100 @@
+
+/*!
+\file
+\brief Multi-Estimator (developped for PMSM)
+
+
+ */
+
+#include "pmsm_ui.h"
+#include <stat/loggers_ui.h>
+#include <estim/KF_ui.h>
+
+#include <itpp/itmex.h>
+#include "../../library/mex/mexlog.h"
+#include "../../library/mex/mexds.h"
+#include "../../library/mex/mexparse.h"
+
+using namespace bdm;
+
+void mexFunction ( int n_output, mxArray *output[], int n_input, const mxArray *input[] ) {
+	// Check the number of inputs and output arguments
+	string fname;
+	if (n_input>0){
+		fname=mxArray2string(input[0]);
+	} else {
+		printf("file name missing");
+	}
+	// ------------------ 
+	
+	UIFile F ( fname.c_str());
+
+	logger* L;
+	DS * DS; 
+	Array<BM*> Es;			// array of estimators
+	int Ndat;				//number of data records
+	int nE;					//number of estimators
+
+	try {
+		UIbuild ( F.lookup ( "logger" ),L );
+		UIbuild ( F.lookup ( "system" ),DS );
+		F.lookupValue ( "experiment.ndat",Ndat );
+		Setting& S=F.lookup ( "estimator" );
+		nE = S.getLength();
+		Es.set_size(nE);
+		for(int i=0;i<nE;i++){		
+			UIbuild (S[i] ,Es(i) );
+		}
+	}
+	catch UICATCH;
+
+	DS->log_add ( *L );
+	string nic="";
+	for (int i=0; i<nE; i++){
+		Es(i)->log_add(*L,nic); // estimate
+	}
+	L->init();
+
+	vec dt=zeros ( DS->_drv()._dsize() );   //data variable
+	Array<datalink*> Dls(nE); 
+	for (int i=0; i<nE; i++){
+		Dls(i)=new datalink( Es(i)->_drv(),DS->_drv() ); //datalink between a datasource and estimator
+	}
+	
+	// Main cycle
+	for ( int tK=1;tK<Ndat;tK++ ) {
+		// Data Source
+		DS->step();							// simulator step
+		DS->getdata ( dt );					// read data
+		DS->logit ( *L );
+		
+		// Estimators
+		for (int i=0; i<nE; i++){
+			Es(i)->bayes ( Dls(i)->pushdown ( dt ) );		// update estimates
+
+			Es(i)->logit (*L);
+		}
+		// Regulators
+		L->step();
+	}
+
+	L->finalize();
+
+	// ------------------ End of routine -----------------------------
+	
+	mexlog* mL=dynamic_cast<mexlog*>(L);
+
+	if (mL) { // user wants output!!
+		if ( n_output<1 ) mexErrMsgTxt ( "Wrong number of output variables!" );
+		output[0] = mL->toCell();
+	}
+	
+	
+	delete L;
+	delete DS;
+	for (int i=0; i<nE; i++){
+		delete Dls(i);
+		delete Es(i);
+	}
+}
+
