Garfield++ 5.0
A toolkit for the detailed simulation of particle detectors based on ionisation measurement in gases and semiconductors
Loading...
Searching...
No Matches
Public Member Functions | Static Public Member Functions | List of all members
GarfieldPhysics Class Reference

#include <GarfieldPhysics.hh>

Public Member Functions

void InitializePhysics ()
 
void DoIt (std::string particleName, double ekin_MeV, double time, double x_cm, double y_cm, double z_cm, double dx, double dy, double dz)
 
void AddParticleName (const std::string particleName, double ekin_min_MeV, double ekin_max_MeV, std::string program)
 
bool FindParticleName (const std::string name, std::string program="garfield")
 
bool FindParticleNameEnergy (std::string name, double ekin_MeV, std::string program="garfield")
 
double GetMinEnergyMeVParticle (std::string name, std::string program="garfield")
 
double GetMaxEnergyMeVParticle (std::string name, std::string program="garfield")
 
void SetIonizationModel (std::string model, bool useDefaults=true)
 
std::string GetIonizationModel ()
 
const std::vector< GarfieldParticle > & GetSecondaryParticles () const
 
void EnableCreateSecondariesInGeant4 (bool flag)
 
bool GetCreateSecondariesInGeant4 () const
 
double GetEnergyDeposit_MeV () const
 
double GetAvalancheSize () const
 
double GetGain () const
 
void Clear ()
 

Static Public Member Functions

static GarfieldPhysics * GetInstance ()
 
static void Dispose ()
 

Detailed Description

Definition at line 76 of file GarfieldPhysics.hh.

Member Function Documentation

◆ AddParticleName()

void GarfieldPhysics::AddParticleName ( const std::string particleName,
double ekin_min_MeV,
double ekin_max_MeV,
std::string program )

Definition at line 108 of file GarfieldPhysics.cc.

110 {
111 if (ekin_min_MeV >= ekin_max_MeV) {
112 std::cout << "Ekin_min=" << ekin_min_MeV
113 << " keV is larger than Ekin_max=" << ekin_max_MeV << " keV"
114 << std::endl;
115 return;
116 }
117
118 if (program == "garfield") {
119 std::cout << "Garfield model (Heed) is applicable for G4Particle "
120 << particleName << " between " << ekin_min_MeV << " MeV and "
121 << ekin_max_MeV << " MeV" << std::endl;
122
123 fMapParticlesEnergyGarfield.insert(std::make_pair(
124 particleName, std::make_pair(ekin_min_MeV, ekin_max_MeV)));
125 } else {
126 std::cout << fIonizationModel << " is applicable for G4Particle "
127 << particleName << " between " << ekin_min_MeV << " MeV and "
128 << ekin_max_MeV << " MeV" << std::endl;
129 fMapParticlesEnergyGeant4.insert(std::make_pair(
130 particleName, std::make_pair(ekin_min_MeV, ekin_max_MeV)));
131 }
132}

Referenced by SetIonizationModel(), and GarfieldMessenger::SetNewValue().

◆ Clear()

void GarfieldPhysics::Clear ( )
inline

Definition at line 110 of file GarfieldPhysics.hh.

110 {
111 fEnergyDeposit = 0;
112 fAvalancheSize = 0;
113 fGain = 0;
114 nsum = 0;
115 }

Referenced by GarfieldEventAction::BeginOfEventAction().

◆ Dispose()

void GarfieldPhysics::Dispose ( )
static

Definition at line 43 of file GarfieldPhysics.cc.

43 {
44 delete fGarfieldPhysics;
45 fGarfieldPhysics = nullptr;
46}

Referenced by main().

◆ DoIt()

void GarfieldPhysics::DoIt ( std::string particleName,
double ekin_MeV,
double time,
double x_cm,
double y_cm,
double z_cm,
double dx,
double dy,
double dz )

Definition at line 238 of file GarfieldPhysics.cc.

240 {
241 G4AnalysisManager* analysisManager = G4AnalysisManager::Instance();
242 fEnergyDeposit = 0;
243 fSecondaryParticles.clear();
244
245 Garfield::AvalancheMC drift;
246 drift.SetSensor(fSensor);
247 drift.SetDistanceSteps(1.e-4);
248 Garfield::AvalancheMicroscopic avalanche;
249 avalanche.SetSensor(fSensor);
250
251 // Wire radius [cm]
252 constexpr double rWire = 25.e-4;
253 // Outer radius of the tube [cm]
254 constexpr double rTube = 1.45;
255 // Half-length of the tube [cm]
256 constexpr double lTube = 10.;
257
258 double eKin_eV = ekin_MeV * 1e+6;
259
260 fEnergyDeposit = 0;
261 if (fIonizationModel != "Heed" || particleName == "gamma") {
262 // Number of electrons produced
263 int nc = 0;
264 if (particleName == "gamma") {
265 fTrackHeed->TransportPhoton(x_cm, y_cm, z_cm, time, eKin_eV, dx, dy, dz,
266 nc);
267 } else {
268 fTrackHeed->TransportDeltaElectron(x_cm, y_cm, z_cm, time, eKin_eV, dx,
269 dy, dz, nc);
270 fEnergyDeposit = eKin_eV;
271 }
272
273 for (int cl = 0; cl < nc; cl++) {
274 double xe, ye, ze, te;
275 double ee, dxe, dye, dze;
276 fTrackHeed->GetElectron(cl, xe, ye, ze, te, ee, dxe, dye, dze);
277 if (fabs(ze) > lTube || sqrt(xe * xe + ye * ye) > rTube) continue;
278 nsum++;
279 if (particleName == "gamma") {
280 fEnergyDeposit += fTrackHeed->GetW();
281 }
282 analysisManager->FillH3(1, ze * 10, xe * 10, ye * 10);
283 if (createSecondariesInGeant4) {
284 double newTime = te;
285 if (newTime < time) {
286 newTime += time;
287 }
288 fSecondaryParticles.emplace_back(GarfieldParticle(
289 "e-", ee, newTime, xe, ye, ze, dxe, dye, dze));
290 }
291
292 drift.DriftElectron(xe, ye, ze, te);
293
294 double xe1, ye1, ze1, te1;
295 double xe2, ye2, ze2, te2;
296
297 int status;
298 drift.GetElectronEndpoint(0, xe1, ye1, ze1, te1, xe2, ye2, ze2, te2,
299 status);
300
301 if (0 < xe2 && xe2 < rWire) {
302 xe2 += 2 * rWire;
303 } else if (0 > xe2 && xe2 > -rWire) {
304 xe2 += -2 * rWire;
305 }
306 if (0 < ye2 && ye2 < rWire) {
307 ye2 += 2 * rWire;
308 } else if (0 > ye2 && ye2 > -rWire) {
309 ye2 += -2 * rWire;
310 }
311
312 double e2 = 0.1;
313 avalanche.AvalancheElectron(xe2, ye2, ze2, te2, e2, 0, 0, 0);
314
315 int ne = 0, ni = 0;
316 avalanche.GetAvalancheSize(ne, ni);
317 fAvalancheSize += ne;
318 }
319 } else {
320 fTrackHeed->SetParticle(particleName);
321 fTrackHeed->SetKineticEnergy(eKin_eV);
322 fTrackHeed->NewTrack(x_cm, y_cm, z_cm, time, dx, dy, dz);
323 for (const auto& cluster : fTrackHeed->GetClusters()) {
324 if (fabs(cluster.z) > lTube) continue;
325 if (sqrt(cluster.x * cluster.x + cluster.y * cluster.y) >= rTube) {
326 continue;
327 }
328 nsum += cluster.electrons.size();
329 fEnergyDeposit += cluster.energy;
330 for (const auto& electron : cluster.electrons) {
331 if (fabs(electron.z) > lTube) continue;
332 if (sqrt(electron.x * electron.x + electron.y * electron.y) >= rTube) {
333 continue;
334 }
335 analysisManager->FillH3(1, electron.z * 10, electron.x * 10, electron.y * 10);
336 if (createSecondariesInGeant4) {
337 double newTime = electron.t;
338 if (newTime < time) {
339 newTime += time;
340 }
341 fSecondaryParticles.emplace_back(GarfieldParticle(
342 "e-", electron.e, newTime, electron.x, electron.y, electron.z,
343 electron.dx, electron.dy, electron.dz));
344 }
345
346 drift.DriftElectron(electron.x, electron.y, electron.z, electron.t);
347
348 double xe1, ye1, ze1, te1;
349 double xe2, ye2, ze2, te2;
350
351 int status;
352 drift.GetElectronEndpoint(0, xe1, ye1, ze1, te1, xe2, ye2, ze2,
353 te2, status);
354
355 if (0 < xe2 && xe2 < rWire) {
356 xe2 += 2 * rWire;
357 } else if (0 > xe2 && xe2 > -rWire) {
358 xe2 += -2 * rWire;
359 }
360 if (0 < ye2 && ye2 < rWire) {
361 ye2 += 2 * rWire;
362 } else if (0 > ye2 && ye2 > -rWire) {
363 ye2 += -2 * rWire;
364 }
365
366 double e2 = 0.1;
367 avalanche.AvalancheElectron(xe2, ye2, ze2, te2, e2, 0, 0, 0);
368
369 int ne = 0, ni = 0;
370 avalanche.GetAvalancheSize(ne, ni);
371 fAvalancheSize += ne;
372 }
373 }
374 }
375 fGain = fAvalancheSize / nsum;
376}
Garfield::AvalancheMC::SetSensor
void SetSensor(Sensor *s)
Set the sensor.
Definition AvalancheMC.cc:68
Garfield::AvalancheMC::GetElectronEndpoint
void GetElectronEndpoint(const size_t i, double &x0, double &y0, double &z0, double &t0, double &x1, double &y1, double &z1, double &t1, int &status) const
Definition AvalancheMC.cc:192
Garfield::AvalancheMC::DriftElectron
bool DriftElectron(const double x, const double y, const double z, const double t)
Simulate the drift line of an electron from a given starting point.
Definition AvalancheMC.cc:212
Garfield::AvalancheMC::SetDistanceSteps
void SetDistanceSteps(const double d=0.001)
Use fixed distance steps (default 10 um).
Definition AvalancheMC.cc:101
Garfield::AvalancheMicroscopic::GetAvalancheSize
void GetAvalancheSize(int &ne, int &ni) const
Return the number of electrons and ions in the avalanche.
Definition AvalancheMicroscopic.hh:138
Garfield::AvalancheMicroscopic::SetSensor
void SetSensor(Sensor *sensor)
Set the sensor.
Definition AvalancheMicroscopic.cc:132
Garfield::AvalancheMicroscopic::AvalancheElectron
bool AvalancheElectron(const double x, const double y, const double z, const double t, const double e, const double dx=0., const double dy=0., const double dz=0.)
Calculate an avalanche initiated by a given electron.
Definition AvalancheMicroscopic.cc:431
Heed::fabs
DoubleAc fabs(const DoubleAc &f)
Definition DoubleAc.h:615
Heed::sqrt
DoubleAc sqrt(const DoubleAc &f)
Definition DoubleAc.cpp:314

◆ EnableCreateSecondariesInGeant4()

void GarfieldPhysics::EnableCreateSecondariesInGeant4 ( bool flag)
inline

Definition at line 101 of file GarfieldPhysics.hh.

101 {
102 createSecondariesInGeant4 = flag;
103 }

◆ FindParticleName()

bool GarfieldPhysics::FindParticleName ( const std::string name,
std::string program = "garfield" )

Definition at line 134 of file GarfieldPhysics.cc.

134 {
135 if (program == "garfield") {
136 auto it = fMapParticlesEnergyGarfield.find(name);
137 if (it != fMapParticlesEnergyGarfield.end()) return true;
138 } else {
139 auto it = fMapParticlesEnergyGeant4.find(name);
140 if (it != fMapParticlesEnergyGeant4.end()) return true;
141 }
142 return false;
143}

Referenced by GarfieldPhysicsList::AddParameterisation().

◆ FindParticleNameEnergy()

bool GarfieldPhysics::FindParticleNameEnergy ( std::string name,
double ekin_MeV,
std::string program = "garfield" )

Definition at line 145 of file GarfieldPhysics.cc.

146 {
147 if (program == "garfield") {
148 auto it = fMapParticlesEnergyGarfield.find(name);
149 if (it != fMapParticlesEnergyGarfield.end()) {
150 EnergyRange_MeV range = it->second;
151 if (range.first <= ekin_MeV && range.second >= ekin_MeV) {
152 return true;
153 }
154 }
155 } else {
156 auto it = fMapParticlesEnergyGeant4.find(name);
157 if (it != fMapParticlesEnergyGeant4.end()) {
158 EnergyRange_MeV range = it->second;
159 if (range.first <= ekin_MeV && range.second >= ekin_MeV) {
160 return true;
161 }
162 }
163 }
164 return false;
165}
EnergyRange_MeV
std::pair< double, double > EnergyRange_MeV
Definition GarfieldPhysics.hh:42

◆ GetAvalancheSize()

double GarfieldPhysics::GetAvalancheSize ( ) const
inline

Definition at line 108 of file GarfieldPhysics.hh.

108{ return fAvalancheSize; }

Referenced by GarfieldEventAction::EndOfEventAction().

◆ GetCreateSecondariesInGeant4()

bool GarfieldPhysics::GetCreateSecondariesInGeant4 ( ) const
inline

Definition at line 104 of file GarfieldPhysics.hh.

104 {
105 return createSecondariesInGeant4;
106 }

◆ GetEnergyDeposit_MeV()

double GarfieldPhysics::GetEnergyDeposit_MeV ( ) const
inline

Definition at line 107 of file GarfieldPhysics.hh.

107{ return fEnergyDeposit / 1000000; }

◆ GetGain()

double GarfieldPhysics::GetGain ( ) const
inline

Definition at line 109 of file GarfieldPhysics.hh.

109{ return fGain; }

Referenced by GarfieldEventAction::EndOfEventAction().

◆ GetInstance()

GarfieldPhysics * GarfieldPhysics::GetInstance ( )
static

Definition at line 38 of file GarfieldPhysics.cc.

38 {
39 if (!fGarfieldPhysics) fGarfieldPhysics = new GarfieldPhysics();
40 return fGarfieldPhysics;
41}

Referenced by GarfieldPhysicsList::AddParameterisation(), GarfieldEventAction::BeginOfEventAction(), GarfieldEventAction::EndOfEventAction(), GarfieldG4FastSimulationModel::GarfieldG4FastSimulationModel(), GarfieldG4FastSimulationModel::GarfieldG4FastSimulationModel(), and GarfieldMessenger::SetNewValue().

◆ GetIonizationModel()

std::string GarfieldPhysics::GetIonizationModel ( )

Definition at line 57 of file GarfieldPhysics.cc.

57{ return fIonizationModel; }

Referenced by GarfieldPhysicsList::AddParameterisation().

◆ GetMaxEnergyMeVParticle()

double GarfieldPhysics::GetMaxEnergyMeVParticle ( std::string name,
std::string program = "garfield" )

Definition at line 185 of file GarfieldPhysics.cc.

186 {
187 if (program == "garfield") {
188 auto it = fMapParticlesEnergyGarfield.find(name);
189 if (it != fMapParticlesEnergyGarfield.end()) {
190 EnergyRange_MeV range = it->second;
191 return range.second;
192 }
193 } else {
194 auto it = fMapParticlesEnergyGeant4.find(name);
195 if (it != fMapParticlesEnergyGeant4.end()) {
196 EnergyRange_MeV range = it->second;
197 return range.second;
198 }
199 }
200 return -1;
201}

Referenced by GarfieldPhysicsList::AddParameterisation().

◆ GetMinEnergyMeVParticle()

double GarfieldPhysics::GetMinEnergyMeVParticle ( std::string name,
std::string program = "garfield" )

Definition at line 167 of file GarfieldPhysics.cc.

168 {
169 if (program == "garfield") {
170 auto it = fMapParticlesEnergyGarfield.find(name);
171 if (it != fMapParticlesEnergyGarfield.end()) {
172 EnergyRange_MeV range = it->second;
173 return range.first;
174 }
175 } else {
176 auto it = fMapParticlesEnergyGeant4.find(name);
177 if (it != fMapParticlesEnergyGeant4.end()) {
178 EnergyRange_MeV range = it->second;
179 return range.first;
180 }
181 }
182 return -1;
183}

Referenced by GarfieldPhysicsList::AddParameterisation().

◆ GetSecondaryParticles()

const std::vector< GarfieldParticle > & GarfieldPhysics::GetSecondaryParticles ( ) const
inline

Definition at line 98 of file GarfieldPhysics.hh.

98 {
99 return fSecondaryParticles;
100 }

◆ InitializePhysics()

void GarfieldPhysics::InitializePhysics ( )

Definition at line 203 of file GarfieldPhysics.cc.

203 {
204 // Define the gas mixture.
205 fMediumMagboltz = new Garfield::MediumMagboltz();
206 fMediumMagboltz->SetComposition("ar", 70., "co2", 30.);
207 fMediumMagboltz->SetTemperature(293.15);
208 fMediumMagboltz->SetPressure(760.);
209 fMediumMagboltz->Initialise(true);
210 // Set the Penning transfer efficiency.
211 const double rPenning = 0.57;
212 const double lambdaPenning = 0.;
213 fMediumMagboltz->EnablePenningTransfer(rPenning, lambdaPenning, "ar");
214 fMediumMagboltz->LoadGasFile("ar_70_co2_30_1000mbar.gas");
215
216 fComponentAnalyticField = new Garfield::ComponentAnalyticField();
217 fComponentAnalyticField->SetMedium(fMediumMagboltz);
218 // Wire radius [cm]
219 constexpr double rWire = 25.e-4;
220 // Tube radius [cm]
221 constexpr double rTube = 1.451;
222 // Voltages
223 constexpr double vWire = 1000.;
224 constexpr double vTube = 0.;
225 // Add the wire in the center.
226 fComponentAnalyticField->AddWire(0., 0., 2 * rWire, vWire, "w");
227 // Add the tube.
228 fComponentAnalyticField->AddTube(rTube, vTube, 0, "t");
229
230 fSensor = new Garfield::Sensor();
231 fSensor->AddComponent(fComponentAnalyticField);
232
233 fTrackHeed = new Garfield::TrackHeed();
234 fTrackHeed->SetSensor(fSensor);
235 fTrackHeed->EnableDeltaElectronTransport();
236}

◆ SetIonizationModel()

void GarfieldPhysics::SetIonizationModel ( std::string model,
bool useDefaults = true )

Definition at line 59 of file GarfieldPhysics.cc.

59 {
60 if (model != "PAIPhot" && model != "PAI" && model != "Heed") {
61 std::cout << "Unknown ionization model " << model << std::endl;
62 std::cout << "Using PAIPhot as default model!" << std::endl;
63 model = "PAIPhot";
64 }
65 fIonizationModel = model;
66
67 if (fIonizationModel == "PAIPhot" || fIonizationModel == "PAI") {
68 if (useDefaults) {
69 // Particle types and energies for which the G4FastSimulationModel with
70 // Garfield++ is valid
71 this->AddParticleName("e-", 1e-6, 1e-3, "garfield");
72 this->AddParticleName("gamma", 1e-6, 1e+8, "garfield");
73
74 // Particle types and energies for which the PAI or PAIPhot model is valid
75 this->AddParticleName("e-", 0, 1e+8, "geant4");
76 this->AddParticleName("e+", 0, 1e+8, "geant4");
77 this->AddParticleName("mu-", 0, 1e+8, "geant4");
78 this->AddParticleName("mu+", 0, 1e+8, "geant4");
79 this->AddParticleName("proton", 0, 1e+8, "geant4");
80 this->AddParticleName("pi+", 0, 1e+8, "geant4");
81 this->AddParticleName("pi-", 0, 1e+8, "geant4");
82 this->AddParticleName("alpha", 0, 1e+8, "geant4");
83 this->AddParticleName("He3", 0, 1e+8, "geant4");
84 this->AddParticleName("GenericIon", 0, 1e+8, "geant4");
85 }
86
87 } else if (fIonizationModel == "Heed") {
88 if (useDefaults) {
89 // Particle types and energies for which the G4FastSimulationModel with
90 // Garfield++ is valid
91 this->AddParticleName("gamma", 1e-6, 1e+8, "garfield");
92 this->AddParticleName("e-", 6e-2, 1e+7, "garfield");
93 this->AddParticleName("e+", 6e-2, 1e+7, "garfield");
94 this->AddParticleName("mu-", 1e+1, 1e+8, "garfield");
95 this->AddParticleName("mu+", 1e+1, 1e+8, "garfield");
96 this->AddParticleName("pi-", 2e+1, 1e+8, "garfield");
97 this->AddParticleName("pi+", 2e+1, 1e+8, "garfield");
98 this->AddParticleName("kaon-", 1e+1, 1e+8, "garfield");
99 this->AddParticleName("kaon+", 1e+1, 1e+8, "garfield");
100 this->AddParticleName("proton", 9.e+1, 1e+8, "garfield");
101 this->AddParticleName("anti_proton", 9.e+1, 1e+8, "garfield");
102 this->AddParticleName("deuteron", 2.e+2, 1e+8, "garfield");
103 this->AddParticleName("alpha", 4.e+2, 1e+8, "garfield");
104 }
105 }
106}
GarfieldPhysics::AddParticleName
void AddParticleName(const std::string particleName, double ekin_min_MeV, double ekin_max_MeV, std::string program)
Definition GarfieldPhysics.cc:108

Referenced by GarfieldMessenger::SetNewValue().

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