Geant4-11
Public Member Functions | Protected Attributes | Private Attributes
RunAction Class Reference

#include <RunAction.hh>

Inheritance diagram for RunAction:
G4UserRunAction

Public Member Functions

void BeginOfRunAction (const G4Run *)
 
void CriticalEnergy ()
 
void EndOfRunAction (const G4Run *)
 
virtual G4RunGenerateRun ()
 
void GetCuts ()
 
G4bool IsMaster () const
 
 RunAction (DetectorConstruction *, PrimaryGeneratorAction *)
 
virtual void SetMaster (G4bool val=true)
 
 ~RunAction ()
 

Protected Attributes

G4bool isMaster = true
 

Private Attributes

DetectorConstructiondetector
 
G4double energyCut [3]
 
PrimaryGeneratorActionprimary
 
G4double rangeCut [3]
 

Detailed Description

Definition at line 44 of file RunAction.hh.

Constructor & Destructor Documentation

◆ RunAction()

RunAction::RunAction ( DetectorConstruction det,
PrimaryGeneratorAction kin 
)

Definition at line 45 of file RunAction.cc.

46:detector(det), primary(kin)
47{ }
PrimaryGeneratorAction * primary
Definition: RunAction.hh:59
DetectorConstruction * detector
Definition: RunAction.hh:58

◆ ~RunAction()

RunAction::~RunAction ( )

Definition at line 51 of file RunAction.cc.

52{ }

Member Function Documentation

◆ BeginOfRunAction()

void RunAction::BeginOfRunAction ( const G4Run )
virtual

Reimplemented from G4UserRunAction.

Definition at line 56 of file RunAction.cc.

57{
58 //set precision for printing
59 G4int prec = G4cout.precision(6);
60
61 // get particle
64 G4String partName = particle->GetParticleName();
65 G4double charge = particle->GetPDGCharge();
67
68 // get material
70 G4String matName = material->GetName();
71 G4double density = material->GetDensity();
72 G4double radl = material->GetRadlen();
73
74 G4cout << "\n " << partName << " ("
75 << G4BestUnit(energy,"Energy") << ") in "
76 << material->GetName() << " (density: "
77 << G4BestUnit(density,"Volumic Mass") << "; radiation length: "
78 << G4BestUnit(radl, "Length") << ")" << G4endl;
79
80 // get cuts
81 GetCuts();
82 if (charge != 0.) {
83 G4cout << "\n Range cuts : \t gamma "
84 << std::setw(8) << G4BestUnit(rangeCut[0],"Length")
85 << "\t e- " << std::setw(8) << G4BestUnit(rangeCut[1],"Length");
86 G4cout << "\n Energy cuts : \t gamma "
87 << std::setw(8) << G4BestUnit(energyCut[0],"Energy")
88 << "\t e- " << std::setw(8) << G4BestUnit(energyCut[1],"Energy")
89 << G4endl;
90 }
91
92 // get processList and extract EM processes (but not MultipleScattering)
93 G4ProcessVector* plist = particle->GetProcessManager()->GetProcessList();
94 G4String procName;
95 G4double cut;
96 std::vector<G4String> emName;
97 std::vector<G4double> enerCut;
98 size_t length = plist->size();
99 for (size_t j=0; j<length; j++) {
100 procName = (*plist)[j]->GetProcessName();
101 cut = energyCut[1];
102 if ((procName == "eBrem")||(procName == "muBrems")) cut = energyCut[0];
103 if (((*plist)[j]->GetProcessType() == fElectromagnetic) &&
104 (procName != "msc")) {
105 emName.push_back(procName);
106 enerCut.push_back(cut);
107 }
108 }
109
110 // print list of processes
111 G4cout << "\n processes : ";
112 for (size_t j=0; j<emName.size();j++)
113 G4cout << "\t" << std::setw(13) << emName[j] << "\t";
114 G4cout << "\t" << std::setw(13) <<"total";
115
116 //instanciate EmCalculator
117 G4EmCalculator emCal;
118 // emCal.SetVerbose(2);
119
120 //compute cross section per atom (only for single material)
121 if (material->GetNumberOfElements() == 1) {
122 G4double Z = material->GetZ();
123 G4double A = material->GetA();
124
125 std::vector<G4double> sigma0;
126 G4double sig, sigtot = 0.;
127
128 for (size_t j=0; j<emName.size();j++) {
130 (energy,particle,emName[j],Z,A,enerCut[j]);
131 sigtot += sig;
132 sigma0.push_back(sig);
133 }
134 sigma0.push_back(sigtot);
135
136 G4cout << "\n \n cross section per atom : ";
137 for (size_t j=0; j<sigma0.size();j++) {
138 G4cout << "\t" << std::setw(13) << G4BestUnit(sigma0[j], "Surface");
139 }
140 G4cout << G4endl;
141 }
142
143 //get cross section per volume
144 std::vector<G4double> sigma1;
145 std::vector<G4double> sigma2;
146 G4double Sig, Sigtot = 0.;
147
148 for (size_t j=0; j<emName.size();j++) {
149 Sig = emCal.GetCrossSectionPerVolume(energy,particle,emName[j],material);
150 if (Sig == 0.) Sig = emCal.ComputeCrossSectionPerVolume
151 (energy,particle,emName[j],material,enerCut[j]);
152 Sigtot += Sig;
153 sigma1.push_back(Sig);
154 sigma2.push_back(Sig/density);
155 }
156 sigma1.push_back(Sigtot);
157 sigma2.push_back(Sigtot/density);
158
159 //print cross sections
160 G4cout << "\n \n cross section per volume : ";
161 for (size_t j=0; j<sigma1.size();j++) {
162 G4cout << "\t" << std::setw(13) << sigma1[j]*cm << " cm^-1";
163 }
164
165 G4cout << "\n cross section per mass : ";
166 for (size_t j=0; j<sigma2.size();j++) {
167 G4cout << "\t" << std::setw(13) << G4BestUnit(sigma2[j], "Surface/Mass");
168 }
169
170 //print mean free path
171
173
174 G4cout << "\n \n mean free path : ";
175 for (size_t j=0; j<sigma1.size();j++) {
176 lambda = DBL_MAX;
177 if (sigma1[j] > 0.) lambda = 1/sigma1[j];
178 G4cout << "\t" << std::setw(13) << G4BestUnit( lambda, "Length");
179 }
180
181 //mean free path (g/cm2)
182 G4cout << "\n (g/cm2) : ";
183 for (size_t j=0; j<sigma2.size();j++) {
184 lambda = DBL_MAX;
185 if (sigma2[j] > 0.) lambda = 1/sigma2[j];
186 G4cout << "\t" << std::setw(13) << G4BestUnit( lambda, "Mass/Surface");
187 }
188 G4cout << G4endl;
189
190 if (charge == 0.) {
191 G4cout.precision(prec);
192 G4cout << "\n-------------------------------------------------------------\n"
193 << G4endl;
194 return;
195 }
196
197 //get stopping power
198 std::vector<G4double> dedx1;
199 std::vector<G4double> dedx2;
200 G4double dedx, dedxtot = 0.;
201
202 for (size_t j=0; j<emName.size();j++) {
203 dedx = emCal.ComputeDEDX(energy,particle,emName[j],material,enerCut[j]);
204 dedx1.push_back(dedx);
205 dedx2.push_back(dedx/density);
206 }
207 dedxtot = emCal.GetDEDX(energy,particle,material);
208 dedx1.push_back(dedxtot);
209 dedx2.push_back(dedxtot/density);
210
211 //print stopping power
212 G4cout << "\n \n restricted dE/dx : ";
213 for (size_t j=0; j<sigma1.size();j++) {
214 G4cout << "\t" << std::setw(13) << G4BestUnit(dedx1[j],"Energy/Length");
215 }
216
217 G4cout << "\n (MeV/g/cm2) : ";
218 for (size_t j=0; j<sigma2.size();j++) {
219 G4cout << "\t" << std::setw(13) << G4BestUnit(dedx2[j],"Energy*Surface/Mass");
220 }
221
222 //get range from restricted dedx
223 G4double range1 = emCal.GetRangeFromRestricteDEDX(energy,particle,material);
224 G4double range2 = range1*density;
225
226 //get range from full dedx
227 G4double Range1 = emCal.GetCSDARange(energy,particle,material);
228 G4double Range2 = Range1*density;
229
230 //print range
231 G4cout << "\n \n range from restrict dE/dx: "
232 << "\t" << std::setw(8) << G4BestUnit(range1,"Length")
233 << " (" << std::setw(8) << G4BestUnit(range2,"Mass/Surface") << ")";
234
235 G4cout << "\n range from full dE/dx : "
236 << "\t" << std::setw(8) << G4BestUnit(Range1,"Length")
237 << " (" << std::setw(8) << G4BestUnit(Range2,"Mass/Surface") << ")";
238
239 //get transport mean free path (for multiple scattering)
240 G4double MSmfp1 = emCal.GetMeanFreePath(energy,particle,"msc",material);
241 G4double MSmfp2 = MSmfp1*density;
242
243 //print transport mean free path
244 G4cout << "\n \n transport mean free path : "
245 << "\t" << std::setw(8) << G4BestUnit(MSmfp1,"Length")
246 << " (" << std::setw(8) << G4BestUnit(MSmfp2,"Mass/Surface") << ")";
247
248 if (particle == G4Electron::Electron()) CriticalEnergy();
249
250 G4cout << "\n-------------------------------------------------------------\n";
251 G4cout << G4endl;
252
253 // reset default precision
254 G4cout.precision(prec);
255}
@ fElectromagnetic
static constexpr double cm
Definition: G4SIunits.hh:99
#define G4BestUnit(a, b)
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85
const G4int Z[17]
const G4double A[17]
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
static G4Electron * Electron()
Definition: G4Electron.cc:93
G4double GetDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, const G4Region *r=nullptr)
G4double GetCSDARange(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, const G4Region *r=nullptr)
G4double ComputeDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, G4double cut=DBL_MAX)
G4double GetRangeFromRestricteDEDX(G4double kinEnergy, const G4ParticleDefinition *, const G4Material *, const G4Region *r=nullptr)
G4double ComputeCrossSectionPerVolume(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, G4double cut=0.0)
G4double ComputeCrossSectionPerAtom(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, G4double Z, G4double A, G4double cut=0.0)
G4double GetCrossSectionPerVolume(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, const G4Region *r=nullptr)
G4double GetMeanFreePath(G4double kinEnergy, const G4ParticleDefinition *, const G4String &processName, const G4Material *, const G4Region *r=nullptr)
G4ProcessManager * GetProcessManager() const
G4double GetPDGCharge() const
const G4String & GetParticleName() const
G4ParticleDefinition * GetParticleDefinition() const
G4double GetParticleEnergy() const
G4ProcessVector * GetProcessList() const
std::size_t size() const
G4ParticleGun * GetParticleGun()
G4double rangeCut[3]
Definition: RunAction.hh:60
void GetCuts()
Definition: RunAction.cc:268
void CriticalEnergy()
Definition: RunAction.cc:299
G4double energyCut[3]
Definition: RunAction.hh:61
static const double prec
Definition: RanecuEngine.cc:61
G4double energy(const ThreeVector &p, const G4double m)
string material
Definition: eplot.py:19
#define DBL_MAX
Definition: templates.hh:62

References A, cm, G4EmCalculator::ComputeCrossSectionPerAtom(), G4EmCalculator::ComputeCrossSectionPerVolume(), G4EmCalculator::ComputeDEDX(), CriticalEnergy(), DBL_MAX, detector, G4Electron::Electron(), G4INCL::KinematicsUtils::energy(), energyCut, fElectromagnetic, G4BestUnit, G4cout, G4endl, G4EmCalculator::GetCrossSectionPerVolume(), G4EmCalculator::GetCSDARange(), GetCuts(), G4EmCalculator::GetDEDX(), DetectorConstruction::GetMaterial(), G4EmCalculator::GetMeanFreePath(), G4ParticleGun::GetParticleDefinition(), G4ParticleGun::GetParticleEnergy(), PrimaryGeneratorAction::GetParticleGun(), G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetPDGCharge(), G4ProcessManager::GetProcessList(), G4ParticleDefinition::GetProcessManager(), G4EmCalculator::GetRangeFromRestricteDEDX(), G4InuclParticleNames::lambda, eplot::material, CLHEP::prec, primary, rangeCut, G4ProcessVector::size(), and Z.

◆ CriticalEnergy()

void RunAction::CriticalEnergy ( )

Definition at line 299 of file RunAction.cc.

300{
301 // compute e- critical energy (Rossi definition) and Moliere radius.
302 // Review of Particle Physics - Eur. Phys. J. C3 (1998) page 147
303 //
304 G4EmCalculator emCal;
305
307 const G4double radl = material->GetRadlen();
308 G4double ekin = 5*MeV;
309 G4double deioni;
310 G4double err = 1., errmax = 0.001;
311 G4int iter = 0 , itermax = 10;
312 while (err > errmax && iter < itermax) {
313 iter++;
314 deioni = radl*
315 emCal.ComputeDEDX(ekin,G4Electron::Electron(),"eIoni",material);
316 err = std::abs(deioni - ekin)/ekin;
317 ekin = deioni;
318 }
319 G4cout << "\n \n critical energy (Rossi) : "
320 << "\t" << std::setw(8) << G4BestUnit(ekin,"Energy");
321
322 //Pdg formula (only for single material)
323 G4double pdga[2] = { 610*MeV, 710*MeV };
324 G4double pdgb[2] = { 1.24, 0.92 };
325 G4double EcPdg = 0.;
326
327 if (material->GetNumberOfElements() == 1) {
328 G4int istat = 0;
329 if (material->GetState() == kStateGas) istat = 1;
330 G4double Zeff = material->GetZ() + pdgb[istat];
331 EcPdg = pdga[istat]/Zeff;
332 G4cout << "\t\t\t (from Pdg formula : "
333 << std::setw(8) << G4BestUnit(EcPdg,"Energy") << ")";
334 }
335
336 const G4double Es = 21.2052*MeV;
337 G4double rMolier1 = Es/ekin, rMolier2 = rMolier1*radl;
338 G4cout << "\n Moliere radius : "
339 << "\t" << std::setw(8) << rMolier1 << " X0 "
340 << "= " << std::setw(8) << G4BestUnit(rMolier2,"Length");
341
342 if (material->GetNumberOfElements() == 1) {
343 G4double rMPdg = radl*Es/EcPdg;
344 G4cout << "\t (from Pdg formula : "
345 << std::setw(8) << G4BestUnit(rMPdg,"Length") << ")";
346 }
347}
@ kStateGas
Definition: G4Material.hh:111
static constexpr double MeV
Definition: G4SIunits.hh:200

References G4EmCalculator::ComputeDEDX(), detector, G4Electron::Electron(), G4BestUnit, G4cout, DetectorConstruction::GetMaterial(), kStateGas, eplot::material, and MeV.

Referenced by BeginOfRunAction().

◆ EndOfRunAction()

void RunAction::EndOfRunAction ( const G4Run )
virtual

Reimplemented from G4UserRunAction.

Definition at line 259 of file RunAction.cc.

260{ }

◆ GenerateRun()

G4Run * G4UserRunAction::GenerateRun ( )
virtualinherited

◆ GetCuts()

void RunAction::GetCuts ( )

Definition at line 268 of file RunAction.cc.

269{
270 G4ProductionCutsTable* theCoupleTable =
272
273 size_t numOfCouples = theCoupleTable->GetTableSize();
274 const G4MaterialCutsCouple* couple = 0;
275 G4int index = 0;
276 for (size_t i=0; i<numOfCouples; i++) {
277 couple = theCoupleTable->GetMaterialCutsCouple(i);
278 if (couple->GetMaterial() == detector->GetMaterial()) {index = i; break;}
279 }
280
281 rangeCut[0] =
282 (*(theCoupleTable->GetRangeCutsVector(idxG4GammaCut)))[index];
283 rangeCut[1] =
284 (*(theCoupleTable->GetRangeCutsVector(idxG4ElectronCut)))[index];
285 rangeCut[2] =
286 (*(theCoupleTable->GetRangeCutsVector(idxG4PositronCut)))[index];
287
288 energyCut[0] =
289 (*(theCoupleTable->GetEnergyCutsVector(idxG4GammaCut)))[index];
290 energyCut[1] =
291 (*(theCoupleTable->GetEnergyCutsVector(idxG4ElectronCut)))[index];
292 energyCut[2] =
293 (*(theCoupleTable->GetEnergyCutsVector(idxG4PositronCut)))[index];
294
295}
@ idxG4ElectronCut
@ idxG4GammaCut
@ idxG4PositronCut
const G4Material * GetMaterial() const
const std::vector< G4double > * GetRangeCutsVector(std::size_t pcIdx) const
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
std::size_t GetTableSize() const
const std::vector< G4double > * GetEnergyCutsVector(std::size_t pcIdx) const
static G4ProductionCutsTable * GetProductionCutsTable()

References detector, energyCut, G4ProductionCutsTable::GetEnergyCutsVector(), DetectorConstruction::GetMaterial(), G4MaterialCutsCouple::GetMaterial(), G4ProductionCutsTable::GetMaterialCutsCouple(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetRangeCutsVector(), G4ProductionCutsTable::GetTableSize(), idxG4ElectronCut, idxG4GammaCut, idxG4PositronCut, and rangeCut.

Referenced by BeginOfRunAction().

◆ IsMaster()

G4bool G4UserRunAction::IsMaster ( ) const
inlineinherited

Definition at line 64 of file G4UserRunAction.hh.

64{ return isMaster; }

References G4UserRunAction::isMaster.

◆ SetMaster()

virtual void G4UserRunAction::SetMaster ( G4bool  val = true)
inlinevirtualinherited

Field Documentation

◆ detector

DetectorConstruction* RunAction::detector
private

Definition at line 58 of file RunAction.hh.

Referenced by BeginOfRunAction(), CriticalEnergy(), and GetCuts().

◆ energyCut

G4double RunAction::energyCut[3]
private

Definition at line 61 of file RunAction.hh.

Referenced by BeginOfRunAction(), and GetCuts().

◆ isMaster

G4bool G4UserRunAction::isMaster = true
protectedinherited

Definition at line 68 of file G4UserRunAction.hh.

Referenced by G4UserRunAction::IsMaster(), and G4UserRunAction::SetMaster().

◆ primary

PrimaryGeneratorAction* RunAction::primary
private

Definition at line 59 of file RunAction.hh.

Referenced by BeginOfRunAction().

◆ rangeCut

G4double RunAction::rangeCut[3]
private

Definition at line 60 of file RunAction.hh.

Referenced by BeginOfRunAction(), and GetCuts().


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