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

#include <Par02FastSimModelEMCal.hh>

Inheritance diagram for Par02FastSimModelEMCal:
Collaboration diagram for Par02FastSimModelEMCal:

Public Member Functions

 Par02FastSimModelEMCal (G4String aModelName, G4Region *aEnvelope, Par02DetectorParametrisation::Parametrisation aParamType)
 
 Par02FastSimModelEMCal (G4String aModelName, G4Region *aEnvelope)
 
 Par02FastSimModelEMCal (G4String aModelName)
 
 ~Par02FastSimModelEMCal ()
 
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 electromagnetic calorimeters.

The fast simulation model describes what should be done instead of a normal tracking. Instead of the ordinary tracking, a particle deposits its energy at the entrance to the electromagnetic calorimeter and its value is smeared (by Par02Smearer::SmearMomentum()). Based on G4 examples/extended/parametrisations/Par01/include/Par01EMShowerModel.hh .

Author
Anna Zaborowska

Definition at line 47 of file Par02FastSimModelEMCal.hh.

Constructor & Destructor Documentation

Par02FastSimModelEMCal::Par02FastSimModelEMCal ( 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 51 of file Par02FastSimModelEMCal.cc.

52  :
53  G4VFastSimulationModel( aModelName, aEnvelope ), fCalculateParametrisation(),
54  fParametrisation( aType ) {}
G4VFastSimulationModel(const G4String &aName)
Par02FastSimModelEMCal::Par02FastSimModelEMCal ( 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 58 of file Par02FastSimModelEMCal.cc.

59  :
60  G4VFastSimulationModel( aModelName, aEnvelope ), fCalculateParametrisation(),
61  fParametrisation( Par02DetectorParametrisation::eCMS ) {}
G4VFastSimulationModel(const G4String &aName)
Par02FastSimModelEMCal::Par02FastSimModelEMCal ( G4String  aModelName)

A constructor.

Parameters
aModelNameA name of the fast simulation model.

Definition at line 65 of file Par02FastSimModelEMCal.cc.

65  :
66  G4VFastSimulationModel( aModelName ), fCalculateParametrisation(),
67  fParametrisation( Par02DetectorParametrisation::eCMS ) {}
G4VFastSimulationModel(const G4String &aName)
Par02FastSimModelEMCal::~Par02FastSimModelEMCal ( )

Definition at line 71 of file Par02FastSimModelEMCal.cc.

71 {}

Member Function Documentation

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

Smears the energy deposit and saves it, together with the position of the deposit, the electromagnetic calorimeter resolution and efficiency to the Par02PrimaryParticleInformation.

Parameters
aFastTrackA track.
aFastStepA step.

Implements G4VFastSimulationModel.

Definition at line 91 of file Par02FastSimModelEMCal.cc.

92  {
93  //G4cout << " ________EMCal model triggered _________" << G4endl;
94 
95  // Kill the parameterised particle at the entrance of the electromagnetic calorimeter
96  aFastStep.KillPrimaryTrack();
97  aFastStep.ProposePrimaryTrackPathLength( 0.0 );
98  G4double Edep = aFastTrack.GetPrimaryTrack()->GetKineticEnergy();
99 
100  // Consider only primary tracks (do nothing else for secondary e-, e+, gammas)
101  G4ThreeVector Pos = aFastTrack.GetPrimaryTrack()->GetPosition();
102  if ( ! aFastTrack.GetPrimaryTrack()->GetParentID() ) {
105  if ( info->GetDoSmearing() ) {
106  // Smearing according to the electromagnetic calorimeter resolution
107  G4ThreeVector Porg = aFastTrack.GetPrimaryTrack()->GetMomentum();
108  G4double res = fCalculateParametrisation->GetResolution(
109  Par02DetectorParametrisation::eEMCAL, fParametrisation, Porg.mag() );
110  G4double eff = fCalculateParametrisation->GetEfficiency(
111  Par02DetectorParametrisation::eEMCAL, fParametrisation, Porg.mag() );
112  G4double Esm;
113  Esm = std::abs( Par02Smearer::Instance()->
114  SmearEnergy( aFastTrack.GetPrimaryTrack(), res ) );
115  Par02Output::Instance()->FillHistogram( 1, (Esm/MeV) / (Edep/MeV) );
116  // Setting the values of Pos, Esm, res and eff
117  ( (Par02PrimaryParticleInformation*) ( const_cast< G4PrimaryParticle* >
118  ( aFastTrack.GetPrimaryTrack()->GetDynamicParticle()->GetPrimaryParticle() )->
119  GetUserInformation() ) )->SetEMCalPosition( Pos );
120  ( (Par02PrimaryParticleInformation*) ( const_cast< G4PrimaryParticle* >
121  ( aFastTrack.GetPrimaryTrack()->GetDynamicParticle()->GetPrimaryParticle() )->
122  GetUserInformation() ) )->SetEMCalEnergy( Esm );
123  ( (Par02PrimaryParticleInformation*) ( const_cast< G4PrimaryParticle* >
124  ( aFastTrack.GetPrimaryTrack()->GetDynamicParticle()->GetPrimaryParticle() )->
125  GetUserInformation() ) )->SetEMCalResolution( res );
126  ( (Par02PrimaryParticleInformation*) ( const_cast< G4PrimaryParticle* >
127  ( aFastTrack.GetPrimaryTrack()->GetDynamicParticle()->GetPrimaryParticle() )->
128  GetUserInformation() ) )->SetEMCalEfficiency( eff );
129  // The (smeared) energy of the particle is deposited in the step
130  // (which corresponds to the entrance of the electromagnetic calorimeter)
131  aFastStep.ProposeTotalEnergyDeposited( Esm );
132  } else {
133  // No smearing: simply setting the value of Edep
134  ( (Par02PrimaryParticleInformation*) ( const_cast< G4PrimaryParticle* >
135  ( aFastTrack.GetPrimaryTrack()->GetDynamicParticle()->GetPrimaryParticle() )->
136  GetUserInformation() ) )->SetEMCalEnergy( Edep );
137  // The (initial) energy of the particle is deposited in the step
138  // (which corresponds to the entrance of the electromagnetic calorimeter)
139  aFastStep.ProposeTotalEnergyDeposited( Edep );
140  }
141  }
142 }
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
ush Pos
Definition: deflate.h:89
const G4ThreeVector & GetPosition() const
void FillHistogram(G4int HNo, G4double value) const
Definition: Par02Output.cc:219
G4double GetKineticEnergy() const
void ProposePrimaryTrackPathLength(G4double)
static Par02Output * Instance()
Definition: Par02Output.cc:60
G4double GetResolution(Detector aDetector, Parametrisation aParametrisation, G4double aMomentum)
G4PrimaryParticle * GetPrimaryParticle() const
G4ThreeVector GetMomentum() const
static G4EventManager * GetEventManager()
void ProposeTotalEnergyDeposited(G4double anEnergyPart)
G4double GetEfficiency(Detector aDetector, Parametrisation aParametrisation, G4double aMomentum)
static constexpr double MeV
Definition: G4SIunits.hh:214
void KillPrimaryTrack()
Definition: G4FastStep.cc:88
double G4double
Definition: G4Types.hh:76
double mag() const
G4bool GetDoSmearing()
Gets the flag indicating if smearing should be done.
G4VUserEventInformation * GetUserInformation()

Here is the call graph for this function:

G4bool Par02FastSimModelEMCal::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 75 of file Par02FastSimModelEMCal.cc.

76  {
77  // Applicable for electrons, positrons, and gammas
78  return &aParticleType == G4Electron::Definition() ||
79  &aParticleType == G4Positron::Definition() ||
80  &aParticleType == G4Gamma::Definition();
81 }
static G4Electron * Definition()
Definition: G4Electron.cc:49
static G4Positron * Definition()
Definition: G4Positron.cc:49
static G4Gamma * Definition()
Definition: G4Gamma.cc:49

Here is the call graph for this function:

G4bool Par02FastSimModelEMCal::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 85 of file Par02FastSimModelEMCal.cc.

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

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