CMS 3D CMS Logo

List of all members | Public Member Functions | Private Types | Private Member Functions | Private Attributes
Geant4ePropagator Class Reference

#include <Geant4ePropagator.h>

Inheritance diagram for Geant4ePropagator:
Propagator

Public Member Functions

Geant4ePropagatorclone () const override
 
template<>
bool configureAnyPropagation (G4ErrorMode &mode, Plane const &pDest, GlobalPoint const &cmsInitPos, GlobalVector const &cmsInitMom) const
 
template<>
bool configureAnyPropagation (G4ErrorMode &mode, Cylinder const &pDest, GlobalPoint const &cmsInitPos, GlobalVector const &cmsInitMom) const
 
 Geant4ePropagator (const MagneticField *field=nullptr, std::string particleName="mu", PropagationDirection dir=alongMomentum, double plimit=1.0)
 
template<>
std::string getSurfaceType (Cylinder const &c) const
 
template<>
std::string getSurfaceType (Plane const &c) const
 
const MagneticFieldmagneticField () const override
 
std::pair< TrajectoryStateOnSurface, double > propagateWithPath (const FreeTrajectoryState &, const Plane &) const override
 
std::pair< TrajectoryStateOnSurface, double > propagateWithPath (const FreeTrajectoryState &, const Cylinder &) const override
 
std::pair< TrajectoryStateOnSurface, double > propagateWithPath (const TrajectoryStateOnSurface &, const Plane &) const override
 
std::pair< TrajectoryStateOnSurface, double > propagateWithPath (const TrajectoryStateOnSurface &, const Cylinder &) const override
 
template<>
Geant4ePropagator::ErrorTargetPair transformToG4SurfaceTarget (const Plane &pDest, bool moveTargetToEndOfSurface) const
 
template<>
Geant4ePropagator::ErrorTargetPair transformToG4SurfaceTarget (const Cylinder &pDest, bool moveTargetToEndOfSurface) const
 
 ~Geant4ePropagator () override
 
- Public Member Functions inherited from Propagator
template<typename STA , typename SUR >
TrajectoryStateOnSurface propagate (STA const &state, SUR const &surface) const
 
virtual FreeTrajectoryState propagate (const FreeTrajectoryState &ftsStart, const GlobalPoint &pDest) const final
 
virtual FreeTrajectoryState propagate (const FreeTrajectoryState &ftsStart, const GlobalPoint &pDest1, const GlobalPoint &pDest2) const final
 
virtual FreeTrajectoryState propagate (const FreeTrajectoryState &ftsStart, const reco::BeamSpot &beamSpot) const final
 
virtual std::pair< TrajectoryStateOnSurface, double > propagateWithPath (const FreeTrajectoryState &, const Surface &) const final
 
virtual std::pair< TrajectoryStateOnSurface, double > propagateWithPath (const TrajectoryStateOnSurface &tsos, const Surface &sur) const final
 
virtual std::pair< FreeTrajectoryState, double > propagateWithPath (const FreeTrajectoryState &ftsStart, const GlobalPoint &pDest) const
 
virtual std::pair< FreeTrajectoryState, double > propagateWithPath (const FreeTrajectoryState &ftsStart, const GlobalPoint &pDest1, const GlobalPoint &pDest2) const
 Propagate to PCA to a line (given by 2 points) given a starting point. More...
 
virtual std::pair< FreeTrajectoryState, double > propagateWithPath (const FreeTrajectoryState &ftsStart, const reco::BeamSpot &beamSpot) const
 Propagate to PCA to a line (given by beamSpot position and slope) given a starting point. More...
 
virtual PropagationDirection propagationDirection () const final
 
 Propagator (PropagationDirection dir=alongMomentum)
 
virtual bool setMaxDirectionChange (float phiMax)
 
virtual void setPropagationDirection (PropagationDirection dir)
 
virtual ~Propagator ()
 

Private Types

typedef std::pair< bool, std::shared_ptr< G4ErrorTarget > > ErrorTargetPair
 
typedef std::pair< TrajectoryStateOnSurface, double > TsosPP
 

Private Member Functions

template<class SurfaceType >
bool configureAnyPropagation (G4ErrorMode &mode, SurfaceType const &pDest, GlobalPoint const &cmsInitPos, GlobalVector const &cmsInitMom) const
 
template<class SurfaceType >
bool configurePropagation (G4ErrorMode &mode, SurfaceType const &pDest, GlobalPoint const &cmsInitPos, GlobalVector const &cmsInitMom) const
 
void debugReportPlaneSetup (GlobalPoint const &posPlane, HepGeom::Point3D< double > const &surfPos, GlobalVector const &normalPlane, HepGeom::Normal3D< double > const &surfNorm, const Plane &pDest) const
 
template<class SurfaceType >
void debugReportTrackState (std::string const &currentContext, GlobalPoint const &cmsInitPos, CLHEP::Hep3Vector const &g4InitPos, GlobalVector const &cmsInitMom, CLHEP::Hep3Vector const &g4InitMom, const SurfaceType &pDest) const
 
void ensureGeant4eIsInitilized (bool forceInit) const
 
std::string generateParticleName (int charge) const
 
template<class SurfaceType >
std::string getSurfaceType (SurfaceType const &surface) const
 
template<class SurfaceType >
std::pair< TrajectoryStateOnSurface, double > propagateGeneric (const FreeTrajectoryState &ftsStart, const SurfaceType &pDest) const
 
template<class SurfaceType >
ErrorTargetPair transformToG4SurfaceTarget (const SurfaceType &pDest, bool moveTargetToEndOfSurface) const
 

Private Attributes

double plimit_
 
const MagneticFieldtheField
 
G4ErrorPropagatorData * theG4eData
 
G4ErrorPropagatorManager * theG4eManager
 
std::string theParticleName
 

Detailed Description

Propagator based on the Geant4e package. Uses the Propagator class in the TrackingTools/GeomPropagators package to define the interface. See that class for more details.

Definition at line 20 of file Geant4ePropagator.h.

Member Typedef Documentation

◆ ErrorTargetPair

typedef std::pair<bool, std::shared_ptr<G4ErrorTarget> > Geant4ePropagator::ErrorTargetPair
private

Definition at line 60 of file Geant4ePropagator.h.

◆ TsosPP

typedef std::pair<TrajectoryStateOnSurface, double> Geant4ePropagator::TsosPP
private

Definition at line 59 of file Geant4ePropagator.h.

Constructor & Destructor Documentation

◆ Geant4ePropagator()

Geant4ePropagator::Geant4ePropagator ( const MagneticField field = nullptr,
std::string  particleName = "mu",
PropagationDirection  dir = alongMomentum,
double  plimit = 1.0 
)

Constructor. Takes as arguments:

  • The magnetic field
  • The particle name whose properties will be used in the propagation. Without the charge, i.e. "mu", "pi", ...
  • The propagation direction. It may be: alongMomentum, oppositeToMomentum

Constructor.

Definition at line 42 of file Geant4ePropagator.cc.

References ensureGeant4eIsInitilized(), and LogDebug.

Referenced by clone().

46  : Propagator(dir),
47  theField(field),
49  theG4eManager(G4ErrorPropagatorManager::GetErrorPropagatorManager()),
50  theG4eData(G4ErrorPropagatorData::GetErrorPropagatorData()),
51  plimit_(plimit) {
52  LogDebug("Geant4e") << "Geant4e Propagator initialized";
53 
54  // has to be called here, doing it later will not load the G4 physics list
55  // properly when using the G4 ES Producer. Reason: unclear
57 }
Propagator(PropagationDirection dir=alongMomentum)
Definition: Propagator.h:46
std::string theParticleName
const MagneticField * theField
G4ErrorPropagatorManager * theG4eManager
void ensureGeant4eIsInitilized(bool forceInit) const
G4ErrorPropagatorData * theG4eData
#define LogDebug(id)

◆ ~Geant4ePropagator()

Geant4ePropagator::~Geant4ePropagator ( )
override

Destructor.

Definition at line 61 of file Geant4ePropagator.cc.

References LogDebug.

61  {
62  LogDebug("Geant4e") << "Geant4ePropagator::~Geant4ePropagator()" << std::endl;
63 
64  // don't close the g4 Geometry here, because the propagator might have been
65  // cloned
66  // but there is only one, globally-shared Geometry
67 }
#define LogDebug(id)

Member Function Documentation

◆ clone()

Geant4ePropagator* Geant4ePropagator::clone ( void  ) const
inlineoverridevirtual

Implements Propagator.

Definition at line 54 of file Geant4ePropagator.h.

References Geant4ePropagator().

54 { return new Geant4ePropagator(*this); }
Geant4ePropagator(const MagneticField *field=nullptr, std::string particleName="mu", PropagationDirection dir=alongMomentum, double plimit=1.0)

◆ configureAnyPropagation() [1/3]

template<class SurfaceType >
bool Geant4ePropagator::configureAnyPropagation ( G4ErrorMode &  mode,
SurfaceType const &  pDest,
GlobalPoint const &  cmsInitPos,
GlobalVector const &  cmsInitMom 
) const
private

Referenced by configurePropagation().

◆ configureAnyPropagation() [2/3]

template<>
bool Geant4ePropagator::configureAnyPropagation ( G4ErrorMode &  mode,
Plane const &  pDest,
GlobalPoint const &  cmsInitPos,
GlobalVector const &  cmsInitMom 
) const

Definition at line 160 of file Geant4ePropagator.cc.

References Plane::localZ(), LogDebug, PV3DBase< T, PVType, FrameType >::mag(), ALCARECOPromptCalibProdSiPixelAli0T_cff::mode, and plimit_.

163  {
164  if (cmsInitMom.mag() < plimit_)
165  return false;
166  if (pDest.localZ(cmsInitPos) * pDest.localZ(cmsInitMom) < 0) {
167  mode = G4ErrorMode_PropForwards;
168  LogDebug("Geant4e") << "G4e - Propagator mode is \'forwards\' indirect "
169  "via the Any direction"
170  << std::endl;
171  } else {
172  mode = G4ErrorMode_PropBackwards;
173  LogDebug("Geant4e") << "G4e - Propagator mode is \'backwards\' indirect "
174  "via the Any direction"
175  << std::endl;
176  }
177 
178  return true;
179 }
#define LogDebug(id)

◆ configureAnyPropagation() [3/3]

template<>
bool Geant4ePropagator::configureAnyPropagation ( G4ErrorMode &  mode,
Cylinder const &  pDest,
GlobalPoint const &  cmsInitPos,
GlobalVector const &  cmsInitMom 
) const

Definition at line 182 of file Geant4ePropagator.cc.

References LogDebug, PV3DBase< T, PVType, FrameType >::mag(), ALCARECOPromptCalibProdSiPixelAli0T_cff::mode, plimit_, SurfaceOrientation::positiveSide, Cylinder::side(), and GloballyPositioned< T >::toLocal().

185  {
186  if (cmsInitMom.mag() < plimit_)
187  return false;
188  //------------------------------------
189  // For cylinder assume outside is backwards, inside is along
190  // General use for particles from collisions
191  LocalPoint lpos = pDest.toLocal(cmsInitPos);
192  Surface::Side theSide = pDest.side(lpos, 0);
193  if (theSide == SurfaceOrientation::positiveSide) { // outside cylinder
194  mode = G4ErrorMode_PropBackwards;
195  LogDebug("Geant4e") << "G4e - Propagator mode is \'backwards\' indirect "
196  "via the Any direction";
197  } else { // inside cylinder
198  mode = G4ErrorMode_PropForwards;
199  LogDebug("Geant4e") << "G4e - Propagator mode is \'forwards\' indirect "
200  "via the Any direction";
201  }
202 
203  return true;
204 }
#define LogDebug(id)

◆ configurePropagation()

template<class SurfaceType >
bool Geant4ePropagator::configurePropagation ( G4ErrorMode &  mode,
SurfaceType const &  pDest,
GlobalPoint const &  cmsInitPos,
GlobalVector const &  cmsInitMom 
) const
private

Definition at line 207 of file Geant4ePropagator.cc.

References alongMomentum, anyDirection, configureAnyPropagation(), LogDebug, PV3DBase< T, PVType, FrameType >::mag(), ALCARECOPromptCalibProdSiPixelAli0T_cff::mode, oppositeToMomentum, plimit_, and Propagator::propagationDirection().

Referenced by propagateGeneric().

210  {
211  if (cmsInitMom.mag() < plimit_)
212  return false;
214  mode = G4ErrorMode_PropBackwards;
215  LogDebug("Geant4e") << "G4e - Propagator mode is \'backwards\' " << std::endl;
216  } else if (propagationDirection() == alongMomentum) {
217  mode = G4ErrorMode_PropForwards;
218  LogDebug("Geant4e") << "G4e - Propagator mode is \'forwards\'" << std::endl;
219  } else if (propagationDirection() == anyDirection) {
220  if (configureAnyPropagation(mode, pDest, cmsInitPos, cmsInitMom) == false)
221  return false;
222  } else {
223  edm::LogError("Geant4e") << "G4e - Unsupported propagation mode";
224  return false;
225  }
226  return true;
227 }
virtual PropagationDirection propagationDirection() const final
Definition: Propagator.h:139
Log< level::Error, false > LogError
bool configureAnyPropagation(G4ErrorMode &mode, SurfaceType const &pDest, GlobalPoint const &cmsInitPos, GlobalVector const &cmsInitMom) const
#define LogDebug(id)

◆ debugReportPlaneSetup()

void Geant4ePropagator::debugReportPlaneSetup ( GlobalPoint const &  posPlane,
HepGeom::Point3D< double > const &  surfPos,
GlobalVector const &  normalPlane,
HepGeom::Normal3D< double > const &  surfNorm,
const Plane pDest 
) const
private

Definition at line 455 of file Geant4ePropagator.cc.

References Geom::Phi< T1, Range >::degrees(), PV3DBase< T, PVType, FrameType >::eta(), Plane::localZ(), LogDebug, PV3DBase< T, PVType, FrameType >::perp(), and PV3DBase< T, PVType, FrameType >::phi().

459  {
460  LogDebug("Geant4e") << "G4e - Destination CMS plane position:" << posPlane << "cm\n"
461  << "G4e - (Ro, eta, phi): (" << posPlane.perp() << " cm, " << posPlane.eta()
462  << ", " << posPlane.phi().degrees() << " deg)\n"
463  << "G4e - Destination G4 plane position: " << surfPos << " mm, Ro = " << surfPos.perp()
464  << " mm";
465  LogDebug("Geant4e") << "G4e - Destination CMS plane normal : " << normalPlane << "\n"
466  << "G4e - Destination G4 plane normal : " << normalPlane;
467  LogDebug("Geant4e") << "G4e - Distance from plane position to plane: " << pDest.localZ(posPlane) << " cm";
468 }
float localZ(const GlobalPoint &gp) const
Definition: Plane.h:45
#define LogDebug(id)

◆ debugReportTrackState()

template<class SurfaceType >
void Geant4ePropagator::debugReportTrackState ( std::string const &  currentContext,
GlobalPoint const &  cmsInitPos,
CLHEP::Hep3Vector const &  g4InitPos,
GlobalVector const &  cmsInitMom,
CLHEP::Hep3Vector const &  g4InitMom,
const SurfaceType &  pDest 
) const
private

Definition at line 471 of file Geant4ePropagator.cc.

References Geom::Phi< T1, Range >::degrees(), PV3DBase< T, PVType, FrameType >::eta(), LogDebug, PV3DBase< T, PVType, FrameType >::perp(), and PV3DBase< T, PVType, FrameType >::phi().

Referenced by propagateGeneric().

476  {
477  LogDebug("Geant4e") << "G4e - Current Context: " << currentContext;
478  LogDebug("Geant4e") << "G4e - CMS point position:" << cmsInitPos << "cm\n"
479  << "G4e - (Ro, eta, phi): (" << cmsInitPos.perp() << " cm, " << cmsInitPos.eta()
480  << ", " << cmsInitPos.phi().degrees() << " deg)\n"
481  << "G4e - G4 point position: " << g4InitPos << " mm, Ro = " << g4InitPos.perp() << " mm";
482  LogDebug("Geant4e") << "G4e - CMS momentum :" << cmsInitMom << "GeV\n"
483  << " pt: " << cmsInitMom.perp() << "G4e - G4 momentum : " << g4InitMom << " MeV";
484 }
#define LogDebug(id)

◆ ensureGeant4eIsInitilized()

void Geant4ePropagator::ensureGeant4eIsInitilized ( bool  forceInit) const
private

Propagate from a free state (e.g. position and momentum in in global cartesian coordinates) to a plane.

Definition at line 77 of file Geant4ePropagator.cc.

References LogDebug, theField, and theG4eManager.

Referenced by Geant4ePropagator().

77  {
78  LogDebug("Geant4ePropagator") << "G4 propagator starts isInitialized, theField: " << theField;
79 
80  auto man = G4RunManagerKernel::GetRunManagerKernel();
81  if (G4StateManager::GetStateManager()->GetCurrentState() == G4State_PreInit) {
82  man->SetVerboseLevel(0);
83  theG4eManager->InitGeant4e();
84 
85  // define 10 mm step limit for propagator
86  G4UImanager::GetUIpointer()->ApplyCommand("/geant4e/limits/stepLength 10.0 mm");
87  }
88  const G4Field *field = G4TransportationManager::GetTransportationManager()->GetFieldManager()->GetDetectorField();
89  if (field == nullptr) {
90  edm::LogError("Geant4e") << "No G4 magnetic field defined";
91  }
92  LogDebug("Geant4ePropagator") << "G4 propagator initialized; field: " << field;
93 }
Log< level::Error, false > LogError
const MagneticField * theField
G4ErrorPropagatorManager * theG4eManager
#define LogDebug(id)

◆ generateParticleName()

std::string Geant4ePropagator::generateParticleName ( int  charge) const
private

Definition at line 144 of file Geant4ePropagator.cc.

References ALCARECOTkAlJpsiMuMu_cff::charge, LogDebug, HiggsValidation_cfi::particleName, AlCaHLTBitMon_QueryRunRegistry::string, and theParticleName.

Referenced by propagateGeneric().

144  {
146 
147  if (charge > 0) {
148  particleName += "+";
149  }
150  if (charge < 0) {
151  particleName += "-";
152  }
153 
154  LogDebug("Geant4e") << "G4e - Particle name: " << particleName;
155 
156  return particleName;
157 }
std::string theParticleName
#define LogDebug(id)

◆ getSurfaceType() [1/3]

template<class SurfaceType >
std::string Geant4ePropagator::getSurfaceType ( SurfaceType const &  surface) const
private

◆ getSurfaceType() [2/3]

template<>
std::string Geant4ePropagator::getSurfaceType ( Cylinder const &  c) const

Definition at line 135 of file Geant4ePropagator.cc.

135  {
136  return "Cylinder";
137 }

◆ getSurfaceType() [3/3]

template<>
std::string Geant4ePropagator::getSurfaceType ( Plane const &  c) const

Definition at line 140 of file Geant4ePropagator.cc.

140  {
141  return "Plane";
142 }

◆ magneticField()

const MagneticField* Geant4ePropagator::magneticField ( ) const
inlineoverridevirtual

Implements Propagator.

Definition at line 56 of file Geant4ePropagator.h.

References theField.

56 { return theField; }
const MagneticField * theField

◆ propagateGeneric()

template<class SurfaceType >
std::pair< TrajectoryStateOnSurface, double > Geant4ePropagator::propagateGeneric ( const FreeTrajectoryState ftsStart,
const SurfaceType &  pDest 
) const
private

Definition at line 230 of file Geant4ePropagator.cc.

References SurfaceSideDefinition::afterSurface, TrackPropagation::algebraicSymMatrix55ToG4ErrorTrajErr(), alongMomentum, SurfaceSideDefinition::beforeSurface, FreeTrajectoryState::charge(), GlobalTrajectoryParameters::charge(), configurePropagation(), FreeTrajectoryState::curvilinearError(), debugReportTrackState(), f, TrackPropagation::g4doubleToCmsDouble(), TrackPropagation::g4ErrorTrajErrToAlgebraicSymMatrix55(), generateParticleName(), TrackPropagation::globalPointToHep3Vector(), TrackPropagation::globalVectorToHep3Vector(), FreeTrajectoryState::hasError(), TrackPropagation::hep3VectorToGlobalVector(), TrackPropagation::hepPoint3DToGlobalPoint(), LogDebug, GlobalTrajectoryParameters::magneticField(), CurvilinearTrajectoryError::matrix(), ALCARECOPromptCalibProdSiPixelAli0T_cff::mode, FreeTrajectoryState::momentum(), oppositeToMomentum, FreeTrajectoryState::parameters(), FreeTrajectoryState::position(), Propagator::propagationDirection(), ALPAKA_ACCELERATOR_NAMESPACE::ecal::reconstruction::internal::barrel::side(), theField, theG4eData, theG4eManager, and transformToG4SurfaceTarget().

Referenced by propagateWithPath().

231  {
233  // Construct the target surface
234  //
235  //* Set the target surface
236 
237  ErrorTargetPair g4eTarget_center = transformToG4SurfaceTarget(pDest, false);
238 
239  // * Get the starting point and direction and convert them to
240  // CLHEP::Hep3Vector
241  // for G4. CMS uses cm and GeV while Geant4 uses mm and MeV
242  GlobalPoint cmsInitPos = ftsStart.position();
243  GlobalVector cmsInitMom = ftsStart.momentum();
244  bool flipped = false;
246  // flip the momentum vector as Geant4 will not do this
247  // on it's own in a backward propagation
248  cmsInitMom = -cmsInitMom;
249  flipped = true;
250  }
251 
252  // Set the mode of propagation according to the propagation direction
253  G4ErrorMode mode = G4ErrorMode_PropForwards;
254  if (!configurePropagation(mode, pDest, cmsInitPos, cmsInitMom))
255  return TsosPP(TrajectoryStateOnSurface(), 0.0f);
256 
257  // re-check propagation direction chosen in case of AnyDirection
258  if (mode == G4ErrorMode_PropBackwards && !flipped)
259  cmsInitMom = -cmsInitMom;
260 
261  CLHEP::Hep3Vector g4InitPos = TrackPropagation::globalPointToHep3Vector(cmsInitPos);
262  CLHEP::Hep3Vector g4InitMom = TrackPropagation::globalVectorToHep3Vector(cmsInitMom * GeV);
263 
264  debugReportTrackState("intitial", cmsInitPos, g4InitPos, cmsInitMom, g4InitMom, pDest);
265 
266  // Set the mode of propagation according to the propagation direction
267  // G4ErrorMode mode = G4ErrorMode_PropForwards;
268 
269  // if (!configurePropagation(mode, pDest, cmsInitPos, cmsInitMom))
270  // return TsosPP(TrajectoryStateOnSurface(), 0.0f);
271 
273  // Set the error and trajectories, and finally propagate
274  //
275  G4ErrorTrajErr g4error(5, 1);
276  if (ftsStart.hasError()) {
278  initErr = ftsStart.curvilinearError();
279  g4error = TrackPropagation::algebraicSymMatrix55ToG4ErrorTrajErr(initErr, ftsStart.charge());
280  LogDebug("Geant4e") << "CMS - Error matrix: " << std::endl << initErr.matrix();
281  } else {
282  LogDebug("Geant4e") << "No error matrix available" << std::endl;
283  return TsosPP(TrajectoryStateOnSurface(), 0.0f);
284  }
285 
286  LogDebug("Geant4e") << "G4e - Error matrix: " << std::endl << g4error;
287 
288  // in CMSSW, the state errors are deflated when performing the backward
289  // propagation
290  if (mode == G4ErrorMode_PropForwards) {
291  G4ErrorPropagatorData::GetErrorPropagatorData()->SetStage(G4ErrorStage_Inflation);
292  } else if (mode == G4ErrorMode_PropBackwards) {
293  G4ErrorPropagatorData::GetErrorPropagatorData()->SetStage(G4ErrorStage_Deflation);
294  }
295 
296  G4ErrorFreeTrajState g4eTrajState(generateParticleName(ftsStart.charge()), g4InitPos, g4InitMom, g4error);
297  LogDebug("Geant4e") << "G4e - Traj. State: " << (g4eTrajState);
298 
300  // Propagate
301  int iterations = 0;
302  double finalPathLength = 0;
303 
304  HepGeom::Point3D<double> finalRecoPos;
305 
306  G4ErrorPropagatorData::GetErrorPropagatorData()->SetMode(mode);
307 
308  theG4eData->SetTarget(g4eTarget_center.second.get());
309  LogDebug("Geant4e") << "Running Propagation to the RECO surface" << std::endl;
310 
311  theG4eManager->InitTrackPropagation();
312 
313  // re-initialize navigator to avoid mismatches and/or segfaults
314  theG4eManager->GetErrorPropagationNavigator()->LocateGlobalPointAndSetup(
315  g4InitPos, &g4InitMom, /*pRelativeSearch = */ false, /*ignoreDirection = */ false);
316 
317  bool continuePropagation = true;
318  while (continuePropagation) {
319  iterations++;
320  LogDebug("Geant4e") << std::endl << "step count " << iterations << " step length " << finalPathLength;
321 
322  // re-initialize navigator to avoid mismatches and/or segfaults
323  theG4eManager->GetErrorPropagationNavigator()->LocateGlobalPointWithinVolume(g4eTrajState.GetPosition());
324 
325  const int ierr = theG4eManager->PropagateOneStep(&g4eTrajState, mode);
326 
327  if (ierr != 0) {
328  // propagation failed, return invalid track state
329  return TsosPP(TrajectoryStateOnSurface(), 0.0f);
330  }
331 
332  const float thisPathLength = TrackPropagation::g4doubleToCmsDouble(g4eTrajState.GetG4Track()->GetStepLength());
333 
334  LogDebug("Geant4e") << "step Length was " << thisPathLength << " cm, current global position: "
335  << TrackPropagation::hepPoint3DToGlobalPoint(g4eTrajState.GetPosition()) << std::endl;
336 
337  finalPathLength += thisPathLength;
338 
339  // if (std::fabs(finalPathLength) > 10000.0f)
340  if (std::fabs(finalPathLength) > 200.0f) {
341  LogDebug("Geant4e") << "ERROR: Quitting propagation: path length mega large" << std::endl;
342  theG4eManager->GetPropagator()->InvokePostUserTrackingAction(g4eTrajState.GetG4Track());
343  continuePropagation = false;
344  LogDebug("Geant4e") << "WARNING: Quitting propagation: max path length "
345  "exceeded, returning invalid state"
346  << std::endl;
347 
348  // reached maximum path length, bail out
349  return TsosPP(TrajectoryStateOnSurface(), 0.0f);
350  }
351 
352  if (theG4eManager->GetPropagator()->CheckIfLastStep(g4eTrajState.GetG4Track())) {
353  theG4eManager->GetPropagator()->InvokePostUserTrackingAction(g4eTrajState.GetG4Track());
354  continuePropagation = false;
355  }
356  }
357 
358  // CMSSW Tracking convention, backward propagations have negative path length
360  finalPathLength = -finalPathLength;
361 
362  // store the correct location for the hit on the RECO surface
363  LogDebug("Geant4e") << "Position on the RECO surface" << g4eTrajState.GetPosition() << std::endl;
364  finalRecoPos = g4eTrajState.GetPosition();
365 
366  theG4eManager->EventTermination();
367 
368  LogDebug("Geant4e") << "Final position of the Track :" << g4eTrajState.GetPosition() << std::endl;
369 
371  // Retrieve the state in the end from Geant4e, convert them to CMS vectors
372  // and points, and build global trajectory parameters.
373  // CMS uses cm and GeV while Geant4 uses mm and MeV
374  //
375  const HepGeom::Vector3D<double> momEnd = g4eTrajState.GetMomentum();
376 
377  // use the hit on the the RECO plane as the final position to be d'accor with
378  // the RecHit measurements
379  const GlobalPoint posEndGV = TrackPropagation::hepPoint3DToGlobalPoint(finalRecoPos);
381 
382  debugReportTrackState("final", posEndGV, finalRecoPos, momEndGV, momEnd, pDest);
383 
384  // Get the error covariance matrix from Geant4e. It comes in curvilinear
385  // coordinates so use the appropiate CMS class
386  G4ErrorTrajErr g4errorEnd = g4eTrajState.GetError();
387 
388  CurvilinearTrajectoryError curvError(
390 
391  if (mode == G4ErrorMode_PropBackwards) {
393  posEndGV, momEndGV, ftsStart.parameters().charge(), &ftsStart.parameters().magneticField());
394 
395  // flip the momentum direction because it has been flipped before running
396  // G4's backwards prop
397  momEndGV = -momEndGV;
398  }
399 
400  LogDebug("Geant4e") << "G4e - Error matrix after propagation: " << std::endl << g4errorEnd;
401 
402  LogDebug("Geant4e") << "CMS - Error matrix after propagation: " << std::endl << curvError.matrix();
403 
404  GlobalTrajectoryParameters tParsDest(posEndGV, momEndGV, ftsStart.charge(), theField);
405 
407 
410 
411  return TsosPP(TrajectoryStateOnSurface(tParsDest, curvError, pDest, side), finalPathLength);
412 }
std::pair< TrajectoryStateOnSurface, double > TsosPP
const CurvilinearTrajectoryError & curvilinearError() const
void debugReportTrackState(std::string const &currentContext, GlobalPoint const &cmsInitPos, CLHEP::Hep3Vector const &g4InitPos, GlobalVector const &cmsInitMom, CLHEP::Hep3Vector const &g4InitMom, const SurfaceType &pDest) const
double g4doubleToCmsDouble(const G4double &d)
virtual PropagationDirection propagationDirection() const final
Definition: Propagator.h:139
CLHEP::Hep3Vector globalVectorToHep3Vector(const GlobalVector &p)
const MagneticField * theField
const GlobalTrajectoryParameters & parameters() const
G4ErrorPropagatorManager * theG4eManager
GlobalPoint position() const
AlgebraicSymMatrix55 g4ErrorTrajErrToAlgebraicSymMatrix55(const G4ErrorTrajErr &e, const int q)
G4ErrorPropagatorData * theG4eData
TrackCharge charge() const
GlobalPoint hepPoint3DToGlobalPoint(const HepGeom::Point3D< double > &r)
std::pair< bool, std::shared_ptr< G4ErrorTarget > > ErrorTargetPair
GlobalVector momentum() const
double f[11][100]
const AlgebraicSymMatrix55 & matrix() const
const MagneticField & magneticField() const
GlobalVector hep3VectorToGlobalVector(const CLHEP::Hep3Vector &p)
std::string generateParticleName(int charge) const
bool configurePropagation(G4ErrorMode &mode, SurfaceType const &pDest, GlobalPoint const &cmsInitPos, GlobalVector const &cmsInitMom) const
G4ErrorTrajErr algebraicSymMatrix55ToG4ErrorTrajErr(const AlgebraicSymMatrix55 &e, const int q)
ErrorTargetPair transformToG4SurfaceTarget(const SurfaceType &pDest, bool moveTargetToEndOfSurface) const
CLHEP::Hep3Vector globalPointToHep3Vector(const GlobalPoint &r)
#define LogDebug(id)

◆ propagateWithPath() [1/4]

std::pair< TrajectoryStateOnSurface, double > Geant4ePropagator::propagateWithPath ( const FreeTrajectoryState ftsStart,
const Plane pDest 
) const
overridevirtual

The methods propagateWithPath() are identical to the corresponding methods propagate() in what concerns the resulting TrajectoryStateOnSurface, but they provide in addition the exact path length along the trajectory. All of these method calls are internally mapped to

The methods propagateWithPath() are identical to the corresponding methods propagate() in what concerns the resulting TrajectoryStateOnSurface, but they provide in addition the exact path length along the trajectory.

Implements Propagator.

Definition at line 424 of file Geant4ePropagator.cc.

References propagateGeneric().

425  {
426  // Finally build the pair<...> that needs to be returned where the second
427  // parameter is the exact path length. Currently calculated with a stepping
428  // action that adds up the length of every step
429  return propagateGeneric(ftsStart, pDest);
430 }
std::pair< TrajectoryStateOnSurface, double > propagateGeneric(const FreeTrajectoryState &ftsStart, const SurfaceType &pDest) const

◆ propagateWithPath() [2/4]

std::pair< TrajectoryStateOnSurface, double > Geant4ePropagator::propagateWithPath ( const FreeTrajectoryState ftsStart,
const Cylinder cDest 
) const
overridevirtual

Implements Propagator.

Definition at line 432 of file Geant4ePropagator.cc.

References propagateGeneric().

433  {
434  // Finally build the pair<...> that needs to be returned where the second
435  // parameter is the exact path length.
436  return propagateGeneric(ftsStart, cDest);
437 }
std::pair< TrajectoryStateOnSurface, double > propagateGeneric(const FreeTrajectoryState &ftsStart, const SurfaceType &pDest) const

◆ propagateWithPath() [3/4]

std::pair< TrajectoryStateOnSurface, double > Geant4ePropagator::propagateWithPath ( const TrajectoryStateOnSurface tsosStart,
const Plane pDest 
) const
overridevirtual

Reimplemented from Propagator.

Definition at line 439 of file Geant4ePropagator.cc.

References TrajectoryStateOnSurface::freeState(), and propagateGeneric().

440  {
441  // Finally build the pair<...> that needs to be returned where the second
442  // parameter is the exact path length.
443  const FreeTrajectoryState ftsStart = *tsosStart.freeState();
444  return propagateGeneric(ftsStart, pDest);
445 }
FreeTrajectoryState const * freeState(bool withErrors=true) const
std::pair< TrajectoryStateOnSurface, double > propagateGeneric(const FreeTrajectoryState &ftsStart, const SurfaceType &pDest) const

◆ propagateWithPath() [4/4]

std::pair< TrajectoryStateOnSurface, double > Geant4ePropagator::propagateWithPath ( const TrajectoryStateOnSurface tsosStart,
const Cylinder cDest 
) const
overridevirtual

Reimplemented from Propagator.

Definition at line 447 of file Geant4ePropagator.cc.

References TrajectoryStateOnSurface::freeState(), and propagateGeneric().

448  {
449  const FreeTrajectoryState ftsStart = *tsosStart.freeState();
450  // Finally build the pair<...> that needs to be returned where the second
451  // parameter is the exact path length.
452  return propagateGeneric(ftsStart, cDest);
453 }
FreeTrajectoryState const * freeState(bool withErrors=true) const
std::pair< TrajectoryStateOnSurface, double > propagateGeneric(const FreeTrajectoryState &ftsStart, const SurfaceType &pDest) const

◆ transformToG4SurfaceTarget() [1/3]

template<class SurfaceType >
ErrorTargetPair Geant4ePropagator::transformToG4SurfaceTarget ( const SurfaceType &  pDest,
bool  moveTargetToEndOfSurface 
) const
private

Referenced by propagateGeneric().

◆ transformToG4SurfaceTarget() [2/3]

template<>
Geant4ePropagator::ErrorTargetPair Geant4ePropagator::transformToG4SurfaceTarget ( const Plane pDest,
bool  moveTargetToEndOfSurface 
) const

Definition at line 96 of file Geant4ePropagator.cc.

References TrackPropagation::globalPointToHepPoint3D(), TrackPropagation::globalVectorToHepNormal3D(), Surface::toGlobal(), and Vector3DBase< T, FrameTag >::unit().

97  {
98  //* Get position and normal (orientation) of the destination plane
99  GlobalPoint posPlane = pDest.toGlobal(LocalPoint(0, 0, 0));
100  GlobalVector normalPlane = pDest.toGlobal(LocalVector(0, 0, 1.));
101  normalPlane = normalPlane.unit();
102 
103  //* Transform this into HepGeom::Point3D<double> and
104  // HepGeom::Normal3D<double> that define a plane for
105  // Geant4e.
106  // CMS uses cm and GeV while Geant4 uses mm and MeV
107  HepGeom::Point3D<double> surfPos = TrackPropagation::globalPointToHepPoint3D(posPlane);
108  HepGeom::Normal3D<double> surfNorm = TrackPropagation::globalVectorToHepNormal3D(normalPlane);
109 
110  //* Set the target surface
111  return ErrorTargetPair(false, std::make_shared<G4ErrorPlaneSurfaceTarget>(surfNorm, surfPos));
112 }
Local3DVector LocalVector
Definition: LocalVector.h:12
Point3DBase< Scalar, LocalTag > LocalPoint
Definition: Definitions.h:30
HepGeom::Normal3D< double > globalVectorToHepNormal3D(const GlobalVector &p)
std::pair< bool, std::shared_ptr< G4ErrorTarget > > ErrorTargetPair
GlobalPoint toGlobal(const Point2DBase< Scalar, LocalTag > lp) const
Definition: Surface.h:79
HepGeom::Point3D< double > globalPointToHepPoint3D(const GlobalPoint &r)
Vector3DBase unit() const
Definition: Vector3DBase.h:54

◆ transformToG4SurfaceTarget() [3/3]

template<>
Geant4ePropagator::ErrorTargetPair Geant4ePropagator::transformToG4SurfaceTarget ( const Cylinder pDest,
bool  moveTargetToEndOfSurface 
) const

Definition at line 115 of file Geant4ePropagator.cc.

References TrackPropagation::globalPointToHep3Vector(), LogDebug, GloballyPositioned< T >::position(), Cylinder::radius(), idealTransformation::rotation, GloballyPositioned< T >::rotation(), and TrackPropagation::tkRotationFToHepRotation().

116  {
117  // Get Cylinder parameters.
118  // CMS uses cm and GeV while Geant4 uses mm and MeV.
119  // - Radius
120  G4float radCyl = pDest.radius() * cm;
121  // - Position: PositionType & GlobalPoint are Basic3DPoint<float,GlobalTag>
122  G4ThreeVector posCyl = TrackPropagation::globalPointToHep3Vector(pDest.position());
123  // - Rotation: Type in CMSSW is RotationType == TkRotation<T>, T=float
124  G4RotationMatrix rotCyl = TrackPropagation::tkRotationFToHepRotation(pDest.rotation());
125 
126  // DEBUG
128  LogDebug("Geant4e") << "G4e - TkRotation" << rotation;
129  LogDebug("Geant4e") << "G4e - G4Rotation" << rotCyl << "mm";
130 
131  return ErrorTargetPair(!moveTargetToEndOfSurface, std::make_shared<G4ErrorCylSurfaceTarget>(radCyl, posCyl, rotCyl));
132 }
CLHEP::HepRotation tkRotationFToHepRotation(const TkRotation< float > &tkr)
std::pair< bool, std::shared_ptr< G4ErrorTarget > > ErrorTargetPair
const PositionType & position() const
Scalar radius() const
Radius of the cylinder.
Definition: Cylinder.h:64
const RotationType & rotation() const
CLHEP::Hep3Vector globalPointToHep3Vector(const GlobalPoint &r)
#define LogDebug(id)

Member Data Documentation

◆ plimit_

double Geant4ePropagator::plimit_
private

Definition at line 71 of file Geant4ePropagator.h.

Referenced by configureAnyPropagation(), and configurePropagation().

◆ theField

const MagneticField* Geant4ePropagator::theField
private

Definition at line 63 of file Geant4ePropagator.h.

Referenced by ensureGeant4eIsInitilized(), magneticField(), and propagateGeneric().

◆ theG4eData

G4ErrorPropagatorData* Geant4ePropagator::theG4eData
private

Definition at line 70 of file Geant4ePropagator.h.

Referenced by propagateGeneric().

◆ theG4eManager

G4ErrorPropagatorManager* Geant4ePropagator::theG4eManager
private

Definition at line 69 of file Geant4ePropagator.h.

Referenced by ensureGeant4eIsInitilized(), and propagateGeneric().

◆ theParticleName

std::string Geant4ePropagator::theParticleName
private

Definition at line 66 of file Geant4ePropagator.h.

Referenced by generateParticleName().