Eigen  3.4.90 (git rev a4098ac676528a83cfb73d4d26ce1b42ec05f47c)
Transform.h
1// This file is part of Eigen, a lightweight C++ template library
2// for linear algebra.
3//
4// Copyright (C) 2008 Gael Guennebaud <gael.guennebaud@inria.fr>
5// Copyright (C) 2009 Benoit Jacob <jacob.benoit.1@gmail.com>
6// Copyright (C) 2010 Hauke Heibel <hauke.heibel@gmail.com>
7//
8// This Source Code Form is subject to the terms of the Mozilla
9// Public License v. 2.0. If a copy of the MPL was not distributed
10// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
11
12#ifndef EIGEN_TRANSFORM_H
13#define EIGEN_TRANSFORM_H
14
15#include "./InternalHeaderCheck.h"
16
17namespace Eigen {
18
19namespace internal {
20
21template<typename Transform>
22struct transform_traits
23{
24 enum
25 {
26 Dim = Transform::Dim,
27 HDim = Transform::HDim,
28 Mode = Transform::Mode,
29 IsProjective = (int(Mode)==int(Projective))
30 };
31};
32
33template< typename TransformType,
34 typename MatrixType,
35 int Case = transform_traits<TransformType>::IsProjective ? 0
36 : int(MatrixType::RowsAtCompileTime) == int(transform_traits<TransformType>::HDim) ? 1
37 : 2,
38 int RhsCols = MatrixType::ColsAtCompileTime>
39struct transform_right_product_impl;
40
41template< typename Other,
42 int Mode,
43 int Options,
44 int Dim,
45 int HDim,
46 int OtherRows=Other::RowsAtCompileTime,
47 int OtherCols=Other::ColsAtCompileTime>
48struct transform_left_product_impl;
49
50template< typename Lhs,
51 typename Rhs,
52 bool AnyProjective =
53 transform_traits<Lhs>::IsProjective ||
54 transform_traits<Rhs>::IsProjective>
55struct transform_transform_product_impl;
56
57template< typename Other,
58 int Mode,
59 int Options,
60 int Dim,
61 int HDim,
62 int OtherRows=Other::RowsAtCompileTime,
63 int OtherCols=Other::ColsAtCompileTime>
64struct transform_construct_from_matrix;
65
66template<typename TransformType> struct transform_take_affine_part;
67
68template<typename Scalar_, int Dim_, int _Mode, int Options_>
69struct traits<Transform<Scalar_,Dim_,_Mode,Options_> >
70{
71 typedef Scalar_ Scalar;
72 typedef Eigen::Index StorageIndex;
73 typedef Dense StorageKind;
74 enum {
75 Dim1 = Dim_==Dynamic ? Dim_ : Dim_ + 1,
76 RowsAtCompileTime = _Mode==Projective ? Dim1 : Dim_,
77 ColsAtCompileTime = Dim1,
78 MaxRowsAtCompileTime = RowsAtCompileTime,
79 MaxColsAtCompileTime = ColsAtCompileTime,
80 Flags = 0
81 };
82};
83
84template<int Mode> struct transform_make_affine;
85
86} // end namespace internal
87
205template<typename Scalar_, int Dim_, int _Mode, int Options_>
207{
208public:
210 enum {
211 Mode = _Mode,
212 Options = Options_,
213 Dim = Dim_,
214 HDim = Dim_+1,
215 Rows = int(Mode)==(AffineCompact) ? Dim : HDim
216 };
218 typedef Scalar_ Scalar;
219 typedef Eigen::Index StorageIndex;
228 typedef Block<MatrixType,Dim,Dim,int(Mode)==(AffineCompact) && (int(Options)&RowMajor)==0> LinearPart;
230 typedef const Block<ConstMatrixType,Dim,Dim,int(Mode)==(AffineCompact) && (int(Options)&RowMajor)==0> ConstLinearPart;
232 typedef typename internal::conditional<int(Mode)==int(AffineCompact),
233 MatrixType&,
236 typedef typename internal::conditional<int(Mode)==int(AffineCompact),
237 const MatrixType&,
242 typedef Block<MatrixType,Dim,1,!(internal::traits<MatrixType>::Flags & RowMajorBit)> TranslationPart;
244 typedef const Block<ConstMatrixType,Dim,1,!(internal::traits<MatrixType>::Flags & RowMajorBit)> ConstTranslationPart;
247
248 // this intermediate enum is needed to avoid an ICE with gcc 3.4 and 4.0
249 enum { TransformTimeDiagonalMode = ((Mode==int(Isometry))?Affine:int(Mode)) };
252
253protected:
254
255 MatrixType m_matrix;
256
257public:
258
261 EIGEN_DEVICE_FUNC inline Transform()
262 {
263 check_template_params();
264 internal::transform_make_affine<(int(Mode)==Affine || int(Mode)==Isometry) ? Affine : AffineCompact>::run(m_matrix);
265 }
266
267 EIGEN_DEVICE_FUNC inline explicit Transform(const TranslationType& t)
268 {
269 check_template_params();
270 *this = t;
271 }
272 EIGEN_DEVICE_FUNC inline explicit Transform(const UniformScaling<Scalar>& s)
273 {
274 check_template_params();
275 *this = s;
276 }
277 template<typename Derived>
278 EIGEN_DEVICE_FUNC inline explicit Transform(const RotationBase<Derived, Dim>& r)
279 {
280 check_template_params();
281 *this = r;
282 }
283
284 typedef internal::transform_take_affine_part<Transform> take_affine_part;
285
287 template<typename OtherDerived>
288 EIGEN_DEVICE_FUNC inline explicit Transform(const EigenBase<OtherDerived>& other)
289 {
290 EIGEN_STATIC_ASSERT((internal::is_same<Scalar,typename OtherDerived::Scalar>::value),
291 YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY);
292
293 check_template_params();
294 internal::transform_construct_from_matrix<OtherDerived,Mode,Options,Dim,HDim>::run(this, other.derived());
295 }
296
298 template<typename OtherDerived>
299 EIGEN_DEVICE_FUNC inline Transform& operator=(const EigenBase<OtherDerived>& other)
300 {
301 EIGEN_STATIC_ASSERT((internal::is_same<Scalar,typename OtherDerived::Scalar>::value),
302 YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY);
303
304 internal::transform_construct_from_matrix<OtherDerived,Mode,Options,Dim,HDim>::run(this, other.derived());
305 return *this;
306 }
307
308 template<int OtherOptions>
309 EIGEN_DEVICE_FUNC inline Transform(const Transform<Scalar,Dim,Mode,OtherOptions>& other)
310 {
311 check_template_params();
312 // only the options change, we can directly copy the matrices
313 m_matrix = other.matrix();
314 }
315
316 template<int OtherMode,int OtherOptions>
317 EIGEN_DEVICE_FUNC inline Transform(const Transform<Scalar,Dim,OtherMode,OtherOptions>& other)
318 {
319 check_template_params();
320 // prevent conversions as:
321 // Affine | AffineCompact | Isometry = Projective
322 EIGEN_STATIC_ASSERT(internal::check_implication(OtherMode==int(Projective), Mode==int(Projective)),
323 YOU_PERFORMED_AN_INVALID_TRANSFORMATION_CONVERSION)
324
325 // prevent conversions as:
326 // Isometry = Affine | AffineCompact
327 EIGEN_STATIC_ASSERT(internal::check_implication(OtherMode==int(Affine)||OtherMode==int(AffineCompact), Mode!=int(Isometry)),
328 YOU_PERFORMED_AN_INVALID_TRANSFORMATION_CONVERSION)
329
330 enum { ModeIsAffineCompact = Mode == int(AffineCompact),
331 OtherModeIsAffineCompact = OtherMode == int(AffineCompact)
332 };
333
334 if(EIGEN_CONST_CONDITIONAL(ModeIsAffineCompact == OtherModeIsAffineCompact))
335 {
336 // We need the block expression because the code is compiled for all
337 // combinations of transformations and will trigger a compile time error
338 // if one tries to assign the matrices directly
339 m_matrix.template block<Dim,Dim+1>(0,0) = other.matrix().template block<Dim,Dim+1>(0,0);
340 makeAffine();
341 }
342 else if(EIGEN_CONST_CONDITIONAL(OtherModeIsAffineCompact))
343 {
344 typedef typename Transform<Scalar,Dim,OtherMode,OtherOptions>::MatrixType OtherMatrixType;
345 internal::transform_construct_from_matrix<OtherMatrixType,Mode,Options,Dim,HDim>::run(this, other.matrix());
346 }
347 else
348 {
349 // here we know that Mode == AffineCompact and OtherMode != AffineCompact.
350 // if OtherMode were Projective, the static assert above would already have caught it.
351 // So the only possibility is that OtherMode == Affine
352 linear() = other.linear();
353 translation() = other.translation();
354 }
355 }
356
357 template<typename OtherDerived>
358 EIGEN_DEVICE_FUNC Transform(const ReturnByValue<OtherDerived>& other)
359 {
360 check_template_params();
361 other.evalTo(*this);
362 }
363
364 template<typename OtherDerived>
365 EIGEN_DEVICE_FUNC Transform& operator=(const ReturnByValue<OtherDerived>& other)
366 {
367 other.evalTo(*this);
368 return *this;
369 }
370
371 #ifdef EIGEN_QT_SUPPORT
372 #if (QT_VERSION < QT_VERSION_CHECK(6, 0, 0))
373 inline Transform(const QMatrix& other);
374 inline Transform& operator=(const QMatrix& other);
375 inline QMatrix toQMatrix(void) const;
376 #endif
377 inline Transform(const QTransform& other);
378 inline Transform& operator=(const QTransform& other);
379 inline QTransform toQTransform(void) const;
380 #endif
381
382 EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT { return int(Mode)==int(Projective) ? m_matrix.cols() : (m_matrix.cols()-1); }
383 EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index cols() const EIGEN_NOEXCEPT { return m_matrix.cols(); }
384
387 EIGEN_DEVICE_FUNC inline Scalar operator() (Index row, Index col) const { return m_matrix(row,col); }
390 EIGEN_DEVICE_FUNC inline Scalar& operator() (Index row, Index col) { return m_matrix(row,col); }
391
393 EIGEN_DEVICE_FUNC inline const MatrixType& matrix() const { return m_matrix; }
395 EIGEN_DEVICE_FUNC inline MatrixType& matrix() { return m_matrix; }
396
398 EIGEN_DEVICE_FUNC inline ConstLinearPart linear() const { return ConstLinearPart(m_matrix,0,0); }
400 EIGEN_DEVICE_FUNC inline LinearPart linear() { return LinearPart(m_matrix,0,0); }
401
403 EIGEN_DEVICE_FUNC inline ConstAffinePart affine() const { return take_affine_part::run(m_matrix); }
405 EIGEN_DEVICE_FUNC inline AffinePart affine() { return take_affine_part::run(m_matrix); }
406
408 EIGEN_DEVICE_FUNC inline ConstTranslationPart translation() const { return ConstTranslationPart(m_matrix,0,Dim); }
410 EIGEN_DEVICE_FUNC inline TranslationPart translation() { return TranslationPart(m_matrix,0,Dim); }
411
436 // note: this function is defined here because some compilers cannot find the respective declaration
437 template<typename OtherDerived>
438 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const typename internal::transform_right_product_impl<Transform, OtherDerived>::ResultType
440 { return internal::transform_right_product_impl<Transform, OtherDerived>::run(*this,other.derived()); }
441
449 template<typename OtherDerived> friend
450 EIGEN_DEVICE_FUNC inline const typename internal::transform_left_product_impl<OtherDerived,Mode,Options,Dim_,Dim_+1>::ResultType
452 { return internal::transform_left_product_impl<OtherDerived,Mode,Options,Dim,HDim>::run(a.derived(),b); }
453
460 template<typename DiagonalDerived>
461 EIGEN_DEVICE_FUNC inline const TransformTimeDiagonalReturnType
462 operator * (const DiagonalBase<DiagonalDerived> &b) const
463 {
465 res.linearExt() *= b;
466 return res;
467 }
468
475 template<typename DiagonalDerived>
476 EIGEN_DEVICE_FUNC friend inline TransformTimeDiagonalReturnType
477 operator * (const DiagonalBase<DiagonalDerived> &a, const Transform &b)
478 {
480 res.linear().noalias() = a*b.linear();
481 res.translation().noalias() = a*b.translation();
482 if (EIGEN_CONST_CONDITIONAL(Mode!=int(AffineCompact)))
483 res.matrix().row(Dim) = b.matrix().row(Dim);
484 return res;
485 }
486
487 template<typename OtherDerived>
488 EIGEN_DEVICE_FUNC inline Transform& operator*=(const EigenBase<OtherDerived>& other) { return *this = *this * other; }
489
491 EIGEN_DEVICE_FUNC inline const Transform operator * (const Transform& other) const
492 {
493 return internal::transform_transform_product_impl<Transform,Transform>::run(*this,other);
494 }
495
496 #if EIGEN_COMP_ICC
497private:
498 // this intermediate structure permits to workaround a bug in ICC 11:
499 // error: template instantiation resulted in unexpected function type of "Eigen::Transform<double, 3, 32, 0>
500 // (const Eigen::Transform<double, 3, 2, 0> &) const"
501 // (the meaning of a name may have changed since the template declaration -- the type of the template is:
502 // "Eigen::internal::transform_transform_product_impl<Eigen::Transform<double, 3, 32, 0>,
503 // Eigen::Transform<double, 3, Mode, Options>, <expression>>::ResultType (const Eigen::Transform<double, 3, Mode, Options> &) const")
504 //
505 template<int OtherMode,int OtherOptions> struct icc_11_workaround
506 {
507 typedef internal::transform_transform_product_impl<Transform,Transform<Scalar,Dim,OtherMode,OtherOptions> > ProductType;
508 typedef typename ProductType::ResultType ResultType;
509 };
510
511public:
513 template<int OtherMode,int OtherOptions>
514 inline typename icc_11_workaround<OtherMode,OtherOptions>::ResultType
515 operator * (const Transform<Scalar,Dim,OtherMode,OtherOptions>& other) const
516 {
517 typedef typename icc_11_workaround<OtherMode,OtherOptions>::ProductType ProductType;
518 return ProductType::run(*this,other);
519 }
520 #else
522 template<int OtherMode,int OtherOptions>
523 EIGEN_DEVICE_FUNC inline typename internal::transform_transform_product_impl<Transform,Transform<Scalar,Dim,OtherMode,OtherOptions> >::ResultType
525 {
526 return internal::transform_transform_product_impl<Transform,Transform<Scalar,Dim,OtherMode,OtherOptions> >::run(*this,other);
527 }
528 #endif
529
531 EIGEN_DEVICE_FUNC void setIdentity() { m_matrix.setIdentity(); }
532
537 EIGEN_DEVICE_FUNC static const Transform Identity()
538 {
540 }
541
542 template<typename OtherDerived>
543 EIGEN_DEVICE_FUNC
544 inline Transform& scale(const MatrixBase<OtherDerived> &other);
545
546 template<typename OtherDerived>
547 EIGEN_DEVICE_FUNC
548 inline Transform& prescale(const MatrixBase<OtherDerived> &other);
549
550 EIGEN_DEVICE_FUNC inline Transform& scale(const Scalar& s);
551 EIGEN_DEVICE_FUNC inline Transform& prescale(const Scalar& s);
552
553 template<typename OtherDerived>
554 EIGEN_DEVICE_FUNC
555 inline Transform& translate(const MatrixBase<OtherDerived> &other);
556
557 template<typename OtherDerived>
558 EIGEN_DEVICE_FUNC
559 inline Transform& pretranslate(const MatrixBase<OtherDerived> &other);
560
561 template<typename RotationType>
562 EIGEN_DEVICE_FUNC
563 inline Transform& rotate(const RotationType& rotation);
564
565 template<typename RotationType>
566 EIGEN_DEVICE_FUNC
567 inline Transform& prerotate(const RotationType& rotation);
568
569 EIGEN_DEVICE_FUNC Transform& shear(const Scalar& sx, const Scalar& sy);
570 EIGEN_DEVICE_FUNC Transform& preshear(const Scalar& sx, const Scalar& sy);
571
572 EIGEN_DEVICE_FUNC inline Transform& operator=(const TranslationType& t);
573
574 EIGEN_DEVICE_FUNC
575 inline Transform& operator*=(const TranslationType& t) { return translate(t.vector()); }
576
577 EIGEN_DEVICE_FUNC inline Transform operator*(const TranslationType& t) const;
578
579 EIGEN_DEVICE_FUNC
580 inline Transform& operator=(const UniformScaling<Scalar>& t);
581
582 EIGEN_DEVICE_FUNC
583 inline Transform& operator*=(const UniformScaling<Scalar>& s) { return scale(s.factor()); }
584
585 EIGEN_DEVICE_FUNC
586 inline TransformTimeDiagonalReturnType operator*(const UniformScaling<Scalar>& s) const
587 {
589 res.scale(s.factor());
590 return res;
591 }
592
593 EIGEN_DEVICE_FUNC
594 inline Transform& operator*=(const DiagonalMatrix<Scalar,Dim>& s) { linearExt() *= s; return *this; }
595
596 template<typename Derived>
597 EIGEN_DEVICE_FUNC inline Transform& operator=(const RotationBase<Derived,Dim>& r);
598 template<typename Derived>
599 EIGEN_DEVICE_FUNC inline Transform& operator*=(const RotationBase<Derived,Dim>& r) { return rotate(r.toRotationMatrix()); }
600 template<typename Derived>
601 EIGEN_DEVICE_FUNC inline Transform operator*(const RotationBase<Derived,Dim>& r) const;
602
603 typedef typename internal::conditional<int(Mode)==Isometry,ConstLinearPart,const LinearMatrixType>::type RotationReturnType;
604 EIGEN_DEVICE_FUNC RotationReturnType rotation() const;
605
606 template<typename RotationMatrixType, typename ScalingMatrixType>
607 EIGEN_DEVICE_FUNC
608 void computeRotationScaling(RotationMatrixType *rotation, ScalingMatrixType *scaling) const;
609 template<typename ScalingMatrixType, typename RotationMatrixType>
610 EIGEN_DEVICE_FUNC
611 void computeScalingRotation(ScalingMatrixType *scaling, RotationMatrixType *rotation) const;
612
613 template<typename PositionDerived, typename OrientationType, typename ScaleDerived>
614 EIGEN_DEVICE_FUNC
615 Transform& fromPositionOrientationScale(const MatrixBase<PositionDerived> &position,
616 const OrientationType& orientation, const MatrixBase<ScaleDerived> &scale);
617
618 EIGEN_DEVICE_FUNC
619 inline Transform inverse(TransformTraits traits = (TransformTraits)Mode) const;
620
622 EIGEN_DEVICE_FUNC const Scalar* data() const { return m_matrix.data(); }
624 EIGEN_DEVICE_FUNC Scalar* data() { return m_matrix.data(); }
625
631 template<typename NewScalarType>
632 EIGEN_DEVICE_FUNC inline typename internal::cast_return_type<Transform,Transform<NewScalarType,Dim,Mode,Options> >::type cast() const
633 { return typename internal::cast_return_type<Transform,Transform<NewScalarType,Dim,Mode,Options> >::type(*this); }
634
636 template<typename OtherScalarType>
637 EIGEN_DEVICE_FUNC inline explicit Transform(const Transform<OtherScalarType,Dim,Mode,Options>& other)
638 {
639 check_template_params();
640 m_matrix = other.matrix().template cast<Scalar>();
641 }
642
647 EIGEN_DEVICE_FUNC bool isApprox(const Transform& other, const typename NumTraits<Scalar>::Real& prec = NumTraits<Scalar>::dummy_precision()) const
648 { return m_matrix.isApprox(other.m_matrix, prec); }
649
652 EIGEN_DEVICE_FUNC void makeAffine()
653 {
654 internal::transform_make_affine<int(Mode)>::run(m_matrix);
655 }
656
661 EIGEN_DEVICE_FUNC inline Block<MatrixType,int(Mode)==int(Projective)?HDim:Dim,Dim> linearExt()
662 { return m_matrix.template block<int(Mode)==int(Projective)?HDim:Dim,Dim>(0,0); }
667 EIGEN_DEVICE_FUNC inline const Block<MatrixType,int(Mode)==int(Projective)?HDim:Dim,Dim> linearExt() const
668 { return m_matrix.template block<int(Mode)==int(Projective)?HDim:Dim,Dim>(0,0); }
669
674 EIGEN_DEVICE_FUNC inline Block<MatrixType,int(Mode)==int(Projective)?HDim:Dim,1> translationExt()
675 { return m_matrix.template block<int(Mode)==int(Projective)?HDim:Dim,1>(0,Dim); }
680 EIGEN_DEVICE_FUNC inline const Block<MatrixType,int(Mode)==int(Projective)?HDim:Dim,1> translationExt() const
681 { return m_matrix.template block<int(Mode)==int(Projective)?HDim:Dim,1>(0,Dim); }
682
683
684 #ifdef EIGEN_TRANSFORM_PLUGIN
685 #include EIGEN_TRANSFORM_PLUGIN
686 #endif
687
688protected:
689 #ifndef EIGEN_PARSED_BY_DOXYGEN
690 EIGEN_DEVICE_FUNC static EIGEN_STRONG_INLINE void check_template_params()
691 {
692 EIGEN_STATIC_ASSERT((Options & (DontAlign|RowMajor)) == Options, INVALID_MATRIX_TEMPLATE_PARAMETERS)
693 }
694 #endif
695
696};
697
699typedef Transform<float,2,Isometry> Isometry2f;
701typedef Transform<float,3,Isometry> Isometry3f;
703typedef Transform<double,2,Isometry> Isometry2d;
705typedef Transform<double,3,Isometry> Isometry3d;
706
708typedef Transform<float,2,Affine> Affine2f;
710typedef Transform<float,3,Affine> Affine3f;
712typedef Transform<double,2,Affine> Affine2d;
714typedef Transform<double,3,Affine> Affine3d;
715
717typedef Transform<float,2,AffineCompact> AffineCompact2f;
719typedef Transform<float,3,AffineCompact> AffineCompact3f;
721typedef Transform<double,2,AffineCompact> AffineCompact2d;
723typedef Transform<double,3,AffineCompact> AffineCompact3d;
724
726typedef Transform<float,2,Projective> Projective2f;
728typedef Transform<float,3,Projective> Projective3f;
730typedef Transform<double,2,Projective> Projective2d;
732typedef Transform<double,3,Projective> Projective3d;
733
734/**************************
735*** Optional QT support ***
736**************************/
737
738#ifdef EIGEN_QT_SUPPORT
739
740#if (QT_VERSION < QT_VERSION_CHECK(6, 0, 0))
745template<typename Scalar, int Dim, int Mode,int Options>
747{
748 check_template_params();
749 *this = other;
750}
751
756template<typename Scalar, int Dim, int Mode,int Options>
757Transform<Scalar,Dim,Mode,Options>& Transform<Scalar,Dim,Mode,Options>::operator=(const QMatrix& other)
758{
759 EIGEN_STATIC_ASSERT(Dim==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
760 if (EIGEN_CONST_CONDITIONAL(Mode == int(AffineCompact)))
761 m_matrix << other.m11(), other.m21(), other.dx(),
762 other.m12(), other.m22(), other.dy();
763 else
764 m_matrix << other.m11(), other.m21(), other.dx(),
765 other.m12(), other.m22(), other.dy(),
766 0, 0, 1;
767 return *this;
768}
769
776template<typename Scalar, int Dim, int Mode, int Options>
777QMatrix Transform<Scalar,Dim,Mode,Options>::toQMatrix(void) const
778{
779 check_template_params();
780 EIGEN_STATIC_ASSERT(Dim==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
781 return QMatrix(m_matrix.coeff(0,0), m_matrix.coeff(1,0),
782 m_matrix.coeff(0,1), m_matrix.coeff(1,1),
783 m_matrix.coeff(0,2), m_matrix.coeff(1,2));
784}
785#endif
786
791template<typename Scalar, int Dim, int Mode,int Options>
793{
794 check_template_params();
795 *this = other;
796}
797
802template<typename Scalar, int Dim, int Mode, int Options>
804{
805 check_template_params();
806 EIGEN_STATIC_ASSERT(Dim==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
807 if (EIGEN_CONST_CONDITIONAL(Mode == int(AffineCompact)))
808 m_matrix << other.m11(), other.m21(), other.dx(),
809 other.m12(), other.m22(), other.dy();
810 else
811 m_matrix << other.m11(), other.m21(), other.dx(),
812 other.m12(), other.m22(), other.dy(),
813 other.m13(), other.m23(), other.m33();
814 return *this;
815}
816
821template<typename Scalar, int Dim, int Mode, int Options>
823{
824 EIGEN_STATIC_ASSERT(Dim==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
825 if (EIGEN_CONST_CONDITIONAL(Mode == int(AffineCompact)))
826 return QTransform(m_matrix.coeff(0,0), m_matrix.coeff(1,0),
827 m_matrix.coeff(0,1), m_matrix.coeff(1,1),
828 m_matrix.coeff(0,2), m_matrix.coeff(1,2));
829 else
830 return QTransform(m_matrix.coeff(0,0), m_matrix.coeff(1,0), m_matrix.coeff(2,0),
831 m_matrix.coeff(0,1), m_matrix.coeff(1,1), m_matrix.coeff(2,1),
832 m_matrix.coeff(0,2), m_matrix.coeff(1,2), m_matrix.coeff(2,2));
833}
834#endif
835
836/*********************
837*** Procedural API ***
838*********************/
839
844template<typename Scalar, int Dim, int Mode, int Options>
845template<typename OtherDerived>
846EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
848{
849 EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,int(Dim))
850 EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
851 linearExt().noalias() = (linearExt() * other.asDiagonal());
852 return *this;
853}
854
859template<typename Scalar, int Dim, int Mode, int Options>
861{
862 EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
863 linearExt() *= s;
864 return *this;
865}
866
871template<typename Scalar, int Dim, int Mode, int Options>
872template<typename OtherDerived>
873EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
875{
876 EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,int(Dim))
877 EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
878 affine().noalias() = (other.asDiagonal() * affine());
879 return *this;
880}
881
886template<typename Scalar, int Dim, int Mode, int Options>
888{
889 EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
890 m_matrix.template topRows<Dim>() *= s;
891 return *this;
892}
893
898template<typename Scalar, int Dim, int Mode, int Options>
899template<typename OtherDerived>
900EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
902{
903 EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,int(Dim))
904 translationExt() += linearExt() * other;
905 return *this;
906}
907
912template<typename Scalar, int Dim, int Mode, int Options>
913template<typename OtherDerived>
914EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
916{
917 EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,int(Dim))
918 if(EIGEN_CONST_CONDITIONAL(int(Mode)==int(Projective)))
919 affine() += other * m_matrix.row(Dim);
920 else
921 translation() += other;
922 return *this;
923}
924
942template<typename Scalar, int Dim, int Mode, int Options>
943template<typename RotationType>
944EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
946{
947 linearExt() *= internal::toRotationMatrix<Scalar,Dim>(rotation);
948 return *this;
949}
950
958template<typename Scalar, int Dim, int Mode, int Options>
959template<typename RotationType>
960EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
962{
963 m_matrix.template block<Dim,HDim>(0,0) = internal::toRotationMatrix<Scalar,Dim>(rotation)
964 * m_matrix.template block<Dim,HDim>(0,0);
965 return *this;
966}
967
973template<typename Scalar, int Dim, int Mode, int Options>
974EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
976{
977 EIGEN_STATIC_ASSERT(int(Dim)==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
978 EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
979 VectorType tmp = linear().col(0)*sy + linear().col(1);
980 linear() << linear().col(0) + linear().col(1)*sx, tmp;
981 return *this;
982}
983
989template<typename Scalar, int Dim, int Mode, int Options>
990EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
992{
993 EIGEN_STATIC_ASSERT(int(Dim)==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
994 EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
995 m_matrix.template block<Dim,HDim>(0,0) = LinearMatrixType(1, sx, sy, 1) * m_matrix.template block<Dim,HDim>(0,0);
996 return *this;
997}
998
999/******************************************************
1000*** Scaling, Translation and Rotation compatibility ***
1001******************************************************/
1002
1003template<typename Scalar, int Dim, int Mode, int Options>
1004EIGEN_DEVICE_FUNC inline Transform<Scalar,Dim,Mode,Options>& Transform<Scalar,Dim,Mode,Options>::operator=(const TranslationType& t)
1005{
1006 linear().setIdentity();
1007 translation() = t.vector();
1008 makeAffine();
1009 return *this;
1010}
1011
1012template<typename Scalar, int Dim, int Mode, int Options>
1013EIGEN_DEVICE_FUNC inline Transform<Scalar,Dim,Mode,Options> Transform<Scalar,Dim,Mode,Options>::operator*(const TranslationType& t) const
1014{
1015 Transform res = *this;
1016 res.translate(t.vector());
1017 return res;
1018}
1019
1020template<typename Scalar, int Dim, int Mode, int Options>
1021EIGEN_DEVICE_FUNC inline Transform<Scalar,Dim,Mode,Options>& Transform<Scalar,Dim,Mode,Options>::operator=(const UniformScaling<Scalar>& s)
1022{
1023 m_matrix.setZero();
1024 linear().diagonal().fill(s.factor());
1025 makeAffine();
1026 return *this;
1027}
1028
1029template<typename Scalar, int Dim, int Mode, int Options>
1030template<typename Derived>
1031EIGEN_DEVICE_FUNC inline Transform<Scalar,Dim,Mode,Options>& Transform<Scalar,Dim,Mode,Options>::operator=(const RotationBase<Derived,Dim>& r)
1032{
1033 linear() = internal::toRotationMatrix<Scalar,Dim>(r);
1034 translation().setZero();
1035 makeAffine();
1036 return *this;
1037}
1038
1039template<typename Scalar, int Dim, int Mode, int Options>
1040template<typename Derived>
1041EIGEN_DEVICE_FUNC inline Transform<Scalar,Dim,Mode,Options> Transform<Scalar,Dim,Mode,Options>::operator*(const RotationBase<Derived,Dim>& r) const
1042{
1043 Transform res = *this;
1044 res.rotate(r.derived());
1045 return res;
1046}
1047
1048/************************
1049*** Special functions ***
1050************************/
1051
1052namespace internal {
1053template<int Mode> struct transform_rotation_impl {
1054 template<typename TransformType>
1055 EIGEN_DEVICE_FUNC static inline
1056 const typename TransformType::LinearMatrixType run(const TransformType& t)
1057 {
1058 typedef typename TransformType::LinearMatrixType LinearMatrixType;
1059 LinearMatrixType result;
1060 t.computeRotationScaling(&result, (LinearMatrixType*)0);
1061 return result;
1062 }
1063};
1064template<> struct transform_rotation_impl<Isometry> {
1065 template<typename TransformType>
1066 EIGEN_DEVICE_FUNC static inline
1067 typename TransformType::ConstLinearPart run(const TransformType& t)
1068 {
1069 return t.linear();
1070 }
1071};
1072}
1083template<typename Scalar, int Dim, int Mode, int Options>
1084EIGEN_DEVICE_FUNC
1085typename Transform<Scalar,Dim,Mode,Options>::RotationReturnType
1087{
1088 return internal::transform_rotation_impl<Mode>::run(*this);
1089}
1090
1091
1103template<typename Scalar, int Dim, int Mode, int Options>
1104template<typename RotationMatrixType, typename ScalingMatrixType>
1105EIGEN_DEVICE_FUNC void Transform<Scalar,Dim,Mode,Options>::computeRotationScaling(RotationMatrixType *rotation, ScalingMatrixType *scaling) const
1106{
1107 // Note that JacobiSVD is faster than BDCSVD for small matrices.
1109
1110 Scalar x = (svd.matrixU() * svd.matrixV().adjoint()).determinant() < Scalar(0) ? Scalar(-1) : Scalar(1); // so x has absolute value 1
1111 VectorType sv(svd.singularValues());
1112 sv.coeffRef(Dim-1) *= x;
1113 if(scaling) *scaling = svd.matrixV() * sv.asDiagonal() * svd.matrixV().adjoint();
1114 if(rotation)
1115 {
1116 LinearMatrixType m(svd.matrixU());
1117 m.col(Dim-1) *= x;
1118 *rotation = m * svd.matrixV().adjoint();
1119 }
1120}
1121
1133template<typename Scalar, int Dim, int Mode, int Options>
1134template<typename ScalingMatrixType, typename RotationMatrixType>
1135EIGEN_DEVICE_FUNC void Transform<Scalar,Dim,Mode,Options>::computeScalingRotation(ScalingMatrixType *scaling, RotationMatrixType *rotation) const
1136{
1137 // Note that JacobiSVD is faster than BDCSVD for small matrices.
1139
1140 Scalar x = (svd.matrixU() * svd.matrixV().adjoint()).determinant() < Scalar(0) ? Scalar(-1) : Scalar(1); // so x has absolute value 1
1141 VectorType sv(svd.singularValues());
1142 sv.coeffRef(Dim-1) *= x;
1143 if(scaling) *scaling = svd.matrixU() * sv.asDiagonal() * svd.matrixU().adjoint();
1144 if(rotation)
1145 {
1146 LinearMatrixType m(svd.matrixU());
1147 m.col(Dim-1) *= x;
1148 *rotation = m * svd.matrixV().adjoint();
1149 }
1150}
1151
1155template<typename Scalar, int Dim, int Mode, int Options>
1156template<typename PositionDerived, typename OrientationType, typename ScaleDerived>
1157EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
1159 const OrientationType& orientation, const MatrixBase<ScaleDerived> &scale)
1160{
1161 linear() = internal::toRotationMatrix<Scalar,Dim>(orientation);
1162 linear() *= scale.asDiagonal();
1163 translation() = position;
1164 makeAffine();
1165 return *this;
1166}
1167
1168namespace internal {
1169
1170template<int Mode>
1171struct transform_make_affine
1172{
1173 template<typename MatrixType>
1174 EIGEN_DEVICE_FUNC static void run(MatrixType &mat)
1175 {
1176 static const int Dim = MatrixType::ColsAtCompileTime-1;
1177 mat.template block<1,Dim>(Dim,0).setZero();
1178 mat.coeffRef(Dim,Dim) = typename MatrixType::Scalar(1);
1179 }
1180};
1181
1182template<>
1183struct transform_make_affine<AffineCompact>
1184{
1185 template<typename MatrixType> EIGEN_DEVICE_FUNC static void run(MatrixType &) { }
1186};
1187
1188// selector needed to avoid taking the inverse of a 3x4 matrix
1189template<typename TransformType, int Mode=TransformType::Mode>
1190struct projective_transform_inverse
1191{
1192 EIGEN_DEVICE_FUNC static inline void run(const TransformType&, TransformType&)
1193 {}
1194};
1195
1196template<typename TransformType>
1197struct projective_transform_inverse<TransformType, Projective>
1198{
1199 EIGEN_DEVICE_FUNC static inline void run(const TransformType& m, TransformType& res)
1200 {
1201 res.matrix() = m.matrix().inverse();
1202 }
1203};
1204
1205} // end namespace internal
1206
1207
1228template<typename Scalar, int Dim, int Mode, int Options>
1229EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>
1231{
1232 Transform res;
1233 if (hint == Projective)
1234 {
1235 internal::projective_transform_inverse<Transform>::run(*this, res);
1236 }
1237 else
1238 {
1239 if (hint == Isometry)
1240 {
1241 res.matrix().template topLeftCorner<Dim,Dim>() = linear().transpose();
1242 }
1243 else if(hint&Affine)
1244 {
1245 res.matrix().template topLeftCorner<Dim,Dim>() = linear().inverse();
1246 }
1247 else
1248 {
1249 eigen_assert(false && "Invalid transform traits in Transform::Inverse");
1250 }
1251 // translation and remaining parts
1252 res.matrix().template topRightCorner<Dim,1>()
1253 = - res.matrix().template topLeftCorner<Dim,Dim>() * translation();
1254 res.makeAffine(); // we do need this, because in the beginning res is uninitialized
1255 }
1256 return res;
1257}
1258
1259namespace internal {
1260
1261/*****************************************************
1262*** Specializations of take affine part ***
1263*****************************************************/
1264
1265template<typename TransformType> struct transform_take_affine_part {
1266 typedef typename TransformType::MatrixType MatrixType;
1267 typedef typename TransformType::AffinePart AffinePart;
1268 typedef typename TransformType::ConstAffinePart ConstAffinePart;
1269 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE AffinePart run(MatrixType& m)
1270 { return m.template block<TransformType::Dim,TransformType::HDim>(0,0); }
1271 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE ConstAffinePart run(const MatrixType& m)
1272 { return m.template block<TransformType::Dim,TransformType::HDim>(0,0); }
1273};
1274
1275template<typename Scalar, int Dim, int Options>
1276struct transform_take_affine_part<Transform<Scalar,Dim,AffineCompact, Options> > {
1278 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE MatrixType& run(MatrixType& m) { return m; }
1279 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const MatrixType& run(const MatrixType& m) { return m; }
1280};
1281
1282/*****************************************************
1283*** Specializations of construct from matrix ***
1284*****************************************************/
1285
1286template<typename Other, int Mode, int Options, int Dim, int HDim>
1287struct transform_construct_from_matrix<Other, Mode,Options,Dim,HDim, Dim,Dim>
1288{
1289 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Transform<typename Other::Scalar,Dim,Mode,Options> *transform, const Other& other)
1290 {
1291 transform->linear() = other;
1292 transform->translation().setZero();
1293 transform->makeAffine();
1294 }
1295};
1296
1297template<typename Other, int Mode, int Options, int Dim, int HDim>
1298struct transform_construct_from_matrix<Other, Mode,Options,Dim,HDim, Dim,HDim>
1299{
1300 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Transform<typename Other::Scalar,Dim,Mode,Options> *transform, const Other& other)
1301 {
1302 transform->affine() = other;
1303 transform->makeAffine();
1304 }
1305};
1306
1307template<typename Other, int Mode, int Options, int Dim, int HDim>
1308struct transform_construct_from_matrix<Other, Mode,Options,Dim,HDim, HDim,HDim>
1309{
1310 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Transform<typename Other::Scalar,Dim,Mode,Options> *transform, const Other& other)
1311 { transform->matrix() = other; }
1312};
1313
1314template<typename Other, int Options, int Dim, int HDim>
1315struct transform_construct_from_matrix<Other, AffineCompact,Options,Dim,HDim, HDim,HDim>
1316{
1317 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE void run(Transform<typename Other::Scalar,Dim,AffineCompact,Options> *transform, const Other& other)
1318 { transform->matrix() = other.template block<Dim,HDim>(0,0); }
1319};
1320
1321/**********************************************************
1322*** Specializations of operator* with rhs EigenBase ***
1323**********************************************************/
1324
1325template<int LhsMode,int RhsMode>
1326struct transform_product_result
1327{
1328 enum
1329 {
1330 Mode =
1331 (LhsMode == (int)Projective || RhsMode == (int)Projective ) ? Projective :
1332 (LhsMode == (int)Affine || RhsMode == (int)Affine ) ? Affine :
1333 (LhsMode == (int)AffineCompact || RhsMode == (int)AffineCompact ) ? AffineCompact :
1334 (LhsMode == (int)Isometry || RhsMode == (int)Isometry ) ? Isometry : Projective
1335 };
1336};
1337
1338template< typename TransformType, typename MatrixType, int RhsCols>
1339struct transform_right_product_impl< TransformType, MatrixType, 0, RhsCols>
1340{
1341 typedef typename MatrixType::PlainObject ResultType;
1342
1343 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE ResultType run(const TransformType& T, const MatrixType& other)
1344 {
1345 return T.matrix() * other;
1346 }
1347};
1348
1349template< typename TransformType, typename MatrixType, int RhsCols>
1350struct transform_right_product_impl< TransformType, MatrixType, 1, RhsCols>
1351{
1352 enum {
1353 Dim = TransformType::Dim,
1354 HDim = TransformType::HDim,
1355 OtherRows = MatrixType::RowsAtCompileTime,
1356 OtherCols = MatrixType::ColsAtCompileTime
1357 };
1358
1359 typedef typename MatrixType::PlainObject ResultType;
1360
1361 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE ResultType run(const TransformType& T, const MatrixType& other)
1362 {
1363 EIGEN_STATIC_ASSERT(OtherRows==HDim, YOU_MIXED_MATRICES_OF_DIFFERENT_SIZES);
1364
1365 typedef Block<ResultType, Dim, OtherCols, int(MatrixType::RowsAtCompileTime)==Dim> TopLeftLhs;
1366
1367 ResultType res(other.rows(),other.cols());
1368 TopLeftLhs(res, 0, 0, Dim, other.cols()).noalias() = T.affine() * other;
1369 res.row(OtherRows-1) = other.row(OtherRows-1);
1370
1371 return res;
1372 }
1373};
1374
1375template< typename TransformType, typename MatrixType, int RhsCols>
1376struct transform_right_product_impl< TransformType, MatrixType, 2, RhsCols>
1377{
1378 enum {
1379 Dim = TransformType::Dim,
1380 HDim = TransformType::HDim,
1381 OtherRows = MatrixType::RowsAtCompileTime,
1382 OtherCols = MatrixType::ColsAtCompileTime
1383 };
1384
1385 typedef typename MatrixType::PlainObject ResultType;
1386
1387 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE ResultType run(const TransformType& T, const MatrixType& other)
1388 {
1389 EIGEN_STATIC_ASSERT(OtherRows==Dim, YOU_MIXED_MATRICES_OF_DIFFERENT_SIZES);
1390
1391 typedef Block<ResultType, Dim, OtherCols, true> TopLeftLhs;
1392 ResultType res(Replicate<typename TransformType::ConstTranslationPart, 1, OtherCols>(T.translation(),1,other.cols()));
1393 TopLeftLhs(res, 0, 0, Dim, other.cols()).noalias() += T.linear() * other;
1394
1395 return res;
1396 }
1397};
1398
1399template< typename TransformType, typename MatrixType >
1400struct transform_right_product_impl< TransformType, MatrixType, 2, 1> // rhs is a vector of size Dim
1401{
1402 typedef typename TransformType::MatrixType TransformMatrix;
1403 enum {
1404 Dim = TransformType::Dim,
1405 HDim = TransformType::HDim,
1406 OtherRows = MatrixType::RowsAtCompileTime,
1407 WorkingRows = plain_enum_min(TransformMatrix::RowsAtCompileTime, HDim)
1408 };
1409
1410 typedef typename MatrixType::PlainObject ResultType;
1411
1412 static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE ResultType run(const TransformType& T, const MatrixType& other)
1413 {
1414 EIGEN_STATIC_ASSERT(OtherRows==Dim, YOU_MIXED_MATRICES_OF_DIFFERENT_SIZES);
1415
1416 Matrix<typename ResultType::Scalar, Dim+1, 1> rhs;
1417 rhs.template head<Dim>() = other; rhs[Dim] = typename ResultType::Scalar(1);
1418 Matrix<typename ResultType::Scalar, WorkingRows, 1> res(T.matrix() * rhs);
1419 return res.template head<Dim>();
1420 }
1421};
1422
1423/**********************************************************
1424*** Specializations of operator* with lhs EigenBase ***
1425**********************************************************/
1426
1427// generic HDim x HDim matrix * T => Projective
1428template<typename Other,int Mode, int Options, int Dim, int HDim>
1429struct transform_left_product_impl<Other,Mode,Options,Dim,HDim, HDim,HDim>
1430{
1431 typedef Transform<typename Other::Scalar,Dim,Mode,Options> TransformType;
1432 typedef typename TransformType::MatrixType MatrixType;
1433 typedef Transform<typename Other::Scalar,Dim,Projective,Options> ResultType;
1434 static ResultType run(const Other& other,const TransformType& tr)
1435 { return ResultType(other * tr.matrix()); }
1436};
1437
1438// generic HDim x HDim matrix * AffineCompact => Projective
1439template<typename Other, int Options, int Dim, int HDim>
1440struct transform_left_product_impl<Other,AffineCompact,Options,Dim,HDim, HDim,HDim>
1441{
1442 typedef Transform<typename Other::Scalar,Dim,AffineCompact,Options> TransformType;
1443 typedef typename TransformType::MatrixType MatrixType;
1444 typedef Transform<typename Other::Scalar,Dim,Projective,Options> ResultType;
1445 static ResultType run(const Other& other,const TransformType& tr)
1446 {
1447 ResultType res;
1448 res.matrix().noalias() = other.template block<HDim,Dim>(0,0) * tr.matrix();
1449 res.matrix().col(Dim) += other.col(Dim);
1450 return res;
1451 }
1452};
1453
1454// affine matrix * T
1455template<typename Other,int Mode, int Options, int Dim, int HDim>
1456struct transform_left_product_impl<Other,Mode,Options,Dim,HDim, Dim,HDim>
1457{
1458 typedef Transform<typename Other::Scalar,Dim,Mode,Options> TransformType;
1459 typedef typename TransformType::MatrixType MatrixType;
1460 typedef TransformType ResultType;
1461 static ResultType run(const Other& other,const TransformType& tr)
1462 {
1463 ResultType res;
1464 res.affine().noalias() = other * tr.matrix();
1465 res.matrix().row(Dim) = tr.matrix().row(Dim);
1466 return res;
1467 }
1468};
1469
1470// affine matrix * AffineCompact
1471template<typename Other, int Options, int Dim, int HDim>
1472struct transform_left_product_impl<Other,AffineCompact,Options,Dim,HDim, Dim,HDim>
1473{
1474 typedef Transform<typename Other::Scalar,Dim,AffineCompact,Options> TransformType;
1475 typedef typename TransformType::MatrixType MatrixType;
1476 typedef TransformType ResultType;
1477 static ResultType run(const Other& other,const TransformType& tr)
1478 {
1479 ResultType res;
1480 res.matrix().noalias() = other.template block<Dim,Dim>(0,0) * tr.matrix();
1481 res.translation() += other.col(Dim);
1482 return res;
1483 }
1484};
1485
1486// linear matrix * T
1487template<typename Other,int Mode, int Options, int Dim, int HDim>
1488struct transform_left_product_impl<Other,Mode,Options,Dim,HDim, Dim,Dim>
1489{
1490 typedef Transform<typename Other::Scalar,Dim,Mode,Options> TransformType;
1491 typedef typename TransformType::MatrixType MatrixType;
1492 typedef TransformType ResultType;
1493 static ResultType run(const Other& other, const TransformType& tr)
1494 {
1495 TransformType res;
1496 if(Mode!=int(AffineCompact))
1497 res.matrix().row(Dim) = tr.matrix().row(Dim);
1498 res.matrix().template topRows<Dim>().noalias()
1499 = other * tr.matrix().template topRows<Dim>();
1500 return res;
1501 }
1502};
1503
1504/**********************************************************
1505*** Specializations of operator* with another Transform ***
1506**********************************************************/
1507
1508template<typename Scalar, int Dim, int LhsMode, int LhsOptions, int RhsMode, int RhsOptions>
1509struct transform_transform_product_impl<Transform<Scalar,Dim,LhsMode,LhsOptions>,Transform<Scalar,Dim,RhsMode,RhsOptions>,false >
1510{
1511 enum { ResultMode = transform_product_result<LhsMode,RhsMode>::Mode };
1512 typedef Transform<Scalar,Dim,LhsMode,LhsOptions> Lhs;
1513 typedef Transform<Scalar,Dim,RhsMode,RhsOptions> Rhs;
1514 typedef Transform<Scalar,Dim,ResultMode,LhsOptions> ResultType;
1515 static ResultType run(const Lhs& lhs, const Rhs& rhs)
1516 {
1517 ResultType res;
1518 res.linear() = lhs.linear() * rhs.linear();
1519 res.translation() = lhs.linear() * rhs.translation() + lhs.translation();
1520 res.makeAffine();
1521 return res;
1522 }
1523};
1524
1525template<typename Scalar, int Dim, int LhsMode, int LhsOptions, int RhsMode, int RhsOptions>
1526struct transform_transform_product_impl<Transform<Scalar,Dim,LhsMode,LhsOptions>,Transform<Scalar,Dim,RhsMode,RhsOptions>,true >
1527{
1528 typedef Transform<Scalar,Dim,LhsMode,LhsOptions> Lhs;
1529 typedef Transform<Scalar,Dim,RhsMode,RhsOptions> Rhs;
1530 typedef Transform<Scalar,Dim,Projective> ResultType;
1531 static ResultType run(const Lhs& lhs, const Rhs& rhs)
1532 {
1533 return ResultType( lhs.matrix() * rhs.matrix() );
1534 }
1535};
1536
1537template<typename Scalar, int Dim, int LhsOptions, int RhsOptions>
1538struct transform_transform_product_impl<Transform<Scalar,Dim,AffineCompact,LhsOptions>,Transform<Scalar,Dim,Projective,RhsOptions>,true >
1539{
1540 typedef Transform<Scalar,Dim,AffineCompact,LhsOptions> Lhs;
1541 typedef Transform<Scalar,Dim,Projective,RhsOptions> Rhs;
1542 typedef Transform<Scalar,Dim,Projective> ResultType;
1543 static ResultType run(const Lhs& lhs, const Rhs& rhs)
1544 {
1545 ResultType res;
1546 res.matrix().template topRows<Dim>() = lhs.matrix() * rhs.matrix();
1547 res.matrix().row(Dim) = rhs.matrix().row(Dim);
1548 return res;
1549 }
1550};
1551
1552template<typename Scalar, int Dim, int LhsOptions, int RhsOptions>
1553struct transform_transform_product_impl<Transform<Scalar,Dim,Projective,LhsOptions>,Transform<Scalar,Dim,AffineCompact,RhsOptions>,true >
1554{
1555 typedef Transform<Scalar,Dim,Projective,LhsOptions> Lhs;
1556 typedef Transform<Scalar,Dim,AffineCompact,RhsOptions> Rhs;
1557 typedef Transform<Scalar,Dim,Projective> ResultType;
1558 static ResultType run(const Lhs& lhs, const Rhs& rhs)
1559 {
1560 ResultType res(lhs.matrix().template leftCols<Dim>() * rhs.matrix());
1561 res.matrix().col(Dim) += lhs.matrix().col(Dim);
1562 return res;
1563 }
1564};
1565
1566} // end namespace internal
1567
1568} // end namespace Eigen
1569
1570#endif // EIGEN_TRANSFORM_H
Expression of a fixed-size or dynamic-size block.
Definition: Block.h:107
TransposeReturnType transpose()
Definition: Transpose.h:184
bool isApprox(const DenseBase< OtherDerived > &other, const RealScalar &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: Fuzzy.h:105
Two-sided Jacobi SVD decomposition of a rectangular matrix.
Definition: JacobiSVD.h:492
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:52
const DiagonalWrapper< const Derived > asDiagonal() const
Definition: DiagonalMatrix.h:325
Derived & setIdentity()
Definition: CwiseNullaryOp.h:875
const Inverse< Derived > inverse() const
Definition: InverseImpl.h:350
static const IdentityReturnType Identity()
Definition: CwiseNullaryOp.h:801
The matrix class, also used for vectors and row-vectors.
Definition: Matrix.h:182
Scalar & coeffRef(Index rowId, Index colId)
Definition: PlainObjectBase.h:187
const Scalar & coeff(Index rowId, Index colId) const
Definition: PlainObjectBase.h:164
Derived & setZero(Index size)
Definition: CwiseNullaryOp.h:564
const Scalar * data() const
Definition: PlainObjectBase.h:259
const SingularValuesType & singularValues() const
Definition: SVDBase.h:131
const MatrixUType & matrixU() const
Definition: SVDBase.h:103
const MatrixVType & matrixV() const
Definition: SVDBase.h:119
Represents an homogeneous transformation in a N dimensional space.
Definition: Transform.h:207
Matrix< Scalar, Dim, 1 > VectorType
Definition: Transform.h:240
const MatrixType ConstMatrixType
Definition: Transform.h:224
Transform()
Definition: Transform.h:261
internal::conditional< int(Mode)==int(AffineCompact), constMatrixType &, constBlock< constMatrixType, Dim, HDim > >::type ConstAffinePart
Definition: Transform.h:238
internal::conditional< int(Mode)==int(AffineCompact), MatrixType &, Block< MatrixType, Dim, HDim > >::type AffinePart
Definition: Transform.h:234
void makeAffine()
Definition: Transform.h:652
Block< MatrixType, Dim, Dim, int(Mode)==(AffineCompact) &&(int(Options)&RowMajor)==0 > LinearPart
Definition: Transform.h:228
ConstTranslationPart translation() const
Definition: Transform.h:408
Transform & preshear(const Scalar &sx, const Scalar &sy)
Definition: Transform.h:991
void computeScalingRotation(ScalingMatrixType *scaling, RotationMatrixType *rotation) const
Definition: Transform.h:1135
MatrixType & matrix()
Definition: Transform.h:395
AffinePart affine()
Definition: Transform.h:405
bool isApprox(const Transform &other, const typename NumTraits< Scalar >::Real &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: Transform.h:647
Scalar * data()
Definition: Transform.h:624
Block< MatrixType, Dim, 1,!(internal::traits< MatrixType >::Flags &RowMajorBit)> TranslationPart
Definition: Transform.h:242
Transform(const EigenBase< OtherDerived > &other)
Definition: Transform.h:288
Eigen::Index Index
Definition: Transform.h:220
LinearPart linear()
Definition: Transform.h:400
Transform< Scalar, Dim, TransformTimeDiagonalMode > TransformTimeDiagonalReturnType
Definition: Transform.h:251
const Scalar * data() const
Definition: Transform.h:622
const internal::transform_right_product_impl< Transform, OtherDerived >::ResultType operator*(const EigenBase< OtherDerived > &other) const
Definition: Transform.h:439
internal::make_proper_matrix_type< Scalar, Rows, HDim, Options >::type MatrixType
Definition: Transform.h:222
ConstAffinePart affine() const
Definition: Transform.h:403
void computeRotationScaling(RotationMatrixType *rotation, ScalingMatrixType *scaling) const
Definition: Transform.h:1105
Transform(const Transform< OtherScalarType, Dim, Mode, Options > &other)
Definition: Transform.h:637
const Block< ConstMatrixType, Dim, Dim, int(Mode)==(AffineCompact) &&(int(Options)&RowMajor)==0 > ConstLinearPart
Definition: Transform.h:230
Translation< Scalar, Dim > TranslationType
Definition: Transform.h:246
Transform & operator=(const EigenBase< OtherDerived > &other)
Definition: Transform.h:299
Scalar operator()(Index row, Index col) const
Definition: Transform.h:387
internal::cast_return_type< Transform, Transform< NewScalarType, Dim, Mode, Options > >::type cast() const
Definition: Transform.h:632
Matrix< Scalar, Dim, Dim, Options > LinearMatrixType
Definition: Transform.h:226
Transform inverse(TransformTraits traits=(TransformTraits) Mode) const
Definition: Transform.h:1230
TranslationPart translation()
Definition: Transform.h:410
void setIdentity()
Definition: Transform.h:531
RotationReturnType rotation() const
Definition: Transform.h:1086
Scalar_ Scalar
Definition: Transform.h:216
static const Transform Identity()
Returns an identity transformation.
Definition: Transform.h:537
EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE_FIXED_SIZE(Scalar_, Dim_==Dynamic ? Dynamic :(Dim_+1) *(Dim_+1)) enum
Definition: Transform.h:209
Transform & shear(const Scalar &sx, const Scalar &sy)
Definition: Transform.h:975
QTransform toQTransform(void) const
Definition: Transform.h:822
const MatrixType & matrix() const
Definition: Transform.h:393
ConstLinearPart linear() const
Definition: Transform.h:398
const Block< ConstMatrixType, Dim, 1,!(internal::traits< MatrixType >::Flags &RowMajorBit)> ConstTranslationPart
Definition: Transform.h:244
Represents a translation transformation.
Definition: Translation.h:33
TransformTraits
Definition: Constants.h:457
@ DontAlign
Definition: Constants.h:327
@ RowMajor
Definition: Constants.h:323
@ ComputeFullV
Definition: Constants.h:399
@ ComputeFullU
Definition: Constants.h:395
@ Affine
Definition: Constants.h:462
@ Projective
Definition: Constants.h:466
@ AffineCompact
Definition: Constants.h:464
@ Isometry
Definition: Constants.h:459
const unsigned int RowMajorBit
Definition: Constants.h:68
Namespace containing all symbols from the Eigen library.
Definition: B01_Experimental.dox:1
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:59
const int Dynamic
Definition: Constants.h:24
Definition: EigenBase.h:32
Derived & derived()
Definition: EigenBase.h:48
Holds information about the various numeric (i.e. scalar) types allowed by Eigen.
Definition: NumTraits.h:235