27 template<
size_t nparts,
size_t dimx,
size_t dimy,
typename resamp_t,
typename float_t,
bool debug=false>
33 using ssv = Eigen::Matrix<float_t,dimx,1>;
35 using osv = Eigen::Matrix<float_t,dimy,1>;
37 using Mat = Eigen::Matrix<float_t,Eigen::Dynamic,Eigen::Dynamic>;
53 APF(
const unsigned int &rs=1);
80 void filter(
const osv &data,
const std::vector<std::function<
const Mat(
const ssv&)> >& fs = std::vector<std::function<
const Mat(
const ssv&)> >());
130 virtual float_t
logGEv (
const osv &yt,
const ssv &xt) = 0;
162 template<
size_t nparts,
size_t dimx,
size_t dimy,
typename resamp_t,
typename float_t,
bool debug>
172 template<
size_t nparts,
size_t dimx,
size_t dimy,
typename resamp_t,
typename float_t,
bool debug>
176 template<
size_t nparts,
size_t dimx,
size_t dimy,
typename resamp_t,
typename float_t,
bool debug>
186 float_t m3(-std::numeric_limits<float_t>::infinity());
187 float_t m2(-std::numeric_limits<float_t>::infinity());
188 for(
size_t ii = 0; ii < nparts; ++ii)
195 ssv xtm1 = oldPartics[ii];
196 logFirstStageUnNormWeights[ii] +=
logGEv(data,
propMu(xtm1));
199 if(logFirstStageUnNormWeights[ii] > m2)
200 m2 = logFirstStageUnNormWeights[ii];
203 if constexpr(debug) {
204 std::cout <<
"time: " <<
m_now 205 <<
", first stage log unnorm weight: " << logFirstStageUnNormWeights[ii]
215 float_t m1(-std::numeric_limits<float_t>::infinity());
216 float_t first_cll_sum(0.0);
217 float_t second_cll_sum(0.0);
218 float_t third_cll_sum(0.0);
221 for(
size_t ii = 0; ii < nparts; ++ii)
224 second_cll_sum += std::exp( logFirstStageUnNormWeights[ii] - m2 );
228 xtm1k = oldPartics[myKs[ii]];
234 std::cout <<
"time: " <<
m_now 235 <<
", transposed sample: " <<
m_particles[ii].transpose()
245 for(
size_t p = 0; p < nparts; ++p)
247 m_logLastCondLike = m1 + std::log(first_cll_sum) + m2 + std::log(second_cll_sum) - 2*m3 - 2*std::log(third_cll_sum);
257 unsigned int rows = testOutput.rows();
258 unsigned int cols = testOutput.cols();
259 Mat numer = Mat::Zero(rows,cols);
262 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
269 std::cout <<
"transposed expectation " << fId <<
"; " <<
m_expectations[fId] <<
"\n";
283 float_t max(-std::numeric_limits<float_t>::infinity());
284 for(
size_t ii = 0; ii < nparts; ++ii)
293 if constexpr(debug) {
294 std::cout <<
"time: " <<
m_now 296 <<
", transposed sample: " <<
m_particles[ii].transpose()
307 for(
size_t i = 0; i < nparts; ++i){
310 m_logLastCondLike = - std::log( static_cast<float_t>(nparts) ) + max + std::log(sumExp);
318 unsigned int rows = testOutput.rows();
319 unsigned int cols = testOutput.cols();
320 Mat numer = Mat::Zero(rows,cols);
322 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
329 std::cout <<
"transposed expectation " << fId <<
"; " <<
m_expectations[fId] <<
"\n";
345 template<
size_t nparts,
size_t dimx,
size_t dimy,
typename resamp_t,
typename float_t,
bool debug>
352 template<
size_t nparts,
size_t dimx,
size_t dimy,
typename resamp_t,
typename float_t,
bool debug>
virtual ssv fSamp(const ssv &xtm1)=0
Samples from f.
rvsamp::k_gen< nparts, float_t > m_kGen
k generator object (default ctor'd)
Definition: auxiliary_pf.h:153
A base-class for Auxiliary Particle Filtering. Filtering only, no smoothing.
Definition: auxiliary_pf.h:28
std::array< float_t, nparts > m_logUnNormWeights
particle unnormalized weights
Definition: auxiliary_pf.h:138
std::vector< Mat > getExpectations() const
return all stored expectations (taken with respect to $p(x_t|y_{1:t})$
Definition: auxiliary_pf.h:353
std::array< unsigned int, N > sample(const std::array< float_t, N > &logWts)
sample N times from (0,1,...N-1)
Definition: rv_samp.h:842
unsigned int m_rs
the resampling schedule
Definition: auxiliary_pf.h:147
virtual float_t logGEv(const osv &yt, const ssv &xt)=0
Evaluates the log of g.
unsigned int m_now
curren time
Definition: auxiliary_pf.h:141
all rv samplers must inherit from this.
float_t m_logLastCondLike
log p(y_t|y_{1:t-1}) or log p(y1)
Definition: auxiliary_pf.h:144
virtual float_t logMuEv(const ssv &x1)=0
Evaluates the log of mu.
resamp_t m_resampler
resampler object (default ctor'd)
Definition: auxiliary_pf.h:150
virtual ssv q1Samp(const osv &y1)=0
Samples from q1.
virtual ~APF()
The (virtual) destructor.
Definition: auxiliary_pf.h:173
virtual ssv propMu(const ssv &xtm1)=0
Evaluates the proposal distribution taking a Eigen::Matrix<float_t,dimx,1> from the previous time's s...
static constexpr unsigned int num_particles
Definition: auxiliary_pf.h:45
std::array< unsigned int, nparts > arrayUInt
Definition: auxiliary_pf.h:43
Eigen::Matrix< float_t, Eigen::Dynamic, Eigen::Dynamic > Mat
Definition: auxiliary_pf.h:37
std::vector< Mat > m_expectations
expectations E[h(x_t) | y_{1:t}] for user defined "h"s
Definition: auxiliary_pf.h:156
All non Rao-Blackwellized particle filters inherit from this.
APF(const unsigned int &rs=1)
The constructor.
Definition: auxiliary_pf.h:163
std::array< ssv, nparts > m_particles
particle samples
Definition: auxiliary_pf.h:135
Eigen::Matrix< float_t, dimx, 1 > ssv
Definition: auxiliary_pf.h:33
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)...
Definition: auxiliary_pf.h:177
float_t getLogCondLike() const
Get the latest log conditional likelihood.
Definition: auxiliary_pf.h:346
std::array< ssv, nparts > arrayVec
Definition: auxiliary_pf.h:41
Eigen::Matrix< float_t, dimy, 1 > osv
Definition: auxiliary_pf.h:35
virtual float_t logQ1Ev(const ssv &x1, const osv &y1)=0
Evaluates the log of q1.
std::array< float_t, nparts > arrayfloat_t
Definition: auxiliary_pf.h:39