Helios++
Helios software for LiDAR simulations
TemporalDesignMatrix.h
1 #ifndef _FLUXIONUM_TEMPORAL_DESIGN_MATRIX_H_
2 
3 #include <fluxionum/DesignMatrix.h>
4 #include <fluxionum/FluxionumTypes.h>
5 
6 #include <memory>
7 #include <unordered_map>
8 
9 namespace fluxionum {
10 
11 template <typename TimeType, typename VarType> class DiffDesignMatrix;
12 
13 using std::string;
14 using std::shared_ptr;
15 using std::make_shared;
16 using std::unordered_map;
17 
73 template <typename TimeType, typename VarType>
74 class TemporalDesignMatrix : public DesignMatrix<VarType> {
75 protected:
76  // *** USING *** //
77  // *************** //
79 
80  // *** ATTRIBUTES *** //
81  // ******************** //
87  arma::Col<TimeType> t;
93  string timeName;
94 
95 public:
96  // *** STATIC METHODS *** //
97  // ************************ //
107  static inline arma::Mat<VarType> extractNonTimeMatrix(
108  arma::Mat<VarType> const &X, size_t const timeColumnIndex
109  ){
110  arma::Mat<VarType> Y(X);
111  Y.shed_col(timeColumnIndex);
112  return Y;
113  }
123  static inline vector<string> extractNonTimeNames(
124  vector<string> const &names, size_t const timeColumnIndex
125  ){
126  vector<string> newNames;
127  size_t numNames = names.size();
128  for(size_t i = 0 ; i < numNames ; ++i){
129  if(i==timeColumnIndex) continue;
130  newNames.push_back(names[i]);
131  }
132  return newNames;
133  }
142  static inline arma::Col<TimeType> extractTimeVector(
143  arma::Mat<VarType> const &X, size_t const timeColumnIndex
144  ){
145  return arma::Col<TimeType>(X.col(timeColumnIndex));
146  }
147 
148 
149  // *** CONSTRUCTION / DESTRUCTION *** //
150  // ************************************ //
161  DesignMatrix<VarType> const &designMatrix,
162  size_t const timeColumnIndex,
163  string const timeName="time",
164  vector<string> const &columnNames=vector<string>(0)
165  ) :
166  DesignMatrix<VarType>(
168  designMatrix.getX(), timeColumnIndex
169  ),
170  designMatrix.hasColumnNames() ?
172  designMatrix.getColumnNames(),
173  timeColumnIndex
174  ) : columnNames
175  ),
176  t(extractTimeVector(designMatrix.getX(), timeColumnIndex)),
177  timeName(
178  designMatrix.hasColumnNames() ?
179  designMatrix.getColumnName(timeColumnIndex) :
180  timeName
181  )
182  {}
191  DesignMatrix<VarType> const &designMatrix,
192  arma::Col<TimeType> const &timeVector,
193  string const timeName="time"
194  ) :
195  DesignMatrix<VarType>(designMatrix),
196  t(timeVector),
198  {}
207  arma::Mat<VarType> const &X,
208  size_t const timeColumnIndex,
209  string const timeName="time",
210  vector<string> const &columnNames=vector<string>(0)
211  ) :
213  X, timeColumnIndex
214  ), columnNames),
215  t(extractTimeVector(X, timeColumnIndex)),
217  {}
226  arma::Mat<VarType> const &X,
227  arma::Col<TimeType> const &t,
228  string const timeName="time",
229  vector<string> const &columnNames=vector<string>(0)
230  ) :
231  DesignMatrix<VarType>(X, columnNames),
232  t(t),
234  {}
245  string const &path,
246  string const &sep=",",
247  string const &timeName="time"
248  ){
250  std::unordered_map<string, string> kv;
251  DesignMatrix<VarType> const dm = reader.read(&kv);
252  size_t const tCol = (size_t) std::strtoul(
253  kv.at("TIME_COLUMN").c_str(), nullptr, 10
254  );
256  dm,
257  tCol,
258  dm.hasColumnNames() ? dm.getColumnName(tCol) : timeName
259  );
260  }
261  virtual ~TemporalDesignMatrix() = default;
262 
263 
264  // *** OPERATORS *** //
265  // ******************* //
274  inline TimeType& operator[] (size_t const i) {return t.at(i);}
275 
276 
277  // *** METHODS *** //
278  // ***************** //
293  DiffDesignMatrixType diffType = \
294  DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES,
295  bool const sort=true
296  ) const;
304  shared_ptr<DiffDesignMatrix<TimeType, VarType>> toDiffDesignMatrixPointer(
305  DiffDesignMatrixType diffType = \
306  DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES,
307  bool const sort=true
308  ) const;
309 
315  void mergeInPlace(DesignMatrix<VarType> const &dm) override{
317  TemporalDesignMatrix const &tdm = static_cast<
319  >(dm);
320  t.insert_rows(t.n_rows, tdm.getTimeVector());
321  }
327  inline void dropRows(arma::uvec const &indices) override{
329  t.shed_rows(indices);
330  }
331 
343  virtual void shiftTime(TimeType const x){
344  t += x;
345  }
351  virtual void sortByTime(){
352  arma::uvec const tSort = arma::sort_index(getTimeVector());
353  t = t(tSort);
354  X = X.rows(tSort);
355  }
356 
391  virtual size_t slopeFilter(VarType const tau){
392  // Obtain DX
393  arma::Mat<VarType> DX = arma::diff(this->getX());
394  DX.each_col() /= arma::diff(getTimeVector());
395  DX = arma::diff(DX);
396  // Filter by DX
397  size_t const m = DX.n_rows;
398  arma::Row<VarType> S(DX.n_cols, arma::fill::zeros);
399  vector<unsigned long long> remove;
400  for(size_t i = 0 ; i < m ; ++i){
401  S += DX.row(i);
402  VarType const maxDiff = arma::max(arma::abs(S));
403  if(maxDiff > tau) S = arma::zeros<arma::Row<VarType>>(DX.n_cols);
404  else remove.push_back(i+1);
405  }
406  dropRows(remove);
407  return remove.size();
408  }
409 
410 
411  // *** GETTERs and SETTERs *** //
412  // ***************************** //
418  inline arma::Col<TimeType> const & getTimeVector() const {return t;}
424  inline string const & getTimeName() const {return timeName;}
430  inline void setTimeName(string const &timeName)
431  {this->timeName = timeName;}
432 };
433 
434 
435 }
436 
437 #define _FLUXIONUM_TEMPORAL_DESIGN_MATRIX_H_
438 #include <fluxionum/TemporalDesignMatrix.tpp>
439 #endif
string const & getColumnName(size_t const j) const
Obtain the name of the -th column.
Definition: AbstractDesignMatrix.h:111
vector< string > columnNames
The column names for the DesignMatrix. It can be either an empty vector when no column names are spec...
Definition: AbstractDesignMatrix.h:35
bool hasColumnNames() const
Check whether there are available column names for the AbstractDesignMatrix (true) or not (false)
Definition: AbstractDesignMatrix.h:104
vector< string > const & getColumnNames() const
Obtain a constant/read reference to the column names.
Definition: AbstractDesignMatrix.h:126
This class represents a DesignMatrix .
Definition: DesignMatrix.h:41
virtual void mergeInPlace(DesignMatrix const &dm)
Merge given DesignMatrix into this DesignMatrix.
Definition: DesignMatrix.h:119
void dropRows(vector< long unsigned int > const &indices)
Definition: DesignMatrix.h:179
arma::Mat< VarType > X
The design matrix .
Definition: DesignMatrix.h:52
arma::Mat< VarType > const & getX() const
Obtain a constant/read reference to the matrix.
Definition: DesignMatrix.h:219
The heart of a differential design matrix is the idea that the columns are defining the values of var...
Definition: DiffDesignMatrix.h:169
This class represents a DesignMatrix where each row satisfy that its elements can be placed into a co...
Definition: TemporalDesignMatrix.h:74
arma::Col< TimeType > t
The time values such that is associated with the -th row of the design matrix . It must strictly sat...
Definition: TemporalDesignMatrix.h:87
TemporalDesignMatrix(arma::Mat< VarType > const &X, arma::Col< TimeType > const &t, string const timeName="time", vector< string > const &columnNames=vector< string >(0))
Build a TemporalDesignMatrix from given DesignMatrix and time vector .
Definition: TemporalDesignMatrix.h:225
shared_ptr< DiffDesignMatrix< TimeType, VarType > > toDiffDesignMatrixPointer(DiffDesignMatrixType diffType=DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES, bool const sort=true) const
Like fluxionum::TemporalDesignMatrix::toDiffDesignMatrix but returning a pointer to the object.
TimeType & operator[](size_t const i)
Access to the component of the time vector.
Definition: TemporalDesignMatrix.h:274
TemporalDesignMatrix(string const &path, string const &sep=",", string const &timeName="time")
Build a TemporalDesignMatrix from data in file at given path and specified time column.
Definition: TemporalDesignMatrix.h:244
TemporalDesignMatrix(DesignMatrix< VarType > const &designMatrix, size_t const timeColumnIndex, string const timeName="time", vector< string > const &columnNames=vector< string >(0))
Build a TemporalDesignMatrix from given DesignMatrix and specified time column.
Definition: TemporalDesignMatrix.h:160
virtual size_t slopeFilter(VarType const tau)
Apply a filter based on the first derivative understood as the slope.
Definition: TemporalDesignMatrix.h:391
string const & getTimeName() const
Obtain the name of the time attribute.
Definition: TemporalDesignMatrix.h:424
static arma::Col< TimeType > extractTimeVector(arma::Mat< VarType > const &X, size_t const timeColumnIndex)
Do a copy of the time column from given DesignMatrix .
Definition: TemporalDesignMatrix.h:142
DiffDesignMatrix< TimeType, VarType > toDiffDesignMatrix(DiffDesignMatrixType diffType=DiffDesignMatrixType::FORWARD_FINITE_DIFFERENCES, bool const sort=true) const
Build a DiffDesignMatrix from the TemporalDesignMatrix.
arma::Col< TimeType > const & getTimeVector() const
Obtain a constant/read reference to the time vector .
Definition: TemporalDesignMatrix.h:418
void dropRows(arma::uvec const &indices) override
Extend DesignMatrix::dropRows(arma::uvec const &) method so the time values are also dropped.
Definition: TemporalDesignMatrix.h:327
virtual void sortByTime()
Sort all rows of matrix and all components of vector so , considering the row is associated with ...
Definition: TemporalDesignMatrix.h:351
static arma::Mat< VarType > extractNonTimeMatrix(arma::Mat< VarType > const &X, size_t const timeColumnIndex)
Do a copy of the DesignMatrix with all its columns but removing the time column.
Definition: TemporalDesignMatrix.h:107
string timeName
The name of the time field in the original DesignMatrix.
Definition: TemporalDesignMatrix.h:93
TemporalDesignMatrix(DesignMatrix< VarType > const &designMatrix, arma::Col< TimeType > const &timeVector, string const timeName="time")
Build a TemporalDesignMatrix from given DesignMatrix and time vector.
Definition: TemporalDesignMatrix.h:190
TemporalDesignMatrix(arma::Mat< VarType > const &X, size_t const timeColumnIndex, string const timeName="time", vector< string > const &columnNames=vector< string >(0))
Build a TemporalDesignMatrix from given matrix and specified time column.
Definition: TemporalDesignMatrix.h:206
void setTimeName(string const &timeName)
Set the name of the time attribute.
Definition: TemporalDesignMatrix.h:430
void mergeInPlace(DesignMatrix< VarType > const &dm) override
Extend DesignMatrix::mergeInPlace method so the time values are also merged.
Definition: TemporalDesignMatrix.h:315
virtual void shiftTime(TimeType const x)
Add given time to all time values in vector .
Definition: TemporalDesignMatrix.h:343
static vector< string > extractNonTimeNames(vector< string > const &names, size_t const timeColumnIndex)
Do a copy of the names of the DesignMatrix but discarding the name of the time column.
Definition: TemporalDesignMatrix.h:123
Class to read design matrices.
Definition: DesignMatrixReader.h:44
virtual fluxionum::DesignMatrix< VarType > read(unordered_map< string, string > *keyval=nullptr)
Read the design matrix.