Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DMXPhysicsList.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 // --------------------------------------------------------------
28 // GEANT 4 - Underground Dark Matter Detector Advanced Example
29 //
30 // For information related to this code contact: Alex Howard
31 // e-mail: alexander.howard@cern.ch
32 // --------------------------------------------------------------
33 // Comments
34 //
35 // Underground Advanced
36 // by A. Howard and H. Araujo
37 // (27th November 2001)
38 //
39 // PhysicsList program
40 //
41 // Modified:
42 //
43 // 14-02-03 Fix bugs in msc and hIon instanciation + cut per region
44 //
45 // 05-02-05 AH - changes to G4Decay - added is not short lived protection
46 // and redefined particles to allow non-static creation
47 // i.e. changed construction to G4MesonConstructor, G4BaryonConstructor
48 //
49 // 23-10-09 LP - migrated EM physics from the LowEnergy processes (not supported) to
50 // the new G4Livermore model implementation. Results unchanged.
51 //
52 // --------------------------------------------------------------
53 
54 #include <iomanip>
55 
56 #include "DMXPhysicsList.hh"
57 
58 #include "globals.hh"
59 #include "G4SystemOfUnits.hh"
60 #include "G4ProcessManager.hh"
61 #include "G4ProcessVector.hh"
62 
63 #include "G4ParticleDefinition.hh"
64 #include "G4ParticleWithCuts.hh"
65 #include "G4ParticleTypes.hh"
66 #include "G4ParticleTable.hh"
67 
68 #include "G4ios.hh"
69 #include "G4UserLimits.hh"
70 
71 // Constructor /////////////////////////////////////////////////////////////
73 {
74 
75  defaultCutValue = 1.0*micrometer; //
76  cutForGamma = defaultCutValue;
77  cutForElectron = 1.0*nanometer;
78  cutForPositron = defaultCutValue;
79 
80  VerboseLevel = 1;
81  OpVerbLevel = 0;
82 
83  SetVerboseLevel(VerboseLevel);
84 }
85 
86 
87 // Destructor //////////////////////////////////////////////////////////////
89 {;}
90 
91 
92 // Construct Particles /////////////////////////////////////////////////////
94 {
95 
96  // In this method, static member functions should be called
97  // for all particles which you want to use.
98  // This ensures that objects of these particle types will be
99  // created in the program.
100 
101  ConstructMyBosons();
102  ConstructMyLeptons();
103  ConstructMyHadrons();
104  ConstructMyShortLiveds();
105 
106 }
107 
108 
109 // construct Bosons://///////////////////////////////////////////////////
110 void DMXPhysicsList::ConstructMyBosons()
111 {
112  // pseudo-particles
115 
116  // gamma
118 
119  //OpticalPhotons
121 
122 }
123 
124 
125 // construct Leptons://///////////////////////////////////////////////////
126 void DMXPhysicsList::ConstructMyLeptons()
127 {
128  // leptons
133 
138 }
139 
140 
141 #include "G4MesonConstructor.hh"
142 #include "G4BaryonConstructor.hh"
143 #include "G4IonConstructor.hh"
144 
145 // construct Hadrons://///////////////////////////////////////////////////
146 void DMXPhysicsList::ConstructMyHadrons()
147 {
148  // mesons
149  G4MesonConstructor mConstructor;
150  mConstructor.ConstructParticle();
151 
152  // baryons
153  G4BaryonConstructor bConstructor;
154  bConstructor.ConstructParticle();
155 
156  // ions
157  G4IonConstructor iConstructor;
158  iConstructor.ConstructParticle();
159 
160 }
161 
163 
164 // construct Shortliveds://///////////////////////////////////////////////////
165 void DMXPhysicsList::ConstructMyShortLiveds()
166 {
167  G4ShortLivedConstructor slConstructor;
168  slConstructor.ConstructParticle();
169 }
170 
171 
172 
173 
174 // Construct Processes //////////////////////////////////////////////////////
176 {
177 
179 
180  ConstructEM();
181 
182  ConstructOp();
183 
184  ConstructHad();
185 
187 
188 }
189 
190 
191 // Transportation ///////////////////////////////////////////////////////////
192 #include "DMXMaxTimeCuts.hh"
193 #include "DMXMinEkineCuts.hh"
194 #include "G4StepLimiter.hh"
195 
197 
199 
201  particleIterator->reset();
202  while( (*particleIterator)() ){
203  G4ParticleDefinition* particle = particleIterator->value();
204  G4ProcessManager* pmanager = particle->GetProcessManager();
205  G4String particleName = particle->GetParticleName();
206  // time cuts for ONLY neutrons:
207  if(particleName == "neutron")
208  pmanager->AddDiscreteProcess(new DMXMaxTimeCuts());
209  // Energy cuts to kill charged (embedded in method) particles:
210  pmanager->AddDiscreteProcess(new DMXMinEkineCuts());
211 
212  // Step limit applied to all particles:
213  pmanager->AddProcess(new G4StepLimiter, -1,-1,1);
214 
215  }
216 }
217 
218 
219 // Electromagnetic Processes ////////////////////////////////////////////////
220 // all charged particles
221 
222 // gamma
223 #include "G4PhotoElectricEffect.hh"
225 
226 #include "G4ComptonScattering.hh"
228 
229 #include "G4GammaConversion.hh"
231 
232 #include "G4RayleighScattering.hh"
234 
235 
236 // e-
237 #include "G4eMultipleScattering.hh"
238 
239 #include "G4eIonisation.hh"
241 
242 #include "G4eBremsstrahlung.hh"
244 
245 
246 // e+
247 #include "G4eIonisation.hh"
248 #include "G4eBremsstrahlung.hh"
249 #include "G4eplusAnnihilation.hh"
250 
251 
252 // alpha and GenericIon and deuterons, triton, He3:
253 
254 //muon:
255 #include "G4MuIonisation.hh"
256 #include "G4MuBremsstrahlung.hh"
257 #include "G4MuPairProduction.hh"
258 #include "G4MuonMinusCapture.hh"
259 
260 //OTHERS:
261 #include "G4hIonisation.hh"
262 #include "G4hMultipleScattering.hh"
263 #include "G4hBremsstrahlung.hh"
264 #include "G4ionIonisation.hh"
266 
267 //em process options to allow msc step-limitation to be switched off
268 #include "G4EmParameters.hh"
269 #include "G4VAtomDeexcitation.hh"
270 #include "G4UAtomicDeexcitation.hh"
271 #include "G4LossTableManager.hh"
272 
274 
275  //set a finer grid of the physic tables in order to improve precision
276  //former LowEnergy models have 200 bins up to 100 GeV
278  param->SetMaxEnergy(100*GeV);
279  param->SetNumberOfBinsPerDecade(20);
281  param->SetFluo(true);
282  param->SetPixe(true);
283  param->SetAuger(true);
286  if(!ad) {
288  }
289 
291  particleIterator->reset();
292  while( (*particleIterator)() ){
293  G4ParticleDefinition* particle = particleIterator->value();
294  G4ProcessManager* pmanager = particle->GetProcessManager();
295  G4String particleName = particle->GetParticleName();
296  G4String particleType = particle->GetParticleType();
297  G4double charge = particle->GetPDGCharge();
298 
299  if (particleName == "gamma")
300  {
301  //gamma
302  G4RayleighScattering* theRayleigh = new G4RayleighScattering();
303  pmanager->AddDiscreteProcess(theRayleigh);
304 
305  G4PhotoElectricEffect* thePhotoElectricEffect = new G4PhotoElectricEffect();
306  thePhotoElectricEffect->SetEmModel(new G4LivermorePhotoElectricModel());
307  pmanager->AddDiscreteProcess(thePhotoElectricEffect);
308 
309  G4ComptonScattering* theComptonScattering = new G4ComptonScattering();
310  theComptonScattering->SetEmModel(new G4LivermoreComptonModel());
311  pmanager->AddDiscreteProcess(theComptonScattering);
312 
313  G4GammaConversion* theGammaConversion = new G4GammaConversion();
314  theGammaConversion->SetEmModel(new G4LivermoreGammaConversionModel());
315  pmanager->AddDiscreteProcess(theGammaConversion);
316 
317  }
318  else if (particleName == "e-")
319  {
320  //electron
321  // process ordering: AddProcess(name, at rest, along step, post step)
322  // Multiple scattering
325  pmanager->AddProcess(msc,-1, 1, -1);
326 
327  // Ionisation
328  G4eIonisation* eIonisation = new G4eIonisation();
329  eIonisation->SetEmModel(new G4LivermoreIonisationModel());
330  eIonisation->SetStepFunction(0.2, 100*um); //improved precision in tracking
331  pmanager->AddProcess(eIonisation,-1, 2, 2);
332 
333  // Bremsstrahlung
334  G4eBremsstrahlung* eBremsstrahlung = new G4eBremsstrahlung();
335  eBremsstrahlung->SetEmModel(new G4LivermoreBremsstrahlungModel());
336  pmanager->AddProcess(eBremsstrahlung, -1,-3, 3);
337  }
338  else if (particleName == "e+")
339  {
340  //positron
343  pmanager->AddProcess(msc,-1, 1, 1);
344 
345  // Ionisation
346  G4eIonisation* eIonisation = new G4eIonisation();
347  eIonisation->SetStepFunction(0.2, 100*um); //
348  pmanager->AddProcess(eIonisation, -1, 2, 2);
349 
350  //Bremsstrahlung (use default, no low-energy available)
351  pmanager->AddProcess(new G4eBremsstrahlung(), -1,-1, 3);
352 
353  //Annihilation
354  pmanager->AddProcess(new G4eplusAnnihilation(),0,-1, 4);
355  }
356  else if( particleName == "mu+" ||
357  particleName == "mu-" )
358  {
359  //muon
360  pmanager->AddProcess(new G4eMultipleScattering, -1, 1, 1);
361  pmanager->AddProcess(new G4MuIonisation(), -1, 2, 2);
362  pmanager->AddProcess(new G4MuBremsstrahlung(), -1,-1, 3);
363  pmanager->AddProcess(new G4MuPairProduction(), -1,-1, 4);
364  if( particleName == "mu-" )
365  pmanager->AddProcess(new G4MuonMinusCapture(), 0,-1,-1);
366  }
367  else if (particleName == "proton" ||
368  particleName == "pi+" ||
369  particleName == "pi-")
370  {
371  //multiple scattering
372  pmanager->AddProcess(new G4hMultipleScattering, -1, 1, 1);
373 
374  //ionisation
375  G4hIonisation* hIonisation = new G4hIonisation();
376  hIonisation->SetStepFunction(0.2, 50*um);
377  pmanager->AddProcess(hIonisation, -1, 2, 2);
378 
379  //bremmstrahlung
380  pmanager->AddProcess(new G4hBremsstrahlung, -1,-3, 3);
381  }
382  else if(particleName == "alpha" ||
383  particleName == "deuteron" ||
384  particleName == "triton" ||
385  particleName == "He3")
386  {
387  //multiple scattering
388  pmanager->AddProcess(new G4hMultipleScattering,-1,1,1);
389 
390  //ionisation
391  G4ionIonisation* ionIoni = new G4ionIonisation();
392  ionIoni->SetStepFunction(0.1, 20*um);
393  pmanager->AddProcess(ionIoni, -1, 2, 2);
394  }
395  else if (particleName == "GenericIon")
396  {
397  // OBJECT may be dynamically created as either a GenericIon or nucleus
398  // G4Nucleus exists and therefore has particle type nucleus
399  // genericIon:
400 
401  //multiple scattering
402  pmanager->AddProcess(new G4hMultipleScattering,-1,1,1);
403 
404  //ionisation
405  G4ionIonisation* ionIoni = new G4ionIonisation();
406  ionIoni->SetEmModel(new G4IonParametrisedLossModel());
407  ionIoni->SetStepFunction(0.1, 20*um);
408  pmanager->AddProcess(ionIoni, -1, 2, 2);
409  }
410 
411  else if ((!particle->IsShortLived()) &&
412  (charge != 0.0) &&
413  (particle->GetParticleName() != "chargedgeantino"))
414  {
415  //all others charged particles except geantino
416  G4hMultipleScattering* aMultipleScattering = new G4hMultipleScattering();
417  G4hIonisation* ahadronIon = new G4hIonisation();
418 
419  //multiple scattering
420  pmanager->AddProcess(aMultipleScattering,-1,1,1);
421 
422  //ionisation
423  pmanager->AddProcess(ahadronIon, -1,2,2);
424  }
425 
426  }
427 }
428 
429 
430 // Optical Processes ////////////////////////////////////////////////////////
431 #include "G4Scintillation.hh"
432 #include "G4OpAbsorption.hh"
433 //#include "G4OpRayleigh.hh"
434 #include "G4OpBoundaryProcess.hh"
435 
437 {
438  // default scintillation process
439  G4Scintillation* theScintProcessDef = new G4Scintillation("Scintillation");
440  // theScintProcessDef->DumpPhysicsTable();
441  theScintProcessDef->SetTrackSecondariesFirst(true);
442  theScintProcessDef->SetScintillationYieldFactor(1.0); //
443  theScintProcessDef->SetScintillationExcitationRatio(0.0); //
444  theScintProcessDef->SetVerboseLevel(OpVerbLevel);
445 
446  // scintillation process for alpha:
447  G4Scintillation* theScintProcessAlpha = new G4Scintillation("Scintillation");
448  // theScintProcessNuc->DumpPhysicsTable();
449  theScintProcessAlpha->SetTrackSecondariesFirst(true);
450  theScintProcessAlpha->SetScintillationYieldFactor(1.1);
451  theScintProcessAlpha->SetScintillationExcitationRatio(1.0);
452  theScintProcessAlpha->SetVerboseLevel(OpVerbLevel);
453 
454  // scintillation process for heavy nuclei
455  G4Scintillation* theScintProcessNuc = new G4Scintillation("Scintillation");
456  // theScintProcessNuc->DumpPhysicsTable();
457  theScintProcessNuc->SetTrackSecondariesFirst(true);
458  theScintProcessNuc->SetScintillationYieldFactor(0.2);
459  theScintProcessNuc->SetScintillationExcitationRatio(1.0);
460  theScintProcessNuc->SetVerboseLevel(OpVerbLevel);
461 
462  // optical processes
463  G4OpAbsorption* theAbsorptionProcess = new G4OpAbsorption();
464  // G4OpRayleigh* theRayleighScatteringProcess = new G4OpRayleigh();
465  G4OpBoundaryProcess* theBoundaryProcess = new G4OpBoundaryProcess();
466  // theAbsorptionProcess->DumpPhysicsTable();
467  // theRayleighScatteringProcess->DumpPhysicsTable();
468  theAbsorptionProcess->SetVerboseLevel(OpVerbLevel);
469  // theRayleighScatteringProcess->SetVerboseLevel(OpVerbLevel);
470  theBoundaryProcess->SetVerboseLevel(OpVerbLevel);
471 
473  particleIterator->reset();
474  while( (*particleIterator)() )
475  {
476  G4ParticleDefinition* particle = particleIterator->value();
477  G4ProcessManager* pmanager = particle->GetProcessManager();
478  G4String particleName = particle->GetParticleName();
479  if (theScintProcessDef->IsApplicable(*particle)) {
480  // if(particle->GetPDGMass() > 5.0*GeV)
481  if(particle->GetParticleName() == "GenericIon") {
482  pmanager->AddProcess(theScintProcessNuc); // AtRestDiscrete
483  pmanager->SetProcessOrderingToLast(theScintProcessNuc,idxAtRest);
484  pmanager->SetProcessOrderingToLast(theScintProcessNuc,idxPostStep);
485  }
486  else if(particle->GetParticleName() == "alpha") {
487  pmanager->AddProcess(theScintProcessAlpha);
488  pmanager->SetProcessOrderingToLast(theScintProcessAlpha,idxAtRest);
489  pmanager->SetProcessOrderingToLast(theScintProcessAlpha,idxPostStep);
490  }
491  else {
492  pmanager->AddProcess(theScintProcessDef);
493  pmanager->SetProcessOrderingToLast(theScintProcessDef,idxAtRest);
494  pmanager->SetProcessOrderingToLast(theScintProcessDef,idxPostStep);
495  }
496  }
497 
498  if (particleName == "opticalphoton") {
499  pmanager->AddDiscreteProcess(theAbsorptionProcess);
500  // pmanager->AddDiscreteProcess(theRayleighScatteringProcess);
501  pmanager->AddDiscreteProcess(theBoundaryProcess);
502  }
503  }
504 }
505 
506 
507 // Hadronic processes ////////////////////////////////////////////////////////
508 
509 // Elastic processes:
510 #include "G4HadronElasticProcess.hh"
511 #include "G4ChipsElasticModel.hh"
512 #include "G4ElasticHadrNucleusHE.hh"
513 
514 // Inelastic processes:
528 
529 // High energy FTFP model and Bertini cascade
530 #include "G4FTFModel.hh"
532 #include "G4ExcitedStringDecay.hh"
533 #include "G4PreCompoundModel.hh"
535 #include "G4TheoFSGenerator.hh"
536 #include "G4CascadeInterface.hh"
537 
538 // Cross sections
539 #include "G4VCrossSectionDataSet.hh"
541 
542 #include "G4CrossSectionElastic.hh"
543 #include "G4BGGPionElasticXS.hh"
544 #include "G4AntiNuclElastic.hh"
545 
548 #include "G4CrossSectionPairGG.hh"
552 
553 #include "G4HadronElastic.hh"
554 #include "G4HadronCaptureProcess.hh"
555 
556 // Neutron high-precision models: <20 MeV
557 #include "G4ParticleHPElastic.hh"
559 #include "G4ParticleHPCapture.hh"
561 #include "G4ParticleHPInelastic.hh"
563 
564 // Stopping processes
568 
569 
570 
572 {
573  //Elastic models
574  const G4double elastic_elimitPi = 1.0*GeV;
575 
576  G4HadronElastic* elastic_lhep0 = new G4HadronElastic();
577  G4HadronElastic* elastic_lhep1 = new G4HadronElastic();
578  elastic_lhep1->SetMaxEnergy( elastic_elimitPi );
579  G4ChipsElasticModel* elastic_chip = new G4ChipsElasticModel();
580  G4ElasticHadrNucleusHE* elastic_he = new G4ElasticHadrNucleusHE();
581  elastic_he->SetMinEnergy( elastic_elimitPi );
582 
583 
584  // Inelastic scattering
585  const G4double theFTFMin0 = 0.0*GeV;
586  const G4double theFTFMin1 = 4.0*GeV;
587  const G4double theFTFMax = 100.0*TeV;
588  const G4double theBERTMin0 = 0.0*GeV;
589  const G4double theBERTMin1 = 19.0*MeV;
590  const G4double theBERTMax = 5.0*GeV;
591  const G4double theHPMin = 0.0*GeV;
592  const G4double theHPMax = 20.0*MeV;
593 
594  G4FTFModel * theStringModel = new G4FTFModel;
596  theStringModel->SetFragmentationModel( theStringDecay );
597  G4PreCompoundModel * thePreEquilib = new G4PreCompoundModel( new G4ExcitationHandler );
598  G4GeneratorPrecompoundInterface * theCascade = new G4GeneratorPrecompoundInterface( thePreEquilib );
599 
600  G4TheoFSGenerator * theFTFModel0 = new G4TheoFSGenerator( "FTFP" );
601  theFTFModel0->SetHighEnergyGenerator( theStringModel );
602  theFTFModel0->SetTransport( theCascade );
603  theFTFModel0->SetMinEnergy( theFTFMin0 );
604  theFTFModel0->SetMaxEnergy( theFTFMax );
605 
606  G4TheoFSGenerator * theFTFModel1 = new G4TheoFSGenerator( "FTFP" );
607  theFTFModel1->SetHighEnergyGenerator( theStringModel );
608  theFTFModel1->SetTransport( theCascade );
609  theFTFModel1->SetMinEnergy( theFTFMin1 );
610  theFTFModel1->SetMaxEnergy( theFTFMax );
611 
612  G4CascadeInterface * theBERTModel0 = new G4CascadeInterface;
613  theBERTModel0->SetMinEnergy( theBERTMin0 );
614  theBERTModel0->SetMaxEnergy( theBERTMax );
615 
616  G4CascadeInterface * theBERTModel1 = new G4CascadeInterface;
617  theBERTModel1->SetMinEnergy( theBERTMin1 );
618  theBERTModel1->SetMaxEnergy( theBERTMax );
619 
622  G4ComponentGGNuclNuclXsc * ggNuclNuclXsec = new G4ComponentGGNuclNuclXsc();
623  G4VCrossSectionDataSet * theGGNuclNuclData = new G4CrossSectionInelastic(ggNuclNuclXsec);
624 
626  particleIterator->reset();
627  while ((*particleIterator)())
628  {
629  G4ParticleDefinition* particle = particleIterator->value();
630  G4ProcessManager* pmanager = particle->GetProcessManager();
631  G4String particleName = particle->GetParticleName();
632 
633  if (particleName == "pi+")
634  {
635  // Elastic scattering
636  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
637  theElasticProcess->AddDataSet( new G4BGGPionElasticXS( particle ) );
638  theElasticProcess->RegisterMe( elastic_lhep1 );
639  theElasticProcess->RegisterMe( elastic_he );
640  pmanager->AddDiscreteProcess( theElasticProcess );
641  //Inelastic scattering
642  G4PionPlusInelasticProcess* theInelasticProcess =
643  new G4PionPlusInelasticProcess("inelastic");
644  theInelasticProcess->AddDataSet( thePiData );
645  theInelasticProcess->RegisterMe( theFTFModel1 );
646  theInelasticProcess->RegisterMe( theBERTModel0 );
647  pmanager->AddDiscreteProcess( theInelasticProcess );
648  }
649 
650  else if (particleName == "pi-")
651  {
652  // Elastic scattering
653  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
654  theElasticProcess->AddDataSet( new G4BGGPionElasticXS( particle ) );
655  theElasticProcess->RegisterMe( elastic_lhep1 );
656  theElasticProcess->RegisterMe( elastic_he );
657  pmanager->AddDiscreteProcess( theElasticProcess );
658  //Inelastic scattering
659  G4PionMinusInelasticProcess* theInelasticProcess =
660  new G4PionMinusInelasticProcess("inelastic");
661  theInelasticProcess->AddDataSet( thePiData );
662  theInelasticProcess->RegisterMe( theFTFModel1 );
663  theInelasticProcess->RegisterMe( theBERTModel0 );
664  pmanager->AddDiscreteProcess( theInelasticProcess );
665  //Absorption
667  }
668 
669  else if (particleName == "kaon+")
670  {
671  // Elastic scattering
672  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
673  theElasticProcess->RegisterMe( elastic_lhep0 );
674  pmanager->AddDiscreteProcess( theElasticProcess );
675  // Inelastic scattering
676  G4KaonPlusInelasticProcess* theInelasticProcess =
677  new G4KaonPlusInelasticProcess("inelastic");
678  theInelasticProcess->AddDataSet( G4CrossSectionDataSetRegistry::Instance()->
679  GetCrossSectionDataSet(G4ChipsKaonPlusInelasticXS::Default_Name()));
680  theInelasticProcess->RegisterMe( theFTFModel1 );
681  theInelasticProcess->RegisterMe( theBERTModel0 );
682  pmanager->AddDiscreteProcess( theInelasticProcess );
683  }
684 
685  else if (particleName == "kaon0S")
686  {
687  // Elastic scattering
688  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
689  theElasticProcess->RegisterMe( elastic_lhep0 );
690  pmanager->AddDiscreteProcess( theElasticProcess );
691  // Inelastic scattering
692  G4KaonZeroSInelasticProcess* theInelasticProcess =
693  new G4KaonZeroSInelasticProcess("inelastic");
694  theInelasticProcess->AddDataSet( G4CrossSectionDataSetRegistry::Instance()->
695  GetCrossSectionDataSet(G4ChipsKaonZeroInelasticXS::Default_Name()));
696  theInelasticProcess->RegisterMe( theFTFModel1 );
697  theInelasticProcess->RegisterMe( theBERTModel0 );
698  pmanager->AddDiscreteProcess( theInelasticProcess );
699  }
700 
701  else if (particleName == "kaon0L")
702  {
703  // Elastic scattering
704  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
705  theElasticProcess->RegisterMe( elastic_lhep0 );
706  pmanager->AddDiscreteProcess( theElasticProcess );
707  // Inelastic scattering
708  G4KaonZeroLInelasticProcess* theInelasticProcess =
709  new G4KaonZeroLInelasticProcess("inelastic");
710  theInelasticProcess->AddDataSet( G4CrossSectionDataSetRegistry::Instance()->
711  GetCrossSectionDataSet(G4ChipsKaonZeroInelasticXS::Default_Name()));
712  theInelasticProcess->RegisterMe( theFTFModel1 );
713  theInelasticProcess->RegisterMe( theBERTModel0 );
714  pmanager->AddDiscreteProcess( theInelasticProcess );
715  }
716 
717  else if (particleName == "kaon-")
718  {
719  // Elastic scattering
720  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
721  theElasticProcess->RegisterMe( elastic_lhep0 );
722  pmanager->AddDiscreteProcess( theElasticProcess );
723  // Inelastic scattering
724  G4KaonMinusInelasticProcess* theInelasticProcess =
725  new G4KaonMinusInelasticProcess("inelastic");
726  theInelasticProcess->AddDataSet( G4CrossSectionDataSetRegistry::Instance()->
727  GetCrossSectionDataSet(G4ChipsKaonMinusInelasticXS::Default_Name()));
728  theInelasticProcess->RegisterMe( theFTFModel1 );
729  theInelasticProcess->RegisterMe( theBERTModel0 );
730  pmanager->AddDiscreteProcess( theInelasticProcess );
732  }
733 
734  else if (particleName == "proton")
735  {
736  // Elastic scattering
737  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
739  GetCrossSectionDataSet(G4ChipsProtonElasticXS::Default_Name()));
740  theElasticProcess->RegisterMe( elastic_chip );
741  pmanager->AddDiscreteProcess( theElasticProcess );
742  // Inelastic scattering
743  G4ProtonInelasticProcess* theInelasticProcess =
744  new G4ProtonInelasticProcess("inelastic");
745  theInelasticProcess->AddDataSet( new G4BGGNucleonInelasticXS( G4Proton::Proton() ) );
746  theInelasticProcess->RegisterMe( theFTFModel1 );
747  theInelasticProcess->RegisterMe( theBERTModel0 );
748  pmanager->AddDiscreteProcess( theInelasticProcess );
749  }
750  else if (particleName == "anti_proton")
751  {
752  // Elastic scattering
753  const G4double elastic_elimitAntiNuc = 100.0*MeV;
754  G4AntiNuclElastic* elastic_anuc = new G4AntiNuclElastic();
755  elastic_anuc->SetMinEnergy( elastic_elimitAntiNuc );
756  G4CrossSectionElastic* elastic_anucxs = new G4CrossSectionElastic( elastic_anuc->GetComponentCrossSection() );
757  G4HadronElastic* elastic_lhep2 = new G4HadronElastic();
758  elastic_lhep2->SetMaxEnergy( elastic_elimitAntiNuc );
759  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
760  theElasticProcess->AddDataSet( elastic_anucxs );
761  theElasticProcess->RegisterMe( elastic_lhep2 );
762  theElasticProcess->RegisterMe( elastic_anuc );
763  pmanager->AddDiscreteProcess( theElasticProcess );
764  // Inelastic scattering
765  G4AntiProtonInelasticProcess* theInelasticProcess =
766  new G4AntiProtonInelasticProcess("inelastic");
767  theInelasticProcess->AddDataSet( theAntiNucleonData );
768  theInelasticProcess->RegisterMe( theFTFModel0 );
769  pmanager->AddDiscreteProcess( theInelasticProcess );
770  // Absorption
772  }
773 
774  else if (particleName == "neutron") {
775  // elastic scattering
776  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
777  theElasticProcess->AddDataSet(G4CrossSectionDataSetRegistry::Instance()->GetCrossSectionDataSet(G4ChipsNeutronElasticXS::Default_Name()));
778  G4HadronElastic* elastic_neutronChipsModel = new G4ChipsElasticModel();
779  elastic_neutronChipsModel->SetMinEnergy( 19.0*MeV );
780  theElasticProcess->RegisterMe( elastic_neutronChipsModel );
781  G4ParticleHPElastic * theElasticNeutronHP = new G4ParticleHPElastic;
782  theElasticNeutronHP->SetMinEnergy( theHPMin );
783  theElasticNeutronHP->SetMaxEnergy( theHPMax );
784  theElasticProcess->RegisterMe( theElasticNeutronHP );
785  theElasticProcess->AddDataSet( new G4ParticleHPElasticData );
786  pmanager->AddDiscreteProcess( theElasticProcess );
787  // inelastic scattering
788  G4NeutronInelasticProcess* theInelasticProcess =
789  new G4NeutronInelasticProcess("inelastic");
790  theInelasticProcess->AddDataSet( new G4BGGNucleonInelasticXS( G4Neutron::Neutron() ) );
791  theInelasticProcess->RegisterMe( theFTFModel1 );
792  theInelasticProcess->RegisterMe( theBERTModel1 );
793  G4ParticleHPInelastic * theNeutronInelasticHPModel = new G4ParticleHPInelastic;
794  theNeutronInelasticHPModel->SetMinEnergy( theHPMin );
795  theNeutronInelasticHPModel->SetMaxEnergy( theHPMax );
796  theInelasticProcess->RegisterMe( theNeutronInelasticHPModel );
797  theInelasticProcess->AddDataSet( new G4ParticleHPInelasticData );
798  pmanager->AddDiscreteProcess(theInelasticProcess);
799  // capture
800  G4HadronCaptureProcess* theCaptureProcess =
802  G4ParticleHPCapture * theLENeutronCaptureModel = new G4ParticleHPCapture;
803  theLENeutronCaptureModel->SetMinEnergy(theHPMin);
804  theLENeutronCaptureModel->SetMaxEnergy(theHPMax);
805  theCaptureProcess->RegisterMe(theLENeutronCaptureModel);
806  theCaptureProcess->AddDataSet( new G4ParticleHPCaptureData);
807  pmanager->AddDiscreteProcess(theCaptureProcess);
808 
809  }
810  else if (particleName == "anti_neutron")
811  {
812  // Elastic scattering
813  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
814  theElasticProcess->RegisterMe( elastic_lhep0 );
815  pmanager->AddDiscreteProcess( theElasticProcess );
816  // Inelastic scattering (include annihilation on-fly)
817  G4AntiNeutronInelasticProcess* theInelasticProcess =
818  new G4AntiNeutronInelasticProcess("inelastic");
819  theInelasticProcess->AddDataSet( theAntiNucleonData );
820  theInelasticProcess->RegisterMe( theFTFModel0 );
821  pmanager->AddDiscreteProcess( theInelasticProcess );
822  }
823 
824  else if (particleName == "deuteron")
825  {
826  // Elastic scattering
827  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
828  theElasticProcess->RegisterMe( elastic_lhep0 );
829  pmanager->AddDiscreteProcess( theElasticProcess );
830  // Inelastic scattering
831  G4DeuteronInelasticProcess* theInelasticProcess =
832  new G4DeuteronInelasticProcess("inelastic");
833  theInelasticProcess->AddDataSet( theGGNuclNuclData );
834  theInelasticProcess->RegisterMe( theFTFModel1 );
835  theInelasticProcess->RegisterMe( theBERTModel0 );
836  pmanager->AddDiscreteProcess( theInelasticProcess );
837  }
838 
839  else if (particleName == "triton")
840  {
841  // Elastic scattering
842  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
843  theElasticProcess->RegisterMe( elastic_lhep0 );
844  pmanager->AddDiscreteProcess( theElasticProcess );
845  // Inelastic scattering
846  G4TritonInelasticProcess* theInelasticProcess =
847  new G4TritonInelasticProcess("inelastic");
848  theInelasticProcess->AddDataSet( theGGNuclNuclData );
849  theInelasticProcess->RegisterMe( theFTFModel1 );
850  theInelasticProcess->RegisterMe( theBERTModel0 );
851  pmanager->AddDiscreteProcess( theInelasticProcess );
852  }
853  else if (particleName == "alpha")
854  {
855  // Elastic scattering
856  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
857  theElasticProcess->RegisterMe( elastic_lhep0 );
858  pmanager->AddDiscreteProcess( theElasticProcess );
859  // Inelastic scattering
860  G4AlphaInelasticProcess* theInelasticProcess =
861  new G4AlphaInelasticProcess("inelastic");
862  theInelasticProcess->AddDataSet( theGGNuclNuclData );
863  theInelasticProcess->RegisterMe( theFTFModel1 );
864  theInelasticProcess->RegisterMe( theBERTModel0 );
865  pmanager->AddDiscreteProcess( theInelasticProcess );
866  }
867 
868  }
869 }
870 
871 
872 // Decays ///////////////////////////////////////////////////////////////////
873 #include "G4Decay.hh"
874 #include "G4RadioactiveDecay.hh"
875 #include "G4IonTable.hh"
876 #include "G4Ions.hh"
877 
879 
880  // Add Decay Process
881  G4Decay* theDecayProcess = new G4Decay();
883  particleIterator->reset();
884  while( (*particleIterator)() )
885  {
886  G4ParticleDefinition* particle = particleIterator->value();
887  G4ProcessManager* pmanager = particle->GetProcessManager();
888 
889  if (theDecayProcess->IsApplicable(*particle) && !particle->IsShortLived())
890  {
891  pmanager ->AddProcess(theDecayProcess);
892  // set ordering for PostStepDoIt and AtRestDoIt
893  pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep);
894  pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest);
895  }
896  }
897 
898  // Declare radioactive decay to the GenericIon in the IonTable.
899  const G4IonTable *theIonTable =
901  G4RadioactiveDecay *theRadioactiveDecay = new G4RadioactiveDecay();
902 
903  for (G4int i=0; i<theIonTable->Entries(); i++)
904  {
905  G4String particleName = theIonTable->GetParticle(i)->GetParticleName();
906  G4String particleType = theIonTable->GetParticle(i)->GetParticleType();
907 
908  if (particleName == "GenericIon")
909  {
910  G4ProcessManager* pmanager =
911  theIonTable->GetParticle(i)->GetProcessManager();
912  pmanager->SetVerboseLevel(VerboseLevel);
913  pmanager ->AddProcess(theRadioactiveDecay);
914  pmanager ->SetProcessOrdering(theRadioactiveDecay, idxPostStep);
915  pmanager ->SetProcessOrdering(theRadioactiveDecay, idxAtRest);
916  }
917  }
918 }
919 
920 // Cuts /////////////////////////////////////////////////////////////////////
922 {
923 
924  if (verboseLevel >1)
925  G4cout << "DMXPhysicsList::SetCuts:";
926 
927  if (verboseLevel>0){
928  G4cout << "DMXPhysicsList::SetCuts:";
929  G4cout << "CutLength : "
930  << G4BestUnit(defaultCutValue,"Length") << G4endl;
931  }
932 
933  //special for low energy physics
934  G4double lowlimit=250*eV;
936 
937  // set cut values for gamma at first and for e- second and next for e+,
938  // because some processes for e+/e- need cut values for gamma
939  SetCutValue(cutForGamma, "gamma");
940  SetCutValue(cutForElectron, "e-");
941  SetCutValue(cutForPositron, "e+");
942 
944 }
945 
static G4Geantino * GeantinoDefinition()
Definition: G4Geantino.cc:82
static G4Electron * ElectronDefinition()
Definition: G4Electron.cc:89
static G4MuonPlus * MuonPlusDefinition()
Definition: G4MuonPlus.cc:94
static G4LossTableManager * Instance()
G4int Entries() const
Definition: G4IonTable.cc:1689
void SetCutValue(G4double aCut, const G4String &pname)
void SetEnergyRange(G4double lowedge, G4double highedge)
void SetMscStepLimitType(G4MscStepLimitType val)
virtual void ConstructGeneral()
virtual G4bool IsApplicable(const G4ParticleDefinition &)
Definition: G4Decay.cc:89
void SetFragmentationModel(G4VStringFragmentation *aModel)
void SetAuger(G4bool val)
G4int AddDiscreteProcess(G4VProcess *aProcess, G4int ord=ordDefault)
void SetVerboseLevel(G4int value)
static constexpr double nanometer
Definition: G4SIunits.hh:101
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
static void ConstructParticle()
static void ConstructParticle()
int G4int
Definition: G4Types.hh:78
void SetMaxEnergy(G4double val)
void SetStepFunction(G4double v1, G4double v2, G4bool lock=true)
G4ComponentAntiNuclNuclearXS * GetComponentCrossSection()
const G4String & GetParticleName() const
virtual void ConstructProcess()
void SetHighEnergyGenerator(G4VHighEnergyGenerator *const value)
static void ConstructParticle()
void RegisterMe(G4HadronicInteraction *a)
static constexpr double TeV
Definition: G4SIunits.hh:218
virtual void AddTransportation()
void SetScintillationYieldFactor(const G4double yieldfactor)
void SetMinEnergy(G4double anEnergy)
void SetEmModel(G4VEmModel *, G4int index=1)
void AddDataSet(G4VCrossSectionDataSet *aDataSet)
G4IonTable * GetIonTable() const
void DumpCutValuesTable(G4int flag=1)
G4GLOB_DLL std::ostream G4cout
static constexpr double um
Definition: G4SIunits.hh:113
void SetVerboseLevel(G4int value)
G4int AddProcess(G4VProcess *aProcess, G4int ordAtRestDoIt=ordInActive, G4int ordAlongSteptDoIt=ordInActive, G4int ordPostStepDoIt=ordInActive)
static G4NeutrinoE * NeutrinoEDefinition()
Definition: G4NeutrinoE.cc:80
static const char * Default_Name()
static G4AntiNeutrinoMu * AntiNeutrinoMuDefinition()
static G4CrossSectionDataSetRegistry * Instance()
void SetNumberOfBinsPerDecade(G4int val)
static G4Proton * Proton()
Definition: G4Proton.cc:93
static constexpr double eV
Definition: G4SIunits.hh:215
const G4String & GetParticleType() const
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
void SetProcessOrdering(G4VProcess *aProcess, G4ProcessVectorDoItIndex idDoIt, G4int ordDoIt=ordDefault)
virtual void ConstructHad()
void SetScintillationExcitationRatio(const G4double ratio)
static G4Positron * PositronDefinition()
Definition: G4Positron.cc:89
static G4ProductionCutsTable * GetProductionCutsTable()
void SetPixe(G4bool val)
virtual void SetCuts()
static G4ParticleTable * GetParticleTable()
void SetTrackSecondariesFirst(const G4bool state)
void SetProcessOrderingToLast(G4VProcess *aProcess, G4ProcessVectorDoItIndex idDoIt)
G4ProcessManager * GetProcessManager() const
static G4EmParameters * Instance()
void SetEmModel(G4VEmModel *, G4int index=1)
static G4MuonMinus * MuonMinusDefinition()
Definition: G4MuonMinus.cc:95
static constexpr double GeV
Definition: G4SIunits.hh:217
virtual void ConstructOp()
void SetMaxEnergy(const G4double anEnergy)
std::vector< G4InuclElementaryParticle >::iterator particleIterator
Definition: G4BigBanger.cc:65
static G4ChargedGeantino * ChargedGeantinoDefinition()
G4ParticleTable::G4PTblDicIterator * GetParticleIterator() const
static const char * Default_Name()
#define G4endl
Definition: G4ios.hh:61
void SetTransport(G4VIntraNuclearTransportModel *const value)
static constexpr double MeV
Definition: G4SIunits.hh:214
G4int AddRestProcess(G4VProcess *aProcess, G4int ord=ordDefault)
G4VAtomDeexcitation * AtomDeexcitation()
static G4AntiNeutrinoE * AntiNeutrinoEDefinition()
static G4OpticalPhoton * OpticalPhotonDefinition()
double G4double
Definition: G4Types.hh:76
G4bool IsApplicable(const G4ParticleDefinition &aParticleType) override
G4double GetPDGCharge() const
static G4NeutrinoMu * NeutrinoMuDefinition()
Definition: G4NeutrinoMu.cc:80
static constexpr double micrometer
Definition: G4SIunits.hh:100
void SetAtomDeexcitation(G4VAtomDeexcitation *)
virtual void ConstructParticle()
void SetVerboseLevel(G4int value)
Definition: G4VProcess.hh:437
virtual void ConstructEM()
void SetFluo(G4bool val)
void SetStepLimitType(G4MscStepLimitType val)
static G4Gamma * GammaDefinition()
Definition: G4Gamma.cc:81
G4ParticleDefinition * GetParticle(G4int index) const
Definition: G4IonTable.cc:1643