Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:26:33

0001 // @(#)root/mathcore:$Name: V02-02-29 $:$Id: Transform3DPJ.h,v 1.1 2007/05/27 20:12:30 pjanot Exp $
0002 // Authors: W. Brown, M. Fischler, L. Moneta    2005
0003 
0004 /**********************************************************************
0005  *                                                                    *
0006  * Copyright (c) 2005 , LCG ROOT MathLib Team                         *
0007  *                                                                    *
0008  *                                                                    *
0009  **********************************************************************/
0010 
0011 // Header file for class Transform3DPJ
0012 //
0013 // Created by: Lorenzo Moneta  October 21 2005
0014 //
0015 //
0016 #ifndef ROOT_Math_GenVector_Transform3DPJ
0017 #define ROOT_Math_GenVector_Transform3DPJ 1
0018 
0019 #include "Math/GenVector/Cartesian3D.h"
0020 #include "Math/GenVector/DisplacementVector3D.h"
0021 #include "Math/GenVector/PositionVector3D.h"
0022 #include "Math/GenVector/LorentzVector.h"
0023 #include "Math/GenVector/Rotation3D.h"
0024 #include "Math/GenVector/AxisAnglefwd.h"
0025 #include "Math/GenVector/EulerAnglesfwd.h"
0026 #include "Math/GenVector/Quaternionfwd.h"
0027 #include "Math/GenVector/RotationXfwd.h"
0028 #include "Math/GenVector/RotationYfwd.h"
0029 #include "Math/GenVector/RotationZfwd.h"
0030 #include "Math/GenVector/Plane3D.h"
0031 
0032 #include <iostream>
0033 
0034 //#include "Math/Vector3Dfwd.h"
0035 
0036 namespace ROOT {
0037 
0038   namespace Math {
0039 
0040     using ROOT::Math::Plane3D;
0041 
0042     /** 
0043      Basic 3D Transformation class describing  a rotation and then a translation
0044      The internal data are a rotation data and a 3D vector data and they can be represented 
0045      like a 3x4 matrix
0046      The class has a template parameter the coordinate system tag of the reference system 
0047      to which the transformatioon will be applied. For example for transforming from 
0048      global to local coordinate systems, the transfrom3D has to be instantiated with the 
0049      coordinate of the traget system
0050 
0051      @ingroup GenVector
0052 
0053   */
0054 
0055     class Transform3DPJ {
0056     public:
0057       typedef DisplacementVector3D<Cartesian3D<double>, DefaultCoordinateSystemTag> Vector;
0058       typedef PositionVector3D<Cartesian3D<double>, DefaultCoordinateSystemTag> Point;
0059 
0060       enum ETransform3DMatrixIndex {
0061         kXX = 0,
0062         kXY = 1,
0063         kXZ = 2,
0064         kDX = 3,
0065         kYX = 4,
0066         kYY = 5,
0067         kYZ = 6,
0068         kDY = 7,
0069         kZX = 8,
0070         kZY = 9,
0071         kZZ = 10,
0072         kDZ = 11
0073       };
0074 
0075       /** 
0076     Default constructor (identy rotation) + zero translation
0077     */
0078       Transform3DPJ() { SetIdentity(); }
0079 
0080       /**
0081        Construct given a pair of pointers or iterators defining the
0082        beginning and end of an array of 12 Scalars
0083     */
0084       template <class IT>
0085       Transform3DPJ(IT begin, IT end) {
0086         SetComponents(begin, end);
0087       }
0088 
0089       /**
0090        Construct from a rotation and then a translation described by a Vector 
0091     */
0092       Transform3DPJ(const Rotation3D &r, const Vector &v) { AssignFrom(r, v); }
0093       /**
0094        Construct from a translation and then a rotation (inverse assignment) 
0095     */
0096       Transform3DPJ(const Vector &v, const Rotation3D &r) {
0097         // is equivalent from having first the rotation and then the translation vector rotated
0098         AssignFrom(r, r(v));
0099       }
0100 
0101       /**
0102        Construct from a 3D Rotation only with zero translation
0103     */
0104       explicit Transform3DPJ(const Rotation3D &r) { AssignFrom(r); }
0105       // convenience methods for the other rotations (cannot use templates for conflict with LA)
0106       explicit Transform3DPJ(const AxisAngle &r) { AssignFrom(Rotation3D(r)); }
0107       explicit Transform3DPJ(const EulerAngles &r) { AssignFrom(Rotation3D(r)); }
0108       explicit Transform3DPJ(const Quaternion &r) { AssignFrom(Rotation3D(r)); }
0109       // TO DO: implement direct methods for axial rotations without going through Rotation3D
0110       explicit Transform3DPJ(const RotationX &r) { AssignFrom(Rotation3D(r)); }
0111       explicit Transform3DPJ(const RotationY &r) { AssignFrom(Rotation3D(r)); }
0112       explicit Transform3DPJ(const RotationZ &r) { AssignFrom(Rotation3D(r)); }
0113 
0114       /**
0115        Construct from a translation only, represented by any DisplacementVector3D 
0116        and with an identity rotation
0117     */
0118       template <class CoordSystem, class Tag>
0119       explicit Transform3DPJ(const DisplacementVector3D<CoordSystem, Tag> &v) {
0120         AssignFrom(Vector(v.X(), v.Y(), v.Z()));
0121       }
0122       /**
0123        Construct from a translation only, represented by a Cartesian 3D Vector,  
0124        and with an identity rotation
0125     */
0126       explicit Transform3DPJ(const Vector &v) { AssignFrom(v); }
0127 
0128       //#if !defined(__MAKECINT__) && !defined(G__DICTIONARY)  // this is ambigous with double * , double *
0129       /**
0130        Construct from a rotation (any rotation object)  and then a translation 
0131        (represented by any DisplacementVector)
0132        The requirements on the rotation and vector objects are that they can be transformed in a 
0133        Rotation3D class and in a Vector
0134     */
0135       // to do : change to displacement vector3D
0136       template <class ARotation, class CoordSystem, class Tag>
0137       Transform3DPJ(const ARotation &r, const DisplacementVector3D<CoordSystem, Tag> &v) {
0138         AssignFrom(Rotation3D(r), Vector(v.X(), v.Y(), v.Z()));
0139       }
0140       /**
0141        Construct from a translation (using any type of DisplacementVector ) 
0142        and then a rotation (any rotation object). 
0143        Requirement on the rotation and vector objects are that they can be transformed in a 
0144        Rotation3D class and in a Vector 
0145     */
0146       template <class ARotation, class CoordSystem, class Tag>
0147       Transform3DPJ(const DisplacementVector3D<CoordSystem, Tag> &v, const ARotation &r) {
0148         // is equivalent from having first the rotation and then the translation vector rotated
0149         Rotation3D r3d(r);
0150         AssignFrom(r3d, r3d(Vector(v.X(), v.Y(), v.Z())));
0151       }
0152 
0153       //#endif
0154 
0155       /**
0156        Construct transformation from one coordinate system defined by three 
0157        points (origin + two axis) to 
0158        a new coordinate system defined by other three points (origin + axis) 
0159        @param fr0  point defining origin of original reference system 
0160        @param fr1  point defining first axis of original reference system 
0161        @param fr2  point defining second axis of original reference system 
0162        @param to0  point defining origin of transformed reference system 
0163        @param to1  point defining first axis transformed reference system 
0164        @param to2  point defining second axis transformed reference system 
0165 
0166      */
0167       Transform3DPJ(
0168           const Point &fr0, const Point &fr1, const Point &fr2, const Point &to0, const Point &to1, const Point &to2);
0169 
0170       // use compiler generated copy ctor, copy assignmet and dtor
0171 
0172       /**
0173        Construct from a linear algebra matrix of size at least 3x4,
0174        which must support operator()(i,j) to obtain elements (0,0) thru (2,3).
0175        The 3x3 sub-block is assumed to be the rotation part and the translations vector 
0176        are described by the 4-th column
0177     */
0178       template <class ForeignMatrix>
0179       explicit Transform3DPJ(const ForeignMatrix &m) {
0180         SetComponents(m);
0181       }
0182 
0183       /**
0184        Raw constructor from 12 Scalar components
0185     */
0186       Transform3DPJ(double xx,
0187                     double xy,
0188                     double xz,
0189                     double dx,
0190                     double yx,
0191                     double yy,
0192                     double yz,
0193                     double dy,
0194                     double zx,
0195                     double zy,
0196                     double zz,
0197                     double dz) {
0198         SetComponents(xx, xy, xz, dx, yx, yy, yz, dy, zx, zy, zz, dz);
0199       }
0200 
0201       /**
0202        Construct from a linear algebra matrix of size at least 3x4,
0203        which must support operator()(i,j) to obtain elements (0,0) thru (2,3).
0204        The 3x3 sub-block is assumed to be the rotation part and the translations vector 
0205        are described by the 4-th column
0206     */
0207       template <class ForeignMatrix>
0208       Transform3DPJ &operator=(const ForeignMatrix &m) {
0209         SetComponents(m);
0210         return *this;
0211       }
0212 
0213       // ======== Components ==============
0214 
0215       /**
0216        Set the 12 matrix components given an iterator to the start of
0217        the desired data, and another to the end (12 past start).
0218     */
0219       template <class IT>
0220       void SetComponents(IT begin, IT end) {
0221         for (int i = 0; i < 12; ++i) {
0222           fM[i] = *begin;
0223           ++begin;
0224         }
0225         assert(end == begin);
0226       }
0227 
0228       /**
0229        Get the 12 matrix components into data specified by an iterator begin
0230        and another to the end of the desired data (12 past start).
0231     */
0232       template <class IT>
0233       void GetComponents(IT begin, IT end) const {
0234         for (int i = 0; i < 12; ++i) {
0235           *begin = fM[i];
0236           ++begin;
0237         }
0238         assert(end == begin);
0239       }
0240 
0241       /**
0242        Get the 12 matrix components into data specified by an iterator begin
0243     */
0244       template <class IT>
0245       void GetComponents(IT begin) const {
0246         std::copy(fM, fM + 12, begin);
0247       }
0248 
0249       /**
0250        Set components from a linear algebra matrix of size at least 3x4,
0251        which must support operator()(i,j) to obtain elements (0,0) thru (2,3).
0252        The 3x3 sub-block is assumed to be the rotation part and the translations vector 
0253        are described by the 4-th column
0254     */
0255       template <class ForeignMatrix>
0256       void SetTransformMatrix(const ForeignMatrix &m) {
0257         fM[kXX] = m(0, 0);
0258         fM[kXY] = m(0, 1);
0259         fM[kXZ] = m(0, 2);
0260         fM[kDX] = m(0, 3);
0261         fM[kYX] = m(1, 0);
0262         fM[kYY] = m(1, 1);
0263         fM[kYZ] = m(1, 2);
0264         fM[kDY] = m(1, 3);
0265         fM[kZX] = m(2, 0);
0266         fM[kZY] = m(2, 1);
0267         fM[kZZ] = m(2, 2);
0268         fM[kDZ] = m(2, 3);
0269       }
0270 
0271       /**
0272        Get components into a linear algebra matrix of size at least 3x4,
0273        which must support operator()(i,j) for write access to elements
0274        (0,0) thru (2,3).
0275     */
0276       template <class ForeignMatrix>
0277       void GetTransformMatrix(ForeignMatrix &m) const {
0278         m(0, 0) = fM[kXX];
0279         m(0, 1) = fM[kXY];
0280         m(0, 2) = fM[kXZ];
0281         m(0, 3) = fM[kDX];
0282         m(1, 0) = fM[kYX];
0283         m(1, 1) = fM[kYY];
0284         m(1, 2) = fM[kYZ];
0285         m(1, 3) = fM[kDY];
0286         m(2, 0) = fM[kZX];
0287         m(2, 1) = fM[kZY];
0288         m(2, 2) = fM[kZZ];
0289         m(2, 3) = fM[kDZ];
0290       }
0291 
0292       /**
0293        Set the components from 12 scalars 
0294     */
0295       void SetComponents(double xx,
0296                          double xy,
0297                          double xz,
0298                          double dx,
0299                          double yx,
0300                          double yy,
0301                          double yz,
0302                          double dy,
0303                          double zx,
0304                          double zy,
0305                          double zz,
0306                          double dz) {
0307         fM[kXX] = xx;
0308         fM[kXY] = xy;
0309         fM[kXZ] = xz;
0310         fM[kDX] = dx;
0311         fM[kYX] = yx;
0312         fM[kYY] = yy;
0313         fM[kYZ] = yz;
0314         fM[kDY] = dy;
0315         fM[kZX] = zx;
0316         fM[kZY] = zy;
0317         fM[kZZ] = zz;
0318         fM[kDZ] = dz;
0319       }
0320 
0321       /**
0322        Get the nine components into 12 scalars
0323     */
0324       void GetComponents(double &xx,
0325                          double &xy,
0326                          double &xz,
0327                          double &dx,
0328                          double &yx,
0329                          double &yy,
0330                          double &yz,
0331                          double &dy,
0332                          double &zx,
0333                          double &zy,
0334                          double &zz,
0335                          double &dz) const {
0336         xx = fM[kXX];
0337         xy = fM[kXY];
0338         xz = fM[kXZ];
0339         dx = fM[kDX];
0340         yx = fM[kYX];
0341         yy = fM[kYY];
0342         yz = fM[kYZ];
0343         dy = fM[kDY];
0344         zx = fM[kZX];
0345         zy = fM[kZY];
0346         zz = fM[kZZ];
0347         dz = fM[kDZ];
0348       }
0349 
0350       /**
0351        Get the rotation and translation vector representing the 3D transformation
0352     */
0353       void GetDecomposition(Rotation3D &r, Vector &v) const;
0354 
0355       // operations on points and vectors
0356 
0357       /**
0358        Transformation operation for Position Vector in Cartesian coordinate 
0359     */
0360       Point operator()(const Point &p) const;
0361 
0362       /**
0363        Transformation operation for Displacement Vectors in Cartesian coordinate 
0364        For the Displacement Vectors only the rotation applies - no translations
0365     */
0366       Vector operator()(const Vector &v) const;
0367 
0368       /**
0369        Transformation operation for Position Vector in any coordinate system 
0370     */
0371       template <class CoordSystem>
0372       PositionVector3D<CoordSystem> operator()(const PositionVector3D<CoordSystem> &p) const {
0373         Point xyzNew = operator()(Point(p));
0374         return PositionVector3D<CoordSystem>(xyzNew);
0375       }
0376 
0377       /**
0378        Transformation operation for Displacement Vector in any coordinate system 
0379     */
0380       template <class CoordSystem>
0381       DisplacementVector3D<CoordSystem> operator()(const DisplacementVector3D<CoordSystem> &v) const {
0382         Vector xyzNew = operator()(Vector(v));
0383         return DisplacementVector3D<CoordSystem>(xyzNew);
0384       }
0385 
0386       /**
0387        Transformation operation for points between different coordinate system tags 
0388     */
0389       template <class CoordSystem, class Tag1, class Tag2>
0390       void Transform(const PositionVector3D<CoordSystem, Tag1> &p1, PositionVector3D<CoordSystem, Tag2> &p2) const {
0391         Point xyzNew = operator()(Point(p1.X(), p1.Y(), p1.Z()));
0392         p2.SetXYZ(xyzNew.X(), xyzNew.Y(), xyzNew.Z());
0393       }
0394 
0395       /**
0396        Transformation operation for Displacement Vector of different coordinate systems 
0397     */
0398       template <class CoordSystem, class Tag1, class Tag2>
0399       void Transform(const DisplacementVector3D<CoordSystem, Tag1> &v1,
0400                      DisplacementVector3D<CoordSystem, Tag2> &v2) const {
0401         Vector xyzNew = operator()(Vector(v1.X(), v1.Y(), v1.Z()));
0402         v2.SetXYZ(xyzNew.X(), xyzNew.Y(), xyzNew.Z());
0403       }
0404 
0405       /**
0406        Transformation operation for a Lorentz Vector in any  coordinate system 
0407     */
0408       template <class CoordSystem>
0409       LorentzVector<CoordSystem> operator()(const LorentzVector<CoordSystem> &q) const {
0410         Vector xyzNew = operator()(Vector(q.Vect()));
0411         return LorentzVector<CoordSystem>(xyzNew.X(), xyzNew.Y(), xyzNew.Z(), q.E());
0412       }
0413 
0414       /**
0415        Transformation on a 3D plane
0416     */
0417       Plane3D operator()(const Plane3D &plane) const;
0418 
0419       // skip transformation for arbitrary vectors - not really defined if point or displacement vectors
0420 
0421       // same but with operator *
0422       /**
0423        Transformation operation for Vectors. Apply same rules as operator() 
0424        depending on type of vector. 
0425        Will work only for DisplacementVector3D, PositionVector3D and LorentzVector
0426     */
0427       template <class AVector>
0428       AVector operator*(const AVector &v) const {
0429         return operator()(v);
0430       }
0431 
0432       /**
0433        multiply (combine) with another transformation in place
0434      */
0435       Transform3DPJ &operator*=(const Transform3DPJ &t);
0436 
0437       /**
0438        multiply (combine) two transformations
0439      */
0440       Transform3DPJ operator*(const Transform3DPJ &t) const {
0441         Transform3DPJ tmp(*this);
0442         tmp *= t;
0443         return tmp;
0444       }
0445 
0446       /** 
0447     Invert the transformation in place
0448     */
0449       void Invert();
0450 
0451       /**
0452        Return the inverse of the transformation.
0453     */
0454       Transform3DPJ Inverse() const {
0455         Transform3DPJ t(*this);
0456         t.Invert();
0457         return t;
0458       }
0459 
0460       /**
0461        Equality/inequality operators
0462     */
0463       bool operator==(const Transform3DPJ &rhs) const {
0464         if (fM[0] != rhs.fM[0])
0465           return false;
0466         if (fM[1] != rhs.fM[1])
0467           return false;
0468         if (fM[2] != rhs.fM[2])
0469           return false;
0470         if (fM[3] != rhs.fM[3])
0471           return false;
0472         if (fM[4] != rhs.fM[4])
0473           return false;
0474         if (fM[5] != rhs.fM[5])
0475           return false;
0476         if (fM[6] != rhs.fM[6])
0477           return false;
0478         if (fM[7] != rhs.fM[7])
0479           return false;
0480         if (fM[8] != rhs.fM[8])
0481           return false;
0482         if (fM[9] != rhs.fM[9])
0483           return false;
0484         if (fM[10] != rhs.fM[10])
0485           return false;
0486         if (fM[11] != rhs.fM[11])
0487           return false;
0488         return true;
0489       }
0490 
0491       bool operator!=(const Transform3DPJ &rhs) const { return !operator==(rhs); }
0492 
0493     protected:
0494       /**
0495        make transformation from first a rotation then a translation
0496      */
0497       void AssignFrom(const Rotation3D &r, const Vector &v);
0498 
0499       /**
0500        make transformation from only rotations (zero translation)
0501      */
0502       void AssignFrom(const Rotation3D &r);
0503 
0504       /**
0505        make transformation from only translation (identity rotations)
0506      */
0507       void AssignFrom(const Vector &v);
0508 
0509       /**
0510        Set identity transformation (identity rotation , zero translation)
0511      */
0512       void SetIdentity();
0513 
0514     private:
0515       double fM[12];
0516     };
0517 
0518     // global functions
0519 
0520     // TODO - I/O should be put in the manipulator form
0521 
0522     std::ostream &operator<<(std::ostream &os, const Transform3DPJ &t);
0523 
0524     // need a function Transform = Translation * Rotation ???
0525 
0526   }  // end namespace Math
0527 
0528 }  // end namespace ROOT
0529 
0530 #endif /* MATHCORE_BASIC3DTRANSFORMATION */