Transform.h
Go to the documentation of this file.
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 // Eigen is free software; you can redistribute it and/or
9 // modify it under the terms of the GNU Lesser General Public
10 // License as published by the Free Software Foundation; either
11 // version 3 of the License, or (at your option) any later version.
12 //
13 // Alternatively, you can redistribute it and/or
14 // modify it under the terms of the GNU General Public License as
15 // published by the Free Software Foundation; either version 2 of
16 // the License, or (at your option) any later version.
17 //
18 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
19 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
21 // GNU General Public License for more details.
22 //
23 // You should have received a copy of the GNU Lesser General Public
24 // License and a copy of the GNU General Public License along with
25 // Eigen. If not, see <http://www.gnu.org/licenses/>.
26 
27 #ifndef EIGEN_TRANSFORM_H
28 #define EIGEN_TRANSFORM_H
29 
30 namespace Eigen {
31 
32 namespace internal {
33 
34 template<typename Transform>
35 struct transform_traits
36 {
37  enum
38  {
39  Dim = Transform::Dim,
40  HDim = Transform::HDim,
41  Mode = Transform::Mode,
42  IsProjective = (Mode==Projective)
43  };
44 };
45 
46 template< typename TransformType,
47  typename MatrixType,
48  int Case = transform_traits<TransformType>::IsProjective ? 0
49  : int(MatrixType::RowsAtCompileTime) == int(transform_traits<TransformType>::HDim) ? 1
50  : 2>
51 struct transform_right_product_impl;
52 
53 template< typename Other,
54  int Mode,
55  int Options,
56  int Dim,
57  int HDim,
58  int OtherRows=Other::RowsAtCompileTime,
59  int OtherCols=Other::ColsAtCompileTime>
60 struct transform_left_product_impl;
61 
62 template< typename Lhs,
63  typename Rhs,
64  bool AnyProjective =
65  transform_traits<Lhs>::IsProjective ||
66  transform_traits<Rhs>::IsProjective>
67 struct transform_transform_product_impl;
68 
69 template< typename Other,
70  int Mode,
71  int Options,
72  int Dim,
73  int HDim,
74  int OtherRows=Other::RowsAtCompileTime,
75  int OtherCols=Other::ColsAtCompileTime>
76 struct transform_construct_from_matrix;
77 
78 template<typename TransformType> struct transform_take_affine_part;
79 
80 } // end namespace internal
81 
190 template<typename _Scalar, int _Dim, int _Mode, int _Options>
192 {
193 public:
195  enum {
196  Mode = _Mode,
197  Options = _Options,
198  Dim = _Dim,
199  HDim = _Dim+1,
200  Rows = int(Mode)==(AffineCompact) ? Dim : HDim
201  };
203  typedef _Scalar Scalar;
204  typedef DenseIndex Index;
208  typedef const MatrixType ConstMatrixType;
216  typedef typename internal::conditional<int(Mode)==int(AffineCompact),
217  MatrixType&,
220  typedef typename internal::conditional<int(Mode)==int(AffineCompact),
221  const MatrixType&,
231 
232  // this intermediate enum is needed to avoid an ICE with gcc 3.4 and 4.0
233  enum { TransformTimeDiagonalMode = ((Mode==int(Isometry))?Affine:int(Mode)) };
236 
237 protected:
238 
239  MatrixType m_matrix;
240 
241 public:
242 
245  inline Transform()
246  {
247  check_template_params();
248  if (int(Mode)==Affine)
249  makeAffine();
250  }
251 
252  inline Transform(const Transform& other)
253  {
254  check_template_params();
255  m_matrix = other.m_matrix;
256  }
257 
258  inline explicit Transform(const TranslationType& t)
259  {
260  check_template_params();
261  *this = t;
262  }
263  inline explicit Transform(const UniformScaling<Scalar>& s)
264  {
265  check_template_params();
266  *this = s;
267  }
268  template<typename Derived>
269  inline explicit Transform(const RotationBase<Derived, Dim>& r)
270  {
271  check_template_params();
272  *this = r;
273  }
274 
275  inline Transform& operator=(const Transform& other)
276  { m_matrix = other.m_matrix; return *this; }
277 
278  typedef internal::transform_take_affine_part<Transform> take_affine_part;
279 
281  template<typename OtherDerived>
282  inline explicit Transform(const EigenBase<OtherDerived>& other)
283  {
284  EIGEN_STATIC_ASSERT((internal::is_same<Scalar,typename OtherDerived::Scalar>::value),
285  YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY);
286 
287  check_template_params();
288  internal::transform_construct_from_matrix<OtherDerived,Mode,Options,Dim,HDim>::run(this, other.derived());
289  }
290 
292  template<typename OtherDerived>
293  inline Transform& operator=(const EigenBase<OtherDerived>& other)
294  {
295  EIGEN_STATIC_ASSERT((internal::is_same<Scalar,typename OtherDerived::Scalar>::value),
296  YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY);
297 
298  internal::transform_construct_from_matrix<OtherDerived,Mode,Options,Dim,HDim>::run(this, other.derived());
299  return *this;
300  }
301 
302  template<int OtherOptions>
304  {
305  check_template_params();
306  // only the options change, we can directly copy the matrices
307  m_matrix = other.matrix();
308  }
309 
310  template<int OtherMode,int OtherOptions>
312  {
313  check_template_params();
314  // prevent conversions as:
315  // Affine | AffineCompact | Isometry = Projective
316  EIGEN_STATIC_ASSERT(EIGEN_IMPLIES(OtherMode==int(Projective), Mode==int(Projective)),
317  YOU_PERFORMED_AN_INVALID_TRANSFORMATION_CONVERSION)
318 
319  // prevent conversions as:
320  // Isometry = Affine | AffineCompact
321  EIGEN_STATIC_ASSERT(EIGEN_IMPLIES(OtherMode==int(Affine)||OtherMode==int(AffineCompact), Mode!=int(Isometry)),
322  YOU_PERFORMED_AN_INVALID_TRANSFORMATION_CONVERSION)
323 
324  enum { ModeIsAffineCompact = Mode == int(AffineCompact),
325  OtherModeIsAffineCompact = OtherMode == int(AffineCompact)
326  };
327 
328  if(ModeIsAffineCompact == OtherModeIsAffineCompact)
329  {
330  // We need the block expression because the code is compiled for all
331  // combinations of transformations and will trigger a compile time error
332  // if one tries to assign the matrices directly
333  m_matrix.template block<Dim,Dim+1>(0,0) = other.matrix().template block<Dim,Dim+1>(0,0);
334  makeAffine();
335  }
336  else if(OtherModeIsAffineCompact)
337  {
338  typedef typename Transform<Scalar,Dim,OtherMode,OtherOptions>::MatrixType OtherMatrixType;
339  internal::transform_construct_from_matrix<OtherMatrixType,Mode,Options,Dim,HDim>::run(this, other.matrix());
340  }
341  else
342  {
343  // here we know that Mode == AffineCompact and OtherMode != AffineCompact.
344  // if OtherMode were Projective, the static assert above would already have caught it.
345  // So the only possibility is that OtherMode == Affine
346  linear() = other.linear();
347  translation() = other.translation();
348  }
349  }
350 
351  template<typename OtherDerived>
353  {
354  check_template_params();
355  other.evalTo(*this);
356  }
357 
358  template<typename OtherDerived>
359  Transform& operator=(const ReturnByValue<OtherDerived>& other)
360  {
361  other.evalTo(*this);
362  return *this;
363  }
364 
365  #ifdef EIGEN_QT_SUPPORT
366  inline Transform(const QMatrix& other);
367  inline Transform& operator=(const QMatrix& other);
368  inline QMatrix toQMatrix(void) const;
369  inline Transform(const QTransform& other);
370  inline Transform& operator=(const QTransform& other);
371  inline QTransform toQTransform(void) const;
372  #endif
373 
376  inline Scalar operator() (Index row, Index col) const { return m_matrix(row,col); }
379  inline Scalar& operator() (Index row, Index col) { return m_matrix(row,col); }
380 
382  inline const MatrixType& matrix() const { return m_matrix; }
384  inline MatrixType& matrix() { return m_matrix; }
385 
387  inline ConstLinearPart linear() const { return ConstLinearPart(m_matrix,0,0); }
389  inline LinearPart linear() { return LinearPart(m_matrix,0,0); }
390 
392  inline ConstAffinePart affine() const { return take_affine_part::run(m_matrix); }
394  inline AffinePart affine() { return take_affine_part::run(m_matrix); }
395 
397  inline ConstTranslationPart translation() const { return ConstTranslationPart(m_matrix,0,Dim); }
399  inline TranslationPart translation() { return TranslationPart(m_matrix,0,Dim); }
400 
412  // note: this function is defined here because some compilers cannot find the respective declaration
413  template<typename OtherDerived>
414  EIGEN_STRONG_INLINE const typename internal::transform_right_product_impl<Transform, OtherDerived>::ResultType
416  { return internal::transform_right_product_impl<Transform, OtherDerived>::run(*this,other.derived()); }
417 
425  template<typename OtherDerived> friend
426  inline const typename internal::transform_left_product_impl<OtherDerived,Mode,Options,_Dim,_Dim+1>::ResultType
428  { return internal::transform_left_product_impl<OtherDerived,Mode,Options,Dim,HDim>::run(a.derived(),b); }
429 
436  template<typename DiagonalDerived>
437  inline const TransformTimeDiagonalReturnType
439  {
440  TransformTimeDiagonalReturnType res(*this);
441  res.linear() *= b;
442  return res;
443  }
444 
451  template<typename DiagonalDerived>
452  friend inline TransformTimeDiagonalReturnType
454  {
455  TransformTimeDiagonalReturnType res;
456  res.linear().noalias() = a*b.linear();
457  res.translation().noalias() = a*b.translation();
458  if (Mode!=int(AffineCompact))
459  res.matrix().row(Dim) = b.matrix().row(Dim);
460  return res;
461  }
462 
463  template<typename OtherDerived>
464  inline Transform& operator*=(const EigenBase<OtherDerived>& other) { return *this = *this * other; }
465 
467  inline const Transform operator * (const Transform& other) const
468  {
469  return internal::transform_transform_product_impl<Transform,Transform>::run(*this,other);
470  }
471 
473  template<int OtherMode,int OtherOptions>
474  inline const typename internal::transform_transform_product_impl<
477  {
478  return internal::transform_transform_product_impl<Transform,Transform<Scalar,Dim,OtherMode,OtherOptions> >::run(*this,other);
479  }
480 
482  void setIdentity() { m_matrix.setIdentity(); }
483 
488  static const Transform Identity()
489  {
490  return Transform(MatrixType::Identity());
491  }
492 
493  template<typename OtherDerived>
494  inline Transform& scale(const MatrixBase<OtherDerived> &other);
495 
496  template<typename OtherDerived>
497  inline Transform& prescale(const MatrixBase<OtherDerived> &other);
498 
499  inline Transform& scale(Scalar s);
500  inline Transform& prescale(Scalar s);
501 
502  template<typename OtherDerived>
503  inline Transform& translate(const MatrixBase<OtherDerived> &other);
504 
505  template<typename OtherDerived>
506  inline Transform& pretranslate(const MatrixBase<OtherDerived> &other);
507 
508  template<typename RotationType>
509  inline Transform& rotate(const RotationType& rotation);
510 
511  template<typename RotationType>
512  inline Transform& prerotate(const RotationType& rotation);
513 
514  Transform& shear(Scalar sx, Scalar sy);
515  Transform& preshear(Scalar sx, Scalar sy);
516 
517  inline Transform& operator=(const TranslationType& t);
518  inline Transform& operator*=(const TranslationType& t) { return translate(t.vector()); }
519  inline Transform operator*(const TranslationType& t) const;
520 
521  inline Transform& operator=(const UniformScaling<Scalar>& t);
522  inline Transform& operator*=(const UniformScaling<Scalar>& s) { return scale(s.factor()); }
523  inline Transform operator*(const UniformScaling<Scalar>& s) const;
524 
525  inline Transform& operator*=(const DiagonalMatrix<Scalar,Dim>& s) { linear() *= s; return *this; }
526 
527  template<typename Derived>
528  inline Transform& operator=(const RotationBase<Derived,Dim>& r);
529  template<typename Derived>
530  inline Transform& operator*=(const RotationBase<Derived,Dim>& r) { return rotate(r.toRotationMatrix()); }
531  template<typename Derived>
532  inline Transform operator*(const RotationBase<Derived,Dim>& r) const;
533 
534  const LinearMatrixType rotation() const;
535  template<typename RotationMatrixType, typename ScalingMatrixType>
536  void computeRotationScaling(RotationMatrixType *rotation, ScalingMatrixType *scaling) const;
537  template<typename ScalingMatrixType, typename RotationMatrixType>
538  void computeScalingRotation(ScalingMatrixType *scaling, RotationMatrixType *rotation) const;
539 
540  template<typename PositionDerived, typename OrientationType, typename ScaleDerived>
541  Transform& fromPositionOrientationScale(const MatrixBase<PositionDerived> &position,
542  const OrientationType& orientation, const MatrixBase<ScaleDerived> &scale);
543 
544  inline Transform inverse(TransformTraits traits = (TransformTraits)Mode) const;
545 
547  const Scalar* data() const { return m_matrix.data(); }
549  Scalar* data() { return m_matrix.data(); }
550 
556  template<typename NewScalarType>
557  inline typename internal::cast_return_type<Transform,Transform<NewScalarType,Dim,Mode,Options> >::type cast() const
558  { return typename internal::cast_return_type<Transform,Transform<NewScalarType,Dim,Mode,Options> >::type(*this); }
559 
561  template<typename OtherScalarType>
562  inline explicit Transform(const Transform<OtherScalarType,Dim,Mode,Options>& other)
563  {
564  check_template_params();
565  m_matrix = other.matrix().template cast<Scalar>();
566  }
567 
572  bool isApprox(const Transform& other, typename NumTraits<Scalar>::Real prec = NumTraits<Scalar>::dummy_precision()) const
573  { return m_matrix.isApprox(other.m_matrix, prec); }
574 
577  void makeAffine()
578  {
579  if(int(Mode)!=int(AffineCompact))
580  {
581  matrix().template block<1,Dim>(Dim,0).setZero();
582  matrix().coeffRef(Dim,Dim) = 1;
583  }
584  }
585 
591  { return m_matrix.template block<int(Mode)==int(Projective)?HDim:Dim,Dim>(0,0); }
597  { return m_matrix.template block<int(Mode)==int(Projective)?HDim:Dim,Dim>(0,0); }
598 
604  { return m_matrix.template block<int(Mode)==int(Projective)?HDim:Dim,1>(0,Dim); }
610  { return m_matrix.template block<int(Mode)==int(Projective)?HDim:Dim,1>(0,Dim); }
611 
612 
613  #ifdef EIGEN_TRANSFORM_PLUGIN
614  #include EIGEN_TRANSFORM_PLUGIN
615  #endif
616 
617 protected:
618  #ifndef EIGEN_PARSED_BY_DOXYGEN
619  static EIGEN_STRONG_INLINE void check_template_params()
620  {
621  EIGEN_STATIC_ASSERT((Options & (DontAlign|RowMajor)) == Options, INVALID_MATRIX_TEMPLATE_PARAMETERS)
622  }
623  #endif
624 
625 };
626 
635 
644 
653 
662 
663 /**************************
664 *** Optional QT support ***
665 **************************/
666 
667 #ifdef EIGEN_QT_SUPPORT
668 
672 template<typename Scalar, int Dim, int Mode,int Options>
674 {
675  check_template_params();
676  *this = other;
677 }
678 
683 template<typename Scalar, int Dim, int Mode,int Options>
685 {
686  EIGEN_STATIC_ASSERT(Dim==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
687  m_matrix << other.m11(), other.m21(), other.dx(),
688  other.m12(), other.m22(), other.dy(),
689  0, 0, 1;
690  return *this;
691 }
692 
699 template<typename Scalar, int Dim, int Mode, int Options>
701 {
702  check_template_params();
703  EIGEN_STATIC_ASSERT(Dim==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
704  return QMatrix(m_matrix.coeff(0,0), m_matrix.coeff(1,0),
705  m_matrix.coeff(0,1), m_matrix.coeff(1,1),
706  m_matrix.coeff(0,2), m_matrix.coeff(1,2));
707 }
708 
713 template<typename Scalar, int Dim, int Mode,int Options>
715 {
716  check_template_params();
717  *this = other;
718 }
719 
724 template<typename Scalar, int Dim, int Mode, int Options>
726 {
727  check_template_params();
728  EIGEN_STATIC_ASSERT(Dim==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
729  if (Mode == int(AffineCompact))
730  m_matrix << other.m11(), other.m21(), other.dx(),
731  other.m12(), other.m22(), other.dy();
732  else
733  m_matrix << other.m11(), other.m21(), other.dx(),
734  other.m12(), other.m22(), other.dy(),
735  other.m13(), other.m23(), other.m33();
736  return *this;
737 }
738 
743 template<typename Scalar, int Dim, int Mode, int Options>
745 {
746  EIGEN_STATIC_ASSERT(Dim==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
747  if (Mode == int(AffineCompact))
748  return QTransform(m_matrix.coeff(0,0), m_matrix.coeff(1,0),
749  m_matrix.coeff(0,1), m_matrix.coeff(1,1),
750  m_matrix.coeff(0,2), m_matrix.coeff(1,2));
751  else
752  return QTransform(m_matrix.coeff(0,0), m_matrix.coeff(1,0), m_matrix.coeff(2,0),
753  m_matrix.coeff(0,1), m_matrix.coeff(1,1), m_matrix.coeff(2,1),
754  m_matrix.coeff(0,2), m_matrix.coeff(1,2), m_matrix.coeff(2,2));
755 }
756 #endif
757 
758 /*********************
759 *** Procedural API ***
760 *********************/
761 
766 template<typename Scalar, int Dim, int Mode, int Options>
767 template<typename OtherDerived>
770 {
771  EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,int(Dim))
772  EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
773  linearExt().noalias() = (linearExt() * other.asDiagonal());
774  return *this;
775 }
776 
781 template<typename Scalar, int Dim, int Mode, int Options>
783 {
784  EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
785  linearExt() *= s;
786  return *this;
787 }
788 
793 template<typename Scalar, int Dim, int Mode, int Options>
794 template<typename OtherDerived>
797 {
798  EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,int(Dim))
799  EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
800  m_matrix.template block<Dim,HDim>(0,0).noalias() = (other.asDiagonal() * m_matrix.template block<Dim,HDim>(0,0));
801  return *this;
802 }
803 
808 template<typename Scalar, int Dim, int Mode, int Options>
810 {
811  EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
812  m_matrix.template topRows<Dim>() *= s;
813  return *this;
814 }
815 
820 template<typename Scalar, int Dim, int Mode, int Options>
821 template<typename OtherDerived>
824 {
825  EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,int(Dim))
826  translationExt() += linearExt() * other;
827  return *this;
828 }
829 
834 template<typename Scalar, int Dim, int Mode, int Options>
835 template<typename OtherDerived>
838 {
839  EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(OtherDerived,int(Dim))
840  if(int(Mode)==int(Projective))
841  affine() += other * m_matrix.row(Dim);
842  else
843  translation() += other;
844  return *this;
845 }
846 
864 template<typename Scalar, int Dim, int Mode, int Options>
865 template<typename RotationType>
867 Transform<Scalar,Dim,Mode,Options>::rotate(const RotationType& rotation)
868 {
869  linearExt() *= internal::toRotationMatrix<Scalar,Dim>(rotation);
870  return *this;
871 }
872 
880 template<typename Scalar, int Dim, int Mode, int Options>
881 template<typename RotationType>
884 {
885  m_matrix.template block<Dim,HDim>(0,0) = internal::toRotationMatrix<Scalar,Dim>(rotation)
886  * m_matrix.template block<Dim,HDim>(0,0);
887  return *this;
888 }
889 
895 template<typename Scalar, int Dim, int Mode, int Options>
898 {
899  EIGEN_STATIC_ASSERT(int(Dim)==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
900  EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
901  VectorType tmp = linear().col(0)*sy + linear().col(1);
902  linear() << linear().col(0) + linear().col(1)*sx, tmp;
903  return *this;
904 }
905 
911 template<typename Scalar, int Dim, int Mode, int Options>
914 {
915  EIGEN_STATIC_ASSERT(int(Dim)==2, YOU_MADE_A_PROGRAMMING_MISTAKE)
916  EIGEN_STATIC_ASSERT(Mode!=int(Isometry), THIS_METHOD_IS_ONLY_FOR_SPECIFIC_TRANSFORMATIONS)
917  m_matrix.template block<Dim,HDim>(0,0) = LinearMatrixType(1, sx, sy, 1) * m_matrix.template block<Dim,HDim>(0,0);
918  return *this;
919 }
920 
921 /******************************************************
922 *** Scaling, Translation and Rotation compatibility ***
923 ******************************************************/
924 
925 template<typename Scalar, int Dim, int Mode, int Options>
927 {
928  linear().setIdentity();
929  translation() = t.vector();
930  makeAffine();
931  return *this;
932 }
933 
934 template<typename Scalar, int Dim, int Mode, int Options>
936 {
937  Transform res = *this;
938  res.translate(t.vector());
939  return res;
940 }
941 
942 template<typename Scalar, int Dim, int Mode, int Options>
944 {
945  m_matrix.setZero();
946  linear().diagonal().fill(s.factor());
947  makeAffine();
948  return *this;
949 }
950 
951 template<typename Scalar, int Dim, int Mode, int Options>
953 {
954  Transform res = *this;
955  res.scale(s.factor());
956  return res;
957 }
958 
959 template<typename Scalar, int Dim, int Mode, int Options>
960 template<typename Derived>
962 {
963  linear() = internal::toRotationMatrix<Scalar,Dim>(r);
964  translation().setZero();
965  makeAffine();
966  return *this;
967 }
968 
969 template<typename Scalar, int Dim, int Mode, int Options>
970 template<typename Derived>
972 {
973  Transform res = *this;
974  res.rotate(r.derived());
975  return res;
976 }
977 
978 /************************
979 *** Special functions ***
980 ************************/
981 
989 template<typename Scalar, int Dim, int Mode, int Options>
992 {
993  LinearMatrixType result;
994  computeRotationScaling(&result, (LinearMatrixType*)0);
995  return result;
996 }
997 
998 
1010 template<typename Scalar, int Dim, int Mode, int Options>
1011 template<typename RotationMatrixType, typename ScalingMatrixType>
1012 void Transform<Scalar,Dim,Mode,Options>::computeRotationScaling(RotationMatrixType *rotation, ScalingMatrixType *scaling) const
1013 {
1015 
1016  Scalar x = (svd.matrixU() * svd.matrixV().adjoint()).determinant(); // so x has absolute value 1
1017  VectorType sv(svd.singularValues());
1018  sv.coeffRef(0) *= x;
1019  if(scaling) scaling->lazyAssign(svd.matrixV() * sv.asDiagonal() * svd.matrixV().adjoint());
1020  if(rotation)
1021  {
1022  LinearMatrixType m(svd.matrixU());
1023  m.col(0) /= x;
1024  rotation->lazyAssign(m * svd.matrixV().adjoint());
1025  }
1026 }
1027 
1039 template<typename Scalar, int Dim, int Mode, int Options>
1040 template<typename ScalingMatrixType, typename RotationMatrixType>
1041 void Transform<Scalar,Dim,Mode,Options>::computeScalingRotation(ScalingMatrixType *scaling, RotationMatrixType *rotation) const
1042 {
1044 
1045  Scalar x = (svd.matrixU() * svd.matrixV().adjoint()).determinant(); // so x has absolute value 1
1046  VectorType sv(svd.singularValues());
1047  sv.coeffRef(0) *= x;
1048  if(scaling) scaling->lazyAssign(svd.matrixU() * sv.asDiagonal() * svd.matrixU().adjoint());
1049  if(rotation)
1050  {
1051  LinearMatrixType m(svd.matrixU());
1052  m.col(0) /= x;
1053  rotation->lazyAssign(m * svd.matrixV().adjoint());
1054  }
1055 }
1056 
1060 template<typename Scalar, int Dim, int Mode, int Options>
1061 template<typename PositionDerived, typename OrientationType, typename ScaleDerived>
1064  const OrientationType& orientation, const MatrixBase<ScaleDerived> &scale)
1065 {
1066  linear() = internal::toRotationMatrix<Scalar,Dim>(orientation);
1067  linear() *= scale.asDiagonal();
1068  translation() = position;
1069  makeAffine();
1070  return *this;
1071 }
1072 
1073 namespace internal {
1074 
1075 // selector needed to avoid taking the inverse of a 3x4 matrix
1076 template<typename TransformType, int Mode=TransformType::Mode>
1077 struct projective_transform_inverse
1078 {
1079  static inline void run(const TransformType&, TransformType&)
1080  {}
1081 };
1082 
1083 template<typename TransformType>
1084 struct projective_transform_inverse<TransformType, Projective>
1085 {
1086  static inline void run(const TransformType& m, TransformType& res)
1087  {
1088  res.matrix() = m.matrix().inverse();
1089  }
1090 };
1091 
1092 } // end namespace internal
1093 
1094 
1115 template<typename Scalar, int Dim, int Mode, int Options>
1116 Transform<Scalar,Dim,Mode,Options>
1118 {
1119  Transform res;
1120  if (hint == Projective)
1121  {
1122  internal::projective_transform_inverse<Transform>::run(*this, res);
1123  }
1124  else
1125  {
1126  if (hint == Isometry)
1127  {
1128  res.matrix().template topLeftCorner<Dim,Dim>() = linear().transpose();
1129  }
1130  else if(hint&Affine)
1131  {
1132  res.matrix().template topLeftCorner<Dim,Dim>() = linear().inverse();
1133  }
1134  else
1135  {
1136  eigen_assert(false && "Invalid transform traits in Transform::Inverse");
1137  }
1138  // translation and remaining parts
1139  res.matrix().template topRightCorner<Dim,1>()
1140  = - res.matrix().template topLeftCorner<Dim,Dim>() * translation();
1141  res.makeAffine(); // we do need this, because in the beginning res is uninitialized
1142  }
1143  return res;
1144 }
1145 
1146 namespace internal {
1147 
1148 /*****************************************************
1149 *** Specializations of take affine part ***
1150 *****************************************************/
1151 
1152 template<typename TransformType> struct transform_take_affine_part {
1153  typedef typename TransformType::MatrixType MatrixType;
1154  typedef typename TransformType::AffinePart AffinePart;
1155  typedef typename TransformType::ConstAffinePart ConstAffinePart;
1156  static inline AffinePart run(MatrixType& m)
1157  { return m.template block<TransformType::Dim,TransformType::HDim>(0,0); }
1158  static inline ConstAffinePart run(const MatrixType& m)
1159  { return m.template block<TransformType::Dim,TransformType::HDim>(0,0); }
1160 };
1161 
1162 template<typename Scalar, int Dim, int Options>
1163 struct transform_take_affine_part<Transform<Scalar,Dim,AffineCompact, Options> > {
1164  typedef typename Transform<Scalar,Dim,AffineCompact,Options>::MatrixType MatrixType;
1165  static inline MatrixType& run(MatrixType& m) { return m; }
1166  static inline const MatrixType& run(const MatrixType& m) { return m; }
1167 };
1168 
1169 /*****************************************************
1170 *** Specializations of construct from matrix ***
1171 *****************************************************/
1172 
1173 template<typename Other, int Mode, int Options, int Dim, int HDim>
1174 struct transform_construct_from_matrix<Other, Mode,Options,Dim,HDim, Dim,Dim>
1175 {
1176  static inline void run(Transform<typename Other::Scalar,Dim,Mode,Options> *transform, const Other& other)
1177  {
1178  transform->linear() = other;
1179  transform->translation().setZero();
1180  transform->makeAffine();
1181  }
1182 };
1183 
1184 template<typename Other, int Mode, int Options, int Dim, int HDim>
1185 struct transform_construct_from_matrix<Other, Mode,Options,Dim,HDim, Dim,HDim>
1186 {
1187  static inline void run(Transform<typename Other::Scalar,Dim,Mode,Options> *transform, const Other& other)
1188  {
1189  transform->affine() = other;
1190  transform->makeAffine();
1191  }
1192 };
1193 
1194 template<typename Other, int Mode, int Options, int Dim, int HDim>
1195 struct transform_construct_from_matrix<Other, Mode,Options,Dim,HDim, HDim,HDim>
1196 {
1197  static inline void run(Transform<typename Other::Scalar,Dim,Mode,Options> *transform, const Other& other)
1198  { transform->matrix() = other; }
1199 };
1200 
1201 template<typename Other, int Options, int Dim, int HDim>
1202 struct transform_construct_from_matrix<Other, AffineCompact,Options,Dim,HDim, HDim,HDim>
1203 {
1204  static inline void run(Transform<typename Other::Scalar,Dim,AffineCompact,Options> *transform, const Other& other)
1205  { transform->matrix() = other.template block<Dim,HDim>(0,0); }
1206 };
1207 
1208 /**********************************************************
1209 *** Specializations of operator* with rhs EigenBase ***
1210 **********************************************************/
1211 
1212 template<int LhsMode,int RhsMode>
1213 struct transform_product_result
1214 {
1215  enum
1216  {
1217  Mode =
1218  (LhsMode == (int)Projective || RhsMode == (int)Projective ) ? Projective :
1219  (LhsMode == (int)Affine || RhsMode == (int)Affine ) ? Affine :
1220  (LhsMode == (int)AffineCompact || RhsMode == (int)AffineCompact ) ? AffineCompact :
1221  (LhsMode == (int)Isometry || RhsMode == (int)Isometry ) ? Isometry : Projective
1222  };
1223 };
1224 
1225 template< typename TransformType, typename MatrixType >
1226 struct transform_right_product_impl< TransformType, MatrixType, 0 >
1227 {
1228  typedef typename MatrixType::PlainObject ResultType;
1229 
1230  static EIGEN_STRONG_INLINE ResultType run(const TransformType& T, const MatrixType& other)
1231  {
1232  return T.matrix() * other;
1233  }
1234 };
1235 
1236 template< typename TransformType, typename MatrixType >
1237 struct transform_right_product_impl< TransformType, MatrixType, 1 >
1238 {
1239  enum {
1240  Dim = TransformType::Dim,
1241  HDim = TransformType::HDim,
1242  OtherRows = MatrixType::RowsAtCompileTime,
1243  OtherCols = MatrixType::ColsAtCompileTime
1244  };
1245 
1246  typedef typename MatrixType::PlainObject ResultType;
1247 
1248  static EIGEN_STRONG_INLINE ResultType run(const TransformType& T, const MatrixType& other)
1249  {
1250  EIGEN_STATIC_ASSERT(OtherRows==HDim, YOU_MIXED_MATRICES_OF_DIFFERENT_SIZES);
1251 
1252  typedef Block<ResultType, Dim, OtherCols, int(MatrixType::RowsAtCompileTime)==Dim> TopLeftLhs;
1253 
1254  ResultType res(other.rows(),other.cols());
1255  TopLeftLhs(res, 0, 0, Dim, other.cols()).noalias() = T.affine() * other;
1256  res.row(OtherRows-1) = other.row(OtherRows-1);
1257 
1258  return res;
1259  }
1260 };
1261 
1262 template< typename TransformType, typename MatrixType >
1263 struct transform_right_product_impl< TransformType, MatrixType, 2 >
1264 {
1265  enum {
1266  Dim = TransformType::Dim,
1267  HDim = TransformType::HDim,
1268  OtherRows = MatrixType::RowsAtCompileTime,
1269  OtherCols = MatrixType::ColsAtCompileTime
1270  };
1271 
1272  typedef typename MatrixType::PlainObject ResultType;
1273 
1274  static EIGEN_STRONG_INLINE ResultType run(const TransformType& T, const MatrixType& other)
1275  {
1276  EIGEN_STATIC_ASSERT(OtherRows==Dim, YOU_MIXED_MATRICES_OF_DIFFERENT_SIZES);
1277 
1278  typedef Block<ResultType, Dim, OtherCols, true> TopLeftLhs;
1279  ResultType res(Replicate<typename TransformType::ConstTranslationPart, 1, OtherCols>(T.translation(),1,other.cols()));
1280  TopLeftLhs(res, 0, 0, Dim, other.cols()).noalias() += T.linear() * other;
1281 
1282  return res;
1283  }
1284 };
1285 
1286 /**********************************************************
1287 *** Specializations of operator* with lhs EigenBase ***
1288 **********************************************************/
1289 
1290 // generic HDim x HDim matrix * T => Projective
1291 template<typename Other,int Mode, int Options, int Dim, int HDim>
1292 struct transform_left_product_impl<Other,Mode,Options,Dim,HDim, HDim,HDim>
1293 {
1294  typedef Transform<typename Other::Scalar,Dim,Mode,Options> TransformType;
1295  typedef typename TransformType::MatrixType MatrixType;
1296  typedef Transform<typename Other::Scalar,Dim,Projective,Options> ResultType;
1297  static ResultType run(const Other& other,const TransformType& tr)
1298  { return ResultType(other * tr.matrix()); }
1299 };
1300 
1301 // generic HDim x HDim matrix * AffineCompact => Projective
1302 template<typename Other, int Options, int Dim, int HDim>
1303 struct transform_left_product_impl<Other,AffineCompact,Options,Dim,HDim, HDim,HDim>
1304 {
1305  typedef Transform<typename Other::Scalar,Dim,AffineCompact,Options> TransformType;
1306  typedef typename TransformType::MatrixType MatrixType;
1307  typedef Transform<typename Other::Scalar,Dim,Projective,Options> ResultType;
1308  static ResultType run(const Other& other,const TransformType& tr)
1309  {
1310  ResultType res;
1311  res.matrix().noalias() = other.template block<HDim,Dim>(0,0) * tr.matrix();
1312  res.matrix().col(Dim) += other.col(Dim);
1313  return res;
1314  }
1315 };
1316 
1317 // affine matrix * T
1318 template<typename Other,int Mode, int Options, int Dim, int HDim>
1319 struct transform_left_product_impl<Other,Mode,Options,Dim,HDim, Dim,HDim>
1320 {
1321  typedef Transform<typename Other::Scalar,Dim,Mode,Options> TransformType;
1322  typedef typename TransformType::MatrixType MatrixType;
1323  typedef TransformType ResultType;
1324  static ResultType run(const Other& other,const TransformType& tr)
1325  {
1326  ResultType res;
1327  res.affine().noalias() = other * tr.matrix();
1328  res.matrix().row(Dim) = tr.matrix().row(Dim);
1329  return res;
1330  }
1331 };
1332 
1333 // affine matrix * AffineCompact
1334 template<typename Other, int Options, int Dim, int HDim>
1335 struct transform_left_product_impl<Other,AffineCompact,Options,Dim,HDim, Dim,HDim>
1336 {
1337  typedef Transform<typename Other::Scalar,Dim,AffineCompact,Options> TransformType;
1338  typedef typename TransformType::MatrixType MatrixType;
1339  typedef TransformType ResultType;
1340  static ResultType run(const Other& other,const TransformType& tr)
1341  {
1342  ResultType res;
1343  res.matrix().noalias() = other.template block<Dim,Dim>(0,0) * tr.matrix();
1344  res.translation() += other.col(Dim);
1345  return res;
1346  }
1347 };
1348 
1349 // linear matrix * T
1350 template<typename Other,int Mode, int Options, int Dim, int HDim>
1351 struct transform_left_product_impl<Other,Mode,Options,Dim,HDim, Dim,Dim>
1352 {
1353  typedef Transform<typename Other::Scalar,Dim,Mode,Options> TransformType;
1354  typedef typename TransformType::MatrixType MatrixType;
1355  typedef TransformType ResultType;
1356  static ResultType run(const Other& other, const TransformType& tr)
1357  {
1358  TransformType res;
1359  if(Mode!=int(AffineCompact))
1360  res.matrix().row(Dim) = tr.matrix().row(Dim);
1361  res.matrix().template topRows<Dim>().noalias()
1362  = other * tr.matrix().template topRows<Dim>();
1363  return res;
1364  }
1365 };
1366 
1367 /**********************************************************
1368 *** Specializations of operator* with another Transform ***
1369 **********************************************************/
1370 
1371 template<typename Scalar, int Dim, int LhsMode, int LhsOptions, int RhsMode, int RhsOptions>
1372 struct transform_transform_product_impl<Transform<Scalar,Dim,LhsMode,LhsOptions>,Transform<Scalar,Dim,RhsMode,RhsOptions>,false >
1373 {
1374  enum { ResultMode = transform_product_result<LhsMode,RhsMode>::Mode };
1375  typedef Transform<Scalar,Dim,LhsMode,LhsOptions> Lhs;
1376  typedef Transform<Scalar,Dim,RhsMode,RhsOptions> Rhs;
1377  typedef Transform<Scalar,Dim,ResultMode,LhsOptions> ResultType;
1378  static ResultType run(const Lhs& lhs, const Rhs& rhs)
1379  {
1380  ResultType res;
1381  res.linear() = lhs.linear() * rhs.linear();
1382  res.translation() = lhs.linear() * rhs.translation() + lhs.translation();
1383  res.makeAffine();
1384  return res;
1385  }
1386 };
1387 
1388 template<typename Scalar, int Dim, int LhsMode, int LhsOptions, int RhsMode, int RhsOptions>
1389 struct transform_transform_product_impl<Transform<Scalar,Dim,LhsMode,LhsOptions>,Transform<Scalar,Dim,RhsMode,RhsOptions>,true >
1390 {
1391  typedef Transform<Scalar,Dim,LhsMode,LhsOptions> Lhs;
1392  typedef Transform<Scalar,Dim,RhsMode,RhsOptions> Rhs;
1393  typedef Transform<Scalar,Dim,Projective> ResultType;
1394  static ResultType run(const Lhs& lhs, const Rhs& rhs)
1395  {
1396  return ResultType( lhs.matrix() * rhs.matrix() );
1397  }
1398 };
1399 
1400 template<typename Scalar, int Dim, int LhsOptions, int RhsOptions>
1401 struct transform_transform_product_impl<Transform<Scalar,Dim,AffineCompact,LhsOptions>,Transform<Scalar,Dim,Projective,RhsOptions>,true >
1402 {
1403  typedef Transform<Scalar,Dim,AffineCompact,LhsOptions> Lhs;
1404  typedef Transform<Scalar,Dim,Projective,RhsOptions> Rhs;
1405  typedef Transform<Scalar,Dim,Projective> ResultType;
1406  static ResultType run(const Lhs& lhs, const Rhs& rhs)
1407  {
1408  ResultType res;
1409  res.matrix().template topRows<Dim>() = lhs.matrix() * rhs.matrix();
1410  res.matrix().row(Dim) = rhs.matrix().row(Dim);
1411  return res;
1412  }
1413 };
1414 
1415 template<typename Scalar, int Dim, int LhsOptions, int RhsOptions>
1416 struct transform_transform_product_impl<Transform<Scalar,Dim,Projective,LhsOptions>,Transform<Scalar,Dim,AffineCompact,RhsOptions>,true >
1417 {
1418  typedef Transform<Scalar,Dim,Projective,LhsOptions> Lhs;
1419  typedef Transform<Scalar,Dim,AffineCompact,RhsOptions> Rhs;
1420  typedef Transform<Scalar,Dim,Projective> ResultType;
1421  static ResultType run(const Lhs& lhs, const Rhs& rhs)
1422  {
1423  ResultType res(lhs.matrix().template leftCols<Dim>() * rhs.matrix());
1424  res.matrix().col(Dim) += lhs.matrix().col(Dim);
1425  return res;
1426  }
1427 };
1428 
1429 } // end namespace internal
1430 
1431 } // end namespace Eigen
1432 
1433 #endif // EIGEN_TRANSFORM_H