Geant4_10
RE06PhysicsList.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 //
28 //
29 // $Id: RE06PhysicsList.cc 66526 2012-12-19 13:41:33Z ihrivnac $
30 //
31 
32 #include "RE06PhysicsList.hh"
33 
34 #include "G4ParticleDefinition.hh"
35 #include "G4ProcessManager.hh"
36 #include "G4ProcessVector.hh"
37 #include "G4ParticleTypes.hh"
38 #include "G4ParticleTable.hh"
39 #include "G4ios.hh"
40 #include "G4SystemOfUnits.hh"
41 
42 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
43 
45 {
46  defaultCutValue = 1.0*mm;
47  SetVerboseLevel(1);
48 }
49 
50 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
51 
53 {}
54 
55 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
56 
58 {
59  // In this method, static member functions should be called
60  // for all particles which you want to use.
61  // This ensures that objects of these particle types will be
62  // created in the program.
63 
68 }
69 
70 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
71 
73 {
74  // pseudo-particles
77 
78  // gamma
80 }
81 
82 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
83 
85 {
86  // leptons
91 
96 }
97 
98 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
99 
101 {
102  // mesons
114 }
115 
116 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
117 
119 {
120 // barions
125 }
126 
127 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
128 
130 {
133  ConstructEM();
134 }
135 
136 #include "G4ComptonScattering.hh"
137 #include "G4GammaConversion.hh"
138 #include "G4PhotoElectricEffect.hh"
139 
140 #include "G4eMultipleScattering.hh"
141 #include "G4MuMultipleScattering.hh"
142 #include "G4hMultipleScattering.hh"
143 
144 #include "G4eIonisation.hh"
145 #include "G4eBremsstrahlung.hh"
146 #include "G4eplusAnnihilation.hh"
147 
148 #include "G4MuIonisation.hh"
149 #include "G4MuBremsstrahlung.hh"
150 #include "G4MuPairProduction.hh"
151 
152 #include "G4hIonisation.hh"
153 
154 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
155 
157 {
158  G4bool displacementFlg = true; // for multiple scattering
159  theParticleIterator->reset();
160  while( (*theParticleIterator)() ){
161  G4ParticleDefinition* particle = theParticleIterator->value();
162  G4ProcessManager* pmanager = particle->GetProcessManager();
163  G4String particleName = particle->GetParticleName();
164 
165  if (particleName == "gamma") {
166  // gamma
167  pmanager->AddDiscreteProcess(new G4GammaConversion());
168  pmanager->AddDiscreteProcess(new G4ComptonScattering());
169  pmanager->AddDiscreteProcess(new G4PhotoElectricEffect());
170 
171  } else if (particleName == "e-") {
172  //electron
173  G4eMultipleScattering* theeminusMultipleScattering
174  = new G4eMultipleScattering();
175  theeminusMultipleScattering->SetLateralDisplasmentFlag(displacementFlg);
176  G4VProcess* theeminusIonisation = new G4eIonisation();
177  G4VProcess* theeminusBremsstrahlung = new G4eBremsstrahlung();
178  //
179  // add processes
180  pmanager->AddProcess(theeminusMultipleScattering);
181  pmanager->AddProcess(theeminusIonisation);
182  pmanager->AddProcess(theeminusBremsstrahlung);
183  //
184  // set ordering for AlongStepDoIt
185  pmanager->SetProcessOrdering(theeminusMultipleScattering, idxAlongStep,1);
186  pmanager->SetProcessOrdering(theeminusIonisation, idxAlongStep,2);
187  pmanager->SetProcessOrdering(theeminusBremsstrahlung, idxAlongStep,3);
188  //
189  // set ordering for PostStepDoIt
190  pmanager->SetProcessOrdering(theeminusMultipleScattering, idxPostStep,1);
191  pmanager->SetProcessOrdering(theeminusIonisation, idxPostStep,2);
192  pmanager->SetProcessOrdering(theeminusBremsstrahlung, idxPostStep,3);
193 
194  } else if (particleName == "e+") {
195  //positron
196  G4eMultipleScattering* theeplusMultipleScattering
197  = new G4eMultipleScattering();
198  theeplusMultipleScattering->SetLateralDisplasmentFlag(displacementFlg);
199  G4VProcess* theeplusIonisation = new G4eIonisation();
200  G4VProcess* theeplusBremsstrahlung = new G4eBremsstrahlung();
201  G4VProcess* theeplusAnnihilation = new G4eplusAnnihilation();
202  //
203  // add processes
204  pmanager->AddProcess(theeplusMultipleScattering);
205  pmanager->AddProcess(theeplusIonisation);
206  pmanager->AddProcess(theeplusBremsstrahlung);
207  pmanager->AddProcess(theeplusAnnihilation);
208  //
209  // set ordering for AtRestDoIt
210  pmanager->SetProcessOrderingToFirst(theeplusAnnihilation, idxAtRest);
211  //
212  // set ordering for AlongStepDoIt
213  pmanager->SetProcessOrdering(theeplusMultipleScattering, idxAlongStep,1);
214  pmanager->SetProcessOrdering(theeplusIonisation, idxAlongStep,2);
215  pmanager->SetProcessOrdering(theeplusBremsstrahlung, idxAlongStep,3);
216  //
217  // set ordering for PostStepDoIt
218  pmanager->SetProcessOrdering(theeplusMultipleScattering, idxPostStep,1);
219  pmanager->SetProcessOrdering(theeplusIonisation, idxPostStep,2);
220  pmanager->SetProcessOrdering(theeplusBremsstrahlung, idxPostStep,3);
221  pmanager->SetProcessOrdering(theeplusAnnihilation, idxPostStep,4);
222 
223  } else if( particleName == "mu+" ||
224  particleName == "mu-" ) {
225  //muon
226  G4MuMultipleScattering* aMultipleScattering
227  = new G4MuMultipleScattering();
228  aMultipleScattering->SetLateralDisplasmentFlag(displacementFlg);
229  G4VProcess* aBremsstrahlung = new G4MuBremsstrahlung();
230  G4VProcess* aPairProduction = new G4MuPairProduction();
231  G4VProcess* anIonisation = new G4MuIonisation();
232  //
233  // add processes
234  pmanager->AddProcess(anIonisation);
235  pmanager->AddProcess(aMultipleScattering);
236  pmanager->AddProcess(aBremsstrahlung);
237  pmanager->AddProcess(aPairProduction);
238  //
239  // set ordering for AlongStepDoIt
240  pmanager->SetProcessOrdering(aMultipleScattering, idxAlongStep,1);
241  pmanager->SetProcessOrdering(anIonisation, idxAlongStep,2);
242  pmanager->SetProcessOrdering(aBremsstrahlung, idxAlongStep,3);
243  pmanager->SetProcessOrdering(aPairProduction, idxAlongStep,4);
244 
245  //
246  // set ordering for PostStepDoIt
247  pmanager->SetProcessOrdering(aMultipleScattering, idxPostStep,1);
248  pmanager->SetProcessOrdering(anIonisation, idxPostStep,2);
249  pmanager->SetProcessOrdering(aBremsstrahlung, idxPostStep,3);
250  pmanager->SetProcessOrdering(aPairProduction, idxPostStep,4);
251 
252  } else if ((!particle->IsShortLived()) &&
253  (particle->GetPDGCharge() != 0.0) &&
254  (particle->GetParticleName() != "chargedgeantino")) {
255  // all others charged particles except geantino
256  G4hMultipleScattering* aMultipleScattering = new G4hMultipleScattering();
257  aMultipleScattering->SetLateralDisplasmentFlag(displacementFlg);
258  G4VProcess* anIonisation = new G4hIonisation();
259  //
260  // add processes
261  pmanager->AddProcess(anIonisation);
262  pmanager->AddProcess(aMultipleScattering);
263  //
264  // set ordering for AlongStepDoIt
265  pmanager->SetProcessOrdering(aMultipleScattering, idxAlongStep,1);
266  pmanager->SetProcessOrdering(anIonisation, idxAlongStep,2);
267  //
268  // set ordering for PostStepDoIt
269  pmanager->SetProcessOrdering(aMultipleScattering, idxPostStep,1);
270  pmanager->SetProcessOrdering(anIonisation, idxPostStep,2);
271  }
272  }
273 }
274 
275 #include "G4Decay.hh"
277 
278 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
279 
281 {
282  // Add Decay Process and Parallel world Scoring Process
283  G4Decay* theDecayProcess = new G4Decay();
284  G4ParallelWorldScoringProcess* theParallelWorldScoringProcess
285  = new G4ParallelWorldScoringProcess("ParaWorldScoringProc");
286  theParallelWorldScoringProcess->SetParallelWorld("ParallelScoringWorld");
287 
288  theParticleIterator->reset();
289  while( (*theParticleIterator)() ){
290  G4ParticleDefinition* particle = theParticleIterator->value();
291  G4ProcessManager* pmanager = particle->GetProcessManager();
292  if (theDecayProcess->IsApplicable(*particle)) {
293  pmanager ->AddProcess(theDecayProcess);
294  pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep);
295  pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest);
296  }
297  pmanager->AddProcess(theParallelWorldScoringProcess);
298  pmanager
299  ->SetProcessOrderingToLast(theParallelWorldScoringProcess, idxAtRest);
300  pmanager
301  ->SetProcessOrdering(theParallelWorldScoringProcess, idxAlongStep, 1);
302  pmanager
303  ->SetProcessOrderingToLast(theParallelWorldScoringProcess, idxPostStep);
304  }
305 }
306 
307 #include "G4Region.hh"
308 #include "G4RegionStore.hh"
309 #include "G4ProductionCuts.hh"
310 
311 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
312 
314 {
315  if (verboseLevel >0){
316  G4cout << "RE06PhysicsList::SetCuts: default cut length : "
317  << G4BestUnit(defaultCutValue,"Length") << G4endl;
318  }
319 
320  // These values are used as the default production thresholds
321  // for the world volume.
323 
324 
325  // Production thresholds for detector regions
326  G4String regName[] = {"Calor-A","Calor-B","Calor-C"};
327  G4double fuc = 1.;
328  for(G4int i=0;i<3;i++)
329  {
330  G4Region* reg = G4RegionStore::GetInstance()->GetRegion(regName[i]);
333  reg->SetProductionCuts(cuts);
334  fuc *= 10.;
335  }
336 }
337 
338 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
static G4Geantino * GeantinoDefinition()
Definition: G4Geantino.cc:82
static G4Electron * ElectronDefinition()
Definition: G4Electron.cc:89
static G4KaonPlus * KaonPlusDefinition()
Definition: G4KaonPlus.cc:108
G4Region * GetRegion(const G4String &name, G4bool verbose=true) const
static G4MuonPlus * MuonPlusDefinition()
Definition: G4MuonPlus.cc:94
void SetProcessOrderingToFirst(G4VProcess *aProcess, G4ProcessVectorDoItIndex idDoIt)
void SetLateralDisplasmentFlag(G4bool val)
virtual G4bool IsApplicable(const G4ParticleDefinition &)
Definition: G4Decay.cc:89
static G4KaonMinus * KaonMinusDefinition()
Definition: G4KaonMinus.cc:108
static G4Proton * ProtonDefinition()
Definition: G4Proton.cc:88
G4int AddDiscreteProcess(G4VProcess *aProcess, G4int ord=ordDefault)
void SetProductionCut(G4double cut, G4int index=-1)
static G4KaonZero * KaonZeroDefinition()
Definition: G4KaonZero.cc:99
static G4AntiKaonZero * AntiKaonZeroDefinition()
static G4KaonZeroShort * KaonZeroShortDefinition()
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
static G4AntiProton * AntiProtonDefinition()
Definition: G4AntiProton.cc:88
G4ProcessManager * GetProcessManager() const
int G4int
Definition: G4Types.hh:78
static G4AntiNeutron * AntiNeutronDefinition()
static G4PionZero * PionZeroDefinition()
Definition: G4PionZero.cc:99
const G4String & GetParticleName() const
static G4RegionStore * GetInstance()
G4GLOB_DLL std::ostream G4cout
static G4PionMinus * PionMinusDefinition()
Definition: G4PionMinus.cc:93
bool G4bool
Definition: G4Types.hh:79
static G4KaonZeroLong * KaonZeroLongDefinition()
void SetVerboseLevel(G4int value)
G4int AddProcess(G4VProcess *aProcess, G4int ordAtRestDoIt=ordInActive, G4int ordAlongSteptDoIt=ordInActive, G4int ordPostStepDoIt=ordInActive)
static G4NeutrinoE * NeutrinoEDefinition()
Definition: G4NeutrinoE.cc:80
void SetParallelWorld(G4String parallelWorldName)
static G4AntiNeutrinoMu * AntiNeutrinoMuDefinition()
static G4PionPlus * PionPlusDefinition()
Definition: G4PionPlus.cc:93
void SetProcessOrdering(G4VProcess *aProcess, G4ProcessVectorDoItIndex idDoIt, G4int ordDoIt=ordDefault)
virtual ~RE06PhysicsList()
Definition of the RE06PhysicsList class.
virtual void ConstructProcess()
static G4Positron * PositronDefinition()
Definition: G4Positron.cc:89
static G4EtaPrime * EtaPrimeDefinition()
Definition: G4EtaPrime.cc:100
void SetProcessOrderingToLast(G4VProcess *aProcess, G4ProcessVectorDoItIndex idDoIt)
static G4MuonMinus * MuonMinusDefinition()
Definition: G4MuonMinus.cc:95
static G4ChargedGeantino * ChargedGeantinoDefinition()
#define G4endl
Definition: G4ios.hh:61
virtual void SetCuts()
void SetProductionCuts(G4ProductionCuts *cut)
static G4AntiNeutrinoE * AntiNeutrinoEDefinition()
double G4double
Definition: G4Types.hh:76
G4double GetPDGCharge() const
static G4NeutrinoMu * NeutrinoMuDefinition()
Definition: G4NeutrinoMu.cc:80
static G4Neutron * NeutronDefinition()
Definition: G4Neutron.cc:99
#define theParticleIterator
virtual void ConstructParticle()
static G4Eta * EtaDefinition()
Definition: G4Eta.cc:104
static G4Gamma * GammaDefinition()
Definition: G4Gamma.cc:81