Changeset 60

Show
Ignore:
Timestamp:
03/31/08 13:53:36 (17 years ago)
Author:
smidl
Message:

nove rozlozeni Q

Files:
1 added
7 modified

Legend:

Unmodified
Added
Removed
  • TabularUnified BDM.kdevelop

    r45 r60  
    2424    <filelistdirectory>/home/smidl/work/mixpp</filelistdirectory> 
    2525    <run> 
    26       <mainprogram>/home/smidl/work/mixpp/tests/pmsm_sim</mainprogram> 
     26      <mainprogram>/home/smidl/work/mixpp/tests/pmsm_sim2</mainprogram> 
    2727      <directoryradio>custom</directoryradio> 
    2828      <customdirectory>/home/smidl/work/mixpp</customdirectory> 
     
    171171    <tree> 
    172172      <hidepatterns>*.o,*.lo,CVS,*~,cmake*,Makefile</hidepatterns> 
    173       <hidenonprojectfiles>true</hidenonprojectfiles> 
     173      <hidenonprojectfiles>false</hidenonprojectfiles> 
    174174      <showvcsfields>false</showvcsfields> 
    175175    </tree> 
    
          
  • TabularUnified BDM.kdevelop.filelist

    r45 r60  
    2828tests/CMakeLists.txt 
    2929tests/pmsm_sim.cpp 
     30tests/pmsm_sim2.cpp 
    3031tests/pmsm_unkQ.cpp 
    3132tests/pmsm_unkQpf.cpp 
     
    3839tests/testResample.cpp 
    3940tests/testSmp.cpp 
     41tests/testbidiff.cpp 
    
          
  • TabularUnified bdm/estim/libPF.cpp

    r32 r60  
    44 
    55using std::endl; 
    6  
    76 
    87void PF::bayes ( const vec &dt ) { 
    
          
  • TabularUnified bdm/estim/libPF.h

    r33 r60  
    1919 
    2020using namespace itpp; 
     21 
     22//UGLY HACK 
     23extern double PF_SSAT;//used for StrSim:06 test... if length>0 the value is written. 
    2124 
    2225/*! 
     
    6063 
    6164class MPF : public PF { 
    62         BM_T* Bms[1000]; 
     65        BM_T* Bms[10000]; 
    6366 
    6467        //! internal class for MPDF providing composition of eEmp with external components 
     
    105108                // 
    106109 
    107                 if ( n>1000 ) it_error ( "increase 1000 here!" ); 
     110                if ( n>10000 ) it_error ( "increase 10000 here!" ); 
    108111 
    109112                for ( int i=0;i<n;i++ ) { 
     
    126129                for ( int i=0;i<n;i++ ) {Bms[i]->condition ( _samples ( i ) );} 
    127130        } 
     131 
     132//SimStr: 
     133        double SSAT;     
    128134}; 
    129135 
     
    132138        int i; 
    133139        vec lls ( n ); 
     140        vec llsP ( n ); 
    134141        ivec ind; 
    135142        double mlls=-std::numeric_limits<double>::infinity(); 
    136143 
     144        // StrSim:06 
     145        double sumLWL=0.0; 
     146        double sumL2WL=0.0; 
     147        double WL = 0.0; 
     148 
    137149        for ( i=0;i<n;i++ ) { 
    138150                //generate new samples from paramater evolution model; 
    139                 _samples ( i ) = par.samplecond ( _samples ( i ), lls ( i ) ); 
     151                _samples ( i ) = par.samplecond ( _samples ( i ), llsP ( i ) ); 
    140152                Bms[i]->condition ( _samples ( i ) ); 
    141153                Bms[i]->bayes ( dt ); 
     
    144156                if ( lls ( i ) >mlls ) mlls=lls ( i ); //find maximum likelihood (for numerical stability) 
    145157        } 
     158         
     159        if ( true ) { 
     160                for ( i=0;i<n;i++ ) { 
     161                        WL = _w ( i ) *exp ( llsP ( i ) ); //using old weights! 
     162                        sumLWL += exp ( lls ( i ) ) *WL; 
     163                        sumL2WL += exp ( 2*lls ( i ) ) *WL; 
     164                } 
     165                SSAT  = sumL2WL/(sumLWL*sumLWL); 
     166        } 
    146167 
    147168        // compute weights 
     
    150171        } 
    151172 
    152         _w /=sum(_w); //? 
    153  
    154         double eff = 1.0/( _w*_w ); 
     173        _w /=sum ( _w ); //? 
     174 
     175 
     176 
     177        double eff = 1.0/ ( _w*_w ); 
    155178        if ( eff < ( 0.1*n ) ) { 
    156179                ind = est.resample(); 
     
    161184                                //fixme this would require new assignment operator 
    162185                                // *Bms[i] = *Bms[ind ( i ) ]; 
    163                                  
     186 
    164187                                // poor-man's solution: replicate constructor here 
    165188                                // copied from MPF::MPF 
     
    170193                        } 
    171194                }; 
    172                 cout << '.';  
     195                cout << '.'; 
    173196        } 
    174197} 
    
          
  • TabularUnified bdm/stat/libEF.cpp

    r32 r60  
    6060 
    6161vec mgamma::samplecond ( vec &cond, double &ll ) { 
    62         *_beta=k / cond ; 
     62        this->condition(cond ); 
    6363        vec smp = epdf.sample(); 
    6464        ll = epdf.evalpdflog ( smp ); 
    
          
  • TabularUnified bdm/stat/libEF.h

    r51 r60  
    9696        double eval ( const vec &val ) const ; 
    9797        double evalpdflog ( const vec &val ) const; 
    98         vec mean()const {return mu;} 
     98        vec mean() const {return mu;} 
    9999 
    100100//Access methods 
     
    140140        //! Returns poiter to alpha and beta. Potentially dengerous: use with care! 
    141141        void _param ( vec* &a, vec* &b ) {a=&alpha;b=&beta;}; 
    142         vec mean()const {vec pom(alpha); pom/=beta; return pom;} 
     142        vec mean() const {vec pom ( alpha ); pom/=beta; return pom;} 
    143143}; 
    144144/* 
     
    228228*/ 
    229229class mgamma : public mEF { 
     230protected: 
    230231        //! Internal epdf that arise by conditioning on \c rvc 
    231232        egamma epdf; 
     
    247248}; 
    248249 
     250/*! 
     251\brief  Gamma random walk around a fixed point 
     252 
     253Mean value, \f$\mu\f$, of this density is given by a geometric combination of \c rvc and given fixed point, $p$. $k$ is the coefficient of the geometric combimation 
     254\f[ \mu = \mu_{t-1} ^{l} p^{1-l}\f] 
     255 
     256Standard deviation of the random walk is proportional to one $k$-th the mean. 
     257This is achieved by setting \f$\alpha=k\f$ and \f$\beta=k/\mu\f$. 
     258 
     259The standard deviation of the walk is then: \f$\mu/\sqrt(k)\f$. 
     260*/ 
     261class mgamma_fix : public mgamma { 
     262protected: 
     263        double l; 
     264        vec refl; 
     265public: 
     266        //! Constructor 
     267        mgamma_fix ( const RV &rv,const RV &rvc ) : mgamma ( rv,rvc ),refl ( rv.count() ) {}; 
     268        //! Set value of \c k 
     269        void set_parameters ( double k0 , vec ref0, double l0 ) { 
     270                mgamma::set_parameters ( k0 ); 
     271                refl=pow ( ref0,1.0-l0 );l=l0; 
     272        }; 
     273 
     274        void condition ( const vec &val ) {vec mean=elem_mult ( refl,pow ( val,l ) ); *_beta=k/mean;}; 
     275}; 
     276 
    249277//! Switch between various resampling methods. 
    250278enum RESAMPLING_METHOD { MULTINOMIAL = 0, STRATIFIED = 1, SYSTEMATIC = 3 }; 
     
    264292public: 
    265293        //! Default constructor 
    266         eEmp ( const RV &rv0 ,int n0) :epdf ( rv0 ),n(n0),w(n),samples(n) {}; 
     294        eEmp ( const RV &rv0 ,int n0 ) :epdf ( rv0 ),n ( n0 ),w ( n ),samples ( n ) {}; 
    267295        //! Set sample 
    268296        void set_parameters ( const vec &w0, epdf* pdf0 ); 
     
    276304        vec sample() const {it_error ( "Not implemented" );return 0;} 
    277305        //! inherited operation : NOT implemneted 
    278         double evalpdflog(const vec &val) const {it_error ( "Not implemented" );return 0.0;} 
    279         vec mean()const {vec pom=zeros(rv.count());  
    280                 for (int i=0;i<n;i++){pom+=samples(i)*w(i);} 
     306        double evalpdflog ( const vec &val ) const {it_error ( "Not implemented" );return 0.0;} 
     307        vec mean() const { 
     308                vec pom=zeros ( rv.count() ); 
     309                for ( int i=0;i<n;i++ ) {pom+=samples ( i ) *w ( i );} 
    281310                return pom; 
    282311        } 
     
    287316 
    288317template<class sq_T> 
    289 enorm<sq_T>::enorm ( RV &rv ) :eEF(rv), mu ( rv.count() ),R ( rv.count() ),_iR ( rv.count() ),cached ( false ),dim ( rv.count() ) {}; 
     318enorm<sq_T>::enorm ( RV &rv ) :eEF ( rv ), mu ( rv.count() ),R ( rv.count() ),_iR ( rv.count() ),cached ( false ),dim ( rv.count() ) {}; 
    290319 
    291320template<class sq_T> 
     
    320349 
    321350template<class sq_T> 
    322 mat enorm<sq_T>::sample ( int N )const { 
     351mat enorm<sq_T>::sample ( int N ) const { 
    323352        mat X ( dim,N ); 
    324353        vec x ( dim ); 
     
    346375template<class sq_T> 
    347376double enorm<sq_T>::evalpdflog ( const vec &val ) const { 
    348         if ( !cached ) {it_error("this should not happen, see cached");} 
     377        if ( !cached ) {it_error ( "this should not happen, see cached" );} 
    349378 
    350379        // 1.83787706640935 = log(2pi) 
    
          
  • TabularUnified tests/CMakeLists.txt

    r54 r60  
    4646link_directories (${BDM_BINARY_DIR}/simulator_zdenek) 
    4747add_executable (pmsm_sim pmsm_sim.cpp) 
     48add_executable (pmsm_sim2 pmsm_sim2.cpp) 
    4849target_link_libraries (pmsm_sim ${BdmLibs} pmsmsim netcdf_c++) 
     50target_link_libraries (pmsm_sim2 ${BdmLibs} pmsmsim netcdf_c++)