pf
pf_base.h
Go to the documentation of this file.
1 #ifndef PF_BASE_H
2 #define PF_BASE_H
3 
4 #include <map>
5 #include <string>
6 #include <vector>
7 #include <Eigen/Dense>
8 
9 
18 template<typename float_t, size_t dimobs, size_t dimstate>
19 class pf_base {
20 public:
21 
22  /* expose float type to users of this ABCTP */
23  using float_type = float_t;
24 
25  /* observation-sized vector */
26  using osv = Eigen::Matrix<float_t,dimobs,1>;
27 
28  /* state-sized vector */
29  using ssv = Eigen::Matrix<float_t,dimstate,1>;
30 
31  /* state-sized vector */
32  using Mat = Eigen::Matrix<float_t,Eigen::Dynamic,Eigen::Dynamic>;
33 
34  /* a function */
35  using func = std::function<const Mat(const ssv&)>;
36 
37  /* functions */
38  using funcs = std::vector<func>;
39 
40  /* the dimension of each observation vector */
41  static constexpr unsigned int dim_obs = dimobs;
42 
43  /* the dimension of the state vector */
44  static constexpr unsigned int dim_state = dimstate;
45 
46 
50  virtual void filter(const osv &data, const funcs& fs = funcs() ) = 0;
51 
52 
56  virtual float_t getLogCondLike() const = 0;
57  virtual ~pf_base(){};
58 };
59 
60 
70 template<typename float_t, size_t dim_s_state, size_t dim_ns_state, size_t dimobs>
71 class rbpf_base {
72 public:
73 
74  /* observation-sized vector */
75  using osv = Eigen::Matrix<float_t,dimobs,1>;
76 
77  /* sampled-state-size vector */
78  using sssv = Eigen::Matrix<float_t,dim_s_state,1>;
79 
80  /* not-sampled-state-sized vector */
81  using nsssv = Eigen::Matrix<float_t,dim_ns_state,1>;
82 
83  /* matrix */
84  using Mat = Eigen::Matrix<float_t,Eigen::Dynamic,Eigen::Dynamic>;
85 
86  /* a function */
87  using func = std::function<const Mat(const nsssv&, const sssv&)>;
88 
89  /* functions */
90  using funcs = std::vector<func>;
91 
92  /* the size of the sampled state portion */
93  static constexpr unsigned int dim_sampled_state = dim_s_state;
94 
95  /* the size of the not sampled state portion */
96  static constexpr unsigned int dim_not_sampled_state = dim_ns_state;
97 
98  /* the size of the observations */
99  static constexpr unsigned int dim_obs = dimobs;
100 
101 
105  virtual void filter(const osv &data, const funcs& fs = funcs() ) = 0;
106  virtual ~rbpf_base(){};
107 };
108 
109 
116 template<size_t dimx, size_t dimy, typename float_t>
117 class ForwardMod {
118 public:
119 
120  /* state-sized vector */
121  using ssv = Eigen::Matrix<float_t, dimx, 1>;
122 
123  /* observation-sized vector */
124  using osv = Eigen::Matrix<float_t, dimy, 1>;
125 
126  /* a pair of paths (xts first, yts second) */
127  using aPair = std::pair<std::vector<ssv>, std::vector<osv> >;
128 
129 
133  aPair sim_forward(unsigned int T);
134 
139  virtual ssv muSamp() = 0;
140 
141 
146  virtual ssv fSamp (const ssv &xtm1) = 0;
147 
148 
153  virtual osv gSamp (const ssv &xt) = 0;
154 };
155 
156 
157 template<size_t dimx, size_t dimy, typename float_t>
158 auto ForwardMod<dimx,dimy,float_t>::sim_forward(unsigned int T) -> aPair {
159 
160  std::vector<ssv> xs;
161  std::vector<osv> ys;
162 
163  // time 1
164  xs.push_back(this->muSamp());
165  ys.push_back(this->gSamp(xs[0]));
166 
167  // time > 1
168  for(size_t i = 1; i < T; ++i) {
169 
170  auto xt = this->fSamp(xs[i-1]);
171  xs.push_back(xt);
172  ys.push_back(this->gSamp(xt));
173  }
174 
175  return aPair(xs, ys);
176 }
177 
178 #endif // PF_BASE_H
aPair sim_forward(unsigned int T)
simulates forward through time
Definition: pf_base.h:158
Definition: pf_base.h:117
Definition: pf_base.h:19
virtual void filter(const osv &data, const funcs &fs=funcs())=0
the filtering function that must be defined
virtual float_t getLogCondLike() const =0
the getter method that must be defined (for conditional log-likelihood)
Definition: pf_base.h:71