Geant4 10.7.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4EmDNAPhysics_stationary_option6 Class Reference

#include <G4EmDNAPhysics_stationary_option6.hh>

+ Inheritance diagram for G4EmDNAPhysics_stationary_option6:

Public Member Functions

 G4EmDNAPhysics_stationary_option6 (G4int ver=1)
 
 G4EmDNAPhysics_stationary_option6 (G4int ver, const G4String &name)
 
virtual ~G4EmDNAPhysics_stationary_option6 ()
 
virtual void ConstructParticle ()
 
virtual void ConstructProcess ()
 
- Public Member Functions inherited from G4VPhysicsConstructor
 G4VPhysicsConstructor (const G4String &="")
 
 G4VPhysicsConstructor (const G4String &name, G4int physics_type)
 
virtual ~G4VPhysicsConstructor ()
 
virtual void ConstructParticle ()=0
 
virtual void ConstructProcess ()=0
 
void SetPhysicsName (const G4String &="")
 
const G4StringGetPhysicsName () const
 
void SetPhysicsType (G4int)
 
G4int GetPhysicsType () const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
G4int GetInstanceID () const
 
virtual void TerminateWorker ()
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VPhysicsConstructor
static const G4VPCManagerGetSubInstanceManager ()
 
- Protected Types inherited from G4VPhysicsConstructor
using PhysicsBuilder_V = G4VPCData::PhysicsBuilders_V
 
- Protected Member Functions inherited from G4VPhysicsConstructor
G4bool RegisterProcess (G4VProcess *process, G4ParticleDefinition *particle)
 
G4ParticleTable::G4PTblDicIteratorGetParticleIterator () const
 
PhysicsBuilder_V GetBuilders () const
 
void AddBuilder (G4PhysicsBuilderInterface *bld)
 
- Protected Attributes inherited from G4VPhysicsConstructor
G4int verboseLevel
 
G4String namePhysics
 
G4int typePhysics
 
G4ParticleTabletheParticleTable
 
G4int g4vpcInstanceID
 
- Static Protected Attributes inherited from G4VPhysicsConstructor
static G4RUN_DLL G4VPCManager subInstanceManager
 

Detailed Description

Definition at line 34 of file G4EmDNAPhysics_stationary_option6.hh.

Constructor & Destructor Documentation

◆ G4EmDNAPhysics_stationary_option6() [1/2]

G4EmDNAPhysics_stationary_option6::G4EmDNAPhysics_stationary_option6 ( G4int  ver = 1)

Definition at line 89 of file G4EmDNAPhysics_stationary_option6.cc.

90 : G4VPhysicsConstructor("G4EmDNAPhysics_stationary_option6"), verbose(ver)
91{
93 param->SetDefaults();
94 param->SetFluo(true);
95 param->SetAuger(true);
96 param->SetAugerCascade(true);
97 param->SetDeexcitationIgnoreCut(true);
98 param->ActivateDNA();
99
101}
@ bElectromagnetic
static G4EmParameters * Instance()
void SetDeexcitationIgnoreCut(G4bool val)
void SetFluo(G4bool val)
void SetAugerCascade(G4bool val)
void SetAuger(G4bool val)

◆ G4EmDNAPhysics_stationary_option6() [2/2]

G4EmDNAPhysics_stationary_option6::G4EmDNAPhysics_stationary_option6 ( G4int  ver,
const G4String name 
)

Definition at line 105 of file G4EmDNAPhysics_stationary_option6.cc.

107 : G4VPhysicsConstructor("G4EmDNAPhysics_stationary_option6"), verbose(ver)
108{
110 param->SetDefaults();
111 param->SetFluo(true);
112 param->SetAuger(true);
113 param->SetAugerCascade(true);
114 param->SetDeexcitationIgnoreCut(true);
115
117}

◆ ~G4EmDNAPhysics_stationary_option6()

G4EmDNAPhysics_stationary_option6::~G4EmDNAPhysics_stationary_option6 ( )
virtual

Definition at line 121 of file G4EmDNAPhysics_stationary_option6.cc.

122{}

Member Function Documentation

◆ ConstructParticle()

void G4EmDNAPhysics_stationary_option6::ConstructParticle ( )
virtual

Implements G4VPhysicsConstructor.

Definition at line 126 of file G4EmDNAPhysics_stationary_option6.cc.

127{
128// bosons
130
131// leptons
134
135// baryons
137
139
140 G4DNAGenericIonsManager * genericIonsManager;
141 genericIonsManager=G4DNAGenericIonsManager::Instance();
142 genericIonsManager->GetIon("alpha++");
143 genericIonsManager->GetIon("alpha+");
144 genericIonsManager->GetIon("helium");
145 genericIonsManager->GetIon("hydrogen");
146
147}
static G4DNAGenericIonsManager * Instance(void)
G4ParticleDefinition * GetIon(const G4String &name)
static G4Electron * Electron()
Definition: G4Electron.cc:93
static G4Gamma * Gamma()
Definition: G4Gamma.cc:85
static G4GenericIon * GenericIonDefinition()
Definition: G4GenericIon.cc:87
static G4Positron * Positron()
Definition: G4Positron.cc:93
static G4Proton * Proton()
Definition: G4Proton.cc:92

◆ ConstructProcess()

void G4EmDNAPhysics_stationary_option6::ConstructProcess ( )
virtual

Implements G4VPhysicsConstructor.

Definition at line 151 of file G4EmDNAPhysics_stationary_option6.cc.

152{
153 if(verbose > 1) {
154 G4cout << "### " << GetPhysicsName() << " Construct Processes " << G4endl;
155 }
157
158 auto myParticleIterator=GetParticleIterator();
159 myParticleIterator->reset();
160 while( (*myParticleIterator)() )
161 {
162 G4ParticleDefinition* particle = myParticleIterator->value();
163 G4String particleName = particle->GetParticleName();
164
165 if (particleName == "e-") {
166
167 // *** Elastic scattering (two alternative models available) ***
168
169 G4DNAElastic* theDNAElasticProcess =
170 new G4DNAElastic("e-_G4DNAElastic");
171 theDNAElasticProcess->SetEmModel(new G4DNACPA100ElasticModel());
173 (theDNAElasticProcess->EmModel()))->SelectStationary(true);
174 ph->RegisterProcess(theDNAElasticProcess, particle);
175
176 // *** Excitation ***
177
178 G4DNAExcitation* theDNAExcitationProcess =
179 new G4DNAExcitation("e-_G4DNAExcitation");
180 theDNAExcitationProcess->SetEmModel(new G4DNACPA100ExcitationModel());
182 (theDNAExcitationProcess->EmModel()))->SelectStationary(true);
183 ph->RegisterProcess(theDNAExcitationProcess, particle);
184
185 // *** Ionisation ***
186
187 G4DNAIonisation* theDNAIonisationProcess =
188 new G4DNAIonisation("e-_G4DNAIonisation");
189 theDNAIonisationProcess->SetEmModel(new G4DNACPA100IonisationModel());
191 (theDNAIonisationProcess->EmModel()))->SelectStationary(true);
192 ph->RegisterProcess(theDNAIonisationProcess, particle);
193
194 // *** Vibrational excitation ***
195
196 G4DNAVibExcitation* theDNAVibExcitationProcess =
197 new G4DNAVibExcitation("e-_G4DNAVibExcitation");
198 theDNAVibExcitationProcess->SetEmModel(new G4DNASancheExcitationModel());
200 (theDNAVibExcitationProcess->EmModel()))->SelectStationary(true);
201 ph->RegisterProcess(theDNAVibExcitationProcess, particle);
202
203 // *** Attachment ***
204
205 G4DNAAttachment* theDNAAttachmentProcess =
206 new G4DNAAttachment("e-_G4DNAAttachment");
207 theDNAAttachmentProcess->SetEmModel(new G4DNAMeltonAttachmentModel());
209 (theDNAAttachmentProcess->EmModel()))->SelectStationary(true);
210 ph->RegisterProcess(theDNAAttachmentProcess, particle);
211
212 } else if ( particleName == "proton" ) {
213
214 // *** Elastic ***
215
216 G4DNAElastic* theDNAElasticProcess =
217 new G4DNAElastic("proton_G4DNAElastic");
218 theDNAElasticProcess->SetEmModel(new G4DNAIonElasticModel());
220 (theDNAElasticProcess->EmModel()))->SelectStationary(true);
221 ph->RegisterProcess(theDNAElasticProcess, particle);
222
223 // *** Excitation ***
224
225 G4DNAExcitation* theDNAExcitationProcess =
226 new G4DNAExcitation("proton_G4DNAExcitation");
227
228 theDNAExcitationProcess->SetEmModel(
230 theDNAExcitationProcess->SetEmModel(
232
234 (theDNAExcitationProcess->EmModel()))->SetLowEnergyLimit(10*eV);
236 (theDNAExcitationProcess->EmModel()))->SetHighEnergyLimit(500*keV);
238 (theDNAExcitationProcess->EmModel()))->SelectStationary(true);
239
241 (theDNAExcitationProcess->EmModel(1)))->SetLowEnergyLimit(500*keV);
243 (theDNAExcitationProcess->EmModel(1)))->SetHighEnergyLimit(100*MeV);
245 (theDNAExcitationProcess->EmModel(1)))->SelectStationary(true);
246
247 ph->RegisterProcess(theDNAExcitationProcess, particle);
248
249 // *** Ionisation ***
250
251 G4DNAIonisation* theDNAIonisationProcess =
252 new G4DNAIonisation("proton_G4DNAIonisation");
253
254 theDNAIonisationProcess->SetEmModel(
256 theDNAIonisationProcess->SetEmModel(
258
260 (theDNAIonisationProcess->EmModel()))->SetLowEnergyLimit(0*eV);
262 (theDNAIonisationProcess->EmModel()))->SetHighEnergyLimit(500*keV);
264 (theDNAIonisationProcess->EmModel()))->SelectStationary(true);
265
267 (theDNAIonisationProcess->EmModel(1)))->SetLowEnergyLimit(500*keV);
269 (theDNAIonisationProcess->EmModel(1)))->SetHighEnergyLimit(100*MeV);
271 (theDNAIonisationProcess->EmModel(1)))->SelectStationary(true);
272 //
274 (theDNAIonisationProcess->EmModel(1)))->SelectFasterComputation(true);
275 //
276
277 ph->RegisterProcess(theDNAIonisationProcess, particle);
278
279 // *** Charge decrease ***
280
281 G4DNAChargeDecrease* theDNAChargeDecreaseProcess =
282 new G4DNAChargeDecrease("proton_G4DNAChargeDecrease");
283 theDNAChargeDecreaseProcess->SetEmModel(
286 (theDNAChargeDecreaseProcess->EmModel()))->SelectStationary(true);
287 ph->RegisterProcess(theDNAChargeDecreaseProcess, particle);
288
289 } else if ( particleName == "hydrogen" ) {
290
291 // *** Elastic ***
292
293 G4DNAElastic* theDNAElasticProcess =
294 new G4DNAElastic("hydrogen_G4DNAElastic");
295 theDNAElasticProcess->SetEmModel(
298 (theDNAElasticProcess->EmModel()))->SelectStationary(true);
299 ph->RegisterProcess(theDNAElasticProcess, particle);
300
301 // *** Excitation ***
302
303 G4DNAExcitation* theDNAExcitationProcess =
304 new G4DNAExcitation("hydrogen_G4DNAExcitation");
305 theDNAExcitationProcess->SetEmModel(
308 (theDNAExcitationProcess->EmModel()))->SelectStationary(true);
309 ph->RegisterProcess(theDNAExcitationProcess, particle);
310
311 // *** Ionisation ***
312
313 G4DNAIonisation* theDNAIonisationProcess =
314 new G4DNAIonisation("hydrogen_G4DNAIonisation");
315 theDNAIonisationProcess->SetEmModel(
318 theDNAIonisationProcess->EmModel()))->SelectStationary(true);
319 ph->RegisterProcess(theDNAIonisationProcess, particle);
320
321 // *** Charge increase ***
322
323 G4DNAChargeIncrease* theDNAChargeIncreaseProcess =
324 new G4DNAChargeIncrease("hydrogen_G4DNAChargeIncrease");
325 theDNAChargeIncreaseProcess->SetEmModel(
328 theDNAChargeIncreaseProcess->EmModel()))->SelectStationary(true);
329 ph->RegisterProcess(theDNAChargeIncreaseProcess, particle);
330
331 } else if ( particleName == "alpha" ) {
332
333 // *** Elastic ***
334
335 G4DNAElastic* theDNAElasticProcess =
336 new G4DNAElastic("alpha_G4DNAElastic");
337 theDNAElasticProcess->SetEmModel(new G4DNAIonElasticModel());
339 (theDNAElasticProcess->EmModel()))->SelectStationary(true);
340 ph->RegisterProcess(theDNAElasticProcess, particle);
341
342 // *** Excitation ***
343
344 G4DNAExcitation* theDNAExcitationProcess =
345 new G4DNAExcitation("alpha_G4DNAExcitation");
346 theDNAExcitationProcess->SetEmModel(
349 (theDNAExcitationProcess->EmModel()))->SelectStationary(true);
350 ph->RegisterProcess(theDNAExcitationProcess, particle);
351
352 // *** Ionisation ***
353
354 G4DNAIonisation* theDNAIonisationProcess =
355 new G4DNAIonisation("alpha_G4DNAIonisation");
356 theDNAIonisationProcess->SetEmModel(
359 (theDNAIonisationProcess->EmModel()))->SelectStationary(true);
360 ph->RegisterProcess(theDNAIonisationProcess, particle);
361
362 // *** Charge decrease ***
363
364 G4DNAChargeDecrease* theDNAChargeDecreaseProcess =
365 new G4DNAChargeDecrease("alpha_G4DNAChargeDecrease");
366 theDNAChargeDecreaseProcess->SetEmModel(
369 (theDNAChargeDecreaseProcess->EmModel()))->SelectStationary(true);
370 ph->RegisterProcess(theDNAChargeDecreaseProcess, particle);
371
372 } else if ( particleName == "alpha+" ) {
373
374 // *** Elastic ***
375
376 G4DNAElastic* theDNAElasticProcess =
377 new G4DNAElastic("alpha+_G4DNAElastic");
378 theDNAElasticProcess->SetEmModel(
381 (theDNAElasticProcess->EmModel()))->SelectStationary(true);
382 ph->RegisterProcess(theDNAElasticProcess, particle);
383
384 // *** Excitation ***
385
386 G4DNAExcitation* theDNAExcitationProcess =
387 new G4DNAExcitation("alpha+_G4DNAExcitation");
388 theDNAExcitationProcess->SetEmModel(
391 (theDNAExcitationProcess->EmModel()))->SelectStationary(true);
392 ph->RegisterProcess(theDNAExcitationProcess, particle);
393
394 // *** Ionisation ***
395
396 G4DNAIonisation* theDNAIonisationProcess =
397 new G4DNAIonisation("alpha+_G4DNAIonisation");
398 theDNAIonisationProcess->SetEmModel(
401 (theDNAIonisationProcess->EmModel()))->SelectStationary(true);
402 ph->RegisterProcess(theDNAIonisationProcess, particle);
403
404 // *** Charge decrease ***
405
406 G4DNAChargeDecrease* theDNAChargeDecreaseProcess =
407 new G4DNAChargeDecrease("alpha+_G4DNAChargeDecrease");
408 theDNAChargeDecreaseProcess->SetEmModel(
411 (theDNAChargeDecreaseProcess->EmModel()))->SelectStationary(true);
412 ph->RegisterProcess(theDNAChargeDecreaseProcess, particle);
413
414 // *** Charge increase ***
415
416 G4DNAChargeIncrease* theDNAChargeIncreaseProcess =
417 new G4DNAChargeIncrease("alpha+_G4DNAChargeIncrease");
418 theDNAChargeIncreaseProcess->SetEmModel(
421 (theDNAChargeIncreaseProcess->EmModel()))->SelectStationary(true);
422 ph->RegisterProcess(theDNAChargeIncreaseProcess, particle);
423
424 } else if ( particleName == "helium" ) {
425
426 // *** Elastic ***
427
428 G4DNAElastic* theDNAElasticProcess =
429 new G4DNAElastic("helium_G4DNAElastic");
430 theDNAElasticProcess->SetEmModel(new G4DNAIonElasticModel());
432 (theDNAElasticProcess->EmModel()))->SelectStationary(true);
433 ph->RegisterProcess(theDNAElasticProcess, particle);
434
435 // *** Excitation ***
436
437 G4DNAExcitation* theDNAExcitationProcess =
438 new G4DNAExcitation("helium_G4DNAExcitation");
439 theDNAExcitationProcess->SetEmModel(
442 (theDNAExcitationProcess->EmModel()))->SelectStationary(true);
443 ph->RegisterProcess(theDNAExcitationProcess, particle);
444
445 // *** Ionisation ***
446
447 G4DNAIonisation* theDNAIonisationProcess =
448 new G4DNAIonisation("helium_G4DNAIonisation");
449 theDNAIonisationProcess->SetEmModel(
452 (theDNAIonisationProcess->EmModel()))->SelectStationary(true);
453 ph->RegisterProcess(theDNAIonisationProcess, particle);
454
455 // *** Charge increase ***
456
457 G4DNAChargeIncrease* theDNAChargeIncreaseProcess =
458 new G4DNAChargeIncrease("helium_G4DNAChargeIncrease");
459 theDNAChargeIncreaseProcess->SetEmModel(
462 (theDNAChargeIncreaseProcess->EmModel()))->SelectStationary(true);
463 ph->RegisterProcess(theDNAChargeIncreaseProcess, particle);
464
465 } else if ( particleName == "GenericIon" ) {
466
467 // *** Ionisation ***
468
469 G4DNAIonisation* theDNAIonisationProcess =
470 new G4DNAIonisation("GenericIon_G4DNAIonisation");
471 theDNAIonisationProcess->SetEmModel(
474 (theDNAIonisationProcess->EmModel()))->SelectStationary(true);
475 ph->RegisterProcess(theDNAIonisationProcess, particle);
476
477 }
478
479 // Warning : the following particles and processes are needed by EM Physics
480 // builders
481 // They are taken from the default Livermore Physics list
482 // These particles are currently not handled by Geant4-DNA
483
484 // e+
485
486 else if (particleName == "e+") {
487
488 // Identical to G4EmStandardPhysics_stationary
489
492 G4eIonisation* eIoni = new G4eIonisation();
493 eIoni->SetStepFunction(0.2, 100*um);
494
495 ph->RegisterProcess(msc, particle);
496 ph->RegisterProcess(eIoni, particle);
497 ph->RegisterProcess(new G4eBremsstrahlung(), particle);
498 ph->RegisterProcess(new G4eplusAnnihilation(), particle);
499
500 } else if (particleName == "gamma") {
501
502 // photoelectric effect - Livermore model only
503 G4PhotoElectricEffect* thePhotoElectricEffect = new G4PhotoElectricEffect();
504 thePhotoElectricEffect->SetEmModel(new G4LivermorePhotoElectricModel());
505 ph->RegisterProcess(thePhotoElectricEffect, particle);
506
507 // Compton scattering - Livermore model only
508 G4ComptonScattering* theComptonScattering = new G4ComptonScattering();
509 theComptonScattering->SetEmModel(new G4LivermoreComptonModel());
510 ph->RegisterProcess(theComptonScattering, particle);
511
512 // gamma conversion - Livermore model below 80 GeV
513 G4GammaConversion* theGammaConversion = new G4GammaConversion();
514 theGammaConversion->SetEmModel(new G4LivermoreGammaConversionModel());
515 ph->RegisterProcess(theGammaConversion, particle);
516
517 // default Rayleigh scattering is Livermore
518 G4RayleighScattering* theRayleigh = new G4RayleighScattering();
519 ph->RegisterProcess(theRayleigh, particle);
520 }
521
522 // Warning : end of particles and processes are needed by EM Physics build.
523
524 }
525
526 // Deexcitation
527 //
530}
G4DNABornExcitationModel1 G4DNABornExcitationModel
#define G4DNABornIonisationModel
@ fUseDistanceToBoundary
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
void SetAtomDeexcitation(G4VAtomDeexcitation *)
static G4LossTableManager * Instance()
const G4String & GetParticleName() const
G4bool RegisterProcess(G4VProcess *process, G4ParticleDefinition *particle)
static G4PhysicsListHelper * GetPhysicsListHelper()
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

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