gr::MatchBase< PointType, _TransformVisitor, OptExts > Class Template Reference

Abstract class for registration algorithms. More...

#include <matchBase.h>

+ Inheritance diagram for gr::MatchBase< PointType, _TransformVisitor, OptExts >:
+ Collaboration diagram for gr::MatchBase< PointType, _TransformVisitor, OptExts >:

Classes

class  Options
 
struct  PosMutablePoint
 A convenience class used to wrap (any) PointType to allow mutation of position of point samples for internal computations. More...
 

Public Types

using LogLevel = Utils::LogLevel
 
using MatrixType = Eigen::Matrix< Scalar, 4, 4 >
 
using OptionsType = gr::Utils::CRTP< OptExts..., Options >
 
using Scalar = typename PointType::Scalar
 
using TransformVisitor = _TransformVisitor
 
using VectorType = typename PointType::VectorType
 

Public Member Functions

template<typename InputRange1 , typename InputRange2 , template< typename > typename Sampler>
Scalar ComputeTransformation (const InputRange1 &P, const InputRange2 &Q, Eigen::Ref< MatrixType > transformation, const Sampler< PointType > &sampler, TransformVisitor &v)
 Computes an approximation of the best LCP (directional) from Q to P and the rigid transformation that realizes it. The input sets may or may not contain normal information for any point. More...
 
const std::vector< PosMutablePoint > & getFirstSampled () const
 Read access to the sampled clouds used for the registration. More...
 
const std::vector< PosMutablePoint > & getSecondSampled () const
 Read access to the sampled clouds used for the registration. More...
 
EIGEN_MAKE_ALIGNED_OPERATOR_NEW MatchBase (const OptionsType &options, const Utils::Logger &logger)
 
virtual ~MatchBase ()
 

Protected Member Functions

template<typename Coordinates >
bool ComputeRigidTransformation (const Coordinates &ref, const Coordinates &candidate, const Eigen::Matrix< Scalar, 3, 1 > &centroid1, Eigen::Matrix< Scalar, 3, 1 > centroid2, Eigen::Ref< MatrixType > transform, Scalar &rms_, bool computeScale) const
 Computes the best rigid transformation between three corresponding pairs. The transformation is characterized by rotation matrix, translation vector and a center about which we rotate. The set of pairs is potentially being updated by the best permutation of the second set. Returns the RMS of the fit. The method is being called with n points but it applies the fit for only 3 after the best permutation is selected in the second set (see bellow). This is done because the solution for planar points is much simpler. The method is the closed-form solution by Horn: people.csail.mit.edu/bkph/papers/Absolute_Orientation.pdf. More...
 
template<typename InputRange1 , typename InputRange2 , template< typename > class Sampler>
void init (const InputRange1 &P, const InputRange2 &Q, const Sampler< PointType > &sampler)
 Initializes the internal state of the Base class. More...
 
virtual void Initialize ()
 Initializes the data structures and needed values before the match computation. This method is called once the internal state of the Base class as been set. More...
 
template<Utils::LogLevel level, typename... Args>
void Log (Args...args) const
 
Scalar MeanDistance () const
 Computes the mean distance between points in Q and their nearest neighbor. We need this for normalization of the user delta (See the paper) to the "scale" of the set. More...
 
bool SelectRandomTriangle (Scalar max_base_diameter, int &base1, int &base2, int &base3)
 Selects a random triangle in the set P (then we add another point to keep the base as planar as possible). We apply a simple heuristic that works in most practical cases. The idea is to accept maximum distance, computed by the estimated overlap, multiplied by the diameter of P, and try to have a triangle with all three edges close to this distance. Wide triangles helps to make the transformation robust while too large triangles makes the probability of having all points in the inliers small so we try to trade-off. More...
 

Protected Attributes

VectorType centroid_P_ {VectorType::Zero()}
 The centroid of P. More...
 
VectorType centroid_Q_ {VectorType::Zero()}
 The centroid of Q. More...
 
KdTree< Scalarkd_tree_
 KdTree used to compute the LCP. More...
 
const Utils::Loggerlogger_
 
OptionsType options_
 
Scalar P_diameter_ {Scalar( 0 )}
 The diameter of P. More...
 
Scalar P_mean_distance_ {Scalar( 0 )}
 Mean distance between points and their nearest neighbor in the set P. Used to normalize the "delta" which is given in terms of this distance. More...
 
VectorType qcentroid1_ {VectorType::Zero()}
 
VectorType qcentroid2_ {VectorType::Zero()}
 
std::mt19937 randomGenerator_
 
std::vector< PosMutablePointsampled_P_3D_
 Sampled P (3D coordinates). More...
 
std::vector< PosMutablePointsampled_Q_3D_
 Sampled Q (3D coordinates). More...
 
Eigen::Matrix< Scalar, 4, 4 > transform_ {Eigen::Matrix<Scalar, 4, 4>::Identity()}
 The transformation matrix by wich we transform Q to P. More...
 

Static Protected Attributes

static constexpr int kNumberOfDiameterTrials = 1000
 

Detailed Description

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
class gr::MatchBase< PointType, _TransformVisitor, OptExts >

Abstract class for registration algorithms.

Member Typedef Documentation

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
using gr::MatchBase< PointType, _TransformVisitor, OptExts >::LogLevel = Utils::LogLevel
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
using gr::MatchBase< PointType, _TransformVisitor, OptExts >::MatrixType = Eigen::Matrix<Scalar, 4, 4>
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
using gr::MatchBase< PointType, _TransformVisitor, OptExts >::OptionsType = gr::Utils::CRTP < OptExts ... , Options >
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
using gr::MatchBase< PointType, _TransformVisitor, OptExts >::Scalar = typename PointType::Scalar
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
using gr::MatchBase< PointType, _TransformVisitor, OptExts >::TransformVisitor = _TransformVisitor
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
using gr::MatchBase< PointType, _TransformVisitor, OptExts >::VectorType = typename PointType::VectorType

Constructor & Destructor Documentation

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
EIGEN_MAKE_ALIGNED_OPERATOR_NEW gr::MatchBase< PointType, _TransformVisitor, OptExts >::MatchBase ( const OptionsType options,
const Utils::Logger logger 
)
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
virtual gr::MatchBase< PointType, _TransformVisitor, OptExts >::~MatchBase ( )
virtual

Member Function Documentation

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
template<typename Coordinates >
bool gr::MatchBase< PointType, _TransformVisitor, OptExts >::ComputeRigidTransformation ( const Coordinates &  ref,
const Coordinates &  candidate,
const Eigen::Matrix< Scalar, 3, 1 > &  centroid1,
Eigen::Matrix< Scalar, 3, 1 >  centroid2,
Eigen::Ref< MatrixType transform,
Scalar rms_,
bool  computeScale 
) const
protected

Computes the best rigid transformation between three corresponding pairs. The transformation is characterized by rotation matrix, translation vector and a center about which we rotate. The set of pairs is potentially being updated by the best permutation of the second set. Returns the RMS of the fit. The method is being called with n points but it applies the fit for only 3 after the best permutation is selected in the second set (see bellow). This is done because the solution for planar points is much simpler. The method is the closed-form solution by Horn: people.csail.mit.edu/bkph/papers/Absolute_Orientation.pdf.

Template Parameters
CoordinatesStruct with operator[](int i) ->i[0:2]
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
template<typename InputRange1 , typename InputRange2 , template< typename > typename Sampler>
Scalar gr::MatchBase< PointType, _TransformVisitor, OptExts >::ComputeTransformation ( const InputRange1 &  P,
const InputRange2 &  Q,
Eigen::Ref< MatrixType transformation,
const Sampler< PointType > &  sampler,
TransformVisitor v 
)
inline

Computes an approximation of the best LCP (directional) from Q to P and the rigid transformation that realizes it. The input sets may or may not contain normal information for any point.

Parameters
[in]PThe first input set.
[in]QThe second input set.
[out]transformationRigid transformation matrix (4x4) that brings Q to the (approximate) optimal LCP. Initial value is considered as a guess
Returns
the computed LCP measure as a fraction of the size of P ([0..1]). Computes an approximation of the best LCP (directional) from Q to P and the rigid transformation that realizes it. The input sets may or may not contain normal information for any point.
Parameters
[in]PThe first input set.
[in]QThe second input set.
[out]transformationRigid transformation matrix (4x4) that brings Q to the (approximate) optimal LCP. Initial value is considered as a guess
Returns
the computed LCP measure as a fraction of the size of P ([0..1]).
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
const std::vector<PosMutablePoint>& gr::MatchBase< PointType, _TransformVisitor, OptExts >::getFirstSampled ( ) const
inline

Read access to the sampled clouds used for the registration.

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
const std::vector<PosMutablePoint>& gr::MatchBase< PointType, _TransformVisitor, OptExts >::getSecondSampled ( ) const
inline

Read access to the sampled clouds used for the registration.

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
template<typename InputRange1 , typename InputRange2 , template< typename > class Sampler>
void gr::MatchBase< PointType, _TransformVisitor, OptExts >::init ( const InputRange1 &  P,
const InputRange2 &  Q,
const Sampler< PointType > &  sampler 
)
protected

Initializes the internal state of the Base class.

Parameters
PThe first input set.
QThe second input set.
samplerThe sampler used to sample the input sets.
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
virtual void gr::MatchBase< PointType, _TransformVisitor, OptExts >::Initialize ( )
inlineprotectedvirtual

Initializes the data structures and needed values before the match computation. This method is called once the internal state of the Base class as been set.

Reimplemented in gr::Match4pcsBase< _Functor, _PointType, _TransformVisitor, _PairFilteringFunctor, PairFilteringOptions >.

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
template<Utils::LogLevel level, typename... Args>
void gr::MatchBase< PointType, _TransformVisitor, OptExts >::Log ( Args...  args) const
inlineprotected
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
Scalar gr::MatchBase< PointType, _TransformVisitor, OptExts >::MeanDistance ( ) const
protected

Computes the mean distance between points in Q and their nearest neighbor. We need this for normalization of the user delta (See the paper) to the "scale" of the set.

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
bool gr::MatchBase< PointType, _TransformVisitor, OptExts >::SelectRandomTriangle ( Scalar  max_base_diameter,
int base1,
int base2,
int base3 
)
protected

Selects a random triangle in the set P (then we add another point to keep the base as planar as possible). We apply a simple heuristic that works in most practical cases. The idea is to accept maximum distance, computed by the estimated overlap, multiplied by the diameter of P, and try to have a triangle with all three edges close to this distance. Wide triangles helps to make the transformation robust while too large triangles makes the probability of having all points in the inliers small so we try to trade-off.

Parameters
max_base_diameterMaximum size allowed between two points of the base

Member Data Documentation

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
VectorType gr::MatchBase< PointType, _TransformVisitor, OptExts >::centroid_P_ {VectorType::Zero()}
protected

The centroid of P.

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
VectorType gr::MatchBase< PointType, _TransformVisitor, OptExts >::centroid_Q_ {VectorType::Zero()}
protected

The centroid of Q.

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
KdTree<Scalar> gr::MatchBase< PointType, _TransformVisitor, OptExts >::kd_tree_
protected

KdTree used to compute the LCP.

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
constexpr int gr::MatchBase< PointType, _TransformVisitor, OptExts >::kNumberOfDiameterTrials = 1000
staticprotected
Todo:
Rationnalize use and name of this variable
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
const Utils::Logger& gr::MatchBase< PointType, _TransformVisitor, OptExts >::logger_
protected
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
OptionsType gr::MatchBase< PointType, _TransformVisitor, OptExts >::options_
protected
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
Scalar gr::MatchBase< PointType, _TransformVisitor, OptExts >::P_diameter_ {Scalar( 0 )}
protected

The diameter of P.

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
Scalar gr::MatchBase< PointType, _TransformVisitor, OptExts >::P_mean_distance_ {Scalar( 0 )}
protected

Mean distance between points and their nearest neighbor in the set P. Used to normalize the "delta" which is given in terms of this distance.

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
VectorType gr::MatchBase< PointType, _TransformVisitor, OptExts >::qcentroid1_ {VectorType::Zero()}
protected
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
VectorType gr::MatchBase< PointType, _TransformVisitor, OptExts >::qcentroid2_ {VectorType::Zero()}
protected
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
std::mt19937 gr::MatchBase< PointType, _TransformVisitor, OptExts >::randomGenerator_
protected
template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
std::vector<PosMutablePoint> gr::MatchBase< PointType, _TransformVisitor, OptExts >::sampled_P_3D_
protected

Sampled P (3D coordinates).

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
std::vector<PosMutablePoint> gr::MatchBase< PointType, _TransformVisitor, OptExts >::sampled_Q_3D_
protected

Sampled Q (3D coordinates).

template<typename PointType, typename _TransformVisitor = DummyTransformVisitor, template< class, class > class... OptExts>
Eigen::Matrix<Scalar, 4, 4> gr::MatchBase< PointType, _TransformVisitor, OptExts >::transform_ {Eigen::Matrix<Scalar, 4, 4>::Identity()}
protected

The transformation matrix by wich we transform Q to P.


The documentation for this class was generated from the following file:
  • /home/travis/build/STORM-IRIT/OpenGR/src/gr/algorithms/matchBase.h