Geometry » Reference » Geometry module module

Contents

This module provides support for:

#include <Eigen/Geometry>

Modules

module Global aligned box typedefs

Classes

template<typename _Scalar, int _AmbientDim>
class Eigen::AlignedBox
An axis aligned box.
template<typename _Scalar>
class Eigen::AngleAxis
Represents a 3D rotation as a rotation angle around an arbitrary 3D axis.
template<typename MatrixType, int _Direction>
class Eigen::Homogeneous
Expression of one (or a set of) homogeneous vector(s)
template<typename _Scalar, int _AmbientDim, int _Options>
class Eigen::Hyperplane
A hyperplane.
template<typename _Scalar, int _Options>
class Eigen::Map<const Quaternion<_Scalar>, _Options>
Quaternion expression mapping a constant memory buffer.
template<typename _Scalar, int _Options>
class Eigen::Map<Quaternion<_Scalar>, _Options>
Expression of a quaternion from a memory buffer.
template<typename _Scalar, int _AmbientDim, int _Options>
class Eigen::ParametrizedLine
A parametrized line.
template<typename _Scalar, int _Options>
class Eigen::Quaternion
The quaternion class used to represent 3D orientations and rotations.
template<class Derived>
class Eigen::QuaternionBase
Base class for quaternion expressions.
template<typename _Scalar>
class Eigen::Rotation2D
Represents a rotation/orientation in a 2 dimensional space.
class Scaling
Represents a generic uniform scaling transformation.
template<typename _Scalar, int _Dim, int _Mode, int _Options>
class Eigen::Transform
Represents an homogeneous transformation in a N dimensional space.
template<typename _Scalar, int _Dim>
class Eigen::Translation
Represents a translation transformation.

Typedefs

using AlignedScaling2d = DiagonalMatrix<double, 2> deprecated
using AlignedScaling2f = DiagonalMatrix<float, 2> deprecated
using AlignedScaling3d = DiagonalMatrix<double, 3> deprecated
using AlignedScaling3f = DiagonalMatrix<float, 3> deprecated
using AngleAxisd = AngleAxis<double>
using AngleAxisf = AngleAxis<float>
using Quaterniond = Quaternion<double>
using Quaternionf = Quaternion<float>
using QuaternionMapAlignedd = Map<Quaternion<double>, Aligned>
using QuaternionMapAlignedf = Map<Quaternion<float>, Aligned>
using QuaternionMapd = Map<Quaternion<double>, 0>
using QuaternionMapf = Map<Quaternion<float>, 0>
using Rotation2Dd = Rotation2D<double>
using Rotation2Df = Rotation2D<float>

Functions

template<typename OtherDerived>
auto cross(const MatrixBase<OtherDerived>& other) const -> PlainObject
template<typename OtherDerived>
auto cross(const MatrixBase<OtherDerived>& other) const -> const CrossReturnType
template<typename OtherDerived>
auto cross3(const MatrixBase<OtherDerived>& other) const -> PlainObject
auto eulerAngles(Index a0, Index a1, Index a2) const -> Matrix<Scalar, 3, 1>
auto hnormalized() const -> const HNormalizedReturnType
homogeneous normalization
auto hnormalized() const -> const HNormalizedReturnType
column or row-wise homogeneous normalization
auto homogeneous() const -> HomogeneousReturnType
auto homogeneous() const -> HomogeneousReturnType
template<typename Derived, typename Scalar>
operator*(const MatrixBase<Derived>& matrix, const UniformScaling<Scalar>& s)
auto Scaling(float s) -> UniformScaling<float>
auto Scaling(double s) -> UniformScaling<double>
template<typename RealScalar>
auto Scaling(const std::complex<RealScalar>& s) -> UniformScaling<std::complex<RealScalar>>
template<typename Scalar>
auto Scaling(const Scalar& sx, const Scalar& sy) -> DiagonalMatrix<Scalar, 2>
template<typename Scalar>
auto Scaling(const Scalar& sx, const Scalar& sy, const Scalar& sz) -> DiagonalMatrix<Scalar, 3>
template<typename Derived>
auto Scaling(const MatrixBase<Derived>& coeffs) -> const DiagonalWrapper<const Derived>
template<typename Derived, typename OtherDerived>
auto umeyama(const MatrixBase<Derived>& src, const MatrixBase<OtherDerived>& dst, bool with_scaling = true) -> internal::umeyama_transform_matrix_type<Derived, OtherDerived>::type
Returns the transformation between two point sets.
auto unitOrthogonal(void) const -> PlainObject

Typedef documentation

typedef DiagonalMatrix<double, 2> AlignedScaling2d

typedef DiagonalMatrix<float, 2> AlignedScaling2f

typedef DiagonalMatrix<double, 3> AlignedScaling3d

typedef DiagonalMatrix<float, 3> AlignedScaling3f

typedef AngleAxis<double> AngleAxisd

double precision angle-axis type

typedef AngleAxis<float> AngleAxisf

single precision angle-axis type

typedef Quaternion<double> Quaterniond

double precision quaternion type

typedef Quaternion<float> Quaternionf

single precision quaternion type

typedef Map<Quaternion<double>, Aligned> QuaternionMapAlignedd

Map a 16-byte aligned array of double precision scalars as a quaternion

typedef Map<Quaternion<float>, Aligned> QuaternionMapAlignedf

Map a 16-byte aligned array of single precision scalars as a quaternion

typedef Map<Quaternion<double>, 0> QuaternionMapd

Map an unaligned array of double precision scalars as a quaternion

typedef Map<Quaternion<float>, 0> QuaternionMapf

Map an unaligned array of single precision scalars as a quaternion

typedef Rotation2D<double> Rotation2Dd

double precision 2D rotation type

typedef Rotation2D<float> Rotation2Df

single precision 2D rotation type

Function documentation

template<typename OtherDerived>
PlainObject cross(const MatrixBase<OtherDerived>& other) const

Returns the cross product of *this and other

This is defined in the Geometry module. #include <Eigen/Geometry>

Here is a very good explanation of cross-product: http://xkcd.com/199/

With complex numbers, the cross product is implemented as $ (\mathbf{a}+i\mathbf{b}) \times (\mathbf{c}+i\mathbf{d}) = (\mathbf{a} \times \mathbf{c} - \mathbf{b} \times \mathbf{d}) - i(\mathbf{a} \times \mathbf{d} - \mathbf{b} \times \mathbf{c})$

template<typename OtherDerived>
const CrossReturnType cross(const MatrixBase<OtherDerived>& other) const

Returns a matrix expression of the cross product of each column or row of the referenced expression with the other vector.

This is defined in the Geometry module. #include <Eigen/Geometry>

The referenced matrix must have one dimension equal to 3. The result matrix has the same dimensions than the referenced one.

template<typename OtherDerived>
PlainObject cross3(const MatrixBase<OtherDerived>& other) const

Returns the cross product of *this and other using only the x, y, and z coefficients

This is defined in the Geometry module. #include <Eigen/Geometry>

The size of *this and other must be four. This function is especially useful when using 4D vectors instead of 3D ones to get advantage of SSE/AltiVec vectorization.

Matrix<Scalar, 3, 1> eulerAngles(Index a0, Index a1, Index a2) const

Returns the Euler-angles of the rotation matrix *this using the convention defined by the triplet (a0,a1,a2)

This is defined in the Geometry module. #include <Eigen/Geometry>

Each of the three parameters a0,a1,a2 represents the respective rotation axis as an integer in {0,1,2}. For instance, in: Vector3f ea = mat.eulerAngles(2, 0, 2); "2" represents the z axis and "0" the x axis, etc. The returned angles are such that we have the following equality:

mat == AngleAxisf(ea[0], Vector3f::UnitZ())
     * AngleAxisf(ea[1], Vector3f::UnitX())
     * AngleAxisf(ea[2], Vector3f::UnitZ());

This corresponds to the right-multiply conventions (with right hand side frames).

The returned angles are in the ranges [0:pi]x[-pi:pi]x[-pi:pi].

const HNormalizedReturnType hnormalized() const

homogeneous normalization

Returns a vector expression of the N-1 first coefficients of *this divided by that last coefficient.

This is defined in the Geometry module. #include <Eigen/Geometry>

This can be used to convert homogeneous coordinates to affine coordinates.

It is essentially a shortcut for: this->head(this->size()-1)/this->coeff(this->size()-1);

Example:

Vector4d v = Vector4d::Random();
Projective3d P(Matrix4d::Random());
cout << "v                   = " << v.transpose() << "]^T" << endl;
cout << "v.hnormalized()     = " << v.hnormalized().transpose() << "]^T" << endl;
cout << "P*v                 = " << (P*v).transpose() << "]^T" << endl;
cout << "(P*v).hnormalized() = " << (P*v).hnormalized().transpose() << "]^T" << endl;

Output:

v                   =   0.68 -0.211  0.566  0.597]^T
v.hnormalized()     =   1.14 -0.354  0.949]^T
P*v                 = 0.663 -0.16 -0.13  0.91]^T
(P*v).hnormalized() =  0.729 -0.176 -0.143]^T

const HNormalizedReturnType hnormalized() const

column or row-wise homogeneous normalization

Returns an expression of the first N-1 coefficients of each column (or row) of *this divided by the last coefficient of each column (or row).

This is defined in the Geometry module. #include <Eigen/Geometry>

This can be used to convert homogeneous coordinates to affine coordinates.

It is conceptually equivalent to calling MatrixBase::hnormalized() to each column (or row) of *this.

Example:

Matrix4Xd M = Matrix4Xd::Random(4,5);
Projective3d P(Matrix4d::Random());
cout << "The matrix M is:" << endl << M << endl << endl;
cout << "M.colwise().hnormalized():" << endl << M.colwise().hnormalized() << endl << endl;
cout << "P*M:" << endl << P*M << endl << endl;
cout << "(P*M).colwise().hnormalized():" << endl << (P*M).colwise().hnormalized() << endl << endl;

Output:

The matrix M is:
   0.68   0.823  -0.444   -0.27   0.271
 -0.211  -0.605   0.108  0.0268   0.435
  0.566   -0.33 -0.0452   0.904  -0.717
  0.597   0.536   0.258   0.832   0.214

M.colwise().hnormalized():
  1.14   1.53  -1.72 -0.325   1.27
-0.354  -1.13  0.419 0.0322   2.03
 0.949 -0.614 -0.175   1.09  -3.35

P*M:
  0.186  -0.589   0.369    1.33   -1.23
 -0.871  -0.337   0.127  -0.715   0.091
 -0.158 -0.0104   0.312   0.429  -0.478
  0.992   0.777  -0.373   0.468  -0.651

(P*M).colwise().hnormalized():
  0.188  -0.759  -0.989    2.85    1.89
 -0.877  -0.433  -0.342   -1.53   -0.14
  -0.16 -0.0134  -0.837   0.915   0.735

HomogeneousReturnType homogeneous() const

Returns a vector expression that is one longer than the vector argument, with the value 1 symbolically appended as the last coefficient.

This is defined in the Geometry module. #include <Eigen/Geometry>

This can be used to convert affine coordinates to homogeneous coordinates.

This is only for vectors (either row-vectors or column-vectors), i.e. matrices which are known at compile-time to have either one row or one column.

Example:

Vector3d v = Vector3d::Random(), w;
Projective3d P(Matrix4d::Random());
cout << "v                                   = [" << v.transpose() << "]^T" << endl;
cout << "h.homogeneous()                     = [" << v.homogeneous().transpose() << "]^T" << endl;
cout << "(P * v.homogeneous())               = [" << (P * v.homogeneous()).transpose() << "]^T" << endl;
cout << "(P * v.homogeneous()).hnormalized() = [" << (P * v.homogeneous()).eval().hnormalized().transpose() << "]^T" << endl;

Output:

v                                   = [  0.68 -0.211  0.566]^T
h.homogeneous()                     = [  0.68 -0.211  0.566      1]^T
(P * v.homogeneous())               = [  1.27  0.772 0.0154 -0.419]^T
(P * v.homogeneous()).hnormalized() = [  -3.03   -1.84 -0.0367]^T

HomogeneousReturnType homogeneous() const

Returns an expression where the value 1 is symbolically appended as the final coefficient to each column (or row) of the matrix.

This is defined in the Geometry module. #include <Eigen/Geometry>

This can be used to convert affine coordinates to homogeneous coordinates.

Example:

Matrix3Xd M = Matrix3Xd::Random(3,5);
Projective3d P(Matrix4d::Random());
cout << "The matrix M is:" << endl << M << endl << endl;
cout << "M.colwise().homogeneous():" << endl << M.colwise().homogeneous() << endl << endl;
cout << "P * M.colwise().homogeneous():" << endl << P * M.colwise().homogeneous() << endl << endl;
cout << "P * M.colwise().homogeneous().hnormalized(): " << endl << (P * M.colwise().homogeneous()).colwise().hnormalized() << endl << endl;

Output:

The matrix M is:
   0.68   0.597   -0.33   0.108   -0.27
 -0.211   0.823   0.536 -0.0452  0.0268
  0.566  -0.605  -0.444   0.258   0.904

M.colwise().homogeneous():
   0.68   0.597   -0.33   0.108   -0.27
 -0.211   0.823   0.536 -0.0452  0.0268
  0.566  -0.605  -0.444   0.258   0.904
      1       1       1       1       1

P * M.colwise().homogeneous():
0.0832 -0.477  -1.21 -0.545 -0.452
 0.998  0.779  0.695  0.894  0.277
-0.271 -0.608 -0.895 -0.544 -0.874
-0.728 -0.551  0.202  -0.21 -0.469

P * M.colwise().homogeneous().hnormalized(): 
-0.114  0.866     -6    2.6  0.962
 -1.37  -1.41   3.44  -4.27 -0.591
 0.373    1.1  -4.43    2.6   1.86

template<typename Derived, typename Scalar>
operator*(const MatrixBase<Derived>& matrix, const UniformScaling<Scalar>& s)

Concatenates a linear transformation matrix and a uniform scaling

UniformScaling<float> Scaling(float s)

Constructs a uniform scaling from scale factor s

UniformScaling<double> Scaling(double s)

Constructs a uniform scaling from scale factor s

template<typename RealScalar>
UniformScaling<std::complex<RealScalar>> Scaling(const std::complex<RealScalar>& s)

Constructs a uniform scaling from scale factor s

template<typename Scalar>
DiagonalMatrix<Scalar, 2> Scaling(const Scalar& sx, const Scalar& sy)

Constructs a 2D axis aligned scaling

template<typename Scalar>
DiagonalMatrix<Scalar, 3> Scaling(const Scalar& sx, const Scalar& sy, const Scalar& sz)

Constructs a 3D axis aligned scaling

template<typename Derived>
const DiagonalWrapper<const Derived> Scaling(const MatrixBase<Derived>& coeffs)

Constructs an axis aligned scaling expression from vector expression coeffs This is an alias for coeffs.asDiagonal()

template<typename Derived, typename OtherDerived>
internal::umeyama_transform_matrix_type<Derived, OtherDerived>::type umeyama(const MatrixBase<Derived>& src, const MatrixBase<OtherDerived>& dst, bool with_scaling = true)

Returns the transformation between two point sets.

Parameters
src Source points $ \mathbf{x} = \left( x_1, \hdots, x_n \right) $ .
dst Destination points $ \mathbf{y} = \left( y_1, \hdots, y_n \right) $ .
with_scaling Sets $ c=1 $ when false is passed.
Returns

The homogeneous transformation

\begin{align*} T = \begin{bmatrix} c\mathbf{R} & \mathbf{t} \\ \mathbf{0} & 1 \end{bmatrix} \end{align*}

minimizing the resudiual above. This transformation is always returned as an Eigen::Matrix.

This is defined in the Geometry module. #include <Eigen/Geometry>

The algorithm is based on: "Least-squares estimation of transformation parameters between two point patterns", Shinji Umeyama, PAMI 1991, DOI: 10.1109/34.88573

It estimates parameters $ c, \mathbf{R}, $ and $ \mathbf{t} $ such that

\begin{align*} \frac{1}{n} \sum_{i=1}^n \vert\vert y_i - (c\mathbf{R}x_i + \mathbf{t}) \vert\vert_2^2 \end{align*}

is minimized.

The algorithm is based on the analysis of the covariance matrix $ \Sigma_{\mathbf{x}\mathbf{y}} \in \mathbb{R}^{d \times d} $ of the input point sets $ \mathbf{x} $ and $ \mathbf{y} $ where $d$ is corresponding to the dimension (which is typically small). The analysis is involving the SVD having a complexity of $O(d^3)$ though the actual computational effort lies in the covariance matrix computation which has an asymptotic lower bound of $O(dm)$ when the input point sets have dimension $d \times m$ .

Currently the method is working only for floating point matrices.

PlainObject unitOrthogonal(void) const

Returns a unit vector which is orthogonal to *this

This is defined in the Geometry module. #include <Eigen/Geometry>

The size of *this must be at least 2. If the size is exactly 2, then the returned vector is a counter clock wise rotation of *this, i.e., (-y,x).normalized().