Geant4 10.7.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4EmDNAPhysics_stationary.cc
Go to the documentation of this file.
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
27
28#include "G4SystemOfUnits.hh"
29
31
32// *** Processes and models for Geant4-DNA
33
34#include "G4DNAElastic.hh"
38
39#include "G4DNAExcitation.hh"
40#include "G4DNAAttachment.hh"
41#include "G4DNAVibExcitation.hh"
42#include "G4DNAIonisation.hh"
45
46// particles
47
48#include "G4Electron.hh"
49#include "G4Proton.hh"
50#include "G4GenericIon.hh"
51
52// Warning : the following is needed in order to use EM Physics builders
53// e+
54#include "G4Positron.hh"
56#include "G4eIonisation.hh"
57#include "G4eBremsstrahlung.hh"
59// gamma
60#include "G4Gamma.hh"
65#include "G4GammaConversion.hh"
69
70#include "G4EmParameters.hh"
71// end of warning
72
73#include "G4LossTableManager.hh"
76#include "G4BuilderType.hh"
77
78// factory
80//
82
83
84//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
85
87 : G4VPhysicsConstructor("G4EmDNAPhysics_stationary"), verbose(ver)
88{
90 param->SetDefaults();
91 param->SetFluo(true);
92 param->SetAuger(true);
93 param->SetAugerCascade(true);
94 param->SetDeexcitationIgnoreCut(true);
95 param->ActivateDNA();
96
98}
99
100//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
101
103 : G4VPhysicsConstructor("G4EmDNAPhysics_stationary"), verbose(ver)
104{
106 param->SetDefaults();
107 param->SetFluo(true);
108 param->SetAuger(true);
109 param->SetAugerCascade(true);
110 param->SetDeexcitationIgnoreCut(true);
111
113}
114
115//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
116
118{}
119
120//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
121
123{
124// bosons
126
127// leptons
130
131// baryons
133
135
136 G4DNAGenericIonsManager * genericIonsManager;
137 genericIonsManager=G4DNAGenericIonsManager::Instance();
138 genericIonsManager->GetIon("alpha++");
139 genericIonsManager->GetIon("alpha+");
140 genericIonsManager->GetIon("helium");
141 genericIonsManager->GetIon("hydrogen");
142
143}
144
145//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
146
148{
149 if(verbose > 1) {
150 G4cout << "### " << GetPhysicsName() << " Construct Processes " << G4endl;
151 }
153
154 auto myParticleIterator=GetParticleIterator();
155 myParticleIterator->reset();
156 while( (*myParticleIterator)() )
157 {
158 G4ParticleDefinition* particle = myParticleIterator->value();
159 G4String particleName = particle->GetParticleName();
160
161 if (particleName == "e-") {
162
163 // *** Elastic scattering (two alternative models available) ***
164
165 G4DNAElastic* theDNAElasticProcess = new G4DNAElastic("e-_G4DNAElastic");
166 theDNAElasticProcess->SetEmModel(new G4DNAChampionElasticModel());
167
168 // or alternative model
169 //theDNAElasticProcess
170 //->SetEmModel(new G4DNAScreenedRutherfordElasticModel());
171
172 ph->RegisterProcess(theDNAElasticProcess, particle);
173
174 // *** Excitation ***
175
176 G4DNAExcitation* theDNAExcitationProcess =
177 new G4DNAExcitation("e-_G4DNAExcitation");
179 theDNAExcitationProcess->SetEmModel(modB);
180 modB->SelectStationary(true);
181 ph->RegisterProcess(theDNAExcitationProcess, particle);
182
183 // *** Ionisation ***
184
185 G4DNAIonisation* theDNAIonisationProcess =
186 new G4DNAIonisation("e-_G4DNAIonisation");
188 theDNAIonisationProcess->SetEmModel(modI);
189 modI->SelectStationary(true);
190 ph->RegisterProcess(theDNAIonisationProcess, particle);
191
192 // *** Vibrational excitation ***
193
194 G4DNAVibExcitation* theDNAVibExcitationProcess =
195 new G4DNAVibExcitation("e-_G4DNAVibExcitation");
197 theDNAVibExcitationProcess->SetEmModel(modS);
198 modS->SelectStationary(true);
199 ph->RegisterProcess(theDNAVibExcitationProcess, particle);
200
201 // *** Attachment ***
202
203 G4DNAAttachment* theDNAAttachmentProcess =
204 new G4DNAAttachment("e-_G4DNAAttachment");
206 theDNAAttachmentProcess->SetEmModel(modM);
207 modM->SelectStationary(true);
208 ph->RegisterProcess(theDNAAttachmentProcess, particle);
209
210 } else if ( particleName == "proton" ) {
211
212 // *** Elastic ***
213
214 G4DNAElastic* theDNAElasticProcess =
215 new G4DNAElastic("proton_G4DNAElastic");
216 theDNAElasticProcess->SetEmModel(new G4DNAIonElasticModel());
217 ((G4DNAIonElasticModel*)(theDNAElasticProcess->EmModel()))
218 ->SelectStationary(true);
219 ph->RegisterProcess(theDNAElasticProcess, particle);
220
221 // *** Excitation ***
222
223 G4DNAExcitation* theDNAExcitationProcess =
224 new G4DNAExcitation("proton_G4DNAExcitation");
225
226 theDNAExcitationProcess
228 theDNAExcitationProcess
230
232 (theDNAExcitationProcess->EmModel(0)))->SetLowEnergyLimit(10*eV);
234 (theDNAExcitationProcess->EmModel(0)))->SetHighEnergyLimit(500*keV);
236 (theDNAExcitationProcess->EmModel(0)))->SelectStationary(true);
237
239 (theDNAExcitationProcess->EmModel(1)))->SetLowEnergyLimit(500*keV);
241 (theDNAExcitationProcess->EmModel(1)))->SetHighEnergyLimit(100*MeV);
243 (theDNAExcitationProcess->EmModel(1)))->SelectStationary(true);
244
245 ph->RegisterProcess(theDNAExcitationProcess, particle);
246
247 // *** Ionisation ***
248
249 G4DNAIonisation* theDNAIonisationProcess =
250 new G4DNAIonisation("proton_G4DNAIonisation");
251
252 theDNAIonisationProcess->SetEmModel(new G4DNARuddIonisationModel);
253 theDNAIonisationProcess->SetEmModel(new G4DNABornIonisationModel);
254
255 ((G4DNARuddIonisationModel*)(theDNAIonisationProcess->EmModel(0)))
256 ->SetLowEnergyLimit(0*eV);
257 ((G4DNARuddIonisationModel*)(theDNAIonisationProcess->EmModel(0)))
258 ->SetHighEnergyLimit(500*keV);
259 ((G4DNARuddIonisationModel*)(theDNAIonisationProcess->EmModel(0)))
260 ->SelectStationary(true);
261
262 ((G4DNABornIonisationModel*)(theDNAIonisationProcess->EmModel(1)))
263 ->SetLowEnergyLimit(500*keV);
264 ((G4DNABornIonisationModel*)(theDNAIonisationProcess->EmModel(1)))
265 ->SetHighEnergyLimit(100*MeV);
266 ((G4DNABornIonisationModel*)(theDNAIonisationProcess->EmModel(1)))
267 ->SelectStationary(true);
268
269 ph->RegisterProcess(theDNAIonisationProcess, particle);
270
271 // *** Charge decrease ***
272
273 G4DNAChargeDecrease* theDNAChargeDecreaseProcess =
274 new G4DNAChargeDecrease("proton_G4DNAChargeDecrease");
275 theDNAChargeDecreaseProcess
278 (theDNAChargeDecreaseProcess->EmModel()))->SelectStationary(true);
279 ph->RegisterProcess(theDNAChargeDecreaseProcess, particle);
280
281 } else if ( particleName == "hydrogen" ) {
282
283 // *** Elastic ***
284
285 G4DNAElastic* theDNAElasticProcess =
286 new G4DNAElastic("hydrogen_G4DNAElastic");
287 theDNAElasticProcess->SetEmModel(new G4DNAIonElasticModel());
288 ((G4DNAIonElasticModel*)(theDNAElasticProcess->EmModel()))
289 ->SelectStationary(true);
290 ph->RegisterProcess(theDNAElasticProcess, particle);
291
292 // *** Excitation ***
293
294 G4DNAExcitation* theDNAExcitationProcess =
295 new G4DNAExcitation("hydrogen_G4DNAExcitation");
296 theDNAExcitationProcess
298 ((G4DNAMillerGreenExcitationModel*)(theDNAExcitationProcess->EmModel()))
299 ->SelectStationary(true);
300 ph->RegisterProcess(theDNAExcitationProcess, particle);
301
302 // *** Ionisation ***
303
304 G4DNAIonisation* theDNAIonisationProcess =
305 new G4DNAIonisation("hydrogen_G4DNAIonisation");
306 theDNAIonisationProcess->SetEmModel(new G4DNARuddIonisationModel());
307 ((G4DNARuddIonisationModel*)(theDNAIonisationProcess->EmModel()))
308 ->SelectStationary(true);
309 ph->RegisterProcess(theDNAIonisationProcess, particle);
310
311 // *** Charge increase ***
312
313 G4DNAChargeIncrease* theDNAChargeIncreaseProcess =
314 new G4DNAChargeIncrease("hydrogen_G4DNAChargeIncrease");
315 theDNAChargeIncreaseProcess
318 (theDNAChargeIncreaseProcess->EmModel()))
319 ->SelectStationary(true);
320 ph->RegisterProcess(theDNAChargeIncreaseProcess, particle);
321
322 } else if ( particleName == "alpha" ) {
323
324 // *** Elastic ***
325
326 G4DNAElastic* theDNAElasticProcess =
327 new G4DNAElastic("alpha_G4DNAElastic");
328 theDNAElasticProcess->SetEmModel(new G4DNAIonElasticModel());
329 ((G4DNAIonElasticModel*)(theDNAElasticProcess->EmModel()))
330 ->SelectStationary(true);
331 ph->RegisterProcess(theDNAElasticProcess, particle);
332
333 // *** Excitation ***
334
335 G4DNAExcitation* theDNAExcitationProcess =
336 new G4DNAExcitation("alpha_G4DNAExcitation");
337 theDNAExcitationProcess->SetEmModel
339 ((G4DNAMillerGreenExcitationModel*)(theDNAExcitationProcess->EmModel()))
340 ->SelectStationary(true);
341 ph->RegisterProcess(theDNAExcitationProcess, particle);
342
343 // *** Ionisation ***
344
345 G4DNAIonisation* theDNAIonisationProcess =
346 new G4DNAIonisation("alpha_G4DNAIonisation");
347 theDNAIonisationProcess->SetEmModel(new G4DNARuddIonisationModel());
348 ((G4DNARuddIonisationModel*)(theDNAIonisationProcess->EmModel()))
349 ->SelectStationary(true);
350 ph->RegisterProcess(theDNAIonisationProcess, particle);
351
352 // *** Charge decrease ***
353
354 G4DNAChargeDecrease* theDNAChargeDecreaseProcess =
355 new G4DNAChargeDecrease("alpha_G4DNAChargeDecrease");
356 theDNAChargeDecreaseProcess->SetEmModel
359 (theDNAChargeDecreaseProcess->EmModel()))
360 ->SelectStationary(true);
361 ph->RegisterProcess(theDNAChargeDecreaseProcess, particle);
362
363 } else if ( particleName == "alpha+" ) {
364
365 // *** Elastic ***
366
367 G4DNAElastic* theDNAElasticProcess =
368 new G4DNAElastic("alpha+_G4DNAElastic");
369 theDNAElasticProcess->SetEmModel(new G4DNAIonElasticModel());
370 ((G4DNAIonElasticModel*)(theDNAElasticProcess->EmModel()))
371 ->SelectStationary(true);
372 ph->RegisterProcess(theDNAElasticProcess, particle);
373
374 // *** Excitation ***
375
376 G4DNAExcitation* theDNAExcitationProcess =
377 new G4DNAExcitation("alpha+_G4DNAExcitation");
378 theDNAExcitationProcess->SetEmModel
381 (theDNAExcitationProcess->EmModel()))->SelectStationary(true);
382 ph->RegisterProcess(theDNAExcitationProcess, particle);
383
384 // *** Ionisation ***
385
386 G4DNAIonisation* theDNAIonisationProcess =
387 new G4DNAIonisation("alpha+_G4DNAIonisation");
388 theDNAIonisationProcess->SetEmModel(new G4DNARuddIonisationModel());
389 ((G4DNARuddIonisationModel*)(theDNAIonisationProcess->EmModel()))
390 ->SelectStationary(true);
391 ph->RegisterProcess(theDNAIonisationProcess, particle);
392
393 // *** Charge decrease ***
394
395 G4DNAChargeDecrease* theDNAChargeDecreaseProcess =
396 new G4DNAChargeDecrease("alpha+_G4DNAChargeDecrease");
397 theDNAChargeDecreaseProcess->SetEmModel
400 (theDNAChargeDecreaseProcess->EmModel()))->SelectStationary(true);
401 ph->RegisterProcess(theDNAChargeDecreaseProcess, particle);
402
403 // *** Charge increase ***
404
405 G4DNAChargeIncrease* theDNAChargeIncreaseProcess =
406 new G4DNAChargeIncrease("alpha+_G4DNAChargeIncrease");
407 theDNAChargeIncreaseProcess->SetEmModel
410 (theDNAChargeIncreaseProcess->EmModel()))->SelectStationary(true);
411 ph->RegisterProcess(theDNAChargeIncreaseProcess, particle);
412
413 } else if ( particleName == "helium" ) {
414
415 // *** Elastic ***
416
417 G4DNAElastic* theDNAElasticProcess =
418 new G4DNAElastic("helium_G4DNAElastic");
419 theDNAElasticProcess->SetEmModel
420 (new G4DNAIonElasticModel());
421 ((G4DNAIonElasticModel*)(theDNAElasticProcess->EmModel()))
422 ->SelectStationary(true);
423 ph->RegisterProcess(theDNAElasticProcess, particle);
424
425 // *** Excitation ***
426
427 G4DNAExcitation* theDNAExcitationProcess =
428 new G4DNAExcitation("helium_G4DNAExcitation");
429 theDNAExcitationProcess->SetEmModel
431 ((G4DNAMillerGreenExcitationModel*)(theDNAExcitationProcess->EmModel()))
432 ->SelectStationary(true);
433 ph->RegisterProcess(theDNAExcitationProcess, particle);
434
435 // *** Ionisation ***
436
437 G4DNAIonisation* theDNAIonisationProcess =
438 new G4DNAIonisation("helium_G4DNAIonisation");
439 theDNAIonisationProcess->SetEmModel
441 ((G4DNARuddIonisationModel*)(theDNAIonisationProcess->EmModel()))
442 ->SelectStationary(true);
443 ph->RegisterProcess(theDNAIonisationProcess, particle);
444
445 // *** Charge increase ***
446
447 G4DNAChargeIncrease* theDNAChargeIncreaseProcess =
448 new G4DNAChargeIncrease("helium_G4DNAChargeIncrease");
449 theDNAChargeIncreaseProcess->SetEmModel
452 (theDNAChargeIncreaseProcess->EmModel()))
453 ->SelectStationary(true);
454 ph->RegisterProcess(theDNAChargeIncreaseProcess, particle);
455
456 } else if ( particleName == "GenericIon" ) {
457
458 // *** Ionisation ***
459
460 G4DNAIonisation* theDNAIonisationProcess =
461 new G4DNAIonisation("GenericIon_G4DNAIonisation");
462 theDNAIonisationProcess->SetEmModel(
464 ((G4DNARuddIonisationExtendedModel*)(theDNAIonisationProcess->EmModel()))
465 ->SelectStationary(true);
466 ph->RegisterProcess(theDNAIonisationProcess, particle);
467
468 }
469
470 // Warning : the following particles and processes are needed by EM Physics
471 // builders
472 // They are taken from the default Livermore Physics list
473 // These particles are currently not handled by Geant4-DNA
474
475 // e+
476
477 else if (particleName == "e+") {
478
479 // Identical to G4EmStandardPhysics_stationary
480
483 G4eIonisation* eIoni = new G4eIonisation();
484 eIoni->SetStepFunction(0.2, 100*um);
485
486 ph->RegisterProcess(msc, particle);
487 ph->RegisterProcess(eIoni, particle);
488 ph->RegisterProcess(new G4eBremsstrahlung(), particle);
489 ph->RegisterProcess(new G4eplusAnnihilation(), particle);
490
491 } else if (particleName == "gamma") {
492 // photoelectric effect - Livermore model only
493 G4PhotoElectricEffect* thePhotoElectricEffect = new G4PhotoElectricEffect();
494 thePhotoElectricEffect->SetEmModel(new G4LivermorePhotoElectricModel());
495 ph->RegisterProcess(thePhotoElectricEffect, particle);
496
497 // Compton scattering - Livermore model only
498 G4ComptonScattering* theComptonScattering = new G4ComptonScattering();
499 theComptonScattering->SetEmModel(new G4LivermoreComptonModel());
500 ph->RegisterProcess(theComptonScattering, particle);
501
502 // gamma conversion - Livermore model below 80 GeV
503 G4GammaConversion* theGammaConversion = new G4GammaConversion();
504 theGammaConversion->SetEmModel(new G4LivermoreGammaConversionModel());
505 ph->RegisterProcess(theGammaConversion, particle);
506
507 // default Rayleigh scattering is Livermore
508 G4RayleighScattering* theRayleigh = new G4RayleighScattering();
509 ph->RegisterProcess(theRayleigh, particle);
510 }
511
512 // Warning : end of particles and processes are needed by EM Physics build.
513
514 }
515
516 // Deexcitation
517 //
520}
521
@ bElectromagnetic
G4DNABornExcitationModel1 G4DNABornExcitationModel
#define G4DNABornIonisationModel
@ fUseDistanceToBoundary
#define G4_DECLARE_PHYSCONSTR_FACTORY(physics_constructor)
int G4int
Definition: G4Types.hh:85
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
static G4DNAGenericIonsManager * Instance(void)
G4ParticleDefinition * GetIon(const G4String &name)
static G4Electron * Electron()
Definition: G4Electron.cc:93
static G4EmParameters * Instance()
void SetDeexcitationIgnoreCut(G4bool val)
void SetFluo(G4bool val)
void SetAugerCascade(G4bool val)
void SetAuger(G4bool val)
static G4Gamma * Gamma()
Definition: G4Gamma.cc:85
static G4GenericIon * GenericIonDefinition()
Definition: G4GenericIon.cc:87
void SetAtomDeexcitation(G4VAtomDeexcitation *)
static G4LossTableManager * Instance()
const G4String & GetParticleName() const
G4bool RegisterProcess(G4VProcess *process, G4ParticleDefinition *particle)
static G4PhysicsListHelper * GetPhysicsListHelper()
static G4Positron * Positron()
Definition: G4Positron.cc:93
static G4Proton * Proton()
Definition: G4Proton.cc:92
G4VEmModel * EmModel(size_t index=0) const
void SetEmModel(G4VEmModel *, G4int index=0)
void SetStepFunction(G4double v1, G4double v2)
void SetStepLimitType(G4MscStepLimitType val)
G4ParticleTable::G4PTblDicIterator * GetParticleIterator() const
const G4String & GetPhysicsName() const