Rao-Blackwellized/Marginal Particle Filter with inner Kalman Filter objectss.
More...
#include <rbpf.h>
|
using | sssv = Eigen::Matrix< float_t, dimss, 1 > |
|
using | nsssv = Eigen::Matrix< float_t, dimnss, 1 > |
|
using | osv = Eigen::Matrix< float_t, dimy, 1 > |
|
using | Mat = Eigen::Matrix< float_t, Eigen::Dynamic, Eigen::Dynamic > |
|
using | nsssMat = Eigen::Matrix< float_t, dimnss, dimnss > |
|
using | arrayMod = std::array< kalman< dimnss, dimy, 0, float_t >, nparts > |
|
using | arrayVec = std::array< sssv, nparts > |
|
using | arrayfloat_t = std::array< float_t, nparts > |
|
using | osv = Eigen::Matrix< float_t, dimobs, 1 > |
|
using | sssv = Eigen::Matrix< float_t, dim_s_state, 1 > |
|
using | nsssv = Eigen::Matrix< float_t, dim_ns_state, 1 > |
|
using | Mat = Eigen::Matrix< float_t, Eigen::Dynamic, Eigen::Dynamic > |
|
using | func = std::function< const Mat(const nsssv &, const sssv &)> |
|
using | funcs = std::vector< func > |
|
|
| rbpf_kalman (const unsigned int &resamp_sched=1) |
| The constructor. More...
|
|
void | filter (const osv &data, const std::vector< std::function< const Mat(const nsssv &x1t, const sssv &x2t)> > &fs=std::vector< std::function< const Mat(const nsssv &x1t, const sssv &x2t)> >()) |
| Filter! More...
|
|
float_t | getLogCondLike () const |
| Get the latest log conditional likelihood. More...
|
|
std::vector< Mat > | getExpectations () const |
| Get the latest filtered expectation E[h(x_1t, x_2t) | y_{1:t}]. More...
|
|
virtual float_t | logMuEv (const sssv &x21)=0 |
| Evaluates the first time state density. More...
|
|
virtual sssv | q1Samp (const osv &y1)=0 |
| Sample from the first time's proposal distribution. More...
|
|
virtual nsssv | initKalmanMean (const sssv &x21)=0 |
| Provides the initial mean vector for each Kalman filter object. More...
|
|
virtual nsssMat | initKalmanVar (const sssv &x21)=0 |
| Provides the initial covariance matrix for each Kalman filter object. More...
|
|
virtual sssv | qSamp (const sssv &x2tm1, const osv &yt)=0 |
| Samples the time t second component. More...
|
|
virtual float_t | logQ1Ev (const sssv &x21, const osv &y1)=0 |
| Evaluates the proposal density of the second state component at time 1. More...
|
|
virtual float_t | logFEv (const sssv &x2t, const sssv &x2tm1)=0 |
| Evaluates the state transition density for the second state component. More...
|
|
virtual float_t | logQEv (const sssv &x2t, const sssv &x2tm1, const osv &yt)=0 |
| Evaluates the proposal density at time t > 1. More...
|
|
virtual void | updateKalman (kalman< dimnss, dimy, 0, float_t > &kMod, const osv &yt, const sssv &x2t)=0 |
| How to update your inner Kalman filter object at each time. More...
|
|
virtual void | filter (const osv &data, const funcs &fs=funcs())=0 |
| the filtering function that must be defined
|
|
|
static constexpr unsigned int | dim_sampled_state |
|
static constexpr unsigned int | dim_not_sampled_state |
|
static constexpr unsigned int | dim_obs |
|
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t, typename float_t>
class rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >
Rao-Blackwellized/Marginal Particle Filter with inner Kalman Filter objectss.
- Author
- t
◆ arrayfloat_t
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
◆ arrayMod
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::arrayMod = std::array<kalman<dimnss,dimy,0,float_t>,nparts> |
◆ arrayVec
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
◆ Mat
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::Mat = Eigen::Matrix<float_t,Eigen::Dynamic,Eigen::Dynamic> |
◆ nsssMat
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::nsssMat = Eigen::Matrix<float_t,dimnss,dimnss> |
"not sampled state size matrix"
◆ nsssv
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::nsssv = Eigen::Matrix<float_t,dimnss,1> |
"not sampled state size vector"
◆ osv
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::osv = Eigen::Matrix<float_t,dimy,1> |
"observation size vector"
◆ sssv
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::sssv = Eigen::Matrix<float_t,dimss,1> |
"sampled state size vector"
◆ rbpf_kalman()
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::rbpf_kalman |
( |
const unsigned int & |
resamp_sched = 1 | ) |
|
The constructor.
- Parameters
-
resamp_sched | how often you want to resample (e.g once every resamp_sched time points) |
◆ filter()
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
void rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::filter |
( |
const osv & |
data, |
|
|
const std::vector< std::function< const Mat(const nsssv &x1t, const sssv &x2t)> > & |
fs = std::vector<std::function<const Mat(const nsssv &x1t, const sssv &x2t)> >() |
|
) |
| |
Filter!
The workhorse function
- Parameters
-
data | the most recent observable portion of the time series. |
fs | a vector of functions computing E[h(x_1t, x_2t^i)| x_2t^i,y_1:t]. to be averaged to yield E[h(x_1t, x_2t)|,y_1:t] |
◆ getExpectations()
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
auto rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::getExpectations |
( |
| ) |
const |
Get the latest filtered expectation E[h(x_1t, x_2t) | y_{1:t}].
Get the expectations you're keeping track of.
- Returns
- a vector of Mats
◆ getLogCondLike()
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
float_t rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::getLogCondLike |
( |
| ) |
const |
Get the latest log conditional likelihood.
- Returns
- the latest log conditional likelihood.
◆ initKalmanMean()
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual nsssv rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::initKalmanMean |
( |
const sssv & |
x21 | ) |
|
|
pure virtual |
Provides the initial mean vector for each Kalman filter object.
provides the initial mean vector for each Kalman filter object.
- Parameters
-
x21 | the second state componenent at time 1. |
- Returns
- a nsssv representing the unconditional mean.
◆ initKalmanVar()
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual nsssMat rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::initKalmanVar |
( |
const sssv & |
x21 | ) |
|
|
pure virtual |
Provides the initial covariance matrix for each Kalman filter object.
provides the initial covariance matrix for each Kalman filter object.
- Parameters
-
x21 | the second state component at time 1. |
- Returns
- a covariance matrix.
◆ logFEv()
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual float_t rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::logFEv |
( |
const sssv & |
x2t, |
|
|
const sssv & |
x2tm1 |
|
) |
| |
|
pure virtual |
Evaluates the state transition density for the second state component.
Evaluates the state transition density for the second state component.
- Parameters
-
x2t | the current second state component. |
x2tm1 | the previous second state component. |
- Returns
- a float_t evaluation.
◆ logMuEv()
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual float_t rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::logMuEv |
( |
const sssv & |
x21 | ) |
|
|
pure virtual |
Evaluates the first time state density.
evaluates log mu(x21).
- Parameters
-
x21 | component two at time 1 |
- Returns
- a float_t evaluation
◆ logQ1Ev()
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual float_t rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::logQ1Ev |
( |
const sssv & |
x21, |
|
|
const osv & |
y1 |
|
) |
| |
|
pure virtual |
Evaluates the proposal density of the second state component at time 1.
Evaluates the proposal density of the second state component at time 1.
- Parameters
-
x21 | the second state component at time 1 you sampled. |
y1 | time 1 observation. |
- Returns
- a float_t evaluation of the density.
◆ logQEv()
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual float_t rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::logQEv |
( |
const sssv & |
x2t, |
|
|
const sssv & |
x2tm1, |
|
|
const osv & |
yt |
|
) |
| |
|
pure virtual |
Evaluates the proposal density at time t > 1.
Evaluates the proposal density at time t > 1.
- Parameters
-
x2t | the current second state component. |
x2tm1 | the previous second state component. |
yt | the current time series observation. |
- Returns
- a float_t evaluation.
◆ q1Samp()
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual sssv rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::q1Samp |
( |
const osv & |
y1 | ) |
|
|
pure virtual |
Sample from the first time's proposal distribution.
samples the second component of the state at time 1.
- Parameters
-
- Returns
- a Vec sample for x21.
◆ qSamp()
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual sssv rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::qSamp |
( |
const sssv & |
x2tm1, |
|
|
const osv & |
yt |
|
) |
| |
|
pure virtual |
Samples the time t second component.
Samples the time t second component.
- Parameters
-
x2tm1 | the previous time's second state component. |
yt | the current observation. |
- Returns
- a sssv sample of the second state component at the current time.
◆ updateKalman()
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual void rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::updateKalman |
( |
kalman< dimnss, dimy, 0, float_t > & |
kMod, |
|
|
const osv & |
yt, |
|
|
const sssv & |
x2t |
|
) |
| |
|
pure virtual |
How to update your inner Kalman filter object at each time.
How to update your inner Kalman filter object at each time.
- Parameters
-
kMod | a Kalman filter object describing the conditional closed-form model. |
yt | the current time series observation. |
x2t | the current second state component. |
◆ m_expectations
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
std::vector<Mat> rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::m_expectations |
|
private |
◆ m_lastLogCondLike
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
float_t rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::m_lastLogCondLike |
|
private |
log p(y_t|y_{1:t-1}) or log p(y1)
◆ m_logUnNormWeights
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
the array of the (log of) unnormalized weights
◆ m_now
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
unsigned int rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::m_now |
|
private |
◆ m_p_innerMods
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
the array of inner Kalman filter objects
◆ m_p_samps
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
the array of particle samples
◆ m_resampler
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
resamp_t rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::m_resampler |
|
private |
◆ m_rs
template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
unsigned int rbpf_kalman< nparts, dimnss, dimss, dimy, resamp_t, float_t >::m_rs |
|
private |
The documentation for this class was generated from the following file: