ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4EmPenelopePhysics.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4EmPenelopePhysics.cc
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 
27 #include "G4EmPenelopePhysics.hh"
28 #include "G4ParticleDefinition.hh"
29 #include "G4SystemOfUnits.hh"
30 #include "G4ParticleTable.hh"
31 
32 // *** Processes and models
33 
34 // gamma
35 
36 #include "G4PhotoElectricEffect.hh"
37 #include "G4PenelopePhotoElectricModel.hh"
38 
39 #include "G4ComptonScattering.hh"
40 #include "G4PenelopeComptonModel.hh"
41 
42 #include "G4GammaConversion.hh"
43 #include "G4PenelopeGammaConversionModel.hh"
44 
45 #include "G4RayleighScattering.hh"
46 #include "G4PenelopeRayleighModel.hh"
47 
48 #include "G4PEEffectFluoModel.hh"
49 #include "G4KleinNishinaModel.hh"
50 
51 // e- and e+
52 
53 #include "G4eMultipleScattering.hh"
54 #include "G4UniversalFluctuation.hh"
55 
56 #include "G4eIonisation.hh"
57 #include "G4PenelopeIonisationModel.hh"
58 
59 #include "G4eBremsstrahlung.hh"
60 #include "G4PenelopeBremsstrahlungModel.hh"
61 
62 // e+ only
63 
64 #include "G4eplusAnnihilation.hh"
65 #include "G4PenelopeAnnihilationModel.hh"
66 
67 // mu
68 
69 #include "G4MuMultipleScattering.hh"
70 #include "G4MuIonisation.hh"
71 #include "G4MuBremsstrahlung.hh"
72 #include "G4MuPairProduction.hh"
73 
74 #include "G4MuBremsstrahlungModel.hh"
75 #include "G4MuPairProductionModel.hh"
76 #include "G4hBremsstrahlungModel.hh"
77 #include "G4hPairProductionModel.hh"
78 
79 // hadrons
80 
81 #include "G4hMultipleScattering.hh"
82 #include "G4MscStepLimitType.hh"
83 
84 #include "G4hBremsstrahlung.hh"
85 #include "G4hPairProduction.hh"
86 #include "G4ePairProduction.hh"
87 
88 #include "G4hIonisation.hh"
89 #include "G4ionIonisation.hh"
90 #include "G4alphaIonisation.hh"
91 #include "G4IonParametrisedLossModel.hh"
92 #include "G4NuclearStopping.hh"
93 
94 // msc models
95 #include "G4UrbanMscModel.hh"
96 #include "G4GoudsmitSaundersonMscModel.hh"
97 #include "G4WentzelVIModel.hh"
98 #include "G4CoulombScattering.hh"
99 #include "G4eCoulombScatteringModel.hh"
100 //
101 
102 #include "G4LossTableManager.hh"
103 #include "G4VAtomDeexcitation.hh"
104 #include "G4UAtomicDeexcitation.hh"
105 
106 // particles
107 
108 #include "G4Gamma.hh"
109 #include "G4Electron.hh"
110 #include "G4Positron.hh"
111 #include "G4MuonPlus.hh"
112 #include "G4MuonMinus.hh"
113 #include "G4PionPlus.hh"
114 #include "G4PionMinus.hh"
115 #include "G4KaonPlus.hh"
116 #include "G4KaonMinus.hh"
117 #include "G4Proton.hh"
118 #include "G4AntiProton.hh"
119 #include "G4Deuteron.hh"
120 #include "G4Triton.hh"
121 #include "G4He3.hh"
122 #include "G4Alpha.hh"
123 #include "G4GenericIon.hh"
124 
125 //
126 #include "G4PhysicsListHelper.hh"
127 #include "G4BuilderType.hh"
128 #include "G4EmModelActivator.hh"
129 
130 // factory
131 #include "G4PhysicsConstructorFactory.hh"
132 //
133 G4_DECLARE_PHYSCONSTR_FACTORY(G4EmPenelopePhysics);
134 
135 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
136 
137 G4EmPenelopePhysics::G4EmPenelopePhysics(G4int ver, const G4String&)
138  : G4VPhysicsConstructor("G4EmPenelope"), verbose(ver)
139 {
140  G4EmParameters* param = G4EmParameters::Instance();
141  param->SetDefaults();
142  param->SetVerbose(verbose);
143  param->SetMinEnergy(100*eV);
144  param->SetLowestElectronEnergy(100*eV);
145  param->SetNumberOfBinsPerDecade(20);
146  param->SetStepFunction(0.2, 10*um);
147  param->SetStepFunctionMuHad(0.2, 50*um);
148  param->SetUseMottCorrection(true);
149  param->SetMscStepLimitType(fUseSafetyPlus);
150  param->SetMscSkin(3);
151  param->SetMscRangeFactor(0.08);
152  param->SetMuHadLateralDisplacement(true);
153  param->SetFluo(true);
154  param->SetMaxNIELEnergy(1*MeV);
155  param->SetPIXEElectronCrossSectionModel("Penelope");
156  SetPhysicsType(bElectromagnetic);
157 }
158 
159 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
160 
161 G4EmPenelopePhysics::~G4EmPenelopePhysics()
162 {}
163 
164 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
165 
166 void G4EmPenelopePhysics::ConstructParticle()
167 {
168  // gamma
169  G4Gamma::Gamma();
170 
171  // leptons
172  G4Electron::Electron();
173  G4Positron::Positron();
174  G4MuonPlus::MuonPlus();
175  G4MuonMinus::MuonMinus();
176 
177  // mesons
178  G4PionPlus::PionPlusDefinition();
179  G4PionMinus::PionMinusDefinition();
180  G4KaonPlus::KaonPlusDefinition();
181  G4KaonMinus::KaonMinusDefinition();
182 
183  // baryons
184  G4Proton::Proton();
185  G4AntiProton::AntiProton();
186 
187  // ions
188  G4Deuteron::Deuteron();
189  G4Triton::Triton();
190  G4He3::He3();
191  G4Alpha::Alpha();
192  G4GenericIon::GenericIonDefinition();
193 }
194 
195 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
196 
197 void G4EmPenelopePhysics::ConstructProcess()
198 {
199  if(verbose > 1) {
200  G4cout << "### " << GetPhysicsName() << " Construct Processes " << G4endl;
201  }
202  G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
203  G4LossTableManager* man = G4LossTableManager::Instance();
204 
205  // muon & hadron bremsstrahlung and pair production
206  G4MuBremsstrahlung* mub = new G4MuBremsstrahlung();
207  G4MuPairProduction* mup = new G4MuPairProduction();
208  G4hBremsstrahlung* pib = new G4hBremsstrahlung();
209  G4hPairProduction* pip = new G4hPairProduction();
210  G4hBremsstrahlung* kb = new G4hBremsstrahlung();
211  G4hPairProduction* kp = new G4hPairProduction();
212  G4hBremsstrahlung* pb = new G4hBremsstrahlung();
213  G4hPairProduction* pp = new G4hPairProduction();
214  G4ePairProduction* ee = new G4ePairProduction();
215 
216  // muon & hadron multiple scattering
217  G4MuMultipleScattering* mumsc = new G4MuMultipleScattering();
218  mumsc->SetEmModel(new G4WentzelVIModel());
219  G4CoulombScattering* muss = new G4CoulombScattering();
220 
221  G4hMultipleScattering* pimsc = new G4hMultipleScattering();
222  pimsc->SetEmModel(new G4WentzelVIModel());
223  G4CoulombScattering* piss = new G4CoulombScattering();
224 
225  G4hMultipleScattering* kmsc = new G4hMultipleScattering();
226  kmsc->SetEmModel(new G4WentzelVIModel());
227  G4CoulombScattering* kss = new G4CoulombScattering();
228 
229  G4hMultipleScattering* hmsc = new G4hMultipleScattering("ionmsc");
230 
231  // high energy limit for e+- scattering models
232  G4double highEnergyLimit = G4EmParameters::Instance()->MscEnergyLimit();
233 
234  // nuclear stopping
235  G4double nielEnergyLimit = G4EmParameters::Instance()->MaxNIELEnergy();
236  G4NuclearStopping* pnuc = new G4NuclearStopping();
237  pnuc->SetMaxKinEnergy(nielEnergyLimit);
238 
239  //Applicability range for Penelope models
240  //for higher energies, the Standard models are used
241  G4double PenelopeHighEnergyLimit = 1.0*GeV;
242 
243  // Add Penelope EM Processes
244  G4ParticleTable* table = G4ParticleTable::GetParticleTable();
245  for(const auto& particleName : partList.PartNames()) {
246  G4ParticleDefinition* particle = table->FindParticle(particleName);
247  if (!particle) { continue; }
248  if (particleName == "gamma") {
249 
250  //Photo-electric effect
251  G4PhotoElectricEffect* thePhotoElectricEffect = new G4PhotoElectricEffect();
252  thePhotoElectricEffect->SetEmModel(new G4PEEffectFluoModel());
253  G4PenelopePhotoElectricModel* thePEPenelopeModel = new
254  G4PenelopePhotoElectricModel();
255  thePEPenelopeModel->SetHighEnergyLimit(PenelopeHighEnergyLimit);
256  thePhotoElectricEffect->AddEmModel(0, thePEPenelopeModel);
257 
258  //Compton scattering
259  G4ComptonScattering* theComptonScattering = new G4ComptonScattering();
260  G4PenelopeComptonModel* theComptonPenelopeModel =
261  new G4PenelopeComptonModel();
262  theComptonScattering->SetEmModel(new G4KleinNishinaModel());
263  theComptonPenelopeModel->SetHighEnergyLimit(PenelopeHighEnergyLimit);
264  theComptonScattering->AddEmModel(0, theComptonPenelopeModel);
265 
266  //Gamma conversion
267  G4GammaConversion* theGammaConversion = new G4GammaConversion();
268  G4PenelopeGammaConversionModel* theGCPenelopeModel =
269  new G4PenelopeGammaConversionModel();
270  theGCPenelopeModel->SetHighEnergyLimit(PenelopeHighEnergyLimit);
271  theGammaConversion->AddEmModel(0, theGCPenelopeModel);
272 
273  //Rayleigh scattering
274  G4RayleighScattering* theRayleigh = new G4RayleighScattering();
275  G4PenelopeRayleighModel* theRayleighPenelopeModel =
276  new G4PenelopeRayleighModel();
277  theRayleigh->SetEmModel(theRayleighPenelopeModel);
278 
279  ph->RegisterProcess(thePhotoElectricEffect, particle);
280  ph->RegisterProcess(theComptonScattering, particle);
281  ph->RegisterProcess(theGammaConversion, particle);
282  ph->RegisterProcess(theRayleigh, particle);
283 
284  } else if (particleName == "e-") {
285 
286  // multiple scattering
287  G4eMultipleScattering* msc = new G4eMultipleScattering;
288  G4GoudsmitSaundersonMscModel* msc1 = new G4GoudsmitSaundersonMscModel();
289  G4WentzelVIModel* msc2 = new G4WentzelVIModel();
290  msc1->SetHighEnergyLimit(highEnergyLimit);
291  msc2->SetLowEnergyLimit(highEnergyLimit);
292  msc->SetEmModel(msc1);
293  msc->SetEmModel(msc2);
294 
295  G4eCoulombScatteringModel* ssm = new G4eCoulombScatteringModel();
296  G4CoulombScattering* ss = new G4CoulombScattering();
297  ss->SetEmModel(ssm);
298  ss->SetMinKinEnergy(highEnergyLimit);
299  ssm->SetLowEnergyLimit(highEnergyLimit);
300  ssm->SetActivationLowEnergyLimit(highEnergyLimit);
301 
302  //Ionisation
303  G4eIonisation* eioni = new G4eIonisation();
304  G4PenelopeIonisationModel* theIoniPenelope =
305  new G4PenelopeIonisationModel();
306  theIoniPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
307  eioni->AddEmModel(0,theIoniPenelope,new G4UniversalFluctuation());
308 
309  //Bremsstrahlung
310  G4eBremsstrahlung* brem = new G4eBremsstrahlung();
311  G4PenelopeBremsstrahlungModel* theBremPenelope = new
312  G4PenelopeBremsstrahlungModel();
313  theBremPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
314  brem->SetEmModel(theBremPenelope);
315 
316  // register processes
317  ph->RegisterProcess(msc, particle);
318  ph->RegisterProcess(eioni, particle);
319  ph->RegisterProcess(brem, particle);
320  ph->RegisterProcess(ee, particle);
321  ph->RegisterProcess(ss, particle);
322 
323  } else if (particleName == "e+") {
324 
325  // multiple scattering
326  G4eMultipleScattering* msc = new G4eMultipleScattering;
327  G4GoudsmitSaundersonMscModel* msc1 = new G4GoudsmitSaundersonMscModel();
328  G4WentzelVIModel* msc2 = new G4WentzelVIModel();
329  msc1->SetHighEnergyLimit(highEnergyLimit);
330  msc2->SetLowEnergyLimit(highEnergyLimit);
331  msc->SetEmModel(msc1);
332  msc->SetEmModel(msc2);
333 
334  G4eCoulombScatteringModel* ssm = new G4eCoulombScatteringModel();
335  G4CoulombScattering* ss = new G4CoulombScattering();
336  ss->SetEmModel(ssm);
337  ss->SetMinKinEnergy(highEnergyLimit);
338  ssm->SetLowEnergyLimit(highEnergyLimit);
339  ssm->SetActivationLowEnergyLimit(highEnergyLimit);
340 
341  //Ionisation
342  G4eIonisation* eioni = new G4eIonisation();
343  G4PenelopeIonisationModel* theIoniPenelope =
344  new G4PenelopeIonisationModel();
345  theIoniPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
346  eioni->AddEmModel(0,theIoniPenelope,new G4UniversalFluctuation());
347 
348  //Bremsstrahlung
349  G4eBremsstrahlung* brem = new G4eBremsstrahlung();
350  G4PenelopeBremsstrahlungModel* theBremPenelope = new
351  G4PenelopeBremsstrahlungModel();
352  theBremPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
353  brem->SetEmModel(theBremPenelope);
354 
355  //Annihilation
356  G4eplusAnnihilation* anni = new G4eplusAnnihilation();
357  G4PenelopeAnnihilationModel* theAnnPenelope = new
358  G4PenelopeAnnihilationModel();
359  theAnnPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
360  anni->AddEmModel(0, theAnnPenelope);
361 
362  // register processes
363  ph->RegisterProcess(msc, particle);
364  ph->RegisterProcess(eioni, particle);
365  ph->RegisterProcess(brem, particle);
366  ph->RegisterProcess(ee, particle);
367  ph->RegisterProcess(anni, particle);
368  ph->RegisterProcess(ss, particle);
369 
370  } else if (particleName == "mu+" ||
371  particleName == "mu-" ) {
372 
373  G4MuIonisation* muIoni = new G4MuIonisation();
374 
375  ph->RegisterProcess(mumsc, particle);
376  ph->RegisterProcess(muIoni, particle);
377  ph->RegisterProcess(mub, particle);
378  ph->RegisterProcess(mup, particle);
379  ph->RegisterProcess(muss, particle);
380 
381  } else if (particleName == "alpha" ||
382  particleName == "He3" ) {
383 
384  G4hMultipleScattering* msc = new G4hMultipleScattering();
385  G4ionIonisation* ionIoni = new G4ionIonisation();
386 
387  ph->RegisterProcess(msc, particle);
388  ph->RegisterProcess(ionIoni, particle);
389  ph->RegisterProcess(pnuc, particle);
390 
391  } else if (particleName == "GenericIon") {
392 
393  G4ionIonisation* ionIoni = new G4ionIonisation();
394  ionIoni->SetEmModel(new G4IonParametrisedLossModel());
395  ionIoni->SetStepFunction(0.1, 1*um);
396 
397  ph->RegisterProcess(hmsc, particle);
398  ph->RegisterProcess(ionIoni, particle);
399  ph->RegisterProcess(pnuc, particle);
400 
401  } else if (particleName == "pi+" ||
402  particleName == "pi-" ) {
403 
404  G4hIonisation* hIoni = new G4hIonisation();
405 
406  ph->RegisterProcess(pimsc, particle);
407  ph->RegisterProcess(hIoni, particle);
408  ph->RegisterProcess(pib, particle);
409  ph->RegisterProcess(pip, particle);
410  ph->RegisterProcess(piss, particle);
411 
412  } else if (particleName == "kaon+" ||
413  particleName == "kaon-" ) {
414 
415  G4hIonisation* hIoni = new G4hIonisation();
416 
417  ph->RegisterProcess(kmsc, particle);
418  ph->RegisterProcess(hIoni, particle);
419  ph->RegisterProcess(kb, particle);
420  ph->RegisterProcess(kp, particle);
421  ph->RegisterProcess(kss, particle);
422 
423  } else if (particleName == "proton" ||
424  particleName == "anti_proton") {
425 
426  G4hMultipleScattering* pmsc = new G4hMultipleScattering();
427  pmsc->SetEmModel(new G4WentzelVIModel());
428  G4hIonisation* hIoni = new G4hIonisation();
429 
430  ph->RegisterProcess(pmsc, particle);
431  ph->RegisterProcess(hIoni, particle);
432  ph->RegisterProcess(pb, particle);
433  ph->RegisterProcess(pp, particle);
434  ph->RegisterProcess(new G4CoulombScattering(), particle);
435 
436  } else if (particleName == "B+" ||
437  particleName == "B-" ||
438  particleName == "D+" ||
439  particleName == "D-" ||
440  particleName == "Ds+" ||
441  particleName == "Ds-" ||
442  particleName == "anti_He3" ||
443  particleName == "anti_alpha" ||
444  particleName == "anti_deuteron" ||
445  particleName == "anti_lambda_c+" ||
446  particleName == "anti_omega-" ||
447  particleName == "anti_sigma_c+" ||
448  particleName == "anti_sigma_c++" ||
449  particleName == "anti_sigma+" ||
450  particleName == "anti_sigma-" ||
451  particleName == "anti_triton" ||
452  particleName == "anti_xi_c+" ||
453  particleName == "anti_xi-" ||
454  particleName == "deuteron" ||
455  particleName == "lambda_c+" ||
456  particleName == "omega-" ||
457  particleName == "sigma_c+" ||
458  particleName == "sigma_c++" ||
459  particleName == "sigma+" ||
460  particleName == "sigma-" ||
461  particleName == "tau+" ||
462  particleName == "tau-" ||
463  particleName == "triton" ||
464  particleName == "xi_c+" ||
465  particleName == "xi-" ) {
466 
467  ph->RegisterProcess(hmsc, particle);
468  ph->RegisterProcess(new G4hIonisation(), particle);
469  ph->RegisterProcess(pnuc, particle);
470  }
471  }
472 
473  // Nuclear stopping
474  pnuc->SetMaxKinEnergy(MeV);
475 
476  // Deexcitation
477  //
478  G4VAtomDeexcitation* deexcitation = new G4UAtomicDeexcitation();
479  man->SetAtomDeexcitation(deexcitation);
480 
481  G4EmModelActivator mact(GetPhysicsName());
482 }
483 
484 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......