Changeset 388

Timestamp:

06/19/09 11:43:48 (15 years ago)

Author:

smidl

Message:

mergers

Files:

• applications/bdmtoolbox/mex/merger_mx.cpp (modified) (2 diffs)
• library/bdm/base/bdmbase.h (modified) (1 diff)
• library/bdm/base/user_info.h (modified) (1 diff)
• library/bdm/mex/mex_parser.h (modified) (1 diff)
• library/bdm/stat/exp_family.h (modified) (1 diff)
• library/bdm/stat/merger.h (modified) (8 diffs)

applications/bdmtoolbox/mex/merger_mx.cpp

@@ -10 +10 @@
 {
     // Check the number of inputs and output arguments
-        if(n_input!=1) mexErrMsgTxt("Configuration structure expected in the form:\n"  
-                "S.sources = { };                          % cell of pdfs (epdfs or mpdfs) to be merged,\n"
-                "S.merger= struct('class','merger_*');     % object to be used for merging,\n\n"
+        if(n_input!=3) mexErrMsgTxt("Usage:\n"  
+                "result=merger(sources, support, merger)\n"
+                "  sources= { struct('class','epdf'),... };  % cell of pdfs (epdfs or mpdfs) to be merged,\n"
+                "  support= struct(\n"
+                "           grid    = {[dim1_start,dim1_end], [dim2_start, dim2_end]...}  %support boundary \n"   
+                "           nbins   = [bins_in_dim1, bins_in_dim2,...]                    %fixed \n"   
+                "         === OR ==\n"
+                "           pdf     = struct('class','epdf'); % pdf to draw samples from\n"
+                "           nsamples= 100;                    % number of samples\n"
+                "           );\n"
+                "        If all elements are present,  (grid,nbins) is used;\n"
+                "  merger = struct('class','merger_*');       % object to be used for merging,\n\n"
                 "see documentation of classes epdf, mpdf, merger_base and their offsprings in BDM.");
 
         // LOAD CONFIG
-        UImxArray F (input[0]);
+        UImxArray Src;
+        Src.addList(input[0],"Sources");
+        
+        UImxArray Sup (input[1]);
+        UImxArray Mer (input[2]);
+        // Sources
+        Array<mpdf*> Sources;
+        //abuse Mer to store sources
 
-        Array<mpdf*> Sources;
-
-        Setting& _Sources=F.lookup("sources");
+        Setting& _Sources=Src.lookup("Sources");
         int Slen=_Sources.getLength();
         Sources.set_size(Slen);
@@ -35 +49 @@
 
         }
+
+        merger_base* Merger=UI::build<merger_base>(Mer.getRoot());
+
+        // Support
+        Setting & _Supp=Sup.getRoot();
         
-        merger_base* Merger=UI::build<merger_base>(F.getRoot(),"merger");
-
+        Array<vec> bounds (0);
+        UI::get (bounds, _Supp, "grid");
+        ivec nbins(0);
+        UI::get (nbins, _Supp, "nbins");
+        
+        epdf *g0=UI::build<epdf> (_Supp, "pdf");
+        int npoints=100;
+        _Supp.lookupValue("nsamples",npoints);
+        
 // COMPUTE RESULTS
-
+        Merger->set_sources(Sources,true); // takes care of deletion of sources
+        if (bounds.length() > 0 && nbins.length()>0) {
+                Merger->set_support (bounds,nbins);
+        } else {
+                if (g0) {
+                        Merger->set_support (*g0,npoints);
+                        delete g0;
+                }
+        }
         Merger->merge();
         
+        // Save results
         if (n_output>0){
                 mxArray* tmp = mxCreateStructMatrix(1,1,0,NULL);

library/bdm/base/bdmbase.h

@@ -581 +581 @@
 class compositepdf {
 protected:
-  //!Number of mpdfs in the composite
-  int n;
   //! Elements of composition
   Array<mpdf*> mpdfs;
-public:
-  compositepdf(Array<mpdf*> A0) : n(A0.length()), mpdfs(A0) {};
+  bool owning_mpdfs;
+public:
+        compositepdf():mpdfs(0){};
+        compositepdf(Array<mpdf*> A0, bool own=false){set_elements(A0,own);};
+        void set_elements(Array<mpdf*> A0, bool own=false) {mpdfs=A0;owning_mpdfs=own;};
   //! find common rv, flag \param checkoverlap modifies whether overlaps are acceptable
   RV getrv(bool checkoverlap = false);
   //! common rvc of all mpdfs is written to rvc
   void setrvc(const RV &rv, RV &rvc);
+  ~compositepdf(){if (owning_mpdfs) for(int i=0;i<mpdfs.length();i++){delete mpdfs(i);}};
 };
 

library/bdm/base/user_info.h

@@ -235 +235 @@
                 return instance;
         }
+        //! VS: addition for root elements
+        template<class T> static T* build( const Setting &element )
+        {
+                T* instance;
+                from_setting<T>( instance,  element );
+                return instance;
+        }
 
         template<class T> static T* build( const Setting &element, const string &name )

library/bdm/mex/mex_parser.h

@@ -17 +17 @@
                 setAutoConvert(true);
         }
-
+        UImxArray() : Config() {                
+                setAutoConvert(true);
+        }       
+        void addList(const mxArray *mxarray, const char* name){
+                Setting & setting = this->getRoot(); //setting is a group
+                Setting & child = setting.add(name,Setting::TypeList);
+                fillList(child,mxarray);
+        }
 private:
         void storeNumeric(Setting &setting, const mxArray *value, string key = "") {

library/bdm/stat/exp_family.h

@@ -155 +155 @@
 
         };
+        UIREGISTER(enorm<chmat>);
+        UIREGISTER(enorm<ldmat>);
+        UIREGISTER(enorm<fsqmat>);
+
 
         /*!

library/bdm/stat/merger.h

@@ -75 +75 @@
 
                 //!Empty constructor
-                merger_base () : compositepdf() {DBG=false;dbg_file=NULL;};
+                merger_base () : compositepdf() {DBG = false;dbg_file = NULL;};
                 //!Constructor from sources
-                merger_base (const Array<mpdf*> &S) {set_sources (S);};
-                //! Function setting the main internal structures 
-                void set_sources (const Array<mpdf*> &Sources) {
-                        compositepdf::set_elements (Sources);
+                merger_base (const Array<mpdf*> &S, bool own=false) {set_sources (S,own);};
+                //! Function setting the main internal structures
+                void set_sources (const Array<mpdf*> &Sources, bool own) {
+                        compositepdf::set_elements (Sources,own);
                         //set sizes
                         dls.set_size (Sources.length());
@@ -109 +109 @@
                         };
                 }
+                //! Rectangular support  each vector of XYZ specifies (begining-end) interval for each dimension. Same number of points, \c dimsize, in each dimension. 
+                void set_support (const Array<vec> &XYZ, const int dimsize) {
+                        set_support(XYZ,dimsize*ones_i(XYZ.length()));
+                }
+                //! Rectangular support  each vector of XYZ specifies (begining-end) interval for each dimension. \c gridsize specifies number of points is each dimension.
+                void set_support (const Array<vec> &XYZ, const ivec &gridsize) {
+                        int dim = XYZ.length();  //check with internal dim!!
+                        Npoints = prod (gridsize); 
+                        eSmp.set_parameters (Npoints, false);
+                        Array<vec> &samples = eSmp._samples();
+                        eSmp._w() = ones (Npoints) / Npoints; //unifrom size of bins
+                        //set samples
+                        ivec ind = zeros_i (dim);      //indeces of dimensions in for cycle;
+                        vec smpi (dim);            // ith sample
+                        vec steps =zeros(dim);            // ith sample
+                        // first corner
+                        for (int j = 0; j < dim; j++) {
+                                smpi (j) = XYZ (j) (0); /* beginning of the interval*/ 
+                                it_assert(gridsize(j)!=0.0,"Zeros in gridsize!");
+                                steps (j) = (smpi(j) - XYZ(j)(1))/gridsize(j);
+                        }
+                        // fill samples
+                        int act_dim=0; //active dimension
+                        for (int i = 0; i < Npoints; i++) {
+                                // copy 
+                                samples(i) = smpi; 
+                                // go through all dimensions
+                                for (int j = 0;j < dim;j++) {
+                                        if (ind (j) == gridsize (j) - 1) { //j-th index is full
+                                                ind (j) = 0; //shift back
+                                                smpi(j) = XYZ(j)(0);
+                                                
+                                                ind (j + 1) ++; //increase the next dimension;
+                                                smpi(j+1) += steps(j+1);
+                                                
+                                                if (ind (j + 1) < gridsize (j + 1) - 1) break;
+                                        } else {
+                                                ind (j) ++; 
+                                                smpi(j) +=steps(j);
+                                                break;
+                                        }
+                                }
+                        }
+                }
                 //! set debug file
-                void set_debug_file (const string fname) { 
-                        if (DBG) delete dbg_file; 
-                        dbg_file = new it_file (fname); 
+                void set_debug_file (const string fname) {
+                        if (DBG) delete dbg_file;
+                        dbg_file = new it_file (fname);
                         if (dbg_file) DBG = true;
                 }
                 //! Set internal parameters used in approximation
-                void set_method (MERGER_METHOD MTH, double beta0=0.0) {
-                        METHOD = MTH; 
+                void set_method (MERGER_METHOD MTH, double beta0 = 0.0) {
+                        METHOD = MTH;
                         beta = beta0;
                 }
                 //! Set support points from a pdf by drawing N samples
-                void set_support(const epdf &overall, int N){
-                        it_assert_debug ( rv.equal ( overall._rv() ),"Incompatible parameter overall!" );
-                        eSmp.set_statistics(&overall,N);
-                        Npoints=N;
-                }
-                
+                void set_support (const epdf &overall, int N) {
+                        eSmp.set_statistics (&overall, N);
+                        Npoints = N;
+                }
+
                 //! Destructor
                 virtual ~merger_base() {
@@ -136 +179 @@
                 };
                 //!@}
-                
+
                 //! \name Mathematical operations
                 //!@{
-                
+
                 //!Merge given sources in given points
                 void merge () {
-                        if (eSmp._w().length() ==0) {it_error("Empty support points use set_support" );}
+                        validate();
+
                         //check if sources overlap:
                         bool OK = true;
@@ -164 +208 @@
                                 eSmp._w() = wtmp / sum (wtmp);
                         } else {
-                                it_error("Sources are not compatible - use merger_mix");
-                        }
-                        ;
+                                it_error ("Sources are not compatible - use merger_mix");
+                        }
                 };
 
@@ -172 +215 @@
                 //! Merge log-likelihood values in points using method specified by parameter METHOD
                 vec merge_points (mat &lW);
-                
-                
+
+
                 //! sample from merged density
 //! weight w is a
@@ -216 +259 @@
                 //! Access function
                 eEmp& _Smp() {return eSmp;}
-                
+
+                //! load from setting
+                void from_setting (const Setting& set) {
+                        // get support
+                        // find which method to use
+                        string meth_str;
+                        UI::get<string> (meth_str, set, "method");
+                        if (!strcmp (meth_str.c_str(), "arithmetic"))
+                                set_method (ARITHMETIC);
+                        else {
+                                if (!strcmp (meth_str.c_str(), "geometric"))
+                                        set_method (GEOMETRIC);
+                                else if (!strcmp (meth_str.c_str(), "lognormal")) {
+                                        set_method (GEOMETRIC);
+                                        set.lookupValue( "beta",beta);
+                                }
+                        }
+                        validate();
+                }
+
+                void validate() {
+                        it_assert (eSmp._w().length() > 0, "Empty support, use set_support().");
+                        it_assert (dim == eSmp._samples() (0).length(), "Support points and rv are not compatible!");
+                        it_assert (isnamed(),"mergers must be named");
+                }
                 //!@}
 };
+UIREGISTER(merger_base);
 
 class merger_mix : public merger_base
@@ -234 +302 @@
                 //!@{
                 merger_mix () {};
-                merger_mix (const Array<mpdf*> &S) {set_sources(S);};
+                merger_mix (const Array<mpdf*> &S,bool own=false) {set_sources (S,own);};
                 //! Set sources and prepare all internal structures
-                void set_sources (const Array<mpdf*> &S) {
-                        merger_base::set_sources(S);
+                void set_sources (const Array<mpdf*> &S, bool own) {
+                        merger_base::set_sources (S,own);
                         Nsources = S.length();
                 }
                 //! Set internal parameters used in approximation
-                void set_parameters (int Ncoms0=10, double effss_coef0 = 0.5) {
+                void set_parameters (int Ncoms0 = 10, double effss_coef0 = 0.5) {
                         Ncoms = Ncoms0;
                         effss_coef = effss_coef0;
                 }
                 //!@}
-                
+
                 //! \name Mathematical operations
                 //!@{
-                
+
                 //!Merge values using mixture approximation
                 void merge ();
@@ -269 +337 @@
                 //! Access function
                 emix* proposal() {emix* tmp = Mix.epredictor(); tmp->set_rv (rv); return tmp;}
-        //! @}
+                //! from_settings
+                void from_settings(const Setting& set){
+                        merger_base::from_setting(set);
+                        set.lookupValue("ncoms",Ncoms);
+                }
                 
+                //! @}
+
 };
+UIREGISTER(merger_mix);
 
 }