Garfield++ 5.0
A toolkit for the detailed simulation of particle detectors based on ionisation measurement in gases and semiconductors
Loading...
Searching...
No Matches
GarfieldG4FastSimulationModel Class Reference

#include <GarfieldG4FastSimulationModel.hh>

+ Inheritance diagram for GarfieldG4FastSimulationModel:

Public Member Functions

 GarfieldG4FastSimulationModel (G4String, G4Region *)
 
 GarfieldG4FastSimulationModel (G4String)
 
 ~GarfieldG4FastSimulationModel ()
 
void SetPhysics (GarfieldPhysics *fGarfieldPhysics)
 
void WriteGeometryToGDML (G4VPhysicalVolume *physicalVolume)
 
virtual G4bool IsApplicable (const G4ParticleDefinition &)
 
virtual G4bool ModelTrigger (const G4FastTrack &)
 
virtual void DoIt (const G4FastTrack &, G4FastStep &)
 

Detailed Description

Definition at line 39 of file GarfieldG4FastSimulationModel.hh.

Constructor & Destructor Documentation

◆ GarfieldG4FastSimulationModel() [1/2]

GarfieldG4FastSimulationModel::GarfieldG4FastSimulationModel ( G4String modelName,
G4Region * envelope )

Definition at line 39 of file GarfieldG4FastSimulationModel.cc.

41 : G4VFastSimulationModel(modelName, envelope) {
42 fGarfieldPhysics = GarfieldPhysics::GetInstance();
43 fGarfieldPhysics->InitializePhysics();
44}
static GarfieldPhysics * GetInstance()

◆ GarfieldG4FastSimulationModel() [2/2]

GarfieldG4FastSimulationModel::GarfieldG4FastSimulationModel ( G4String modelName)

Definition at line 46 of file GarfieldG4FastSimulationModel.cc.

47 : G4VFastSimulationModel(modelName) {
48 fGarfieldPhysics = GarfieldPhysics::GetInstance();
49 fGarfieldPhysics->InitializePhysics();
50}

◆ ~GarfieldG4FastSimulationModel()

GarfieldG4FastSimulationModel::~GarfieldG4FastSimulationModel ( )

Definition at line 52 of file GarfieldG4FastSimulationModel.cc.

52{}

Member Function Documentation

◆ DoIt()

void GarfieldG4FastSimulationModel::DoIt ( const G4FastTrack & fastTrack,
G4FastStep & fastStep )
virtual

Definition at line 83 of file GarfieldG4FastSimulationModel.cc.

84 {
85
86 G4ThreeVector localdir = fastTrack.GetPrimaryTrackLocalDirection();
87 G4ThreeVector localpos = fastTrack.GetPrimaryTrackLocalPosition();
88
89 double ekin_MeV = fastTrack.GetPrimaryTrack()->GetKineticEnergy() / MeV;
90 double globalTime = fastTrack.GetPrimaryTrack()->GetGlobalTime();
91
92 G4String particleName =
93 fastTrack.GetPrimaryTrack()->GetParticleDefinition()->GetParticleName();
94
95 fastStep.KillPrimaryTrack();
96 fastStep.SetPrimaryTrackPathLength(0.0);
97
98 if (particleName == "kaon+") {
99 particleName = "K+";
100 } else if (particleName == "kaon-") {
101 particleName = "K-";
102 } else if (particleName == "anti_proton") {
103 particleName = "anti-proton";
104 }
105
106 fGarfieldPhysics->DoIt(
107 particleName, ekin_MeV, globalTime, localpos.x() / CLHEP::cm,
108 localpos.y() / CLHEP::cm, localpos.z() / CLHEP::cm,
109 localdir.x(), localdir.y(), localdir.z());
110
111 fastStep.SetTotalEnergyDeposited(fGarfieldPhysics->GetEnergyDeposit_MeV());
112
113 if (!fGarfieldPhysics->GetCreateSecondariesInGeant4()) return;
114 const auto& secondaryParticles = fGarfieldPhysics->GetSecondaryParticles();
115
116 if (secondaryParticles.empty()) return;
117 fastStep.SetNumberOfSecondaryTracks(secondaryParticles.size());
118
119 G4double totalEnergySecondaries_MeV = 0.;
120
121 for (const auto& sp : secondaryParticles) {
122 G4double eKin_MeV = sp.getEkin_MeV();
123 G4double time = sp.getTime();
124 G4ThreeVector momentumDirection(sp.getDX(), sp.getDY(), sp.getDZ());
125 G4ThreeVector position(sp.getX_mm(), sp.getY_mm(), sp.getZ_mm());
126 if (sp.getParticleName() == "e-") {
127 G4DynamicParticle particle(G4Electron::ElectronDefinition(),
128 momentumDirection, eKin_MeV);
129 fastStep.CreateSecondaryTrack(particle, position, time, true);
130 totalEnergySecondaries_MeV += eKin_MeV;
131 } else if (sp.getParticleName() == "gamma") {
132 G4DynamicParticle particle(G4Gamma::GammaDefinition(),
133 momentumDirection, eKin_MeV);
134 fastStep.CreateSecondaryTrack(particle, position, time, true);
135 totalEnergySecondaries_MeV += eKin_MeV;
136 }
137 }
138
139}

◆ IsApplicable()

G4bool GarfieldG4FastSimulationModel::IsApplicable ( const G4ParticleDefinition & particleType)
virtual

Definition at line 62 of file GarfieldG4FastSimulationModel.cc.

63 {
64 G4String particleName = particleType.GetParticleName();
65 if (fGarfieldPhysics->FindParticleName(particleName, "garfield")) {
66 return true;
67 }
68 return false;
69}

◆ ModelTrigger()

G4bool GarfieldG4FastSimulationModel::ModelTrigger ( const G4FastTrack & fastTrack)
virtual

Definition at line 71 of file GarfieldG4FastSimulationModel.cc.

72 {
73 double ekin_MeV = fastTrack.GetPrimaryTrack()->GetKineticEnergy() / MeV;
74 G4String particleName =
75 fastTrack.GetPrimaryTrack()->GetParticleDefinition()->GetParticleName();
76 if (fGarfieldPhysics->FindParticleNameEnergy(particleName, ekin_MeV,
77 "garfield")) {
78 return true;
79 }
80 return false;
81}

◆ SetPhysics()

void GarfieldG4FastSimulationModel::SetPhysics ( GarfieldPhysics * fGarfieldPhysics)

◆ WriteGeometryToGDML()

void GarfieldG4FastSimulationModel::WriteGeometryToGDML ( G4VPhysicalVolume * physicalVolume)

Definition at line 54 of file GarfieldG4FastSimulationModel.cc.

55 {
56 G4GDMLParser* parser = new G4GDMLParser();
57 remove("garfieldGeometry.gdml");
58 parser->Write("garfieldGeometry.gdml", physicalVolume, false);
59 delete parser;
60}

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