BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSEmStandardPhysicsOp4Channelling.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// $Id: G4EmStandardPhysics_option4_channeling.cc 99938 2016-10-12 08:06:52Z gcosmo $
27//
28//---------------------------------------------------------------------------
29//
30// ClassName: G4EmStandardPhysics_option4_channeling
31//
32// Author: V.Ivanchenko 28.09.2012 from Option3 physics constructor
33//
34// Modified:
35//
36//----------------------------------------------------------------------------
37//
38#include "G4Version.hh"
39// just exclude whole physics list for < 10.4 as it's only used with channelling that's in 10.4
40#if G4VERSION_NUMBER > 1039
41
42
43
44#include "BDSEmStandardPhysicsOp4Channelling.hh"
45
46#include "G4SystemOfUnits.hh"
47#include "G4ParticleDefinition.hh"
48#include "G4LossTableManager.hh"
49#include "G4EmParameters.hh"
50
51#include "G4ComptonScattering.hh"
52#include "G4GammaConversion.hh"
53#include "G4PhotoElectricEffect.hh"
54#include "G4RayleighScattering.hh"
55#include "G4PEEffectFluoModel.hh"
56#include "G4KleinNishinaModel.hh"
57#include "G4LowEPComptonModel.hh"
58#include "G4PenelopeGammaConversionModel.hh"
59#include "G4LivermorePhotoElectricModel.hh"
60
61#include "G4eMultipleScattering.hh"
62#include "G4MuMultipleScattering.hh"
63#include "G4hMultipleScattering.hh"
64#include "G4MscStepLimitType.hh"
65#include "G4UrbanMscModel.hh"
66#include "G4DummyModel.hh"
67#include "G4WentzelVIModel.hh"
68#include "G4CoulombScattering.hh"
69#include "G4eCoulombScatteringModel.hh"
70
71#include "G4eIonisation.hh"
72#include "G4eBremsstrahlung.hh"
73#include "G4Generator2BS.hh"
74#include "G4Generator2BN.hh"
75#include "G4SeltzerBergerModel.hh"
76#include "G4PenelopeIonisationModel.hh"
77#include "G4UniversalFluctuation.hh"
78
79#include "G4eplusAnnihilation.hh"
80#include "G4UAtomicDeexcitation.hh"
81
82#include "G4MuIonisation.hh"
83#include "G4MuBremsstrahlung.hh"
84#include "G4MuPairProduction.hh"
85#include "G4hBremsstrahlung.hh"
86#include "G4hPairProduction.hh"
87#include "G4ePairProduction.hh"
88
89#include "G4MuBremsstrahlungModel.hh"
90#include "G4MuPairProductionModel.hh"
91#include "G4hBremsstrahlungModel.hh"
92#include "G4hPairProductionModel.hh"
93
94#include "G4hIonisation.hh"
95#include "G4ionIonisation.hh"
96#include "G4IonParametrisedLossModel.hh"
97#include "G4NuclearStopping.hh"
98
99#include "G4Gamma.hh"
100#include "G4Electron.hh"
101#include "G4Positron.hh"
102#include "G4MuonPlus.hh"
103#include "G4MuonMinus.hh"
104#include "G4PionPlus.hh"
105#include "G4PionMinus.hh"
106#include "G4KaonPlus.hh"
107#include "G4KaonMinus.hh"
108#include "G4Proton.hh"
109#include "G4AntiProton.hh"
110#include "G4Deuteron.hh"
111#include "G4Triton.hh"
112#include "G4He3.hh"
113#include "G4Alpha.hh"
114#include "G4GenericIon.hh"
115
116#include "G4PhysicsListHelper.hh"
117#include "G4BuilderType.hh"
118#include "G4EmModelActivator.hh"
119
120// factory
121#include "G4PhysicsConstructorFactory.hh"
122//
123//G4_DECLARE_PHYSCONSTR_FACTORY(BDSEmStandardPhysicsOp4Channelling);
124
125//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
126
127BDSEmStandardPhysicsOp4Channelling::BDSEmStandardPhysicsOp4Channelling(G4int ver,
128 const G4String&)
129 : G4VPhysicsConstructor("G4EmStandard_opt4_channeling"), verbose(ver)
130{
131 G4EmParameters* param = G4EmParameters::Instance();
132 param->SetDefaults();
133 param->SetVerbose(verbose);
134 param->SetMinEnergy(100*eV);
135 param->SetMaxEnergy(10*TeV);
136 param->SetLowestElectronEnergy(10*eV);
137 param->SetNumberOfBinsPerDecade(20);
138 param->ActivateAngularGeneratorForIonisation(true);
139 param->SetMscThetaLimit(0.0);
140 param->SetFluo(true);
141 param->SetAuger(true);
142 param->SetPixe(true);
143 SetPhysicsType(bElectromagnetic);
144}
145
146//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
147
148BDSEmStandardPhysicsOp4Channelling::~BDSEmStandardPhysicsOp4Channelling()
149{;}
150
151void BDSEmStandardPhysicsOp4Channelling::ConstructParticle()
152{
153 // gamma
154 G4Gamma::Gamma();
155
156 // leptons
157 G4Electron::Electron();
158 G4Positron::Positron();
159 G4MuonPlus::MuonPlus();
160 G4MuonMinus::MuonMinus();
161
162 // mesons
163 G4PionPlus::PionPlusDefinition();
164 G4PionMinus::PionMinusDefinition();
165 G4KaonPlus::KaonPlusDefinition();
166 G4KaonMinus::KaonMinusDefinition();
167
168 // barions
169 G4Proton::Proton();
170 G4AntiProton::AntiProton();
171
172 // ions
173 G4Deuteron::Deuteron();
174 G4Triton::Triton();
175 G4He3::He3();
176 G4Alpha::Alpha();
177 G4GenericIon::GenericIonDefinition();
178}
179
180//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
181
182void BDSEmStandardPhysicsOp4Channelling::ConstructProcess()
183{
184 if(verbose > 1) {
185 G4cout << "### " << GetPhysicsName() << " Construct Processes " << G4endl;
186 }
187 G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
188
189 // muon & hadron bremsstrahlung and pair production
190 G4MuBremsstrahlung* mub = new G4MuBremsstrahlung();
191 G4MuPairProduction* mup = new G4MuPairProduction();
192 G4hBremsstrahlung* pib = new G4hBremsstrahlung();
193 G4hPairProduction* pip = new G4hPairProduction();
194 G4hBremsstrahlung* kb = new G4hBremsstrahlung();
195 G4hPairProduction* kp = new G4hPairProduction();
196 G4hBremsstrahlung* pb = new G4hBremsstrahlung();
197 G4hPairProduction* pp = new G4hPairProduction();
198 G4ePairProduction* ee = new G4ePairProduction();
199
200 // muon & hadron multiple scattering
201 G4CoulombScattering* muss = new G4CoulombScattering();
202 G4CoulombScattering* piss = new G4CoulombScattering();
203 G4CoulombScattering* kss = new G4CoulombScattering();
204 G4CoulombScattering* pss = new G4CoulombScattering();
205
206 // energy limits for e+- scattering models
207 //G4double highEnergyLimit = 100*MeV;
208 // energy limits for e+- ionisation models
209 G4double penEnergyLimit = 1*MeV;
210
211 // nuclear stopping
212 G4NuclearStopping* pnuc = new G4NuclearStopping();
213
214 // Add standard EM Processes
215#if G4VERSION_NUMBER > 1029
216 auto aParticleIterator=GetParticleIterator();
217#endif
218 aParticleIterator->reset();
219 while( (*aParticleIterator)() ){
220 G4ParticleDefinition* particle = aParticleIterator->value();
221 G4String particleName = particle->GetParticleName();
222
223 if (particleName == "gamma") {
224
225 // Photoelectric
226 G4PhotoElectricEffect* pe = new G4PhotoElectricEffect();
227 G4VEmModel* theLivermorePEModel = new G4LivermorePhotoElectricModel();
228 pe->SetEmModel(theLivermorePEModel,1);
229 ph->RegisterProcess(pe, particle);
230
231 // Compton scattering
232 G4ComptonScattering* cs = new G4ComptonScattering;
233 cs->SetEmModel(new G4KleinNishinaModel(),1);
234 G4VEmModel* theLowEPComptonModel = new G4LowEPComptonModel();
235 theLowEPComptonModel->SetHighEnergyLimit(20*MeV);
236 cs->AddEmModel(0, theLowEPComptonModel);
237 ph->RegisterProcess(cs, particle);
238
239 // Gamma conversion
240 G4GammaConversion* gc = new G4GammaConversion();
241 G4VEmModel* thePenelopeGCModel = new G4PenelopeGammaConversionModel();
242 thePenelopeGCModel->SetHighEnergyLimit(1*GeV);
243 gc->SetEmModel(thePenelopeGCModel,1);
244 ph->RegisterProcess(gc, particle);
245
246 // Rayleigh scattering
247 ph->RegisterProcess(new G4RayleighScattering(), particle);
248
249 } else if (particleName == "e-") {
250
251 // multiple scattering
252 G4eCoulombScatteringModel* ssm = new G4eCoulombScatteringModel();
253 G4CoulombScattering* ss = new G4CoulombScattering();
254 ss->SetEmModel(ssm, 1);
255
256 // ionisation
257 G4eIonisation* eIoni = new G4eIonisation();
258 eIoni->SetStepFunction(0.2, 100*um);
259 G4PenelopeIonisationModel* pen = new G4PenelopeIonisationModel();
260 pen->SetHighEnergyLimit(penEnergyLimit);
261 eIoni->AddEmModel(0, pen, new G4UniversalFluctuation());
262
263 // bremsstrahlung
264 G4eBremsstrahlung* brem = new G4eBremsstrahlung();
265 G4SeltzerBergerModel* br1 = new G4SeltzerBergerModel();
266 G4eBremsstrahlungRelModel* br2 = new G4eBremsstrahlungRelModel();
267 br1->SetAngularDistribution(new G4Generator2BS());
268 br2->SetAngularDistribution(new G4Generator2BS());
269 brem->SetEmModel(br1,1);
270 brem->SetEmModel(br2,2);
271 br2->SetLowEnergyLimit(GeV);
272
273 // register processes
274 ph->RegisterProcess(eIoni, particle);
275 ph->RegisterProcess(brem, particle);
276#if G4VERSION_NUMBER > 1029
277 ph->RegisterProcess(ee, particle);
278#endif
279 ph->RegisterProcess(ss, particle);
280
281 } else if (particleName == "e+") {
282
283 // multiple scattering
284 G4eCoulombScatteringModel* ssm = new G4eCoulombScatteringModel();
285 G4CoulombScattering* ss = new G4CoulombScattering();
286 ss->SetEmModel(ssm, 1);
287
288 // ionisation
289 G4eIonisation* eIoni = new G4eIonisation();
290 eIoni->SetStepFunction(0.2, 100*um);
291 G4PenelopeIonisationModel* pen = new G4PenelopeIonisationModel();
292 pen->SetHighEnergyLimit(penEnergyLimit);
293 eIoni->AddEmModel(0, pen, new G4UniversalFluctuation());
294
295 // bremsstrahlung
296 G4eBremsstrahlung* brem = new G4eBremsstrahlung();
297 G4SeltzerBergerModel* br1 = new G4SeltzerBergerModel();
298 G4eBremsstrahlungRelModel* br2 = new G4eBremsstrahlungRelModel();
299 br1->SetAngularDistribution(new G4Generator2BS());
300 br2->SetAngularDistribution(new G4Generator2BS());
301 brem->SetEmModel(br1,1);
302 brem->SetEmModel(br2,2);
303 br2->SetLowEnergyLimit(GeV);
304
305 // register processes
306 ph->RegisterProcess(eIoni, particle);
307 ph->RegisterProcess(brem, particle);
308#if G4VERSION_NUMBER > 1029
309 ph->RegisterProcess(ee, particle);
310#endif
311 ph->RegisterProcess(new G4eplusAnnihilation(), particle);
312 ph->RegisterProcess(ss, particle);
313
314 } else if (particleName == "mu+" ||
315 particleName == "mu-" ) {
316
317 G4MuIonisation* muIoni = new G4MuIonisation();
318 muIoni->SetStepFunction(0.2, 50*um);
319
320 ph->RegisterProcess(muIoni, particle);
321 ph->RegisterProcess(mub, particle);
322 ph->RegisterProcess(mup, particle);
323 ph->RegisterProcess(muss, particle);
324
325 } else if (particleName == "alpha" ||
326 particleName == "He3") {
327
328 G4ionIonisation* ionIoni = new G4ionIonisation();
329 ionIoni->SetStepFunction(0.1, 10*um);
330
331 ph->RegisterProcess(new G4CoulombScattering(), particle);
332 ph->RegisterProcess(ionIoni, particle);
333 //ph->RegisterProcess(pnuc, particle);
334
335 } else if (particleName == "GenericIon") {
336
337 G4ionIonisation* ionIoni = new G4ionIonisation();
338 ionIoni->SetEmModel(new G4IonParametrisedLossModel());
339 ionIoni->SetStepFunction(0.1, 1*um);
340
341 ph->RegisterProcess(new G4CoulombScattering(), particle);
342 ph->RegisterProcess(ionIoni, particle);
343 //ph->RegisterProcess(pnuc, particle);
344
345 } else if (particleName == "pi+" ||
346 particleName == "pi-" ) {
347
348 //G4hMultipleScattering* pimsc = new G4hMultipleScattering();
349 G4hIonisation* hIoni = new G4hIonisation();
350 hIoni->SetStepFunction(0.2, 50*um);
351
352 ph->RegisterProcess(hIoni, particle);
353 ph->RegisterProcess(pib, particle);
354 ph->RegisterProcess(pip, particle);
355 ph->RegisterProcess(piss, particle);
356
357 } else if (particleName == "kaon+" ||
358 particleName == "kaon-" ) {
359
360 //G4hMultipleScattering* kmsc = new G4hMultipleScattering();
361 G4hIonisation* hIoni = new G4hIonisation();
362 hIoni->SetStepFunction(0.2, 50*um);
363
364 ph->RegisterProcess(hIoni, particle);
365 ph->RegisterProcess(kb, particle);
366 ph->RegisterProcess(kp, particle);
367 ph->RegisterProcess(kss, particle);
368
369 } else if (particleName == "proton" ||
370 particleName == "anti_proton") {
371
372 G4hIonisation* hIoni = new G4hIonisation();
373 hIoni->SetStepFunction(0.1, 20*um);
374
375 ph->RegisterProcess(hIoni, particle);
376 ph->RegisterProcess(pb, particle);
377 ph->RegisterProcess(pp, particle);
378 ph->RegisterProcess(pss, particle);
379 //ph->RegisterProcess(pnuc, particle);
380
381 } else if (particleName == "B+" ||
382 particleName == "B-" ||
383 particleName == "D+" ||
384 particleName == "D-" ||
385 particleName == "Ds+" ||
386 particleName == "Ds-" ||
387 particleName == "anti_He3" ||
388 particleName == "anti_alpha" ||
389 particleName == "anti_deuteron" ||
390 particleName == "anti_lambda_c+" ||
391 particleName == "anti_omega-" ||
392 particleName == "anti_sigma_c+" ||
393 particleName == "anti_sigma_c++" ||
394 particleName == "anti_sigma+" ||
395 particleName == "anti_sigma-" ||
396 particleName == "anti_triton" ||
397 particleName == "anti_xi_c+" ||
398 particleName == "anti_xi-" ||
399 particleName == "deuteron" ||
400 particleName == "lambda_c+" ||
401 particleName == "omega-" ||
402 particleName == "sigma_c+" ||
403 particleName == "sigma_c++" ||
404 particleName == "sigma+" ||
405 particleName == "sigma-" ||
406 particleName == "tau+" ||
407 particleName == "tau-" ||
408 particleName == "triton" ||
409 particleName == "xi_c+" ||
410 particleName == "xi-" ) {
411
412 ph->RegisterProcess(new G4CoulombScattering(), particle);
413 ph->RegisterProcess(new G4hIonisation(), particle);
414 //ph->RegisterProcess(pnuc, particle);
415 }
416 }
417
418 // Nuclear stopping
419 pnuc->SetMaxKinEnergy(MeV);
420
421 // Deexcitation
422 G4VAtomDeexcitation* de = new G4UAtomicDeexcitation();
423 G4LossTableManager::Instance()->SetAtomDeexcitation(de);
424#if G4VERSION_NUMBER > 1029
425 G4EmModelActivator mact(GetPhysicsName());
426#else
427 G4EmModelActivator mact;
428#endif
429}
430
431#endif
432
433//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......