Changeset 271 for bdm/stat/libDS.h

Show
Ignore:
Timestamp:
02/16/09 10:03:13 (15 years ago)
Author:
smidl
Message:

Next major revision

Files:
1 modified

Legend:

Unmodified
Added
Removed

  • bdm/stat/libDS.h

    r270 r271  
    2020 
    2121namespace bdm { 
    22         /*! 
    23         * \brief Memory storage of off-line data column-wise 
     22/*! 
     23* \brief Memory storage of off-line data column-wise 
    2424 
    25         The data are stored in an internal matrix \c Data . Each column of Data corresponds to one discrete time observation \f$t\f$. Access to this matrix is via indices \c rowid and \c delays. 
     25The data are stored in an internal matrix \c Data . Each column of Data corresponds to one discrete time observation \f$t\f$. Access to this matrix is via indices \c rowid and \c delays. 
    2626 
    27         The data can be loaded from a file. 
    28         */ 
    29         class MemDS : public DS { 
    30                 //! internal matrix of data 
    31                 mat Data; 
    32                 //! active column in the Data matrix 
    33                 int time; 
    34                 //!  vector of rows that are presented in Dt 
    35                 ivec rowid; 
    36                 //! vector of delays that are presented in Dt 
    37                 ivec delays; 
     27The data can be loaded from a file. 
     28*/ 
     29class MemDS : public DS { 
     30        //! internal matrix of data 
     31        mat Data; 
     32        //! active column in the Data matrix 
     33        int time; 
     34        //!  vector of rows that are presented in Dt 
     35        ivec rowid; 
     36        //! vector of delays that are presented in Dt 
     37        ivec delays; 
    3838 
    39         public: 
    40                 void getdata ( vec &dt ); 
    41                 void getdata ( vec &dt, const ivec &indeces ); 
    42                 void set_rvs ( RV &drv, RV &urv ); 
    43                 void write ( vec &ut ) {it_error ( "MemDS::write is not supported" );} 
    44                 void write ( vec &ut,ivec &indices ) {it_error ( "MemDS::write is not supported" );} 
    45                 void step(); 
    46                 //!Default constructor 
    47                 MemDS ( mat &Dat, ivec &rowid, ivec &delays ); 
     39public: 
     40        void getdata ( vec &dt ); 
     41        void getdata ( vec &dt, const ivec &indeces ); 
     42        void set_rvs ( RV &drv, RV &urv ); 
     43        void write ( vec &ut ) {it_error ( "MemDS::write is not supported" );} 
     44        void write ( vec &ut,ivec &indices ) {it_error ( "MemDS::write is not supported" );} 
     45        void step(); 
     46        //!Default constructor 
     47        MemDS ( mat &Dat, ivec &rowid, ivec &delays ); 
     48}; 
     49 
     50/*! 
     51\brief Generator of ARX data 
     52 
     53*/ 
     54class ArxDS : public DS { 
     55protected: 
     56        //! Rv of the regressor 
     57        RV Rrv; 
     58        //! History, ordered as \f$[y_t, u_t, y_{t-1 }, u_{t-1}, \ldots]\f$ 
     59        vec H; 
     60        //! (future) input 
     61        vec U; 
     62        //! temporary variable for regressor 
     63        vec rgr; 
     64        //! data link: H -> rgr 
     65        datalink rgrlnk; 
     66        //! model of Y - linear Gaussian 
     67        mlnorm<chmat> model; 
     68        //! options 
     69        bool opt_L_theta; 
     70        //! loggers 
     71        int L_theta; 
     72        int L_R; 
     73        int dt_size; 
     74public: 
     75        void getdata ( vec &dt ) { 
     76                //it_assert_debug ( dt.length() ==Drv.count(),"ArxDS" ); 
     77                dt=H; 
     78        }; 
     79        void getdata ( vec &dt, const ivec &indices ) { 
     80                it_assert_debug ( dt.length() ==indices.length(),"ArxDS" ); 
     81                dt=H ( indices ); 
     82        }; 
     83        void write ( vec &ut ) { 
     84                //it_assert_debug ( ut.length() ==Urv.count(),"ArxDS" ); 
     85                U=ut; 
     86        }; 
     87        void write ( vec &ut, const ivec &indices ) { 
     88                it_assert_debug ( ut.length() ==indices.length(),"ArxDS" ); 
     89                set_subvector ( U, indices,ut ); 
     90        }; 
     91        void step(); 
     92        //!Default constructor 
     93        ArxDS ( ) {}; 
     94        //! Set parameters of the internal model, H is maximum time delay 
     95        void set_parameters ( const mat &Th0, const vec mu0, const chmat &sqR0 ) 
     96        { model.set_parameters ( Th0, mu0, sqR0 );}; 
     97        //! Set 
     98        void set_drv(RV &yrv, RV &urv, RV &rrv){ 
     99                Rrv = rrv; 
     100                Urv = urv; 
     101                dt_size = yrv._dsize()+urv._dsize(); 
     102                 
     103                RV drv = concat(yrv,urv); 
     104                Drv = drv; 
     105                int td = rrv.mint(); 
     106                H.set_size(drv._dsize()*(-td+1)); 
     107                U.set_size(Urv._dsize()); 
     108                for (int i=-1;i>=td;i--){ 
     109                        drv.t(-1); 
     110                        Drv.add(drv); //shift u1 
     111                } 
     112                rgrlnk.set_connection(rrv,Drv); 
     113                 
     114                dtsize = Drv._dsize(); 
     115                utsize = Urv._dsize(); 
     116        } 
     117        //! set options from a string 
     118        void set_options ( const string &s ) { 
     119                opt_L_theta= ( s.find ( "L_theta" ) !=string::npos ); 
     120        }; 
     121        virtual void log_add ( logger &L ) { 
     122                //DS::log_add ( L ); too long!! 
     123                L_dt=L.add ( Drv(0,dt_size),"" ); 
     124                L_ut=L.add ( Urv,"" ); 
     125 
     126                mat &A =model._A(); 
     127                mat R =model._R(); 
     128                if ( opt_L_theta ) {L_theta=L.add ( RV ( "{th }", vec_1 ( A.rows() *A.cols() ) ),"t" );} 
     129                if ( opt_L_theta ) {L_R=L.add ( RV ( "{R }", vec_1 ( R.rows() *R.cols() ) ),"r" );} 
     130        } 
     131        virtual void logit ( logger &L ) { 
     132                //DS::logit ( L ); 
     133                L.logit( L_dt, H.left(dt_size)); 
     134                L.logit(L_ut, U); 
     135                 
     136                mat &A =model._A(); 
     137                mat R =model._R(); 
     138                if ( opt_L_theta ) {L.logit ( L_theta,vec ( A._data(), A.rows() *A.cols() ) );}; 
     139                if ( opt_L_theta ) {L.logit ( L_R, vec ( R._data(), R.rows() *R.rows() ) );}; 
     140        } 
     141 
     142}; 
     143 
     144class stateDS : public DS { 
     145protected: 
     146        //!conditional pdf of the state evolution \f$ f(x_t|x_{t-1}) \f$ 
     147        mpdf* IM; 
     148        //!conditional pdf of the observations \f$ f(d_t|x_t) \f$ 
     149        mpdf* OM; 
     150        //! result storage 
     151        vec dt; 
     152        //! state storage 
     153        vec xt; 
     154        //! input storage 
     155        vec ut; 
     156        //! Logger 
     157        int L_xt; 
     158public: 
     159        void getdata ( vec &dt0 ) {dt0=dt;} 
     160        void getdata ( vec &dt0, const ivec &indeces ) {dt0=dt ( indeces );} 
     161 
     162        stateDS ( mpdf* IM0, mpdf* OM0, int usize ) :DS ( ),IM ( IM0 ),OM ( OM0 ), 
     163                        dt ( OM0->dimension() ), xt ( IM0->dimension() ), ut ( usize ) {} 
     164        ~stateDS() {delete IM; delete OM;} 
     165        virtual void step() { 
     166                xt=IM->samplecond ( concat ( xt,ut ) ); 
     167                dt=OM->samplecond ( concat ( xt,ut ) ); 
    48168        }; 
    49169 
    50         /*! 
    51         \brief Generator of ARX data 
     170        virtual void log_add ( logger &L ) { 
     171                DS::log_add ( L ); 
     172                L_xt=L.add ( IM->_rv(),"true" ); 
     173        } 
     174        virtual void logit ( logger &L ) { 
     175                DS::logit ( L ); 
     176                L.logit ( L_xt,xt ); 
     177        } 
    52178 
    53         */ 
    54         class ArxDS : public DS { 
    55         protected: 
    56                 //! Rv of the regressor 
    57                 RV Rrv; 
    58                 //! Rv of the history (full regressor) 
    59                 RV Hrv; 
    60                 //! History, ordered as \f$[y_t, u_t, y_{t-1 }, u_{t-1}, \ldots]\f$ 
    61                 vec H; 
    62                 //! (future) input 
    63                 vec U; 
    64                 //! temporary variable for regressor 
    65                 vec rgr; 
    66                 //! data link: H -> rgr 
    67                 datalink *rgrlnk; 
    68                 //! model of Y - linear Gaussian 
    69                 mlnorm<chmat> model; 
    70                 //! options 
    71                 bool opt_L_theta; 
    72                 //! loggers 
    73                 int L_theta; 
    74                 int L_R; 
    75         public: 
    76                 void getdata ( vec &dt ) { 
    77                         //it_assert_debug ( dt.length() ==Drv.count(),"ArxDS" ); 
    78                         dt=H.left ( Urv._dsize() +Drv._dsize() ); 
    79                 }; 
    80                 void getdata ( vec &dt, const ivec &indices ) { 
    81                         it_assert_debug ( dt.length() ==indices.length(),"ArxDS" ); 
    82                         dt=H ( indices ); 
    83                 }; 
    84                 void write ( vec &ut ) { 
    85                         //it_assert_debug ( ut.length() ==Urv.count(),"ArxDS" ); 
    86                         U=ut; 
    87                 }; 
    88                 void write ( vec &ut, const ivec &indices ) { 
    89                         it_assert_debug ( ut.length() ==indices.length(),"ArxDS" ); 
    90                         set_subvector ( U, indices,ut ); 
    91                 }; 
    92                 void step(); 
    93                 //!Default constructor 
    94                 ArxDS ( ){}; 
    95                 //! Set parameters of the internal model 
    96                 void set_parameters ( const mat &Th0, const vec mu0, const chmat &sqR0 ) 
    97                 { model.set_parameters ( Th0, mu0, sqR0 ); }; 
    98                 //! set options from a string 
    99                 void set_options ( const string &s ) { 
    100                         opt_L_theta= ( s.find ( "L_theta" ) !=string::npos ); 
    101                 }; 
    102                 virtual void log_add ( logger &L ) { 
    103                         DS::log_add ( L ); 
    104                         mat &A =model._A(); 
    105                         mat R =model._R(); 
    106                         if ( opt_L_theta ) {L_theta=L.add ( RV ( "{theta }", vec_1 ( A.rows() *A.cols() ) ),"t" );} 
    107                         if ( opt_L_theta ) {L_R=L.add ( RV ( "{R }", vec_1 ( R.rows() *R.cols() ) ),"r" );} 
    108                 } 
    109                 virtual void logit ( logger &L ) { 
    110                         DS::logit ( L ); 
    111                         mat &A =model._A(); 
    112                         mat R =model._R(); 
    113                         if ( opt_L_theta ) {L.logit ( L_theta,vec ( A._data(), A.rows() *A.cols() ) );}; 
    114                         if ( opt_L_theta ) {L.logit ( L_R, vec ( R._data(), R.rows() *R.rows() ) );}; 
    115                 } 
    116  
    117         }; 
    118  
    119         class ARXDS : public ArxDS { 
    120         public: 
    121                 ARXDS ( ) : ArxDS ( ) {} 
    122  
    123                 void getdata ( vec &dt ) {dt=H;} 
    124                 void getdata ( vec &dt, const ivec &indeces ) {dt=H ( indeces );} 
    125                 virtual RV _drv() const {return Hrv;} 
    126  
    127         }; 
    128  
    129         class stateDS : public DS { 
    130         protected: 
    131                 //!conditional pdf of the state evolution \f$ f(x_t|x_{t-1}) \f$ 
    132                 mpdf* IM; 
    133                 //!conditional pdf of the observations \f$ f(d_t|x_t) \f$ 
    134                 mpdf* OM; 
    135                 //! result storage 
    136                 vec dt; 
    137                 //! state storage 
    138                 vec xt; 
    139                 //! input storage 
    140                 vec ut; 
    141                 //! Logger 
    142                 int L_xt; 
    143         public: 
    144                 void getdata ( vec &dt0 ) {dt0=dt;} 
    145                 void getdata ( vec &dt0, const ivec &indeces ) {dt0=dt ( indeces );} 
    146  
    147                 stateDS ( mpdf* IM0, mpdf* OM0, int usize ) :DS ( ),IM ( IM0 ),OM ( OM0 ), 
    148                                 dt ( OM0->dimension() ), xt ( IM0->dimension() ), ut ( usize) {} 
    149                 ~stateDS() {delete IM; delete OM;} 
    150                 virtual void step() { 
    151                         xt=IM->samplecond(concat ( xt,ut )); 
    152                         dt=OM->samplecond(concat ( xt,ut )); 
    153                 }; 
    154                  
    155                 virtual void log_add ( logger &L ) { 
    156                         DS::log_add ( L ); 
    157                         L_xt=L.add(IM->_rv(),"true"); 
    158                 } 
    159                 virtual void logit ( logger &L ) { 
    160                         DS::logit ( L ); 
    161                         L.logit ( L_xt,xt); 
    162                 } 
    163  
    164         }; 
     179}; 
    165180 
    166181}; //namespace