Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Par02FastSimModelTracker Class Reference

#include <Par02FastSimModelTracker.hh>

Inheritance diagram for Par02FastSimModelTracker:
Collaboration diagram for Par02FastSimModelTracker:

Public Member Functions

 Par02FastSimModelTracker (G4String aModelName, G4Region *aEnvelope, Par02DetectorParametrisation::Parametrisation aParamType)
 
 Par02FastSimModelTracker (G4String aModelName, G4Region *aEnvelope)
 
 Par02FastSimModelTracker (G4String aModelName)
 
 ~Par02FastSimModelTracker ()
 
virtual G4bool IsApplicable (const G4ParticleDefinition &aParticle)
 
virtual G4bool ModelTrigger (const G4FastTrack &aFastTrack)
 
virtual void DoIt (const G4FastTrack &aFastTrack, G4FastStep &aFastStep)
 
- Public Member Functions inherited from G4VFastSimulationModel
 G4VFastSimulationModel (const G4String &aName)
 
 G4VFastSimulationModel (const G4String &aName, G4Envelope *, G4bool IsUnique=FALSE)
 
virtual ~G4VFastSimulationModel ()
 
virtual G4bool AtRestModelTrigger (const G4FastTrack &)
 
virtual void AtRestDoIt (const G4FastTrack &, G4FastStep &)
 
const G4String GetName () const
 
G4bool operator== (const G4VFastSimulationModel &) const
 

Detailed Description

Shortcut to the ordinary tracking for tracking detectors.

The fast simulation model describes what should be done instead of a normal tracking. Instead of the ordinary tracking, a particle momentum at the entrance of the tracking detector is smeared (by Par02Smearer::SmearMomentum()) and the particle is placed at the tracking detector exit, at the place it would reach without the change of its momentum. Based on G4 examples/extended/parametrisations/Par01/include/Par01EMShowerModel.hh .

Author
Anna Zaborowska

Definition at line 50 of file Par02FastSimModelTracker.hh.

Constructor & Destructor Documentation

Par02FastSimModelTracker::Par02FastSimModelTracker ( G4String  aModelName,
G4Region aEnvelope,
Par02DetectorParametrisation::Parametrisation  aParamType 
)

A constructor.

Parameters
aModelNameA name of the fast simulation model.
aEnvelopeA region where the model can take over the ordinary tracking.
aParamTypeA parametrisation type.

Definition at line 55 of file Par02FastSimModelTracker.cc.

56  :
57  G4VFastSimulationModel( aModelName, aEnvelope ), fCalculateParametrisation(),
58  fParametrisation( aType ) {}
G4VFastSimulationModel(const G4String &aName)
Par02FastSimModelTracker::Par02FastSimModelTracker ( G4String  aModelName,
G4Region aEnvelope 
)

A constructor.

Parameters
aModelNameA name of the fast simulation model.
aEnvelopeA region where the model can take over the ordinary tracking.

Definition at line 62 of file Par02FastSimModelTracker.cc.

63  :
64  G4VFastSimulationModel( aModelName, aEnvelope ), fCalculateParametrisation(),
65  fParametrisation( Par02DetectorParametrisation::eCMS ) {}
G4VFastSimulationModel(const G4String &aName)
Par02FastSimModelTracker::Par02FastSimModelTracker ( G4String  aModelName)

A constructor.

Parameters
aModelNameA name of the fast simulation model.

Definition at line 69 of file Par02FastSimModelTracker.cc.

69  :
70  G4VFastSimulationModel( aModelName ), fCalculateParametrisation(),
71  fParametrisation( Par02DetectorParametrisation::eCMS ) {}
G4VFastSimulationModel(const G4String &aName)
Par02FastSimModelTracker::~Par02FastSimModelTracker ( )

Definition at line 75 of file Par02FastSimModelTracker.cc.

75 {}

Member Function Documentation

void Par02FastSimModelTracker::DoIt ( const G4FastTrack aFastTrack,
G4FastStep aFastStep 
)
virtual

Calculates the final position (at the outer boundary of the tracking detector) of a particle with the momentum at the entrance of the tracking detector. Smears the particle momentum and saves it, together with the tracking detector resolution and efficiency to the Par02PrimaryParticleInformation.

Parameters
aFastTrackA track.
aFastStepA step.

Implements G4VFastSimulationModel.

Definition at line 92 of file Par02FastSimModelTracker.cc.

93  {
94 
95  G4cout << " ________Tracker model triggered _________" << G4endl;
96 
97  // Calculate the final position (at the outer boundary of the tracking detector)
98  // of the particle with the momentum at the entrance of the tracking detector.
99 
100  G4Track track = * aFastTrack.GetPrimaryTrack();
101  G4FieldTrack aFieldTrack( '0' );
102  G4FieldTrackUpdator::Update( &aFieldTrack, &track );
103 
104  G4double retSafety = -1.0;
105  ELimited retStepLimited;
106  G4FieldTrack endTrack( 'a' );
107  G4double currentMinimumStep = 10.0*m; // Temporary: change that to sth connected
108  // to particle momentum.
109  G4PathFinder* fPathFinder = G4PathFinder::GetInstance();
110  /*G4double lengthAlongCurve = */
111  fPathFinder->ComputeStep( aFieldTrack,
112  currentMinimumStep,
113  0,
114  aFastTrack.GetPrimaryTrack()->GetCurrentStepNumber(),
115  retSafety,
116  retStepLimited,
117  endTrack,
118  aFastTrack.GetPrimaryTrack()->GetVolume() );
119 
120  // Place the particle at the tracking detector exit
121  // (at the place it would reach without the change of its momentum).
122  aFastStep.ProposePrimaryTrackFinalPosition( endTrack.GetPosition() );
123 
124  // Consider only primary tracks (do nothing else for secondary charged particles)
125  G4ThreeVector Porg = aFastTrack.GetPrimaryTrack()->GetMomentum();
126  if ( ! aFastTrack.GetPrimaryTrack()->GetParentID() ) {
129  if ( info->GetDoSmearing() ) {
130  // Smearing according to the tracking detector resolution
131  G4double res = fCalculateParametrisation->
133  fParametrisation, Porg.mag() );
134  G4double eff = fCalculateParametrisation->
136  fParametrisation, Porg.mag() );
137  G4ThreeVector Psm;
138  Psm = Par02Smearer::Instance()->
139  SmearMomentum( aFastTrack.GetPrimaryTrack(), res );
140  Par02Output::Instance()->FillHistogram( 0, ((Psm.mag()/MeV) / (Porg.mag()/MeV)) );
141  // Setting the values of Psm, res and eff
142  ( (Par02PrimaryParticleInformation*) ( const_cast< G4PrimaryParticle* >
143  ( aFastTrack.GetPrimaryTrack()->GetDynamicParticle()->GetPrimaryParticle() )->
144  GetUserInformation() ) )->SetTrackerMomentum( Psm );
145  ( (Par02PrimaryParticleInformation*) ( const_cast< G4PrimaryParticle* >
146  ( aFastTrack.GetPrimaryTrack()->GetDynamicParticle()->GetPrimaryParticle() )->
147  GetUserInformation() ) )->SetTrackerResolution( res );
148  ( (Par02PrimaryParticleInformation*) ( const_cast< G4PrimaryParticle* >
149  ( aFastTrack.GetPrimaryTrack()->GetDynamicParticle()->GetPrimaryParticle() )->
150  GetUserInformation() ) )->SetTrackerEfficiency( eff );
151  } else {
152  // No smearing: simply setting the value of Porg
153  ( (Par02PrimaryParticleInformation*) ( const_cast< G4PrimaryParticle* >
154  ( aFastTrack.GetPrimaryTrack()->GetDynamicParticle()->GetPrimaryParticle() )->
155  GetUserInformation() ) )->SetTrackerMomentum( Porg );
156  }
157  }
158 }
static G4PathFinder * GetInstance()
Definition: G4PathFinder.cc:57
const XML_Char XML_Encoding * info
Definition: expat.h:530
G4int GetParentID() const
const G4Track * GetPrimaryTrack() const
Definition: G4FastTrack.hh:208
static Par02Smearer * Instance()
Definition: Par02Smearer.cc:59
const G4DynamicParticle * GetDynamicParticle() const
ELimited
void FillHistogram(G4int HNo, G4double value) const
Definition: Par02Output.cc:219
void ProposePrimaryTrackFinalPosition(const G4ThreeVector &, G4bool localCoordinates=true)
Definition: G4FastStep.cc:99
G4GLOB_DLL std::ostream G4cout
G4int GetCurrentStepNumber() const
static constexpr double m
Definition: G4SIunits.hh:129
G4double ComputeStep(const G4FieldTrack &pFieldTrack, G4double pCurrentProposedStepLength, G4int navigatorId, G4int stepNo, G4double &pNewSafety, ELimited &limitedStep, G4FieldTrack &EndState, G4VPhysicalVolume *currentVolume)
static Par02Output * Instance()
Definition: Par02Output.cc:60
G4PrimaryParticle * GetPrimaryParticle() const
G4ThreeVector GetMomentum() const
static G4EventManager * GetEventManager()
G4VPhysicalVolume * GetVolume() const
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
double G4double
Definition: G4Types.hh:76
double mag() const
G4bool GetDoSmearing()
Gets the flag indicating if smearing should be done.
static void Update(G4FieldTrack *, const G4Track *)
G4VUserEventInformation * GetUserInformation()

Here is the call graph for this function:

G4bool Par02FastSimModelTracker::IsApplicable ( const G4ParticleDefinition aParticle)
virtual

Checks if this model should be applied to this particle type.

Parameters
aParticleA particle definition (type).

Implements G4VFastSimulationModel.

Definition at line 79 of file Par02FastSimModelTracker.cc.

80  {
81  return aParticleType.GetPDGCharge() != 0; // Applicable for all charged particles
82 }

Here is the call graph for this function:

G4bool Par02FastSimModelTracker::ModelTrigger ( const G4FastTrack aFastTrack)
virtual

Checks if the model should be applied taking into account the kinematics of a track.

Parameters
aFastTrackA track.

Implements G4VFastSimulationModel.

Definition at line 86 of file Par02FastSimModelTracker.cc.

86  {
87  return true; // No kinematical restrictions to apply the parametrisation
88 }

The documentation for this class was generated from the following files: