Index: /library/bdm/estim/particles.cpp
===================================================================
--- /library/bdm/estim/particles.cpp (revision 384)
+++ /library/bdm/estim/particles.cpp (revision 461)
@@ -14,5 +14,5 @@
 		//generate new samples from paramater evolution model;
 		_samples ( i ) = par->samplecond ( _samples ( i ) );
-		lls ( i ) = par->_e()->evallog ( _samples ( i ) );
+		lls ( i ) = par->e()->evallog ( _samples ( i ) );
 		lls ( i ) *= obs->evallogcond ( dt,_samples ( i ) );
 
Index: /library/bdm/estim/particles.h
===================================================================
--- /library/bdm/estim/particles.h (revision 384)
+++ /library/bdm/estim/particles.h (revision 461)
@@ -216,5 +216,5 @@
 		//generate new samples from paramater evolution model;
 		_samples ( i ) = par->samplecond ( _samples ( i ) );
-		llsP ( i ) = par->_e()->evallog ( _samples ( i ) );
+		llsP ( i ) = par->e()->evallog ( _samples ( i ) );
 		BMs ( i )->condition ( _samples ( i ) );
 		BMs ( i )->bayes ( dt );
Index: /library/bdm/shared_ptr.h
===================================================================
--- /library/bdm/shared_ptr.h (revision 420)
+++ /library/bdm/shared_ptr.h (revision 461)
@@ -28,4 +28,6 @@
 template <typename T>
 class shared_ptr {
+    template<class U> friend class shared_ptr;
+
 private:
     T *payload;
@@ -51,5 +53,15 @@
     //! If other is empty, constructs an empty shared_ptr; otherwise,
     //! constructs a shared_ptr that shares ownership with other.
-    shared_ptr(const shared_ptr &other):
+    shared_ptr(const shared_ptr<T> &other):
+	payload(other.payload),
+	refCnt(other.refCnt)
+    {
+	add_ref();
+    }
+      
+    //! If other is empty, constructs an empty shared_ptr; otherwise,
+    //! constructs a shared_ptr that shares ownership with other.
+    template<typename U>
+    shared_ptr(const shared_ptr<U> &other):
 	payload(other.payload),
 	refCnt(other.refCnt)
Index: /library/bdm/base/bdmbase.h
===================================================================
--- /library/bdm/base/bdmbase.h (revision 436)
+++ /library/bdm/base/bdmbase.h (revision 461)
@@ -18,6 +18,6 @@
 #include "../itpp_ext.h"
 #include "../bdmroot.h"
+#include "../shared_ptr.h"
 #include "user_info.h"
-
 
 using namespace libconfig;
@@ -350,86 +350,68 @@
   //! random variable in condition
   RV rvc;
+
+private:
   //! pointer to internal epdf
-  epdf* ep;
+  shared_ptr<epdf> shep;
+
 public:
   //! \name Constructors
   //! @{
 
-  mpdf() : dimc(0), rvc() {};
-  //! copy constructor does not set pointer \c ep - has to be done in offsprings!
-  mpdf(const mpdf &m) : dimc(m.dimc), rvc(m.rvc) {};
-  //!@}
-
-  //! \name Matematical operations
-  //!@{
-
-  //! Returns a sample from the density conditioned on \c cond, \f$x \sim epdf(rv|cond)\f$. \param cond is numeric value of \c rv
-  virtual vec samplecond(const vec &cond) {
-    this->condition(cond);
-    vec temp = ep->sample();
-    return temp;
-  };
-  //! Returns \param N samples from the density conditioned on \c cond, \f$x \sim epdf(rv|cond)\f$. \param cond is numeric value of \c rv
-  virtual mat samplecond_m(const vec &cond, int N) {
-    this->condition(cond);
-    mat temp(ep->dimension(), N);
-    vec smp(ep->dimension());
-    for (int i = 0; i < N; i++) {smp = ep->sample() ; temp.set_col(i, smp);}
-    return temp;
-  };
-  //! Update \c ep so that it represents this mpdf conditioned on \c rvc = cond
-  virtual void condition(const vec &cond) {it_error("Not implemented");};
+    mpdf():dimc(0), rvc() { }
+
+    mpdf(const mpdf &m):dimc(m.dimc), rvc(m.rvc), shep(m.shep) { }
+    //!@}
+
+    //! \name Matematical operations
+    //!@{
+
+    //! Returns a sample from the density conditioned on \c cond, \f$x \sim epdf(rv|cond)\f$. \param cond is numeric value of \c rv
+    virtual vec samplecond(const vec &cond);
+
+    //! Returns \param N samples from the density conditioned on \c cond, \f$x \sim epdf(rv|cond)\f$. \param cond is numeric value of \c rv
+    virtual mat samplecond_m(const vec &cond, int N);
+
+    //! Update \c ep so that it represents this mpdf conditioned on \c rvc = cond
+    virtual void condition(const vec &cond) {it_error("Not implemented");};
 
   //! Shortcut for conditioning and evaluation of the internal epdf. In some cases,  this operation can be implemented efficiently.
-  virtual double evallogcond(const vec &dt, const vec &cond) {
-    double tmp;
-    this->condition(cond);
-    tmp = ep->evallog(dt);
- //   it_assert_debug(std::isfinite(tmp), "Infinite value");
-    return tmp;
-  };
-
-  //! Matrix version of evallogcond
-  virtual vec evallogcond_m(const mat &Dt, const vec &cond) {this->condition(cond); return ep->evallog_m(Dt);};
-  //! Array<vec> version of evallogcond
-  virtual vec evallogcond_m(const Array<vec> &Dt, const vec &cond) {this->condition(cond); return ep->evallog_m(Dt);};
-
-  //! \name Access to attributes
-  //! @{
-
-  RV _rv() {return ep->_rv();}
-  RV _rvc() {it_assert_debug(isnamed(), ""); return rvc;}
-  int dimension() {return ep->dimension();}
-  int dimensionc() {return dimc;}
-  epdf& _epdf() {return *ep;}
-  epdf* _e() {return ep;}
-  //! Load from structure with elements:
-  //!  \code
-  //! { rv = {class="RV", names=(...),}; // RV describing meaning of random variable
-  //!   rvc= {class="RV", names=(...),}; // RV describing meaning of random variable in condition
-  //!   // elements of offsprings
-  //! }
-  //! \endcode
-  //!@}
-  void from_setting(const Setting &set){
-	  if (set.exists("rv")){
-		  RV* r = UI::build<RV>(set,"rv");
-		  set_rv(*r); 
-		  delete r;
-	  }
-	  if (set.exists("rvc")){
-		  RV* r = UI::build<RV>(set,"rvc");
-		  set_rvc(*r); 
-		  delete r;
-	  }
-  }
-  //!@}
-
-  //! \name Connection to other objects
-  //!@{
-  void set_rvc(const RV &rvc0) {rvc = rvc0;}
-  void set_rv(const RV &rv0) {ep->set_rv(rv0);}
-  bool isnamed() {return (ep->isnamed()) && (dimc == rvc._dsize());}
-  //!@}
+    virtual double evallogcond(const vec &dt, const vec &cond);
+
+    //! Matrix version of evallogcond
+    virtual vec evallogcond_m(const mat &Dt, const vec &cond);
+
+    //! Array<vec> version of evallogcond
+    virtual vec evallogcond_m(const Array<vec> &Dt, const vec &cond);
+
+    //! \name Access to attributes
+    //! @{
+
+    RV _rv() { return shep->_rv(); }
+    RV _rvc() { it_assert_debug(isnamed(), ""); return rvc; }
+    int dimension() { return shep->dimension(); }
+    int dimensionc() { return dimc; }
+
+    epdf *e() { return shep.get(); }
+
+    void set_ep(shared_ptr<epdf> ep) { shep = ep; }
+
+    //! Load from structure with elements:
+    //!  \code
+    //! { rv = {class="RV", names=(...),}; // RV describing meaning of random variable
+    //!   rvc= {class="RV", names=(...),}; // RV describing meaning of random variable in condition
+    //!   // elements of offsprings
+    //! }
+    //! \endcode
+    //!@}
+    void from_setting(const Setting &set);
+    //!@}
+
+    //! \name Connection to other objects
+    //!@{
+    void set_rvc(const RV &rvc0) { rvc = rvc0; }
+    void set_rv(const RV &rv0) { shep->set_rv(rv0); }
+    bool isnamed() { return (shep->isnamed()) && (dimc == rvc._dsize()); }
+    //!@}
 };
 
@@ -628,12 +610,12 @@
 */
 class mepdf : public mpdf {
-	bool owning_ep; // flag trigers deleting ep by destructor
 public:
 	//!Default constructor
-	mepdf(){};
-	mepdf ( epdf* em, bool owning_ep0=false ) :mpdf ( ) {ep= em ;owning_ep=owning_ep0;dimc=0;};
-	mepdf (const epdf* em ) :mpdf ( ) {ep=const_cast<epdf*>( em );};
+	mepdf() { }
+
+	mepdf(shared_ptr<epdf> em) :mpdf ( ) { set_ep(em); dimc=0; }
+
 	void condition ( const vec &cond ) {}
-	~mepdf(){if (owning_ep) delete ep;}
+
   //! Load from structure with elements:
   //!  \code
@@ -644,7 +626,6 @@
   //!@}
 	void from_setting(const Setting &set){
-		epdf* e = UI::build<epdf>(set,"epdf");
-		ep=	e; 
-		owning_ep=true;
+		epdf *e = UI::build<epdf>(set, "epdf");
+		set_ep(e);
 	}
 };
Index: /library/bdm/base/bdmbase.cpp
===================================================================
--- /library/bdm/base/bdmbase.cpp (revision 436)
+++ /library/bdm/base/bdmbase.cpp (revision 461)
@@ -135,6 +135,55 @@
 	for ( int i = 0;i < N;i++ ) X.set_col ( i, this->sample() );
 	return X;
-};
-
+}
+
+vec mpdf::samplecond(const vec &cond) {
+    condition(cond);
+    vec temp = shep->sample();
+    return temp;
+}
+
+mat mpdf::samplecond_m(const vec &cond, int N) {
+    condition(cond);
+    mat temp(shep->dimension(), N);
+    vec smp(shep->dimension());
+    for (int i = 0; i < N; i++) {
+        smp = shep->sample();
+	temp.set_col(i, smp);
+    }
+
+    return temp;
+}
+
+double mpdf::evallogcond(const vec &dt, const vec &cond) {
+    double tmp;
+    condition(cond);
+    tmp = shep->evallog(dt);
+    // it_assert_debug(std::isfinite(tmp), "Infinite value");
+    return tmp;
+}
+
+vec mpdf::evallogcond_m(const mat &Dt, const vec &cond) {
+    condition(cond);
+    return shep->evallog_m(Dt);
+}
+
+vec mpdf::evallogcond_m(const Array<vec> &Dt, const vec &cond) {
+    condition(cond);
+    return shep->evallog_m(Dt);
+}
+
+void mpdf::from_setting(const Setting &set){
+    if (set.exists("rv")) {
+        RV *r = UI::build<RV>(set, "rv");
+	set_rv(*r);
+	delete r;
+    }
+
+    if (set.exists("rvc")) {
+        RV *r = UI::build<RV>(set, "rvc");
+	set_rvc(*r); 
+	delete r;
+    }
+}
 
 std::ostream &operator<< ( std::ostream &os, const RV &rv ) {
Index: /library/bdm/stat/merger.cpp
===================================================================
--- /library/bdm/stat/merger.cpp (revision 423)
+++ /library/bdm/stat/merger.cpp (revision 461)
@@ -156,5 +156,5 @@
 				{
 					// no need for conditioning or marginalization
-					lw_src = mpdfs ( i )->_epdf().evallog_m ( Smp  );
+					lw_src = mpdfs ( i )->e()->evallog_m ( Smp  );
 				}
 				else
@@ -202,5 +202,5 @@
 					{
 						mpdfs ( i )->condition ( dls ( i )->get_cond ( Smp ( k ) ) );
-						lw_src ( k ) += mpdfs ( i )->_epdf().evallog ( dls ( i )->pushdown ( Smp ( k ) ) );
+						lw_src ( k ) += mpdfs ( i )->e()->evallog ( dls ( i )->pushdown ( Smp ( k ) ) );
 					}
 
Index: /library/bdm/stat/exp_family.h
===================================================================
--- /library/bdm/stat/exp_family.h (revision 450)
+++ /library/bdm/stat/exp_family.h (revision 461)
@@ -15,4 +15,5 @@
 
 
+#include "../shared_ptr.h"
 #include "../base/bdmbase.h"
 #include "../math/chmat.h"
@@ -530,5 +531,5 @@
 		protected:
 			//! Internal epdf that arise by conditioning on \c rvc
-			enorm<sq_T> epdf;
+			shared_ptr<enorm<sq_T> > iepdf;
 			mat A;
 			vec mu_const;
@@ -537,9 +538,10 @@
 			//! \name Constructors
 			//!@{
-			mlnorm ( ) :mEF (),epdf ( ),A ( ),_mu ( epdf._mu() ) {ep =&epdf; };
-			mlnorm ( const  mat &A, const vec &mu0, const sq_T &R ) :epdf ( ),_mu ( epdf._mu() )
-			{
-				ep =&epdf; set_parameters ( A,mu0,R );
-			};
+			mlnorm():iepdf(new enorm<sq_T>()), _mu(iepdf->_mu()) { set_ep(iepdf); };
+			mlnorm (const mat &A, const vec &mu0, const sq_T &R ) :iepdf(new enorm<sq_T>()), _mu(iepdf->_mu())
+			{
+				set_ep(iepdf); set_parameters ( A,mu0,R );
+			}
+
 			//! Set \c A and \c R
 			void set_parameters ( const  mat &A, const vec &mu0, const sq_T &R );
@@ -553,5 +555,5 @@
 			mat& _A() {return A;}
 			//!access function
-			mat _R() {return epdf._R().to_mat();}
+			mat _R() { return iepdf->_R().to_mat(); }
 
 			template<class sq_M>
@@ -577,12 +579,12 @@
 		protected:
 			//! Internal epdf that arise by conditioning on \c rvc
-			enorm<sq_T> epdf;
+			shared_ptr<enorm<sq_T> > iepdf;
 			vec &mu;
 			fnc* g;
 		public:
 			//!default constructor
-			mgnorm() :mu ( epdf._mu() ) {ep=&epdf;}
+			mgnorm():iepdf(new enorm<sq_T>()), mu(iepdf->_mu()) { set_ep(iepdf); }
 			//!set mean function
-			void set_parameters ( fnc* g0, const sq_T &R0 ) {g=g0; epdf.set_parameters ( zeros ( g->dimension() ), R0 );}
+			void set_parameters ( fnc* g0, const sq_T &R0 ) {g=g0; iepdf->set_parameters ( zeros ( g->dimension() ), R0 );}
 			void condition ( const vec &cond ) {mu=g->eval ( cond );};
 
@@ -662,5 +664,5 @@
 		public:
 			mlstudent ( ) :mlnorm<ldmat> (),
-					Lambda (),	_R ( epdf._R() ) {}
+					Lambda (),	_R ( iepdf->_R() ) {}
 			void set_parameters ( const mat &A0, const vec &mu0, const ldmat &R0, const ldmat& Lambda0 )
 			{
@@ -668,5 +670,5 @@
 				it_assert_debug ( R0.rows() ==A0.rows(),"" );
 
-				epdf.set_parameters ( mu0,Lambda ); //
+				iepdf->set_parameters ( mu0,Lambda ); //
 				A = A0;
 				mu_const = mu0;
@@ -705,15 +707,22 @@
 		protected:
 			//! Internal epdf that arise by conditioning on \c rvc
-			egamma epdf;
+			shared_ptr<egamma> iepdf;
+
 			//! Constant \f$k\f$
 			double k;
-			//! cache of epdf.beta
+
+			//! cache of iepdf.beta
 			vec &_beta;
 
 		public:
 			//! Constructor
-			mgamma ( ) : mEF ( ), epdf (), _beta ( epdf._beta() ) {ep=&epdf;};
+			mgamma():iepdf(new egamma()), k(0),
+			    _beta(iepdf->_beta()) {
+			    set_ep(iepdf);
+			}
+
 			//! Set value of \c k
-			void set_parameters ( double k, const vec &beta0 );
+			void set_parameters(double k, const vec &beta0);
+
 			void condition ( const vec &val ) {_beta=k/val;};
 			//! Load from structure with elements:
@@ -749,10 +758,13 @@
 		protected:
 			//! Internal epdf that arise by conditioning on \c rvc
-			eigamma epdf;
+			shared_ptr<eigamma> iepdf;
+
 			//! Constant \f$k\f$
 			double k;
-			//! cache of epdf.alpha
+
+			//! cache of iepdf.alpha
 			vec &_alpha;
-			//! cache of epdf.beta
+
+			//! cache of iepdf.beta
 			vec &_beta;
 
@@ -760,6 +772,17 @@
 			//! \name Constructors
 			//!@{
-			migamma ( ) : mEF (), epdf ( ), _alpha ( epdf._alpha() ), _beta ( epdf._beta() ) {ep=&epdf;};
-			migamma ( const migamma &m ) : mEF (), epdf ( m.epdf ), _alpha ( epdf._alpha() ), _beta ( epdf._beta() ) {ep=&epdf;};
+			migamma():iepdf(new eigamma()),
+			    k(0),
+			    _alpha(iepdf->_alpha()),
+			    _beta(iepdf->_beta()) {
+			    set_ep(iepdf);
+			}
+
+			migamma(const migamma &m):iepdf(m.iepdf),
+			    k(0),
+			    _alpha(iepdf->_alpha()),
+			    _beta(iepdf->_beta()) {
+			    set_ep(iepdf);
+			}
 			//!@}
 
@@ -768,5 +791,5 @@
 			{
 				k=k0;
-				epdf.set_parameters ( ( 1.0/ ( k*k ) +2.0 ) *ones ( len ) /*alpha*/, ones ( len ) /*beta*/ );
+				iepdf->set_parameters ( ( 1.0/ ( k*k ) +2.0 ) *ones ( len ) /*alpha*/, ones ( len ) /*beta*/ );
 				dimc = dimension();
 			};
@@ -907,17 +930,24 @@
 	{
 		protected:
-			elognorm eno;
+			shared_ptr<elognorm> eno;
+
 			//! parameter 1/2*sigma^2
 			double sig2;
+
 			//! access
 			vec &mu;
 		public:
 			//! Constructor
-			mlognorm ( ) : eno (), mu ( eno._mu() ) {ep=&eno;};
+			mlognorm():eno(new elognorm()),
+			    sig2(0),
+			    mu(eno->_mu()) {
+			    set_ep(eno);
+			}
+
 			//! Set value of \c k
 			void set_parameters ( int size, double k )
 			{
 				sig2 = 0.5*log ( k*k+1 );
-				eno.set_parameters ( zeros ( size ),2*sig2*eye ( size ) );
+				eno->set_parameters ( zeros ( size ),2*sig2*eye ( size ) );
 
 				dimc = size;
@@ -1053,5 +1083,5 @@
 	{
 		protected:
-			eiWishartCh eiW;
+			shared_ptr<eiWishartCh> eiW;
 			//!square root of \f$ \nu-p-1 \f$ - needed for computation of \f$ \Psi \f$ from conditions
 			double sqd;
@@ -1060,5 +1090,11 @@
 			double l;
 			int p;
-		public:
+
+		public:
+			rwiWishartCh():eiW(new eiWishartCh()),
+			    sqd(0), l(0), p(0) {
+			    set_ep(eiW);
+			}
+
 			void set_parameters ( int p0, double k, vec ref0, double l0  )
 			{
@@ -1069,7 +1105,6 @@
 				refl=pow(ref0,1-l);
 				
-				eiW.set_parameters ( eye ( p ),delta );
-				ep=&eiW;
-				dimc=eiW.dimension();
+				eiW->set_parameters ( eye ( p ),delta );
+				dimc=eiW->dimension();
 			}
 			void condition ( const vec &c ) {
@@ -1081,5 +1116,5 @@
 				}
 
-				eiW._setY ( sqd*z );
+				eiW->_setY ( sqd*z );
 			}
 	};
@@ -1239,8 +1274,8 @@
 		it_assert_debug ( A0.rows() ==R0.rows(),"" );
 
-		epdf.set_parameters ( zeros ( A0.rows() ),R0 );
+		iepdf->set_parameters(zeros(A0.rows()), R0);
 		A = A0;
 		mu_const = mu0;
-		dimc=A0.cols();
+		dimc = A0.cols();
 	}
 
@@ -1332,5 +1367,5 @@
 		os << "A:"<< ml.A<<endl;
 		os << "mu:"<< ml.mu_const<<endl;
-		os << "R:" << ml.epdf._R().to_mat() <<endl;
+		os << "R:" << ml.iepdf->_R().to_mat() <<endl;
 		return os;
 	};
Index: /library/bdm/stat/emix.h
===================================================================
--- /library/bdm/stat/emix.h (revision 395)
+++ /library/bdm/stat/emix.h (revision 461)
@@ -16,4 +16,5 @@
 #define LOG2  0.69314718055995	
 
+#include "../shared_ptr.h"
 #include "exp_family.h"
 
@@ -51,7 +52,7 @@
 		rvc = nom0->_rv().subt ( rv );
 		dimc = rvc._dsize();
-		ep = new epdf;
-		ep->set_parameters(rv._dsize());
-		ep->set_rv(rv);
+		set_ep(shared_ptr<epdf>(new epdf));
+		e()->set_parameters(rv._dsize());
+		e()->set_rv(rv);
 		
 		//prepare data structures
@@ -66,5 +67,5 @@
 	double evallogcond ( const vec &val, const vec &cond ) {
 		double tmp;
-		vec nom_val ( ep->dimension() + dimc );
+		vec nom_val ( e()->dimension() + dimc );
 		dl.pushup_cond ( nom_val,val,cond );
 		tmp = exp ( nom->evallog ( nom_val ) - den->evallog ( cond ) );
@@ -228,11 +229,17 @@
 	//! Data link for each mpdfs
 	Array<datalink_m2m*> dls;
+
 	//! dummy ep
-	epdf dummy;
+	shared_ptr<epdf> dummy;
+
 public:
 	/*!\brief Constructor from list of mFacs,
 	*/
-	mprod (){};
-	mprod (Array<mpdf*> mFacs ){set_elements( mFacs );};
+	mprod():dummy(new epdf()) { }
+	mprod (Array<mpdf*> mFacs):
+	    dummy(new epdf()) {
+	    set_elements(mFacs);
+	}
+
 	void set_elements(Array<mpdf*> mFacs , bool own=false) {
 		
@@ -241,9 +248,9 @@
 		epdfs.set_size(mFacs.length());
 				
-		ep=&dummy;
+		set_ep(dummy);
 		RV rv=getrv ( true );
 		set_rv ( rv );
-		dummy.set_parameters ( rv._dsize() );
-		setrvc ( ep->_rv(),rvc );
+		dummy->set_parameters(rv._dsize());
+		setrvc(e()->_rv(), rvc);
 		// rv and rvc established = > we can link them with mpdfs
 		for ( int i = 0;i < mpdfs.length(); i++ ) {
@@ -253,5 +260,5 @@
 
 		for ( int i=0; i<mpdfs.length(); i++ ) {
-			epdfs ( i ) =& ( mpdfs ( i )->_epdf() );
+			epdfs(i) = mpdfs(i)->e();
 		}
 	};
@@ -291,6 +298,6 @@
 	//TODO smarter...
 	vec samplecond ( const vec &cond ) {
-		//! Ugly hack to help to discover if mpfs are not in proper order. Correct solution = check that explicitely.
-		vec smp= std::numeric_limits<double>::infinity() * ones ( ep->dimension() );
+	    //! Ugly hack to help to discover if mpfs are not in proper order. Correct solution = check that explicitely.
+	    vec smp= std::numeric_limits<double>::infinity() * ones ( e()->dimension() );
 		vec smpi;
 		// Hard assumption here!!! We are going backwards, to assure that samples that are needed from smp are already generated!
@@ -417,9 +424,14 @@
 	//! Component (epdfs)
 	Array<mpdf*> Coms;
+
 	//!Internal epdf
-	emix Epdf;
+	shared_ptr<emix> iepdf;
+
 public:
 	//!Default constructor
-	mmix ( ) : mpdf ( ), Epdf () {ep = &Epdf;};
+	mmix():iepdf(new emix()) {
+	    set_ep(iepdf);
+	}
+
 	//! Set weights \c w and components \c R
 	void set_parameters ( const vec &w, const Array<mpdf*> &Coms ) {
@@ -427,8 +439,9 @@
 
 		for ( int i = 0;i < Coms.length();i++ ) {
-			Eps ( i ) = & ( Coms ( i )->_epdf() );
-		}
-		Epdf.set_parameters ( w, Eps );
-	};
+			Eps(i) = Coms(i)->e();
+		}
+
+		iepdf->set_parameters(w, Eps);
+	}
 
 	void condition ( const vec &cond ) {
Index: /library/bdm/stat/exp_family.cpp
===================================================================
--- /library/bdm/stat/exp_family.cpp (revision 404)
+++ /library/bdm/stat/exp_family.cpp (revision 461)
@@ -211,11 +211,9 @@
 }
 
-//MGamma
-void mgamma::set_parameters ( double k0, const vec &beta0 ) {
-	k=k0;
-	ep = &epdf;
-	epdf.set_parameters ( k*ones ( beta0.length()),beta0 );
-	dimc = ep->dimension();
-};
+void mgamma::set_parameters(double k0, const vec &beta0) {
+	k = k0;
+	iepdf->set_parameters(k * ones(beta0.length()), beta0);
+	dimc = e()->dimension();
+}
 
 ivec eEmp::resample (RESAMPLING_METHOD method) {
Index: /library/tests/emix_test.cpp
===================================================================
--- /library/tests/emix_test.cpp (revision 454)
+++ /library/tests/emix_test.cpp (revision 461)
@@ -1,4 +1,6 @@
+#include "shared_ptr.h"
 #include "stat/exp_family.h"
 #include "stat/emix.h"
+
 using namespace bdm;
 
@@ -80,5 +82,5 @@
 	cout << "original empirical variance: " <<Er <<endl;
 
-	epdf* Mg = (epdf*)M2.marginal ( y );
+	shared_ptr<epdf> Mg(dynamic_cast<epdf *>(M2.marginal(y)));
 	mpdf* Cn = (mpdf*)M2.condition ( x );
 
Index: /library/tests/testSmp.cpp
===================================================================
--- /library/tests/testSmp.cpp (revision 450)
+++ /library/tests/testSmp.cpp (revision 461)
@@ -94,5 +94,5 @@
 	
 	Smp = mMix.samplecond_m(mu0,N);
-	disp(mMix._epdf().mean(),zeros(2),Smp);
+	disp(mMix.e()->mean(), zeros(2), Smp);
 
 	cout << "======= EProd ======== " << endl;
Index: /library/tests/epdf_harness.cpp
===================================================================
--- /library/tests/epdf_harness.cpp (revision 458)
+++ /library/tests/epdf_harness.cpp (revision 461)
@@ -80,8 +80,8 @@
 	mat smp = hepdf->sample_m(nsamples);
 	int n = smp.cols();
-	vec Emu = smp * ones(n) / n;
-	mat Er = (smp*smp.transpose())/n - outer_product(Emu,Emu);
-	CHECK_CLOSE_EX(mean, Emu, tolerance);
-	CHECK_CLOSE_EX(R, Er, tolerance);
+	vec emu = smp * ones(n) / n;
+	mat er = (smp * smp.T()) / n - outer_product(emu, emu);
+	CHECK_CLOSE_EX(mean, emu, tolerance);
+	CHECK_CLOSE_EX(R, er, tolerance);
     }
 }
Index: /library/tests/test_kalman_QR.cpp
===================================================================
--- /library/tests/test_kalman_QR.cpp (revision 386)
+++ /library/tests/test_kalman_QR.cpp (revision 461)
@@ -64,5 +64,5 @@
 	KF_QR.set_parameters(&evolQR,&evolQR,100);
 	evolQR.condition("1 1 1");
-	KF_QR.set_statistics(evolQR._e(), &KF);
+	KF_QR.set_statistics(evolQR.e(), &KF);
 	const epdf& mpost=KF_QR.posterior();
 	const epdf& mposttr=KFtr.posterior();
