Changeset 286 for bdm/stat/emix.h

Timestamp:

03/05/09 14:03:35 (15 years ago)

Author:

smidl

Message:

make mpdm work again

Files:

• bdm/stat/emix.h (modified) (12 diffs)

bdm/stat/emix.h

@@ -18 +18 @@
 //#include <std>
 
-namespace bdm{
+namespace bdm {
 
 //this comes first because it is used inside emix!
@@ -47 +47 @@
         //!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 ( rv,rvc,nom0->_rv() ) {
-                ep->set_rv(_rv());
-                set_rvc(nom0->_rv().subt ( rv ) );
-                if ( copy ) {
-//                      nom = nom0->_copy_();
-                        it_error ( "todo" );
-                        destroynom=true;
-                }
-                else {
-                        nom = nom0;
-                        destroynom = false;
-                }
-                it_assert_debug ( rvc.length() >0,"Makes no sense to use this object!" );
+        mratio ( const epdf* nom0, const RV &rv, bool copy=false ) :mpdf ( ), dl ( ) {
+                // adjust rv and rvc
+                rvc = nom0->_rv().subt ( rv );
+                dimc = rvc._dsize();
+                ep = new epdf;
+                ep->set_parameters(rv._dsize());
+                ep->set_rv(rv);
+                
+                //prepare data structures
+                if ( copy ) {it_error ( "todo" ); destroynom=true; }
+                else { nom = nom0; destroynom = false; }
+                it_assert_debug ( rvc.length() >0,"Makes no sense to use this object!" );               
+                
+                // build denominator
                 den = nom->marginal ( rvc );
+                dl.set_connection(rv,rvc,nom0->_rv());
         };
         double evallogcond ( const vec &val, const vec &cond ) {
@@ -67 +69 @@
                 dl.pushup_cond ( nom_val,val,cond );
                 tmp = exp ( nom->evallog ( nom_val ) - den->evallog ( cond ) );
-                it_assert_debug(std::isfinite(tmp),"Infinite value");
+                it_assert_debug ( std::isfinite ( tmp ),"Infinite value" );
                 return tmp;
         }
@@ -96 +98 @@
 public:
         //!Default constructor
-        emix (  ) : epdf ( ) {};
+        emix ( ) : epdf ( ) {};
         //! Set weights \c w and components \c Coms
         //!By default Coms are copied inside. Parameter \c copy can be set to false if Coms live externally. Use method ownComs() if Coms should be destroyed by the destructor.
@@ -103 +105 @@
         vec sample() const;
         vec mean() const {
-                int i; vec mu = zeros ( dim);
+                int i; vec mu = zeros ( dim );
                 for ( i = 0;i < w.length();i++ ) {mu += w ( i ) * Coms ( i )->mean(); }
                 return mu;
@@ -109 +111 @@
         vec variance() const {
                 //non-central moment
-                vec mom2 = zeros(dim);
-                for ( int i = 0;i < w.length();i++ ) {mom2 += w ( i ) * pow(Coms ( i )->mean(),2); }
+                vec mom2 = zeros ( dim );
+                for ( int i = 0;i < w.length();i++ ) {mom2 += w ( i ) * pow ( Coms ( i )->mean(),2 ); }
                 //central moment
-                return mom2-pow(mean(),2);
+                return mom2-pow ( mean(),2 );
         }
         double evallog ( const vec &val ) const {
@@ -118 +120 @@
                 double sum = 0.0;
                 for ( i = 0;i < w.length();i++ ) {sum += w ( i ) * exp ( Coms ( i )->evallog ( val ) );}
-                if (sum==0.0){sum=std::numeric_limits<double>::epsilon();}
+                if ( sum==0.0 ) {sum=std::numeric_limits<double>::epsilon();}
                 double tmp=log ( sum );
-                it_assert_debug(std::isfinite(tmp),"Infinite");
+                it_assert_debug ( std::isfinite ( tmp ),"Infinite" );
                 return tmp;
         };
@@ -151 +153 @@
         //!access function
         epdf* _Coms ( int i ) {return Coms ( i );}
+        void set_rv(const RV &rv){
+                epdf::set_rv(rv);
+                for(int i=0;i<Coms.length();i++){Coms(i)->set_rv(rv);}
+        }
 };
 
@@ -174 +180 @@
         mprod ( Array<mpdf*> mFacs ) : compositepdf ( mFacs ), mpdf (), epdfs ( n ), dls ( n ) {
                 ep=&dummy;
-                RV rv=getrv(true);
-                set_rv(rv);dummy.set_parameters(rv._dsize());
+                RV rv=getrv ( true );
+                set_rv ( rv );dummy.set_parameters ( rv._dsize() );
                 setrvc ( ep->_rv(),rvc );
                 // rv and rvc established = > we can link them with mpdfs
                 for ( int i = 0;i < n;i++ ) {
-                        dls ( i ) = new datalink_m2m ( mpdfs ( i )->_rv(), mpdfs ( i )->_rvc(), _rv(), _rvc() );
+                        dls ( i ) = new datalink_m2m;
+                        dls(i)->set_connection( mpdfs ( i )->_rv(), mpdfs ( i )->_rvc(), _rv(), _rvc() );
                 }
 
@@ -237 +244 @@
         Array<datalink*> dls;
 public:
-        eprod ( const Array<const epdf*> epdfs0 ) : epdf (  ),epdfs ( epdfs0 ),dls ( epdfs.length() ) {
+        eprod () : epdfs ( 0 ),dls ( 0 ) {};
+        void set_parameters ( const Array<const epdf*> &epdfs0, bool named=true ) {
+                epdfs=epdfs0;//.set_length ( epdfs0.length() );
+                dls.set_length ( epdfs.length() );
+
                 bool independent=true;
-                for ( int i=0;i<epdfs.length();i++ ) {
-                        independent=rv.add ( epdfs ( i )->_rv() );
-                        it_assert_debug ( independent==true, "eprod:: given components are not independent ." );
-                }
-                for ( int i=0;i<epdfs.length();i++ ) {
-                        dls ( i ) = new datalink ( epdfs ( i )->_rv() , rv );
+                if ( named ) {
+                        for ( int i=0;i<epdfs.length();i++ ) {
+                                independent=rv.add ( epdfs ( i )->_rv() );
+                                it_assert_debug ( independent==true, "eprod:: given components are not independent." );
+                        }
+                        dim=rv._dsize();
+                }
+                else {
+                        dim =0; for ( int i=0;i<epdfs.length();i++ ) {
+                                dim+=epdfs ( i )->dimension();
+                        }
+                }
+                //
+                int cumdim=0;
+                int dimi=0;
+                int i;
+                for ( i=0;i<epdfs.length();i++ ) {
+                        dls ( i ) = new datalink;
+                        if ( named ) {dls ( i )->set_connection ( epdfs ( i )->_rv() , rv );}
+                        else {
+                                dimi = epdfs ( i )->dimension();
+                                dls ( i )->set_connection ( dimi, dim, linspace ( cumdim,cumdim+dimi-1 ) );
+                                cumdim+=dimi;
+                        }
                 }
         }
@@ -260 +289 @@
                 for ( int i=0;i<epdfs.length();i++ ) {
                         vec pom = epdfs ( i )->mean();
-                        dls ( i )->pushup ( tmp, pow(pom,2) );
-                }
-                return tmp-pow(mean(),2);
+                        dls ( i )->pushup ( tmp, pow ( pom,2 ) );
+                }
+                return tmp-pow ( mean(),2 );
         }
         vec sample() const {
@@ -277 +306 @@
                         tmp+=epdfs ( i )->evallog ( dls ( i )->pushdown ( val ) );
                 }
-                it_assert_debug(std::isfinite(tmp),"Infinite");
+                it_assert_debug ( std::isfinite ( tmp ),"Infinite" );
                 return tmp;
         }