Eigen  3.3.0
 
Loading...
Searching...
No Matches
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
15namespace Eigen {
16
17namespace internal {
18
19template<typename Transform>
20struct transform_traits
21{
22 enum
23 {
24 Dim = Transform::Dim,
25 HDim = Transform::HDim,
26 Mode = Transform::Mode,
27 IsProjective = (int(Mode)==int(Projective))
28 };
29};
30
31template< typename TransformType,
32 typename MatrixType,
33 int Case = transform_traits<TransformType>::IsProjective ? 0
34 : int(MatrixType::RowsAtCompileTime) == int(transform_traits<TransformType>::HDim) ? 1
35 : 2,
36 int RhsCols = MatrixType::ColsAtCompileTime>
37struct transform_right_product_impl;
38
39template< typename Other,
40 int Mode,
41 int Options,
42 int Dim,
43 int HDim,
44 int OtherRows=Other::RowsAtCompileTime,
45 int OtherCols=Other::ColsAtCompileTime>
46struct transform_left_product_impl;
47
48template< typename Lhs,
49 typename Rhs,
50 bool AnyProjective =
51 transform_traits<Lhs>::IsProjective ||
52 transform_traits<Rhs>::IsProjective>
53struct transform_transform_product_impl;
54
55template< typename Other,
56 int Mode,
57 int Options,
58 int Dim,
59 int HDim,
60 int OtherRows=Other::RowsAtCompileTime,
61 int OtherCols=Other::ColsAtCompileTime>
62struct transform_construct_from_matrix;
63
64template<typename TransformType> struct transform_take_affine_part;
65
66template<typename _Scalar, int _Dim, int _Mode, int _Options>
67struct traits<Transform<_Scalar,_Dim,_Mode,_Options> >
68{
69 typedef _Scalar Scalar;
70 typedef Eigen::Index StorageIndex;
71 typedef Dense StorageKind;
72 enum {
73 Dim1 = _Dim==Dynamic ? _Dim : _Dim + 1,
74 RowsAtCompileTime = _Mode==Projective ? Dim1 : _Dim,
75 ColsAtCompileTime = Dim1,
76 MaxRowsAtCompileTime = RowsAtCompileTime,
77 MaxColsAtCompileTime = ColsAtCompileTime,
78 Flags = 0
79 };
80};
81
82template<int Mode> struct transform_make_affine;
83
84} // end namespace internal
85
200template<typename _Scalar, int _Dim, int _Mode, int _Options>
202{
203public:
205 enum {
206 Mode = _Mode,
207 Options = _Options,
208 Dim = _Dim,
209 HDim = _Dim+1,
210 Rows = int(Mode)==(AffineCompact) ? Dim : HDim
211 };
213 typedef _Scalar Scalar;
214 typedef Eigen::Index StorageIndex;
223 typedef Block<MatrixType,Dim,Dim,int(Mode)==(AffineCompact) && (Options&RowMajor)==0> LinearPart;
225 typedef const Block<ConstMatrixType,Dim,Dim,int(Mode)==(AffineCompact) && (Options&RowMajor)==0> ConstLinearPart;
227 typedef typename internal::conditional<int(Mode)==int(AffineCompact),
228 MatrixType&,
231 typedef typename internal::conditional<int(Mode)==int(AffineCompact),
232 const MatrixType&,
237 typedef Block<MatrixType,Dim,1,!(internal::traits<MatrixType>::Flags & RowMajorBit)> TranslationPart;
239 typedef const Block<ConstMatrixType,Dim,1,!(internal::traits<MatrixType>::Flags & RowMajorBit)> ConstTranslationPart;
242
243 // this intermediate enum is needed to avoid an ICE with gcc 3.4 and 4.0
244 enum { TransformTimeDiagonalMode = ((Mode==int(Isometry))?Affine:int(Mode)) };
247
248protected:
249
250 MatrixType m_matrix;
251
252public:
253
256 EIGEN_DEVICE_FUNC inline Transform()
257 {
258 check_template_params();
259 internal::transform_make_affine<(int(Mode)==Affine) ? Affine : AffineCompact>::run(m_matrix);
260 }
261
262 EIGEN_DEVICE_FUNC inline Transform(const Transform& other)
263 {
264 check_template_params();
265 m_matrix = other.m_matrix;
266 }
267
268 EIGEN_DEVICE_FUNC inline explicit Transform(const TranslationType& t)
269 {
270 check_template_params();
271 *this = t;
272 }
273 EIGEN_DEVICE_FUNC inline explicit Transform(const UniformScaling<Scalar>& s)
274 {
275 check_template_params();
276 *this = s;
277 }
278 template<typename Derived>
279 EIGEN_DEVICE_FUNC inline explicit Transform(const RotationBase<Derived, Dim>& r)
280 {
281 check_template_params();
282 *this = r;
283 }
284
285 EIGEN_DEVICE_FUNC inline Transform& operator=(const Transform& other)
286 { m_matrix = other.m_matrix; return *this; }
287
288 typedef internal::transform_take_affine_part<Transform> take_affine_part;
289
291 template<typename OtherDerived>
292 EIGEN_DEVICE_FUNC inline explicit Transform(const EigenBase<OtherDerived>& other)
293 {
294 EIGEN_STATIC_ASSERT((internal::is_same<Scalar,typename OtherDerived::Scalar>::value),
295 YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY);
296
297 check_template_params();
298 internal::transform_construct_from_matrix<OtherDerived,Mode,Options,Dim,HDim>::run(this, other.derived());
299 }
300
302 template<typename OtherDerived>
303 EIGEN_DEVICE_FUNC inline Transform& operator=(const EigenBase<OtherDerived>& other)
304 {
305 EIGEN_STATIC_ASSERT((internal::is_same<Scalar,typename OtherDerived::Scalar>::value),
306 YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY);
307
308 internal::transform_construct_from_matrix<OtherDerived,Mode,Options,Dim,HDim>::run(this, other.derived());
309 return *this;
310 }
311
312 template<int OtherOptions>
313 EIGEN_DEVICE_FUNC inline Transform(const Transform<Scalar,Dim,Mode,OtherOptions>& other)
314 {
315 check_template_params();
316 // only the options change, we can directly copy the matrices
317 m_matrix = other.matrix();
318 }
319
320 template<int OtherMode,int OtherOptions>
321 EIGEN_DEVICE_FUNC inline Transform(const Transform<Scalar,Dim,OtherMode,OtherOptions>& other)
322 {
323 check_template_params();
324 // prevent conversions as:
325 // Affine | AffineCompact | Isometry = Projective
326 EIGEN_STATIC_ASSERT(EIGEN_IMPLIES(OtherMode==int(Projective), Mode==int(Projective)),
327 YOU_PERFORMED_AN_INVALID_TRANSFORMATION_CONVERSION)
328
329 // prevent conversions as:
330 // Isometry = Affine | AffineCompact
331 EIGEN_STATIC_ASSERT(EIGEN_IMPLIES(OtherMode==int(Affine)||OtherMode==int(AffineCompact), Mode!=int(Isometry)),
332 YOU_PERFORMED_AN_INVALID_TRANSFORMATION_CONVERSION)
333
334 enum { ModeIsAffineCompact = Mode == int(AffineCompact),
335 OtherModeIsAffineCompact = OtherMode == int(AffineCompact)
336 };
337
338 if(ModeIsAffineCompact == OtherModeIsAffineCompact)
339 {
340 // We need the block expression because the code is compiled for all
341 // combinations of transformations and will trigger a compile time error
342 // if one tries to assign the matrices directly
343 m_matrix.template block<Dim,Dim+1>(0,0) = other.matrix().template block<Dim,Dim+1>(0,0);
344 makeAffine();
345 }
346 else if(OtherModeIsAffineCompact)
347 {
348 typedef typename Transform<Scalar,Dim,OtherMode,OtherOptions>::MatrixType OtherMatrixType;
349 internal::transform_construct_from_matrix<OtherMatrixType,Mode,Options,Dim,HDim>::run(this, other.matrix());
350 }
351 else
352 {
353 // here we know that Mode == AffineCompact and OtherMode != AffineCompact.
354 // if OtherMode were Projective, the static assert above would already have caught it.
355 // So the only possibility is that OtherMode == Affine
356 linear() = other.linear();
357 translation() = other.translation();
358 }
359 }
360
361 template<typename OtherDerived>
362 EIGEN_DEVICE_FUNC Transform(const ReturnByValue<OtherDerived>& other)
363 {
364 check_template_params();
365 other.evalTo(*this);
366 }
367
368 template<typename OtherDerived>
369 EIGEN_DEVICE_FUNC Transform& operator=(const ReturnByValue<OtherDerived>& other)
370 {
371 other.evalTo(*this);
372 return *this;
373 }
374
375 #ifdef EIGEN_QT_SUPPORT
376 inline Transform(const QMatrix& other);
377 inline Transform& operator=(const QMatrix& other);
378 inline QMatrix toQMatrix(void) const;
379 inline Transform(const QTransform& other);
380 inline Transform& operator=(const QTransform& other);
381 inline QTransform toQTransform(void) const;
382 #endif
383
384 EIGEN_DEVICE_FUNC Index rows() const { return int(Mode)==int(Projective) ? m_matrix.cols() : (m_matrix.cols()-1); }
385 EIGEN_DEVICE_FUNC Index cols() const { return m_matrix.cols(); }
386
389 EIGEN_DEVICE_FUNC inline Scalar operator() (Index row, Index col) const { return m_matrix(row,col); }
392 EIGEN_DEVICE_FUNC inline Scalar& operator() (Index row, Index col) { return m_matrix(row,col); }
393
395 EIGEN_DEVICE_FUNC inline const MatrixType& matrix() const { return m_matrix; }
397 EIGEN_DEVICE_FUNC inline MatrixType& matrix() { return m_matrix; }
398
400 EIGEN_DEVICE_FUNC inline ConstLinearPart linear() const { return ConstLinearPart(m_matrix,0,0); }
402 EIGEN_DEVICE_FUNC inline LinearPart linear() { return LinearPart(m_matrix,0,0); }
403
405 EIGEN_DEVICE_FUNC inline ConstAffinePart affine() const { return take_affine_part::run(m_matrix); }
407 EIGEN_DEVICE_FUNC inline AffinePart affine() { return take_affine_part::run(m_matrix); }
408
410 EIGEN_DEVICE_FUNC inline ConstTranslationPart translation() const { return ConstTranslationPart(m_matrix,0,Dim); }
412 EIGEN_DEVICE_FUNC inline TranslationPart translation() { return TranslationPart(m_matrix,0,Dim); }
413
438 // note: this function is defined here because some compilers cannot find the respective declaration
439 template<typename OtherDerived>
440 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const typename internal::transform_right_product_impl<Transform, OtherDerived>::ResultType
442 { return internal::transform_right_product_impl<Transform, OtherDerived>::run(*this,other.derived()); }
443
451 template<typename OtherDerived> friend
452 EIGEN_DEVICE_FUNC inline const typename internal::transform_left_product_impl<OtherDerived,Mode,Options,_Dim,_Dim+1>::ResultType
454 { return internal::transform_left_product_impl<OtherDerived,Mode,Options,Dim,HDim>::run(a.derived(),b); }
455
462 template<typename DiagonalDerived>
463 EIGEN_DEVICE_FUNC inline const TransformTimeDiagonalReturnType
464 operator * (const DiagonalBase<DiagonalDerived> &b) const
465 {
467 res.linearExt() *= b;
468 return res;
469 }
470
477 template<typename DiagonalDerived>
478 EIGEN_DEVICE_FUNC friend inline TransformTimeDiagonalReturnType
479 operator * (const DiagonalBase<DiagonalDerived> &a, const Transform &b)
480 {
482 res.linear().noalias() = a*b.linear();
483 res.translation().noalias() = a*b.translation();
484 if (Mode!=int(AffineCompact))
485 res.matrix().row(Dim) = b.matrix().row(Dim);
486 return res;
487 }
488
489 template<typename OtherDerived>
490 EIGEN_DEVICE_FUNC inline Transform& operator*=(const EigenBase<OtherDerived>& other) { return *this = *this * other; }
491
493 EIGEN_DEVICE_FUNC inline const Transform operator * (const Transform& other) const
494 {
495 return internal::transform_transform_product_impl<Transform,Transform>::run(*this,other);
496 }
497
498 #if EIGEN_COMP_ICC
499private:
500 // this intermediate structure permits to workaround a bug in ICC 11:
501 // error: template instantiation resulted in unexpected function type of "Eigen::Transform<double, 3, 32, 0>
502 // (const Eigen::Transform<double, 3, 2, 0> &) const"
503 // (the meaning of a name may have changed since the template declaration -- the type of the template is:
504 // "Eigen::internal::transform_transform_product_impl<Eigen::Transform<double, 3, 32, 0>,
505 // Eigen::Transform<double, 3, Mode, Options>, <expression>>::ResultType (const Eigen::Transform<double, 3, Mode, Options> &) const")
506 //
507 template<int OtherMode,int OtherOptions> struct icc_11_workaround
508 {
509 typedef internal::transform_transform_product_impl<Transform,Transform<Scalar,Dim,OtherMode,OtherOptions> > ProductType;
510 typedef typename ProductType::ResultType ResultType;
511 };
512
513public:
515 template<int OtherMode,int OtherOptions>
516 inline typename icc_11_workaround<OtherMode,OtherOptions>::ResultType
517 operator * (const Transform<Scalar,Dim,OtherMode,OtherOptions>& other) const
518 {
519 typedef typename icc_11_workaround<OtherMode,OtherOptions>::ProductType ProductType;
520 return ProductType::run(*this,other);
521 }
522 #else
524 template<int OtherMode,int OtherOptions>
525 EIGEN_DEVICE_FUNC inline typename internal::transform_transform_product_impl<Transform,Transform<Scalar,Dim,OtherMode,OtherOptions> >::ResultType
527 {
528 return internal::transform_transform_product_impl<Transform,Transform<Scalar,Dim,OtherMode,OtherOptions> >::run(*this,other);
529 }
530 #endif
531
533 EIGEN_DEVICE_FUNC void setIdentity() { m_matrix.setIdentity(); }
534
539 EIGEN_DEVICE_FUNC static const Transform Identity()
540 {
541 return Transform(MatrixType::Identity());
542 }
543
544 template<typename OtherDerived>
545 EIGEN_DEVICE_FUNC
546 inline Transform& scale(const MatrixBase<OtherDerived> &other);
547
548 template<typename OtherDerived>
549 EIGEN_DEVICE_FUNC
550 inline Transform& prescale(const MatrixBase<OtherDerived> &other);
551
552 EIGEN_DEVICE_FUNC inline Transform& scale(const Scalar& s);
553 EIGEN_DEVICE_FUNC inline Transform& prescale(const Scalar& s);
554
555 template<typename OtherDerived>
556 EIGEN_DEVICE_FUNC
557 inline Transform& translate(const MatrixBase<OtherDerived> &other);
558
559 template<typename OtherDerived>
560 EIGEN_DEVICE_FUNC
561 inline Transform& pretranslate(const MatrixBase<OtherDerived> &other);
562
563 template<typename RotationType>
564 EIGEN_DEVICE_FUNC
565 inline Transform& rotate(const RotationType& rotation);
566
567 template<typename RotationType>
568 EIGEN_DEVICE_FUNC
569 inline Transform& prerotate(const RotationType& rotation);
570
571 EIGEN_DEVICE_FUNC Transform& shear(const Scalar& sx, const Scalar& sy);
572 EIGEN_DEVICE_FUNC Transform& preshear(const Scalar& sx, const Scalar& sy);
573
574 EIGEN_DEVICE_FUNC inline Transform& operator=(const TranslationType& t);
575
576 EIGEN_DEVICE_FUNC
577 inline Transform& operator*=(const TranslationType& t) { return translate(t.vector()); }
578
579 EIGEN_DEVICE_FUNC inline Transform operator*(const TranslationType& t) const;
580
581 EIGEN_DEVICE_FUNC
582 inline Transform& operator=(const UniformScaling<Scalar>& t);
583
584 EIGEN_DEVICE_FUNC
585 inline Transform& operator*=(const UniformScaling<Scalar>& s) { return scale(s.factor()); }
586
587 EIGEN_DEVICE_FUNC
588 inline TransformTimeDiagonalReturnType operator*(const UniformScaling<Scalar>& s) const
589 {
591 res.scale(s.factor());
592 return res;
593 }
594
595 EIGEN_DEVICE_FUNC
596 inline Transform& operator*=(const DiagonalMatrix<Scalar,Dim>& s) { linearExt() *= s; return *this; }
597
598 template<typename Derived>
599 EIGEN_DEVICE_FUNC inline Transform& operator=(const RotationBase<Derived,Dim>& r);
600 template<typename Derived>
601 EIGEN_DEVICE_FUNC inline Transform& operator*=(const RotationBase<Derived,Dim>& r) { return rotate(r.toRotationMatrix()); }
602 template<typename Derived>
603 EIGEN_DEVICE_FUNC inline Transform operator*(const RotationBase<Derived,Dim>& r) const;
604
605 EIGEN_DEVICE_FUNC const LinearMatrixType rotation() const;
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
743template<typename Scalar, int Dim, int Mode,int Options>
745{
746 check_template_params();
747 *this = other;
748}
749
754template<typename Scalar, int Dim, int Mode,int Options>
756{
757 EIGEN_STATIC_ASSERT(Dim==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
758 if (Mode == int(AffineCompact))
759 m_matrix << other.m11(), other.m21(), other.dx(),
760 other.m12(), other.m22(), other.dy();
761 else
762 m_matrix << other.m11(), other.m21(), other.dx(),
763 other.m12(), other.m22(), other.dy(),
764 0, 0, 1;
765 return *this;
766}
767
774template<typename Scalar, int Dim, int Mode, int Options>
776{
777 check_template_params();
778 EIGEN_STATIC_ASSERT(Dim==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
779 return QMatrix(m_matrix.coeff(0,0), m_matrix.coeff(1,0),
780 m_matrix.coeff(0,1), m_matrix.coeff(1,1),
781 m_matrix.coeff(0,2), m_matrix.coeff(1,2));
782}
783
788template<typename Scalar, int Dim, int Mode,int Options>
790{
791 check_template_params();
792 *this = other;
793}
794
799template<typename Scalar, int Dim, int Mode, int Options>
801{
802 check_template_params();
803 EIGEN_STATIC_ASSERT(Dim==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
804 if (Mode == int(AffineCompact))
805 m_matrix << other.m11(), other.m21(), other.dx(),
806 other.m12(), other.m22(), other.dy();
807 else
808 m_matrix << other.m11(), other.m21(), other.dx(),
809 other.m12(), other.m22(), other.dy(),
810 other.m13(), other.m23(), other.m33();
811 return *this;
812}
813
818template<typename Scalar, int Dim, int Mode, int Options>
820{
821 EIGEN_STATIC_ASSERT(Dim==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
822 if (Mode == int(AffineCompact))
823 return QTransform(m_matrix.coeff(0,0), m_matrix.coeff(1,0),
824 m_matrix.coeff(0,1), m_matrix.coeff(1,1),
825 m_matrix.coeff(0,2), m_matrix.coeff(1,2));
826 else
827 return QTransform(m_matrix.coeff(0,0), m_matrix.coeff(1,0), m_matrix.coeff(2,0),
828 m_matrix.coeff(0,1), m_matrix.coeff(1,1), m_matrix.coeff(2,1),
829 m_matrix.coeff(0,2), m_matrix.coeff(1,2), m_matrix.coeff(2,2));
830}
831#endif
832
833/*********************
834*** Procedural API ***
835*********************/
836
841template<typename Scalar, int Dim, int Mode, int Options>
842template<typename OtherDerived>
843EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
845{
846 EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,int(Dim))
847 EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
848 linearExt().noalias() = (linearExt() * other.asDiagonal());
849 return *this;
850}
851
856template<typename Scalar, int Dim, int Mode, int Options>
858{
859 EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
860 linearExt() *= s;
861 return *this;
862}
863
868template<typename Scalar, int Dim, int Mode, int Options>
869template<typename OtherDerived>
870EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
872{
873 EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,int(Dim))
874 EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
875 affine().noalias() = (other.asDiagonal() * affine());
876 return *this;
877}
878
883template<typename Scalar, int Dim, int Mode, int Options>
885{
886 EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
887 m_matrix.template topRows<Dim>() *= s;
888 return *this;
889}
890
895template<typename Scalar, int Dim, int Mode, int Options>
896template<typename OtherDerived>
897EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
899{
900 EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,int(Dim))
901 translationExt() += linearExt() * other;
902 return *this;
903}
904
909template<typename Scalar, int Dim, int Mode, int Options>
910template<typename OtherDerived>
911EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
913{
914 EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,int(Dim))
915 if(int(Mode)==int(Projective))
916 affine() += other * m_matrix.row(Dim);
917 else
918 translation() += other;
919 return *this;
920}
921
939template<typename Scalar, int Dim, int Mode, int Options>
940template<typename RotationType>
941EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
943{
944 linearExt() *= internal::toRotationMatrix<Scalar,Dim>(rotation);
945 return *this;
946}
947
955template<typename Scalar, int Dim, int Mode, int Options>
956template<typename RotationType>
957EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
959{
960 m_matrix.template block<Dim,HDim>(0,0) = internal::toRotationMatrix<Scalar,Dim>(rotation)
961 * m_matrix.template block<Dim,HDim>(0,0);
962 return *this;
963}
964
970template<typename Scalar, int Dim, int Mode, int Options>
971EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
973{
974 EIGEN_STATIC_ASSERT(int(Dim)==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
975 EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
976 VectorType tmp = linear().col(0)*sy + linear().col(1);
977 linear() << linear().col(0) + linear().col(1)*sx, tmp;
978 return *this;
979}
980
986template<typename Scalar, int Dim, int Mode, int Options>
987EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
989{
990 EIGEN_STATIC_ASSERT(int(Dim)==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
991 EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
992 m_matrix.template block<Dim,HDim>(0,0) = LinearMatrixType(1, sx, sy, 1) * m_matrix.template block<Dim,HDim>(0,0);
993 return *this;
994}
995
996/******************************************************
997*** Scaling, Translation and Rotation compatibility ***
998******************************************************/
999
1000template<typename Scalar, int Dim, int Mode, int Options>
1001EIGEN_DEVICE_FUNC inline Transform<Scalar,Dim,Mode,Options>& Transform<Scalar,Dim,Mode,Options>::operator=(const TranslationType& t)
1002{
1003 linear().setIdentity();
1004 translation() = t.vector();
1005 makeAffine();
1006 return *this;
1007}
1008
1009template<typename Scalar, int Dim, int Mode, int Options>
1010EIGEN_DEVICE_FUNC inline Transform<Scalar,Dim,Mode,Options> Transform<Scalar,Dim,Mode,Options>::operator*(const TranslationType& t) const
1011{
1012 Transform res = *this;
1013 res.translate(t.vector());
1014 return res;
1015}
1016
1017template<typename Scalar, int Dim, int Mode, int Options>
1018EIGEN_DEVICE_FUNC inline Transform<Scalar,Dim,Mode,Options>& Transform<Scalar,Dim,Mode,Options>::operator=(const UniformScaling<Scalar>& s)
1019{
1020 m_matrix.setZero();
1021 linear().diagonal().fill(s.factor());
1022 makeAffine();
1023 return *this;
1024}
1025
1026template<typename Scalar, int Dim, int Mode, int Options>
1027template<typename Derived>
1028EIGEN_DEVICE_FUNC inline Transform<Scalar,Dim,Mode,Options>& Transform<Scalar,Dim,Mode,Options>::operator=(const RotationBase<Derived,Dim>& r)
1029{
1030 linear() = internal::toRotationMatrix<Scalar,Dim>(r);
1031 translation().setZero();
1032 makeAffine();
1033 return *this;
1034}
1035
1036template<typename Scalar, int Dim, int Mode, int Options>
1037template<typename Derived>
1038EIGEN_DEVICE_FUNC inline Transform<Scalar,Dim,Mode,Options> Transform<Scalar,Dim,Mode,Options>::operator*(const RotationBase<Derived,Dim>& r) const
1039{
1040 Transform res = *this;
1041 res.rotate(r.derived());
1042 return res;
1043}
1044
1045/************************
1046*** Special functions ***
1047************************/
1048
1056template<typename Scalar, int Dim, int Mode, int Options>
1057EIGEN_DEVICE_FUNC const typename Transform<Scalar,Dim,Mode,Options>::LinearMatrixType
1059{
1060 LinearMatrixType result;
1061 computeRotationScaling(&result, (LinearMatrixType*)0);
1062 return result;
1063}
1064
1065
1077template<typename Scalar, int Dim, int Mode, int Options>
1078template<typename RotationMatrixType, typename ScalingMatrixType>
1079EIGEN_DEVICE_FUNC void Transform<Scalar,Dim,Mode,Options>::computeRotationScaling(RotationMatrixType *rotation, ScalingMatrixType *scaling) const
1080{
1082
1083 Scalar x = (svd.matrixU() * svd.matrixV().adjoint()).determinant(); // so x has absolute value 1
1084 VectorType sv(svd.singularValues());
1085 sv.coeffRef(0) *= x;
1086 if(scaling) scaling->lazyAssign(svd.matrixV() * sv.asDiagonal() * svd.matrixV().adjoint());
1087 if(rotation)
1088 {
1089 LinearMatrixType m(svd.matrixU());
1090 m.col(0) /= x;
1091 rotation->lazyAssign(m * svd.matrixV().adjoint());
1092 }
1093}
1094
1106template<typename Scalar, int Dim, int Mode, int Options>
1107template<typename ScalingMatrixType, typename RotationMatrixType>
1108EIGEN_DEVICE_FUNC void Transform<Scalar,Dim,Mode,Options>::computeScalingRotation(ScalingMatrixType *scaling, RotationMatrixType *rotation) const
1109{
1111
1112 Scalar x = (svd.matrixU() * svd.matrixV().adjoint()).determinant(); // so x has absolute value 1
1113 VectorType sv(svd.singularValues());
1114 sv.coeffRef(0) *= x;
1115 if(scaling) scaling->lazyAssign(svd.matrixU() * sv.asDiagonal() * svd.matrixU().adjoint());
1116 if(rotation)
1117 {
1118 LinearMatrixType m(svd.matrixU());
1119 m.col(0) /= x;
1120 rotation->lazyAssign(m * svd.matrixV().adjoint());
1121 }
1122}
1123
1127template<typename Scalar, int Dim, int Mode, int Options>
1128template<typename PositionDerived, typename OrientationType, typename ScaleDerived>
1129EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>&
1131 const OrientationType& orientation, const MatrixBase<ScaleDerived> &scale)
1132{
1133 linear() = internal::toRotationMatrix<Scalar,Dim>(orientation);
1134 linear() *= scale.asDiagonal();
1135 translation() = position;
1136 makeAffine();
1137 return *this;
1138}
1139
1140namespace internal {
1141
1142template<int Mode>
1143struct transform_make_affine
1144{
1145 template<typename MatrixType>
1146 EIGEN_DEVICE_FUNC static void run(MatrixType &mat)
1147 {
1148 static const int Dim = MatrixType::ColsAtCompileTime-1;
1149 mat.template block<1,Dim>(Dim,0).setZero();
1150 mat.coeffRef(Dim,Dim) = typename MatrixType::Scalar(1);
1151 }
1152};
1153
1154template<>
1155struct transform_make_affine<AffineCompact>
1156{
1157 template<typename MatrixType> EIGEN_DEVICE_FUNC static void run(MatrixType &) { }
1158};
1159
1160// selector needed to avoid taking the inverse of a 3x4 matrix
1161template<typename TransformType, int Mode=TransformType::Mode>
1162struct projective_transform_inverse
1163{
1164 EIGEN_DEVICE_FUNC static inline void run(const TransformType&, TransformType&)
1165 {}
1166};
1167
1168template<typename TransformType>
1169struct projective_transform_inverse<TransformType, Projective>
1170{
1171 EIGEN_DEVICE_FUNC static inline void run(const TransformType& m, TransformType& res)
1172 {
1173 res.matrix() = m.matrix().inverse();
1174 }
1175};
1176
1177} // end namespace internal
1178
1179
1200template<typename Scalar, int Dim, int Mode, int Options>
1201EIGEN_DEVICE_FUNC Transform<Scalar,Dim,Mode,Options>
1203{
1204 Transform res;
1205 if (hint == Projective)
1206 {
1207 internal::projective_transform_inverse<Transform>::run(*this, res);
1208 }
1209 else
1210 {
1211 if (hint == Isometry)
1212 {
1213 res.matrix().template topLeftCorner<Dim,Dim>() = linear().transpose();
1214 }
1215 else if(hint&Affine)
1216 {
1217 res.matrix().template topLeftCorner<Dim,Dim>() = linear().inverse();
1218 }
1219 else
1220 {
1221 eigen_assert(false && "Invalid transform traits in Transform::Inverse");
1222 }
1223 // translation and remaining parts
1224 res.matrix().template topRightCorner<Dim,1>()
1225 = - res.matrix().template topLeftCorner<Dim,Dim>() * translation();
1226 res.makeAffine(); // we do need this, because in the beginning res is uninitialized
1227 }
1228 return res;
1229}
1230
1231namespace internal {
1232
1233/*****************************************************
1234*** Specializations of take affine part ***
1235*****************************************************/
1236
1237template<typename TransformType> struct transform_take_affine_part {
1238 typedef typename TransformType::MatrixType MatrixType;
1239 typedef typename TransformType::AffinePart AffinePart;
1240 typedef typename TransformType::ConstAffinePart ConstAffinePart;
1241 static inline AffinePart run(MatrixType& m)
1242 { return m.template block<TransformType::Dim,TransformType::HDim>(0,0); }
1243 static inline ConstAffinePart run(const MatrixType& m)
1244 { return m.template block<TransformType::Dim,TransformType::HDim>(0,0); }
1245};
1246
1247template<typename Scalar, int Dim, int Options>
1248struct transform_take_affine_part<Transform<Scalar,Dim,AffineCompact, Options> > {
1250 static inline MatrixType& run(MatrixType& m) { return m; }
1251 static inline const MatrixType& run(const MatrixType& m) { return m; }
1252};
1253
1254/*****************************************************
1255*** Specializations of construct from matrix ***
1256*****************************************************/
1257
1258template<typename Other, int Mode, int Options, int Dim, int HDim>
1259struct transform_construct_from_matrix<Other, Mode,Options,Dim,HDim, Dim,Dim>
1260{
1261 static inline void run(Transform<typename Other::Scalar,Dim,Mode,Options> *transform, const Other& other)
1262 {
1263 transform->linear() = other;
1264 transform->translation().setZero();
1265 transform->makeAffine();
1266 }
1267};
1268
1269template<typename Other, int Mode, int Options, int Dim, int HDim>
1270struct transform_construct_from_matrix<Other, Mode,Options,Dim,HDim, Dim,HDim>
1271{
1272 static inline void run(Transform<typename Other::Scalar,Dim,Mode,Options> *transform, const Other& other)
1273 {
1274 transform->affine() = other;
1275 transform->makeAffine();
1276 }
1277};
1278
1279template<typename Other, int Mode, int Options, int Dim, int HDim>
1280struct transform_construct_from_matrix<Other, Mode,Options,Dim,HDim, HDim,HDim>
1281{
1282 static inline void run(Transform<typename Other::Scalar,Dim,Mode,Options> *transform, const Other& other)
1283 { transform->matrix() = other; }
1284};
1285
1286template<typename Other, int Options, int Dim, int HDim>
1287struct transform_construct_from_matrix<Other, AffineCompact,Options,Dim,HDim, HDim,HDim>
1288{
1289 static inline void run(Transform<typename Other::Scalar,Dim,AffineCompact,Options> *transform, const Other& other)
1290 { transform->matrix() = other.template block<Dim,HDim>(0,0); }
1291};
1292
1293/**********************************************************
1294*** Specializations of operator* with rhs EigenBase ***
1295**********************************************************/
1296
1297template<int LhsMode,int RhsMode>
1298struct transform_product_result
1299{
1300 enum
1301 {
1302 Mode =
1303 (LhsMode == (int)Projective || RhsMode == (int)Projective ) ? Projective :
1304 (LhsMode == (int)Affine || RhsMode == (int)Affine ) ? Affine :
1305 (LhsMode == (int)AffineCompact || RhsMode == (int)AffineCompact ) ? AffineCompact :
1306 (LhsMode == (int)Isometry || RhsMode == (int)Isometry ) ? Isometry : Projective
1307 };
1308};
1309
1310template< typename TransformType, typename MatrixType, int RhsCols>
1311struct transform_right_product_impl< TransformType, MatrixType, 0, RhsCols>
1312{
1313 typedef typename MatrixType::PlainObject ResultType;
1314
1315 static EIGEN_STRONG_INLINE ResultType run(const TransformType& T, const MatrixType& other)
1316 {
1317 return T.matrix() * other;
1318 }
1319};
1320
1321template< typename TransformType, typename MatrixType, int RhsCols>
1322struct transform_right_product_impl< TransformType, MatrixType, 1, RhsCols>
1323{
1324 enum {
1325 Dim = TransformType::Dim,
1326 HDim = TransformType::HDim,
1327 OtherRows = MatrixType::RowsAtCompileTime,
1328 OtherCols = MatrixType::ColsAtCompileTime
1329 };
1330
1331 typedef typename MatrixType::PlainObject ResultType;
1332
1333 static EIGEN_STRONG_INLINE ResultType run(const TransformType& T, const MatrixType& other)
1334 {
1335 EIGEN_STATIC_ASSERT(OtherRows==HDim, YOU_MIXED_MATRICES_OF_DIFFERENT_SIZES);
1336
1337 typedef Block<ResultType, Dim, OtherCols, int(MatrixType::RowsAtCompileTime)==Dim> TopLeftLhs;
1338
1339 ResultType res(other.rows(),other.cols());
1340 TopLeftLhs(res, 0, 0, Dim, other.cols()).noalias() = T.affine() * other;
1341 res.row(OtherRows-1) = other.row(OtherRows-1);
1342
1343 return res;
1344 }
1345};
1346
1347template< typename TransformType, typename MatrixType, int RhsCols>
1348struct transform_right_product_impl< TransformType, MatrixType, 2, RhsCols>
1349{
1350 enum {
1351 Dim = TransformType::Dim,
1352 HDim = TransformType::HDim,
1353 OtherRows = MatrixType::RowsAtCompileTime,
1354 OtherCols = MatrixType::ColsAtCompileTime
1355 };
1356
1357 typedef typename MatrixType::PlainObject ResultType;
1358
1359 static EIGEN_STRONG_INLINE ResultType run(const TransformType& T, const MatrixType& other)
1360 {
1361 EIGEN_STATIC_ASSERT(OtherRows==Dim, YOU_MIXED_MATRICES_OF_DIFFERENT_SIZES);
1362
1363 typedef Block<ResultType, Dim, OtherCols, true> TopLeftLhs;
1364 ResultType res(Replicate<typename TransformType::ConstTranslationPart, 1, OtherCols>(T.translation(),1,other.cols()));
1365 TopLeftLhs(res, 0, 0, Dim, other.cols()).noalias() += T.linear() * other;
1366
1367 return res;
1368 }
1369};
1370
1371template< typename TransformType, typename MatrixType >
1372struct transform_right_product_impl< TransformType, MatrixType, 2, 1> // rhs is a vector of size Dim
1373{
1374 typedef typename TransformType::MatrixType TransformMatrix;
1375 enum {
1376 Dim = TransformType::Dim,
1377 HDim = TransformType::HDim,
1378 OtherRows = MatrixType::RowsAtCompileTime,
1379 WorkingRows = EIGEN_PLAIN_ENUM_MIN(TransformMatrix::RowsAtCompileTime,HDim)
1380 };
1381
1382 typedef typename MatrixType::PlainObject ResultType;
1383
1384 static EIGEN_STRONG_INLINE ResultType run(const TransformType& T, const MatrixType& other)
1385 {
1386 EIGEN_STATIC_ASSERT(OtherRows==Dim, YOU_MIXED_MATRICES_OF_DIFFERENT_SIZES);
1387
1388 Matrix<typename ResultType::Scalar, Dim+1, 1> rhs;
1389 rhs.template head<Dim>() = other; rhs[Dim] = typename ResultType::Scalar(1);
1390 Matrix<typename ResultType::Scalar, WorkingRows, 1> res(T.matrix() * rhs);
1391 return res.template head<Dim>();
1392 }
1393};
1394
1395/**********************************************************
1396*** Specializations of operator* with lhs EigenBase ***
1397**********************************************************/
1398
1399// generic HDim x HDim matrix * T => Projective
1400template<typename Other,int Mode, int Options, int Dim, int HDim>
1401struct transform_left_product_impl<Other,Mode,Options,Dim,HDim, HDim,HDim>
1402{
1403 typedef Transform<typename Other::Scalar,Dim,Mode,Options> TransformType;
1404 typedef typename TransformType::MatrixType MatrixType;
1405 typedef Transform<typename Other::Scalar,Dim,Projective,Options> ResultType;
1406 static ResultType run(const Other& other,const TransformType& tr)
1407 { return ResultType(other * tr.matrix()); }
1408};
1409
1410// generic HDim x HDim matrix * AffineCompact => Projective
1411template<typename Other, int Options, int Dim, int HDim>
1412struct transform_left_product_impl<Other,AffineCompact,Options,Dim,HDim, HDim,HDim>
1413{
1414 typedef Transform<typename Other::Scalar,Dim,AffineCompact,Options> TransformType;
1415 typedef typename TransformType::MatrixType MatrixType;
1416 typedef Transform<typename Other::Scalar,Dim,Projective,Options> ResultType;
1417 static ResultType run(const Other& other,const TransformType& tr)
1418 {
1419 ResultType res;
1420 res.matrix().noalias() = other.template block<HDim,Dim>(0,0) * tr.matrix();
1421 res.matrix().col(Dim) += other.col(Dim);
1422 return res;
1423 }
1424};
1425
1426// affine matrix * T
1427template<typename Other,int Mode, int Options, int Dim, int HDim>
1428struct transform_left_product_impl<Other,Mode,Options,Dim,HDim, Dim,HDim>
1429{
1430 typedef Transform<typename Other::Scalar,Dim,Mode,Options> TransformType;
1431 typedef typename TransformType::MatrixType MatrixType;
1432 typedef TransformType ResultType;
1433 static ResultType run(const Other& other,const TransformType& tr)
1434 {
1435 ResultType res;
1436 res.affine().noalias() = other * tr.matrix();
1437 res.matrix().row(Dim) = tr.matrix().row(Dim);
1438 return res;
1439 }
1440};
1441
1442// affine matrix * AffineCompact
1443template<typename Other, int Options, int Dim, int HDim>
1444struct transform_left_product_impl<Other,AffineCompact,Options,Dim,HDim, Dim,HDim>
1445{
1446 typedef Transform<typename Other::Scalar,Dim,AffineCompact,Options> TransformType;
1447 typedef typename TransformType::MatrixType MatrixType;
1448 typedef TransformType ResultType;
1449 static ResultType run(const Other& other,const TransformType& tr)
1450 {
1451 ResultType res;
1452 res.matrix().noalias() = other.template block<Dim,Dim>(0,0) * tr.matrix();
1453 res.translation() += other.col(Dim);
1454 return res;
1455 }
1456};
1457
1458// linear matrix * T
1459template<typename Other,int Mode, int Options, int Dim, int HDim>
1460struct transform_left_product_impl<Other,Mode,Options,Dim,HDim, Dim,Dim>
1461{
1462 typedef Transform<typename Other::Scalar,Dim,Mode,Options> TransformType;
1463 typedef typename TransformType::MatrixType MatrixType;
1464 typedef TransformType ResultType;
1465 static ResultType run(const Other& other, const TransformType& tr)
1466 {
1467 TransformType res;
1468 if(Mode!=int(AffineCompact))
1469 res.matrix().row(Dim) = tr.matrix().row(Dim);
1470 res.matrix().template topRows<Dim>().noalias()
1471 = other * tr.matrix().template topRows<Dim>();
1472 return res;
1473 }
1474};
1475
1476/**********************************************************
1477*** Specializations of operator* with another Transform ***
1478**********************************************************/
1479
1480template<typename Scalar, int Dim, int LhsMode, int LhsOptions, int RhsMode, int RhsOptions>
1481struct transform_transform_product_impl<Transform<Scalar,Dim,LhsMode,LhsOptions>,Transform<Scalar,Dim,RhsMode,RhsOptions>,false >
1482{
1483 enum { ResultMode = transform_product_result<LhsMode,RhsMode>::Mode };
1484 typedef Transform<Scalar,Dim,LhsMode,LhsOptions> Lhs;
1485 typedef Transform<Scalar,Dim,RhsMode,RhsOptions> Rhs;
1486 typedef Transform<Scalar,Dim,ResultMode,LhsOptions> ResultType;
1487 static ResultType run(const Lhs& lhs, const Rhs& rhs)
1488 {
1489 ResultType res;
1490 res.linear() = lhs.linear() * rhs.linear();
1491 res.translation() = lhs.linear() * rhs.translation() + lhs.translation();
1492 res.makeAffine();
1493 return res;
1494 }
1495};
1496
1497template<typename Scalar, int Dim, int LhsMode, int LhsOptions, int RhsMode, int RhsOptions>
1498struct transform_transform_product_impl<Transform<Scalar,Dim,LhsMode,LhsOptions>,Transform<Scalar,Dim,RhsMode,RhsOptions>,true >
1499{
1500 typedef Transform<Scalar,Dim,LhsMode,LhsOptions> Lhs;
1501 typedef Transform<Scalar,Dim,RhsMode,RhsOptions> Rhs;
1502 typedef Transform<Scalar,Dim,Projective> ResultType;
1503 static ResultType run(const Lhs& lhs, const Rhs& rhs)
1504 {
1505 return ResultType( lhs.matrix() * rhs.matrix() );
1506 }
1507};
1508
1509template<typename Scalar, int Dim, int LhsOptions, int RhsOptions>
1510struct transform_transform_product_impl<Transform<Scalar,Dim,AffineCompact,LhsOptions>,Transform<Scalar,Dim,Projective,RhsOptions>,true >
1511{
1512 typedef Transform<Scalar,Dim,AffineCompact,LhsOptions> Lhs;
1513 typedef Transform<Scalar,Dim,Projective,RhsOptions> Rhs;
1514 typedef Transform<Scalar,Dim,Projective> ResultType;
1515 static ResultType run(const Lhs& lhs, const Rhs& rhs)
1516 {
1517 ResultType res;
1518 res.matrix().template topRows<Dim>() = lhs.matrix() * rhs.matrix();
1519 res.matrix().row(Dim) = rhs.matrix().row(Dim);
1520 return res;
1521 }
1522};
1523
1524template<typename Scalar, int Dim, int LhsOptions, int RhsOptions>
1525struct transform_transform_product_impl<Transform<Scalar,Dim,Projective,LhsOptions>,Transform<Scalar,Dim,AffineCompact,RhsOptions>,true >
1526{
1527 typedef Transform<Scalar,Dim,Projective,LhsOptions> Lhs;
1528 typedef Transform<Scalar,Dim,AffineCompact,RhsOptions> Rhs;
1529 typedef Transform<Scalar,Dim,Projective> ResultType;
1530 static ResultType run(const Lhs& lhs, const Rhs& rhs)
1531 {
1532 ResultType res(lhs.matrix().template leftCols<Dim>() * rhs.matrix());
1533 res.matrix().col(Dim) += lhs.matrix().col(Dim);
1534 return res;
1535 }
1536};
1537
1538} // end namespace internal
1539
1540} // end namespace Eigen
1541
1542#endif // EIGEN_TRANSFORM_H
Expression of a fixed-size or dynamic-size block.
Definition: Block.h:105
Two-sided Jacobi SVD decomposition of a rectangular matrix.
Definition: JacobiSVD.h:485
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:50
const DiagonalWrapper< const Derived > asDiagonal() const
Definition: DiagonalMatrix.h:277
The matrix class, also used for vectors and row-vectors.
Definition: Matrix.h:180
Scalar & coeffRef(Index rowId, Index colId)
Definition: PlainObjectBase.h:177
Derived & lazyAssign(const DenseBase< OtherDerived > &other)
Definition: PlainObjectBase.h:459
const Scalar & coeff(Index rowId, Index colId) const
Definition: PlainObjectBase.h:154
Derived & setZero(Index size)
Definition: CwiseNullaryOp.h:515
const Scalar * data() const
Definition: PlainObjectBase.h:249
const SingularValuesType & singularValues() const
Definition: SVDBase.h:111
const MatrixUType & matrixU() const
Definition: SVDBase.h:83
const MatrixVType & matrixV() const
Definition: SVDBase.h:99
Represents an homogeneous transformation in a N dimensional space.
Definition: Transform.h:202
const MatrixType & matrix() const
Definition: Transform.h:395
Block< MatrixType, Dim, Dim, int(Mode)==(AffineCompact) &&(Options &RowMajor)==0 > LinearPart
Definition: Transform.h:223
const LinearMatrixType rotation() const
Definition: Transform.h:1058
Transform & preshear(const Scalar &sx, const Scalar &sy)
Definition: Transform.h:988
Transform()
Definition: Transform.h:256
Transform & operator=(const EigenBase< OtherDerived > &other)
Definition: Transform.h:303
void computeScalingRotation(ScalingMatrixType *scaling, RotationMatrixType *rotation) const
Definition: Transform.h:1108
Transform(const Transform< OtherScalarType, Dim, Mode, Options > &other)
Definition: Transform.h:637
MatrixType & matrix()
Definition: Transform.h:397
EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF_VECTORIZABLE_FIXED_SIZE(_Scalar, _Dim==Dynamic ? Dynamic :(_Dim+1) *(_Dim+1)) enum
Definition: Transform.h:204
internal::make_proper_matrix_type< Scalar, Rows, HDim, Options >::type MatrixType
Definition: Transform.h:217
const Block< ConstMatrixType, Dim, 1,!(internal::traits< MatrixType >::Flags &RowMajorBit)> ConstTranslationPart
Definition: Transform.h:239
Matrix< Scalar, Dim, Dim, Options > LinearMatrixType
Definition: Transform.h:221
Eigen::Index Index
Definition: Transform.h:215
_Scalar Scalar
Definition: Transform.h:213
Transform< Scalar, Dim, TransformTimeDiagonalMode > TransformTimeDiagonalReturnType
Definition: Transform.h:246
static const Transform Identity()
Returns an identity transformation.
Definition: Transform.h:539
Scalar * data()
Definition: Transform.h:624
internal::cast_return_type< Transform, Transform< NewScalarType, Dim, Mode, Options > >::type cast() const
Definition: Transform.h:632
internal::conditional< int(Mode)==int(AffineCompact), MatrixType &, Block< MatrixType, Dim, HDim > >::type AffinePart
Definition: Transform.h:229
Translation< Scalar, Dim > TranslationType
Definition: Transform.h:241
internal::conditional< int(Mode)==int(AffineCompact), constMatrixType &, constBlock< constMatrixType, Dim, HDim > >::type ConstAffinePart
Definition: Transform.h:233
void computeRotationScaling(RotationMatrixType *rotation, ScalingMatrixType *scaling) const
Definition: Transform.h:1079
ConstAffinePart affine() const
Definition: Transform.h:405
Scalar operator()(Index row, Index col) const
Definition: Transform.h:389
AffinePart affine()
Definition: Transform.h:407
friend const internal::transform_left_product_impl< OtherDerived, Mode, Options, _Dim, _Dim+1 >::ResultType operator*(const EigenBase< OtherDerived > &a, const Transform &b)
Definition: Transform.h:453
TranslationPart translation()
Definition: Transform.h:412
Block< MatrixType, Dim, 1,!(internal::traits< MatrixType >::Flags &RowMajorBit)> TranslationPart
Definition: Transform.h:237
Matrix< Scalar, Dim, 1 > VectorType
Definition: Transform.h:235
Transform inverse(TransformTraits traits=(TransformTraits) Mode) const
Definition: Transform.h:1202
bool isApprox(const Transform &other, const typename NumTraits< Scalar >::Real &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: Transform.h:647
void makeAffine()
Definition: Transform.h:652
const Block< ConstMatrixType, Dim, Dim, int(Mode)==(AffineCompact) &&(Options &RowMajor)==0 > ConstLinearPart
Definition: Transform.h:225
void setIdentity()
Definition: Transform.h:533
LinearPart linear()
Definition: Transform.h:402
QMatrix toQMatrix(void) const
Definition: Transform.h:775
ConstTranslationPart translation() const
Definition: Transform.h:410
const Scalar * data() const
Definition: Transform.h:622
Transform & shear(const Scalar &sx, const Scalar &sy)
Definition: Transform.h:972
QTransform toQTransform(void) const
Definition: Transform.h:819
Transform(const EigenBase< OtherDerived > &other)
Definition: Transform.h:292
const MatrixType ConstMatrixType
Definition: Transform.h:219
ConstLinearPart linear() const
Definition: Transform.h:400
Represents a translation transformation.
Definition: Translation.h:31
TransformTraits
Definition: Constants.h:445
@ DontAlign
Definition: Constants.h:326
@ RowMajor
Definition: Constants.h:322
@ ComputeFullV
Definition: Constants.h:387
@ ComputeFullU
Definition: Constants.h:383
@ Affine
Definition: Constants.h:450
@ Projective
Definition: Constants.h:454
@ AffineCompact
Definition: Constants.h:452
@ Isometry
Definition: Constants.h:447
const unsigned int RowMajorBit
Definition: Constants.h:61
Namespace containing all symbols from the Eigen library.
Definition: Core:287
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:33
const int Dynamic
Definition: Constants.h:21
Definition: EigenBase.h:29
Derived & derived()
Definition: EigenBase.h:44
Holds information about the various numeric (i.e. scalar) types allowed by Eigen.
Definition: NumTraits.h:151