| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // 26 //
27 // ------------------------------------------- 27 // --------------------------------------------------------------
28 // GEANT 4 - Underground Dark Matter Detecto 28 // GEANT 4 - Underground Dark Matter Detector Advanced Example
29 // 29 //
30 // For information related to this code c 30 // For information related to this code contact: Alex Howard
31 // e-mail: alexander.howard@cern.ch 31 // e-mail: alexander.howard@cern.ch
32 // ------------------------------------------- 32 // --------------------------------------------------------------
33 // Comments 33 // Comments
34 // 34 //
35 // Underground Advanced 35 // Underground Advanced
36 // by A. Howard and H. Araujo 36 // by A. Howard and H. Araujo
37 // (27th November 2001) 37 // (27th November 2001)
38 // 38 //
39 // PhysicsList program 39 // PhysicsList program
40 // 40 //
41 // Modified: 41 // Modified:
42 // 42 //
43 // 14-02-03 Fix bugs in msc and hIon instancia 43 // 14-02-03 Fix bugs in msc and hIon instanciation + cut per region
44 // 44 //
45 // 05-02-05 AH - changes to G4Decay - added is 45 // 05-02-05 AH - changes to G4Decay - added is not short lived protection
46 // and redefined particles to allow n 46 // and redefined particles to allow non-static creation
47 // i.e. changed construction to G4Mes 47 // i.e. changed construction to G4MesonConstructor, G4BaryonConstructor
48 // <<
49 // 23-10-09 LP - migrated EM physics from the <<
50 // the new G4Livermore model implemen <<
51 // 48 //
52 // ------------------------------------------- 49 // --------------------------------------------------------------
53 50
54 #include <iomanip> <<
55 <<
56 #include "DMXPhysicsList.hh" 51 #include "DMXPhysicsList.hh"
57 52
58 #include "globals.hh" 53 #include "globals.hh"
59 #include "G4SystemOfUnits.hh" <<
60 #include "G4ProcessManager.hh" 54 #include "G4ProcessManager.hh"
61 #include "G4ProcessVector.hh" 55 #include "G4ProcessVector.hh"
62 56
63 #include "G4ParticleDefinition.hh" 57 #include "G4ParticleDefinition.hh"
64 #include "G4ParticleWithCuts.hh" 58 #include "G4ParticleWithCuts.hh"
65 #include "G4ParticleTypes.hh" 59 #include "G4ParticleTypes.hh"
66 #include "G4ParticleTable.hh" 60 #include "G4ParticleTable.hh"
67 61
68 #include "G4ios.hh" 62 #include "G4ios.hh"
69 #include "G4UserLimits.hh" << 63 #include <iomanip>
70 <<
71 #include "G4MesonConstructor.hh" <<
72 #include "G4BaryonConstructor.hh" <<
73 #include "G4IonConstructor.hh" <<
74 #include "G4ShortLivedConstructor.hh" <<
75 <<
76 #include "DMXMaxTimeCuts.hh" <<
77 #include "DMXMinEkineCuts.hh" <<
78 #include "G4StepLimiter.hh" <<
79 <<
80 // gamma <<
81 #include "G4PhotoElectricEffect.hh" <<
82 #include "G4LivermorePhotoElectricModel.hh" <<
83 <<
84 #include "G4ComptonScattering.hh" <<
85 #include "G4LivermoreComptonModel.hh" <<
86 <<
87 #include "G4GammaConversion.hh" <<
88 #include "G4BetheHeitler5DModel.hh" <<
89 <<
90 #include "G4RayleighScattering.hh" <<
91 #include "G4LivermoreRayleighModel.hh" <<
92 <<
93 // e- <<
94 #include "G4eMultipleScattering.hh" <<
95 <<
96 #include "G4eIonisation.hh" <<
97 #include "G4LivermoreIonisationModel.hh" <<
98 <<
99 #include "G4eBremsstrahlung.hh" <<
100 #include "G4UniversalFluctuation.hh" <<
101 <<
102 // e+ <<
103 #include "G4eIonisation.hh" <<
104 #include "G4eBremsstrahlung.hh" <<
105 #include "G4eplusAnnihilation.hh" <<
106 <<
107 // alpha and GenericIon and deuterons, triton, <<
108 //muon: <<
109 #include "G4MuIonisation.hh" <<
110 #include "G4MuBremsstrahlung.hh" <<
111 #include "G4MuPairProduction.hh" <<
112 #include "G4MuMultipleScattering.hh" <<
113 #include "G4MuonMinusCapture.hh" <<
114 <<
115 //OTHERS: <<
116 #include "G4hIonisation.hh" <<
117 #include "G4hMultipleScattering.hh" <<
118 #include "G4hBremsstrahlung.hh" <<
119 #include "G4ionIonisation.hh" <<
120 #include "G4IonParametrisedLossModel.hh" <<
121 #include "G4NuclearStopping.hh" <<
122 <<
123 //em process options to allow msc step-limitat <<
124 #include "G4EmParameters.hh" <<
125 #include "G4VAtomDeexcitation.hh" <<
126 #include "G4UAtomicDeexcitation.hh" <<
127 #include "G4LossTableManager.hh" <<
128 <<
129 #include "G4Scintillation.hh" <<
130 #include "G4OpAbsorption.hh" <<
131 #include "G4OpBoundaryProcess.hh" <<
132 #include "G4OpticalParameters.hh" <<
133 <<
134 // Elastic processes: <<
135 #include "G4HadronElasticProcess.hh" <<
136 #include "G4ChipsElasticModel.hh" <<
137 #include "G4ElasticHadrNucleusHE.hh" <<
138 <<
139 // Inelastic processes: <<
140 #include "G4HadronInelasticProcess.hh" <<
141 <<
142 // High energy FTFP model and Bertini cascade <<
143 #include "G4FTFModel.hh" <<
144 #include "G4LundStringFragmentation.hh" <<
145 #include "G4ExcitedStringDecay.hh" <<
146 #include "G4PreCompoundModel.hh" <<
147 #include "G4GeneratorPrecompoundInterface.hh" <<
148 #include "G4TheoFSGenerator.hh" <<
149 #include "G4CascadeInterface.hh" <<
150 <<
151 // Cross sections <<
152 #include "G4VCrossSectionDataSet.hh" <<
153 #include "G4CrossSectionDataSetRegistry.hh" <<
154 <<
155 #include "G4CrossSectionElastic.hh" <<
156 #include "G4CrossSectionInelastic.hh" <<
157 #include "G4BGGPionElasticXS.hh" <<
158 #include "G4BGGPionInelasticXS.hh" <<
159 #include "G4AntiNuclElastic.hh" <<
160 <<
161 #include "G4CrossSectionInelastic.hh" <<
162 #include "G4BGGNucleonInelasticXS.hh" <<
163 #include "G4BGGNucleonElasticXS.hh" <<
164 #include "G4NeutronInelasticXS.hh" <<
165 #include "G4NeutronElasticXS.hh" <<
166 #include "G4ComponentAntiNuclNuclearXS.hh" <<
167 #include "G4ComponentGGNuclNuclXsc.hh" <<
168 #include "G4ComponentGGHadronNucleusXsc.hh" <<
169 <<
170 #include "G4HadronElastic.hh" <<
171 #include "G4NeutronCaptureProcess.hh" <<
172 <<
173 // Neutron high-precision models: <20 MeV <<
174 #include "G4ParticleHPElastic.hh" <<
175 #include "G4ParticleHPElasticData.hh" <<
176 #include "G4ParticleHPCapture.hh" <<
177 #include "G4ParticleHPCaptureData.hh" <<
178 #include "G4ParticleHPInelastic.hh" <<
179 #include "G4ParticleHPInelasticData.hh" <<
180 <<
181 // Stopping processes <<
182 #include "G4HadronStoppingProcess.hh" <<
183 #include "G4HadronicAbsorptionBertini.hh" <<
184 #include "G4HadronicAbsorptionFritiof.hh" <<
185 64
186 #include "G4HadronicParameters.hh" << 65 #include "G4UserLimits.hh"
187 66
188 #include "G4Decay.hh" <<
189 #include "G4RadioactiveDecay.hh" <<
190 #include "G4PhysicsListHelper.hh" <<
191 #include "G4NuclideTable.hh" <<
192 #include "G4NuclearLevelData.hh" <<
193 67
194 // Constructor /////////////////////////////// 68 // Constructor /////////////////////////////////////////////////////////////
195 DMXPhysicsList::DMXPhysicsList() : G4VUserPhys 69 DMXPhysicsList::DMXPhysicsList() : G4VUserPhysicsList()
196 { 70 {
197 71
198 defaultCutValue = 1.0*micrometer; // 72 defaultCutValue = 1.0*micrometer; //
199 cutForGamma = defaultCutValue; 73 cutForGamma = defaultCutValue;
200 cutForElectron = 1.0*nanometer; 74 cutForElectron = 1.0*nanometer;
201 cutForPositron = defaultCutValue; 75 cutForPositron = defaultCutValue;
202 76
203 VerboseLevel = 1; 77 VerboseLevel = 1;
204 OpVerbLevel = 0; 78 OpVerbLevel = 0;
205 79
206 //set a finer grid of the physic tables in o <<
207 //former LowEnergy models have 200 bins up t <<
208 G4EmParameters* param = G4EmParameters::Inst <<
209 param->SetMaxEnergy(100*GeV); <<
210 param->SetNumberOfBinsPerDecade(20); <<
211 param->SetMscStepLimitType(fMinimal); <<
212 param->SetFluo(true); <<
213 param->SetPixe(true); <<
214 param->SetAuger(true); <<
215 <<
216 G4EmParameters::Instance()->AddPhysics("Worl <<
217 G4DeexPrecoParameters* deex = G4NuclearLevel <<
218 deex->SetStoreICLevelData(true); <<
219 deex->SetMaxLifeTime(G4NuclideTable::GetInst <<
220 /std::log(2.)); <<
221 SetVerboseLevel(VerboseLevel); 80 SetVerboseLevel(VerboseLevel);
222 } 81 }
223 82
>> 83
224 // Destructor //////////////////////////////// 84 // Destructor //////////////////////////////////////////////////////////////
225 DMXPhysicsList::~DMXPhysicsList() 85 DMXPhysicsList::~DMXPhysicsList()
226 {;} 86 {;}
227 87
>> 88
228 // Construct Particles /////////////////////// 89 // Construct Particles /////////////////////////////////////////////////////
229 void DMXPhysicsList::ConstructParticle() 90 void DMXPhysicsList::ConstructParticle()
230 { 91 {
>> 92
231 // In this method, static member functions s 93 // In this method, static member functions should be called
232 // for all particles which you want to use. 94 // for all particles which you want to use.
233 // This ensures that objects of these partic 95 // This ensures that objects of these particle types will be
234 // created in the program. 96 // created in the program.
235 97
236 ConstructMyBosons(); 98 ConstructMyBosons();
237 ConstructMyLeptons(); 99 ConstructMyLeptons();
238 ConstructMyHadrons(); 100 ConstructMyHadrons();
239 ConstructMyShortLiveds(); 101 ConstructMyShortLiveds();
>> 102
240 } 103 }
241 104
>> 105
242 // construct Bosons:////////////////////////// 106 // construct Bosons://///////////////////////////////////////////////////
243 void DMXPhysicsList::ConstructMyBosons() 107 void DMXPhysicsList::ConstructMyBosons()
244 { 108 {
245 // pseudo-particles 109 // pseudo-particles
246 G4Geantino::GeantinoDefinition(); 110 G4Geantino::GeantinoDefinition();
247 G4ChargedGeantino::ChargedGeantinoDefinition 111 G4ChargedGeantino::ChargedGeantinoDefinition();
248 112
249 // gamma 113 // gamma
250 G4Gamma::GammaDefinition(); 114 G4Gamma::GammaDefinition();
251 115
252 //OpticalPhotons 116 //OpticalPhotons
253 G4OpticalPhoton::OpticalPhotonDefinition(); 117 G4OpticalPhoton::OpticalPhotonDefinition();
254 118
255 } 119 }
256 120
>> 121
257 // construct Leptons:///////////////////////// 122 // construct Leptons://///////////////////////////////////////////////////
258 void DMXPhysicsList::ConstructMyLeptons() 123 void DMXPhysicsList::ConstructMyLeptons()
259 { 124 {
260 // leptons 125 // leptons
261 G4Electron::ElectronDefinition(); 126 G4Electron::ElectronDefinition();
262 G4Positron::PositronDefinition(); 127 G4Positron::PositronDefinition();
263 G4MuonPlus::MuonPlusDefinition(); 128 G4MuonPlus::MuonPlusDefinition();
264 G4MuonMinus::MuonMinusDefinition(); 129 G4MuonMinus::MuonMinusDefinition();
265 130
266 G4NeutrinoE::NeutrinoEDefinition(); 131 G4NeutrinoE::NeutrinoEDefinition();
267 G4AntiNeutrinoE::AntiNeutrinoEDefinition(); 132 G4AntiNeutrinoE::AntiNeutrinoEDefinition();
268 G4NeutrinoMu::NeutrinoMuDefinition(); 133 G4NeutrinoMu::NeutrinoMuDefinition();
269 G4AntiNeutrinoMu::AntiNeutrinoMuDefinition() 134 G4AntiNeutrinoMu::AntiNeutrinoMuDefinition();
270 } 135 }
271 136
>> 137
>> 138 #include "G4MesonConstructor.hh"
>> 139 #include "G4BaryonConstructor.hh"
>> 140 #include "G4IonConstructor.hh"
>> 141
272 // construct Hadrons:///////////////////////// 142 // construct Hadrons://///////////////////////////////////////////////////
273 void DMXPhysicsList::ConstructMyHadrons() 143 void DMXPhysicsList::ConstructMyHadrons()
274 { 144 {
275 // mesons 145 // mesons
276 G4MesonConstructor mConstructor; 146 G4MesonConstructor mConstructor;
277 mConstructor.ConstructParticle(); 147 mConstructor.ConstructParticle();
278 148
279 // baryons 149 // baryons
280 G4BaryonConstructor bConstructor; 150 G4BaryonConstructor bConstructor;
281 bConstructor.ConstructParticle(); 151 bConstructor.ConstructParticle();
282 152
283 // ions 153 // ions
284 G4IonConstructor iConstructor; 154 G4IonConstructor iConstructor;
285 iConstructor.ConstructParticle(); 155 iConstructor.ConstructParticle();
>> 156
286 } 157 }
287 158
>> 159
288 // construct Shortliveds:///////////////////// 160 // construct Shortliveds://///////////////////////////////////////////////////
289 void DMXPhysicsList::ConstructMyShortLiveds() 161 void DMXPhysicsList::ConstructMyShortLiveds()
290 { 162 {
291 G4ShortLivedConstructor slConstructor; << 163 // ShortLiveds
292 slConstructor.ConstructParticle(); << 164 ;
293 } 165 }
294 166
>> 167
>> 168
>> 169
295 // Construct Processes /////////////////////// 170 // Construct Processes //////////////////////////////////////////////////////
296 void DMXPhysicsList::ConstructProcess() 171 void DMXPhysicsList::ConstructProcess()
297 { 172 {
>> 173
298 AddTransportation(); 174 AddTransportation();
299 175
300 ConstructEM(); 176 ConstructEM();
301 177
302 ConstructOp(); 178 ConstructOp();
303 179
304 ConstructHad(); 180 ConstructHad();
305 181
306 ConstructGeneral(); 182 ConstructGeneral();
>> 183
307 } 184 }
308 185
>> 186
309 // Transportation //////////////////////////// 187 // Transportation ///////////////////////////////////////////////////////////
>> 188 #include "DMXMaxTimeCuts.hh"
>> 189 #include "DMXMinEkineCuts.hh"
>> 190 #include "G4StepLimiter.hh"
>> 191
310 void DMXPhysicsList::AddTransportation() { 192 void DMXPhysicsList::AddTransportation() {
311 193
312 G4VUserPhysicsList::AddTransportation(); 194 G4VUserPhysicsList::AddTransportation();
313 195
314 auto particleIterator=GetParticleIterator(); << 196 theParticleIterator->reset();
315 particleIterator->reset(); << 197 while( (*theParticleIterator)() ){
316 while( (*particleIterator)() ){ << 198 G4ParticleDefinition* particle = theParticleIterator->value();
317 G4ParticleDefinition* particle = particleI <<
318 G4ProcessManager* pmanager = particle->Get 199 G4ProcessManager* pmanager = particle->GetProcessManager();
319 G4String particleName = particle->GetParti 200 G4String particleName = particle->GetParticleName();
320 // time cuts for ONLY neutrons: 201 // time cuts for ONLY neutrons:
321 if(particleName == "neutron") 202 if(particleName == "neutron")
322 pmanager->AddDiscreteProcess(new DMXMaxT 203 pmanager->AddDiscreteProcess(new DMXMaxTimeCuts());
323 // Energy cuts to kill charged (embedded i 204 // Energy cuts to kill charged (embedded in method) particles:
324 pmanager->AddDiscreteProcess(new DMXMinEki 205 pmanager->AddDiscreteProcess(new DMXMinEkineCuts());
325 206
326 // Step limit applied to all particles: 207 // Step limit applied to all particles:
327 pmanager->AddDiscreteProcess(new G4StepLim << 208 pmanager->AddProcess(new G4StepLimiter, -1,-1,1);
>> 209
328 } 210 }
329 } 211 }
330 212
>> 213
331 // Electromagnetic Processes ///////////////// 214 // Electromagnetic Processes ////////////////////////////////////////////////
332 // all charged particles 215 // all charged particles
>> 216 #include "G4MultipleScattering.hh"
>> 217
>> 218 // gamma
>> 219 #include "G4LowEnergyRayleigh.hh"
>> 220 #include "G4LowEnergyPhotoElectric.hh"
>> 221 #include "G4LowEnergyCompton.hh"
>> 222 #include "G4LowEnergyGammaConversion.hh"
>> 223
>> 224
>> 225 // e-
>> 226 #include "G4LowEnergyIonisation.hh"
>> 227 #include "G4LowEnergyBremsstrahlung.hh"
>> 228
>> 229 // e+
>> 230 #include "G4eIonisation.hh"
>> 231 #include "G4eBremsstrahlung.hh"
>> 232 #include "G4eplusAnnihilation.hh"
>> 233
>> 234
>> 235 // alpha and GenericIon and deuterons, triton, He3:
>> 236 //hIonisation #include "G4hLowEnergyIonisation.hh" -> moved to G4hIonisation
>> 237 #include "G4EnergyLossTables.hh"
>> 238 // hLowEnergyIonisation uses Ziegler 1988 as the default
>> 239
>> 240
>> 241 //muon:
>> 242 #include "G4MuIonisation.hh"
>> 243 #include "G4MuBremsstrahlung.hh"
>> 244 #include "G4MuPairProduction.hh"
>> 245 #include "G4MuonMinusCaptureAtRest.hh"
>> 246
>> 247 //OTHERS:
>> 248 #include "G4hIonisation.hh" // standard hadron ionisation
>> 249
>> 250 //em process options to allow msc step-limitation to be switched off
>> 251 #include "G4EmProcessOptions.hh"
>> 252
333 void DMXPhysicsList::ConstructEM() { 253 void DMXPhysicsList::ConstructEM() {
334 <<
335 G4LossTableManager* man = G4LossTableManager <<
336 man->SetAtomDeexcitation(new G4UAtomicDeexci <<
337 254
338 G4EmParameters* em_params = G4EmParameters:: << 255 // processes:
>> 256
>> 257 G4LowEnergyPhotoElectric* lowePhot = new G4LowEnergyPhotoElectric();
>> 258 G4LowEnergyIonisation* loweIon = new G4LowEnergyIonisation();
>> 259 G4LowEnergyBremsstrahlung* loweBrem = new G4LowEnergyBremsstrahlung();
>> 260
>> 261 // note LowEIon uses proton as basis for its data-base, therefore
>> 262 // cannot specify different LowEnergyIonisation models for different
>> 263 // particles, but can change model globally for Ion, Alpha and Proton.
>> 264
>> 265
>> 266 //fluorescence apply specific cut for fluorescence from photons, electrons
>> 267 //and bremsstrahlung photons:
>> 268 G4double fluorcut = 250*eV;
>> 269 lowePhot->SetCutForLowEnSecPhotons(fluorcut);
>> 270 loweIon->SetCutForLowEnSecPhotons(fluorcut);
>> 271 loweBrem->SetCutForLowEnSecPhotons(fluorcut);
339 272
340 auto particleIterator=GetParticleIterator(); << 273 // setting tables explicitly for electronic stopping power
341 particleIterator->reset(); << 274 // ahadronLowEIon->SetElectronicStoppingPowerModel
342 while( (*particleIterator)() ){ << 275 // (G4GenericIon::GenericIonDefinition(), "ICRU_R49p") ;
343 G4ParticleDefinition* particle = particleI << 276 // ahadronLowEIon->SetElectronicStoppingPowerModel
>> 277 // (G4Proton::ProtonDefinition(), "ICRU_R49p") ;
>> 278
>> 279 // Switch off the Barkas and Bloch corrections
>> 280 // ahadronLowEIon->SetBarkasOff();
>> 281
>> 282
>> 283 theParticleIterator->reset();
>> 284 while( (*theParticleIterator)() ){
>> 285 G4ParticleDefinition* particle = theParticleIterator->value();
344 G4ProcessManager* pmanager = particle->Get 286 G4ProcessManager* pmanager = particle->GetProcessManager();
345 G4String particleName = particle->GetParti 287 G4String particleName = particle->GetParticleName();
346 G4String particleType = particle->GetParti 288 G4String particleType = particle->GetParticleType();
347 G4double charge = particle->GetPDGCharge() 289 G4double charge = particle->GetPDGCharge();
348 290
349 if (particleName == "gamma") 291 if (particleName == "gamma")
350 { 292 {
351 //gamma 293 //gamma
352 G4RayleighScattering* theRayleigh = new G4Ra << 294 pmanager->AddDiscreteProcess(new G4LowEnergyRayleigh());
353 pmanager->AddDiscreteProcess(theRayleigh); << 295 pmanager->AddDiscreteProcess(lowePhot);
354 << 296 pmanager->AddDiscreteProcess(new G4LowEnergyCompton());
355 G4PhotoElectricEffect* thePhotoElectricEffec << 297 pmanager->AddDiscreteProcess(new G4LowEnergyGammaConversion());
356 thePhotoElectricEffect->SetEmModel(new G4Liv <<
357 pmanager->AddDiscreteProcess(thePhotoElectri <<
358 <<
359 G4ComptonScattering* theComptonScattering = <<
360 theComptonScattering->SetEmModel(new G4Liver <<
361 pmanager->AddDiscreteProcess(theComptonScatt <<
362 <<
363 G4GammaConversion* theGammaConversion = new <<
364 theGammaConversion->SetEmModel(new G4BetheHe <<
365 pmanager->AddDiscreteProcess(theGammaConvers <<
366 <<
367 } 298 }
368 else if (particleName == "e-") 299 else if (particleName == "e-")
369 { 300 {
370 //electron 301 //electron
371 // process ordering: AddProcess(name, at res 302 // process ordering: AddProcess(name, at rest, along step, post step)
372 // Multiple scattering << 303 // -1 = not implemented, then ordering
373 G4eMultipleScattering* msc = new G4eMultiple << 304 G4MultipleScattering* aMultipleScattering = new G4MultipleScattering();
374 em_params->SetMscStepLimitType(fUseDistanceT << 305 pmanager->AddProcess(aMultipleScattering, -1, 1, 1);
375 pmanager->AddProcess(msc,-1, 1, -1); << 306 pmanager->AddProcess(loweIon, -1, 2, 2);
376 << 307 pmanager->AddProcess(loweBrem, -1,-1, 3);
377 // Ionisation <<
378 G4eIonisation* eIonisation = new G4eIonisati <<
379 G4VEmModel* theIoniLiv = new G4Livermo <<
380 theIoniLiv->SetHighEnergyLimit(0.1*MeV <<
381 eIonisation->AddEmModel(0, theIoniLiv, <<
382 em_params->SetStepFunction(0.2, 100*um); //i <<
383 pmanager->AddProcess(eIonisation,-1, 2, 1); <<
384 <<
385 // Bremsstrahlung <<
386 G4eBremsstrahlung* eBremsstrahlung = new G4e <<
387 pmanager->AddProcess(eBremsstrahlung, -1,-3, <<
388 } 308 }
389 else if (particleName == "e+") 309 else if (particleName == "e+")
390 { 310 {
391 //positron << 311 //positron
392 G4eMultipleScattering* msc = new G4eMultiple << 312 G4MultipleScattering* aMultipleScattering = new G4MultipleScattering();
393 msc->SetStepLimitType(fUseDistanceToBoundary << 313 pmanager->AddProcess(aMultipleScattering, -1, 1, 1);
394 pmanager->AddProcess(msc,-1, 1, -1); << 314 pmanager->AddProcess(new G4eIonisation(), -1, 2, 2);
395 << 315 pmanager->AddProcess(new G4eBremsstrahlung(), -1,-1, 3);
396 // Ionisation << 316 pmanager->AddProcess(new G4eplusAnnihilation(),0,-1, 4);
397 G4eIonisation* eIonisation = new G4eIonisati <<
398 // eIonisation->SetStepFunction(0.2, 100*um) <<
399 pmanager->AddProcess(eIonisation, <<
400 <<
401 //Bremsstrahlung (use default, no low-energy <<
402 pmanager->AddProcess(new G4eBremsstrahlung() <<
403 <<
404 //Annihilation <<
405 pmanager->AddProcess(new G4eplusAnnihilation <<
406 } 317 }
407 else if( particleName == "mu+" || 318 else if( particleName == "mu+" ||
408 particleName == "mu-" ) 319 particleName == "mu-" )
409 { 320 {
410 //muon 321 //muon
411 pmanager->AddProcess(new G4MuMultipleScatter << 322 G4MultipleScattering* aMultipleScattering = new G4MultipleScattering();
412 pmanager->AddProcess(new G4MuIonisation(), << 323 pmanager->AddProcess(aMultipleScattering, -1, 1, 1);
413 pmanager->AddProcess(new G4MuBremsstrahlung( << 324 pmanager->AddProcess(new G4MuIonisation(), -1, 2, 2);
414 pmanager->AddProcess(new G4MuPairProduction( << 325 pmanager->AddProcess(new G4MuBremsstrahlung(), -1,-1, 3);
>> 326 pmanager->AddProcess(new G4MuPairProduction(), -1,-1, 4);
415 if( particleName == "mu-" ) 327 if( particleName == "mu-" )
416 pmanager->AddProcess(new G4MuonMinusCaptur << 328 pmanager->AddProcess(new G4MuonMinusCaptureAtRest(), 0,-1,-1);
417 } 329 }
418 else if (particleName == "proton" || << 330 else if (particleName == "proton" ||
419 particleName == "pi+" || << 331 particleName == "alpha" ||
420 particleName == "pi-") <<
421 { <<
422 //multiple scattering <<
423 pmanager->AddProcess(new G4hMultipleScatteri <<
424 <<
425 //ionisation <<
426 G4hIonisation* hIonisation = new G4hIonisati <<
427 em_params->SetStepFunctionMuHad(0.2, 50*um); <<
428 pmanager->AddProcess(hIonisation, <<
429 <<
430 //bremmstrahlung <<
431 pmanager->AddProcess(new G4hBremsstrahlung, <<
432 } <<
433 else if(particleName == "alpha" || <<
434 particleName == "deuteron" || 332 particleName == "deuteron" ||
435 particleName == "triton" || 333 particleName == "triton" ||
436 particleName == "He3") << 334 particleName == "He3" ||
437 { << 335 particleName == "GenericIon" ||
438 //multiple scattering << 336 (particleType == "nucleus" && charge != 0))
439 pmanager->AddProcess(new G4hMultipleScatteri <<
440 <<
441 //ionisation <<
442 G4ionIonisation* ionIoni = new G4ionIonisati <<
443 em_params->SetStepFunctionLightIons(0.1, 1*C <<
444 pmanager->AddProcess(ionIoni, <<
445 pmanager->AddProcess(new G4NuclearStopping() <<
446 } <<
447 else if (particleName == "GenericIon") <<
448 { 337 {
449 // OBJECT may be dynamically created as eith 338 // OBJECT may be dynamically created as either a GenericIon or nucleus
450 // G4Nucleus exists and therefore has partic 339 // G4Nucleus exists and therefore has particle type nucleus
451 // genericIon: 340 // genericIon:
452 << 341 G4MultipleScattering* aMultipleScattering = new G4MultipleScattering();
453 //multiple scattering << 342 //hIonisation G4hLowEnergyIonisation* ahadronLowEIon = new G4hLowEnergyIonisation();
454 pmanager->AddProcess(new G4hMultipleScatteri << 343 G4hIonisation* ahadronIon = new G4hIonisation();
455 << 344 pmanager->AddProcess(aMultipleScattering,-1,1,1);
456 //ionisation << 345 //hIonisation pmanager->AddProcess(ahadronLowEIon,-1,2,2);
457 G4ionIonisation* ionIoni = new G4ionIonisati << 346 pmanager->AddProcess(ahadronIon,-1,2,2);
458 em_params->SetStepFunctionIons(0.1, 1*CLHEP: << 347 // ahadronLowEIon->SetNuclearStoppingOff() ;
459 pmanager->AddProcess(ionIoni, << 348 // ahadronLowEIon->SetNuclearStoppingPowerModel("ICRU_R49") ;
460 pmanager->AddProcess(new G4NuclearStopping() << 349 // ahadronLowEIon->SetNuclearStoppingOn() ;
>> 350
>> 351 //fluorescence switch off for hadrons (for now) PIXE:
>> 352 //hIonisation ahadronLowEIon->SetFluorescence(false);
461 } 353 }
462 <<
463 else if ((!particle->IsShortLived()) && 354 else if ((!particle->IsShortLived()) &&
464 (charge != 0.0) && 355 (charge != 0.0) &&
465 (particle->GetParticleName() != "charge 356 (particle->GetParticleName() != "chargedgeantino"))
466 { 357 {
467 //all others charged particles except geanti 358 //all others charged particles except geantino
468 G4hMultipleScattering* aMultipleScatte << 359 G4MultipleScattering* aMultipleScattering = new G4MultipleScattering();
>> 360 //hIonisation G4hLowEnergyIonisation* ahadronLowEIon = new G4hLowEnergyIonisation();
469 G4hIonisation* ahadronIon = new G4hIon 361 G4hIonisation* ahadronIon = new G4hIonisation();
470 << 362 pmanager->AddProcess(aMultipleScattering,-1,1,1);
471 //multiple scattering << 363 //hIonisation pmanager->AddProcess(ahadronLowEIon, -1,2,2);
472 pmanager->AddProcess(aMultipleScattering,-1, << 364 pmanager->AddProcess(ahadronIon, -1,2,2);
473 << 365 // pmanager->AddProcess(new G4hIonisation(), -1,2,2);
474 //ionisation <<
475 pmanager->AddProcess(ahadronIon, -1,2, <<
476 } 366 }
>> 367
477 } 368 }
>> 369
>> 370 // turn off msc step-limitation - especially as electron cut 1nm
>> 371 G4EmProcessOptions opt;
>> 372 // opt.SetMscStepLimitation(false);
>> 373 opt.SetMscStepLimitation(fMinimal);
>> 374
478 } 375 }
479 376
>> 377
480 // Optical Processes ///////////////////////// 378 // Optical Processes ////////////////////////////////////////////////////////
>> 379 #include "G4Scintillation.hh"
>> 380 #include "G4OpAbsorption.hh"
>> 381 //#include "G4OpRayleigh.hh"
>> 382 #include "G4OpBoundaryProcess.hh"
>> 383
481 void DMXPhysicsList::ConstructOp() 384 void DMXPhysicsList::ConstructOp()
482 { 385 {
483 G4OpticalParameters* opParams = G4OpticalPar << 386 // default scintillation process
484 G4Scintillation* theScintProcessDef = new G4 387 G4Scintillation* theScintProcessDef = new G4Scintillation("Scintillation");
485 opParams->SetScintTrackSecondariesFirst(true << 388 // theScintProcessDef->DumpPhysicsTable();
486 opParams->SetScintByParticleType(true); << 389 theScintProcessDef->SetTrackSecondariesFirst(true);
>> 390 theScintProcessDef->SetScintillationYieldFactor(1.0); //
>> 391 theScintProcessDef->SetScintillationExcitationRatio(0.0); //
>> 392 theScintProcessDef->SetVerboseLevel(OpVerbLevel);
>> 393
>> 394 // scintillation process for alpha:
>> 395 G4Scintillation* theScintProcessAlpha = new G4Scintillation("Scintillation");
>> 396 // theScintProcessNuc->DumpPhysicsTable();
>> 397 theScintProcessAlpha->SetTrackSecondariesFirst(true);
>> 398 theScintProcessAlpha->SetScintillationYieldFactor(1.1);
>> 399 theScintProcessAlpha->SetScintillationExcitationRatio(1.0);
>> 400 theScintProcessAlpha->SetVerboseLevel(OpVerbLevel);
>> 401
>> 402 // scintillation process for heavy nuclei
>> 403 G4Scintillation* theScintProcessNuc = new G4Scintillation("Scintillation");
>> 404 // theScintProcessNuc->DumpPhysicsTable();
>> 405 theScintProcessNuc->SetTrackSecondariesFirst(true);
>> 406 theScintProcessNuc->SetScintillationYieldFactor(0.2);
>> 407 theScintProcessNuc->SetScintillationExcitationRatio(1.0);
>> 408 theScintProcessNuc->SetVerboseLevel(OpVerbLevel);
487 409
488 // optical processes 410 // optical processes
489 G4OpAbsorption* theAbsorptionProcess = new G 411 G4OpAbsorption* theAbsorptionProcess = new G4OpAbsorption();
>> 412 // G4OpRayleigh* theRayleighScatteringProcess = new G4OpRayleigh();
490 G4OpBoundaryProcess* theBoundaryProcess = ne 413 G4OpBoundaryProcess* theBoundaryProcess = new G4OpBoundaryProcess();
>> 414 // theAbsorptionProcess->DumpPhysicsTable();
>> 415 // theRayleighScatteringProcess->DumpPhysicsTable();
>> 416 theAbsorptionProcess->SetVerboseLevel(OpVerbLevel);
>> 417 // theRayleighScatteringProcess->SetVerboseLevel(OpVerbLevel);
>> 418 theBoundaryProcess->SetVerboseLevel(OpVerbLevel);
>> 419 G4OpticalSurfaceModel themodel = unified;
>> 420 theBoundaryProcess->SetModel(themodel);
491 421
492 auto particleIterator=GetParticleIterator(); << 422 theParticleIterator->reset();
493 particleIterator->reset(); << 423 while( (*theParticleIterator)() )
494 while( (*particleIterator)() ) <<
495 { 424 {
496 G4ParticleDefinition* particle = particl << 425 G4ParticleDefinition* particle = theParticleIterator->value();
497 G4ProcessManager* pmanager = particle->G 426 G4ProcessManager* pmanager = particle->GetProcessManager();
498 G4String particleName = particle->GetPar 427 G4String particleName = particle->GetParticleName();
499 if (theScintProcessDef->IsApplicable(*pa 428 if (theScintProcessDef->IsApplicable(*particle)) {
500 pmanager->AddProcess(theScintProcessDe << 429 // if(particle->GetPDGMass() > 5.0*GeV)
501 pmanager->SetProcessOrderingToLast(the << 430 if(particle->GetParticleName() == "GenericIon") {
502 pmanager->SetProcessOrderingToLast(the << 431 pmanager->AddProcess(theScintProcessNuc); // AtRestDiscrete
>> 432 pmanager->SetProcessOrderingToLast(theScintProcessNuc,idxAtRest);
>> 433 pmanager->SetProcessOrderingToLast(theScintProcessNuc,idxPostStep);
>> 434 }
>> 435 else if(particle->GetParticleName() == "alpha") {
>> 436 pmanager->AddProcess(theScintProcessAlpha);
>> 437 pmanager->SetProcessOrderingToLast(theScintProcessAlpha,idxAtRest);
>> 438 pmanager->SetProcessOrderingToLast(theScintProcessAlpha,idxPostStep);
>> 439 }
>> 440 else {
>> 441 pmanager->AddProcess(theScintProcessDef);
>> 442 pmanager->SetProcessOrderingToLast(theScintProcessDef,idxAtRest);
>> 443 pmanager->SetProcessOrderingToLast(theScintProcessDef,idxPostStep);
>> 444 }
503 } 445 }
504 446
505 if (particleName == "opticalphoton") { 447 if (particleName == "opticalphoton") {
506 pmanager->AddDiscreteProcess(theAbsorptionPr 448 pmanager->AddDiscreteProcess(theAbsorptionProcess);
>> 449 // pmanager->AddDiscreteProcess(theRayleighScatteringProcess);
507 pmanager->AddDiscreteProcess(theBoundaryProc 450 pmanager->AddDiscreteProcess(theBoundaryProcess);
508 } 451 }
509 } 452 }
510 } 453 }
511 454
>> 455
512 // Hadronic processes //////////////////////// 456 // Hadronic processes ////////////////////////////////////////////////////////
513 457
>> 458 // Elastic processes:
>> 459 #include "G4HadronElasticProcess.hh"
>> 460
>> 461 // Inelastic processes:
>> 462 #include "G4PionPlusInelasticProcess.hh"
>> 463 #include "G4PionMinusInelasticProcess.hh"
>> 464 #include "G4KaonPlusInelasticProcess.hh"
>> 465 #include "G4KaonZeroSInelasticProcess.hh"
>> 466 #include "G4KaonZeroLInelasticProcess.hh"
>> 467 #include "G4KaonMinusInelasticProcess.hh"
>> 468 #include "G4ProtonInelasticProcess.hh"
>> 469 #include "G4AntiProtonInelasticProcess.hh"
>> 470 #include "G4NeutronInelasticProcess.hh"
>> 471 #include "G4AntiNeutronInelasticProcess.hh"
>> 472 #include "G4DeuteronInelasticProcess.hh"
>> 473 #include "G4TritonInelasticProcess.hh"
>> 474 #include "G4AlphaInelasticProcess.hh"
>> 475
>> 476 // Low-energy Models: < 20GeV
>> 477 #include "G4LElastic.hh"
>> 478 #include "G4LEPionPlusInelastic.hh"
>> 479 #include "G4LEPionMinusInelastic.hh"
>> 480 #include "G4LEKaonPlusInelastic.hh"
>> 481 #include "G4LEKaonZeroSInelastic.hh"
>> 482 #include "G4LEKaonZeroLInelastic.hh"
>> 483 #include "G4LEKaonMinusInelastic.hh"
>> 484 #include "G4LEProtonInelastic.hh"
>> 485 #include "G4LEAntiProtonInelastic.hh"
>> 486 #include "G4LENeutronInelastic.hh"
>> 487 #include "G4LEAntiNeutronInelastic.hh"
>> 488 #include "G4LEDeuteronInelastic.hh"
>> 489 #include "G4LETritonInelastic.hh"
>> 490 #include "G4LEAlphaInelastic.hh"
>> 491 #include "G4HadronCaptureProcess.hh"
>> 492 // High-energy Models: >20 GeV
>> 493 #include "G4HEPionPlusInelastic.hh"
>> 494 #include "G4HEPionMinusInelastic.hh"
>> 495 #include "G4HEKaonPlusInelastic.hh"
>> 496 #include "G4HEKaonZeroInelastic.hh"
>> 497 #include "G4HEKaonZeroInelastic.hh"
>> 498 #include "G4HEKaonMinusInelastic.hh"
>> 499 #include "G4HEProtonInelastic.hh"
>> 500 #include "G4HEAntiProtonInelastic.hh"
>> 501 #include "G4HENeutronInelastic.hh"
>> 502 #include "G4HEAntiNeutronInelastic.hh"
>> 503
>> 504 // Neutron high-precision models: <20 MeV
>> 505 #include "G4NeutronHPElastic.hh"
>> 506 #include "G4NeutronHPElasticData.hh"
>> 507 #include "G4NeutronHPCapture.hh"
>> 508 #include "G4NeutronHPCaptureData.hh"
>> 509 #include "G4NeutronHPInelastic.hh"
>> 510 #include "G4NeutronHPInelasticData.hh"
>> 511 #include "G4LCapture.hh"
>> 512
>> 513 // Stopping processes
>> 514 #include "G4PiMinusAbsorptionAtRest.hh"
>> 515 #include "G4KaonMinusAbsorptionAtRest.hh"
>> 516 #include "G4AntiProtonAnnihilationAtRest.hh"
>> 517 #include "G4AntiNeutronAnnihilationAtRest.hh"
>> 518
>> 519
>> 520 // ConstructHad()
>> 521 // Makes discrete physics processes for the hadrons, at present limited
>> 522 // to those particles with GHEISHA interactions (INTRC > 0).
>> 523 // The processes are: Elastic scattering and Inelastic scattering.
>> 524 // F.W.Jones 09-JUL-1998
514 void DMXPhysicsList::ConstructHad() 525 void DMXPhysicsList::ConstructHad()
515 { 526 {
516 //Elastic models << 527 G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
517 G4HadronElastic* elastic_lhep0 = new G4Hadro << 528 G4LElastic* theElasticModel = new G4LElastic;
518 G4ChipsElasticModel* elastic_chip = new G4Ch << 529 theElasticProcess->RegisterMe(theElasticModel);
519 G4ElasticHadrNucleusHE* elastic_he = new G4E <<
520 530
521 // Inelastic scattering << 531 theParticleIterator->reset();
522 const G4double theFTFMin0 = 0.0*GeV; << 532 while ((*theParticleIterator)())
523 const G4double theFTFMin1 = 3.0*GeV; <<
524 const G4double theFTFMax = G4HadronicParamet <<
525 const G4double theBERTMin0 = 0.0*GeV; <<
526 const G4double theBERTMin1 = 19.0*MeV; <<
527 const G4double theBERTMax = 6.0*GeV; <<
528 const G4double theHPMin = 0.0*GeV; <<
529 const G4double theHPMax = 20.0*MeV; <<
530 <<
531 G4FTFModel * theStringModel = new G4FTFModel <<
532 G4ExcitedStringDecay * theStringDecay = new <<
533 theStringModel->SetFragmentationModel( theSt <<
534 G4PreCompoundModel * thePreEquilib = new G4P <<
535 G4GeneratorPrecompoundInterface * theCascade <<
536 <<
537 G4TheoFSGenerator * theFTFModel0 = new G4The <<
538 theFTFModel0->SetHighEnergyGenerator( theStr <<
539 theFTFModel0->SetTransport( theCascade ); <<
540 theFTFModel0->SetMinEnergy( theFTFMin0 ); <<
541 theFTFModel0->SetMaxEnergy( theFTFMax ); <<
542 <<
543 G4TheoFSGenerator * theFTFModel1 = new G4The <<
544 theFTFModel1->SetHighEnergyGenerator( theStr <<
545 theFTFModel1->SetTransport( theCascade ); <<
546 theFTFModel1->SetMinEnergy( theFTFMin1 ); <<
547 theFTFModel1->SetMaxEnergy( theFTFMax ); <<
548 <<
549 G4CascadeInterface * theBERTModel0 = new G4C <<
550 theBERTModel0->SetMinEnergy( theBERTMin0 ); <<
551 theBERTModel0->SetMaxEnergy( theBERTMax ); <<
552 <<
553 G4CascadeInterface * theBERTModel1 = new G4C <<
554 theBERTModel1->SetMinEnergy( theBERTMin1 ); <<
555 theBERTModel1->SetMaxEnergy( theBERTMax ); <<
556 <<
557 G4VCrossSectionDataSet * theAntiNucleonData <<
558 G4ComponentGGNuclNuclXsc * ggNuclNuclXsec = <<
559 G4VCrossSectionDataSet * theGGNuclNuclData = <<
560 G4VCrossSectionDataSet * theGGNNEl = new G4C <<
561 G4ComponentGGHadronNucleusXsc * ggHNXsec = n <<
562 G4VCrossSectionDataSet * theGGHNEl = new G4C <<
563 G4VCrossSectionDataSet * theGGHNInel = new G <<
564 <<
565 auto particleIterator=GetParticleIterator(); <<
566 particleIterator->reset(); <<
567 while ((*particleIterator)()) <<
568 { 533 {
569 G4ParticleDefinition* particle = particl << 534 G4ParticleDefinition* particle = theParticleIterator->value();
570 G4ProcessManager* pmanager = particle->G 535 G4ProcessManager* pmanager = particle->GetProcessManager();
571 G4String particleName = particle->GetPar 536 G4String particleName = particle->GetParticleName();
572 537
573 if (particleName == "pi+") 538 if (particleName == "pi+")
574 { 539 {
575 // Elastic scattering << 540 pmanager->AddDiscreteProcess(theElasticProcess);
576 G4HadronElasticProcess* theElasticPr << 541 G4PionPlusInelasticProcess* theInelasticProcess =
577 theElasticProcess->AddDataSet( new G << 542 new G4PionPlusInelasticProcess("inelastic");
578 theElasticProcess->RegisterMe( elast << 543 G4LEPionPlusInelastic* theLEInelasticModel =
579 pmanager->AddDiscreteProcess( theElasticPr << 544 new G4LEPionPlusInelastic;
580 //Inelastic scattering << 545 theInelasticProcess->RegisterMe(theLEInelasticModel);
581 G4HadronInelasticProcess* theInelasticProc << 546 G4HEPionPlusInelastic* theHEInelasticModel =
582 new G4HadronInelasticProcess( "inelastic << 547 new G4HEPionPlusInelastic;
583 theInelasticProcess->AddDataSet( new G4BGG << 548 theInelasticProcess->RegisterMe(theHEInelasticModel);
584 theInelasticProcess->RegisterMe( theFTFMod << 549 pmanager->AddDiscreteProcess(theInelasticProcess);
585 theInelasticProcess->RegisterMe( the <<
586 pmanager->AddDiscreteProcess( theInelastic <<
587 } 550 }
588 551
589 else if (particleName == "pi-") 552 else if (particleName == "pi-")
590 { 553 {
591 // Elastic scattering << 554 pmanager->AddDiscreteProcess(theElasticProcess);
592 G4HadronElasticProcess* theElasticPr << 555 G4PionMinusInelasticProcess* theInelasticProcess =
593 theElasticProcess->AddDataSet( new G << 556 new G4PionMinusInelasticProcess("inelastic");
594 theElasticProcess->RegisterMe( elast << 557 G4LEPionMinusInelastic* theLEInelasticModel =
595 pmanager->AddDiscreteProcess( theElasticPr << 558 new G4LEPionMinusInelastic;
596 //Inelastic scattering << 559 theInelasticProcess->RegisterMe(theLEInelasticModel);
597 G4HadronInelasticProcess* theInelasticProc << 560 G4HEPionMinusInelastic* theHEInelasticModel =
598 new G4HadronInelasticProcess( "inelastic << 561 new G4HEPionMinusInelastic;
599 theInelasticProcess->AddDataSet( new G4BGG << 562 theInelasticProcess->RegisterMe(theHEInelasticModel);
600 theInelasticProcess->RegisterMe( theFTFMod << 563 pmanager->AddDiscreteProcess(theInelasticProcess);
601 theInelasticProcess->RegisterMe( the << 564 G4String prcNam;
602 pmanager->AddDiscreteProcess( theInelastic << 565 pmanager->AddRestProcess(new G4PiMinusAbsorptionAtRest, ordDefault);
603 //Absorption <<
604 pmanager->AddRestProcess(new G4HadronicAbs <<
605 } 566 }
>> 567
606 else if (particleName == "kaon+") 568 else if (particleName == "kaon+")
607 { 569 {
608 // Elastic scattering << 570 pmanager->AddDiscreteProcess(theElasticProcess);
609 G4HadronElasticProcess* theElasticPr << 571 G4KaonPlusInelasticProcess* theInelasticProcess =
610 theElasticProcess->AddDataSet( theGGHNEl ) << 572 new G4KaonPlusInelasticProcess("inelastic");
611 theElasticProcess->RegisterMe( elast << 573 G4LEKaonPlusInelastic* theLEInelasticModel =
612 pmanager->AddDiscreteProcess( theElasticPr << 574 new G4LEKaonPlusInelastic;
613 // Inelastic scattering << 575 theInelasticProcess->RegisterMe(theLEInelasticModel);
614 G4HadronInelasticProcess* theInelasticProc << 576 G4HEKaonPlusInelastic* theHEInelasticModel =
615 new G4HadronInelasticProcess( "inelastic << 577 new G4HEKaonPlusInelastic;
616 theInelasticProcess->AddDataSet( theGGHNIn << 578 theInelasticProcess->RegisterMe(theHEInelasticModel);
617 theInelasticProcess->RegisterMe( theFTFMod << 579 pmanager->AddDiscreteProcess(theInelasticProcess);
618 theInelasticProcess->RegisterMe( the << 580 }
619 pmanager->AddDiscreteProcess( theInelastic << 581
620 } <<
621 else if (particleName == "kaon0S") 582 else if (particleName == "kaon0S")
622 { 583 {
623 // Elastic scattering << 584 pmanager->AddDiscreteProcess(theElasticProcess);
624 G4HadronElasticProcess* theElasticPr << 585 G4KaonZeroSInelasticProcess* theInelasticProcess =
625 theElasticProcess->AddDataSet( theGGHNEl ) << 586 new G4KaonZeroSInelasticProcess("inelastic");
626 theElasticProcess->RegisterMe( elast << 587 G4LEKaonZeroSInelastic* theLEInelasticModel =
627 pmanager->AddDiscreteProcess( theElasticPr << 588 new G4LEKaonZeroSInelastic;
628 // Inelastic scattering << 589 theInelasticProcess->RegisterMe(theLEInelasticModel);
629 G4HadronInelasticProcess* theInelasticProc << 590 G4HEKaonZeroInelastic* theHEInelasticModel =
630 new G4HadronInelasticProcess( "inelastic << 591 new G4HEKaonZeroInelastic;
631 theInelasticProcess->AddDataSet( theGGHNIn << 592 theInelasticProcess->RegisterMe(theHEInelasticModel);
632 theInelasticProcess->RegisterMe( theFTFMod << 593 pmanager->AddDiscreteProcess(theInelasticProcess);
633 theInelasticProcess->RegisterMe( the <<
634 pmanager->AddDiscreteProcess( theInelastic <<
635 } 594 }
636 595
637 else if (particleName == "kaon0L") 596 else if (particleName == "kaon0L")
638 { 597 {
639 // Elastic scattering << 598 pmanager->AddDiscreteProcess(theElasticProcess);
640 G4HadronElasticProcess* theElasticPr << 599 G4KaonZeroLInelasticProcess* theInelasticProcess =
641 theElasticProcess->AddDataSet( theGGHNEl ) << 600 new G4KaonZeroLInelasticProcess("inelastic");
642 theElasticProcess->RegisterMe( elast << 601 G4LEKaonZeroLInelastic* theLEInelasticModel =
643 pmanager->AddDiscreteProcess( theElasticPr << 602 new G4LEKaonZeroLInelastic;
644 // Inelastic scattering << 603 theInelasticProcess->RegisterMe(theLEInelasticModel);
645 G4HadronInelasticProcess* theInelasticProc << 604 G4HEKaonZeroInelastic* theHEInelasticModel =
646 new G4HadronInelasticProcess( "inelastic << 605 new G4HEKaonZeroInelastic;
647 theInelasticProcess->AddDataSet( theGGHNIn << 606 theInelasticProcess->RegisterMe(theHEInelasticModel);
648 theInelasticProcess->RegisterMe( theFTFMod << 607 pmanager->AddDiscreteProcess(theInelasticProcess);
649 theInelasticProcess->RegisterMe( the <<
650 pmanager->AddDiscreteProcess( theInelastic <<
651 } 608 }
652 609
653 else if (particleName == "kaon-") 610 else if (particleName == "kaon-")
654 { 611 {
655 // Elastic scattering << 612 pmanager->AddDiscreteProcess(theElasticProcess);
656 G4HadronElasticProcess* theElasticPr << 613 G4KaonMinusInelasticProcess* theInelasticProcess =
657 theElasticProcess->AddDataSet( theGGHNEl ) << 614 new G4KaonMinusInelasticProcess("inelastic");
658 theElasticProcess->RegisterMe( elast << 615 G4LEKaonMinusInelastic* theLEInelasticModel =
659 pmanager->AddDiscreteProcess( theElasticPr << 616 new G4LEKaonMinusInelastic;
660 // Inelastic scattering << 617 theInelasticProcess->RegisterMe(theLEInelasticModel);
661 G4HadronInelasticProcess* theInelasticProc << 618 G4HEKaonMinusInelastic* theHEInelasticModel =
662 new G4HadronInelasticProcess( "inelastic << 619 new G4HEKaonMinusInelastic;
663 theInelasticProcess->AddDataSet( the << 620 theInelasticProcess->RegisterMe(theHEInelasticModel);
664 theInelasticProcess->RegisterMe( theFTFMod << 621 pmanager->AddDiscreteProcess(theInelasticProcess);
665 theInelasticProcess->RegisterMe( the << 622 pmanager->AddRestProcess(new G4KaonMinusAbsorptionAtRest, ordDefault);
666 pmanager->AddDiscreteProcess( theInelastic <<
667 pmanager->AddRestProcess(new G4HadronicAbs <<
668 } 623 }
669 624
670 else if (particleName == "proton") 625 else if (particleName == "proton")
671 { 626 {
672 // Elastic scattering << 627 pmanager->AddDiscreteProcess(theElasticProcess);
673 G4HadronElasticProcess* theElasticPr << 628 G4ProtonInelasticProcess* theInelasticProcess =
674 theElasticProcess->AddDataSet( new G << 629 new G4ProtonInelasticProcess("inelastic");
675 theElasticProcess->RegisterMe( elast << 630 G4LEProtonInelastic* theLEInelasticModel = new G4LEProtonInelastic;
676 pmanager->AddDiscreteProcess( theElasticPr << 631 theInelasticProcess->RegisterMe(theLEInelasticModel);
677 // Inelastic scattering << 632 G4HEProtonInelastic* theHEInelasticModel = new G4HEProtonInelastic;
678 G4HadronInelasticProcess* theInelasticProc << 633 theInelasticProcess->RegisterMe(theHEInelasticModel);
679 new G4HadronInelasticProcess( "inelastic << 634 pmanager->AddDiscreteProcess(theInelasticProcess);
680 theInelasticProcess->AddDataSet( new G4BGG <<
681 theInelasticProcess->RegisterMe( theFTFMod <<
682 theInelasticProcess->RegisterMe( the <<
683 pmanager->AddDiscreteProcess( theInelastic <<
684 } 635 }
>> 636
685 else if (particleName == "anti_proton") 637 else if (particleName == "anti_proton")
686 { 638 {
687 // Elastic scattering << 639 pmanager->AddDiscreteProcess(theElasticProcess);
688 const G4double elastic_elimitAntiNuc << 640 G4AntiProtonInelasticProcess* theInelasticProcess =
689 G4AntiNuclElastic* elastic_anuc = ne << 641 new G4AntiProtonInelasticProcess("inelastic");
690 elastic_anuc->SetMinEnergy( elastic_ << 642 G4LEAntiProtonInelastic* theLEInelasticModel =
691 G4CrossSectionElastic* elastic_anucx << 643 new G4LEAntiProtonInelastic;
692 G4HadronElastic* elastic_lhep2 = new << 644 theInelasticProcess->RegisterMe(theLEInelasticModel);
693 elastic_lhep2->SetMaxEnergy( elastic << 645 G4HEAntiProtonInelastic* theHEInelasticModel =
694 G4HadronElasticProcess* theElasticPr << 646 new G4HEAntiProtonInelastic;
695 theElasticProcess->AddDataSet( elast << 647 theInelasticProcess->RegisterMe(theHEInelasticModel);
696 theElasticProcess->RegisterMe( elast << 648 pmanager->AddDiscreteProcess(theInelasticProcess);
697 theElasticProcess->RegisterMe( elast <<
698 pmanager->AddDiscreteProcess( theElasticPr <<
699 // Inelastic scattering <<
700 G4HadronInelasticProcess* theInelasticProc <<
701 new G4HadronInelasticProcess( "inelastic <<
702 theInelasticProcess->AddDataSet( theAntiNu <<
703 theInelasticProcess->RegisterMe( theFTFMod <<
704 pmanager->AddDiscreteProcess( theInelastic <<
705 // Absorption <<
706 pmanager->AddRestProcess(new G4HadronicAbs <<
707 } 649 }
>> 650
708 else if (particleName == "neutron") { 651 else if (particleName == "neutron") {
709 // elastic scattering 652 // elastic scattering
710 G4HadronElasticProcess* theElasticProcess = << 653 G4HadronElasticProcess* theNeutronElasticProcess =
711 theElasticProcess->AddDataSet(new G4Ne << 654 new G4HadronElasticProcess;
712 G4HadronElastic* elastic_neutronChipsM << 655 G4LElastic* theElasticModel1 = new G4LElastic;
713 elastic_neutronChipsModel->SetMinEnergy( 19. << 656 G4NeutronHPElastic * theElasticNeutron = new G4NeutronHPElastic;
714 theElasticProcess->RegisterMe( elastic << 657 theNeutronElasticProcess->RegisterMe(theElasticModel1);
715 G4ParticleHPElastic * theElasticNeutronHP = << 658 theElasticModel1->SetMinEnergy(19*MeV);
716 theElasticNeutronHP->SetMinEnergy( the << 659 theNeutronElasticProcess->RegisterMe(theElasticNeutron);
717 theElasticNeutronHP->SetMaxEnergy( the << 660 G4NeutronHPElasticData * theNeutronData = new G4NeutronHPElasticData;
718 theElasticProcess->RegisterMe( theElasticNeu << 661 theNeutronElasticProcess->AddDataSet(theNeutronData);
719 theElasticProcess->AddDataSet( new G4Particl << 662 pmanager->AddDiscreteProcess(theNeutronElasticProcess);
720 pmanager->AddDiscreteProcess( theElasticProc << 663 // inelastic scattering
721 // inelastic scattering << 664 G4NeutronInelasticProcess* theInelasticProcess =
722 G4HadronInelasticProcess* theInelasticProces << 665 new G4NeutronInelasticProcess("inelastic");
723 new G4HadronInelasticProcess( "inelastic", << 666 G4LENeutronInelastic* theInelasticModel = new G4LENeutronInelastic;
724 theInelasticProcess->AddDataSet( new G4Neutr << 667 theInelasticModel->SetMinEnergy(19*MeV);
725 theInelasticProcess->RegisterMe( theFTFModel << 668 theInelasticProcess->RegisterMe(theInelasticModel);
726 theInelasticProcess->RegisterMe( theBE << 669 G4NeutronHPInelastic * theLENeutronInelasticModel =
727 G4ParticleHPInelastic * theNeutronInelasticH << 670 new G4NeutronHPInelastic;
728 theNeutronInelasticHPModel->SetMinEner << 671 theInelasticProcess->RegisterMe(theLENeutronInelasticModel);
729 theNeutronInelasticHPModel->SetMaxEner << 672 G4NeutronHPInelasticData * theNeutronData1 =
730 theInelasticProcess->RegisterMe( theNeutronI << 673 new G4NeutronHPInelasticData;
731 theInelasticProcess->AddDataSet( new G4Parti << 674 theInelasticProcess->AddDataSet(theNeutronData1);
732 pmanager->AddDiscreteProcess(theInelasticPro 675 pmanager->AddDiscreteProcess(theInelasticProcess);
733 // capture 676 // capture
734 G4NeutronCaptureProcess* theCaptureProcess = << 677 G4HadronCaptureProcess* theCaptureProcess =
735 new G4NeutronCaptureProcess; << 678 new G4HadronCaptureProcess;
736 G4ParticleHPCapture * theLENeutronCaptureMod << 679 G4LCapture* theCaptureModel = new G4LCapture;
737 theLENeutronCaptureModel->SetMinEnergy(theHP << 680 theCaptureModel->SetMinEnergy(19*MeV);
738 theLENeutronCaptureModel->SetMaxEnergy(theHP << 681 theCaptureProcess->RegisterMe(theCaptureModel);
>> 682 G4NeutronHPCapture * theLENeutronCaptureModel = new G4NeutronHPCapture;
739 theCaptureProcess->RegisterMe(theLENeutronCa 683 theCaptureProcess->RegisterMe(theLENeutronCaptureModel);
740 theCaptureProcess->AddDataSet( new G4Particl << 684 G4NeutronHPCaptureData * theNeutronData3 = new G4NeutronHPCaptureData;
>> 685 theCaptureProcess->AddDataSet(theNeutronData3);
741 pmanager->AddDiscreteProcess(theCaptureProce 686 pmanager->AddDiscreteProcess(theCaptureProcess);
>> 687 // G4ProcessManager* pmanager = G4Neutron::Neutron->GetProcessManager();
>> 688 // pmanager->AddProcess(new G4UserSpecialCuts(),-1,-1,1);
742 } 689 }
743 else if (particleName == "anti_neutron") 690 else if (particleName == "anti_neutron")
744 { 691 {
745 // Elastic scattering << 692 pmanager->AddDiscreteProcess(theElasticProcess);
746 G4HadronElasticProcess* theElasticPr << 693 G4AntiNeutronInelasticProcess* theInelasticProcess =
747 theElasticProcess->AddDataSet( theGGHNEl ) << 694 new G4AntiNeutronInelasticProcess("inelastic");
748 theElasticProcess->RegisterMe( elast << 695 G4LEAntiNeutronInelastic* theLEInelasticModel =
749 pmanager->AddDiscreteProcess( theElasticPr << 696 new G4LEAntiNeutronInelastic;
750 // Inelastic scattering (include ann << 697 theInelasticProcess->RegisterMe(theLEInelasticModel);
751 G4HadronInelasticProcess* theInelasticProc << 698 G4HEAntiNeutronInelastic* theHEInelasticModel =
752 new G4HadronInelasticProcess( "inelastic << 699 new G4HEAntiNeutronInelastic;
753 theInelasticProcess->AddDataSet( theAntiNu << 700 theInelasticProcess->RegisterMe(theHEInelasticModel);
754 theInelasticProcess->RegisterMe( theFTFMod << 701 pmanager->AddDiscreteProcess(theInelasticProcess);
755 pmanager->AddDiscreteProcess( theInelastic <<
756 } 702 }
>> 703
757 else if (particleName == "deuteron") 704 else if (particleName == "deuteron")
758 { 705 {
759 // Elastic scattering << 706 pmanager->AddDiscreteProcess(theElasticProcess);
760 G4HadronElasticProcess* theElasticPr << 707 G4DeuteronInelasticProcess* theInelasticProcess =
761 theElasticProcess->AddDataSet( theGGNNEl ) << 708 new G4DeuteronInelasticProcess("inelastic");
762 theElasticProcess->RegisterMe( elast << 709 G4LEDeuteronInelastic* theLEInelasticModel =
763 pmanager->AddDiscreteProcess( theElasticPr << 710 new G4LEDeuteronInelastic;
764 // Inelastic scattering << 711 theInelasticProcess->RegisterMe(theLEInelasticModel);
765 G4HadronInelasticProcess* theInelasticProc << 712 pmanager->AddDiscreteProcess(theInelasticProcess);
766 new G4HadronInelasticProcess( "inelastic <<
767 theInelasticProcess->AddDataSet( theGGNucl <<
768 theInelasticProcess->RegisterMe( theFTFMod <<
769 theInelasticProcess->RegisterMe( the <<
770 pmanager->AddDiscreteProcess( theInelastic <<
771 } 713 }
>> 714
772 else if (particleName == "triton") 715 else if (particleName == "triton")
773 { 716 {
774 // Elastic scattering << 717 pmanager->AddDiscreteProcess(theElasticProcess);
775 G4HadronElasticProcess* theElasticPr << 718 G4TritonInelasticProcess* theInelasticProcess =
776 theElasticProcess->AddDataSet( theGGNNEl ) << 719 new G4TritonInelasticProcess("inelastic");
777 theElasticProcess->RegisterMe( elast << 720 G4LETritonInelastic* theLEInelasticModel =
778 pmanager->AddDiscreteProcess( theElasticPr << 721 new G4LETritonInelastic;
779 // Inelastic scattering << 722 theInelasticProcess->RegisterMe(theLEInelasticModel);
780 G4HadronInelasticProcess* theInelasticProc << 723 pmanager->AddDiscreteProcess(theInelasticProcess);
781 new G4HadronInelasticProcess( "inelastic <<
782 theInelasticProcess->AddDataSet( theGGNucl <<
783 theInelasticProcess->RegisterMe( theFTFMod <<
784 theInelasticProcess->RegisterMe( the <<
785 pmanager->AddDiscreteProcess( theInelastic <<
786 } 724 }
>> 725
787 else if (particleName == "alpha") 726 else if (particleName == "alpha")
788 { 727 {
789 // Elastic scattering << 728 pmanager->AddDiscreteProcess(theElasticProcess);
790 G4HadronElasticProcess* theElasticPr << 729 G4AlphaInelasticProcess* theInelasticProcess =
791 theElasticProcess->AddDataSet( theGGNNEl ) << 730 new G4AlphaInelasticProcess("inelastic");
792 theElasticProcess->RegisterMe( elast << 731 G4LEAlphaInelastic* theLEInelasticModel =
793 pmanager->AddDiscreteProcess( theElasticPr << 732 new G4LEAlphaInelastic;
794 // Inelastic scattering << 733 theInelasticProcess->RegisterMe(theLEInelasticModel);
795 G4HadronInelasticProcess* theInelasticProc << 734 pmanager->AddDiscreteProcess(theInelasticProcess);
796 new G4HadronInelasticProcess( "inelastic <<
797 theInelasticProcess->AddDataSet( the <<
798 theInelasticProcess->RegisterMe( theFTFMod <<
799 theInelasticProcess->RegisterMe( the <<
800 pmanager->AddDiscreteProcess( theInelastic <<
801 } 735 }
>> 736
802 } 737 }
803 } 738 }
804 739
>> 740
805 // Decays //////////////////////////////////// 741 // Decays ///////////////////////////////////////////////////////////////////
>> 742 #include "G4Decay.hh"
>> 743 #include "G4RadioactiveDecay.hh"
>> 744 #include "G4IonTable.hh"
>> 745 #include "G4Ions.hh"
>> 746
806 void DMXPhysicsList::ConstructGeneral() { 747 void DMXPhysicsList::ConstructGeneral() {
807 748
808 // Add Decay Process 749 // Add Decay Process
809 G4Decay* theDecayProcess = new G4Decay(); 750 G4Decay* theDecayProcess = new G4Decay();
810 auto particleIterator=GetParticleIterator(); << 751 theParticleIterator->reset();
811 particleIterator->reset(); << 752 while( (*theParticleIterator)() )
812 while( (*particleIterator)() ) <<
813 { 753 {
814 G4ParticleDefinition* particle = particl << 754 G4ParticleDefinition* particle = theParticleIterator->value();
815 G4ProcessManager* pmanager = particle->G 755 G4ProcessManager* pmanager = particle->GetProcessManager();
816 756
817 if (theDecayProcess->IsApplicable(*parti 757 if (theDecayProcess->IsApplicable(*particle) && !particle->IsShortLived())
818 { 758 {
819 pmanager ->AddProcess(theDecayProcess); 759 pmanager ->AddProcess(theDecayProcess);
820 // set ordering for PostStepDoIt and AtRes 760 // set ordering for PostStepDoIt and AtRestDoIt
821 pmanager ->SetProcessOrdering(theDecayProc 761 pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep);
822 pmanager ->SetProcessOrdering(theDecayProc 762 pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest);
823 } 763 }
824 } 764 }
825 765
826 // Declare radioactive decay to the GenericI 766 // Declare radioactive decay to the GenericIon in the IonTable.
827 G4LossTableManager* man = G4LossTableManager << 767 const G4IonTable *theIonTable =
828 G4VAtomDeexcitation* ad = man->AtomDeexcitat << 768 G4ParticleTable::GetParticleTable()->GetIonTable();
829 if(!ad) { << 769 G4RadioactiveDecay *theRadioactiveDecay = new G4RadioactiveDecay();
830 G4EmParameters::Instance()->SetAugerCascad <<
831 ad = new G4UAtomicDeexcitation(); <<
832 man->SetAtomDeexcitation(ad); <<
833 ad->InitialiseAtomicDeexcitation(); <<
834 } <<
835 770
836 G4PhysicsListHelper::GetPhysicsListHelper()- << 771 for (G4int i=0; i<theIonTable->Entries(); i++)
837 RegisterProcess(new G4RadioactiveDecay(), << 772 {
>> 773 G4String particleName = theIonTable->GetParticle(i)->GetParticleName();
>> 774 G4String particleType = theIonTable->GetParticle(i)->GetParticleType();
>> 775
>> 776 if (particleName == "GenericIon")
>> 777 {
>> 778 G4ProcessManager* pmanager =
>> 779 theIonTable->GetParticle(i)->GetProcessManager();
>> 780 pmanager->SetVerboseLevel(VerboseLevel);
>> 781 pmanager ->AddProcess(theRadioactiveDecay);
>> 782 pmanager ->SetProcessOrdering(theRadioactiveDecay, idxPostStep);
>> 783 pmanager ->SetProcessOrdering(theRadioactiveDecay, idxAtRest);
>> 784 }
>> 785 }
838 } 786 }
839 787
840 // Cuts ////////////////////////////////////// 788 // Cuts /////////////////////////////////////////////////////////////////////
841 void DMXPhysicsList::SetCuts() 789 void DMXPhysicsList::SetCuts()
842 { 790 {
843 791
844 if (verboseLevel >1) 792 if (verboseLevel >1)
845 G4cout << "DMXPhysicsList::SetCuts:"; 793 G4cout << "DMXPhysicsList::SetCuts:";
846 794
847 if (verboseLevel>0){ 795 if (verboseLevel>0){
848 G4cout << "DMXPhysicsList::SetCuts:"; 796 G4cout << "DMXPhysicsList::SetCuts:";
849 G4cout << "CutLength : " 797 G4cout << "CutLength : "
850 << G4BestUnit(defaultCutValue,"Length") < 798 << G4BestUnit(defaultCutValue,"Length") << G4endl;
851 } 799 }
852 800
853 //special for low energy physics 801 //special for low energy physics
854 G4double lowlimit=250*eV; 802 G4double lowlimit=250*eV;
855 G4ProductionCutsTable::GetProductionCutsTabl 803 G4ProductionCutsTable::GetProductionCutsTable()->SetEnergyRange(lowlimit,100.*GeV);
856 804
857 // set cut values for gamma at first and for 805 // set cut values for gamma at first and for e- second and next for e+,
858 // because some processes for e+/e- need cut 806 // because some processes for e+/e- need cut values for gamma
859 SetCutValue(cutForGamma, "gamma"); 807 SetCutValue(cutForGamma, "gamma");
860 SetCutValue(cutForElectron, "e-"); 808 SetCutValue(cutForElectron, "e-");
861 SetCutValue(cutForPositron, "e+"); 809 SetCutValue(cutForPositron, "e+");
862 810
863 if (verboseLevel>0) DumpCutValuesTable(); 811 if (verboseLevel>0) DumpCutValuesTable();
864 } 812 }
865 813
866 814