Geant4 10.7.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4EmPenelopePhysics.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
29#include "G4SystemOfUnits.hh"
30
31// *** Processes and models
32
33// gamma
34
37
40
41#include "G4GammaConversion.hh"
43
46
49
50// e- and e+
51
54
55#include "G4eIonisation.hh"
57
58#include "G4eBremsstrahlung.hh"
60
61// e+ only
62
65
66// hadrons
67
69#include "G4MscStepLimitType.hh"
70
71#include "G4ePairProduction.hh"
72
73#include "G4hIonisation.hh"
74#include "G4ionIonisation.hh"
76#include "G4NuclearStopping.hh"
77
78// msc models
79#include "G4UrbanMscModel.hh"
81#include "G4WentzelVIModel.hh"
84//
85
86#include "G4EmBuilder.hh"
87
88// particles
89
90#include "G4Gamma.hh"
91#include "G4Electron.hh"
92#include "G4Positron.hh"
93#include "G4GenericIon.hh"
94
95//
97#include "G4BuilderType.hh"
98#include "G4EmModelActivator.hh"
99
100// factory
102//
104
105//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
106
108 : G4VPhysicsConstructor("G4EmPenelope"), verbose(ver)
109{
111 param->SetDefaults();
112 param->SetVerbose(verbose);
113 param->SetMinEnergy(100*CLHEP::eV);
114 param->SetLowestElectronEnergy(100*CLHEP::eV);
115 param->SetNumberOfBinsPerDecade(20);
117 param->SetStepFunction(0.2, 10*CLHEP::um);
118 param->SetStepFunctionMuHad(0.1, 50*CLHEP::um);
119 param->SetStepFunctionLightIons(0.1, 20*CLHEP::um);
120 param->SetStepFunctionIons(0.1, 1*CLHEP::um);
121 param->SetUseMottCorrection(true);
123 param->SetMscSkin(3);
124 param->SetMscRangeFactor(0.08);
125 param->SetMuHadLateralDisplacement(true);
126 param->SetFluo(true);
127 param->SetMaxNIELEnergy(1*CLHEP::MeV);
128 param->SetPIXEElectronCrossSectionModel("Penelope");
130}
131
132//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
133
135{}
136
137//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
138
140{
141 // minimal set of particles for EM physics
143}
144
145//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
146
148{
149 if(verbose > 1) {
150 G4cout << "### " << GetPhysicsName() << " Construct Processes " << G4endl;
151 }
154
155 // processes used by several particles
157 G4hMultipleScattering* hmsc = new G4hMultipleScattering("ionmsc");
158
159 // high energy limit for e+- scattering models
160 G4double highEnergyLimit = G4EmParameters::Instance()->MscEnergyLimit();
161
162 // nuclear stopping
163 G4double nielEnergyLimit = G4EmParameters::Instance()->MaxNIELEnergy();
165 pnuc->SetMaxKinEnergy(nielEnergyLimit);
166
167 //Applicability range for Penelope models
168 //for higher energies, the Standard models are used
169 G4double PenelopeHighEnergyLimit = 1.0*CLHEP::GeV;
170
171 // Add gamma EM Processes
173
174 //Photo-electric effect
175 G4PhotoElectricEffect* thePhotoElectricEffect = new G4PhotoElectricEffect();
176 thePhotoElectricEffect->SetEmModel(new G4PEEffectFluoModel());
178 thePEPenelopeModel->SetHighEnergyLimit(PenelopeHighEnergyLimit);
179 thePhotoElectricEffect->AddEmModel(0, thePEPenelopeModel);
180
181 //Compton scattering
182 G4ComptonScattering* theComptonScattering = new G4ComptonScattering();
183 G4PenelopeComptonModel* theComptonPenelopeModel = new G4PenelopeComptonModel();
184 theComptonScattering->SetEmModel(new G4KleinNishinaModel());
185 theComptonPenelopeModel->SetHighEnergyLimit(PenelopeHighEnergyLimit);
186 theComptonScattering->AddEmModel(0, theComptonPenelopeModel);
187
188 //Gamma conversion
189 G4GammaConversion* theGammaConversion = new G4GammaConversion();
191 theGCPenelopeModel->SetHighEnergyLimit(PenelopeHighEnergyLimit);
192 theGammaConversion->AddEmModel(0, theGCPenelopeModel);
193
194 //Rayleigh scattering
195 G4RayleighScattering* theRayleigh = new G4RayleighScattering();
196 G4PenelopeRayleighModel* theRayleighPenelopeModel = new G4PenelopeRayleighModel();
197 theRayleigh->SetEmModel(theRayleighPenelopeModel);
198
199 ph->RegisterProcess(thePhotoElectricEffect, particle);
200 ph->RegisterProcess(theComptonScattering, particle);
201 ph->RegisterProcess(theGammaConversion, particle);
202 ph->RegisterProcess(theRayleigh, particle);
203
204 // e-
205 particle = G4Electron::Electron();
206
207 // multiple scattering
211 msc1->SetHighEnergyLimit(highEnergyLimit);
212 msc2->SetLowEnergyLimit(highEnergyLimit);
213 msc->SetEmModel(msc1);
214 msc->SetEmModel(msc2);
215
218 ss->SetEmModel(ssm);
219 ss->SetMinKinEnergy(highEnergyLimit);
220 ssm->SetLowEnergyLimit(highEnergyLimit);
221 ssm->SetActivationLowEnergyLimit(highEnergyLimit);
222
223 //Ionisation
224 G4eIonisation* eioni = new G4eIonisation();
226 theIoniPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
227 eioni->AddEmModel(0, theIoniPenelope, new G4UniversalFluctuation());
228
229 //Bremsstrahlung
232 theBremPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
233 brem->SetEmModel(theBremPenelope);
234
235 // register processes
236 ph->RegisterProcess(msc, particle);
237 ph->RegisterProcess(eioni, particle);
238 ph->RegisterProcess(brem, particle);
239 ph->RegisterProcess(ee, particle);
240 ph->RegisterProcess(ss, particle);
241
242 // e+
243 particle = G4Positron::Positron();
244
245 // multiple scattering
246 msc = new G4eMultipleScattering;
247 msc1 = new G4GoudsmitSaundersonMscModel();
248 msc2 = new G4WentzelVIModel();
249 msc1->SetHighEnergyLimit(highEnergyLimit);
250 msc2->SetLowEnergyLimit(highEnergyLimit);
251 msc->SetEmModel(msc1);
252 msc->SetEmModel(msc2);
253
254 ssm = new G4eCoulombScatteringModel();
255 ss = new G4CoulombScattering();
256 ss->SetEmModel(ssm);
257 ss->SetMinKinEnergy(highEnergyLimit);
258 ssm->SetLowEnergyLimit(highEnergyLimit);
259 ssm->SetActivationLowEnergyLimit(highEnergyLimit);
260
261 //Ionisation
262 eioni = new G4eIonisation();
263 theIoniPenelope = new G4PenelopeIonisationModel();
264 theIoniPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
265 eioni->AddEmModel(0,theIoniPenelope, new G4UniversalFluctuation());
266
267 //Bremsstrahlung
268 brem = new G4eBremsstrahlung();
269 theBremPenelope = new G4PenelopeBremsstrahlungModel();
270 theBremPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
271 brem->SetEmModel(theBremPenelope);
272
273 //Annihilation
276 theAnnPenelope->SetHighEnergyLimit(PenelopeHighEnergyLimit);
277 anni->AddEmModel(0, theAnnPenelope);
278
279 // register processes
280 ph->RegisterProcess(msc, particle);
281 ph->RegisterProcess(eioni, particle);
282 ph->RegisterProcess(brem, particle);
283 ph->RegisterProcess(ee, particle);
284 ph->RegisterProcess(anni, particle);
285 ph->RegisterProcess(ss, particle);
286
287 // generic ion
288 particle = G4GenericIon::GenericIon();
289 G4ionIonisation* ionIoni = new G4ionIonisation();
291 ph->RegisterProcess(hmsc, particle);
292 ph->RegisterProcess(ionIoni, particle);
293 ph->RegisterProcess(pnuc, particle);
294
295 // muons, hadrons, ions
297
298 // extra configuration
300}
301
302//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
@ bElectromagnetic
@ fUseSafetyPlus
#define G4_DECLARE_PHYSCONSTR_FACTORY(physics_constructor)
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
static G4Electron * Electron()
Definition: G4Electron.cc:93
static void ConstructCharged(G4hMultipleScattering *hmsc, G4NuclearStopping *nucStopping, G4bool isWVI=true)
Definition: G4EmBuilder.cc:170
static void ConstructMinimalEmSet()
Definition: G4EmBuilder.cc:234
static void PrepareEMPhysics()
Definition: G4EmBuilder.cc:259
void SetMinEnergy(G4double val)
void SetLowestElectronEnergy(G4double val)
void SetStepFunctionLightIons(G4double v1, G4double v2)
void SetNumberOfBinsPerDecade(G4int val)
static G4EmParameters * Instance()
G4double MaxNIELEnergy() const
G4double MscEnergyLimit() const
void SetPIXEElectronCrossSectionModel(const G4String &)
void SetMuHadLateralDisplacement(G4bool val)
void ActivateAngularGeneratorForIonisation(G4bool val)
void SetStepFunction(G4double v1, G4double v2)
void SetFluo(G4bool val)
void SetStepFunctionMuHad(G4double v1, G4double v2)
void SetVerbose(G4int val)
void SetMscSkin(G4double val)
void SetMaxNIELEnergy(G4double val)
void SetStepFunctionIons(G4double v1, G4double v2)
void SetMscStepLimitType(G4MscStepLimitType val)
void SetUseMottCorrection(G4bool val)
void SetMscRangeFactor(G4double val)
void ConstructParticle() override
G4EmPenelopePhysics(G4int ver=1, const G4String &name="")
void ConstructProcess() override
static G4Gamma * Gamma()
Definition: G4Gamma.cc:85
static G4GenericIon * GenericIon()
Definition: G4GenericIon.cc:92
G4bool RegisterProcess(G4VProcess *process, G4ParticleDefinition *particle)
static G4PhysicsListHelper * GetPhysicsListHelper()
static G4Positron * Positron()
Definition: G4Positron.cc:93
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:757
void SetActivationLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:778
void SetLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:764
void SetMinKinEnergy(G4double e)
void AddEmModel(G4int, G4VEmModel *, const G4Region *region=nullptr)
void SetEmModel(G4VEmModel *, G4int index=0)
void SetMaxKinEnergy(G4double e)
void SetEmModel(G4VEmModel *, G4int index=0)
void AddEmModel(G4int, G4VEmModel *, G4VEmFluctuationModel *fluc=0, const G4Region *region=nullptr)
void SetEmModel(G4VMscModel *, size_t index=0)
const G4String & GetPhysicsName() const