Show
Ignore:
Timestamp:
07/02/09 22:16:23 (15 years ago)
Author:
smidl
Message:

arx example for estimator workig

Files:
1 modified

Legend:

Unmodified
Added
Removed

  • applications/bdmtoolbox/mex/estimator.cpp

    r394 r412  
    1 #include <itpp/itmex.h> 
     1/*! 
     2\file 
     3\brief Application Estimator 
     4 
     5The general task of estimation is defined on the following scheme: 
     6\dot 
     7digraph estimation{ 
     8        node [shape=box]; 
     9        {rank="same"; "Data Source"; "Bayesian Model"} 
     10        "Data Source" -> "Bayesian Model" [label="data"]; 
     11        "Bayesian Model" -> "Result Logger" [label="estimated\n statistics"];  
     12        "Data Source" -> "Result Logger" [label="Simulated\n data"]; 
     13} 
     14\enddot 
     15 
     16Here,  
     17\li Data Source is an object (class DS) providing sequential data, \f$ [d_1, d_2, \ldots d_t] \f$. 
     18\li Bayesian Model is an object (class BM) performing Bayesian filtering, 
     19\li Result Logger is an object (class logger) dedicated to storing important data from the experiment. 
     20 
     21\section  cmd Command-line usage 
     22Execute command: 
     23\code 
     24$> estimator config_file.cfg 
     25\endcode 
     26 
     27Full description of the experiment is in the file config_file.cfg which is expected to have the following structure: 
     28\code 
     29system = {type = "DS_offspring", ...};      // definition of a data source 
     30estimator = {type = "BM_offspring", ...};   // definition of an estimator 
     31logger = {type = "logger_type",...};        // definition of a logger 
     32experiment = {ndat = 11000; };              // definition of number of data records 
     33\endcode 
     34 
     35The above description must be specialized to specific classes. See, \subpage arx_ui how to do it for estimation of an ARX model. 
     36 
     37\section ex Matlab usage 
     38Execute command: 
     39\code 
     40>> estimator('config_file.cfg'); 
     41\endcode 
     42when using loggers storing results on hard drives, and 
     43\code 
     44>> Res=estimator('config_file.cfg'); 
     45\endcode 
     46when using logger of the type \c "mex_logger". The results will be stored in structure \c M. 
     47 
     48 */ 
    249 
    350#include "estim/arx.h" 
    4 #include "stat/datasources.h" 
    5 #include "stat/loggers.h" 
    6 #include "mex_logger.h" 
     51#include "base/datasources.h" 
     52#include "base/loggers.h" 
    753 
    854//#include "mex_datasource.h" 
    … …  
    1056using namespace bdm; 
    1157 
     58#ifdef MEX 
     59#include <itpp/itmex.h> 
     60#include "mex/mex_logger.h" 
     61#include "mex/mex_parser.h" 
     62 
    1263void mexFunction ( int n_output, mxArray *output[], int n_input, const mxArray *input[] ) { 
    1364        // Check the number of inputs and output arguments 
    14         string fname; 
    15         if (n_input>0){ 
    16                 fname=mxArray2string(input[0]); 
     65        if(n_input<3) mexErrMsgTxt("Usage:\n"   
     66                "result=estimator(system, estimators, experiment, logger)\n" 
     67                "  system     = struct('class','datasource',...);  % Estimated system\n" 
     68                "  estimators = {struct('class','estimator',...),  % Estimators\n" 
     69                "                struct('class','estimator',...),...} \n" 
     70                "  experiment = struct('ndat',100);                % number of data in experiment \n" 
     71                "  === optional ===" 
     72                "  logger     = struct('class','mexlogger');       % How to store results, default=mexlog, i.e. matlab structure\n\n" 
     73                "see documentation of classes datasource, BM, and mexlogger and their offsprings in BDM."); 
     74         
     75        UImxArray Cfg; 
     76        try{ 
     77        Cfg.addGroup(input[0],"system"); 
     78        Cfg.addList(input[1],"estimators"); 
     79        Cfg.addGroup(input[2],"experiment"); 
     80        if (n_input>3){ 
     81                Cfg.addGroup(input[3],"logger"); 
     82        }else{ 
     83                // define logger as mexlog 
     84                Setting &S=Cfg.getRoot(); 
     85                S.add("logger",Setting::TypeGroup); 
     86                S["logger"].add("class",Setting::TypeString); 
     87                S["logger"]["class"]="mexlog"; 
     88                S["logger"].add("maxlen",Setting::TypeInt); 
     89                int maxlen; 
     90                S["experiment"].lookupValue("ndat",maxlen); 
     91                S["logger"]["maxlen"]=maxlen; 
     92        } 
     93        } catch(SettingException e){it_error("error: "+string(e.getPath()));} 
     94 
     95        //DBG 
     96        Cfg.writeFile("estimator.cfg"); 
     97 
     98#else 
     99int main ( int argc, char* argv[] ) { 
     100        const char *fname; 
     101        if ( argc>1 ) { 
     102                fname = argv[1]; 
    17103        } else { 
    18                 fname = "UIArxDS_test.cfg"; 
     104                fname="estimator.cfg"; 
    19105        } 
    20         // ------------------  
     106        UIFile Cfg ( fname ); 
     107#endif 
    21108         
    22         printf("name: %s", fname.c_str()); 
    23         UIFile F ( fname.c_str()); 
    24  
    25         logger* L = UI::build<logger>( F, "logger"); 
    26         ArxDS * DS = UI::build<ArxDS>( F, "system" ); 
    27         BM* E = UI::build<BM>( F, "estimator" ); 
    28         int Ndat = F.lookupValue ( "experiment.ndat",Ndat ); 
     109        logger* L = UI::build<logger>( Cfg, "logger"); 
     110        ArxDS * DS = UI::build<ArxDS>( Cfg, "system" ); 
     111        Array<BM*> Es;  UI::get(Es,Cfg, "estimators" ); 
     112        int Ndat; 
     113        Cfg.lookupValue ( "experiment.ndat",Ndat ); 
    29114 
    30115        DS->log_add ( *L ); 
    31         int L_est= L->add ( E->posterior()._rv(), "est" ); // estimate 
    32         int L_lb = L->add ( E->posterior()._rv(), "lb" ); // lower bound 
    33         int L_ub = L->add ( E->posterior()._rv(), "ub" ); // upper bound 
     116        string Ename; 
     117        Setting &S=Cfg; 
     118        for (int i=0; i<Es.length(); i++){ 
     119                try{ 
     120                        UI::get(Ename, S["estimators"][i], "name"); 
     121                } catch (UIException e){ 
     122                        Ename="Est"+num2str(i); 
     123                } 
     124                 
     125                Es(i)->log_add(*L,Ename); // estimate 
     126        } 
    34127        L->init(); 
    35128 
    36129        vec dt=zeros ( DS->_drv()._dsize() );   //data variable 
    37         datalink dl ( E->_drv(),DS->_drv() ); //datalink between a datasource and estimator 
    38  
     130        Array<datalink*> Dls(Es.length());  
     131        for (int i=0; i<Es.length(); i++){ 
     132                Dls(i)=new datalink( Es(i)->_drv(),DS->_drv() ); //datalink between a datasource and estimator 
     133        } 
     134         
    39135        for ( int tK=1;tK<Ndat;tK++ ) { 
    40136                DS->step();                                                     // simulator step 
    41137                DS->getdata ( dt );                                     // read data 
    42                 E->bayes ( dl.pushdown ( dt ) );                // update estimates 
    43  
    44138                DS->logit ( *L ); 
    45                 L->logit ( L_est, E->posterior().mean() ); 
    46                 L->logit ( L_lb,  E->posterior().mean()-2*sqrt ( E->posterior().variance() ) ); 
    47                 L->logit ( L_ub,  E->posterior().mean() +2*sqrt ( E->posterior().variance() ) ); 
    48  
     139                 
     140                for (int i=0; i<Es.length(); i++){ 
     141                        Es(i)->bayes ( Dls(i)->pushdown ( dt ) );               // update estimates 
     142                        Es(i)->logit (*L); 
     143                } 
    49144                L->step(); 
    50145        } 
    51146 
    52147        L->finalize(); 
     148        // ------------------ End of routine ----------------------------- 
    53149 
    54         // ------------------ End of routine ----------------------------- 
    55          
    56         mex_logger* mL=dynamic_cast<mex_logger*>(L); 
     150#ifdef MEX 
     151        mexlog* mL=dynamic_cast<mexlog*>(L); 
    57152 
    58153        if (mL) { // user wants output!! 
    … …  
    60155                output[0] = mL->toCell(); 
    61156        } 
    62          
     157#endif 
    63158        /////// 
    64159        delete L; 
    65160        delete DS; 
    66         delete E; 
     161        for (int i;i<Es.length();i++){delete Es(i);delete Dls(i);} 
    67162}