Helios++
Helios software for LiDAR simulations
fluxionum::DiffDesignMatrix< TimeType, VarType > Class Template Reference

The heart of a differential design matrix is the idea that the columns are defining the values of variables differentiated over time. The \(m-1\) rows are samples of the differential behavior of the function \(x_j(t_i)\) for the \(j\)-th column as variable. For then, it is possible to interpolate the following vector of derivatives. More...

#include <DiffDesignMatrix.h>

Inheritance diagram for fluxionum::DiffDesignMatrix< TimeType, VarType >:
Collaboration diagram for fluxionum::DiffDesignMatrix< TimeType, VarType >:

Public Member Functions

 DiffDesignMatrix (vector< string > const &columnNames=vector< string >(0), string const timeName="time", DiffDesignMatrixType diffType=DiffDesignMatrixType::UNKNOWN)
 Build a DiffDesignMatrix with no data. More...
 
 DiffDesignMatrix (arma::Col< TimeType > const &t, arma::Mat< VarType > const &A, string const &timeName="time", vector< string > const &columnNames=vector< string >(0), DiffDesignMatrixType diffType=DiffDesignMatrixType::UNKNOWN)
 Build a DiffDesignMatrix from given armadillo column vector of time values and given armadillo matrix of attributes. More...
 
 DiffDesignMatrix (string const &path, string const &timeName="time", DiffDesignMatrixType diffType=DiffDesignMatrixType::UNKNOWN)
 Build a DiffDesignMatrix from data in file at given path. More...
 
 DiffDesignMatrix (TemporalDesignMatrix< TimeType, VarType > const &tdm, DiffDesignMatrixType const diffType, bool const sort=true)
 Build a DiffDesignMatrix from given TemporalDesignMatrix. More...
 
VarType & operator() (size_t const i, size_t const j) override
 Access to the \(a_{ij}\) component of the \(A\) design matrix. More...
 
TimeType & operator[] (size_t const i)
 Access to the \(t_{i}\) component of the \(\vec{t}\) sorted time vector. More...
 
virtual void shiftTime (TimeType const x)
 Add given time \(x\) to all time values in vector \(\vec{t}\). More...
 
size_t getNumRows () const override
 Obtain the number of rows of the DiffDesignMatrix \(A\). More...
 
size_t getNumColumns () const override
 Obtain the number of columns of the DiffDesignMatrix \(A\). More...
 
size_t getNumElements () const override
 Obtain the number of elements of the DiffDesignMatrix \(A\). More...
 
DiffDesignMatrixType getDiffType () const
 Obtain the type of differential of the DiffDesignMatrix. More...
 
string const & getTimeName () const
 Obtain the name of the time attribute. More...
 
arma::Col< TimeType > const & getTimeVector () const
 Obtain the sorted series of time values DiffDesignMatrix::t. More...
 
arma::Mat< VarType > const & getA () const
 Obtain a constant/read reference to the \(A\) matrix. More...
 
- Public Member Functions inherited from fluxionum::AbstractDesignMatrix< VarType >
 AbstractDesignMatrix (vector< string > const &columnNames=vector< string >(0))
 Default constructor for the AbstractDesignMatrix. More...
 
VarType & operator() (size_t const i, string const columnName)
 Like the AbstractDesignMatrix::operator()(size_t const, size_t const) method but specifying the column by its name instead of its index. More...
 
bool hasColumnNames () const
 Check whether there are available column names for the AbstractDesignMatrix (true) or not (false) More...
 
string const & getColumnName (size_t const j) const
 Obtain the name of the \(j\)-th column. More...
 
void setColumnName (size_t const j, string const &columnName)
 Set the name of the \(j\)-th column. More...
 
vector< string > const & getColumnNames () const
 Obtain a constant/read reference to the column names. More...
 
void setColumnNames (vector< string > const &columnNames)
 Obtain a constant/read reference to the column names. More...
 

Protected Attributes

enum DiffDesignMatrixType diffType
 Specify the type of differential of the DiffDesignMatrix. More...
 
TimeType ta
 The starting value for differential time interval, \(t_a\).
 
TimeType tb
 The end value for the differential time interval, \(t_b\).
 
arma::Col< TimeType > t
 The sorted series of time values so \(\forall i, t_{i+1} > t_{i}\) is strictly satisfied and the \(i\)-th time value corresponds to the approximated \(\overrightarrow{\frac{dx_i}{dt}}\) row vector in the \(\frac{DX}{DT} \in \mathbb{R}^{(m-1) \times n}\) matrix. It is used to build the \(T\) matrix. More...
 
arma::Mat< VarType > A
 The \(\frac{DX}{DT}\) matrix, whether it comes from forward finite differences or from central finite differences.
 
string timeName
 The name of the time field. More...
 
- Protected Attributes inherited from fluxionum::AbstractDesignMatrix< VarType >
vector< string > columnNames
 The column names for the DesignMatrix. It can be either an empty vector when no column names are specified or a vector with as many names as columns (in the same order)
 

Additional Inherited Members

- Protected Member Functions inherited from fluxionum::AbstractDesignMatrix< VarType >
size_t translateColumnNameToIndex (string const &columnName) const
 Find the corresponding column name for given index. More...
 

Detailed Description

template<typename TimeType, typename VarType>
class fluxionum::DiffDesignMatrix< TimeType, VarType >

The heart of a differential design matrix is the idea that the columns are defining the values of variables differentiated over time. The \(m-1\) rows are samples of the differential behavior of the function \(x_j(t_i)\) for the \(j\)-th column as variable. For then, it is possible to interpolate the following vector of derivatives.

Author
Alberto M. Esmoris Pena
Version
1.0

\[ \overrightarrow{\frac{dx}{dt}} = \left[\begin{array}{c} \frac{dx_1}{dt} \\ \vdots \\ \frac{dx_n}{dt} \end{array}\right] \]

using each \(i\)-th row to fit each \(j\)-th variable

To explain in detail the fundamentals of a DiffDesignMatrix, assume a restricted time domain such that \(t \in [t_a, t_b]\). The values defining the extremes of the interval correspond to the minimum and the maximum time values in the TemporalDesignMatrix that was used to generate the DiffDesignMatrix. There is also a time vector available, such that for the \(m\) points (matrix rows) there are also \(m\) associated time values \(\vec{t} = \left(t_1, \ldots, t_m\right)\). Thus, for each \(t_i\) there exists a vector \(\vec{x_i'} \in \mathbb{R}^{n}\) which defines the values that the attributes took at time \(t_i\). From now on, it will be assumed that \(\forall i, t_{i+1} > t_{i}\). It is, the indices of time values and their associated vectors are sorted in time with no repetitions.

Now, lets define the \(\frac{DX}{DT} \in \mathbb{R}^{(m-1) \times n}\) matrix as a forward finite differences based matrix:

\[ \frac{DX}{DT} = \frac{\Delta X}{\Delta T} = \left[\begin{array}{ccc} \frac{x_{21}' - x_{11}'}{t_2 - t_1} & \ldots & \frac{x_{2n}' - x_{1n}'}{t_2 - t_1} \\ \vdots & \ddots & \vdots \\ \frac{x_{m1}' - x_{m-1,1}'}{t_m - t_{m-1}} & \ldots & \frac{x_{mn}' - x_{m-1,n}'}{t_m - t_{m-1}} \end{array}\right] \]

The former matrix can also be defined as a central finite differences based matrix, as will be later explained. However, no matter how the matrix is defined, it will be noted as matrix \(A \in \mathbb{R}^{l \times n}\) from now on, for the sake of simplicity. This matrix is expected to define accurate samples of the derivative for each of the variables.

\[ A = \left[\begin{array}{ccc} a_{11} & \ldots & a_{1n} \\ \vdots & \ddots & \vdots \\ a_{m1} & \ldots & a_{mn} \end{array}\right] \]

An example of fitting is developed based on the Moore-Penrose pseudoinverse of a matrix. For this purpose, assume \(l\) time samples of \(n\) variables each, for which it is possible to build an approximation by means of a polynomial series. Then, it must be that there are \(n\) functions, such that for any \(j\)-th variable there is a function \(f_j(t_i) = \sum_{k=0}^{k_j} \omega_{jk}t_{i}^{k}\approx a_{ij}\) that approximates \(\frac{dx_j}{dt}\). Since there are \(n\) polynomial series of order \(k_j\) each, it is known that there are \(n\) coefficient vectors \(\vec{\omega_{j}} \in \mathbb{R}^{k_j + 1}\). Also, there is an implicit matrix \(T_j \in \mathbb{R}^{l \times (k_j+1)}\) for each variable as follows:

\[ T_j = \left[\begin{array}{ccc} t_1^0 & \dots & t_1^k \\ \vdots & \ddots & \vdots \\ t_l^0 & \dots & t_l^k \end{array}\right] \]

Let \(\vec{a_j} = \left(a_{1j}, \ldots, a_{lj}\right)\) be the \(j\)-th column vector of matrix \(A_{l \times n}\) and thus the system \(T_j \vec{\omega_j} = \vec{a_j}\) arises. The polynomial coefficients can then be found simply by solving them with Moore-Penrose pseudoinverse such that \(\vec{\omega_j} = (T_j^T T_j)^{-1}T_j^T \vec{a_j} = T_j^\dagger \vec{a_j}\). This way the derivatives for each of the \(n\) variables are approximated each by its own function:

\[ \forall j, f_j(t) = \sum_{k=0}^{k_j}{\omega_{jk}t^k} \approx \frac{d}{dt}x_j(t) \]

It was stated before that the \(A \in \mathbb{R}^{l \times n}\) matrix was generated from forward finite differences so \(l=m-1\) and thus \(A = \frac{DX}{DT} = \frac{\Delta X}{\Delta T} \in \mathbb{R}^{m-1 \times n}\). However, it is also possible to define \(A = \frac{DX}{DT} = \frac{\delta X}{\delta T} \in \mathbb{R}^{m-2 \times n}\) so it comes from central finite differences. Central finite differences might lead to a more accurate approximation of the derivatives, at the expenses of losing \(2\) data points instead of \(1\). The first step to compute central finite differences is to define a new time vector such that \(\vec{\tau} = \left(\frac{t_1+t_2}{2}, \ldots, \frac{t_{m-1}+t_m}{2}\right) = \left(\tau_1, \ldots, \tau_{m-1}\right) \in \mathbb{R}^{m-1}\). Thus, for central finite differences the \(T \in \mathbb{R}^{(m-2) \times (k_{j}+1)}\) matrix has one less row and its values come from the first \(m-2\) values of vector \(\vec{\tau} \in \mathbb{R}^{m-1}\) instead of the first \(m-1\) values of the time vector \(\vec{t} \in \mathbb{R}^{m}\). The second step is to compute the matrix \(\frac{\delta X}{\delta T}\) as:

\[ \frac{DX}{DT} = \frac{\delta X}{\delta T} = \left[\begin{array}{ccc} \frac{x_{31}' - x_{11}'}{2(\tau_2 - \tau_1)} & \ldots & \frac{x_{3n}' - x_{1n}'}{2(\tau_2 - \tau_1)} \\ \vdots & \ddots & \vdots \\ \frac{x_{m1}' - x_{m-2,1}'}{2(\tau_{m-1} - \tau_{m-2})} & \ldots & \frac{x_{mn}' - x_{m-2,n}'}{2(\tau_{m-1} - \tau_{m-2})} \end{array}\right] = \left[\begin{array}{ccc} \frac{x_{31}' - x_{11}'}{t_3 - t_1} & \ldots & \frac{x_{3n}' - x_{1n}'}{t_3 - t_1} \\ \vdots & \ddots & \vdots \\ \frac{x_{m1}' - x_{m-2,1}'}{t_m - t_{m-2}} & \ldots & \frac{x_{mn}' - x_{m-2,n}'}{t_m - t_{m-2}} \end{array}\right] \]

The fitting example for central finite differences would be exactly as described before, but now the approximation would be based on \(\vec{\tau}\) and \(A=\frac{\delta X}{\delta T}\) instead of \(\vec{t}\) and \(A=\frac{\Delta X}{\Delta T}\). It is, the function \(f_j(t)\) is fitted from \(f_j(\tau_{i}) = \sum_{k=0}^{k_j}{\omega_{jk}\tau_i^k} \approx a_{ij}\).

Template Parameters
TimeTypeThe time's domain
VarTypeThe non time's domain (basic DesignMatrix domain)
See also
fluxionum::AbstractDesignMatrix
fluxionum::TemporalDesignMatrix
fluxionum::DiffDesignMatrixType

Constructor & Destructor Documentation

◆ DiffDesignMatrix() [1/4]

template<typename TimeType , typename VarType >
fluxionum::DiffDesignMatrix< TimeType, VarType >::DiffDesignMatrix ( vector< string > const &  columnNames = vector<string>(0),
string const  timeName = "time",
DiffDesignMatrixType  diffType = DiffDesignMatrixType::UNKNOWN 
)
inline

Build a DiffDesignMatrix with no data.

Parameters
columnNamesEither the name for each column or an empty vector if there are no names

◆ DiffDesignMatrix() [2/4]

template<typename TimeType , typename VarType >
fluxionum::DiffDesignMatrix< TimeType, VarType >::DiffDesignMatrix ( arma::Col< TimeType > const &  t,
arma::Mat< VarType > const &  A,
string const &  timeName = "time",
vector< string > const &  columnNames = vector<string>(0),
DiffDesignMatrixType  diffType = DiffDesignMatrixType::UNKNOWN 
)
inline

Build a DiffDesignMatrix from given armadillo column vector of time values and given armadillo matrix of attributes.

Parameters
tArmadillo column vector of time values to build DiffDesignMatrix from
AArmadillo matrix of attributes to build DiffDesignMatrix from
columnNamesEither the name for each column or an empty vector if there are no names

◆ DiffDesignMatrix() [3/4]

template<typename TimeType , typename VarType >
fluxionum::DiffDesignMatrix< TimeType, VarType >::DiffDesignMatrix ( string const &  path,
string const &  timeName = "time",
DiffDesignMatrixType  diffType = DiffDesignMatrixType::UNKNOWN 
)
inline

Build a DiffDesignMatrix from data in file at given path.

Parameters
pathPath to the file containing the data for the DiffDesignMatrix
columnNamesThe default column names to be used if no column names are read from file at given path

◆ DiffDesignMatrix() [4/4]

template<typename TimeType , typename VarType >
fluxionum::DiffDesignMatrix< TimeType, VarType >::DiffDesignMatrix ( TemporalDesignMatrix< TimeType, VarType > const &  tdm,
DiffDesignMatrixType const  diffType,
bool const  sort = true 
)
inline

Build a DiffDesignMatrix from given TemporalDesignMatrix.

Parameters
tdmThe TemporalDesignMatrix to be used to build the DiffDesignMatrix
sortWhen true, the DiffDesignMatrix construction will lead to an expensive sort operation to assure the finite differences computation makes sense. When false, the expensive sort computation will be skipped as it is expected that the TemporalDesignMatrix is already sorted
diffTypeThe type of differencing method
See also
TemporalDesignMatrix::toDiffDesignMatrix

Member Function Documentation

◆ getA()

template<typename TimeType , typename VarType >
arma::Mat<VarType> const& fluxionum::DiffDesignMatrix< TimeType, VarType >::getA ( ) const
inline

Obtain a constant/read reference to the \(A\) matrix.

Returns
Constant/read reference to the \(A\) matrix
See also
fluxionum::DiffDesignMatrix::A

◆ getDiffType()

template<typename TimeType , typename VarType >
DiffDesignMatrixType fluxionum::DiffDesignMatrix< TimeType, VarType >::getDiffType ( ) const
inline

Obtain the type of differential of the DiffDesignMatrix.

Returns
The type of tdifferential of the DiffDesignMatrix
See also
fluxionum::DiffDesignMatrix::diffType

◆ getNumColumns()

template<typename TimeType , typename VarType >
size_t fluxionum::DiffDesignMatrix< TimeType, VarType >::getNumColumns ( ) const
inlineoverridevirtual

Obtain the number of columns of the DiffDesignMatrix \(A\).

Returns
The number of columns of the DiffDesignMatrix \(A\)
See also
fluxionum::DiffDesignMatrix::A
fluxionum::AbstractDesignMatrix::getNumColumns

Implements fluxionum::AbstractDesignMatrix< VarType >.

◆ getNumElements()

template<typename TimeType , typename VarType >
size_t fluxionum::DiffDesignMatrix< TimeType, VarType >::getNumElements ( ) const
inlineoverridevirtual

Obtain the number of elements of the DiffDesignMatrix \(A\).

Returns
The number of elements of the DiffDesignMatrix \(A\)
See also
fluxionum::DiffDesignMatrix::A
fluxionum::AbstractDesignMatrix::getNumElements

Implements fluxionum::AbstractDesignMatrix< VarType >.

◆ getNumRows()

template<typename TimeType , typename VarType >
size_t fluxionum::DiffDesignMatrix< TimeType, VarType >::getNumRows ( ) const
inlineoverridevirtual

Obtain the number of rows of the DiffDesignMatrix \(A\).

Returns
The number of rows of the DiffDesignMatrix \(A\)
See also
fluxionum::DiffDesignMatrix::A
fluxionum::AbstractDesignMatrix::getNumRows

Implements fluxionum::AbstractDesignMatrix< VarType >.

◆ getTimeName()

template<typename TimeType , typename VarType >
string const& fluxionum::DiffDesignMatrix< TimeType, VarType >::getTimeName ( ) const
inline

Obtain the name of the time attribute.

Returns
The name of the time attribute
See also
fluxionum::TemporalDesignMatrix::timeName

◆ getTimeVector()

template<typename TimeType , typename VarType >
arma::Col<TimeType> const& fluxionum::DiffDesignMatrix< TimeType, VarType >::getTimeVector ( ) const
inline

Obtain the sorted series of time values DiffDesignMatrix::t.

Returns
The sorted series of time values DiffDesignMatrix::t
See also
fluxionum::DiffDesignMatrix::t

◆ operator()()

template<typename TimeType , typename VarType >
VarType& fluxionum::DiffDesignMatrix< TimeType, VarType >::operator() ( size_t const  i,
size_t const  j 
)
inlineoverridevirtual

Access to the \(a_{ij}\) component of the \(A\) design matrix.

Parameters
iThe row of the element being accessed
jThe column of the element being accessed
Returns
Reference to the element at \(i\)-th row and \(j\)-th column
See also
fluxionum::DiffDesignMatrix::A
AbstractDesignMatrix::operator()(size_t const, size_t const)

Implements fluxionum::AbstractDesignMatrix< VarType >.

◆ operator[]()

template<typename TimeType , typename VarType >
TimeType& fluxionum::DiffDesignMatrix< TimeType, VarType >::operator[] ( size_t const  i)
inline

Access to the \(t_{i}\) component of the \(\vec{t}\) sorted time vector.

Parameters
iThe index of the time being accessed (its position in the time vector)
Returns
Reference to the \(t_{i}\) component of the sorted time vector \(\vec{t}\)

◆ shiftTime()

template<typename TimeType , typename VarType >
virtual void fluxionum::DiffDesignMatrix< TimeType, VarType >::shiftTime ( TimeType const  x)
inlinevirtual

Add given time \(x\) to all time values in vector \(\vec{t}\).

\vec{t'} = \vec{t} + x = \left[\begin{array}{c} t_1 + x \ \vdots t_m + x \end{array}\right.

Parameters
xThe time shift for each time component

Member Data Documentation

◆ diffType

template<typename TimeType , typename VarType >
enum DiffDesignMatrixType fluxionum::DiffDesignMatrix< TimeType, VarType >::diffType
protected

Specify the type of differential of the DiffDesignMatrix.

See also
fluxionum::DiffDesignMatrixType

◆ t

template<typename TimeType , typename VarType >
arma::Col<TimeType> fluxionum::DiffDesignMatrix< TimeType, VarType >::t
protected

The sorted series of time values so \(\forall i, t_{i+1} > t_{i}\) is strictly satisfied and the \(i\)-th time value corresponds to the approximated \(\overrightarrow{\frac{dx_i}{dt}}\) row vector in the \(\frac{DX}{DT} \in \mathbb{R}^{(m-1) \times n}\) matrix. It is used to build the \(T\) matrix.

Notice the \(\vec{t}\) vector is in fact the \(\vec{\tau}\) vector when the DiffDesignMatrix is built from central finite differences

WARNING!: The DiffDesignMatrix::t vector must not be confused with TemporalDesignMatrix::t vector. For the case of forward finite differences the DiffDesignMatrix::t vector is \(\vec{t} = \left(t_1, \ldots, t_{m-1}\right) \in \mathbb{R}^{m-1}\). It is, all values from DiffDesignMatrix::t but the last one. In the case of central finite differences the DiffDesignMatrix::t vector is \(\vec{t} = \left(\frac{t_1+t_2}{2}, \ldots, \frac{t_{m-2}-t_{m-1}}{2}\right) \in \mathbb{R}^{m-2}\). It is, all values from the \(\vec{\tau}\) vector but the last one.

◆ timeName

template<typename TimeType , typename VarType >
string fluxionum::DiffDesignMatrix< TimeType, VarType >::timeName
protected

The name of the time field.

By default, it is "time"


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