1 #ifndef BOOTSTRAP_FILTER_H 2 #define BOOTSTRAP_FILTER_H 22 template<
size_t nparts,
size_t dimx,
size_t dimy,
typename resamp_t,
typename float_t,
bool debug=false>
28 using ssv = Eigen::Matrix<float_t, dimx, 1>;
30 using osv = Eigen::Matrix<float_t, dimy, 1>;
32 using Mat = Eigen::Matrix<float_t, Eigen::Dynamic, Eigen::Dynamic>;
45 BSFilter(
const unsigned int &rs = 1);
67 void filter(
const osv &data,
const std::vector<std::function<
const Mat(
const ssv&)> >& fs = std::vector<std::function<
const Mat(
const ssv&)> >());
108 virtual float_t
logGEv (const
osv &yt, const
ssv &xt ) = 0;
143 template<
size_t nparts,
size_t dimx,
size_t dimy, typename resamp_t, typename float_t,
bool debug>
146 , m_logLastCondLike(0.0)
150 std::fill(m_logUnNormWeights.begin(), m_logUnNormWeights.end(), 0.0);
154 template<
size_t nparts,
size_t dimx,
size_t dimy,
typename resamp_t,
typename float_t,
bool debug>
158 template<
size_t nparts,
size_t dimx,
size_t dimy,
typename resamp_t,
typename float_t,
bool debug>
167 float_t maxOldLogUnNormWts(-std::numeric_limits<float_t>::infinity());
169 for(
size_t ii = 0; ii < nparts; ++ii)
189 float_t sumExp1(0.0);
190 float_t sumExp2(0.0);
191 for(
size_t i = 0; i < nparts; ++i){
193 sumExp2 += std::exp(oldLogUnNormWts[i] - maxOldLogUnNormWts);
195 m_logLastCondLike = maxNumer + std::log(sumExp1) - maxOldLogUnNormWts - std::log(sumExp2);
202 unsigned int rows = testOutput.rows();
203 unsigned int cols = testOutput.cols();
204 Mat numer = Mat::Zero(rows,cols);
205 float_t weightNormConst (0.0);
206 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
214 std::cout <<
"transposed expectation " << fId <<
": " <<
m_expectations[fId].transpose() <<
"\n";
229 for(
size_t ii = 0; ii < nparts; ++ii)
246 for(
size_t i = 0; i < nparts; ++i){
258 unsigned int rows = testOutput.rows();
259 unsigned int cols = testOutput.cols();
260 Mat numer = Mat::Zero(rows,cols);
261 float_t weightNormConst (0.0);
262 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
270 std::cout <<
"transposed expectation " << fId <<
": " <<
m_expectations[fId].transpose() <<
"\n";
287 template<
size_t nparts,
size_t dimx,
size_t dimy,
typename resamp_t,
typename float_t,
bool debug>
294 template<
size_t nparts,
size_t dimx,
size_t dimy,
typename resamp_t,
typename float_t,
bool debug>
302 #endif // BOOTSTRAP_FILTER_H virtual float_t logMuEv(const ssv &x1)=0
Calculate muEv or logmuEv.
Eigen::Matrix< float_t, dimx, 1 > ssv
Definition: bootstrap_filter.h:28
Eigen::Matrix< float_t, dimy, 1 > osv
Definition: bootstrap_filter.h:30
virtual float_t logQ1Ev(const ssv &x1, const osv &y1)=0
Calculate q1Ev or log q1Ev.
std::array< ssv, nparts > arrayStates
Definition: bootstrap_filter.h:34
resamp_t m_resampler
resampler object
Definition: bootstrap_filter.h:133
virtual ssv fSamp(const ssv &xtm1)=0
Sample from the state transition distribution.
std::vector< Mat > m_expectations
expectations E[h(x_t) | y_{1:t}] for user defined "h"s
Definition: bootstrap_filter.h:136
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...
Definition: bootstrap_filter.h:159
BSFilter(const unsigned int &rs=1)
The constructor.
Definition: bootstrap_filter.h:144
auto getExpectations() const -> std::vector< Mat >
return all stored expectations (taken with respect to $p(x_t|y_{1:t})$
Definition: bootstrap_filter.h:295
virtual ~BSFilter()
The (virtual) destructor.
Definition: bootstrap_filter.h:155
Eigen::Matrix< float_t, Eigen::Dynamic, Eigen::Dynamic > Mat
Definition: bootstrap_filter.h:32
virtual ssv q1Samp(const osv &y1)=0
Samples from time 1 proposal.
static constexpr unsigned int num_particles
Definition: bootstrap_filter.h:38
All non Rao-Blackwellized particle filters inherit from this.
arrayFloat m_logUnNormWeights
particle unnormalized weights
Definition: bootstrap_filter.h:124
float_t getLogCondLike() const
Returns the most recent (log-) conditiona likelihood.
Definition: bootstrap_filter.h:288
A base class for the bootstrap particle filter.
Definition: bootstrap_filter.h:23
float_t m_logLastCondLike
log p(y_t|y_{1:t-1}) or log p(y1)
Definition: bootstrap_filter.h:130
unsigned int m_resampSched
resampling schedule (e.g. resample every __ time points)
Definition: bootstrap_filter.h:139
std::array< float_t, nparts > arrayFloat
Definition: bootstrap_filter.h:36
virtual float_t logGEv(const osv &yt, const ssv &xt)=0
Calculate gEv or logGEv.
arrayStates m_particles
particle samples
Definition: bootstrap_filter.h:121
unsigned int m_now
time point
Definition: bootstrap_filter.h:127