Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-01-18 09:54:36

0001 // -*- C++ -*-
0002 // $Id: Transform3D.h,v 1.5 2010/06/16 16:21:27 garren Exp $
0003 // ---------------------------------------------------------------------------
0004 //
0005 // This file is a part of the CLHEP - a Class Library for High Energy Physics.
0006 //
0007 // Hep geometrical 3D Transformation class
0008 //
0009 // Author: Evgeni Chernyaev <Evgueni.Tcherniaev@cern.ch>
0010 //
0011 //          ******************************************
0012 //          *                                        *
0013 //          *               Transform                *
0014 //          *               /  / \  \                *
0015 //          *       --------  /   \  --------        *
0016 //          *      /         /     \         \       *
0017 //          *   Rotate Translate  Reflect   Scale    *
0018 //          *    / | \    / | \    / | \    / | \    *
0019 //          *   X  Y  Z  X  Y  Z  X  Y  Z  X  Y  Z   *
0020 //          *                                        *
0021 //          ******************************************
0022 //
0023 // Identity transformation:
0024 //   Transform3D::Identity   - global identity transformation;
0025 //   any constructor without parameters, e.g. Transform3D();
0026 //   m.setIdentity()            - set "m" to identity;
0027 //
0028 // General transformations:
0029 //   Transform3D(m,v)         - transformation given by Rotation "m"
0030 //                              and CLHEP::Hep3Vector "v";
0031 //   Transform3D(a0,a1,a2, b0,b1,b2) - transformation given by initial
0032 //                               and transformed positions of three points;
0033 // Rotations:
0034 //   Rotate3D(m)              - rotation given by CLHEP::HepRotation "m";
0035 //   Rotate3D(ang,v)          - rotation through the angle "ang" around
0036 //                              vector "v";
0037 //   Rotate3D(ang,p1,p2)      - rotation through the angle "ang"
0038 //                              counterclockwise around the axis given by
0039 //                              two points p1->p2;
0040 //   Rotate3D(a1,a2, b1,b2)   - rotation around the origin defined by initial
0041 //                              and transformed positions of two points;
0042 //   RotateX3D(ang)           - rotation around X-axis;
0043 //   RotateY3D(ang)           - rotation around Y-axis;
0044 //   RotateZ3D(ang)           - rotation around Z-axis;
0045 //
0046 // Translations:
0047 //   Translate3D(v)           - translation given by CLHEP::Hep3Vector "v";
0048 //   Translate3D(dx,dy,dz)    - translation on vector (dx,dy,dz);
0049 //   TraslateX3D(dx)          - translation along X-axis;
0050 //   TraslateY3D(dy)          - translation along Y-axis;
0051 //   TraslateZ3D(dz)          - translation along Z-axis;
0052 //
0053 // Reflections:
0054 //   Reflect3D(a,b,c,d)       - reflection in the plane a*x+b*y+c*z+d=0;
0055 //   Reflect3D(normal,p)      - reflection in the plane going through "p"
0056 //                              and whose normal is equal to "normal";
0057 //   ReflectX3D(a)            - reflect X in the plane x=a (default a=0);
0058 //   ReflectY3D(a)            - reflect Y in the plane y=a (default a=0);
0059 //   ReflectZ3D(a)            - reflect Z in the plane z=a (default a=0);
0060 //
0061 // Scalings:
0062 //   Scale3D(sx,sy,sz)        - general scaling with factors "sx","sy","sz"
0063 //                                 along X, Y and Z;
0064 //   Scale3D(s)               - scaling with constant factor "s" along all
0065 //                                 directions;
0066 //   ScaleX3D(sx)             - scale X;
0067 //   ScaleY3D(sy)             - scale Y;
0068 //   ScaleZ3D(sz)             - scale Z;
0069 //
0070 // Inverse transformation:
0071 //   m.inverse() or           - returns inverse transformation;
0072 //
0073 // Compound transformation:
0074 //   m3 = m2 * m1             - it is relatively slow in comparison with
0075 //                              transformation of a vector. Use parenthesis
0076 //                              to avoid this operation (see example below);
0077 // Transformation of point:
0078 //   p2 = m * p1
0079 //
0080 // Transformation of vector:
0081 //   v2 = m * v1
0082 //
0083 // Transformation of normal:
0084 //   n2 = m * n1
0085 //
0086 // The following table explains how different transformations affect
0087 // point, vector and normal. "+" means affect, "-" means do not affect,
0088 // "*" meas affect but in different way than "+"
0089 //
0090 //                     Point  Vector  Normal
0091 //      -------------+-------+-------+-------
0092 //       Rotation    !   +   !   +   !   +
0093 //       Translation !   +   !   -   !   -
0094 //       Reflection  !   +   !   +   !   *
0095 //       Scaling     !   +   !   +   !   *
0096 //      -------------+-------+-------+-------
0097 //
0098 // Example of the usage:
0099 //
0100 //   Transform3D m1, m2, m3;
0101 //   HepVector3D    v2, v1(0,0,0);
0102 //
0103 //   m1 = Rotate3D(angle, Vector3D(1,1,1));
0104 //   m2 = Translate3D(dx,dy,dz);
0105 //   m3 = m1.inverse();
0106 //
0107 //   v2 = m3*(m2*(m1*v1));
0108 //
0109 // History:
0110 // 24.09.96 E.Chernyaev - initial version
0111 //
0112 // 26.02.97 E.Chernyaev
0113 // - added global Identity by request of John Allison
0114 //   (to avoid problems with compilation on HP)
0115 // - added getRotation and getTranslation
0116 //
0117 // 29.01.01 E.Chernyaev - added subscripting
0118 // 11.06.01 E.Chernyaev - added getDecomposition
0119 
0120 #ifndef HEP_TRANSFROM3D_H
0121 #define HEP_TRANSFROM3D_H
0122 
0123 #include "CLHEP/Geometry/defs.h"
0124 #include "CLHEP/Vector/ThreeVector.h"
0125 
0126 namespace HepGeom {
0127 
0128   template<class T> class Point3D;
0129   template<class T> class Vector3D;
0130   template<class T> class Normal3D;
0131 
0132   class Translate3D;
0133   class Rotate3D;
0134   class Scale3D;
0135 
0136   /**
0137    * Class for transformation of 3D geometrical objects.
0138    * It allows different translations, rotations, scalings and reflections.
0139    * Several specialized classes are derived from it:
0140    *
0141    * TranslateX3D, TranslateY3D, TranslateZ3D, Translate3D,<br>
0142    * RotateX3D,    RotateY3D,    RotateZ3D,    Rotate3D,   <br>
0143    * ScaleX3D,     ScaleY3D,     ScaleZ3D,     Scale3D,    <br>
0144    * ReflectX3D,   ReflectY3D,   ReflectZ3D,   Reflect3D.
0145    *
0146    * The idea behind these classes is to provide some additional constructors
0147    * for Transform3D, they normally should not be used as separate classes.
0148    *
0149    * Example:
0150    * @code
0151    *   HepGeom::Transform3D m;
0152    *   m = HepGeom::TranslateX3D(10.*cm);
0153    * @endcode
0154    *
0155    * Remark:
0156    * For the reason that the operator* is left associative, the notation
0157    * @code
0158    *   v2 = m3*(m2*(m1*v1));
0159    * @endcode
0160    * is much more effective then the notation
0161    * @code
0162    *   v2 = m3*m2*m1*v1;
0163    * @endcode
0164    * In the first case three operations Transform3D*Vector3D are executed,
0165    * in the second case two operations Transform3D*Transform3D and one
0166    * Transform3D*Vector3D are performed. Transform3D*Transform3D is
0167    * roughly three times slower than Transform3D*Vector3D.
0168    *
0169    * @author <Evgueni.Tcherniaev@cern.ch>
0170    * @ingroup geometry
0171    */
0172   class Transform3D {
0173   protected:
0174     double xx_, xy_, xz_, dx_,     // 4x3  Transformation Matrix
0175            yx_, yy_, yz_, dy_,
0176            zx_, zy_, zz_, dz_;
0177 
0178     // Protected constructor
0179     Transform3D(double XX, double XY, double XZ, double DX,
0180         double YX, double YY, double YZ, double DY,
0181         double ZX, double ZY, double ZZ, double DZ)
0182       : xx_(XX), xy_(XY), xz_(XZ), dx_(DX),
0183     yx_(YX), yy_(YY), yz_(YZ), dy_(DY),
0184     zx_(ZX), zy_(ZY), zz_(ZZ), dz_(DZ) {}
0185 
0186     // Set transformation matrix
0187     void setTransform(double XX, double XY, double XZ, double DX,
0188               double YX, double YY, double YZ, double DY,
0189               double ZX, double ZY, double ZZ, double DZ) {
0190       xx_ = XX; xy_ = XY; xz_ = XZ; dx_ = DX;
0191       yx_ = YX; yy_ = YY; yz_ = YZ; dy_ = DY;
0192       zx_ = ZX; zy_ = ZY; zz_ = ZZ; dz_ = DZ;
0193     }
0194 
0195   public:
0196     /**
0197      * Global identity transformation. */
0198     static const Transform3D Identity;
0199 
0200     // Helper class for implemention of C-style subscripting r[i][j]
0201     class Transform3D_row {
0202     public:
0203       inline Transform3D_row(const Transform3D &, int);
0204       inline double operator [] (int) const;
0205     private:
0206       const Transform3D & rr;
0207       int ii;
0208     };
0209 
0210     /**
0211      * Default constructor - sets the Identity transformation. */
0212     Transform3D()
0213       : xx_(1), xy_(0), xz_(0), dx_(0),
0214     yx_(0), yy_(1), yz_(0), dy_(0),
0215     zx_(0), zy_(0), zz_(1), dz_(0) {}
0216 
0217     /**
0218      * Constructor: rotation and then translation. */
0219     inline Transform3D(const CLHEP::HepRotation & mt, const CLHEP::Hep3Vector & v);
0220 
0221     /**
0222      * Constructor: transformation of basis (assumed - no reflection). */
0223     Transform3D(const Point3D<double> & fr0,
0224         const Point3D<double> & fr1,
0225         const Point3D<double> & fr2,
0226         const Point3D<double> & to0,
0227         const Point3D<double> & to1,
0228         const Point3D<double> & to2);
0229 
0230     /**
0231      * Copy constructor. */
0232     Transform3D(const Transform3D & mt) = default;
0233 
0234     /**
0235      * Move constructor. */
0236     Transform3D(Transform3D && mt) = default;
0237 
0238     /**
0239      * Destructor. */
0240     ~Transform3D() = default;
0241 
0242     /**
0243      * Assignment. */
0244     Transform3D & operator=(const Transform3D & mt) = default;
0245 
0246     /**
0247      * Move assignment. */
0248     Transform3D & operator=(Transform3D && mt) = default;
0249 
0250     /**
0251      * Returns object of the helper class for C-style subscripting r[i][j] */
0252     inline const Transform3D_row operator [] (int) const;
0253 
0254     /** Fortran-style subscripting: returns (i,j) element of the matrix. */
0255     double operator () (int, int) const;
0256 
0257     /**
0258      * Gets xx-element of the transformation matrix. */
0259     double xx() const { return xx_; }
0260     /**
0261      * Gets xy-element of the transformation matrix. */
0262     double xy() const { return xy_; }
0263     /**
0264      * Gets xz-element of the transformation matrix. */
0265     double xz() const { return xz_; }
0266     /**
0267      * Gets yx-element of the transformation matrix. */
0268     double yx() const { return yx_; }
0269     /**
0270      * Gets yy-element of the transformation matrix. */
0271     double yy() const { return yy_; }
0272     /**
0273      * Gets yz-element of the transformation matrix. */
0274     double yz() const { return yz_; }
0275     /**
0276      * Gets zx-element of the transformation matrix. */
0277     double zx() const { return zx_; }
0278     /**
0279      * Gets zy-element of the transformation matrix. */
0280     double zy() const { return zy_; }
0281     /**
0282      * Gets zz-element of the transformation matrix. */
0283     double zz() const { return zz_; }
0284     /**
0285      * Gets dx-element of the transformation matrix. */
0286     double dx() const { return dx_; }
0287     /**
0288      * Gets dy-element of the transformation matrix. */
0289     double dy() const { return dy_; }
0290     /**
0291      * Gets dz-element of the transformation matrix. */
0292     double dz() const { return dz_; }
0293 
0294     /**
0295      * Sets the Identity transformation. */
0296     void setIdentity() {
0297       xy_= xz_= dx_= yx_= yz_= dy_= zx_= zy_= dz_= 0; xx_= yy_= zz_= 1;
0298     }
0299 
0300     /**
0301      * Returns the inverse transformation. */
0302     Transform3D inverse() const;
0303 
0304     /**
0305      * Transformation by another Transform3D. */
0306     Transform3D operator*(const Transform3D & b) const;
0307 
0308     /**
0309      * Decomposition of general transformation.
0310      * This function gets decomposition of the transformation
0311      * in three consequentive specific transformations: Scale3D,
0312      * then Rotate3D, then Translate3, i.e.
0313      * @code
0314      *   Transform3D = Translate3D * Rotate3D * Scale3D
0315      * @endcode
0316      *
0317      * @param scale       output: scaling transformation;
0318      *                    if there was a reflection, then scale factor for
0319      *                    z-component (scale(2,2)) will be negative.
0320      * @param rotation    output: rotation transformaion.
0321      * @param translation output: translation transformaion.
0322      */
0323     void getDecomposition(Scale3D & scale,
0324               Rotate3D & rotation,
0325               Translate3D & translation) const;
0326 
0327     /**
0328      * Returns true if the difference between corresponding
0329      * matrix elements is less than the tolerance.
0330      */
0331     bool isNear(const Transform3D & t, double tolerance = 2.2E-14 ) const;
0332 
0333     /**
0334      * Extracts the rotation matrix.
0335      * This functions is obsolete - use getDecomposition() instead.
0336      */
0337     inline CLHEP::HepRotation getRotation() const;
0338 
0339     /**
0340      * Extracts the translation vector.
0341      * This functions is obsolete - use getDecomposition() instead.
0342      */
0343     inline CLHEP::Hep3Vector getTranslation() const;
0344 
0345     /**
0346      * Test for equality. */
0347     bool operator == (const Transform3D & transform) const;
0348 
0349     /**
0350      * Test for inequality. */
0351     bool operator != (const Transform3D & transform) const {
0352       return ! operator==(transform);
0353     }
0354   };
0355 
0356   //   R O T A T I O N S
0357 
0358   /**
0359    * Constructs a rotation transformation.
0360    * This class provides additional constructors for Transform3D
0361    * and should not be used as a separate class.
0362    *
0363    * Example of use:
0364    * @code
0365    *   Transform3D m;
0366    *   m = Rotate3D(30.*deg, HepVector3D(1.,1.,1.));
0367    * @endcode
0368    *
0369    * @author <Evgueni.Tcherniaev@cern.ch>
0370    * @ingroup geometry
0371    */
0372   class Rotate3D : public Transform3D {
0373   public:
0374     /**
0375      * Default constructor: sets the Identity transformation. */
0376     Rotate3D() : Transform3D() {}
0377 
0378     /**
0379      * Constructor from CLHEP::HepRotation. */
0380     inline Rotate3D(const CLHEP::HepRotation &mt);
0381 
0382     /**
0383      * Constructor from angle and axis given by two points.
0384      * @param a  angle of rotation
0385      * @param p1 begin point of the axis
0386      * @param p2 end point of the axis
0387      */
0388     Rotate3D(double a,
0389          const Point3D<double> & p1,
0390          const Point3D<double> & p2);
0391 
0392     /**
0393      * Constructor from angle and axis.
0394      * @param a angle of rotation
0395      * @param v axis of rotation
0396      */
0397     inline Rotate3D(double a, const Vector3D<double> & v);
0398 
0399     /**
0400      * Constructor for rotation given by original and rotated position of
0401      * two points. It is assumed that there is no reflection.
0402      * @param fr1 original position of 1st point
0403      * @param fr2 original position of 2nd point
0404      * @param to1 rotated position of 1st point
0405      * @param to2 rotated position of 2nd point
0406      */
0407     inline Rotate3D(const Point3D<double> & fr1,
0408             const Point3D<double> & fr2,
0409             const Point3D<double> & to1,
0410             const Point3D<double> & to2);
0411   };
0412 
0413   /**
0414    * Constructs a rotation around x-axis.
0415    * This class provides additional constructors for Transform3D
0416    * and should not be used as a separate class.
0417    *
0418    * Example of use:
0419    * @code
0420    *   Transform3D m;
0421    *   m = RotateX3D(30.*deg);
0422    * @endcode
0423    *
0424    * @author <Evgueni.Tcherniaev@cern.ch>
0425    * @ingroup geometry
0426    */
0427   class RotateX3D : public Rotate3D {
0428   public:
0429     /**
0430      * Default constructor: sets the Identity transformation. */
0431     RotateX3D() : Rotate3D() {}
0432 
0433     /**
0434      * Constructs a rotation around x-axis by angle a. */
0435     RotateX3D(double a) {
0436       double cosa = std::cos(a), sina = std::sin(a);
0437       setTransform(1,0,0,0,  0,cosa,-sina,0,  0,sina,cosa,0);
0438     }
0439   };
0440 
0441   /**
0442    * Constructs a rotation around y-axis.
0443    * This class provides additional constructors for Transform3D
0444    * and should not be used as a separate class.
0445    *
0446    * Example of use:
0447    * @code
0448    *   Transform3D m;
0449    *   m = RotateY3D(30.*deg);
0450    * @endcode
0451    *
0452    * @author <Evgueni.Tcherniaev@cern.ch>
0453    * @ingroup geometry
0454    */
0455   class RotateY3D : public Rotate3D {
0456   public:
0457     /**
0458      * Default constructor: sets the Identity transformation. */
0459     RotateY3D() : Rotate3D() {}
0460 
0461     /**
0462      * Constructs a rotation around y-axis by angle a. */
0463     RotateY3D(double a) {
0464       double cosa = std::cos(a), sina = std::sin(a);
0465       setTransform(cosa,0,sina,0,  0,1,0,0,  -sina,0,cosa,0);
0466     }
0467   };
0468 
0469   /**
0470    * Constructs a rotation around z-axis.
0471    * This class provides additional constructors for Transform3D
0472    * and should not be used as a separate class.
0473    *
0474    * Example of use:
0475    * @code
0476    *   Transform3D m;
0477    *   m = RotateZ3D(30.*deg);
0478    * @endcode
0479    *
0480    * @author <Evgueni.Tcherniaev@cern.ch>
0481    * @ingroup geometry
0482    */
0483   class RotateZ3D : public Rotate3D {
0484   public:
0485     /**
0486      * Default constructor: sets the Identity transformation. */
0487     RotateZ3D() : Rotate3D() {}
0488 
0489     /**
0490      * Constructs a rotation around z-axis by angle a. */
0491     RotateZ3D(double a) {
0492       double cosa = std::cos(a), sina = std::sin(a);
0493       setTransform(cosa,-sina,0,0,  sina,cosa,0,0,  0,0,1,0);
0494     }
0495   };
0496 
0497   //   T R A N S L A T I O N S
0498 
0499   /**
0500    * Constructs a translation transformation.
0501    * This class provides additional constructors for Transform3D
0502    * and should not be used as a separate class.
0503    *
0504    * Example of use:
0505    * @code
0506    *   Transform3D m;
0507    *   m = Translate3D(10.,20.,30.);
0508    * @endcode
0509    *
0510    * @author <Evgueni.Tcherniaev@cern.ch>
0511    * @ingroup geometry
0512    */
0513   class Translate3D : public Transform3D {
0514   public:
0515     /**
0516      * Default constructor: sets the Identity transformation. */
0517     Translate3D() : Transform3D() {}
0518 
0519     /**
0520      * Constructor from CLHEP::Hep3Vector. */
0521     inline Translate3D(const CLHEP::Hep3Vector &v);
0522 
0523     /**
0524      * Constructor from three numbers. */
0525     Translate3D(double x, double y, double z)
0526       : Transform3D(1,0,0,x, 0,1,0,y, 0,0,1,z) {}
0527   };
0528 
0529   /**
0530    * Constructs a translation along x-axis.
0531    * This class provides additional constructors for Transform3D
0532    * and should not be used as a separate class.
0533    *
0534    * Example of use:
0535    * @code
0536    *   Transform3D m;
0537    *   m = TranslateX3D(10.);
0538    * @endcode
0539    *
0540    * @author <Evgueni.Tcherniaev@cern.ch>
0541    * @ingroup geometry
0542    */
0543   class TranslateX3D : public Translate3D {
0544   public:
0545     /**
0546      * Default constructor: sets the Identity transformation. */
0547     TranslateX3D() : Translate3D() {}
0548 
0549     /**
0550      * Constructor from a number. */
0551     TranslateX3D(double x) : Translate3D(x, 0, 0) {}
0552   };
0553 
0554   /**
0555    * Constructs a translation along y-axis.
0556    * This class provides additional constructors for Transform3D
0557    * and should not be used as a separate class.
0558    *
0559    * Example of use:
0560    * @code
0561    *   Transform3D m;
0562    *   m = TranslateY3D(10.);
0563    * @endcode
0564    *
0565    * @author <Evgueni.Tcherniaev@cern.ch>
0566    * @ingroup geometry
0567    */
0568   class TranslateY3D : public Translate3D {
0569   public:
0570     /**
0571      * Default constructor: sets the Identity transformation. */
0572     TranslateY3D() : Translate3D() {}
0573 
0574     /**
0575      * Constructor from a number. */
0576     TranslateY3D(double y) : Translate3D(0, y, 0) {}
0577   };
0578 
0579   /**
0580    * Constructs a translation along z-axis.
0581    * This class provides additional constructors for Transform3D
0582    * and should not be used as a separate class.
0583    *
0584    * Example of use:
0585    * @code
0586    *   Transform3D m;
0587    *   m = TranslateZ3D(10.);
0588    * @endcode
0589    *
0590    * @author <Evgueni.Tcherniaev@cern.ch>
0591    * @ingroup geometry
0592    */
0593   class TranslateZ3D : public Translate3D {
0594   public:
0595     /**
0596      * Default constructor: sets the Identity transformation. */
0597     TranslateZ3D() : Translate3D() {}
0598 
0599     /**
0600      * Constructor from a number. */
0601     TranslateZ3D(double z) : Translate3D(0, 0, z) {}
0602   };
0603 
0604   //   R E F L E C T I O N S
0605 
0606   /**
0607    * Constructs a reflection transformation.
0608    * This class provides additional constructors for Transform3D
0609    * and should not be used as a separate class.
0610    *
0611    * Example of use:
0612    * @code
0613    *   Transform3D m;
0614    *   m = Reflect3D(1.,1.,1.,0.);
0615    * @endcode
0616    *
0617    * @author <Evgueni.Tcherniaev@cern.ch>
0618    * @ingroup geometry
0619    */
0620   class Reflect3D : public Transform3D {
0621   protected:
0622     Reflect3D(double XX, double XY, double XZ, double DX,
0623               double YX, double YY, double YZ, double DY,
0624               double ZX, double ZY, double ZZ, double DZ)
0625       : Transform3D(XX,XY,XZ,DX, YX,YY,YZ,DY, ZX,ZY,ZZ,DZ) {}
0626 
0627   public:
0628     /**
0629      * Default constructor: sets the Identity transformation. */
0630     Reflect3D() : Transform3D() {}
0631 
0632     /**
0633      * Constructor from four numbers.
0634      * Sets reflection in a plane a*x+b*y+c*z+d=0
0635      */
0636     Reflect3D(double a, double b, double c, double d);
0637 
0638     /**
0639      * Constructor from a plane given by its normal and origin. */
0640     inline Reflect3D(const Normal3D<double> & normal,
0641                      const Point3D<double> & point);
0642   };
0643 
0644   /**
0645    * Constructs reflection in a plane x=const.
0646    * This class provides additional constructors for Transform3D
0647    * and should not be used as a separate class.
0648    *
0649    * Example of use:
0650    * @code
0651    *   Transform3D m;
0652    *   m = ReflectX3D(1.);
0653    * @endcode
0654    *
0655    * @author <Evgueni.Tcherniaev@cern.ch>
0656    * @ingroup geometry
0657    */
0658   class ReflectX3D : public Reflect3D {
0659   public:
0660     /**
0661      * Constructor from a number. */
0662     ReflectX3D(double x=0) : Reflect3D(-1,0,0,x+x, 0,1,0,0, 0,0,1,0) {}
0663   };
0664 
0665   /**
0666    * Constructs reflection in a plane y=const.
0667    * This class provides additional constructors for Transform3D
0668    * and should not be used as a separate class.
0669    *
0670    * Example of use:
0671    * @code
0672    *   Transform3D m;
0673    *   m = ReflectY3D(1.);
0674    * @endcode
0675    *
0676    * @author <Evgueni.Tcherniaev@cern.ch>
0677    * @ingroup geometry
0678    */
0679   class ReflectY3D : public Reflect3D {
0680   public:
0681     /**
0682      * Constructor from a number. */
0683     ReflectY3D(double y=0) : Reflect3D(1,0,0,0, 0,-1,0,y+y, 0,0,1,0) {}
0684   };
0685 
0686   /**
0687    * Constructs reflection in a plane z=const.
0688    * This class provides additional constructors for Transform3D
0689    * and should not be used as a separate class.
0690    *
0691    * Example of use:
0692    * @code
0693    *   Transform3D m;
0694    *   m = ReflectZ3D(1.);
0695    * @endcode
0696    *
0697    * @author <Evgueni.Tcherniaev@cern.ch>
0698    * @ingroup geometry
0699    */
0700   class ReflectZ3D : public Reflect3D {
0701   public:
0702     /**
0703      *  Constructor from a number. */
0704     ReflectZ3D(double z=0) : Reflect3D(1,0,0,0, 0,1,0,0, 0,0,-1,z+z) {}
0705   };
0706 
0707   //   S C A L I N G S
0708 
0709   /**
0710    * Constructs a scaling transformation.
0711    * This class provides additional constructors for Transform3D
0712    * and should not be used as a separate class.
0713    *
0714    * Example of use:
0715    * @code
0716    *   Transform3D m;
0717    *   m = Scale3D(2.);
0718    * @endcode
0719    *
0720    * @author <Evgueni.Tcherniaev@cern.ch>
0721    * @ingroup geometry
0722    */
0723   class Scale3D : public Transform3D {
0724   public:
0725     /**
0726      * Default constructor: sets the Identity transformation. */
0727     Scale3D() : Transform3D() {}
0728 
0729     /**
0730      * Constructor from three numbers - scale factors in different directions.
0731      */
0732     Scale3D(double x, double y, double z)
0733       : Transform3D(x,0,0,0, 0,y,0,0, 0,0,z,0) {}
0734 
0735     /**
0736      * Constructor from a number: sets uniform scaling in all directions. */
0737     Scale3D(double sc)
0738       : Transform3D(sc,0,0,0, 0,sc,0,0, 0,0,sc,0) {}
0739   };
0740 
0741   /**
0742    * Constructs a scaling transformation in x-direction.
0743    * This class provides additional constructors for Transform3D
0744    * and should not be used as a separate class.
0745    *
0746    * Example of use:
0747    * @code
0748    *   Transform3D m;
0749    *   m = ScaleX3D(2.);
0750    * @endcode
0751    *
0752    * @author <Evgueni.Tcherniaev@cern.ch>
0753    * @ingroup geometry
0754    */
0755   class ScaleX3D : public Scale3D {
0756   public:
0757     /**
0758      * Default constructor: sets the Identity transformation. */
0759     ScaleX3D() : Scale3D() {}
0760 
0761     /**
0762      * Constructor from a number (scale factor in x-direction). */
0763     ScaleX3D(double x) : Scale3D(x, 1, 1) {}
0764   };
0765 
0766   /**
0767    * Constructs a scaling transformation in y-direction.
0768    * This class provides additional constructors for Transform3D
0769    * and should not be used as a separate class.
0770    *
0771    * Example of use:
0772    * @code
0773    *   Transform3D m;
0774    *   m = ScaleY3D(2.);
0775    * @endcode
0776    *
0777    * @author <Evgueni.Tcherniaev@cern.ch>
0778    * @ingroup geometry
0779    */
0780   class ScaleY3D : public Scale3D {
0781   public:
0782     /**
0783      * Default constructor: sets the Identity transformation. */
0784     ScaleY3D() : Scale3D() {}
0785 
0786     /**
0787      * Constructor from a number (scale factor in y-direction). */
0788     ScaleY3D(double y) : Scale3D(1, y, 1) {}
0789   };
0790 
0791   /**
0792    * Constructs a scaling transformation in z-direction.
0793    * This class provides additional constructors for Transform3D
0794    * and should not be used as a separate class.
0795    *
0796    * Example of use:
0797    * @code
0798    *   Transform3D m;
0799    *   m = ScaleZ3D(2.);
0800    * @endcode
0801    *
0802    * @author <Evgueni.Tcherniaev@cern.ch>
0803    * @ingroup geometry
0804    */
0805   class ScaleZ3D : public Scale3D {
0806   public:
0807     /**
0808      * Default constructor: sets the Identity transformation. */
0809     ScaleZ3D() : Scale3D() {}
0810     /**
0811      * Constructor from a number (scale factor in z-direction). */
0812     ScaleZ3D(double z) : Scale3D(1, 1, z) {}
0813   };
0814 } /* namespace HepGeom */
0815 
0816 #ifdef ENABLE_BACKWARDS_COMPATIBILITY
0817 //  backwards compatibility will be enabled ONLY in CLHEP 1.9
0818 typedef HepGeom::Transform3D  HepTransform3D;
0819 typedef HepGeom::Rotate3D     HepRotate3D;
0820 typedef HepGeom::RotateX3D    HepRotateX3D;
0821 typedef HepGeom::RotateY3D    HepRotateY3D;
0822 typedef HepGeom::RotateZ3D    HepRotateZ3D;
0823 typedef HepGeom::Translate3D  HepTranslate3D;
0824 typedef HepGeom::TranslateX3D HepTranslateX3D;
0825 typedef HepGeom::TranslateY3D HepTranslateY3D;
0826 typedef HepGeom::TranslateZ3D HepTranslateZ3D;
0827 typedef HepGeom::Reflect3D    HepReflect3D;
0828 typedef HepGeom::ReflectX3D   HepReflectX3D;
0829 typedef HepGeom::ReflectY3D   HepReflectY3D;
0830 typedef HepGeom::ReflectZ3D   HepReflectZ3D;
0831 typedef HepGeom::Scale3D      HepScale3D;
0832 typedef HepGeom::ScaleX3D     HepScaleX3D;
0833 typedef HepGeom::ScaleY3D     HepScaleY3D;
0834 typedef HepGeom::ScaleZ3D     HepScaleZ3D;
0835 #endif
0836 
0837 #include "CLHEP/Geometry/Transform3D.icc"
0838 
0839 #endif /* HEP_TRANSFROM3D_H */