| Geant4 Cross Reference |
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 // $Id: G4EmPenelopePhysicsMI.cc 109526 2018-04-30 07:11:52Z gcosmo $
27 // customized by gpaterno for MI in Rayleigh Scattering, March 2019
28 //
29 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
30
31 #include "G4EmPenelopePhysicsMI.hh"
32
33 #include "G4ParticleDefinition.hh"
34 #include "G4ParticleTable.hh"
35 #include "G4SystemOfUnits.hh"
36
37 // Processes and models
38 // gamma
39 #include "G4ComptonScattering.hh"
40 #include "G4GammaConversion.hh"
41 #include "G4PenelopeComptonModel.hh"
42 #include "G4PenelopeGammaConversionModel.hh"
43 #include "G4PenelopePhotoElectricModel.hh"
44 #include "G4PenelopeRayleighModel.hh"
45 #include "G4PenelopeRayleighModelMI.hh"
46 #include "G4PhotoElectricEffect.hh"
47 #include "G4RayleighScattering.hh"
48
49 // e- and e+
50 #include "G4PenelopeBremsstrahlungModel.hh"
51 #include "G4PenelopeIonisationModel.hh"
52 #include "G4UniversalFluctuation.hh"
53 #include "G4eBremsstrahlung.hh"
54 #include "G4eIonisation.hh"
55 #include "G4eMultipleScattering.hh"
56
57 // e+ only
58 #include "G4PenelopeAnnihilationModel.hh"
59 #include "G4eplusAnnihilation.hh"
60
61 // mu
62 #include "G4MuBremsstrahlung.hh"
63 #include "G4MuBremsstrahlungModel.hh"
64 #include "G4MuIonisation.hh"
65 #include "G4MuMultipleScattering.hh"
66 #include "G4MuPairProduction.hh"
67 #include "G4MuPairProductionModel.hh"
68 #include "G4hBremsstrahlungModel.hh"
69 #include "G4hPairProductionModel.hh"
70
71 // hadrons
72 #include "G4IonParametrisedLossModel.hh"
73 #include "G4MscStepLimitType.hh"
74 #include "G4NuclearStopping.hh"
75 #include "G4hBremsstrahlung.hh"
76 #include "G4hIonisation.hh"
77 #include "G4hMultipleScattering.hh"
78 #include "G4hPairProduction.hh"
79 #include "G4ionIonisation.hh"
80
81 // msc models
82 #include "G4CoulombScattering.hh"
83 #include "G4GoudsmitSaundersonMscModel.hh"
84 #include "G4LossTableManager.hh"
85 #include "G4UAtomicDeexcitation.hh"
86 #include "G4UrbanMscModel.hh"
87 #include "G4VAtomDeexcitation.hh"
88 #include "G4WentzelVIModel.hh"
89 #include "G4eCoulombScatteringModel.hh"
90
91 // particles
92 #include "G4Alpha.hh"
93 #include "G4AntiProton.hh"
94 #include "G4BuilderType.hh"
95 #include "G4Deuteron.hh"
96 #include "G4Electron.hh"
97 #include "G4EmModelActivator.hh"
98 #include "G4Gamma.hh"
99 #include "G4GenericIon.hh"
100 #include "G4He3.hh"
101 #include "G4KaonMinus.hh"
102 #include "G4KaonPlus.hh"
103 #include "G4MuonMinus.hh"
104 #include "G4MuonPlus.hh"
105 #include "G4PhysicsListHelper.hh"
106 #include "G4PionMinus.hh"
107 #include "G4PionPlus.hh"
108 #include "G4Positron.hh"
109 #include "G4Proton.hh"
110 #include "G4Triton.hh"
111
112 // factory
113 #include "G4PhysicsConstructorFactory.hh"
114
115 G4_DECLARE_PHYSCONSTR_FACTORY(G4EmPenelopePhysicsMI);
116
117 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
118
119 G4EmPenelopePhysicsMI::G4EmPenelopePhysicsMI(G4int ver, const G4String&, G4bool UseMIFlag)
120 : G4VPhysicsConstructor("G4EmPenelopeMI"), fVerbose(ver), fUseMIFlag(UseMIFlag)
121 {
122 G4EmParameters* param = G4EmParameters::Instance();
123 param->SetDefaults();
124 param->SetVerbose(fVerbose);
125 param->SetMinEnergy(100 * eV);
126 param->SetLowestElectronEnergy(100 * eV);
127 param->SetNumberOfBinsPerDecade(20);
128 param->SetMscRangeFactor(0.02);
129 param->SetMscStepLimitType(fUseDistanceToBoundary);
130 param->SetMuHadLateralDisplacement(true);
131 param->SetFluo(true);
132 param->SetPIXEElectronCrossSectionModel("Penelope");
133 SetPhysicsType(bElectromagnetic);
134 }
135
136 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
137
138 G4EmPenelopePhysicsMI::~G4EmPenelopePhysicsMI() {}
139
140 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
141
142 void G4EmPenelopePhysicsMI::ConstructParticle()
143 {
144 // gamma
145 G4Gamma::Gamma();
146
147 // leptons
148 G4Electron::Electron();
149 G4Positron::Positron();
150 G4MuonPlus::MuonPlus();
151 G4MuonMinus::MuonMinus();
152
153 // mesons
154 G4PionPlus::PionPlusDefinition();
155 G4PionMinus::PionMinusDefinition();
156 G4KaonPlus::KaonPlusDefinition();
157 G4KaonMinus::KaonMinusDefinition();
158
159 // baryons
160 G4Proton::Proton();
161 G4AntiProton::AntiProton();
162
163 // ions
164 G4Deuteron::Deuteron();
165 G4Triton::Triton();
166 G4He3::He3();
167 G4Alpha::Alpha();
168 G4GenericIon::GenericIonDefinition();
169 }
170
171 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
172
173 void G4EmPenelopePhysicsMI::ConstructProcess()
174 {
175 if (fVerbose > 1) {
176 G4cout << "### " << GetPhysicsName() << " Construct Processes " << G4endl;
177 }
178
179 G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
180
181 // muon & hadron bremsstrahlung and pair production
182 G4MuBremsstrahlung* mub = new G4MuBremsstrahlung();
183 G4MuPairProduction* mup = new G4MuPairProduction();
184 G4hBremsstrahlung* pib = new G4hBremsstrahlung();
185 G4hPairProduction* pip = new G4hPairProduction();
186 G4hBremsstrahlung* kb = new G4hBremsstrahlung();
187 G4hPairProduction* kp = new G4hPairProduction();
188 G4hBremsstrahlung* pb = new G4hBremsstrahlung();
189 G4hPairProduction* pp = new G4hPairProduction();
190
191 // muon & hadron multiple scattering
192 G4MuMultipleScattering* mumsc = new G4MuMultipleScattering();
193 mumsc->SetEmModel(new G4WentzelVIModel());
194 G4hMultipleScattering* hmsc = new G4hMultipleScattering("ionmsc");
195
196 // high energy limit for e+- scattering models
197 G4double highEnergyLimit = G4EmParameters::Instance()->MscEnergyLimit();
198
199 // nuclear stopping
200 G4NuclearStopping* pnuc = new G4NuclearStopping();
201
202 // Applicability range for Penelope models
203 // for higher energies, the Standard models are used
204 G4double PenelopeHighEnergyLimit = 1.0 * GeV;
205
206 // Add Penelope EM Processes
207 G4ParticleTable* table = G4ParticleTable::GetParticleTable();
208 for (const auto& particleName : fPartList.PartNames()) {
209 G4ParticleDefinition* particle = table->FindParticle(particleName);
210 if (!particle) {
211 continue;
212 }
213 if (particleName == "gamma") {
214 // Photo-electric effect
215 G4PhotoElectricEffect* thePhotoElectricEffect = new G4PhotoElectricEffect();
216 G4PenelopePhotoElectricModel* thePEPenelopeModel = new G4PenelopePhotoElectricModel();
217 thePEPenelopeModel->SetHighEnergyLimit(PenelopeHighEnergyLimit);
218 thePhotoElectricEffect->SetEmModel(thePEPenelopeModel);
219 ph->RegisterProcess(thePhotoElectricEffect, particle);
220
221 // Compton scattering
222 G4ComptonScattering* theComptonScattering = new G4ComptonScattering();
223 G4PenelopeComptonModel* theComptonPenelopeModel = new G4PenelopeComptonModel();
224 theComptonPenelopeModel->SetHighEnergyLimit(PenelopeHighEnergyLimit);
225 theComptonScattering->SetEmModel(theComptonPenelopeModel);
226 ph->RegisterProcess(theComptonScattering, particle);
227
228 // Gamma conversion
229 G4GammaConversion* theGammaConversion = new G4GammaConversion();
230 G4PenelopeGammaConversionModel* theGCPenelopeModel = new G4PenelopeGammaConversionModel();
231 theGammaConversion->SetEmModel(theGCPenelopeModel);
232 ph->RegisterProcess(theGammaConversion, particle);
233
234 // Rayleigh scattering (modified by gpaterno)
235 G4RayleighScattering* theRayleigh = new G4RayleighScattering();
236 G4PenelopeRayleighModelMI* theRayleighPenelopeModel = new G4PenelopeRayleighModelMI();
237 theRayleighPenelopeModel->SetVerbosityLevel(1);
238 theRayleighPenelopeModel->SetMIActive(fUseMIFlag);
239 // theRayleighPenelopeModel->SetHighEnergyLimit(PenelopeHighEnergyLimit);
240 theRayleigh->SetEmModel(theRayleighPenelopeModel);
241 ph->RegisterProcess(theRayleigh, particle);
242 }
243 else if (particleName == "e-") {
244 // multiple scattering
245 G4eMultipleScattering* msc = new G4eMultipleScattering;
246 G4UrbanMscModel* msc1 = new G4UrbanMscModel();
247 G4WentzelVIModel* msc2 = new G4WentzelVIModel();
248 msc1->SetHighEnergyLimit(highEnergyLimit);
249 msc2->SetLowEnergyLimit(highEnergyLimit);
250 msc->SetEmModel(msc1);
251 msc->SetEmModel(msc2);
252
253 G4eCoulombScatteringModel* ssm = new G4eCoulombScatteringModel();
254 G4CoulombScattering* ss = new G4CoulombScattering();
255 ss->SetEmModel(ssm);
256 ss->SetMinKinEnergy(highEnergyLimit);
257 ssm->SetLowEnergyLimit(highEnergyLimit);
258 ssm->SetActivationLowEnergyLimit(highEnergyLimit);
259
260 // Ionisation
261 G4eIonisation* eIoni = new G4eIonisation();
262 G4PenelopeIonisationModel* theIoniPenelope = new G4PenelopeIonisationModel();
263 theIoniPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
264 eIoni->AddEmModel(0, theIoniPenelope, new G4UniversalFluctuation());
265 eIoni->SetStepFunction(0.2, 100 * um); //
266
267 // Bremsstrahlung
268 G4eBremsstrahlung* eBrem = new G4eBremsstrahlung();
269 G4PenelopeBremsstrahlungModel* theBremPenelope = new G4PenelopeBremsstrahlungModel();
270 theBremPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
271 eBrem->SetEmModel(theBremPenelope);
272
273 // register processes
274 ph->RegisterProcess(msc, particle);
275 ph->RegisterProcess(eIoni, particle);
276 ph->RegisterProcess(eBrem, particle);
277 ph->RegisterProcess(ss, particle);
278 }
279 else if (particleName == "e+") {
280 // multiple scattering
281 G4eMultipleScattering* msc = new G4eMultipleScattering;
282 G4UrbanMscModel* msc1 = new G4UrbanMscModel();
283 G4WentzelVIModel* msc2 = new G4WentzelVIModel();
284 msc1->SetHighEnergyLimit(highEnergyLimit);
285 msc2->SetLowEnergyLimit(highEnergyLimit);
286 msc->SetEmModel(msc1);
287 msc->SetEmModel(msc2);
288
289 G4eCoulombScatteringModel* ssm = new G4eCoulombScatteringModel();
290 G4CoulombScattering* ss = new G4CoulombScattering();
291 ss->SetEmModel(ssm);
292 ss->SetMinKinEnergy(highEnergyLimit);
293 ssm->SetLowEnergyLimit(highEnergyLimit);
294 ssm->SetActivationLowEnergyLimit(highEnergyLimit);
295
296 // Ionisation
297 G4eIonisation* eIoni = new G4eIonisation();
298 G4PenelopeIonisationModel* theIoniPenelope = new G4PenelopeIonisationModel();
299 theIoniPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
300 eIoni->AddEmModel(0, theIoniPenelope, new G4UniversalFluctuation());
301 eIoni->SetStepFunction(0.2, 100 * um); //
302
303 // Bremsstrahlung
304 G4eBremsstrahlung* eBrem = new G4eBremsstrahlung();
305 G4PenelopeBremsstrahlungModel* theBremPenelope = new G4PenelopeBremsstrahlungModel();
306 theBremPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
307 eBrem->SetEmModel(theBremPenelope);
308
309 // Annihilation
310 G4eplusAnnihilation* eAnni = new G4eplusAnnihilation();
311 G4PenelopeAnnihilationModel* theAnnPenelope = new G4PenelopeAnnihilationModel();
312 theAnnPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
313 eAnni->AddEmModel(0, theAnnPenelope);
314
315 // register processes
316 ph->RegisterProcess(msc, particle);
317 ph->RegisterProcess(eIoni, particle);
318 ph->RegisterProcess(eBrem, particle);
319 ph->RegisterProcess(eAnni, particle);
320 ph->RegisterProcess(ss, particle);
321 }
322 else if (particleName == "mu+" || particleName == "mu-") {
323 G4MuIonisation* muIoni = new G4MuIonisation();
324 muIoni->SetStepFunction(0.2, 50 * um);
325
326 ph->RegisterProcess(mumsc, particle);
327 ph->RegisterProcess(muIoni, particle);
328 ph->RegisterProcess(mub, particle);
329 ph->RegisterProcess(mup, particle);
330 ph->RegisterProcess(new G4CoulombScattering(), particle);
331 }
332 else if (particleName == "alpha" || particleName == "He3") {
333 G4hMultipleScattering* msc = new G4hMultipleScattering();
334 G4ionIonisation* ionIoni = new G4ionIonisation();
335 ionIoni->SetStepFunction(0.1, 10 * um);
336
337 ph->RegisterProcess(msc, particle);
338 ph->RegisterProcess(ionIoni, particle);
339 ph->RegisterProcess(pnuc, particle);
340 }
341 else if (particleName == "GenericIon") {
342 G4ionIonisation* ionIoni = new G4ionIonisation();
343 ionIoni->SetEmModel(new G4IonParametrisedLossModel());
344 ionIoni->SetStepFunction(0.1, 1 * um);
345
346 ph->RegisterProcess(hmsc, particle);
347 ph->RegisterProcess(ionIoni, particle);
348 ph->RegisterProcess(pnuc, particle);
349 }
350 else if (particleName == "pi+" || particleName == "pi-") {
351 G4hMultipleScattering* pimsc = new G4hMultipleScattering();
352 G4hIonisation* hIoni = new G4hIonisation();
353 hIoni->SetStepFunction(0.2, 50 * um);
354
355 ph->RegisterProcess(pimsc, particle);
356 ph->RegisterProcess(hIoni, particle);
357 ph->RegisterProcess(pib, particle);
358 ph->RegisterProcess(pip, particle);
359 }
360 else if (particleName == "kaon+" || particleName == "kaon-") {
361 G4hMultipleScattering* kmsc = new G4hMultipleScattering();
362 G4hIonisation* hIoni = new G4hIonisation();
363 hIoni->SetStepFunction(0.2, 50 * um);
364
365 ph->RegisterProcess(kmsc, particle);
366 ph->RegisterProcess(hIoni, particle);
367 ph->RegisterProcess(kb, particle);
368 ph->RegisterProcess(kp, particle);
369 }
370 else if (particleName == "proton" || particleName == "anti_proton") {
371 G4hMultipleScattering* pmsc = new G4hMultipleScattering();
372 G4hIonisation* hIoni = new G4hIonisation();
373 hIoni->SetStepFunction(0.2, 50 * um);
374
375 ph->RegisterProcess(pmsc, particle);
376 ph->RegisterProcess(hIoni, particle);
377 ph->RegisterProcess(pb, particle);
378 ph->RegisterProcess(pp, particle);
379 ph->RegisterProcess(pnuc, particle);
380 }
381 else if (particleName == "B+" || particleName == "B-" || particleName == "D+"
382 || particleName == "D-" || particleName == "Ds+" || particleName == "Ds-"
383 || particleName == "anti_He3" || particleName == "anti_alpha"
384 || particleName == "anti_deuteron" || particleName == "anti_lambda_c+"
385 || particleName == "anti_omega-" || particleName == "anti_sigma_c+"
386 || particleName == "anti_sigma_c++" || particleName == "anti_sigma+"
387 || particleName == "anti_sigma-" || particleName == "anti_triton"
388 || particleName == "anti_xi_c+" || particleName == "anti_xi-"
389 || particleName == "deuteron" || particleName == "lambda_c+"
390 || particleName == "omega-" || particleName == "sigma_c+"
391 || particleName == "sigma_c++" || particleName == "sigma+" || particleName == "sigma-"
392 || particleName == "tau+" || particleName == "tau-" || particleName == "triton"
393 || particleName == "xi_c+" || particleName == "xi-")
394 {
395 ph->RegisterProcess(hmsc, particle);
396 ph->RegisterProcess(new G4hIonisation(), particle);
397 }
398 }
399
400 // Nuclear stopping
401 pnuc->SetMaxKinEnergy(MeV);
402
403 // Deexcitation
404 G4VAtomDeexcitation* deexcitation = new G4UAtomicDeexcitation();
405 G4LossTableManager::Instance()->SetAtomDeexcitation(deexcitation);
406
407 G4EmModelActivator mact(GetPhysicsName());
408 }
409
410 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
411