Changeset 487

Timestamp:

08/08/09 13:42:18 (15 years ago)

Author:

smidl

Message:

1st step of mpdf redesign - BROKEN compile

Location:

library/bdm

Files:

• base/bdmbase.cpp (modified) (3 diffs)
• base/bdmbase.h (modified) (11 diffs)
• estim/particles.cpp (modified) (1 diff)
• estim/particles.h (modified) (1 diff)
• stat/emix.h (modified) (7 diffs)
• stat/exp_family.cpp (modified) (2 diffs)
• stat/exp_family.h (modified) (25 diffs)
• stat/merger.cpp (modified) (2 diffs)

library/bdm/base/bdmbase.cpp

@@ -140 +140 @@
 }
 
-vec mpdf::samplecond ( const vec &cond ) {
-        condition ( cond );
-        vec temp = shep->sample();
+template<class EPDF>
+vec mpdf_internal<EPDF>::samplecond ( const vec &cond ) {
+        condition ( cond );
+        vec temp = iepdf.sample();
         return temp;
 }
 
-mat mpdf::samplecond_m ( const vec &cond, int N ) {
-        condition ( cond );
-        mat temp ( shep->dimension(), N );
-        vec smp ( shep->dimension() );
+template<class EPDF>
+mat mpdf_internal<EPDF>::samplecond_m ( const vec &cond, int N ) {
+        condition ( cond );
+        mat temp ( dimension(), N );
+        vec smp ( dimension() );
         for ( int i = 0; i < N; i++ ) {
-                smp = shep->sample();
+                smp = iepdf.sample();
                 temp.set_col ( i, smp );
         }
@@ -158 +160 @@
 }
 
-double mpdf::evallogcond ( const vec &dt, const vec &cond ) {
+template<class EPDF>
+double mpdf_internal<EPDF>::evallogcond ( const vec &dt, const vec &cond ) {
         double tmp;
         condition ( cond );
-        tmp = shep->evallog ( dt );
+        tmp = iepdf.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 );
+template<class EPDF>
+vec mpdf_internal<EPDF>::evallogcond_m ( const mat &Dt, const vec &cond ) {
+        condition ( cond );
+        return iepdf.evallog_m ( Dt );
+}
+
+template<class EPDF>
+                vec mpdf_internal<EPDF>::evallogcond_m ( const Array<vec> &Dt, const vec &cond ) {
+        condition ( cond );
+        return iepdf.evallog_m ( Dt );
 }
 
@@ -313 +318 @@
 void mepdf::from_setting ( const Setting &set ) {
         shared_ptr<epdf> e ( UI::build<epdf> ( set, "epdf", UI::compulsory ) );
-        set_ep ( e );
+        ipdf = e;
 }
 

library/bdm/base/bdmbase.h

@@ -399 +399 @@
 //! Conditional probability density, e.g. modeling some dependencies.
 //TODO Samplecond can be generalized
-
 class mpdf : public root {
 protected:
@@ -406 +405 @@
         //! random variable in condition
         RV rvc;
-
 private:
-        //! pointer to internal epdf
-        shared_ptr<epdf> shep;
+        //! internal epdf, used only as cache to avoid virtual calls of \c _rv() and \c _dimension()
+        epdf* ep;
+        
+protected:
+        void set_ep(epdf &iepdf) {
+                ep = &iepdf;
+        }
 
 public:
@@ -415 +418 @@
         //! @{
 
-        mpdf() : dimc ( 0 ), rvc() { }
-
-        mpdf ( const mpdf &m ) : dimc ( m.dimc ), rvc ( m.rvc ), shep ( m.shep ) { }
+        mpdf() : dimc ( 0 ), rvc(), ep(NULL) { }
+
+        mpdf ( const mpdf &m ) : dimc ( m.dimc ), rvc ( m.rvc ), ep ( m.ep ) { }
         //!@}
 
@@ -424 +427 @@
 
         //! 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 );
+        virtual vec samplecond ( const vec &cond ){it_error("Not implemented");return vec(0);};
 
         //! 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" );
-        };
+        virtual mat samplecond_m ( const vec &cond, int N ){it_error("Not implemented");return mat();}
+
 
         //! 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 );
+        virtual double evallogcond ( const vec &dt, const vec &cond ){it_error("Not implemented");return 0.0;}
 
         //! Matrix version of evallogcond
-        virtual vec evallogcond_m ( const mat &Dt, const vec &cond );
+        virtual vec evallogcond_m ( const mat &Dt, const vec &cond ){it_error("Not implemented");return vec(0);}
 
         //! Array<vec> version of evallogcond
-        virtual vec evallogcond_m ( const Array<vec> &Dt, const vec &cond );
+        virtual vec evallogcond_m ( const Array<vec> &Dt, const vec &cond ){it_error("Not implemented");return vec(0);}
 
         //! \name Access to attributes
@@ -447 +446 @@
 
         RV _rv() {
-                return shep->_rv();
+                return ep->_rv();
         }
         RV _rvc() {
@@ -454 +453 @@
         }
         int dimension() {
-                return shep->dimension();
+                return ep->dimension();
         }
         int dimensionc() {
@@ -460 +459 @@
         }
 
-        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
+        //! { class = "mpdf_offspring",
+        //!   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
@@ -485 +477 @@
         }
         void set_rv ( const RV &rv0 ) {
-                shep->set_rv ( rv0 );
+                ep->set_rv ( rv0 );
         }
         bool isnamed() {
-                return ( shep->isnamed() ) && ( dimc == rvc._dsize() );
-        }
-        //!@}
+                return ( ep->isnamed() ) && ( dimc == rvc._dsize() );
+        }
+        //!@}
+};
+
+template <class EPDF>
+class mpdf_internal: public mpdf{
+        protected :
+                EPDF iepdf;
+        public:
+                //! constructor
+                mpdf_internal(): mpdf(),iepdf(){set_ep(iepdf);}
+                //! Update \c iepdf so that it represents this mpdf conditioned on \c rvc = cond
+                //! This function provides convenient reimplementation in offsprings
+                        virtual void condition ( const vec &cond ) {
+                                it_error ( "Not implemented" );
+                        };
+                //!access function to iepdf
+                EPDF& e(){return iepdf;}
+                        
+                //! Reimplements samplecond using \c condition()
+                vec samplecond ( const vec &cond );
+                //! Reimplements evallogcond using \c condition()
+                double evallogcond ( const vec &val, const vec &cond );
+                //! Efficient version of evallogcond for matrices 
+                virtual vec evallogcond_m ( const mat &Dt, const vec &cond );
+                //! Efficient version of evallogcond for Array<vec>
+                virtual vec evallogcond_m ( const Array<vec> &Dt, const vec &cond );
+                //! Efficient version of samplecond 
+                virtual mat samplecond_m ( const vec &cond, int N );
 };
 
@@ -686 +705 @@
 */
 class mepdf : public mpdf {
+
+shared_ptr<epdf> ipdf;
 public:
         //!Default constructor
@@ -691 +712 @@
 
         mepdf ( shared_ptr<epdf> em ) {
-                set_ep ( em );
+                ipdf = em;
+                set_ep(*ipdf.get());
                 dimc = 0;
         }
@@ -701 +723 @@
         //!  \code
         //! { class = "mepdf",
-        //!   epdfs = {class="epdfs",...}
+        //!   epdf = {class="epdf_offspring",...}
         //! }
         //! \endcode

library/bdm/estim/particles.cpp

@@ -13 +13 @@
         for ( i = 0; i < n; i++ ) {
                 //generate new samples from paramater evolution model;
-                _samples ( i ) = par->samplecond ( _samples ( i ) );
-                lls ( i ) = par->e()->evallog ( _samples ( i ) );
+                vec old_smp=_samples ( i );
+                _samples ( i ) = par->samplecond ( old_smp );
+                lls ( i ) = par->evallogcond ( _samples ( i ), old_smp );
                 lls ( i ) *= obs->evallogcond ( dt, _samples ( i ) );
 

library/bdm/estim/particles.h

@@ -250 +250 @@
         for ( i = 0; i < n; i++ ) {
                 //generate new samples from paramater evolution model;
-                _samples ( i ) = par->samplecond ( _samples ( i ) );
-                llsP ( i ) = par->e()->evallog ( _samples ( i ) );
+                vec old_smp=_samples ( i );
+                _samples ( i ) = par->samplecond ( old_smp );
+                llsP ( i ) = par->evallogcond ( _samples ( i ), old_smp );
                 BMs ( i )->condition ( _samples ( i ) );
                 BMs ( i )->bayes ( dt );

library/bdm/stat/emix.h

@@ -45 +45 @@
         //!datalink between conditional and nom
         datalink_m2e dl;
+        //! dummy epdf that stores only rv and dim
+        epdf iepdf;
 public:
         //!Default constructor. By default, the given epdf is not copied!
         //! It is assumed that this function will be used only temporarily.
-        mratio ( const epdf* nom0, const RV &rv, bool copy = false ) : mpdf ( ), dl ( ) {
+        mratio ( const epdf* nom0, const RV &rv, bool copy = false ) : mpdf ( ), dl ( ),iepdf() {
                 // adjust rv and rvc
                 rvc = nom0->_rv().subt ( rv );
                 dimc = rvc._dsize();
-                set_ep ( shared_ptr<epdf> ( new epdf ) );
-                e()->set_parameters ( rv._dsize() );
-                e()->set_rv ( rv );
+                set_ep ( iepdf );
+                iepdf.set_parameters ( rv._dsize() );
+                iepdf.set_rv ( rv );
 
                 //prepare data structures
@@ -72 +74 @@
         double evallogcond ( const vec &val, const vec &cond ) {
                 double tmp;
-                vec nom_val ( e()->dimension() + dimc );
+                vec nom_val ( dimension() + dimc );
                 dl.pushup_cond ( nom_val, val, cond );
                 tmp = exp ( nom->evallog ( nom_val ) - den->evallog ( cond ) );
-                it_assert_debug ( std::isfinite ( tmp ), "Infinite value" );
                 return tmp;
         }
@@ -278 +279 @@
 class mprod: public compositepdf, public mpdf {
 protected:
-        //! pointers to epdfs - shortcut to mpdfs().posterior()
-        Array<epdf*> epdfs;
         //! Data link for each mpdfs
         Array<datalink_m2m*> dls;
 
-        //! dummy ep
-        shared_ptr<epdf> dummy;
+        //! dummy epdf used only as storage for RV and dim
+        epdf iepdf;
 
 public:
         /*!\brief Constructor from list of mFacs,
         */
-        mprod() : dummy ( new epdf() ) { }
+        mprod() : iepdf( ) { }
         mprod ( Array<mpdf*> mFacs ) :
-                        dummy ( new epdf() ) {
+                        iepdf ( ) {
                 set_elements ( mFacs );
         }
@@ -299 +298 @@
                 compositepdf::set_elements ( mFacs, own );
                 dls.set_size ( mFacs.length() );
-                epdfs.set_size ( mFacs.length() );
-
-                set_ep ( dummy );
+
+                set_ep ( iepdf);
                 RV rv = getrv ( true );
                 set_rv ( rv );
-                dummy->set_parameters ( rv._dsize() );
-                setrvc ( e()->_rv(), rvc );
+                iepdf.set_parameters ( rv._dsize() );
+                setrvc (_rv(), rvc );
                 // rv and rvc established = > we can link them with mpdfs
                 for ( int i = 0; i < mpdfs.length(); i++ ) {
@@ -312 +310 @@
                 }
 
-                for ( int i = 0; i < mpdfs.length(); i++ ) {
-                        epdfs ( i ) = mpdfs ( i )->e();
-                }
         };
 
@@ -352 +347 @@
         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 ( e()->dimension() );
+                vec smp = std::numeric_limits<double>::infinity() * ones ( dimension() );
                 vec smpi;
                 // Hard assumption here!!! We are going backwards, to assure that samples that are needed from smp are already generated!
                 for ( int i = ( mpdfs.length() - 1 ); i >= 0; i-- ) {
-                        if ( mpdfs ( i )->dimensionc() ) {
-                                mpdfs ( i )->condition ( dls ( i )->get_cond ( smp , cond ) ); // smp is val here!!
-                        }
-                        smpi = epdfs ( i )->sample();
+                        // generate contribution of this mpdf
+                        smpi = mpdfs(i)->samplecond(dls ( i )->get_cond ( smp , cond ));                        
                         // copy contribution of this pdf into smp
                         dls ( i )->pushup ( smp, smpi );
-                        // add ith likelihood contribution
                 }
                 return smp;
@@ -483 +475 @@
 
 */
-class mmix : public mpdf {
-protected:
-        //! Component (epdfs)
-        Array<mpdf*> Coms;
-
-        //!Internal epdf
-        shared_ptr<emix> iepdf;
-
-public:
-        //!Default constructor
-        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 ) {
-                Array<epdf*> Eps ( Coms.length() );
-
-                for ( int i = 0; i < Coms.length(); i++ ) {
-                        Eps ( i ) = Coms ( i )->e();
-                }
-
-                iepdf->set_parameters ( w, Eps );
-        }
-
-        void condition ( const vec &cond ) {
-                for ( int i = 0; i < Coms.length(); i++ ) {
-                        Coms ( i )->condition ( cond );
-                }
-        };
-};
+// class mmix : public mpdf {
+// protected:
+//      //! Component (epdfs)
+//      Array<mpdf*> Coms;
+// 
+//      //!Internal epdf
+//      emix iepdf;
+// 
+// public:
+//      //!Default constructor
+//      mmix() : iepdf ( ) {
+//              set_ep ( iepdf );
+//      }
+// 
+//      //! Set weights \c w and components \c R
+//      void set_parameters ( const vec &w, const Array<mpdf*> &Coms ) {
+//              Array<epdf*> Eps ( Coms.length() );
+// 
+//              for ( int i = 0; i < Coms.length(); i++ ) {
+//                      Eps ( i ) = Coms ( i )->e();
+//              }
+// 
+//              iepdf->set_parameters ( w, Eps );
+//      }
+// 
+//      void condition ( const vec &cond ) {
+//              for ( int i = 0; i < Coms.length(); i++ ) {
+//                      Coms ( i )->condition ( cond );
+//              }
+//      };
+// };
 
 }

library/bdm/stat/exp_family.cpp

@@ -11 +11 @@
 
 using std::cout;
+
+        ///////////
 
 void BMEF::bayes ( const vec &dt ) {
@@ -214 +216 @@
 void mgamma::set_parameters ( double k0, const vec &beta0 ) {
         k = k0;
-        iepdf->set_parameters ( k * ones ( beta0.length() ), beta0 );
-        dimc = e()->dimension();
+        iepdf.set_parameters ( k * ones ( beta0.length() ), beta0 );
+        dimc = iepdf.dimension();
 }
 

library/bdm/stat/exp_family.h

@@ -55 +55 @@
                                 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 +65 @@
         };
 
-        /*!
-        * \brief Exponential family model.
-
-        * More?...
-        */
-
-        class mEF : public mpdf
-        {
-
-                public:
-                        //! Default constructor
-                        mEF ( ) :mpdf ( ) {};
-        };
 
 //! Estimator for Exponential family
@@ -106 +93 @@
         };
 
-        template<class sq_T>
+        template<class sq_T, template <typename> class TEpdf >
         class mlnorm;
 
@@ -527 +514 @@
          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 +553 @@
                         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 +574 @@
 //! 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 +639 @@
         Perhaps a moment-matching technique?
         */
-        class mlstudent : public mlnorm<ldmat>
+        class mlstudent : public mlnorm<ldmat,enorm>
         {
                 protected:
@@ -651 +646 @@
                         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 +653 @@
                                 it_assert_debug ( R0.rows() ==A0.rows(),"" );
 
-                                iepdf->set_parameters ( mu0,Lambda ); //
+                                iepdf.set_parameters ( mu0,Lambda ); //
                                 A = A0;
                                 mu_const = mu0;
@@ -666 +661 @@
                         void condition ( const vec &cond )
                         {
-                                _mu = A*cond + mu_const;
+                                iepdf._mu() = A*cond + mu_const;
                                 double zeta;
                                 //ugly hack!
@@ -691 +686 @@
         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 +698 @@
                 public:
                         //! Constructor
-                        mgamma():iepdf(new egamma()), k(0),
-                            _beta(iepdf->_beta()) {
-                            set_ep(iepdf);
+                        mgamma():mpdf_internal<egamma>(), k(0),
+                            _beta(iepdf._beta()) {
                         }
 
@@ -742 +734 @@
         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 +749 @@
                         //! \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 +766 @@
                         {
                                 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 +902 @@
         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 +912 @@
                 public:
                         //! Constructor
-                        mlognorm():eno(new elognorm()),
+                        mlognorm():mpdf_internal<elognorm>(),
                             sig2(0),
-                            mu(eno->_mu()) {
-                            set_ep(eno);
+                            mu(iepdf._mu()) {
                         }
 
@@ -937 +921 @@
                         {
                                 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 +1055 @@
         {
                 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 +1064 @@
 
                 public:
-                        rwiWishartCh():eiW(new eiWishartCh()),
+                        rwiWishartCh():eiW(),
                             sqd(0), l(0), p(0) {
                             set_ep(eiW);
@@ -1093 +1077 @@
                                 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 +1088 @@
                                 }
 
-                                eiW->_setY ( sqd*z );
+                                eiW._setY ( sqd*z );
                         }
         };
@@ -1256 +1240 @@
         };
 
-        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 +1266 @@
 // }
 
-        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 +1315 @@
         }
 
-///////////
+        ////
+        ///////
+        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>

library/bdm/stat/merger.cpp

@@ -152 +152 @@
                         if ( mpdfs ( i )->dimension() == dim ) {
                                 // no need for conditioning or marginalization
-                                lw_src = mpdfs ( i )->e()->evallog_m ( Smp );
+                                lw_src = mpdfs ( i )->evallogcond_m ( Smp , vec(0));
                         } else {
                                 // compute likelihood of marginal on the conditional variable
@@ -190 +190 @@
                                 // Compute likelihood of the partial source
                                 for ( int k = 0; k < Npoints; k++ ) {
-                                        mpdfs ( i )->condition ( dls ( i )->get_cond ( Smp ( k ) ) );
-                                        lw_src ( k ) += mpdfs ( i )->e()->evallog ( dls ( i )->pushdown ( Smp ( k ) ) );
+                                        lw_src ( k ) += mpdfs ( i )->evallogcond ( dls ( i )->pushdown ( Smp ( k ) ),  
+                                                         dls ( i )->get_cond ( Smp ( k ) ));
                                 }