Geant4
9.6.p02
Main Page
Related Pages
Modules
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Groups
Pages
geant4_9_6_p02
examples
advanced
microdosimetry
src
PhysicsList.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
//
26
// -------------------------------------------------------------------
27
// $Id$
28
// -------------------------------------------------------------------
29
30
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
31
32
#include "PhysicsList.hh"
33
#include "
G4SystemOfUnits.hh
"
34
35
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
36
37
PhysicsList::PhysicsList
():
G4VUserPhysicsList
()
38
{
39
defaultCutValue
= 1*
micrometer
;
40
cutForGamma =
defaultCutValue
;
41
cutForElectron =
defaultCutValue
;
42
cutForPositron =
defaultCutValue
;
43
cutForProton =
defaultCutValue
;
44
45
SetVerboseLevel
(1);
46
}
47
48
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
49
50
PhysicsList::~PhysicsList
()
51
{}
52
53
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
54
55
void
PhysicsList::ConstructParticle
()
56
{
57
ConstructBosons
();
58
ConstructLeptons
();
59
ConstructBarions
();
60
}
61
62
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
63
64
void
PhysicsList::ConstructBosons
()
65
{
66
// gamma
67
G4Gamma::GammaDefinition
();
68
}
69
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
70
71
void
PhysicsList::ConstructLeptons
()
72
{
73
// leptons
74
G4Electron::ElectronDefinition
();
75
G4Positron::PositronDefinition
();
76
}
77
78
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
79
80
//DNA
81
#include "
G4DNAGenericIonsManager.hh
"
82
//ENDDNA
83
84
void
PhysicsList::ConstructBarions
()
85
{
86
// baryons
87
G4Proton::ProtonDefinition
();
88
G4GenericIon::GenericIonDefinition
();
89
90
// Geant4 DNA new particles
91
G4DNAGenericIonsManager
* genericIonsManager;
92
genericIonsManager=
G4DNAGenericIonsManager::Instance
();
93
genericIonsManager->
GetIon
(
"alpha++"
);
94
genericIonsManager->
GetIon
(
"alpha+"
);
95
genericIonsManager->
GetIon
(
"helium"
);
96
genericIonsManager->
GetIon
(
"hydrogen"
);
97
}
98
99
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
100
101
void
PhysicsList::ConstructProcess
()
102
{
103
AddTransportation
();
104
ConstructEM();
105
ConstructGeneral
();
106
}
107
108
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
109
110
// Geant4-DNA MODELS
111
112
#include "
G4DNAElastic.hh
"
113
#include "
G4DNAChampionElasticModel.hh
"
114
#include "
G4DNAScreenedRutherfordElasticModel.hh
"
115
116
#include "
G4DNAExcitation.hh
"
117
#include "
G4DNAMillerGreenExcitationModel.hh
"
118
#include "
G4DNABornExcitationModel.hh
"
119
120
#include "
G4DNAIonisation.hh
"
121
#include "
G4DNABornIonisationModel.hh
"
122
#include "
G4DNARuddIonisationModel.hh
"
123
124
#include "
G4DNAChargeDecrease.hh
"
125
#include "
G4DNADingfelderChargeDecreaseModel.hh
"
126
127
#include "
G4DNAChargeIncrease.hh
"
128
#include "
G4DNADingfelderChargeIncreaseModel.hh
"
129
130
#include "
G4DNAAttachment.hh
"
131
#include "
G4DNAMeltonAttachmentModel.hh
"
132
133
#include "
G4DNAVibExcitation.hh
"
134
#include "
G4DNASancheExcitationModel.hh
"
135
136
//
137
138
#include "
G4LossTableManager.hh
"
139
#include "
G4EmConfigurator.hh
"
140
#include "
G4VEmModel.hh
"
141
#include "
G4DummyModel.hh
"
142
#include "
G4eIonisation.hh
"
143
#include "
G4hIonisation.hh
"
144
#include "
G4ionIonisation.hh
"
145
#include "
G4eMultipleScattering.hh
"
146
#include "
G4hMultipleScattering.hh
"
147
#include "
G4BraggIonGasModel.hh
"
148
#include "
G4BetheBlochIonGasModel.hh
"
149
#include "
G4UrbanMscModel93.hh
"
150
#include "
G4MollerBhabhaModel.hh
"
151
#include "
G4IonFluctuations.hh
"
152
#include "
G4UniversalFluctuation.hh
"
153
154
#include "
G4ElectronCapture.hh
"
155
156
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
157
158
void
PhysicsList::ConstructEM()
159
{
160
161
theParticleIterator
->
reset
();
162
163
while
( (*
theParticleIterator
)() )
164
{
165
166
G4ParticleDefinition
* particle =
theParticleIterator
->
value
();
167
G4ProcessManager
* pmanager = particle->
GetProcessManager
();
168
G4String
particleName = particle->
GetParticleName
();
169
170
// *********************************
171
// 1) Processes for the World region
172
// *********************************
173
174
if
(particleName ==
"e-"
) {
175
176
// STANDARD msc is active in the world
177
G4eMultipleScattering
* msc =
new
G4eMultipleScattering
();
178
pmanager->
AddProcess
(msc, -1, 1, 1);
179
180
// STANDARD ionisation is active in the world
181
G4eIonisation
* eion =
new
G4eIonisation
();
182
eion->
SetEmModel
(
new
G4MollerBhabhaModel
(), 1);
183
pmanager->
AddProcess
(eion, -1, 2, 2);
184
185
// DNA elastic is not active in the world
186
G4DNAElastic
* theDNAElasticProcess =
new
G4DNAElastic
(
"e-_G4DNAElastic"
);
187
theDNAElasticProcess->
SetEmModel
(
new
G4DummyModel
(),1);
188
pmanager->
AddDiscreteProcess
(theDNAElasticProcess);
189
190
// DNA excitation is not active in the world
191
G4DNAExcitation
* dnaex =
new
G4DNAExcitation
(
"e-_G4DNAExcitation"
);
192
dnaex->
SetEmModel
(
new
G4DummyModel
(),1);
193
pmanager->
AddDiscreteProcess
(dnaex);
194
195
// DNA ionisation is not active in the world
196
G4DNAIonisation
* dnaioni =
new
G4DNAIonisation
(
"e-_G4DNAIonisation"
);
197
dnaioni->
SetEmModel
(
new
G4DummyModel
(),1);
198
pmanager->
AddDiscreteProcess
(dnaioni);
199
200
// DNA attachment is not active in the world
201
G4DNAAttachment
* dnaatt =
new
G4DNAAttachment
(
"e-_G4DNAAttachment"
);
202
dnaatt->
SetEmModel
(
new
G4DummyModel
(),1);
203
pmanager->
AddDiscreteProcess
(dnaatt);
204
205
// DNA vib. excitation is not active in the world
206
G4DNAVibExcitation
* dnavib =
new
G4DNAVibExcitation
(
"e-_G4DNAVibExcitation"
);
207
dnavib->
SetEmModel
(
new
G4DummyModel
(),1);
208
pmanager->
AddDiscreteProcess
(dnavib);
209
210
// THE FOLLOWING PROCESS WILL KILL ALL ELECTRONS BELOW A SELECTED ENERY THRESHOLD
211
// Capture of low-energy e-
212
G4ElectronCapture
* ecap =
new
G4ElectronCapture
(
"Target"
, 5.1*
eV
);
213
pmanager->
AddDiscreteProcess
(ecap);
214
215
}
else
if
( particleName ==
"proton"
) {
216
217
// STANDARD msc is active in the world
218
G4hMultipleScattering
* msc =
new
G4hMultipleScattering
();
219
pmanager->
AddProcess
(msc, -1, 1, 1);
220
221
// STANDARD ionisation is active in the world
222
G4hIonisation
* hion =
new
G4hIonisation
();
223
hion->
SetEmModel
(
new
G4BraggIonGasModel
(), 1);
224
hion->
SetEmModel
(
new
G4BetheBlochIonGasModel
(), 2);
225
pmanager->
AddProcess
(hion, -1, 2, 2);
226
227
// DNA excitation is not active in the world
228
G4DNAExcitation
* dnaex =
new
G4DNAExcitation
(
"proton_G4DNAExcitation"
);
229
dnaex->
SetEmModel
(
new
G4DummyModel
(),1);
230
dnaex->
SetEmModel
(
new
G4DummyModel
(),2);
231
pmanager->
AddDiscreteProcess
(dnaex);
232
233
// DNA ionisation is not active in the world
234
G4DNAIonisation
* dnaioni =
new
G4DNAIonisation
(
"proton_G4DNAIonisation"
);
235
dnaioni->
SetEmModel
(
new
G4DummyModel
(),1);
236
dnaioni->
SetEmModel
(
new
G4DummyModel
(),2);
237
pmanager->
AddDiscreteProcess
(dnaioni);
238
239
// DNA charge decrease is ACTIVE in the world since no corresponding STANDARD process exist
240
pmanager->
AddDiscreteProcess
(
new
G4DNAChargeDecrease
(
"proton_G4DNAChargeDecrease"
));
241
242
}
else
if
( particleName ==
"hydrogen"
) {
243
244
// DNA processes are ACTIVE in the world since no corresponding STANDARD processes exist
245
pmanager->
AddDiscreteProcess
(
new
G4DNAIonisation
(
"hydrogen_G4DNAIonisation"
));
246
pmanager->
AddDiscreteProcess
(
new
G4DNAExcitation
(
"hydrogen_G4DNAExcitation"
));
247
pmanager->
AddDiscreteProcess
(
new
G4DNAChargeIncrease
(
"hydrogen_G4DNAChargeIncrease"
));
248
249
}
else
if
(particleName ==
"GenericIon"
) {
// THIS IS NEEDED FOR STANDARD ALPHA G4ionIonisation PROCESS
250
251
// STANDARD msc is active in the world
252
pmanager->
AddProcess
(
new
G4hMultipleScattering
, -1, 1, 1);
253
254
// STANDARD ionisation is active in the world
255
G4ionIonisation
* hion =
new
G4ionIonisation
();
256
hion->
SetEmModel
(
new
G4BraggIonGasModel
(),1);
257
hion->
SetEmModel
(
new
G4BetheBlochIonGasModel
(), 2);
258
pmanager->
AddProcess
(hion, -1, 2, 2);
259
260
}
else
if
( particleName ==
"alpha"
) {
261
262
// STANDARD msc is active in the world
263
G4hMultipleScattering
* msc =
new
G4hMultipleScattering
();
264
pmanager->
AddProcess
(msc, -1, 1, 1);
265
266
// STANDARD ionisation is active in the world
267
G4ionIonisation
* hion =
new
G4ionIonisation
();
268
hion->
SetEmModel
(
new
G4BraggIonGasModel
(),1);
269
hion->
SetEmModel
(
new
G4BetheBlochIonGasModel
(), 2);
270
pmanager->
AddProcess
(hion, -1, 2, 2);
271
272
// DNA excitation is not active in the world
273
G4DNAExcitation
* dnaex =
new
G4DNAExcitation
(
"alpha_G4DNAExcitation"
);
274
dnaex->
SetEmModel
(
new
G4DummyModel
(),1);
275
pmanager->
AddDiscreteProcess
(dnaex);
276
277
// DNA ionisation is not active in the world
278
G4DNAIonisation
* dnaioni =
new
G4DNAIonisation
(
"alpha_G4DNAIonisation"
);
279
dnaioni->
SetEmModel
(
new
G4DummyModel
(),1);
280
pmanager->
AddDiscreteProcess
(dnaioni);
281
282
// DNA charge decrease is ACTIVE in the world since no corresponding STANDARD process exist
283
pmanager->
AddDiscreteProcess
(
new
G4DNAChargeDecrease
(
"alpha_G4DNAChargeDecrease"
));
284
285
}
else
if
( particleName ==
"alpha+"
) {
286
287
// DNA processes are ACTIVE in the world since no corresponding STANDARD processes exist
288
pmanager->
AddDiscreteProcess
(
new
G4DNAExcitation
(
"alpha+_G4DNAExcitation"
));
289
pmanager->
AddDiscreteProcess
(
new
G4DNAIonisation
(
"alpha+_G4DNAIonisation"
));
290
pmanager->
AddDiscreteProcess
(
new
G4DNAChargeDecrease
(
"alpha+_G4DNAChargeDecrease"
));
291
pmanager->
AddDiscreteProcess
(
new
G4DNAChargeIncrease
(
"alpha+_G4DNAChargeIncrease"
));
292
293
}
else
if
( particleName ==
"helium"
) {
294
295
// DNA processes are ACTIVE in the world since no corresponding STANDARD processes exist
296
pmanager->
AddDiscreteProcess
(
new
G4DNAExcitation
(
"helium_G4DNAExcitation"
));
297
pmanager->
AddDiscreteProcess
(
new
G4DNAIonisation
(
"helium_G4DNAIonisation"
));
298
pmanager->
AddDiscreteProcess
(
new
G4DNAChargeIncrease
(
"helium_G4DNAChargeIncrease"
));
299
300
}
301
}
302
303
// **************************************
304
// 2) Define processes for Target region
305
// **************************************
306
307
// STANDARD EM processes should be inactivated when corresponding DNA processes are used
308
// - STANDARD EM e- processes are inactivated below 1 MeV
309
// - STANDARD EM proton & alpha processes are inactivated below standEnergyLimit
310
G4double
standEnergyLimit = 9.9*
MeV
;
311
//
312
313
G4double
massFactor = 1.0079/4.0026;
314
G4EmConfigurator
* em_config =
G4LossTableManager::Instance
()->
EmConfigurator
();
315
316
G4VEmModel
* mod;
317
318
// *** e-
319
320
// ---> STANDARD EM processes are inactivated below 1 MeV
321
322
mod =
new
G4UrbanMscModel93
();
323
mod->
SetActivationLowEnergyLimit
(1*
MeV
);
324
em_config->
SetExtraEmModel
(
"e-"
,
"msc"
,mod,
"Target"
);
325
326
mod =
new
G4MollerBhabhaModel
();
327
mod->
SetActivationLowEnergyLimit
(0.99*
MeV
);
328
em_config->
SetExtraEmModel
(
"e-"
,
"eIoni"
,mod,
"Target"
,0.0,100*
TeV
,
new
G4UniversalFluctuation
());
329
330
// ---> DNA processes activated
331
332
mod =
new
G4DNAChampionElasticModel
();
333
em_config->
SetExtraEmModel
(
"e-"
,
"e-_G4DNAElastic"
,mod,
"Target"
,0.0,1*
MeV
);
334
335
mod =
new
G4DNABornIonisationModel
();
336
em_config->
SetExtraEmModel
(
"e-"
,
"e-_G4DNAIonisation"
,mod,
"Target"
,11*
eV
,1*
MeV
);
337
338
mod =
new
G4DNABornExcitationModel
();
339
em_config->
SetExtraEmModel
(
"e-"
,
"e-_G4DNAExcitation"
,mod,
"Target"
,9*
eV
,1*
MeV
);
340
341
mod =
new
G4DNAMeltonAttachmentModel
();
342
em_config->
SetExtraEmModel
(
"e-"
,
"e-_G4DNAAttachment"
,mod,
"Target"
,4*
eV
,13*
eV
);
343
344
mod =
new
G4DNASancheExcitationModel
();
345
em_config->
SetExtraEmModel
(
"e-"
,
"e-_G4DNAVibExcitation"
,mod,
"Target"
,2*
eV
,100*
eV
);
346
347
// *** proton
348
349
// ---> STANDARD EM processes inactivated below standEnergyLimit
350
351
// STANDARD msc is still active
352
// Inactivate following STANDARD processes
353
354
mod =
new
G4BraggIonGasModel
();
355
mod->
SetActivationLowEnergyLimit
(standEnergyLimit);
356
em_config->
SetExtraEmModel
(
"proton"
,
"hIoni"
,mod,
"Target"
,0.0,2*
MeV
,
new
G4IonFluctuations
());
357
358
mod =
new
G4BetheBlochIonGasModel
();
359
mod->
SetActivationLowEnergyLimit
(standEnergyLimit);
360
em_config->
SetExtraEmModel
(
"proton"
,
"hIoni"
,mod,
"Target"
,2*
MeV
,100*
TeV
,
new
G4UniversalFluctuation
());
361
362
// ---> DNA processes activated
363
364
mod =
new
G4DNARuddIonisationModel
();
365
em_config->
SetExtraEmModel
(
"proton"
,
"proton_G4DNAIonisation"
,mod,
"Target"
,0.0,0.5*
MeV
);
366
367
mod =
new
G4DNABornIonisationModel
();
368
em_config->
SetExtraEmModel
(
"proton"
,
"proton_G4DNAIonisation"
,mod,
"Target"
,0.5*
MeV
,10*
MeV
);
369
370
mod =
new
G4DNAMillerGreenExcitationModel
();
371
em_config->
SetExtraEmModel
(
"proton"
,
"proton_G4DNAExcitation"
,mod,
"Target"
,10*
eV
,0.5*
MeV
);
372
373
mod =
new
G4DNABornExcitationModel
();
374
em_config->
SetExtraEmModel
(
"proton"
,
"proton_G4DNAExcitation"
,mod,
"Target"
,0.5*
MeV
,10*
MeV
);
375
376
// *** alpha
377
378
// ---> STANDARD EM processes inactivated below standEnergyLimit
379
380
// STANDARD msc is still active
381
// Inactivate following STANDARD processes
382
383
mod =
new
G4BraggIonGasModel
();
384
mod->
SetActivationLowEnergyLimit
(standEnergyLimit);
385
em_config->
SetExtraEmModel
(
"alpha"
,
"ionIoni"
,mod,
"Target"
,0.0,2*
MeV
/massFactor,
new
G4IonFluctuations
());
386
387
mod =
new
G4BetheBlochIonGasModel
();
388
mod->
SetActivationLowEnergyLimit
(standEnergyLimit);
389
em_config->
SetExtraEmModel
(
"alpha"
,
"ionIoni"
,mod,
"Target"
,2*
MeV
/massFactor,100*
TeV
,
new
G4UniversalFluctuation
());
390
391
// ---> DNA processes activated
392
393
mod =
new
G4DNARuddIonisationModel
();
394
em_config->
SetExtraEmModel
(
"alpha"
,
"alpha_G4DNAIonisation"
,mod,
"Target"
,0.0,10*
MeV
);
395
396
mod =
new
G4DNAMillerGreenExcitationModel
();
397
em_config->
SetExtraEmModel
(
"alpha"
,
"alpha_G4DNAExcitation"
,mod,
"Target"
,1*
keV
,10*
MeV
);
398
399
}
400
401
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
402
403
void
PhysicsList::ConstructGeneral
()
404
{ }
405
406
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
407
408
void
PhysicsList::SetCuts
()
409
{
410
if
(
verboseLevel
>0)
411
{
412
G4cout
<<
"PhysicsList::SetCuts:"
;
413
G4cout
<<
"CutLength : "
<<
G4BestUnit
(
defaultCutValue
,
"Length"
) <<
G4endl
;
414
}
415
416
// set cut values for gamma at first and for e- second and next for e+,
417
// because some processes for e+/e- need cut values for gamma
418
SetCutValue
(cutForGamma,
"gamma"
);
419
SetCutValue
(cutForElectron,
"e-"
);
420
SetCutValue
(cutForPositron,
"e+"
);
421
SetCutValue
(cutForProton,
"proton"
);
422
423
if
(
verboseLevel
>0) {
DumpCutValuesTable
(); }
424
}
Generated on Sat May 25 2013 14:31:57 for Geant4 by
1.8.4