Logo ROOT   6.10/09
Reference Guide
Transform3D.h
Go to the documentation of this file.
1 // @(#)root/mathcore:$Id$
2 // Authors: W. Brown, M. Fischler, L. Moneta 2005
3 
4 /**********************************************************************
5  * *
6  * Copyright (c) 2005 , LCG ROOT MathLib Team *
7  * *
8  * *
9  **********************************************************************/
10 
11 // Header file for class Transform3D
12 //
13 // Created by: Lorenzo Moneta October 21 2005
14 //
15 //
16 #ifndef ROOT_Math_GenVector_Transform3D
17 #define ROOT_Math_GenVector_Transform3D 1
18 
19 
20 
22 
24 
26 
28 
29 
37 
38 #include <iostream>
39 #include <type_traits>
40 #include <cmath>
41 
42 //#include "Math/Vector3Dfwd.h"
43 
44 
45 
46 namespace ROOT {
47 
48 namespace Math {
49 
50 namespace Impl {
51 
52 //_________________________________________________________________________________________
53 /**
54  Basic 3D Transformation class describing a rotation and then a translation
55  The internal data are a 3D rotation data (represented as a 3x3 matrix) and a 3D vector data.
56  They are represented and held in this class like a 3x4 matrix (a simple array of 12 numbers).
57 
58  The class can be constructed from any 3D rotation object
59  (ROOT::Math::Rotation3D, ROOT::Math::AxisAngle, ROOT::Math::Quaternion, etc...) and/or
60  a 3D Vector (ROOT::Math::DislacementVector3D or via ROOT::Math::Translation ) representing a Translation.
61  The Transformation is defined by applying first the rotation and then the translation.
62  A transformation defined by applying first a translation and then a rotation is equivalent to the
63  transformation obtained applying first the rotation and then a translation equivalent to the rotated vector.
64  The operator * can be used to obtain directly such transformations, in addition to combine various
65  transformations.
66  Keep in mind that the operator * (like in the case of rotations ) is not commutative.
67  The operator * is used (in addition to operator() ) to apply a transformations on the vector
68  (DisplacementVector3D and LorentzVector classes) and point (PositionVector3D) classes.
69  In the case of Vector objects the transformation only rotates them and does not translate them.
70  Only Point objects are able to be both rotated and translated.
71 
72 
73  @ingroup GenVector
74 
75 */
76 
77 template <typename T = double>
78 class Transform3D {
79 
80 public:
81  typedef T Scalar;
82 
85 
87  kXX = 0, kXY = 1, kXZ = 2, kDX = 3,
88  kYX = 4, kYY = 5, kYZ = 6, kDY = 7,
89  kZX = 8, kZY = 9, kZZ =10, kDZ = 11
90  };
91 
92 
93 
94  /**
95  Default constructor (identy rotation) + zero translation
96  */
98  {
99  SetIdentity();
100  }
101 
102  /**
103  Construct given a pair of pointers or iterators defining the
104  beginning and end of an array of 12 Scalars
105  */
106  template<class IT>
107  Transform3D(IT begin, IT end)
108  {
109  SetComponents(begin,end);
110  }
111 
112  /**
113  Construct from a rotation and then a translation described by a Vector
114  */
115  Transform3D( const Rotation3D & r, const Vector & v)
116  {
117  AssignFrom( r, v );
118  }
119  /**
120  Construct from a rotation and then a translation described by a Translation3D class
121  */
122  Transform3D(const Rotation3D &r, const Translation3D<T> &t) { AssignFrom(r, t.Vect()); }
123 
124  /**
125  Construct from a rotation (any rotation object) and then a translation
126  (represented by any DisplacementVector)
127  The requirements on the rotation and vector objects are that they can be transformed in a
128  Rotation3D class and in a Cartesian3D Vector
129  */
130  template <class ARotation, class CoordSystem, class Tag>
132  {
133  AssignFrom( Rotation3D(r), Vector (v.X(),v.Y(),v.Z()) );
134  }
135 
136  /**
137  Construct from a rotation (any rotation object) and then a translation
138  represented by a Translation3D class
139  The requirements on the rotation is that it can be transformed in a
140  Rotation3D class
141  */
142  template <class ARotation>
143  Transform3D(const ARotation &r, const Translation3D<T> &t)
144  {
145  AssignFrom( Rotation3D(r), t.Vect() );
146  }
147 
148 
149 #ifdef OLD_VERSION
150  /**
151  Construct from a translation and then a rotation (inverse assignment)
152  */
153  Transform3D( const Vector & v, const Rotation3D & r)
154  {
155  // is equivalent from having first the rotation and then the translation vector rotated
156  AssignFrom( r, r(v) );
157  }
158 #endif
159 
160  /**
161  Construct from a 3D Rotation only with zero translation
162  */
163  explicit Transform3D( const Rotation3D & r) {
164  AssignFrom(r);
165  }
166 
167  // convenience methods for constructing a Transform3D from all the 3D rotations classes
168  // (cannot use templates for conflict with LA)
169 
170  explicit Transform3D( const AxisAngle & r) {
172  }
173  explicit Transform3D( const EulerAngles & r) {
175  }
176  explicit Transform3D( const Quaternion & r) {
178  }
179  explicit Transform3D( const RotationZYX & r) {
181  }
182 
183  // Constructors from axial rotations
184  // TO DO: implement direct methods for axial rotations without going through Rotation3D
185  explicit Transform3D( const RotationX & r) {
187  }
188  explicit Transform3D( const RotationY & r) {
190  }
191  explicit Transform3D( const RotationZ & r) {
193  }
194 
195  /**
196  Construct from a translation only, represented by any DisplacementVector3D
197  and with an identity rotation
198  */
199  template<class CoordSystem, class Tag>
201  AssignFrom(Vector(v.X(),v.Y(),v.Z()));
202  }
203  /**
204  Construct from a translation only, represented by a Cartesian 3D Vector,
205  and with an identity rotation
206  */
207  explicit Transform3D( const Vector & v) {
208  AssignFrom(v);
209  }
210  /**
211  Construct from a translation only, represented by a Translation3D class
212  and with an identity rotation
213  */
214  explicit Transform3D(const Translation3D<T> &t) { AssignFrom(t.Vect()); }
215 
216  //#if !defined(__MAKECINT__) && !defined(G__DICTIONARY) // this is ambigous with double * , double *
217 
218 
219 #ifdef OLD_VERSION
220  /**
221  Construct from a translation (using any type of DisplacementVector )
222  and then a rotation (any rotation object).
223  Requirement on the rotation and vector objects are that they can be transformed in a
224  Rotation3D class and in a Vector
225  */
226  template <class ARotation, class CoordSystem, class Tag>
227  Transform3D(const DisplacementVector3D<CoordSystem,Tag> & v , const ARotation & r)
228  {
229  // is equivalent from having first the rotation and then the translation vector rotated
230  Rotation3D r3d(r);
231  AssignFrom( r3d, r3d( Vector(v.X(),v.Y(),v.Z()) ) );
232  }
233 #endif
234 
235 public:
236  /**
237  Construct transformation from one coordinate system defined by three
238  points (origin + two axis) to
239  a new coordinate system defined by other three points (origin + axis)
240  Scalar version.
241  @param fr0 point defining origin of original reference system
242  @param fr1 point defining first axis of original reference system
243  @param fr2 point defining second axis of original reference system
244  @param to0 point defining origin of transformed reference system
245  @param to1 point defining first axis transformed reference system
246  @param to2 point defining second axis transformed reference system
247  */
248  template <typename SCALAR = T, typename std::enable_if<std::is_arithmetic<SCALAR>::value>::type * = nullptr>
249  Transform3D(const Point &fr0, const Point &fr1, const Point &fr2, const Point &to0, const Point &to1,
250  const Point &to2)
251  {
252  // takes impl. from CLHEP ( E.Chernyaev). To be checked
253 
254  Vector x1 = (fr1 - fr0).Unit();
255  Vector y1 = (fr2 - fr0).Unit();
256  Vector x2 = (to1 - to0).Unit();
257  Vector y2 = (to2 - to0).Unit();
258 
259  // C H E C K A N G L E S
260 
261  const T cos1 = x1.Dot(y1);
262  const T cos2 = x2.Dot(y2);
263 
264  if (std::fabs(T(1) - cos1) <= T(0.000001) || std::fabs(T(1) - cos2) <= T(0.000001)) {
265  std::cerr << "Transform3D: Error : zero angle between axes" << std::endl;
266  SetIdentity();
267  } else {
268  if (std::fabs(cos1 - cos2) > T(0.000001)) {
269  std::cerr << "Transform3D: Warning: angles between axes are not equal" << std::endl;
270  }
271 
272  // F I N D R O T A T I O N M A T R I X
273 
274  Vector z1 = (x1.Cross(y1)).Unit();
275  y1 = z1.Cross(x1);
276 
277  Vector z2 = (x2.Cross(y2)).Unit();
278  y2 = z2.Cross(x2);
279 
280  T x1x = x1.x();
281  T x1y = x1.y();
282  T x1z = x1.z();
283  T y1x = y1.x();
284  T y1y = y1.y();
285  T y1z = y1.z();
286  T z1x = z1.x();
287  T z1y = z1.y();
288  T z1z = z1.z();
289 
290  T x2x = x2.x();
291  T x2y = x2.y();
292  T x2z = x2.z();
293  T y2x = y2.x();
294  T y2y = y2.y();
295  T y2z = y2.z();
296  T z2x = z2.x();
297  T z2y = z2.y();
298  T z2z = z2.z();
299 
300  T detxx = (y1y * z1z - z1y * y1z);
301  T detxy = -(y1x * z1z - z1x * y1z);
302  T detxz = (y1x * z1y - z1x * y1y);
303  T detyx = -(x1y * z1z - z1y * x1z);
304  T detyy = (x1x * z1z - z1x * x1z);
305  T detyz = -(x1x * z1y - z1x * x1y);
306  T detzx = (x1y * y1z - y1y * x1z);
307  T detzy = -(x1x * y1z - y1x * x1z);
308  T detzz = (x1x * y1y - y1x * x1y);
309 
310  T txx = x2x * detxx + y2x * detyx + z2x * detzx;
311  T txy = x2x * detxy + y2x * detyy + z2x * detzy;
312  T txz = x2x * detxz + y2x * detyz + z2x * detzz;
313  T tyx = x2y * detxx + y2y * detyx + z2y * detzx;
314  T tyy = x2y * detxy + y2y * detyy + z2y * detzy;
315  T tyz = x2y * detxz + y2y * detyz + z2y * detzz;
316  T tzx = x2z * detxx + y2z * detyx + z2z * detzx;
317  T tzy = x2z * detxy + y2z * detyy + z2z * detzy;
318  T tzz = x2z * detxz + y2z * detyz + z2z * detzz;
319 
320  // S E T T R A N S F O R M A T I O N
321 
322  T dx1 = fr0.x(), dy1 = fr0.y(), dz1 = fr0.z();
323  T dx2 = to0.x(), dy2 = to0.y(), dz2 = to0.z();
324 
325  SetComponents(txx, txy, txz, dx2 - txx * dx1 - txy * dy1 - txz * dz1, tyx, tyy, tyz,
326  dy2 - tyx * dx1 - tyy * dy1 - tyz * dz1, tzx, tzy, tzz, dz2 - tzx * dx1 - tzy * dy1 - tzz * dz1);
327  }
328  }
329 
330  /**
331  Construct transformation from one coordinate system defined by three
332  points (origin + two axis) to
333  a new coordinate system defined by other three points (origin + axis)
334  Vectorised version.
335  @param fr0 point defining origin of original reference system
336  @param fr1 point defining first axis of original reference system
337  @param fr2 point defining second axis of original reference system
338  @param to0 point defining origin of transformed reference system
339  @param to1 point defining first axis transformed reference system
340  @param to2 point defining second axis transformed reference system
341  */
342  template <typename SCALAR = T, typename std::enable_if<!std::is_arithmetic<SCALAR>::value>::type * = nullptr>
343  Transform3D(const Point &fr0, const Point &fr1, const Point &fr2, const Point &to0, const Point &to1,
344  const Point &to2)
345  {
346  // takes impl. from CLHEP ( E.Chernyaev). To be checked
347 
348  Vector x1 = (fr1 - fr0).Unit();
349  Vector y1 = (fr2 - fr0).Unit();
350  Vector x2 = (to1 - to0).Unit();
351  Vector y2 = (to2 - to0).Unit();
352 
353  // C H E C K A N G L E S
354 
355  const T cos1 = x1.Dot(y1);
356  const T cos2 = x2.Dot(y2);
357 
358  const auto m1 = (abs(T(1) - cos1) <= T(0.000001) || abs(T(1) - cos2) <= T(0.000001));
359 
360  const auto m2 = (abs(cos1 - cos2) > T(0.000001));
361  if (any_of(m2)) {
362  std::cerr << "Transform3D: Warning: angles between axes are not equal" << std::endl;
363  }
364 
365  // F I N D R O T A T I O N M A T R I X
366 
367  Vector z1 = (x1.Cross(y1)).Unit();
368  y1 = z1.Cross(x1);
369 
370  Vector z2 = (x2.Cross(y2)).Unit();
371  y2 = z2.Cross(x2);
372 
373  T x1x = x1.x();
374  T x1y = x1.y();
375  T x1z = x1.z();
376  T y1x = y1.x();
377  T y1y = y1.y();
378  T y1z = y1.z();
379  T z1x = z1.x();
380  T z1y = z1.y();
381  T z1z = z1.z();
382 
383  T x2x = x2.x();
384  T x2y = x2.y();
385  T x2z = x2.z();
386  T y2x = y2.x();
387  T y2y = y2.y();
388  T y2z = y2.z();
389  T z2x = z2.x();
390  T z2y = z2.y();
391  T z2z = z2.z();
392 
393  T detxx = (y1y * z1z - z1y * y1z);
394  T detxy = -(y1x * z1z - z1x * y1z);
395  T detxz = (y1x * z1y - z1x * y1y);
396  T detyx = -(x1y * z1z - z1y * x1z);
397  T detyy = (x1x * z1z - z1x * x1z);
398  T detyz = -(x1x * z1y - z1x * x1y);
399  T detzx = (x1y * y1z - y1y * x1z);
400  T detzy = -(x1x * y1z - y1x * x1z);
401  T detzz = (x1x * y1y - y1x * x1y);
402 
403  T txx = x2x * detxx + y2x * detyx + z2x * detzx;
404  T txy = x2x * detxy + y2x * detyy + z2x * detzy;
405  T txz = x2x * detxz + y2x * detyz + z2x * detzz;
406  T tyx = x2y * detxx + y2y * detyx + z2y * detzx;
407  T tyy = x2y * detxy + y2y * detyy + z2y * detzy;
408  T tyz = x2y * detxz + y2y * detyz + z2y * detzz;
409  T tzx = x2z * detxx + y2z * detyx + z2z * detzx;
410  T tzy = x2z * detxy + y2z * detyy + z2z * detzy;
411  T tzz = x2z * detxz + y2z * detyz + z2z * detzz;
412 
413  // S E T T R A N S F O R M A T I O N
414 
415  T dx1 = fr0.x(), dy1 = fr0.y(), dz1 = fr0.z();
416  T dx2 = to0.x(), dy2 = to0.y(), dz2 = to0.z();
417 
418  SetComponents(txx, txy, txz, dx2 - txx * dx1 - txy * dy1 - txz * dz1, tyx, tyy, tyz,
419  dy2 - tyx * dx1 - tyy * dy1 - tyz * dz1, tzx, tzy, tzz, dz2 - tzx * dx1 - tzy * dy1 - tzz * dz1);
420 
421  if (any_of(m1)) {
422  std::cerr << "Transform3D: Error : zero angle between axes" << std::endl;
423  SetIdentity(m1);
424  }
425  }
426 
427  // use compiler generated copy ctor, copy assignmet and dtor
428 
429  /**
430  Construct from a linear algebra matrix of size at least 3x4,
431  which must support operator()(i,j) to obtain elements (0,0) thru (2,3).
432  The 3x3 sub-block is assumed to be the rotation part and the translations vector
433  are described by the 4-th column
434  */
435  template<class ForeignMatrix>
436  explicit Transform3D(const ForeignMatrix & m) {
437  SetComponents(m);
438  }
439 
440  /**
441  Raw constructor from 12 Scalar components
442  */
443  Transform3D(T xx, T xy, T xz, T dx, T yx, T yy, T yz, T dy, T zx, T zy, T zz, T dz)
444  {
445  SetComponents (xx, xy, xz, dx, yx, yy, yz, dy, zx, zy, zz, dz);
446  }
447 
448 
449  /**
450  Construct from a linear algebra matrix of size at least 3x4,
451  which must support operator()(i,j) to obtain elements (0,0) thru (2,3).
452  The 3x3 sub-block is assumed to be the rotation part and the translations vector
453  are described by the 4-th column
454  */
455  template <class ForeignMatrix>
456  Transform3D<T> &operator=(const ForeignMatrix &m)
457  {
458  SetComponents(m);
459  return *this;
460  }
461 
462 
463  // ======== Components ==============
464 
465 
466  /**
467  Set the 12 matrix components given an iterator to the start of
468  the desired data, and another to the end (12 past start).
469  */
470  template<class IT>
471 #ifndef NDEBUG
472  void SetComponents(IT begin, IT end) {
473 #else
474  void SetComponents(IT begin, IT ) {
475 #endif
476  for (int i = 0; i <12; ++i) {
477  fM[i] = *begin;
478  ++begin;
479  }
480  assert (end==begin);
481  }
482 
483  /**
484  Get the 12 matrix components into data specified by an iterator begin
485  and another to the end of the desired data (12 past start).
486  */
487  template<class IT>
488 #ifndef NDEBUG
489  void GetComponents(IT begin, IT end) const {
490 #else
491  void GetComponents(IT begin, IT ) const {
492 #endif
493  for (int i = 0; i <12; ++i) {
494  *begin = fM[i];
495  ++begin;
496  }
497  assert (end==begin);
498  }
499 
500  /**
501  Get the 12 matrix components into data specified by an iterator begin
502  */
503  template<class IT>
504  void GetComponents(IT begin) const {
505  std::copy(fM, fM + 12, begin);
506  }
507 
508  /**
509  Set components from a linear algebra matrix of size at least 3x4,
510  which must support operator()(i,j) to obtain elements (0,0) thru (2,3).
511  The 3x3 sub-block is assumed to be the rotation part and the translations vector
512  are described by the 4-th column
513  */
514  template<class ForeignMatrix>
515  void
516  SetTransformMatrix (const ForeignMatrix & m) {
517  fM[kXX]=m(0,0); fM[kXY]=m(0,1); fM[kXZ]=m(0,2); fM[kDX]=m(0,3);
518  fM[kYX]=m(1,0); fM[kYY]=m(1,1); fM[kYZ]=m(1,2); fM[kDY]=m(1,3);
519  fM[kZX]=m(2,0); fM[kZY]=m(2,1); fM[kZZ]=m(2,2); fM[kDZ]=m(2,3);
520  }
521 
522  /**
523  Get components into a linear algebra matrix of size at least 3x4,
524  which must support operator()(i,j) for write access to elements
525  (0,0) thru (2,3).
526  */
527  template<class ForeignMatrix>
528  void
529  GetTransformMatrix (ForeignMatrix & m) const {
530  m(0,0)=fM[kXX]; m(0,1)=fM[kXY]; m(0,2)=fM[kXZ]; m(0,3)=fM[kDX];
531  m(1,0)=fM[kYX]; m(1,1)=fM[kYY]; m(1,2)=fM[kYZ]; m(1,3)=fM[kDY];
532  m(2,0)=fM[kZX]; m(2,1)=fM[kZY]; m(2,2)=fM[kZZ]; m(2,3)=fM[kDZ];
533  }
534 
535 
536  /**
537  Set the components from 12 scalars
538  */
539  void SetComponents(T xx, T xy, T xz, T dx, T yx, T yy, T yz, T dy, T zx, T zy, T zz, T dz)
540  {
541  fM[kXX]=xx; fM[kXY]=xy; fM[kXZ]=xz; fM[kDX]=dx;
542  fM[kYX]=yx; fM[kYY]=yy; fM[kYZ]=yz; fM[kDY]=dy;
543  fM[kZX]=zx; fM[kZY]=zy; fM[kZZ]=zz; fM[kDZ]=dz;
544  }
545 
546  /**
547  Get the components into 12 scalars
548  */
549  void GetComponents(T &xx, T &xy, T &xz, T &dx, T &yx, T &yy, T &yz, T &dy, T &zx, T &zy, T &zz, T &dz) const
550  {
551  xx=fM[kXX]; xy=fM[kXY]; xz=fM[kXZ]; dx=fM[kDX];
552  yx=fM[kYX]; yy=fM[kYY]; yz=fM[kYZ]; dy=fM[kDY];
553  zx=fM[kZX]; zy=fM[kZY]; zz=fM[kZZ]; dz=fM[kDZ];
554  }
555 
556 
557  /**
558  Get the rotation and translation vector representing the 3D transformation
559  in any rotation and any vector (the Translation class could also be used)
560  */
561  template<class AnyRotation, class V>
562  void GetDecomposition(AnyRotation &r, V &v) const {
563  GetRotation(r);
564  GetTranslation(v);
565  }
566 
567 
568  /**
569  Get the rotation and translation vector representing the 3D transformation
570  */
571  void GetDecomposition(Rotation3D &r, Vector &v) const {
572  GetRotation(r);
573  GetTranslation(v);
574  }
575 
576  /**
577  Get the 3D rotation representing the 3D transformation
578  */
580  return Rotation3D( fM[kXX], fM[kXY], fM[kXZ],
581  fM[kYX], fM[kYY], fM[kYZ],
582  fM[kZX], fM[kZY], fM[kZZ] );
583  }
584 
585  /**
586  Get the rotation representing the 3D transformation
587  */
588  template <class AnyRotation>
589  AnyRotation Rotation() const {
590  return AnyRotation(Rotation3D(fM[kXX], fM[kXY], fM[kXZ], fM[kYX], fM[kYY], fM[kYZ], fM[kZX], fM[kZY], fM[kZZ]));
591  }
592 
593  /**
594  Get the rotation (any type) representing the 3D transformation
595  */
596  template <class AnyRotation>
597  void GetRotation(AnyRotation &r) const {
598  r = Rotation();
599  }
600 
601  /**
602  Get the translation representing the 3D transformation in a Cartesian vector
603  */
604  Translation3D<T> Translation() const { return Translation3D<T>(fM[kDX], fM[kDY], fM[kDZ]); }
605 
606  /**
607  Get the translation representing the 3D transformation in any vector
608  which implements the SetXYZ method
609  */
610  template <class AnyVector>
611  void GetTranslation(AnyVector &v) const {
612  v.SetXYZ(fM[kDX], fM[kDY], fM[kDZ]);
613  }
614 
615 
616 
617  // operations on points and vectors
618 
619  /**
620  Transformation operation for Position Vector in Cartesian coordinate
621  For a Position Vector first a rotation and then a translation is applied
622  */
623  Point operator() (const Point & p) const {
624  return Point ( fM[kXX]*p.X() + fM[kXY]*p.Y() + fM[kXZ]*p.Z() + fM[kDX],
625  fM[kYX]*p.X() + fM[kYY]*p.Y() + fM[kYZ]*p.Z() + fM[kDY],
626  fM[kZX]*p.X() + fM[kZY]*p.Y() + fM[kZZ]*p.Z() + fM[kDZ] );
627  }
628 
629 
630  /**
631  Transformation operation for Displacement Vectors in Cartesian coordinate
632  For the Displacement Vectors only the rotation applies - no translations
633  */
634  Vector operator() (const Vector & v) const {
635  return Vector( fM[kXX]*v.X() + fM[kXY]*v.Y() + fM[kXZ]*v.Z() ,
636  fM[kYX]*v.X() + fM[kYY]*v.Y() + fM[kYZ]*v.Z() ,
637  fM[kZX]*v.X() + fM[kZY]*v.Y() + fM[kZZ]*v.Z() );
638  }
639 
640 
641  /**
642  Transformation operation for Position Vector in any coordinate system
643  */
644  template <class CoordSystem>
646  {
648  }
649  /**
650  Transformation operation for Position Vector in any coordinate system
651  */
652  template <class CoordSystem>
654  {
655  return operator()(v);
656  }
657 
658  /**
659  Transformation operation for Displacement Vector in any coordinate system
660  */
661  template<class CoordSystem >
664  }
665  /**
666  Transformation operation for Displacement Vector in any coordinate system
667  */
668  template <class CoordSystem>
670  {
671  return operator()(v);
672  }
673 
674  /**
675  Directly apply the inverse affine transformation on vectors.
676  Avoids having to calculate the inverse as an intermediate result.
677  This is possible since the inverse of a rotation is its transpose.
678  */
679  Vector ApplyInverse(const Vector &v) const
680  {
681  return Vector(fM[kXX] * v.X() + fM[kYX] * v.Y() + fM[kZX] * v.Z(),
682  fM[kXY] * v.X() + fM[kYY] * v.Y() + fM[kZY] * v.Z(),
683  fM[kXZ] * v.X() + fM[kYZ] * v.Y() + fM[kZZ] * v.Z());
684  }
685 
686  /**
687  Directly apply the inverse affine transformation on points
688  (first inverse translation then inverse rotation).
689  Avoids having to calculate the inverse as an intermediate result.
690  This is possible since the inverse of a rotation is its transpose.
691  */
692  Point ApplyInverse(const Point &p) const
693  {
694  Point tmp(p.X() - fM[kDX], p.Y() - fM[kDY], p.Z() - fM[kDZ]);
695  return Point(fM[kXX] * tmp.X() + fM[kYX] * tmp.Y() + fM[kZX] * tmp.Z(),
696  fM[kXY] * tmp.X() + fM[kYY] * tmp.Y() + fM[kZY] * tmp.Z(),
697  fM[kXZ] * tmp.X() + fM[kYZ] * tmp.Y() + fM[kZZ] * tmp.Z());
698  }
699 
700  /**
701  Directly apply the inverse affine transformation on an arbitrary
702  coordinate-system point.
703  Involves casting to Point(p) type.
704  */
705  template <class CoordSystem>
707  {
709  }
710 
711  /**
712  Directly apply the inverse affine transformation on an arbitrary
713  coordinate-system vector.
714  Involves casting to Vector(p) type.
715  */
716  template <class CoordSystem>
718  {
720  }
721 
722  /**
723  Transformation operation for points between different coordinate system tags
724  */
725  template <class CoordSystem, class Tag1, class Tag2>
727  {
728  const Point xyzNew = operator()(Point(p1.X(), p1.Y(), p1.Z()));
729  p2.SetXYZ( xyzNew.X(), xyzNew.Y(), xyzNew.Z() );
730  }
731 
732 
733  /**
734  Transformation operation for Displacement Vector of different coordinate systems
735  */
736  template <class CoordSystem, class Tag1, class Tag2>
738  {
739  const Vector xyzNew = operator()(Vector(v1.X(), v1.Y(), v1.Z()));
740  v2.SetXYZ( xyzNew.X(), xyzNew.Y(), xyzNew.Z() );
741  }
742 
743  /**
744  Transformation operation for a Lorentz Vector in any coordinate system
745  */
746  template <class CoordSystem >
748  const Vector xyzNew = operator()(Vector(q.Vect()));
749  return LorentzVector<CoordSystem>(xyzNew.X(), xyzNew.Y(), xyzNew.Z(), q.E());
750  }
751  /**
752  Transformation operation for a Lorentz Vector in any coordinate system
753  */
754  template <class CoordSystem>
756  {
757  return operator()(q);
758  }
759 
760  /**
761  Transformation on a 3D plane
762  */
763  template <typename TYPE>
765  {
766  // transformations on a 3D plane
767  const auto n = plane.Normal();
768  // take a point on the plane. Use origin projection on the plane
769  // ( -ad, -bd, -cd) if (a**2 + b**2 + c**2 ) = 1
770  const auto d = plane.HesseDistance();
771  Point p(-d * n.X(), -d * n.Y(), -d * n.Z());
772  return Plane3D<TYPE>(operator()(n), operator()(p));
773  }
774 
775  /// Multiplication operator for 3D plane
776  template <typename TYPE>
778  {
779  return operator()(plane);
780  }
781 
782  // skip transformation for arbitrary vectors - not really defined if point or displacement vectors
783 
784  /**
785  multiply (combine) with another transformation in place
786  */
787  inline Transform3D<T> &operator*=(const Transform3D<T> &t);
788 
789  /**
790  multiply (combine) two transformations
791  */
792  inline Transform3D<T> operator*(const Transform3D<T> &t) const;
793 
794  /**
795  Invert the transformation in place (scalar)
796  */
797  template <typename SCALAR = T, typename std::enable_if<std::is_arithmetic<SCALAR>::value>::type * = nullptr>
798  void Invert()
799  {
800  //
801  // Name: Transform3D::inverse Date: 24.09.96
802  // Author: E.Chernyaev (IHEP/Protvino) Revised:
803  //
804  // Function: Find inverse affine transformation.
805 
806  T detxx = fM[kYY] * fM[kZZ] - fM[kYZ] * fM[kZY];
807  T detxy = fM[kYX] * fM[kZZ] - fM[kYZ] * fM[kZX];
808  T detxz = fM[kYX] * fM[kZY] - fM[kYY] * fM[kZX];
809  T det = fM[kXX] * detxx - fM[kXY] * detxy + fM[kXZ] * detxz;
810  if (det == T(0)) {
811  std::cerr << "Transform3D::inverse error: zero determinant" << std::endl;
812  return;
813  }
814  det = T(1) / det;
815  detxx *= det;
816  detxy *= det;
817  detxz *= det;
818  T detyx = (fM[kXY] * fM[kZZ] - fM[kXZ] * fM[kZY]) * det;
819  T detyy = (fM[kXX] * fM[kZZ] - fM[kXZ] * fM[kZX]) * det;
820  T detyz = (fM[kXX] * fM[kZY] - fM[kXY] * fM[kZX]) * det;
821  T detzx = (fM[kXY] * fM[kYZ] - fM[kXZ] * fM[kYY]) * det;
822  T detzy = (fM[kXX] * fM[kYZ] - fM[kXZ] * fM[kYX]) * det;
823  T detzz = (fM[kXX] * fM[kYY] - fM[kXY] * fM[kYX]) * det;
824  SetComponents(detxx, -detyx, detzx, -detxx * fM[kDX] + detyx * fM[kDY] - detzx * fM[kDZ], -detxy, detyy, -detzy,
825  detxy * fM[kDX] - detyy * fM[kDY] + detzy * fM[kDZ], detxz, -detyz, detzz,
826  -detxz * fM[kDX] + detyz * fM[kDY] - detzz * fM[kDZ]);
827  }
828 
829  /**
830  Invert the transformation in place (vectorised)
831  */
832  template <typename SCALAR = T, typename std::enable_if<!std::is_arithmetic<SCALAR>::value>::type * = nullptr>
833  void Invert()
834  {
835  //
836  // Name: Transform3D::inverse Date: 24.09.96
837  // Author: E.Chernyaev (IHEP/Protvino) Revised:
838  //
839  // Function: Find inverse affine transformation.
840 
841  T detxx = fM[kYY] * fM[kZZ] - fM[kYZ] * fM[kZY];
842  T detxy = fM[kYX] * fM[kZZ] - fM[kYZ] * fM[kZX];
843  T detxz = fM[kYX] * fM[kZY] - fM[kYY] * fM[kZX];
844  T det = fM[kXX] * detxx - fM[kXY] * detxy + fM[kXZ] * detxz;
845  const auto detZmask = (det == T(0));
846  if (any_of(detZmask)) {
847  std::cerr << "Transform3D::inverse error: zero determinant" << std::endl;
848  det(detZmask) = T(1);
849  }
850  det = T(1) / det;
851  detxx *= det;
852  detxy *= det;
853  detxz *= det;
854  T detyx = (fM[kXY] * fM[kZZ] - fM[kXZ] * fM[kZY]) * det;
855  T detyy = (fM[kXX] * fM[kZZ] - fM[kXZ] * fM[kZX]) * det;
856  T detyz = (fM[kXX] * fM[kZY] - fM[kXY] * fM[kZX]) * det;
857  T detzx = (fM[kXY] * fM[kYZ] - fM[kXZ] * fM[kYY]) * det;
858  T detzy = (fM[kXX] * fM[kYZ] - fM[kXZ] * fM[kYX]) * det;
859  T detzz = (fM[kXX] * fM[kYY] - fM[kXY] * fM[kYX]) * det;
860  // Set det=0 cases to 0
861  if (any_of(detZmask)) {
862  detxx(detZmask) = T(0);
863  detxy(detZmask) = T(0);
864  detxz(detZmask) = T(0);
865  detyx(detZmask) = T(0);
866  detyy(detZmask) = T(0);
867  detyz(detZmask) = T(0);
868  detzx(detZmask) = T(0);
869  detzy(detZmask) = T(0);
870  detzz(detZmask) = T(0);
871  }
872  // set final components
873  SetComponents(detxx, -detyx, detzx, -detxx * fM[kDX] + detyx * fM[kDY] - detzx * fM[kDZ], -detxy, detyy, -detzy,
874  detxy * fM[kDX] - detyy * fM[kDY] + detzy * fM[kDZ], detxz, -detyz, detzz,
875  -detxz * fM[kDX] + detyz * fM[kDY] - detzz * fM[kDZ]);
876  }
877 
878  /**
879  Return the inverse of the transformation.
880  */
882  {
883  Transform3D<T> t(*this);
884  t.Invert();
885  return t;
886  }
887 
888  /**
889  Equality operator. Check equality for each element
890  To do: use T tolerance
891  */
892  bool operator==(const Transform3D<T> &rhs) const
893  {
894  return (fM[0] == rhs.fM[0] && fM[1] == rhs.fM[1] && fM[2] == rhs.fM[2] && fM[3] == rhs.fM[3] &&
895  fM[4] == rhs.fM[4] && fM[5] == rhs.fM[5] && fM[6] == rhs.fM[6] && fM[7] == rhs.fM[7] &&
896  fM[8] == rhs.fM[8] && fM[9] == rhs.fM[9] && fM[10] == rhs.fM[10] && fM[11] == rhs.fM[11]);
897  }
898 
899  /**
900  Inequality operator. Check equality for each element
901  To do: use T tolerance
902  */
903  bool operator!=(const Transform3D<T> &rhs) const { return !operator==(rhs); }
904 
905 protected:
906 
907  /**
908  make transformation from first a rotation then a translation
909  */
910  void AssignFrom(const Rotation3D &r, const Vector &v)
911  {
912  // assignment from rotation + translation
913 
914  T rotData[9];
915  r.GetComponents(rotData, rotData + 9);
916  // first raw
917  for (int i = 0; i < 3; ++i) fM[i] = rotData[i];
918  // second raw
919  for (int i = 0; i < 3; ++i) fM[kYX + i] = rotData[3 + i];
920  // third raw
921  for (int i = 0; i < 3; ++i) fM[kZX + i] = rotData[6 + i];
922 
923  // translation data
924  T vecData[3];
925  v.GetCoordinates(vecData, vecData + 3);
926  fM[kDX] = vecData[0];
927  fM[kDY] = vecData[1];
928  fM[kDZ] = vecData[2];
929  }
930 
931  /**
932  make transformation from only rotations (zero translation)
933  */
934  void AssignFrom(const Rotation3D &r)
935  {
936  // assign from only a rotation (null translation)
937  T rotData[9];
938  r.GetComponents(rotData, rotData + 9);
939  for (int i = 0; i < 3; ++i) {
940  for (int j = 0; j < 3; ++j) fM[4 * i + j] = rotData[3 * i + j];
941  // empty vector data
942  fM[4 * i + 3] = T(0);
943  }
944  }
945 
946  /**
947  make transformation from only translation (identity rotations)
948  */
949  void AssignFrom(const Vector &v)
950  {
951  // assign from a translation only (identity rotations)
952  fM[kXX] = T(1);
953  fM[kXY] = T(0);
954  fM[kXZ] = T(0);
955  fM[kDX] = v.X();
956  fM[kYX] = T(0);
957  fM[kYY] = T(1);
958  fM[kYZ] = T(0);
959  fM[kDY] = v.Y();
960  fM[kZX] = T(0);
961  fM[kZY] = T(0);
962  fM[kZZ] = T(1);
963  fM[kDZ] = v.Z();
964  }
965 
966  /**
967  Set identity transformation (identity rotation , zero translation)
968  */
969  void SetIdentity()
970  {
971  // set identity ( identity rotation and zero translation)
972  fM[kXX] = T(1);
973  fM[kXY] = T(0);
974  fM[kXZ] = T(0);
975  fM[kDX] = T(0);
976  fM[kYX] = T(0);
977  fM[kYY] = T(1);
978  fM[kYZ] = T(0);
979  fM[kDY] = T(0);
980  fM[kZX] = T(0);
981  fM[kZY] = T(0);
982  fM[kZZ] = T(1);
983  fM[kDZ] = T(0);
984  }
985 
986  /**
987  Set identity transformation (identity rotation , zero translation)
988  vectorised version that sets using a mask
989  */
990  template <typename SCALAR = T, typename std::enable_if<!std::is_arithmetic<SCALAR>::value>::type * = nullptr>
991  void SetIdentity(const typename SCALAR::mask_type m)
992  {
993  // set identity ( identity rotation and zero translation)
994  fM[kXX](m) = T(1);
995  fM[kXY](m) = T(0);
996  fM[kXZ](m) = T(0);
997  fM[kDX](m) = T(0);
998  fM[kYX](m) = T(0);
999  fM[kYY](m) = T(1);
1000  fM[kYZ](m) = T(0);
1001  fM[kDY](m) = T(0);
1002  fM[kZX](m) = T(0);
1003  fM[kZY](m) = T(0);
1004  fM[kZZ](m) = T(1);
1005  fM[kDZ](m) = T(0);
1006  }
1007 
1008 private:
1009  T fM[12]; // transformation elements (3x4 matrix)
1010 };
1011 
1012 
1013 
1014 
1015 // inline functions (combination of transformations)
1016 
1017 template <class T>
1019 {
1020  // combination of transformations
1021 
1022  SetComponents(fM[kXX]*t.fM[kXX]+fM[kXY]*t.fM[kYX]+fM[kXZ]*t.fM[kZX],
1023  fM[kXX]*t.fM[kXY]+fM[kXY]*t.fM[kYY]+fM[kXZ]*t.fM[kZY],
1024  fM[kXX]*t.fM[kXZ]+fM[kXY]*t.fM[kYZ]+fM[kXZ]*t.fM[kZZ],
1025  fM[kXX]*t.fM[kDX]+fM[kXY]*t.fM[kDY]+fM[kXZ]*t.fM[kDZ]+fM[kDX],
1026 
1027  fM[kYX]*t.fM[kXX]+fM[kYY]*t.fM[kYX]+fM[kYZ]*t.fM[kZX],
1028  fM[kYX]*t.fM[kXY]+fM[kYY]*t.fM[kYY]+fM[kYZ]*t.fM[kZY],
1029  fM[kYX]*t.fM[kXZ]+fM[kYY]*t.fM[kYZ]+fM[kYZ]*t.fM[kZZ],
1030  fM[kYX]*t.fM[kDX]+fM[kYY]*t.fM[kDY]+fM[kYZ]*t.fM[kDZ]+fM[kDY],
1031 
1032  fM[kZX]*t.fM[kXX]+fM[kZY]*t.fM[kYX]+fM[kZZ]*t.fM[kZX],
1033  fM[kZX]*t.fM[kXY]+fM[kZY]*t.fM[kYY]+fM[kZZ]*t.fM[kZY],
1034  fM[kZX]*t.fM[kXZ]+fM[kZY]*t.fM[kYZ]+fM[kZZ]*t.fM[kZZ],
1035  fM[kZX]*t.fM[kDX]+fM[kZY]*t.fM[kDY]+fM[kZZ]*t.fM[kDZ]+fM[kDZ]);
1036 
1037  return *this;
1038 }
1039 
1040 template <class T>
1042 {
1043  // combination of transformations
1044 
1045  return Transform3D<T>(fM[kXX] * t.fM[kXX] + fM[kXY] * t.fM[kYX] + fM[kXZ] * t.fM[kZX],
1046  fM[kXX] * t.fM[kXY] + fM[kXY] * t.fM[kYY] + fM[kXZ] * t.fM[kZY],
1047  fM[kXX] * t.fM[kXZ] + fM[kXY] * t.fM[kYZ] + fM[kXZ] * t.fM[kZZ],
1048  fM[kXX] * t.fM[kDX] + fM[kXY] * t.fM[kDY] + fM[kXZ] * t.fM[kDZ] + fM[kDX],
1049 
1050  fM[kYX] * t.fM[kXX] + fM[kYY] * t.fM[kYX] + fM[kYZ] * t.fM[kZX],
1051  fM[kYX] * t.fM[kXY] + fM[kYY] * t.fM[kYY] + fM[kYZ] * t.fM[kZY],
1052  fM[kYX] * t.fM[kXZ] + fM[kYY] * t.fM[kYZ] + fM[kYZ] * t.fM[kZZ],
1053  fM[kYX] * t.fM[kDX] + fM[kYY] * t.fM[kDY] + fM[kYZ] * t.fM[kDZ] + fM[kDY],
1054 
1055  fM[kZX] * t.fM[kXX] + fM[kZY] * t.fM[kYX] + fM[kZZ] * t.fM[kZX],
1056  fM[kZX] * t.fM[kXY] + fM[kZY] * t.fM[kYY] + fM[kZZ] * t.fM[kZY],
1057  fM[kZX] * t.fM[kXZ] + fM[kZY] * t.fM[kYZ] + fM[kZZ] * t.fM[kZZ],
1058  fM[kZX] * t.fM[kDX] + fM[kZY] * t.fM[kDY] + fM[kZZ] * t.fM[kDZ] + fM[kDZ]);
1059 }
1060 
1061 
1062 
1063 
1064 //--- global functions resulting in Transform3D -------
1065 
1066 
1067 // ------ combination of a translation (first) and a rotation ------
1068 
1069 
1070 /**
1071  combine a translation and a rotation to give a transform3d
1072  First the translation then the rotation
1073  */
1074 template <class T>
1076 {
1077  return Transform3D<T>(r, r(t.Vect()));
1078 }
1079 template <class T>
1081 {
1082  Rotation3D r3(r);
1083  return Transform3D<T>(r3, r3(t.Vect()));
1084 }
1085 template <class T>
1087 {
1088  Rotation3D r3(r);
1089  return Transform3D<T>(r3, r3(t.Vect()));
1090 }
1091 template <class T>
1093 {
1094  Rotation3D r3(r);
1095  return Transform3D<T>(r3, r3(t.Vect()));
1096 }
1097 template <class T>
1099 {
1100  Rotation3D r3(r);
1101  return Transform3D<T>(r3, r3(t.Vect()));
1102 }
1103 template <class T>
1105 {
1106  Rotation3D r3(r);
1107  return Transform3D<T>(r3, r3(t.Vect()));
1108 }
1109 template <class T>
1111 {
1112  Rotation3D r3(r);
1113  return Transform3D<T>(r3, r3(t.Vect()));
1114 }
1115 template <class T>
1117 {
1118  Rotation3D r3(r);
1119  return Transform3D<T>(r3, r3(t.Vect()));
1120 }
1121 
1122 // ------ combination of a rotation (first) and then a translation ------
1123 
1124 /**
1125  combine a rotation and a translation to give a transform3d
1126  First a rotation then the translation
1127  */
1128 template <class T>
1130 {
1131  return Transform3D<T>(r, t.Vect());
1132 }
1133 template <class T>
1135 {
1136  return Transform3D<T>(Rotation3D(r), t.Vect());
1137 }
1138 template <class T>
1140 {
1141  return Transform3D<T>(Rotation3D(r), t.Vect());
1142 }
1143 template <class T>
1145 {
1146  return Transform3D<T>(Rotation3D(r), t.Vect());
1147 }
1148 template <class T>
1150 {
1151  return Transform3D<T>(Rotation3D(r), t.Vect());
1152 }
1153 template <class T>
1155 {
1156  return Transform3D<T>(Rotation3D(r), t.Vect());
1157 }
1158 template <class T>
1160 {
1161  return Transform3D<T>(Rotation3D(r), t.Vect());
1162 }
1163 template <class T>
1165 {
1166  return Transform3D<T>(Rotation3D(r), t.Vect());
1167 }
1168 
1169 // ------ combination of a Transform3D and a pure translation------
1170 
1171 /**
1172  combine a transformation and a translation to give a transform3d
1173  First the translation then the transform3D
1174  */
1175 template <class T>
1177 {
1178  Rotation3D r = t.Rotation();
1179  return Transform3D<T>(r, r(d.Vect()) + t.Translation().Vect());
1180 }
1181 
1182 /**
1183  combine a translation and a transformation to give a transform3d
1184  First the transformation then the translation
1185  */
1186 template <class T>
1188 {
1189  return Transform3D<T>(t.Rotation(), t.Translation().Vect() + d.Vect());
1190 }
1191 
1192 // ------ combination of a Transform3D and any rotation------
1193 
1194 
1195 /**
1196  combine a transformation and a rotation to give a transform3d
1197  First the rotation then the transform3D
1198  */
1199 template <class T>
1201 {
1202  return Transform3D<T>(t.Rotation() * r, t.Translation());
1203 }
1204 template <class T>
1206 {
1207  return Transform3D<T>(t.Rotation() * r, t.Translation());
1208 }
1209 template <class T>
1211 {
1212  return Transform3D<T>(t.Rotation() * r, t.Translation());
1213 }
1214 template <class T>
1216 {
1217  return Transform3D<T>(t.Rotation() * r, t.Translation());
1218 }
1219 template <class T>
1221 {
1222  return Transform3D<T>(t.Rotation() * r, t.Translation());
1223 }
1224 template <class T>
1226 {
1227  return Transform3D<T>(t.Rotation() * r, t.Translation());
1228 }
1229 template <class T>
1231 {
1232  return Transform3D<T>(t.Rotation() * r, t.Translation());
1233 }
1234 template <class T>
1236 {
1237  return Transform3D<T>(t.Rotation() * r, t.Translation());
1238 }
1239 
1240 
1241 
1242 /**
1243  combine a rotation and a transformation to give a transform3d
1244  First the transformation then the rotation
1245  */
1246 template <class T>
1248 {
1249  return Transform3D<T>(r * t.Rotation(), r * t.Translation().Vect());
1250 }
1251 template <class T>
1253 {
1254  Rotation3D r3d(r);
1255  return Transform3D<T>(r3d * t.Rotation(), r3d * t.Translation().Vect());
1256 }
1257 template <class T>
1259 {
1260  Rotation3D r3d(r);
1261  return Transform3D<T>(r3d * t.Rotation(), r3d * t.Translation().Vect());
1262 }
1263 template <class T>
1265 {
1266  Rotation3D r3d(r);
1267  return Transform3D<T>(r3d * t.Rotation(), r3d * t.Translation().Vect());
1268 }
1269 template <class T>
1271 {
1272  Rotation3D r3d(r);
1273  return Transform3D<T>(r3d * t.Rotation(), r3d * t.Translation().Vect());
1274 }
1275 template <class T>
1277 {
1278  Rotation3D r3d(r);
1279  return Transform3D<T>(r3d * t.Rotation(), r3d * t.Translation().Vect());
1280 }
1281 template <class T>
1283 {
1284  Rotation3D r3d(r);
1285  return Transform3D<T>(r3d * t.Rotation(), r3d * t.Translation().Vect());
1286 }
1287 template <class T>
1289 {
1290  Rotation3D r3d(r);
1291  return Transform3D<T>(r3d * t.Rotation(), r3d * t.Translation().Vect());
1292 }
1293 
1294 
1295 //---I/O functions
1296 // TODO - I/O should be put in the manipulator form
1297 
1298 /**
1299  print the 12 components of the Transform3D
1300  */
1301 template <class T>
1302 std::ostream &operator<<(std::ostream &os, const Transform3D<T> &t)
1303 {
1304  // TODO - this will need changing for machine-readable issues
1305  // and even the human readable form needs formatting improvements
1306 
1307  T m[12];
1308  t.GetComponents(m, m + 12);
1309  os << "\n" << m[0] << " " << m[1] << " " << m[2] << " " << m[3];
1310  os << "\n" << m[4] << " " << m[5] << " " << m[6] << " " << m[7];
1311  os << "\n" << m[8] << " " << m[9] << " " << m[10] << " " << m[11] << "\n";
1312  return os;
1313 }
1314 
1315 } // end namespace Impl
1316 
1317 // typedefs for double and float versions
1320 
1321 } // end namespace Math
1322 
1323 } // end namespace ROOT
1324 
1325 
1326 #endif /* ROOT_Math_GenVector_Transform3D */
void AssignFrom(const Rotation3D &r, const Vector &v)
make transformation from first a rotation then a translation
Definition: Transform3D.h:910
Class describing a generic LorentzVector in the 4D space-time, using the specified coordinate system ...
Definition: LorentzVector.h:48
Transform3D(IT begin, IT end)
Construct given a pair of pointers or iterators defining the beginning and end of an array of 12 Scal...
Definition: Transform3D.h:107
bool operator==(const Transform3D< T > &rhs) const
Equality operator.
Definition: Transform3D.h:892
Point ApplyInverse(const Point &p) const
Directly apply the inverse affine transformation on points (first inverse translation then inverse ro...
Definition: Transform3D.h:692
void AssignFrom(const Rotation3D &r)
make transformation from only rotations (zero translation)
Definition: Transform3D.h:934
Scalar HesseDistance() const
Return the Hesse Distance (distance from the origin) of the plane or the d coefficient expressed in n...
Definition: Plane3D.h:168
LorentzVector< CoordSystem > operator*(const LorentzVector< CoordSystem > &q) const
Transformation operation for a Lorentz Vector in any coordinate system.
Definition: Transform3D.h:755
DisplacementVector3D< CoordSystem, Tag > & SetXYZ(Scalar a, Scalar b, Scalar c)
set the values of the vector from the cartesian components (x,y,z) (if the vector is held in polar or...
Namespace for new ROOT classes and functions.
Definition: StringConv.hxx:21
Transform3D(const Vector &v)
Construct from a translation only, represented by a Cartesian 3D Vector, and with an identity rotatio...
Definition: Transform3D.h:207
void GetComponents(ForeignVector &v1, ForeignVector &v2, ForeignVector &v3) const
Get components into three vectors which will be the (orthonormal) columns of the rotation matrix...
Definition: Rotation3D.h:249
void SetIdentity(const typename SCALAR::mask_type m)
Set identity transformation (identity rotation , zero translation) vectorised version that sets using...
Definition: Transform3D.h:991
Plane3D< TYPE > operator*(const Plane3D< TYPE > &plane) const
Multiplication operator for 3D plane.
Definition: Transform3D.h:777
Vector Normal() const
Return normal vector to the plane as Cartesian DisplacementVector.
Definition: Plane3D.h:162
void GetDecomposition(AnyRotation &r, V &v) const
Get the rotation and translation vector representing the 3D transformation in any rotation and any ve...
Definition: Transform3D.h:562
double T(double x)
Definition: ChebyshevPol.h:34
Rotation class representing a 3D rotation about the Z axis by the angle of rotation.
Definition: RotationZ.h:43
void GetComponents(IT begin, IT end) const
Get the 12 matrix components into data specified by an iterator begin and another to the end of the d...
Definition: Transform3D.h:489
void GetComponents(T &xx, T &xy, T &xz, T &dx, T &yx, T &yy, T &yz, T &dy, T &zx, T &zy, T &zz, T &dz) const
Get the components into 12 scalars.
Definition: Transform3D.h:549
Impl::Transform3D< float > Transform3DF
Definition: Transform3D.h:1319
void SetComponents(IT begin, IT end)
Set the 12 matrix components given an iterator to the start of the desired data, and another to the e...
Definition: Transform3D.h:472
Class describing a generic position vector (point) in 3 dimensions.
void GetDecomposition(Rotation3D &r, Vector &v) const
Get the rotation and translation vector representing the 3D transformation.
Definition: Transform3D.h:571
Transform3D(const Point &fr0, const Point &fr1, const Point &fr2, const Point &to0, const Point &to1, const Point &to2)
Construct transformation from one coordinate system defined by three points (origin + two axis) to a ...
Definition: Transform3D.h:249
Rotation class with the (3D) rotation represented by angles describing first a rotation of an angle p...
Definition: RotationZYX.h:59
Translation3D< T > Translation() const
Get the translation representing the 3D transformation in a Cartesian vector.
Definition: Transform3D.h:604
void SetTransformMatrix(const ForeignMatrix &m)
Set components from a linear algebra matrix of size at least 3x4, which must support operator()(i...
Definition: Transform3D.h:516
Rotation class with the (3D) rotation represented by a unit quaternion (u, i, j, k).
Definition: Quaternion.h:47
Transform3D(const DisplacementVector3D< CoordSystem, Tag > &v)
Construct from a translation only, represented by any DisplacementVector3D and with an identity rotat...
Definition: Transform3D.h:200
void Transform(const PositionVector3D< CoordSystem, Tag1 > &p1, PositionVector3D< CoordSystem, Tag2 > &p2) const
Transformation operation for points between different coordinate system tags.
Definition: Transform3D.h:726
void GetComponents(IT begin) const
Get the 12 matrix components into data specified by an iterator begin.
Definition: Transform3D.h:504
Transform3D(const Rotation3D &r)
Construct from a 3D Rotation only with zero translation.
Definition: Transform3D.h:163
Transform3D(const ARotation &r, const DisplacementVector3D< CoordSystem, Tag > &v)
Construct from a rotation (any rotation object) and then a translation (represented by any Displaceme...
Definition: Transform3D.h:131
Transform3D< T > & operator=(const ForeignMatrix &m)
Construct from a linear algebra matrix of size at least 3x4, which must support operator()(i,j) to obtain elements (0,0) thru (2,3).
Definition: Transform3D.h:456
Transform3D(const ARotation &r, const Translation3D< T > &t)
Construct from a rotation (any rotation object) and then a translation represented by a Translation3D...
Definition: Transform3D.h:143
AxisAngle class describing rotation represented with direction axis (3D Vector) and an angle of rotat...
Definition: AxisAngle.h:41
Scalar X() const
Cartesian X, converting if necessary from internal coordinate system.
static const double x2[5]
Transform3D(const Rotation3D &r, const Translation3D< T > &t)
Construct from a rotation and then a translation described by a Translation3D class.
Definition: Transform3D.h:122
void GetTransformMatrix(ForeignMatrix &m) const
Get components into a linear algebra matrix of size at least 3x4, which must support operator()(i...
Definition: Transform3D.h:529
Rotation class representing a 3D rotation about the Y axis by the angle of rotation.
Definition: RotationY.h:43
void AssignFrom(const Vector &v)
make transformation from only translation (identity rotations)
Definition: Transform3D.h:949
unsigned int r3[N_CITIES]
Definition: simanTSP.cxx:323
Scalar Z() const
Cartesian Z, converting if necessary from internal coordinate system.
Scalar Y() const
Cartesian Y, converting if necessary from internal coordinate system.
void SetIdentity()
Set identity transformation (identity rotation , zero translation)
Definition: Transform3D.h:969
static double p2(double t, double a, double b, double c)
Vector ApplyInverse(const Vector &v) const
Directly apply the inverse affine transformation on vectors.
Definition: Transform3D.h:679
Scalar E() const
return 4-th component (time, or energy for a 4-momentum vector)
void Invert()
Invert the transformation in place (scalar)
Definition: Transform3D.h:798
Transform3D(const RotationZ &r)
Definition: Transform3D.h:191
Class describing a 3 dimensional translation.
Definition: Translation3D.h:51
Scalar Z() const
Cartesian Z, converting if necessary from internal coordinate system.
bool operator!=(const Transform3D< T > &rhs) const
Inequality operator.
Definition: Transform3D.h:903
Class describing a generic displacement vector in 3 dimensions.
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
TRandom2 r(17)
SVector< double, 2 > v
Definition: Dict.h:5
XPoint xy[kMAXMK]
Definition: TGX11.cxx:122
void GetTranslation(AnyVector &v) const
Get the translation representing the 3D transformation in any vector which implements the SetXYZ meth...
Definition: Transform3D.h:611
Scalar X() const
Cartesian X, converting if necessary from internal coordinate system.
Rotation class representing a 3D rotation about the X axis by the angle of rotation.
Definition: RotationX.h:43
DisplacementVector3D Cross(const DisplacementVector3D< OtherCoords, Tag > &v) const
Return vector (cross) product of two displacement vectors, as a vector in the coordinate system of th...
Rotation class with the (3D) rotation represented by a 3x3 orthogonal matrix.
Definition: Rotation3D.h:65
const Vector & Vect() const
return a const reference to the underline vector representing the translation
TMarker * m
Definition: textangle.C:8
Transform3D(const ForeignMatrix &m)
Construct from a linear algebra matrix of size at least 3x4, which must support operator()(i,j) to obtain elements (0,0) thru (2,3).
Definition: Transform3D.h:436
static double p1(double t, double a, double b)
void SetComponents(T xx, T xy, T xz, T dx, T yx, T yy, T yz, T dy, T zx, T zy, T zz, T dz)
Set the components from 12 scalars.
Definition: Transform3D.h:539
void GetRotation(AnyRotation &r) const
Get the rotation (any type) representing the 3D transformation.
Definition: Transform3D.h:597
void GetCoordinates(Scalar &a, Scalar &b, Scalar &c) const
get internal data into 3 Scalar numbers
Point operator()(const Point &p) const
Transformation operation for Position Vector in Cartesian coordinate For a Position Vector first a ro...
Definition: Transform3D.h:623
PositionVector3D< CoordSystem > operator()(const PositionVector3D< CoordSystem > &p) const
Transformation operation for Position Vector in any coordinate system.
Definition: Transform3D.h:645
static const double x1[5]
Transform3D(const Quaternion &r)
Definition: Transform3D.h:176
EulerAngles class describing rotation as three angles (Euler Angles).
Definition: EulerAngles.h:43
int type
Definition: TGX11.cxx:120
Transform3D(const RotationY &r)
Definition: Transform3D.h:188
Transform3D(const AxisAngle &r)
Definition: Transform3D.h:170
Transform3D(const RotationX &r)
Definition: Transform3D.h:185
PositionVector3D< CoordSystem > operator*(const PositionVector3D< CoordSystem > &v) const
Transformation operation for Position Vector in any coordinate system.
Definition: Transform3D.h:653
Transform3D(const Translation3D< T > &t)
Construct from a translation only, represented by a Translation3D class and with an identity rotation...
Definition: Transform3D.h:214
Namespace for new Math classes and functions.
Plane3D< TYPE > operator()(const Plane3D< TYPE > &plane) const
Transformation on a 3D plane.
Definition: Transform3D.h:764
Transform3D< T > & operator*=(const Transform3D< T > &t)
multiply (combine) with another transformation in place
Definition: Transform3D.h:1018
Scalar Y() const
Cartesian Y, converting if necessary from internal coordinate system.
PositionVector3D< CoordSystem, Tag > & SetXYZ(Scalar a, Scalar b, Scalar c)
set the values of the vector from the cartesian components (x,y,z) (if the vector is held in polar or...
Transform3D< T > Inverse() const
Return the inverse of the transformation.
Definition: Transform3D.h:881
PositionVector3D< Cartesian3D< T >, DefaultCoordinateSystemTag > Point
Definition: Transform3D.h:84
DisplacementVector3D< Cartesian3D< T >, DefaultCoordinateSystemTag > Vector
Definition: Transform3D.h:83
Transform3D()
Default constructor (identy rotation) + zero translation.
Definition: Transform3D.h:97
Transform3D(const EulerAngles &r)
Definition: Transform3D.h:173
Transform3D(const Rotation3D &r, const Vector &v)
Construct from a rotation and then a translation described by a Vector.
Definition: Transform3D.h:115
DisplacementVector3D< CoordSystem > ApplyInverse(const DisplacementVector3D< CoordSystem > &p) const
Directly apply the inverse affine transformation on an arbitrary coordinate-system vector...
Definition: Transform3D.h:717
Class describing a geometrical plane in 3 dimensions.
Definition: Plane3D.h:51
SVector< T, D > Unit(const SVector< T, D > &rhs)
Unit.
Definition: Functions.h:381
AnyRotation Rotation() const
Get the rotation representing the 3D transformation.
Definition: Transform3D.h:589
DefaultCoordinateSystemTag Default tag for identifying any coordinate system.
Basic 3D Transformation class describing a rotation and then a translation The internal data are a 3D...
Definition: Transform3D.h:78
float * q
Definition: THbookFile.cxx:87
Transform3D(T xx, T xy, T xz, T dx, T yx, T yy, T yz, T dy, T zx, T zy, T zz, T dz)
Raw constructor from 12 Scalar components.
Definition: Transform3D.h:443
const Int_t n
Definition: legend1.C:16
::ROOT::Math::DisplacementVector3D< Cartesian3D< Scalar > > Vect() const
get the spatial components of the Vector in a DisplacementVector based on Cartesian Coordinates ...
void Transform(const DisplacementVector3D< CoordSystem, Tag1 > &v1, DisplacementVector3D< CoordSystem, Tag2 > &v2) const
Transformation operation for Displacement Vector of different coordinate systems. ...
Definition: Transform3D.h:737
DisplacementVector3D< CoordSystem > operator*(const DisplacementVector3D< CoordSystem > &v) const
Transformation operation for Displacement Vector in any coordinate system.
Definition: Transform3D.h:669
Transform3D(const RotationZYX &r)
Definition: Transform3D.h:179
Rotation3D Rotation() const
Get the 3D rotation representing the 3D transformation.
Definition: Transform3D.h:579
PositionVector3D< CoordSystem > ApplyInverse(const PositionVector3D< CoordSystem > &p) const
Directly apply the inverse affine transformation on an arbitrary coordinate-system point...
Definition: Transform3D.h:706
Scalar Dot(const DisplacementVector3D< OtherCoords, Tag > &v) const
Return the scalar (dot) product of two displacement vectors.