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

A base-class for Auxiliary Particle Filtering. Filtering only, no smoothing. More...

#include <auxiliary_pf.h>

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

Public Member Functions

 APF (const unsigned int &rs=1)
 The constructor. More...
 
virtual ~APF ()
 The (virtual) destructor.
 
float_t getLogCondLike () const
 Get the latest log conditional likelihood. More...
 
std::vector< MatgetExpectations () const
 return all stored expectations (taken with respect to $p(x_t|y_{1:t})$ More...
 
void filter (const osv &data, const std::vector< std::function< const Mat(const ssv &)> > &fs=std::vector< std::function< const Mat(const ssv &)> >())
 Use a new datapoint to update the filtering distribution (or smoothing if pathLength > 0). More...
 
virtual float_t logMuEv (const ssv &x1)=0
 Evaluates the log of mu. More...
 
virtual ssv propMu (const ssv &xtm1)=0
 Evaluates the proposal distribution taking a Eigen::Matrix<float_t,dimx,1> from the previous time's state, and returning a state for the current time. More...
 
virtual ssv q1Samp (const osv &y1)=0
 Samples from q1. More...
 
virtual ssv fSamp (const ssv &xtm1)=0
 Samples from f. More...
 
virtual float_t logQ1Ev (const ssv &x1, const osv &y1)=0
 Evaluates the log of q1. More...
 
virtual float_t logGEv (const osv &yt, const ssv &xt)=0
 Evaluates the log of g. 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

std::array< ssv, nparts > m_particles
 particle samples
 
std::array< float_t, nparts > m_logUnNormWeights
 particle unnormalized weights
 
unsigned int m_now
 curren time
 
float_t m_logLastCondLike
 log p(y_t|y_{1:t-1}) or log p(y1)
 
unsigned int m_rs
 the resampling schedule
 
resamp_t m_resampler
 resampler object (default ctor'd)
 
rvsamp::k_gen< nparts, float_t > m_kGen
 k generator object (default ctor'd)
 
std::vector< Matm_expectations
 expectations E[h(x_t) | y_{1:t}] for user defined "h"s
 

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 arrayfloat_t = std::array< float_t, nparts >
 
using arrayVec = std::array< ssv, nparts >
 
using arrayUInt = std::array< unsigned int, 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::APF< nparts, dimx, dimy, resamp_t, float_t, debug >

A base-class for Auxiliary Particle Filtering. Filtering only, no smoothing.

Author
taylor

Member Typedef Documentation

◆ arrayfloat_t

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

type alias for array of float_ts

◆ arrayUInt

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug = false>
using pf::filters::APF< nparts, dimx, dimy, resamp_t, float_t, debug >::arrayUInt = std::array<unsigned int, nparts>
private

type alias for array of unsigned ints

◆ arrayVec

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

type alias for array of state vectors

◆ Mat

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

type alias for linear algebra stuff (dimension of the state ^2)

◆ osv

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

"observation 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::APF< 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

◆ APF()

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

The constructor.

Parameters
rsresampling schedule (e.g. resample every rs time points).

Member Function Documentation

◆ filter()

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

Use a new datapoint to update the filtering distribution (or smoothing if pathLength > 0).

Parameters
dataa Eigen::Matrix<float_t,dimy,1> representing the data
fsa std::vector of callback functions that are used to calculate expectations with respect to the filtering distribution.

◆ fSamp()

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

Samples from f.

Parameters
xtm1a Eigen::Matrix<float_t,dimx,1> representing the previous time's state.
Returns
a Eigen::Matrix<float_t,dimx,1> state sample for the current time.

◆ getExpectations()

template<size_t nparts, size_t dimx, size_t dimy, typename resamp_t , typename float_t , bool debug>
auto pf::filters::APF< nparts, dimx, dimy, resamp_t, float_t, debug >::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 pf::filters::APF< nparts, dimx, dimy, resamp_t, float_t, debug >::getLogCondLike

Get the latest log conditional likelihood.

Returns
a float_t of the most recent conditional likelihood.

◆ 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::APF< nparts, dimx, dimy, resamp_t, float_t, debug >::logGEv ( const osv yt,
const ssv xt 
)
pure virtual

Evaluates the log of g.

Parameters
yta Eigen::Matrix<float_t,dimy,1> representing time t's data observation.
xta Eigen::Matrix<float_t,dimx,1> representing time t's state.
Returns
a float_t 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::APF< nparts, dimx, dimy, resamp_t, float_t, debug >::logMuEv ( const ssv x1)
pure virtual

Evaluates the log of mu.

Parameters
x1a Eigen::Matrix<float_t,dimx,1> representing time 1's state.
Returns
a float_t 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::APF< nparts, dimx, dimy, resamp_t, float_t, debug >::logQ1Ev ( const ssv x1,
const osv y1 
)
pure virtual

Evaluates the log of q1.

Parameters
x1a Eigen::Matrix<float_t,dimx,1> representing time 1's state.
y1a Eigen::Matrix<float_t,dimy,1> representing time 1's data observation.
Returns
a float_t evaluation.

◆ propMu()

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

Evaluates the proposal distribution taking a Eigen::Matrix<float_t,dimx,1> from the previous time's state, and returning a state for the current time.

Parameters
xtm1a Eigen::Matrix<float_t,dimx,1> representing the previous time's state.
Returns
a Eigen::Matrix<float_t,dimx,1> representing a likely current time state, to be used by the observation density.

◆ q1Samp()

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

Samples from q1.

Parameters
y1a Eigen::Matrix<float_t,dimy,1> representing time 1's data point.
Returns
a Eigen::Matrix<float_t,dimx,1> sample for time 1's state.

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