Geometry_Module

Modules

 Global aligned box typedefs
 

Classes

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

Typedefs

typedef AngleAxis< float > Eigen::AngleAxisf
 
typedef AngleAxis< doubleEigen::AngleAxisd
 
typedef Quaternion< float > Eigen::Quaternionf
 
typedef Quaternion< doubleEigen::Quaterniond
 
typedef Map< Quaternion< float >, 0 > Eigen::QuaternionMapf
 
typedef Map< Quaternion< double >, 0 > Eigen::QuaternionMapd
 
typedef Map< Quaternion< float >, AlignedEigen::QuaternionMapAlignedf
 
typedef Map< Quaternion< double >, AlignedEigen::QuaternionMapAlignedd
 
typedef Rotation2D< float > Eigen::Rotation2Df
 
typedef Rotation2D< doubleEigen::Rotation2Dd
 
typedef Transform< float, 2, IsometryEigen::Isometry2f
 
typedef Transform< float, 3, IsometryEigen::Isometry3f
 
typedef Transform< double, 2, IsometryEigen::Isometry2d
 
typedef Transform< double, 3, IsometryEigen::Isometry3d
 
typedef Transform< float, 2, AffineEigen::Affine2f
 
typedef Transform< float, 3, AffineEigen::Affine3f
 
typedef Transform< double, 2, AffineEigen::Affine2d
 
typedef Transform< double, 3, AffineEigen::Affine3d
 
typedef Transform< float, 2, AffineCompactEigen::AffineCompact2f
 
typedef Transform< float, 3, AffineCompactEigen::AffineCompact3f
 
typedef Transform< double, 2, AffineCompactEigen::AffineCompact2d
 
typedef Transform< double, 3, AffineCompactEigen::AffineCompact3d
 
typedef Transform< float, 2, ProjectiveEigen::Projective2f
 
typedef Transform< float, 3, ProjectiveEigen::Projective3f
 
typedef Transform< double, 2, ProjectiveEigen::Projective2d
 
typedef Transform< double, 3, ProjectiveEigen::Projective3d
 

Functions

template<typename Derived , typename OtherDerived >
internal::umeyama_transform_matrix_type< Derived, OtherDerived >::type Eigen::umeyama (const MatrixBase< Derived > &src, const MatrixBase< OtherDerived > &dst, bool with_scaling=true)
 Returns the transformation between two point sets. More...
 
EIGEN_DEVICE_FUNC Matrix< Scalar, 3, 1 > Eigen::MatrixBase< Derived >::canonicalEulerAngles (Index a0, Index a1, Index a2) const
 
EIGEN_DEPRECATED EIGEN_DEVICE_FUNC Matrix< Scalar, 3, 1 > Eigen::MatrixBase< Derived >::eulerAngles (Index a0, Index a1, Index a2) const
 
EIGEN_DEVICE_FUNC HomogeneousReturnType Eigen::MatrixBase< Derived >::homogeneous () const
 
EIGEN_DEVICE_FUNC HomogeneousReturnType Eigen::VectorwiseOp< ExpressionType, Direction >::homogeneous () const
 
EIGEN_DEVICE_FUNC const HNormalizedReturnType Eigen::MatrixBase< Derived >::hnormalized () const
 homogeneous normalization More...
 
EIGEN_DEVICE_FUNC const HNormalizedReturnType Eigen::VectorwiseOp< ExpressionType, Direction >::hnormalized () const
 column or row-wise homogeneous normalization More...
 
template<typename OtherDerived >
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE internal::cross_impl< Derived, OtherDerived >::return_type Eigen::MatrixBase< Derived >::cross (const MatrixBase< OtherDerived > &other) const
 
template<typename OtherDerived >
EIGEN_DEVICE_FUNC PlainObject Eigen::MatrixBase< Derived >::cross3 (const MatrixBase< OtherDerived > &other) const
 
template<typename OtherDerived >
EIGEN_DEVICE_FUNC const CrossReturnType Eigen::VectorwiseOp< ExpressionType, Direction >::cross (const MatrixBase< OtherDerived > &other) const
 
EIGEN_DEVICE_FUNC PlainObject Eigen::MatrixBase< Derived >::unitOrthogonal (void) const
 

Detailed Description

Typedef Documentation

◆ Affine2d

◆ Affine2f

typedef Transform<float, 2, Affine> Eigen::Affine2f

◆ Affine3d

◆ Affine3f

typedef Transform<float, 3, Affine> Eigen::Affine3f

◆ AffineCompact2d

◆ AffineCompact2f

◆ AffineCompact3d

◆ AffineCompact3f

◆ AngleAxisd

double precision angle-axis type

◆ AngleAxisf

typedef AngleAxis<float> Eigen::AngleAxisf

single precision angle-axis type

◆ Isometry2d

◆ Isometry2f

typedef Transform<float, 2, Isometry> Eigen::Isometry2f

◆ Isometry3d

◆ Isometry3f

typedef Transform<float, 3, Isometry> Eigen::Isometry3f

◆ Projective2d

◆ Projective2f

◆ Projective3d

◆ Projective3f

◆ Quaterniond

double precision quaternion type

◆ Quaternionf

single precision quaternion type

◆ QuaternionMapAlignedd

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

◆ QuaternionMapAlignedf

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

◆ QuaternionMapd

Map an unaligned array of double precision scalars as a quaternion

◆ QuaternionMapf

typedef Map<Quaternion<float>, 0> Eigen::QuaternionMapf

Map an unaligned array of single precision scalars as a quaternion

◆ Rotation2Dd

double precision 2D rotation type

◆ Rotation2Df

single precision 2D rotation type

Function Documentation

◆ canonicalEulerAngles()

template<typename Derived >
EIGEN_DEVICE_FUNC Matrix< typename MatrixBase< Derived >::Scalar, 3, 1 > Eigen::MatrixBase< Derived >::canonicalEulerAngles ( Index  a0,
Index  a1,
Index  a2 
) const
inline

\geometry_module

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

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());
AngleAxis< float > AngleAxisf
Definition: AngleAxis.h:165

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

For Tait-Bryan angle configurations (a0 != a2), the returned angles are in the ranges [-pi:pi]x[-pi/2:pi/2]x[-pi:pi]. For proper Euler angle configurations (a0 == a2), the returned angles are in the ranges [-pi:pi]x[0:pi]x[-pi:pi].

The approach used is also described here: https://d3cw3dd2w32x2b.cloudfront.net/wp-content/uploads/2012/07/euler-angles.pdf

See also
class AngleAxis
46  {
47  /* Implemented from Graphics Gems IV */
49 
50  Matrix<Scalar, 3, 1> res;
51 
52  const Index odd = ((a0 + 1) % 3 == a1) ? 0 : 1;
53  const Index i = a0;
54  const Index j = (a0 + 1 + odd) % 3;
55  const Index k = (a0 + 2 - odd) % 3;
56 
57  if (a0 == a2) {
58  // Proper Euler angles (same first and last axis).
59  // The i, j, k indices enable addressing the input matrix as the XYX archetype matrix (see Graphics Gems IV),
60  // where e.g. coeff(k, i) means third column, first row in the XYX archetype matrix:
61  // c2 s2s1 s2c1
62  // s2s3 -c2s1s3 + c1c3 -c2c1s3 - s1c3
63  // -s2c3 c2s1c3 + c1s3 c2c1c3 - s1s3
64 
65  // Note: s2 is always positive.
66  Scalar s2 = numext::hypot(coeff(j, i), coeff(k, i));
67  if (odd) {
68  res[0] = numext::atan2(coeff(j, i), coeff(k, i));
69  // s2 is always positive, so res[1] will be within the canonical [0, pi] range
70  res[1] = numext::atan2(s2, coeff(i, i));
71  } else {
72  // In the !odd case, signs of all three angles are flipped at the very end. To keep the solution within the
73  // canonical range, we flip the solution and make res[1] always negative here (since s2 is always positive,
74  // -atan2(s2, c2) will always be negative). The final flip at the end due to !odd will thus make res[1] positive
75  // and canonical. NB: in the general case, there are two correct solutions, but only one is canonical. For proper
76  // Euler angles, flipping from one solution to the other involves flipping the sign of the second angle res[1] and
77  // adding/subtracting pi to the first and third angles. The addition/subtraction of pi to the first angle res[0]
78  // is handled here by flipping the signs of arguments to atan2, while the calculation of the third angle does not
79  // need special adjustment since it uses the adjusted res[0] as the input and produces a correct result.
80  res[0] = numext::atan2(-coeff(j, i), -coeff(k, i));
81  res[1] = -numext::atan2(s2, coeff(i, i));
82  }
83 
84  // With a=(0,1,0), we have i=0; j=1; k=2, and after computing the first two angles,
85  // we can compute their respective rotation, and apply its inverse to M. Since the result must
86  // be a rotation around x, we have:
87  //
88  // c2 s1.s2 c1.s2 1 0 0
89  // 0 c1 -s1 * M = 0 c3 s3
90  // -s2 s1.c2 c1.c2 0 -s3 c3
91  //
92  // Thus: m11.c1 - m21.s1 = c3 & m12.c1 - m22.s1 = s3
93 
94  Scalar s1 = numext::sin(res[0]);
95  Scalar c1 = numext::cos(res[0]);
96  res[2] = numext::atan2(c1 * coeff(j, k) - s1 * coeff(k, k), c1 * coeff(j, j) - s1 * coeff(k, j));
97  } else {
98  // Tait-Bryan angles (all three axes are different; typically used for yaw-pitch-roll calculations).
99  // The i, j, k indices enable addressing the input matrix as the XYZ archetype matrix (see Graphics Gems IV),
100  // where e.g. coeff(k, i) means third column, first row in the XYZ archetype matrix:
101  // c2c3 s2s1c3 - c1s3 s2c1c3 + s1s3
102  // c2s3 s2s1s3 + c1c3 s2c1s3 - s1c3
103  // -s2 c2s1 c2c1
104 
105  res[0] = numext::atan2(coeff(j, k), coeff(k, k));
106 
107  Scalar c2 = numext::hypot(coeff(i, i), coeff(i, j));
108  // c2 is always positive, so the following atan2 will always return a result in the correct canonical middle angle
109  // range [-pi/2, pi/2]
110  res[1] = numext::atan2(-coeff(i, k), c2);
111 
112  Scalar s1 = numext::sin(res[0]);
113  Scalar c1 = numext::cos(res[0]);
114  res[2] = numext::atan2(s1 * coeff(k, i) - c1 * coeff(j, i), c1 * coeff(j, j) - s1 * coeff(k, j));
115  }
116  if (!odd) {
117  res = -res;
118  }
119 
120  return res;
121 }
int i
Definition: BiCGSTAB_step_by_step.cpp:9
cout<< "Here is the matrix m:"<< endl<< m<< endl;Matrix< ptrdiff_t, 3, 1 > res
Definition: PartialRedux_count.cpp:3
#define EIGEN_STATIC_ASSERT_MATRIX_SPECIFIC_SIZE(TYPE, ROWS, COLS)
Definition: StaticAssert.h:55
SCALAR Scalar
Definition: bench_gemm.cpp:45
if(UPLO(*uplo)==INVALID) info
Definition: level3_impl.h:428
char char char int int * k
Definition: level2_impl.h:374
EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE T cos(const T &x)
Definition: MathFunctions.h:1559
EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE T sin(const T &x)
Definition: MathFunctions.h:1581
EIGEN_DEVICE_FUNC EIGEN_ALWAYS_INLINE T atan2(const T &y, const T &x)
Definition: MathFunctions.h:1689
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:83
std::ptrdiff_t j
Definition: tut_arithmetic_redux_minmax.cpp:2

References Eigen::numext::atan2(), Eigen::numext::cos(), EIGEN_STATIC_ASSERT_MATRIX_SPECIFIC_SIZE, i, j, k, res, and Eigen::numext::sin().

◆ cross() [1/2]

template<typename ExpressionType , int Direction>
template<typename OtherDerived >
EIGEN_DEVICE_FUNC const VectorwiseOp< ExpressionType, Direction >::CrossReturnType Eigen::VectorwiseOp< ExpressionType, Direction >::cross ( const MatrixBase< OtherDerived > &  other) const

\geometry_module

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

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

See also
MatrixBase::cross()
153  {
156  (internal::is_same<Scalar, typename OtherDerived::Scalar>::value),
157  YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY)
158 
159  typename internal::nested_eval<ExpressionType, 2>::type mat(_expression());
160  typename internal::nested_eval<OtherDerived, 2>::type vec(other.derived());
161 
163  if (Direction == Vertical) {
164  eigen_assert(CrossReturnType::RowsAtCompileTime == 3 && "the matrix must have exactly 3 rows");
165  res.row(0) = (mat.row(1) * vec.coeff(2) - mat.row(2) * vec.coeff(1)).conjugate();
166  res.row(1) = (mat.row(2) * vec.coeff(0) - mat.row(0) * vec.coeff(2)).conjugate();
167  res.row(2) = (mat.row(0) * vec.coeff(1) - mat.row(1) * vec.coeff(0)).conjugate();
168  } else {
169  eigen_assert(CrossReturnType::ColsAtCompileTime == 3 && "the matrix must have exactly 3 columns");
170  res.col(0) = (mat.col(1) * vec.coeff(2) - mat.col(2) * vec.coeff(1)).conjugate();
171  res.col(1) = (mat.col(2) * vec.coeff(0) - mat.col(0) * vec.coeff(2)).conjugate();
172  res.col(2) = (mat.col(0) * vec.coeff(1) - mat.col(1) * vec.coeff(0)).conjugate();
173  }
174  return res;
175 }
Direction
An enum that indicates the direction in Cartesian coordinates.
Definition: GeneralDefine.h:56
#define eigen_assert(x)
Definition: Macros.h:910
#define EIGEN_STATIC_ASSERT(X, MSG)
Definition: StaticAssert.h:26
#define EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(TYPE, SIZE)
Definition: StaticAssert.h:50
int rows
Definition: Tutorial_commainit_02.cpp:1
int cols
Definition: Tutorial_commainit_02.cpp:1
Scalar coeff(Index row, Index col) const
Definition: SparseMatrix.h:211
ExpressionType::PlainObject CrossReturnType
Definition: VectorwiseOp.h:653
EIGEN_DEVICE_FUNC const ExpressionType & _expression() const
Definition: VectorwiseOp.h:256
@ Vertical
Definition: Constants.h:266
squared absolute value
Definition: GlobalFunctions.h:87
type
Definition: compute_granudrum_aor.py:141
Definition: Eigen_Colamd.h:49

References cols, eigen_assert, EIGEN_STATIC_ASSERT, EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE, res, rows, and Eigen::Vertical.

◆ cross() [2/2]

template<typename Derived >
template<typename OtherDerived >
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE internal::cross_impl<Derived, OtherDerived>::return_type Eigen::MatrixBase< Derived >::cross ( const MatrixBase< OtherDerived > &  other) const

\geometry_module

Returns
the cross product of *this and other. This is either a scalar for size-2 vectors or a size-3 vector for size-3 vectors.

This method is implemented for two different cases: between vectors of fixed size 2 and between vectors of fixed size 3.

For vectors of size 3, the output is simply the traditional cross product.

For vectors of size 2, the output is a scalar. Given vectors \( v = \begin{bmatrix} v_1 & v_2 \end{bmatrix} \) and \( w = \begin{bmatrix} w_1 & w_2 \end{bmatrix} \), the result is simply \( v\times w = \overline{v_1 w_2 - v_2 w_1} = \text{conj}\left|\begin{smallmatrix} v_1 & w_1 \\ v_2 & w_2 \end{smallmatrix}\right| \); or, to put it differently, it is the third coordinate of the cross product of \( \begin{bmatrix} v_1 & v_2 & v_3 \end{bmatrix} \) and \( \begin{bmatrix} w_1 & w_2 & w_3 \end{bmatrix} \). For real-valued inputs, the result can be interpreted as the signed area of a parallelogram spanned by the two vectors.

Note
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}).\]

This definition preserves the orthogonality condition that \(\mathbf{u} \cdot (\mathbf{u} \times \mathbf{v}) = \mathbf{v} \cdot (\mathbf{u} \times \mathbf{v}) = 0\).
See also
MatrixBase::cross3()
96  {
98 }
static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE return_type run(const MatrixBase< Derived > &first, const MatrixBase< OtherDerived > &second)
Definition: OrthoMethods.h:28

References Eigen::internal::cross_impl< Derived, OtherDerived, Size >::run().

◆ cross3()

template<typename Derived >
template<typename OtherDerived >
EIGEN_DEVICE_FUNC MatrixBase< Derived >::PlainObject Eigen::MatrixBase< Derived >::cross3 ( const MatrixBase< OtherDerived > &  other) const
inline

\geometry_module

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

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.

See also
MatrixBase::cross()
128  {
131 
132  typedef typename internal::nested_eval<Derived, 2>::type DerivedNested;
133  typedef typename internal::nested_eval<OtherDerived, 2>::type OtherDerivedNested;
134  DerivedNested lhs(derived());
135  OtherDerivedNested rhs(other.derived());
136 
137  return internal::cross3_impl<Architecture::Target, internal::remove_all_t<DerivedNested>,
138  internal::remove_all_t<OtherDerivedNested>>::run(lhs, rhs);
139 }
if n return
Definition: level1_cplx_impl.h:31
@ Target
Definition: Constants.h:495
typename remove_all< T >::type remove_all_t
Definition: Meta.h:142
auto run(Kernel kernel, Args &&... args) -> decltype(kernel(args...))
Definition: gpu_test_helper.h:414

References EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE, and Eigen::run().

◆ eulerAngles()

template<typename Derived >
EIGEN_DEPRECATED EIGEN_DEVICE_FUNC Matrix< typename MatrixBase< Derived >::Scalar, 3, 1 > Eigen::MatrixBase< Derived >::eulerAngles ( Index  a0,
Index  a1,
Index  a2 
) const
inline

\geometry_module

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

NB: The returned angles are in non-canonical ranges [0:pi]x[-pi:pi]x[-pi:pi]. For canonical Tait-Bryan/proper Euler ranges, use canonicalEulerAngles.

See also
MatrixBase::canonicalEulerAngles
class AngleAxis
137  {
138  /* Implemented from Graphics Gems IV */
140 
141  Matrix<Scalar, 3, 1> res;
142 
143  const Index odd = ((a0 + 1) % 3 == a1) ? 0 : 1;
144  const Index i = a0;
145  const Index j = (a0 + 1 + odd) % 3;
146  const Index k = (a0 + 2 - odd) % 3;
147 
148  if (a0 == a2) {
149  res[0] = numext::atan2(coeff(j, i), coeff(k, i));
150  if ((odd && res[0] < Scalar(0)) || ((!odd) && res[0] > Scalar(0))) {
151  if (res[0] > Scalar(0)) {
152  res[0] -= Scalar(EIGEN_PI);
153  } else {
154  res[0] += Scalar(EIGEN_PI);
155  }
156 
157  Scalar s2 = numext::hypot(coeff(j, i), coeff(k, i));
158  res[1] = -numext::atan2(s2, coeff(i, i));
159  } else {
160  Scalar s2 = numext::hypot(coeff(j, i), coeff(k, i));
161  res[1] = numext::atan2(s2, coeff(i, i));
162  }
163 
164  // With a=(0,1,0), we have i=0; j=1; k=2, and after computing the first two angles,
165  // we can compute their respective rotation, and apply its inverse to M. Since the result must
166  // be a rotation around x, we have:
167  //
168  // c2 s1.s2 c1.s2 1 0 0
169  // 0 c1 -s1 * M = 0 c3 s3
170  // -s2 s1.c2 c1.c2 0 -s3 c3
171  //
172  // Thus: m11.c1 - m21.s1 = c3 & m12.c1 - m22.s1 = s3
173 
174  Scalar s1 = numext::sin(res[0]);
175  Scalar c1 = numext::cos(res[0]);
176  res[2] = numext::atan2(c1 * coeff(j, k) - s1 * coeff(k, k), c1 * coeff(j, j) - s1 * coeff(k, j));
177  } else {
178  res[0] = numext::atan2(coeff(j, k), coeff(k, k));
179  Scalar c2 = numext::hypot(coeff(i, i), coeff(i, j));
180  if ((odd && res[0] < Scalar(0)) || ((!odd) && res[0] > Scalar(0))) {
181  if (res[0] > Scalar(0)) {
182  res[0] -= Scalar(EIGEN_PI);
183  } else {
184  res[0] += Scalar(EIGEN_PI);
185  }
186  res[1] = numext::atan2(-coeff(i, k), -c2);
187  } else {
188  res[1] = numext::atan2(-coeff(i, k), c2);
189  }
190  Scalar s1 = numext::sin(res[0]);
191  Scalar c1 = numext::cos(res[0]);
192  res[2] = numext::atan2(s1 * coeff(k, i) - c1 * coeff(j, i), c1 * coeff(j, j) - s1 * coeff(k, j));
193  }
194  if (!odd) {
195  res = -res;
196  }
197 
198  return res;
199 }
#define EIGEN_PI
Definition: MathFunctions.h:16
internal::traits< Derived >::Scalar Scalar
Definition: MatrixBase.h:58

References Eigen::numext::atan2(), Eigen::numext::cos(), EIGEN_PI, EIGEN_STATIC_ASSERT_MATRIX_SPECIFIC_SIZE, i, j, k, res, and Eigen::numext::sin().

◆ hnormalized() [1/2]

template<typename Derived >
EIGEN_DEVICE_FUNC const MatrixBase< Derived >::HNormalizedReturnType Eigen::MatrixBase< Derived >::hnormalized
inline

homogeneous normalization

\geometry_module

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

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);
Scalar Scalar int size
Definition: benchVecAdd.cpp:17

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;
Array< int, Dynamic, 1 > v
Definition: Array_initializer_list_vector_cxx11.cpp:1
EIGEN_DEVICE_FUNC TransposeReturnType transpose()
Definition: Transpose.h:162
Transform< double, 3, Projective > Projective3d
Definition: Transform.h:704
EIGEN_DEVICE_FUNC const HNormalizedReturnType hnormalized() const
homogeneous normalization
Definition: Homogeneous.h:166
double P
Uniform pressure.
Definition: TwenteMeshGluing.cpp:77

Output:

See also
VectorwiseOp::hnormalized()
167  {
169  return ConstStartMinusOne(derived(), 0, 0, ColsAtCompileTime == 1 ? size() - 1 : 1,
170  ColsAtCompileTime == 1 ? 1 : size() - 1) /
171  coeff(size() - 1);
172 }
#define EIGEN_STATIC_ASSERT_VECTOR_ONLY(TYPE)
Definition: StaticAssert.h:36
@ ColsAtCompileTime
Definition: DenseBase.h:102
Block< const Derived, internal::traits< Derived >::ColsAtCompileTime==1 ? SizeMinusOne :1, internal::traits< Derived >::ColsAtCompileTime==1 ? 1 :SizeMinusOne > ConstStartMinusOne
Definition: MatrixBase.h:412

References EIGEN_STATIC_ASSERT_VECTOR_ONLY, and size.

◆ hnormalized() [2/2]

template<typename ExpressionType , int Direction>
EIGEN_DEVICE_FUNC const VectorwiseOp< ExpressionType, Direction >::HNormalizedReturnType Eigen::VectorwiseOp< ExpressionType, Direction >::hnormalized
inline

column or row-wise homogeneous normalization

\geometry_module

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 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;
Matrix< RealScalar, Dynamic, Dynamic > M
Definition: bench_gemm.cpp:50
The matrix class, also used for vectors and row-vectors.
Definition: Eigen/Eigen/src/Core/Matrix.h:186

Output:

See also
MatrixBase::hnormalized()
191  {
192  return HNormalized_Block(_expression(), 0, 0, Direction == Vertical ? _expression().rows() - 1 : _expression().rows(),
194  .cwiseQuotient(Replicate < HNormalized_Factors, Direction == Vertical ? HNormalized_SizeMinusOne : 1,
198  Direction == Horizontal ? _expression().cols() - 1 : 0,
199  Direction == Vertical ? 1 : _expression().rows(),
200  Direction == Horizontal ? 1 : _expression().cols()),
201  Direction == Vertical ? _expression().rows() - 1 : 1,
202  Direction == Horizontal ? _expression().cols() - 1 : 1));
203 }
@ HNormalized_SizeMinusOne
Definition: VectorwiseOp.h:660
Block< const ExpressionType, Direction==Vertical ? 1 :int(internal::traits< ExpressionType >::RowsAtCompileTime), Direction==Horizontal ? 1 :int(internal::traits< ExpressionType >::ColsAtCompileTime)> HNormalized_Factors
Definition: VectorwiseOp.h:671
Block< const ExpressionType, Direction==Vertical ? int(HNormalized_SizeMinusOne) :int(internal::traits< ExpressionType >::RowsAtCompileTime), Direction==Horizontal ? int(HNormalized_SizeMinusOne) :int(internal::traits< ExpressionType >::ColsAtCompileTime)> HNormalized_Block
Definition: VectorwiseOp.h:667
@ Horizontal
Definition: Constants.h:269

References cols, Eigen::Horizontal, rows, and Eigen::Vertical.

◆ homogeneous() [1/2]

template<typename Derived >
EIGEN_DEVICE_FUNC MatrixBase< Derived >::HomogeneousReturnType Eigen::MatrixBase< Derived >::homogeneous
inline

\geometry_module

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

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

\only_for_vectors

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;
RowVector3d w
Definition: Matrix_resize_int.cpp:3
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE EvalReturnType eval() const
Definition: DenseBase.h:381

Output:

See also
VectorwiseOp::homogeneous(), class Homogeneous
126  {
128  return HomogeneousReturnType(derived());
129 }
Homogeneous< Derived, HomogeneousReturnTypeDirection > HomogeneousReturnType
Definition: MatrixBase.h:406

References EIGEN_STATIC_ASSERT_VECTOR_ONLY.

◆ homogeneous() [2/2]

template<typename ExpressionType , int Direction>
EIGEN_DEVICE_FUNC Homogeneous< ExpressionType, Direction > Eigen::VectorwiseOp< ExpressionType, Direction >::homogeneous
inline

\geometry_module

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

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:

See also
MatrixBase::homogeneous(), class Homogeneous
144  {
146 }
Homogeneous< ExpressionType, Direction > HomogeneousReturnType
Definition: VectorwiseOp.h:650

◆ umeyama()

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

Returns the transformation between two point sets.

\geometry_module

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.

Todo:
Should the return type of umeyama() become a Transform?
Parameters
srcSource points \( \mathbf{x} = \left( x_1, \hdots, x_n \right) \).
dstDestination points \( \mathbf{y} = \left( y_1, \hdots, y_n \right) \).
with_scalingSets \( 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 residual above. This transformation is always returned as an Eigen::Matrix.
95  {
96  typedef typename internal::umeyama_transform_matrix_type<Derived, OtherDerived>::type TransformationMatrixType;
98  typedef typename NumTraits<Scalar>::Real RealScalar;
99 
100  EIGEN_STATIC_ASSERT(!NumTraits<Scalar>::IsComplex, NUMERIC_TYPE_MUST_BE_REAL)
102  (internal::is_same<Scalar, typename internal::traits<OtherDerived>::Scalar>::value),
103  YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY)
104 
105  enum { Dimension = internal::min_size_prefer_dynamic(Derived::RowsAtCompileTime, OtherDerived::RowsAtCompileTime) };
106 
107  typedef Matrix<Scalar, Dimension, 1> VectorType;
108  typedef Matrix<Scalar, Dimension, Dimension> MatrixType;
109  typedef typename internal::plain_matrix_type_row_major<Derived>::type RowMajorMatrixType;
110 
111  const Index m = src.rows(); // dimension
112  const Index n = src.cols(); // number of measurements
113 
114  // required for demeaning ...
115  const RealScalar one_over_n = RealScalar(1) / static_cast<RealScalar>(n);
116 
117  // computation of mean
118  const VectorType src_mean = src.rowwise().sum() * one_over_n;
119  const VectorType dst_mean = dst.rowwise().sum() * one_over_n;
120 
121  // demeaning of src and dst points
122  const RowMajorMatrixType src_demean = src.colwise() - src_mean;
123  const RowMajorMatrixType dst_demean = dst.colwise() - dst_mean;
124 
125  // Eq. (38)
126  const MatrixType sigma = one_over_n * dst_demean * src_demean.transpose();
127 
128  JacobiSVD<MatrixType, ComputeFullU | ComputeFullV> svd(sigma);
129 
130  // Initialize the resulting transformation with an identity matrix...
131  TransformationMatrixType Rt = TransformationMatrixType::Identity(m + 1, m + 1);
132 
133  // Eq. (39)
134  VectorType S = VectorType::Ones(m);
135 
136  if (svd.matrixU().determinant() * svd.matrixV().determinant() < 0) {
137  Index tmp = m - 1;
138  S(tmp) = -1;
139  }
140 
141  // Eq. (40) and (43)
142  Rt.block(0, 0, m, m).noalias() = svd.matrixU() * S.asDiagonal() * svd.matrixV().transpose();
143 
144  if (with_scaling) {
145  // Eq. (36)-(37)
146  const Scalar src_var = src_demean.rowwise().squaredNorm().sum() * one_over_n;
147 
148  // Eq. (42)
149  const Scalar c = Scalar(1) / src_var * svd.singularValues().dot(S);
150 
151  // Eq. (41)
152  Rt.col(m).head(m) = dst_mean;
153  Rt.col(m).head(m).noalias() -= c * Rt.topLeftCorner(m, m) * src_mean;
154  Rt.block(0, 0, m, m) *= c;
155  } else {
156  Rt.col(m).head(m) = dst_mean;
157  Rt.col(m).head(m).noalias() -= Rt.topLeftCorner(m, m) * src_mean;
158  }
159 
160  return Rt;
161 }
const unsigned n
Definition: CG3DPackingUnitTest.cpp:11
cout<< "Here is the matrix m:"<< endl<< m<< endl;JacobiSVD< MatrixXf, ComputeThinU|ComputeThinV > svd(m)
NumTraits< Scalar >::Real RealScalar
Definition: bench_gemm.cpp:46
MatrixXf MatrixType
Definition: benchmark-blocking-sizes.cpp:52
int * m
Definition: level2_cplx_impl.h:294
Eigen::Matrix< Scalar, Dynamic, Dynamic, ColMajor > tmp
Definition: level3_impl.h:365
constexpr int min_size_prefer_dynamic(A a, B b)
Definition: Meta.h:668
int c
Definition: calibrate.py:100
int sigma
Definition: calibrate.py:179
@ S
Definition: quadtree.h:62
@ IsComplex
Definition: NumTraits.h:176
T Real
Definition: NumTraits.h:183
Definition: fft_test_shared.h:66

References calibrate::c, Eigen::DenseBase< Derived >::colwise(), EIGEN_STATIC_ASSERT, m, Eigen::internal::min_size_prefer_dynamic(), n, Eigen::DenseBase< Derived >::rowwise(), oomph::QuadTreeNames::S, calibrate::sigma, Eigen::VectorwiseOp< ExpressionType, Direction >::sum(), svd(), tmp, and Eigen::value.

Referenced by run_fixed_size_test(), and run_test().

◆ unitOrthogonal()

template<typename Derived >
EIGEN_DEVICE_FUNC MatrixBase< Derived >::PlainObject Eigen::MatrixBase< Derived >::unitOrthogonal ( void  ) const
inline

\geometry_module

Returns
a unit vector which is orthogonal to *this

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().

See also
cross()
254  {
257 }
static EIGEN_DEVICE_FUNC VectorType run(const Derived &src)
Definition: OrthoMethods.h:185

References EIGEN_STATIC_ASSERT_VECTOR_ONLY, and Eigen::internal::unitOrthogonal_selector< Derived, Size >::run().