pf
pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug > Class Template Referenceabstract

A base class for the bootstrap particle filter. More...

#include <bootstrap_filter.h>

Inheritance diagram for pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug >:
Collaboration diagram for pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug >:

Public Member Functions

 BSFilter (const unsigned int &rs=1)
 The constructor. More...
 
virtual ~BSFilter ()
 The (virtual) destructor.
 
float_t getLogCondLike () const
 Returns the most recent (log-) conditional likelihood. More...
 
void filter (const osv &data, const std::vector< std::function< const Mat(const ssv &)> > &fs=std::vector< std::function< const Mat(const ssv &)> >())
 updates filtering distribution on a new datapoint. Optionally stores expectations of functionals. More...
 
auto getExpectations () const -> std::vector< Mat >
 return all stored expectations (taken with respect to $p(x_t|y_{1:t})$ More...
 
virtual float_t logMuEv (const ssv &x1)=0
 Calculate muEv or logmuEv. More...
 
virtual ssv q1Samp (const osv &y1)=0
 Samples from time 1 proposal. More...
 
virtual float_t logQ1Ev (const ssv &x1, const osv &y1)=0
 Calculate q1Ev or log q1Ev. More...
 
virtual float_t logGEv (const osv &yt, const ssv &xt)=0
 Calculate gEv or logGEv. More...
 
virtual ssv fSamp (const ssv &xtm1)=0
 Sample from the state transition distribution. More...
 
- Public Member Functions inherited from pf::bases::pf_base< float_t, dimy, dimx >
virtual void filter (const obs_sized_vec &data, const func_vec &fs=func_vec())=0
 the filtering function that must be defined More...
 
virtual float_t getLogCondLike () const=0
 the getter method that must be defined (for conditional log-likelihood) More...
 
virtual ~pf_base ()
 virtual destructor
 

Protected Attributes

arrayStates m_particles
 particle samples
 
arrayFloat m_logUnNormWeights
 particle unnormalized weights
 
unsigned int m_now
 time point
 
float_t m_logLastCondLike
 log p(y_t|y_{1:t-1}) or log p(y1)

 
resamp_t m_resampler
 resampler object
 
std::vector< Matm_expectations
 expectations E[h(x_t) | y_{1:t}] for user defined "h"s
 
unsigned int m_resampSched
 resampling schedule (e.g. resample every __ time points)
 

Private Types

using ssv = Eigen::Matrix< float_t, dimx, 1 >
 
using osv = Eigen::Matrix< float_t, dimy, 1 >
 
using Mat = Eigen::Matrix< float_t, Eigen::Dynamic, Eigen::Dynamic >
 
using arrayStates = std::array< ssv, nparts >
 
using arrayFloat = std::array< float_t, nparts >
 

Additional Inherited Members

- Public Types inherited from pf::bases::pf_base< float_t, dimy, dimx >
using float_type = float_t
 
using obs_sized_vec = Eigen::Matrix< float_t, dimobs, 1 >
 
using state_sized_vec = Eigen::Matrix< float_t, dimstate, 1 >
 
using dynamic_matrix = Eigen::Matrix< float_t, Eigen::Dynamic, Eigen::Dynamic >
 
using func = std::function< const dynamic_matrix(const state_sized_vec &)>
 
using func_vec = std::vector< func >
 
- Static Public Attributes inherited from pf::bases::pf_base< float_t, dimy, dimx >
static constexpr unsigned int dim_obs
 
static constexpr unsigned int dim_state
 

Detailed Description

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t, typename float_t, bool debug = false>
class pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug >

A base class for the bootstrap particle filter.

Author
taylor

Member Typedef Documentation

◆ arrayFloat

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug = false>
using pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug >::arrayFloat = std::array<float_t, nparts>
private

type alias for array of floating points

◆ arrayStates

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug = false>
using pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug >::arrayStates = std::array<ssv, nparts>
private

type alias for linear algebra stuff

◆ Mat

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug = false>
using pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug >::Mat = Eigen::Matrix<float_t, Eigen::Dynamic, Eigen::Dynamic>
private

type alias for dynamically sized matrix

◆ osv

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug = false>
using pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug >::osv = Eigen::Matrix<float_t, dimy, 1>
private

"obs size vector" type alias for linear algebra stuff

◆ ssv

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug = false>
using pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug >::ssv = Eigen::Matrix<float_t, dimx, 1>
private

"state size vector" type alias for linear algebra stuff

Constructor & Destructor Documentation

◆ BSFilter()

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug>
BSFilter::BSFilter ( const unsigned int &  rs = 1)

The constructor.

Parameters
rsthe resampling schedule (e.g. every rs time point)

Member Function Documentation

◆ filter()

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug>
void BSFilter::filter ( const osv data,
const std::vector< std::function< const Mat(const ssv &)> > &  fs = std::vector<std::function<const Mat(const ssv&)> >() 
)

updates filtering distribution on a new datapoint. Optionally stores expectations of functionals.

Parameters
datathe most recent data point
fsa vector of functions if you want to calculate expectations.

◆ fSamp()

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug = false>
virtual ssv pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug >::fSamp ( const ssv xtm1)
pure virtual

Sample from the state transition distribution.

Parameters
xtm1is a const Vec& describing the time t-1 state
Returns
the sample as a Vec

◆ getExpectations()

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug>
auto BSFilter::getExpectations

return all stored expectations (taken with respect to $p(x_t|y_{1:t})$

Returns
return a std::vector<Mat> of expectations. How many depends on how many callbacks you gave to

◆ getLogCondLike()

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug>
float_t BSFilter::getLogCondLike

Returns the most recent (log-) conditional likelihood.

Returns
log p(y_t | y_{1:t-1})

◆ logGEv()

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug = false>
virtual float_t pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug >::logGEv ( const osv yt,
const ssv xt 
)
pure virtual

Calculate gEv or logGEv.

Parameters
ytis a const Vec& describing the time t datum
xtis a const Vec& describing the time t state
Returns
the density or log-density evaluation

◆ logMuEv()

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug = false>
virtual float_t pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug >::logMuEv ( const ssv x1)
pure virtual

Calculate muEv or logmuEv.

Parameters
x1is a const Vec& describing the state sample
Returns
the density or log-density evaluation

◆ logQ1Ev()

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug = false>
virtual float_t pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug >::logQ1Ev ( const ssv x1,
const osv y1 
)
pure virtual

Calculate q1Ev or log q1Ev.

Parameters
x1is a const Vec& describing the time 1 state sample
y1is a const Vec& describing the time 1 datum
Returns
the density or log-density evaluation

◆ q1Samp()

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug = false>
virtual ssv pf::filters::BSFilter< nparts, dimx, dimy, resamp_t, float_t, debug >::q1Samp ( const osv y1)
pure virtual

Samples from time 1 proposal.

Parameters
y1is a const Vec& representing the first observed datum
Returns
the sample as a Vec

The documentation for this class was generated from the following file: