Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:31:54

0001 #ifndef TrackPropagation_ConvertFromToCLHEP_h
0002 #define TrackPropagation_ConvertFromToCLHEP_h
0003 
0004 // CLHEP
0005 #include "CLHEP/Geometry/Normal3D.h"
0006 #include "CLHEP/Geometry/Point3D.h"
0007 #include "CLHEP/Units/GlobalSystemOfUnits.h"
0008 #include "CLHEP/Vector/Rotation.h"
0009 #include "CLHEP/Vector/ThreeVector.h"
0010 
0011 // CMS
0012 #include "DataFormats/CLHEP/interface/AlgebraicObjects.h"
0013 #include "DataFormats/GeometrySurface/interface/TkRotation.h"
0014 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0015 #include "DataFormats/GeometryVector/interface/GlobalVector.h"
0016 
0017 // Geant4
0018 #include "G4ErrorFreeTrajState.hh"
0019 
0020 /** Utilities to convert among CLHEP and CMS points and vectors
0021  */
0022 
0023 namespace TrackPropagation {
0024   /**
0025  Convert a CMS GlobalPoint to a CLHEP HepGeom::Point3D<double>
0026  CMS uses cm while Geant4 uses mm. This is taken into account in the
0027  conversion.
0028  */
0029 
0030   inline HepGeom::Point3D<double> globalPointToHepPoint3D(const GlobalPoint &r) {
0031     return HepGeom::Point3D<double>(r.x() * cm, r.y() * cm, r.z() * cm);
0032   }
0033 
0034   /** Convert a CLHEP HepGeom::Point3D<double>  to a CMS GlobalPoint
0035  CMS uses cms while Geant4 uses mm. This is taken into account in the
0036  conversion.
0037  */
0038   inline GlobalPoint hepPoint3DToGlobalPoint(const HepGeom::Point3D<double> &r) {
0039     return GlobalPoint(r.x() / cm, r.y() / cm, r.z() / cm);
0040   }
0041 
0042   /** Convert a G4double representing a scaler measure ( e.g. Track length ) to
0043  the CMS convention, which is using mm
0044  */
0045   inline double g4doubleToCmsDouble(const G4double &d) { return d / cm; }
0046 
0047   /** Convert a CMS GlobalVector to a CLHEP HepGeom::Normal3D<double>
0048  CMS uses GeV while G4 uses MeV
0049  */
0050   inline HepGeom::Normal3D<double> globalVectorToHepNormal3D(const GlobalVector &p) {
0051     return HepGeom::Normal3D<double>(p.x(), p.y(), p.z());
0052   }
0053 
0054   /** Convert a CLHEP HepGeom::Normal3D<double>  to a CMS GlobalVector
0055  CMS uses GeV while G4 uses MeV
0056  */
0057   inline GlobalVector hepNormal3DToGlobalVector(const HepGeom::Normal3D<double> &p) {
0058     return GlobalVector(p.x(), p.y(), p.z());
0059   }
0060 
0061   /** Convert a CMS GlobalVector to a CLHEP CLHEP::Hep3Vector
0062  */
0063   inline CLHEP::Hep3Vector globalVectorToHep3Vector(const GlobalVector &p) {
0064     return CLHEP::Hep3Vector(p.x(), p.y(), p.z());
0065   }
0066 
0067   /** Convert a CLHEP CLHEP::Hep3Vector to a CMS GlobalVector
0068  */
0069   inline GlobalVector hep3VectorToGlobalVector(const CLHEP::Hep3Vector &p) { return GlobalVector(p.x(), p.y(), p.z()); }
0070 
0071   /** Convert a CMS GlobalPoint to a CLHEP CLHEP::Hep3Vector
0072  CMS uses cm while Geant4 uses mm. This is taken into account in the
0073  conversion.
0074  */
0075   inline CLHEP::Hep3Vector globalPointToHep3Vector(const GlobalPoint &r) {
0076     return CLHEP::Hep3Vector(r.x() * cm, r.y() * cm, r.z() * cm);
0077   }
0078 
0079   /** Convert a CLHEP CLHEP::Hep3Vector to a CMS GlobalPoint
0080  CMS uses cm while Geant4 uses mm. This is taken into account in the
0081  conversion.
0082  */
0083   inline GlobalPoint hep3VectorToGlobalPoint(const CLHEP::Hep3Vector &v) {
0084     return GlobalPoint(v.x() / cm, v.y() / cm, v.z() / cm);
0085   }
0086 
0087   /** Convert a CMS TkRotation<float> to a CLHEP
0088  * CLHEP::HepRotation=G4RotationMatrix
0089  */
0090   inline CLHEP::HepRotation tkRotationFToHepRotation(const TkRotation<float> &tkr) {
0091     return CLHEP::HepRotation(CLHEP::Hep3Vector(tkr.xx(), tkr.yx(), tkr.zx()),
0092                               CLHEP::Hep3Vector(tkr.xy(), tkr.yy(), tkr.zy()),
0093                               CLHEP::Hep3Vector(tkr.xz(), tkr.yz(), tkr.zz()));
0094   }
0095 
0096   /** Convert a CLHEP CLHEP::Hep3Vector to a CMS GlobalPoint
0097  */
0098   inline TkRotation<float> hepRotationToTkRotationF(const CLHEP::HepRotation &r) {
0099     return TkRotation<float>(r.xx(), r.xy(), r.xz(), r.yx(), r.yy(), r.yz(), r.zx(), r.zy(), r.zz());
0100   }
0101 
0102   /** Convert a G4 Trajectory Error Matrix to the CMS Algebraic Sym Matrix
0103  CMS uses q/p as first parameter, G4 uses 1/p
0104  */
0105 
0106   inline AlgebraicSymMatrix55 g4ErrorTrajErrToAlgebraicSymMatrix55(const G4ErrorTrajErr &e, const int q) {
0107     assert(q != 0);
0108     // From DataFormats/CLHEP/interface/Migration.h
0109     // typedef ROOT::Math::SMatrix<double,5,5,ROOT::Math::MatRepSym<double,5> >
0110     // AlgebraicSymMatrix55;
0111     AlgebraicSymMatrix55 m55;
0112     for (unsigned int i = 0; i < 5; i++)
0113       for (unsigned int j = 0; j < 5; j++) {
0114         m55(i, j) = e(i + 1, j + 1);
0115         if (i == 0)
0116           m55(i, j) = double(q) * m55(i, j);
0117         if (j == 0)
0118           m55(i, j) = double(q) * m55(i, j);
0119       }
0120     return m55;
0121   }
0122 
0123   /** Convert a CMS Algebraic Sym Matrix (for curv error) to a G4 Trajectory Error
0124  * Matrix
0125  */
0126   inline G4ErrorTrajErr algebraicSymMatrix55ToG4ErrorTrajErr(const AlgebraicSymMatrix55 &e, const int q) {
0127     assert(q != 0);
0128     G4ErrorTrajErr g4err(5, 1);
0129     for (unsigned int i = 0; i < 5; i++)
0130       for (unsigned int j = 0; j < 5; j++) {
0131         g4err(i + 1, j + 1) = e(i, j);
0132 
0133         if (i == 0)
0134           g4err(i + 1, j + 1) = g4err(i + 1, j + 1) * double(q);
0135         if (j == 0)
0136           g4err(i + 1, j + 1) = g4err(i + 1, j + 1) * double(q);
0137       }
0138     return g4err;
0139   }
0140 
0141 }  // namespace TrackPropagation
0142 
0143 #endif