pf
pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug > Class Template Reference

A class template for Kalman filtering. More...

#include <cf_filters.h>

Inheritance diagram for pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >:
Collaboration diagram for pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >:

Public Types

using ssv = Eigen::Matrix< float_t, dimstate, 1 >
 
using osv = Eigen::Matrix< float_t, dimobs, 1 >
 
using isv = Eigen::Matrix< float_t, diminput, 1 >
 
using ssMat = Eigen::Matrix< float_t, dimstate, dimstate >
 
using osMat = Eigen::Matrix< float_t, dimobs, dimobs >
 
using siMat = Eigen::Matrix< float_t, dimstate, diminput >
 
using oiMat = Eigen::Matrix< float_t, dimobs, diminput >
 
using obsStateSizeMat = Eigen::Matrix< float_t, dimobs, dimstate >
 
using stateObsSizeMat = Eigen::Matrix< float_t, dimstate, dimobs >
 
- Public Types inherited from pf::bases::cf_filter< dimstate, dimobs, float_t >
using obs_sized_vec = Eigen::Matrix< float_t, dimobs, 1 >
 
using state_sized_vec = Eigen::Matrix< float_t, dimstate, 1 >
 

Public Member Functions

 kalman ()
 Default constructor. More...
 
 kalman (const ssv &initStateMean, const ssMat &initStateVar)
 Non-default constructor. More...
 
virtual ~kalman ()
 The (virtual) destructor.
 
float_t getLogCondLike () const
 returns the log of the latest conditional likelihood. More...
 
ssv getFiltMean () const
 Get the current filter mean. More...
 
ssMat getFiltVar () const
 Get the current filter variance-covariance matrix. More...
 
osv getPredYMean (const ssMat &stateTrans, const obsStateSizeMat &obsMat, const siMat &stateInptAffector, const oiMat &obsInptAffector, const isv &inputData) const
 get the one-step-ahead point forecast for y More...
 
osMat getPredYVar (const ssMat &stateTrans, const ssMat &cholStateVar, const obsStateSizeMat &obsMat, const osMat &cholObsVar) const
 get the one-step-ahead forecast variance More...
 
void update (const osv &yt, const ssMat &stateTrans, const ssMat &cholStateVar, const siMat &stateInptAffector, const isv &inputData, const obsStateSizeMat &obsMat, const oiMat &obsInptAffector, const osMat &cholObsVar)
 Perform a Kalman filter predict-and-update. More...
 
- Public Member Functions inherited from pf::bases::cf_filter< dimstate, dimobs, float_t >
virtual ~cf_filter ()
 The (virtual) destructor.
 

Private Member Functions

void updatePrior (const ssMat &stateTransMat, const ssMat &cholStateVar, const siMat &stateInptAffector, const isv &inputData)
 Predicts the next state. More...
 
void updatePosterior (const osv &yt, const obsStateSizeMat &obsMat, const oiMat &obsInptAffector, const isv &inputData, const osMat &cholObsVar)
 Turns prediction into new filtering distribution. More...
 

Private Attributes

ssv m_predMean
 predictive state mean
 
ssv m_filtMean
 filter mean
 
ssMat m_predVar
 predictive var matrix
 
ssMat m_filtVar
 filter var matrix
 
float_t m_lastLogCondLike
 latest log conditional likelihood
 
bool m_fresh
 has data been observed?
 
const float_t m_pi
 pi
 

Detailed Description

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t, bool debug = false>
class pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >

A class template for Kalman filtering.

Author
taylor

Member Typedef Documentation

◆ isv

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug = false>
using pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::isv = Eigen::Matrix<float_t,diminput,1>

"input size vector" type alias for linear algebra stuff

◆ obsStateSizeMat

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug = false>
using pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::obsStateSizeMat = Eigen::Matrix<float_t,dimobs,dimstate>

"observation dimension by state dimension -sized matrix"

◆ oiMat

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug = false>
using pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::oiMat = Eigen::Matrix<float_t,dimobs,diminput>

"observation dimension by input dim matrix"

◆ osMat

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug = false>
using pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::osMat = Eigen::Matrix<float_t,dimobs,dimobs>

"observation size matrix" type alias for linear algebra stuff

◆ osv

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug = false>
using pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::osv = Eigen::Matrix<float_t,dimobs,1>

"observation size vector" type alias for linear algebra stuff

◆ siMat

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug = false>
using pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::siMat = Eigen::Matrix<float_t,dimstate,diminput>

"state dim by input dimension matrix"

◆ ssMat

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug = false>
using pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::ssMat = Eigen::Matrix<float_t,dimstate,dimstate>

"state size matrix" type alias for linear algebra stuff

◆ ssv

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug = false>
using pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::ssv = Eigen::Matrix<float_t,dimstate,1>

"state size vector" type alias for linear algebra stuff

◆ stateObsSizeMat

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug = false>
using pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::stateObsSizeMat = Eigen::Matrix<float_t,dimstate,dimobs>

"state dimension by observation dimension matrix

Constructor & Destructor Documentation

◆ kalman() [1/2]

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug>
pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::kalman

Default constructor.

Need ths fir constructing default std::array<>s. Fills all vectors and matrices with zeros.

◆ kalman() [2/2]

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug>
pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::kalman ( const ssv initStateMean,
const ssMat initStateVar 
)

Non-default constructor.

Non-default constructor.

Member Function Documentation

◆ getFiltMean()

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug>
auto pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::getFiltMean

Get the current filter mean.

Returns
E[x_t | y_{1:t}]

◆ getFiltVar()

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug>
auto pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::getFiltVar

Get the current filter variance-covariance matrix.

Returns
V[x_t | y_{1:t}]

◆ getLogCondLike()

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug>
float_t pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::getLogCondLike
virtual

returns the log of the latest conditional likelihood.

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

Implements pf::bases::cf_filter< dimstate, dimobs, float_t >.

◆ getPredYMean()

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug>
auto pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::getPredYMean ( const ssMat stateTrans,
const obsStateSizeMat obsMat,
const siMat stateInptAffector,
const oiMat obsInptAffector,
const isv inputData 
) const

get the one-step-ahead point forecast for y

Returns
E[y_{t+1} | y_{1:t}, params]

◆ getPredYVar()

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug>
auto pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::getPredYVar ( const ssMat stateTrans,
const ssMat cholStateVar,
const obsStateSizeMat obsMat,
const osMat cholObsVar 
) const

get the one-step-ahead forecast variance

Returns
V[y_{t+1} | y_{1:t}, params]

◆ update()

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug>
void pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::update ( const osv yt,
const ssMat stateTrans,
const ssMat cholStateVar,
const siMat stateInptAffector,
const isv inputData,
const obsStateSizeMat obsMat,
const oiMat obsInptAffector,
const osMat cholObsVar 
)

Perform a Kalman filter predict-and-update.

Parameters
ytthe new data point.
stateTransthe transition matrix of the state
cholStateVarthe Cholesky Decomposition of the state noise covariance matrix.
stateInptAffectorthe matrix affecting how input data affects state transition.
inputDataexogenous input data
obsMatthe observation/emission matrix of the observation's conditional (on the state) distn.
obsInptAffectorthe matrix affecting how input data affects the observational distribution.
cholObsVarthe Cholesky Decomposition of the observatio noise covariance matrix.

◆ updatePosterior()

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug>
void pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::updatePosterior ( const osv yt,
const obsStateSizeMat obsMat,
const oiMat obsInptAffector,
const isv inputData,
const osMat cholObsVar 
)
private

Turns prediction into new filtering distribution.

Parameters
yt
obsMat
obsInptAffector
inputData
cholObsVar

◆ updatePrior()

template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t , bool debug>
void pf::filters::kalman< dimstate, dimobs, diminput, float_t, debug >::updatePrior ( const ssMat stateTransMat,
const ssMat cholStateVar,
const siMat stateInptAffector,
const isv inputData 
)
private

Predicts the next state.

Todo:
handle diagonal variance matrices, and ensure symmetricness in other ways
Parameters
stateTransMat
cholStateVar
stateInptAffector
inputData

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