Index: library/bdm/stat/exp_family.h
===================================================================
--- library/bdm/stat/exp_family.h (revision 476)
+++ library/bdm/stat/exp_family.h (revision 487)
@@ -55,5 +55,5 @@
 				return tmp;}
 			//!Evaluate normalized log-probability for many samples
-			virtual vec evallog ( const mat &Val ) const
+			virtual vec evallog_m ( const mat &Val ) const
 			{
 				vec x ( Val.cols() );
@@ -65,17 +65,4 @@
 	};
 
-	/*!
-	* \brief Exponential family model.
-
-	* More?...
-	*/
-
-	class mEF : public mpdf
-	{
-
-		public:
-			//! Default constructor
-			mEF ( ) :mpdf ( ) {};
-	};
 
 //! Estimator for Exponential family
@@ -106,5 +93,5 @@
 	};
 
-	template<class sq_T>
+	template<class sq_T, template <typename> class TEpdf >
 	class mlnorm;
 
@@ -527,27 +514,37 @@
 	 Mean value \f$mu=A*rvc+mu_0\f$.
 	*/
-	template<class sq_T>
-	class mlnorm : public mEF
+	template<class sq_T, template <typename> class TEpdf =enorm>
+	class mlnorm : public mpdf_internal< TEpdf<sq_T> >
 	{
 		protected:
 			//! Internal epdf that arise by conditioning on \c rvc
-			shared_ptr<enorm<sq_T> > iepdf;
 			mat A;
 			vec mu_const;
-			vec& _mu; //cached epdf.mu;
+//			vec& _mu; //cached epdf.mu; !!!!!! WHY NOT?
 		public:
 			//! \name Constructors
 			//!@{
-			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 );
+			mlnorm():mpdf_internal< TEpdf<sq_T> >() {};
+			mlnorm (const mat &A, const vec &mu0, const sq_T &R ) : mpdf_internal< TEpdf<sq_T> >()
+			{
+				set_parameters ( A,mu0,R );
 			}
 
 			//! Set \c A and \c R
-			void set_parameters ( const  mat &A, const vec &mu0, const sq_T &R );
+			void set_parameters ( const  mat &A0, const vec &mu0, const sq_T &R0 ){
+				it_assert_debug ( A0.rows() ==mu0.length(),"" );
+				it_assert_debug ( A0.rows() ==R0.rows(),"" );
+
+				this->iepdf.set_parameters(zeros(A0.rows()), R0);
+				A = A0;
+				mu_const = mu0;
+				this->dimc = A0.cols();
+			}
 			//!@}
 			//! Set value of \c rvc . Result of this operation is stored in \c epdf use function \c _ep to access it.
-			void condition ( const vec &cond );
+			void condition ( const vec &cond ){
+				this->iepdf._mu()= A*cond + mu_const;
+//R is already assigned;
+			}
 
 			//!access function
@@ -556,8 +553,8 @@
 			mat& _A() {return A;}
 			//!access function
-			mat _R() { return iepdf->_R().to_mat(); }
-
-			template<class sq_M>
-			friend std::ostream &operator<< ( std::ostream &os,  mlnorm<sq_M> &ml );
+			mat _R(){ return this->iepdf._R().to_mat(); }
+		
+			template<typename sq_M>
+			friend std::ostream &operator<< ( std::ostream &os,  mlnorm<sq_M,enorm> &ml );
 			
 			void from_setting(const Setting &set){
@@ -577,17 +574,15 @@
 //! Mpdf with general function for mean value
 	template<class sq_T>
-	class mgnorm : public mEF
-	{
-		protected:
-			//! Internal epdf that arise by conditioning on \c rvc
-			shared_ptr<enorm<sq_T> > iepdf;
-			vec &mu;
+	class mgnorm : public mpdf_internal< enorm< sq_T > >
+	{
+		protected:
+//			vec &mu; WHY NOT?
 			fnc* g;
 		public:
 			//!default constructor
-			mgnorm():iepdf(new enorm<sq_T>()), mu(iepdf->_mu()) { set_ep(iepdf); }
+			mgnorm():mpdf_internal<enorm<sq_T> >() { }
 			//!set mean function
-			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 );};
+			inline void set_parameters ( fnc* g0, const sq_T &R0 );
+			inline void condition ( const vec &cond );
 
 
@@ -644,5 +639,5 @@
 	Perhaps a moment-matching technique?
 	*/
-	class mlstudent : public mlnorm<ldmat>
+	class mlstudent : public mlnorm<ldmat,enorm>
 	{
 		protected:
@@ -651,6 +646,6 @@
 			ldmat Re;
 		public:
-			mlstudent ( ) :mlnorm<ldmat> (),
-					Lambda (),	_R ( iepdf->_R() ) {}
+			mlstudent ( ) :mlnorm<ldmat,enorm> (),
+					Lambda (),	_R ( iepdf._R() ) {}
 			void set_parameters ( const mat &A0, const vec &mu0, const ldmat &R0, const ldmat& Lambda0 )
 			{
@@ -658,5 +653,5 @@
 				it_assert_debug ( R0.rows() ==A0.rows(),"" );
 
-				iepdf->set_parameters ( mu0,Lambda ); //
+				iepdf.set_parameters ( mu0,Lambda ); //
 				A = A0;
 				mu_const = mu0;
@@ -666,5 +661,5 @@
 			void condition ( const vec &cond )
 			{
-				_mu = A*cond + mu_const;
+				iepdf._mu() = A*cond + mu_const;
 				double zeta;
 				//ugly hack!
@@ -691,9 +686,7 @@
 	The standard deviation of the walk is then: \f$\mu/\sqrt(k)\f$.
 	*/
-	class mgamma : public mEF
-	{
-		protected:
-			//! Internal epdf that arise by conditioning on \c rvc
-			shared_ptr<egamma> iepdf;
+	class mgamma : public mpdf_internal<egamma>
+	{
+		protected:
 
 			//! Constant \f$k\f$
@@ -705,7 +698,6 @@
 		public:
 			//! Constructor
-			mgamma():iepdf(new egamma()), k(0),
-			    _beta(iepdf->_beta()) {
-			    set_ep(iepdf);
+			mgamma():mpdf_internal<egamma>(), k(0),
+			    _beta(iepdf._beta()) {
 			}
 
@@ -742,10 +734,7 @@
 	The standard deviation of the walk is then: \f$ \mu/\sqrt(k)\f$.
 	 */
-	class migamma : public mEF
-	{
-		protected:
-			//! Internal epdf that arise by conditioning on \c rvc
-			shared_ptr<eigamma> iepdf;
-
+	class migamma : public mpdf_internal<eigamma>
+	{
+		protected:
 			//! Constant \f$k\f$
 			double k;
@@ -760,16 +749,14 @@
 			//! \name Constructors
 			//!@{
-			migamma():iepdf(new eigamma()),
+			migamma():mpdf_internal<eigamma>(),
 			    k(0),
-			    _alpha(iepdf->_alpha()),
-			    _beta(iepdf->_beta()) {
-			    set_ep(iepdf);
-			}
-
-			migamma(const migamma &m):iepdf(m.iepdf),
+			    _alpha(iepdf._alpha()),
+			    _beta(iepdf._beta()) {
+			}
+
+			migamma(const migamma &m):mpdf_internal<eigamma>(),
 			    k(0),
-			    _alpha(iepdf->_alpha()),
-			    _beta(iepdf->_beta()) {
-			    set_ep(iepdf);
+			    _alpha(iepdf._alpha()),
+			    _beta(iepdf._beta()) {
 			}
 			//!@}
@@ -779,5 +766,5 @@
 			{
 				k=k0;
-				iepdf->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();
 			};
@@ -915,9 +902,7 @@
 	The standard deviation of the walk is then: \f$\mu/\sqrt(k)\f$.
 	 */
-	class mlognorm : public mpdf
-	{
-		protected:
-			shared_ptr<elognorm> eno;
-
+	class mlognorm : public mpdf_internal<elognorm>
+	{
+		protected:
 			//! parameter 1/2*sigma^2
 			double sig2;
@@ -927,8 +912,7 @@
 		public:
 			//! Constructor
-			mlognorm():eno(new elognorm()),
+			mlognorm():mpdf_internal<elognorm>(),
 			    sig2(0),
-			    mu(eno->_mu()) {
-			    set_ep(eno);
+			    mu(iepdf._mu()) {
 			}
 
@@ -937,5 +921,5 @@
 			{
 				sig2 = 0.5*log ( k*k+1 );
-				eno->set_parameters ( zeros ( size ),2*sig2*eye ( size ) );
+				iepdf.set_parameters ( zeros ( size ),2*sig2*eye ( size ) );
 
 				dimc = size;
@@ -1071,5 +1055,5 @@
 	{
 		protected:
-			shared_ptr<eiWishartCh> eiW;
+			eiWishartCh eiW;
 			//!square root of \f$ \nu-p-1 \f$ - needed for computation of \f$ \Psi \f$ from conditions
 			double sqd;
@@ -1080,5 +1064,5 @@
 
 		public:
-			rwiWishartCh():eiW(new eiWishartCh()),
+			rwiWishartCh():eiW(),
 			    sqd(0), l(0), p(0) {
 			    set_ep(eiW);
@@ -1093,6 +1077,6 @@
 				refl=pow(ref0,1-l);
 				
-				eiW->set_parameters ( eye ( p ),delta );
-				dimc=eiW->dimension();
+				eiW.set_parameters ( eye ( p ),delta );
+				dimc=eiW.dimension();
 			}
 			void condition ( const vec &c ) {
@@ -1104,5 +1088,5 @@
 				}
 
-				eiW->_setY ( sqd*z );
+				eiW._setY ( sqd*z );
 			}
 	};
@@ -1256,15 +1240,4 @@
 	};
 
-	template<class sq_T>
-	void mlnorm<sq_T>::set_parameters ( const mat &A0, const vec &mu0, const sq_T &R0 )
-	{
-		it_assert_debug ( A0.rows() ==mu0.length(),"" );
-		it_assert_debug ( A0.rows() ==R0.rows(),"" );
-
-		iepdf->set_parameters(zeros(A0.rows()), R0);
-		A = A0;
-		mu_const = mu0;
-		dimc = A0.cols();
-	}
 
 // template<class sq_T>
@@ -1293,10 +1266,4 @@
 // }
 
-	template<class sq_T>
-	void mlnorm<sq_T>::condition ( const vec &cond )
-	{
-		_mu = A*cond + mu_const;
-//R is already assigned;
-	}
 
 	template<class sq_T>
@@ -1348,5 +1315,10 @@
 	}
 
-///////////
+	////
+	///////
+	template<class sq_T>				
+			void mgnorm<sq_T >::set_parameters ( fnc* g0, const sq_T &R0 ) {g=g0; this->iepdf.set_parameters ( zeros ( g->dimension() ), R0 );}
+	template<class sq_T>				
+			void mgnorm<sq_T >::condition ( const vec &cond ) {this->iepdf._mu()=g->eval ( cond );};
 
 	template<class sq_T>
