Public Types | Public Member Functions | Public Attributes | Protected Attributes

TRTK::EstimateProjectiveTransformation2D< T > Class Template Reference

Estimates a 2D projective transformation from two point sets. More...

#include <EstimateProjectiveTransformation2D.hpp>

Inheritance diagram for TRTK::EstimateProjectiveTransformation2D< T >:
Collaboration diagram for TRTK::EstimateProjectiveTransformation2D< T >:

List of all members.

Public Types

typedef super::ArrayXT ArrayXT
 General-purpose array of arbitrary size with value type T.
typedef super::MatrixXT MatrixXT
 Matrix of arbitrary size with value type T.
typedef super::VectorXT VectorXT
 Column vector of arbitrary size with value type T.
typedef super::Vector2T Vector2T
 2D column vector with value type T.
typedef super::Vector3T Vector3T
 3D column vector with value type T.
typedef super::Vector4T Vector4T
 4D column vector with value type T.
typedef super::RowVector2T RowVector2T
 2D row vector with value type T.
typedef super::RowVector3T RowVector3T
 3D row vector with value type T.
typedef super::Matrix2T Matrix2T
 2 x 2 matrix with value type T.
typedef super::Matrix3T Matrix3T
 3 x 3 matrix with value type T.
typedef super::Matrix4T Matrix4T
 4 x 4 matrix with value type T.
enum  Error { NOT_ENOUGH_POINTS, UNEQUAL_NUMBER_OF_POINTS, UNKNOWN_ERROR, WRONG_POINT_SIZE }

Public Member Functions

 EstimateProjectiveTransformation2D ()
 Constructs a new instance of this class.
 EstimateProjectiveTransformation2D (const std::vector< Coordinate< T > > &source_points, const std::vector< Coordinate< T > > &target_points)
 Constructs a new instance of this class.
 EstimateProjectiveTransformation2D (const std::vector< Vector2T, Eigen::aligned_allocator< Vector2T > > &source_points, const std::vector< Vector2T, Eigen::aligned_allocator< Vector2T > > &target_points)
 Constructs a new instance of this class.
 EstimateProjectiveTransformation2D (const std::vector< Vector3T > &source_points, const std::vector< Vector3T > &target_points)
 Constructs a new instance of this class.
virtual ~EstimateProjectiveTransformation2D ()
 Deletes the instance.
void compute ()
 Computes the transformation estimation.
void setMaxIterations (int value)
 Sets the maximum number of iterations within Newton's method.
void setMaxRMS (value_type value)
 Sets the maximum root mean square (RMS) error value.
virtual value_type getRMS () const
 Returns the root mean square (RMS) error of the estimated transformation.
virtual const Matrix3TgetTransformationMatrix () const
void setSourcePoints (const std::vector< Coordinate< T > > &)
 Sets the source points.
void setSourcePoints (const std::vector< Vector2T, Eigen::aligned_allocator< Vector2T > > &)
 Sets the source points.
void setSourcePoints (const std::vector< Vector3T > &)
 Sets the source points.
void setTargetPoints (const std::vector< Coordinate< T > > &)
 Sets the target points.
void setTargetPoints (const std::vector< Vector2T, Eigen::aligned_allocator< Vector2T > > &)
 Sets the target points.
void setTargetPoints (const std::vector< Vector3T > &)
 Sets the target points.

Public Attributes

EIGEN_MAKE_ALIGNED_OPERATOR_NEW
typedef T 
value_type
 Internally used value type (should be a floating point type).

Protected Attributes

MatrixXT m_source_points
MatrixXT m_target_points
Matrix3T m_transformation_matrix

Detailed Description

template<class T>
class TRTK::EstimateProjectiveTransformation2D< T >

Estimates a 2D projective transformation from two point sets.

This class estimates a 2D projective transformation between two point sets by using Newton's method. The point sets must have the same cardinality, and the points must correspond to each other. There must be at least four corresponding point pairs.

The algorithm estimates a transformation matrix as shown below:

\[ \begin{pmatrix} y_1 \\ y_2 \end{pmatrix} = \frac{1}{y'_3} \begin{pmatrix} y'_1 \\ y'_2 \end{pmatrix} \quad \wedge \quad \begin{pmatrix} y'_1 \\ y'_2 \\ y'_3 \end{pmatrix} = \begin{pmatrix} a_{11} & a_{12} & a_{13} \\ a_{21} & a_{22} & a_{23} \\ a_{31} & a_{32} & a_{33} \end{pmatrix} \begin{pmatrix} x_1 \\ x_2 \\ 1 \end{pmatrix} \]

The source points are the set of all \( x \) and the target points are the set of all \( y \).

Note:
  • To obtain stable estimation results, at least six point pairs should be given.
  • Some functions might throw an error object. See the appropriate function for more details.
  • The template parameter T is the internally used value type and must be a floating number type.
See also:
Homogeneous Coordinates

Here is an more elaborate example to see, how to use the class:

 #include <cstdlib>
 #include <iostream>
 #include <vector>

 #include <TRTK/Coordinate.hpp>
 #include <TRTK/Transform2D.hpp>
 #include <TRTK/EstimateProjectiveTransformation2D.hpp>

 using std::cout;
 using std::endl;
 using std::vector;

 using namespace TRTK;


 int main()
 {
     // Construct a projective transformation.

     Transform2D<double> transform;
     transform.a12() = 1;
     transform.a13() = 3;
     transform.a21() = 2;
     transform.a31() = 3;
     transform.a33() = 1;

     // Construct two sets with source and target points, respectively.
     // Add some noise to the target points.

     vector<Coordinate<double> > source_points;
     vector<Coordinate<double> > target_points;

     Coordinate<double> source_point1( 1,  2);
     Coordinate<double> source_point2( 3, -2);
     Coordinate<double> source_point3(-1,  2);
     Coordinate<double> source_point4( 2,  0);
     Coordinate<double> source_point5(-1, -1);
     Coordinate<double> source_point6( 3,  3);

     source_points.push_back(source_point1);
     source_points.push_back(source_point2);
     source_points.push_back(source_point3);
     source_points.push_back(source_point4);
     source_points.push_back(source_point5);
     source_points.push_back(source_point6);

     Coordinate<double> target_point1 = transform * source_point1 + 0.1 * double(std::rand()) / RAND_MAX;
     Coordinate<double> target_point2 = transform * source_point2 + 0.1 * double(std::rand()) / RAND_MAX;
     Coordinate<double> target_point3 = transform * source_point3 + 0.1 * double(std::rand()) / RAND_MAX;
     Coordinate<double> target_point4 = transform * source_point4 + 0.1 * double(std::rand()) / RAND_MAX;
     Coordinate<double> target_point5 = transform * source_point5 + 0.1 * double(std::rand()) / RAND_MAX;
     Coordinate<double> target_point6 = transform * source_point6 + 0.1 * double(std::rand()) / RAND_MAX;

     target_points.push_back(target_point1);
     target_points.push_back(target_point2);
     target_points.push_back(target_point3);
     target_points.push_back(target_point4);
     target_points.push_back(target_point5);
     target_points.push_back(target_point6);

     // Perform the transformation estimation.

     EstimateProjectiveTransformation2D<double> estimateProjectiveTransformation2D(source_points,
                                                                                   target_points);
     estimateProjectiveTransformation2D.compute();

     // Display the results.

     cout.precision(4);
     cout << std::fixed;

     cout << "Original transformation matrix:" << endl << endl
          << transform.getTransformationMatrix() << endl << endl;

     cout << "Estimated transformation matrix:" << endl << endl
          << estimateProjectiveTransformation2D.getTransformationMatrix() << endl << endl;

     // Example of how to use the result.

     Transform2D<double> transform2 = estimateProjectiveTransformation2D.getTransformationMatrix();

     cout << "Source point 1:               " << source_point1 << endl
          << "Target point 1 (original):    " << target_point1 << endl
          << "Target point 1 (transformed): " << transform2 * source_point1 << endl << endl;

     cout << "RMS: " << estimateProjectiveTransformation2D.getRMS() << endl;

     return 0;
 }

Output:

 Original transformation matrix:

 1.0000 1.0000 3.0000
 2.0000 1.0000 0.0000
 3.0000 0.0000 1.0000

 Estimated transformation matrix:

     1.0286     0.7833     2.6236
     1.7938     0.7647     0.2280
     2.5321    -0.0311     1.0000

 Source point 1:               (1.0000, 2.0000)
 Target point 1 (original):    (1.5001, 1.0001)
 Target point 1 (transformed): (1.5039, 1.0234)

 RMS: 0.0458
See also:
EstimateRigidTransformation2D, EstimateProjectiveTransformation2D
Bug:
Depending on the compiler this function may be numerically unstable. It has to be investigated what is the reason for it.
Todo:
Check the estimation algorithm for numerical issues.
Author:
Christoph Haenisch
Version:
0.2.4
Date:
2019-07-12

Definition at line 222 of file EstimateProjectiveTransformation2D.hpp.


Member Enumeration Documentation

template<class T>
enum TRTK::EstimateTransformation::Error [inherited]
Enumerator:
NOT_ENOUGH_POINTS 

More points are required to estimate the transformation.

UNEQUAL_NUMBER_OF_POINTS 

The two point sets do not have the same cardinality.

UNKNOWN_ERROR 

An unknown error occured.

WRONG_POINT_SIZE 

One or more points have a wrong size.

Definition at line 62 of file EstimateTransformation.hpp.


Constructor & Destructor Documentation

Constructs a new instance of this class.

Template Parameters:
Tscalar type of the coordinates

Definition at line 289 of file EstimateProjectiveTransformation2D.hpp.

template<class T >
TRTK::EstimateProjectiveTransformation2D< T >::EstimateProjectiveTransformation2D ( const std::vector< Coordinate< T > > &  source_points,
const std::vector< Coordinate< T > > &  target_points 
)

Constructs a new instance of this class.

Template Parameters:
Tscalar type of the coordinates
Parameters:
[in]source_pointsThe points can either be plain 2D or 3D homogeneous coordinates.
[in]target_pointsThe points can either be plain 2D or 3D homogeneous coordinates.

No transformation estimation is done. To do so, please call compute().

Exceptions:
ErrorObjIf there are any coordinates other than 2D or 3D coordinates, an error object is thrown and its error code is set to WRONG_POINT_SIZE.
Note:
The point sets must have the same cardinality, and the points must correspond to each other.

Definition at line 318 of file EstimateProjectiveTransformation2D.hpp.

template<class T >
TRTK::EstimateProjectiveTransformation2D< T >::EstimateProjectiveTransformation2D ( const std::vector< Vector2T, Eigen::aligned_allocator< Vector2T > > &  source_points,
const std::vector< Vector2T, Eigen::aligned_allocator< Vector2T > > &  target_points 
)

Constructs a new instance of this class.

Template Parameters:
Tscalar type of the coordinates
Parameters:
[in]source_points2D coordinates.
[in]target_points2D coordinates.

No transformation estimation is done. To do so, please call compute().

Note:
The point sets must have the same cardinality, and the points must correspond to each other.

Definition at line 342 of file EstimateProjectiveTransformation2D.hpp.

template<class T >
TRTK::EstimateProjectiveTransformation2D< T >::EstimateProjectiveTransformation2D ( const std::vector< Vector3T > &  source_points,
const std::vector< Vector3T > &  target_points 
)

Constructs a new instance of this class.

Template Parameters:
Tscalar type of the coordinates
Parameters:
[in]source_points3D homogeneous coordinates.
[in]target_points3D homogeneous coordinates.

No transformation estimation is done. To do so, please call compute().

Note:
The point sets must have the same cardinality, and the points must correspond to each other.

Definition at line 366 of file EstimateProjectiveTransformation2D.hpp.

Deletes the instance.

Template Parameters:
Tscalar type of the coordinates

Definition at line 382 of file EstimateProjectiveTransformation2D.hpp.


Member Function Documentation

template<class T >
void TRTK::EstimateProjectiveTransformation2D< T >::compute (  ) [virtual]

Computes the transformation estimation.

Template Parameters:
Tscalar type of the coordinates

This function uses Newton's method (after Simpson) to estimate the sought transformation matrix entries. For a more detailed explanation of the algorithm, please have a look at the source code. (A link should be below.)

The direction of the transformation is from the source points to the target points: \( p_{target} = \cal{T} \; \{ p_{source} \} \).

Exceptions:
ErrorObjIf the point sets do not have the same cardinality, an error object is thrown and its error code is set to UNEQUAL_NUMBER_OF_POINTS.
ErrorObjIf there are not enough points (at least three) to perform the transformation estimation, an error object is thrown and its error code is set to NOT_ENOUGH_POINTS.

Implements TRTK::EstimateTransformation2D< T >.

Definition at line 408 of file EstimateProjectiveTransformation2D.hpp.

template<class T >
EstimateTransformation2D< T >::value_type TRTK::EstimateTransformation2D< T >::getRMS (  ) const [virtual, inherited]

Returns the root mean square (RMS) error of the estimated transformation.

Template Parameters:
Tscalar type of the coordinates

It is assumed, that the computation was done before.

The value type T must provide a function T sqrt(T value) which yields the square root of value.

Exceptions:
ErrorObjIf the point sets do not have the same cardinality, an error object is thrown and its error code is set to UNEQUAL_NUMBER_OF_POINTS.
See also:
compute()

Implements TRTK::EstimateTransformation< T >.

Definition at line 148 of file EstimateTransformation2D.hpp.

template<class T >
const EstimateTransformation2D< T >::Matrix3T & TRTK::EstimateTransformation2D< T >::getTransformationMatrix (  ) const [virtual, inherited]
Template Parameters:
Tscalar type of the coordinates
Returns:
Returns the sought transformation matrix in the form of a homogeneous 3x3 matrix. (This comprises the rotation as well as the translation, for instance.)
See also:
compute()

Definition at line 194 of file EstimateTransformation2D.hpp.

template<class T >
void TRTK::EstimateProjectiveTransformation2D< T >::setMaxIterations ( int  value )

Sets the maximum number of iterations within Newton's method.

Template Parameters:
Tscalar type of the coordinates

Sets the maximum number of iterations within Newton's method which is used to compute the transformation matrix.

The default value is 0.1.

See also:
setMaxRMS()

Definition at line 601 of file EstimateProjectiveTransformation2D.hpp.

template<class T >
void TRTK::EstimateProjectiveTransformation2D< T >::setMaxRMS ( value_type  value )

Sets the maximum root mean square (RMS) error value.

Template Parameters:
Tscalar type of the coordinates

Sets the maximum root mean square (RMS) error value within Newton's method which is used to compute the transformation matrix. That is, the iteration is stopped, if the RMS is below value.

The default value is 10.

See also:
setMaxIterations()

Definition at line 621 of file EstimateProjectiveTransformation2D.hpp.

template<class T >
void TRTK::EstimateTransformation2D< T >::setSourcePoints ( const std::vector< Vector2T, Eigen::aligned_allocator< Vector2T > > &  source_points ) [inherited]

Sets the source points.

Template Parameters:
Tscalar type of the coordinates
Parameters:
[in]source_points2D coordinates.

Definition at line 246 of file EstimateTransformation2D.hpp.

template<class T >
void TRTK::EstimateTransformation2D< T >::setSourcePoints ( const std::vector< Vector3T > &  source_points ) [inherited]

Sets the source points.

Template Parameters:
Tscalar type of the coordinates
Parameters:
[in]source_points3D homogeneous coordinates.

Definition at line 264 of file EstimateTransformation2D.hpp.

template<class T >
void TRTK::EstimateTransformation2D< T >::setSourcePoints ( const std::vector< Coordinate< T > > &  source_points ) [inherited]

Sets the source points.

Template Parameters:
Tscalar type of the coordinates
Parameters:
[in]source_pointsThe points can either be plain 2D or 3D homogeneous coordinates.
Exceptions:
ErrorObjIf there are any coordinates other than 2D or 3D coordinates, an error object is thrown and its error code is set to WRONG_POINT_SIZE.

Definition at line 214 of file EstimateTransformation2D.hpp.

template<class T >
void TRTK::EstimateTransformation2D< T >::setTargetPoints ( const std::vector< Coordinate< T > > &  target_points ) [inherited]

Sets the target points.

Template Parameters:
Tscalar type of the coordinates
Parameters:
[in]target_pointsThe points can either be plain 2D or 3D homogeneous coordinates.
Exceptions:
ErrorObjIf there are any coordinates other than 2D or 3D coordinates, an error object is thrown and its error code is set to WRONG_POINT_SIZE.

Definition at line 292 of file EstimateTransformation2D.hpp.

template<class T >
void TRTK::EstimateTransformation2D< T >::setTargetPoints ( const std::vector< Vector3T > &  target_points ) [inherited]

Sets the target points.

Template Parameters:
Tscalar type of the coordinates
Parameters:
[in]target_points3D homogeneous coordinates.

Definition at line 342 of file EstimateTransformation2D.hpp.

template<class T >
void TRTK::EstimateTransformation2D< T >::setTargetPoints ( const std::vector< Vector2T, Eigen::aligned_allocator< Vector2T > > &  target_points ) [inherited]

Sets the target points.

Template Parameters:
Tscalar type of the coordinates
Parameters:
[in]target_points2D coordinates.

Definition at line 324 of file EstimateTransformation2D.hpp.


The documentation for this class was generated from the following file:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines