CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_6_1_2_SLHC2_patch1/src/DataFormats/GeometrySurface/interface/oldTkRotation.h

Go to the documentation of this file.
00001 #ifndef Geom_oldTkRotation_H
00002 #define Geom_oldTkRotation_H
00003 
00004 #include "DataFormats/GeometryVector/interface/Basic2DVector.h"
00005 #include "DataFormats/GeometryVector/interface/Basic3DVector.h"
00006 #include "DataFormats/GeometryVector/interface/GlobalVector.h"
00007 /*
00008 #include "DataFormats/GeometrySurface/interface/GlobalError.h"
00009 #include "DataFormats/GeometrySurface/interface/LocalError.h"
00010 */
00011 #include <iosfwd>
00012 
00013 template <class T> class TkRotation;
00014 
00015 template <class T>
00016 std::ostream & operator<<( std::ostream& s, const TkRotation<T>& r);
00017 
00018 namespace geometryDetails {
00019   void TkRotationErr1();
00020   void TkRotationErr2();
00021 
00022 }
00023 
00024 
00028 template <class T>
00029 class TkRotation {
00030 public:
00031   typedef Vector3DBase< T, GlobalTag>  GlobalVector;
00032 
00033   TkRotation() : 
00034     R11( 1), R12( 0), R13( 0), 
00035     R21( 0), R22( 1), R23( 0),
00036     R31( 0), R32( 0), R33( 1) {}
00037 
00038   TkRotation( T xx, T xy, T xz, T yx, T yy, T yz, T zx, T zy, T zz) :
00039     R11(xx), R12(xy), R13(xz), 
00040     R21(yx), R22(yy), R23(yz),
00041     R31(zx), R32(zy), R33(zz) {}
00042 
00043   TkRotation( const T* p) :
00044     R11(p[0]), R12(p[1]), R13(p[2]), 
00045     R21(p[3]), R22(p[4]), R23(p[5]),
00046     R31(p[6]), R32(p[7]), R33(p[8]) {}
00047 
00048   TkRotation( const GlobalVector & aX, const GlobalVector & aY)  {
00049 
00050     GlobalVector uX = aX.unit();
00051     GlobalVector uY = aY.unit();
00052     GlobalVector uZ(uX.cross(uY));
00053  
00054     R11 = uX.x(); R12 = uX.y();  R13 = uX.z();
00055     R21 = uY.x(); R22 = uY.y();  R23 = uY.z();
00056     R31 = uZ.x(); R32 = uZ.y();  R33 = uZ.z();
00057 
00058   }
00059 
00064   TkRotation( const GlobalVector & aX, const GlobalVector & aY, 
00065               const GlobalVector & aZ) :
00066     R11( aX.x()), R12( aX.y()), R13( aX.z()), 
00067     R21( aY.x()), R22( aY.y()), R23( aY.z()),
00068     R31( aZ.x()), R32( aZ.y()), R33( aZ.z()) {}
00069     
00070 
00079   TkRotation( const Basic3DVector<T>& axis, T phi) :
00080     R11( cos(phi) ), R12( sin(phi)), R13( 0), 
00081     R21( -sin(phi)), R22( cos(phi)), R23( 0),
00082     R31( 0), R32( 0), R33( 1) {
00083 
00084     //rotation around the z-axis by  phi
00085     TkRotation tmpRotz ( cos(axis.phi()), sin(axis.phi()), 0.,
00086                         -sin(axis.phi()), cos(axis.phi()), 0.,
00087                          0.,              0.,              1. );
00088     //rotation around y-axis by theta
00089     TkRotation tmpRoty ( cos(axis.theta()), 0.,-sin(axis.theta()),
00090                          0.,              1.,              0.,
00091                          sin(axis.theta()), 0., cos(axis.theta()) );
00092     (*this)*=tmpRoty;
00093     (*this)*=tmpRotz;      // =  this * tmpRoty * tmpRotz 
00094    
00095     // (tmpRoty * tmpRotz)^-1 * this * tmpRoty * tmpRotz
00096 
00097     *this = (tmpRoty*tmpRotz).multiplyInverse(*this);
00098 
00099   }
00100   /* this is the same thing...derived from the CLHEP ... it gives the
00101      same results MODULO the sign of the rotation....  but I don't want
00102      that... had 
00103   TkRotation (const Basic3DVector<T>& axis, float phi) {
00104     T cx = axis.x();
00105     T cy = axis.y();
00106     T cz = axis.z();
00107     
00108     T ll = axis.mag();
00109     if (ll == 0) {
00110       geometryDetails::TkRotationErr1();
00111     }else{
00112       
00113       float cosa = cos(phi), sina = sin(phi);
00114       cx /= ll; cy /= ll; cz /= ll;   
00115       
00116        R11 = cosa + (1-cosa)*cx*cx;
00117        R12 =        (1-cosa)*cx*cy - sina*cz;
00118        R13 =        (1-cosa)*cx*cz + sina*cy;
00119       
00120        R21 =        (1-cosa)*cy*cx + sina*cz;
00121        R22 = cosa + (1-cosa)*cy*cy; 
00122        R23 =        (1-cosa)*cy*cz - sina*cx;
00123       
00124        R31 =        (1-cosa)*cz*cx - sina*cy;
00125        R32 =        (1-cosa)*cz*cy + sina*cx;
00126        R33 = cosa + (1-cosa)*cz*cz;
00127       
00128     }
00129     
00130   }
00131   */
00132 
00133   template <typename U>
00134   TkRotation( const TkRotation<U>& a) : 
00135     R11(a.xx()), R12(a.xy()), R13(a.xz()), 
00136     R21(a.yx()), R22(a.yy()), R23(a.yz()),
00137     R31(a.zx()), R32(a.zy()), R33(a.zz()) {}
00138 
00139   TkRotation transposed() const {
00140       return TkRotation( R11, R21, R31, 
00141                          R12, R22, R32,
00142                          R13, R23, R33);
00143   }
00144 
00145   Basic3DVector<T> operator*( const Basic3DVector<T>& v) const {
00146     return Basic3DVector<T>( R11*v.x() + R12*v.y() + R13*v.z(),
00147                              R21*v.x() + R22*v.y() + R23*v.z(),
00148                              R31*v.x() + R32*v.y() + R33*v.z());
00149   }
00150 
00151   Basic3DVector<T> multiplyInverse( const Basic3DVector<T>& v) const {
00152     return Basic3DVector<T>( R11*v.x() + R21*v.y() + R31*v.z(),
00153                              R12*v.x() + R22*v.y() + R32*v.z(),
00154                              R13*v.x() + R23*v.y() + R33*v.z());
00155   }
00156 
00157 #ifndef CMS_NO_TEMPLATE_MEMBERS
00158   template <class Scalar>
00159   Basic3DVector<Scalar> multiplyInverse( const Basic3DVector<Scalar>& v) const {
00160     return Basic3DVector<Scalar>( R11*v.x() + R21*v.y() + R31*v.z(),
00161                                   R12*v.x() + R22*v.y() + R32*v.z(),
00162                                   R13*v.x() + R23*v.y() + R33*v.z());
00163   }
00164 #endif
00165 
00166   Basic3DVector<T> operator*( const Basic2DVector<T>& v) const {
00167     return Basic3DVector<T>( R11*v.x() + R12*v.y(),
00168                              R21*v.x() + R22*v.y(),
00169                              R31*v.x() + R32*v.y());
00170   }
00171   Basic3DVector<T> multiplyInverse( const Basic2DVector<T>& v) const {
00172     return Basic3DVector<T>( R11*v.x() + R21*v.y(),
00173                              R12*v.x() + R22*v.y(),
00174                              R13*v.x() + R23*v.y());
00175   }
00176 
00177   
00178 
00179   TkRotation operator*( const TkRotation& b) const {
00180     return TkRotation(R11*b.R11 + R12*b.R21 + R13*b.R31,
00181                       R11*b.R12 + R12*b.R22 + R13*b.R32,
00182                       R11*b.R13 + R12*b.R23 + R13*b.R33,
00183                       R21*b.R11 + R22*b.R21 + R23*b.R31,
00184                       R21*b.R12 + R22*b.R22 + R23*b.R32,
00185                       R21*b.R13 + R22*b.R23 + R23*b.R33,
00186                       R31*b.R11 + R32*b.R21 + R33*b.R31,
00187                       R31*b.R12 + R32*b.R22 + R33*b.R32,
00188                       R31*b.R13 + R32*b.R23 + R33*b.R33);
00189   }
00190 
00191   TkRotation multiplyInverse( const TkRotation& b) const {
00192     return TkRotation(R11*b.R11 + R21*b.R21 + R31*b.R31,
00193                       R11*b.R12 + R21*b.R22 + R31*b.R32,
00194                       R11*b.R13 + R21*b.R23 + R31*b.R33,
00195                       R12*b.R11 + R22*b.R21 + R32*b.R31,
00196                       R12*b.R12 + R22*b.R22 + R32*b.R32,
00197                       R12*b.R13 + R22*b.R23 + R32*b.R33,
00198                       R13*b.R11 + R23*b.R21 + R33*b.R31,
00199                       R13*b.R12 + R23*b.R22 + R33*b.R32,
00200                       R13*b.R13 + R23*b.R23 + R33*b.R33);
00201   }
00202 
00203   TkRotation& operator*=( const TkRotation& b) {
00204     return *this = operator * (b);
00205   }
00206 
00207   // Note a *= b; <=> a = a * b; while a.transform(b); <=> a = b * a;
00208   TkRotation& transform(const TkRotation& b) {
00209     return *this = b.operator * (*this);
00210   }
00211 
00212   TkRotation & rotateAxes(const Basic3DVector<T>& newX,
00213                           const Basic3DVector<T>& newY,
00214                           const Basic3DVector<T>& newZ) {
00215     T del = 0.001;
00216 
00217     if (
00218         
00219         // the check for right-handedness is not needed since
00220         // we want to change the handedness when it's left in cmsim
00221         //
00222         //       fabs(newZ.x()-w.x()) > del ||
00223         //       fabs(newZ.y()-w.y()) > del ||
00224         //       fabs(newZ.z()-w.z()) > del ||
00225         fabs(newX.mag2()-1.) > del ||
00226         fabs(newY.mag2()-1.) > del || 
00227         fabs(newZ.mag2()-1.) > del ||
00228         fabs(newX.dot(newY)) > del ||
00229         fabs(newY.dot(newZ)) > del ||
00230         fabs(newZ.dot(newX)) > del) {
00231       geometryDetails::TkRotationErr2();
00232       return *this;
00233     } else {
00234       return transform(TkRotation(newX.x(), newY.x(), newZ.x(),
00235                                   newX.y(), newY.y(), newZ.y(),
00236                                   newX.z(), newY.z(), newZ.z()));
00237     }
00238   }
00239 
00240   Basic3DVector<T> x() const { return  Basic3DVector<T>(xx(),xy(),xz());}
00241   Basic3DVector<T> y() const { return  Basic3DVector<T>(yx(),yy(),yz());}
00242   Basic3DVector<T> z() const { return  Basic3DVector<T>(zx(),zy(),zz());}
00243 
00244 
00245   T const &xx() const { return R11;} 
00246   T const &xy() const { return R12;} 
00247   T const &xz() const { return R13;} 
00248   T const &yx() const { return R21;} 
00249   T const &yy() const { return R22;} 
00250   T const &yz() const { return R23;} 
00251   T const &zx() const { return R31;} 
00252   T const &zy() const { return R32;} 
00253   T const &zz() const { return R33;} 
00254 
00255 private:
00256 
00257   T R11, R12, R13;
00258   T R21, R22, R23;
00259   T R31, R32, R33;
00260 };
00261 
00262 
00263 template<>
00264 std::ostream & operator<< <float>( std::ostream& s, const TkRotation<float>& r);
00265 
00266 template<>
00267 std::ostream & operator<< <double>( std::ostream& s, const TkRotation<double>& r);
00268 
00269 
00270 template <class T, class U>
00271 inline Basic3DVector<U> operator*( const TkRotation<T>& r, const Basic3DVector<U>& v) {
00272   return Basic3DVector<U>( r.xx()*v.x() + r.xy()*v.y() + r.xz()*v.z(),
00273                            r.yx()*v.x() + r.yy()*v.y() + r.yz()*v.z(),
00274                            r.zx()*v.x() + r.zy()*v.y() + r.zz()*v.z());
00275 }
00276 
00277 template <class T, class U>
00278 inline TkRotation<typename PreciseFloatType<T,U>::Type>
00279 operator*( const TkRotation<T>& a, const TkRotation<U>& b) {
00280     typedef TkRotation<typename PreciseFloatType<T,U>::Type> RT;
00281     return RT( a.xx()*b.xx() + a.xy()*b.yx() + a.xz()*b.zx(),
00282                a.xx()*b.xy() + a.xy()*b.yy() + a.xz()*b.zy(),
00283                a.xx()*b.xz() + a.xy()*b.yz() + a.xz()*b.zz(),
00284                a.yx()*b.xx() + a.yy()*b.yx() + a.yz()*b.zx(),
00285                a.yx()*b.xy() + a.yy()*b.yy() + a.yz()*b.zy(),
00286                a.yx()*b.xz() + a.yy()*b.yz() + a.yz()*b.zz(),
00287                a.zx()*b.xx() + a.zy()*b.yx() + a.zz()*b.zx(),
00288                a.zx()*b.xy() + a.zy()*b.yy() + a.zz()*b.zy(),
00289                a.zx()*b.xz() + a.zy()*b.yz() + a.zz()*b.zz());
00290 }
00291 
00292 #endif
00293 
00294 
00295 
00296