Back to home page

EIC code displayed by LXR

 
 

    


Warning, file /geant4/examples/advanced/underground_physics/src/DMXPhysicsList.cc was not indexed or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 //
0027 // --------------------------------------------------------------
0028 //   GEANT 4 - Underground Dark Matter Detector Advanced Example
0029 //
0030 //      For information related to this code contact: Alex Howard
0031 //      e-mail: alexander.howard@cern.ch
0032 // --------------------------------------------------------------
0033 // Comments
0034 //
0035 //                  Underground Advanced
0036 //               by A. Howard and H. Araujo 
0037 //                    (27th November 2001)
0038 //
0039 // PhysicsList program
0040 //
0041 // Modified:
0042 //
0043 // 14-02-03 Fix bugs in msc and hIon instanciation + cut per region
0044 //
0045 // 05-02-05 AH - changes to G4Decay - added is not short lived protection
0046 //          and redefined particles to allow non-static creation
0047 //          i.e. changed construction to G4MesonConstructor, G4BaryonConstructor
0048 // 
0049 // 23-10-09 LP - migrated EM physics from the LowEnergy processes (not supported) to 
0050 //          the new G4Livermore model implementation. Results unchanged.
0051 //
0052 // --------------------------------------------------------------
0053 
0054 #include <iomanip>  
0055 
0056 #include "DMXPhysicsList.hh"
0057 
0058 #include "globals.hh"
0059 #include "G4SystemOfUnits.hh"
0060 #include "G4ProcessManager.hh"
0061 #include "G4ProcessVector.hh"
0062 
0063 #include "G4ParticleDefinition.hh"
0064 #include "G4ParticleWithCuts.hh"
0065 #include "G4ParticleTypes.hh"
0066 #include "G4ParticleTable.hh"
0067 
0068 #include "G4ios.hh"
0069 #include "G4UserLimits.hh"
0070 
0071 #include "G4MesonConstructor.hh"
0072 #include "G4BaryonConstructor.hh"
0073 #include "G4IonConstructor.hh"
0074 #include "G4ShortLivedConstructor.hh"
0075 
0076 #include "DMXMaxTimeCuts.hh"
0077 #include "DMXMinEkineCuts.hh"
0078 #include "G4StepLimiter.hh"
0079 
0080 // gamma
0081 #include "G4PhotoElectricEffect.hh"
0082 #include "G4LivermorePhotoElectricModel.hh"
0083 
0084 #include "G4ComptonScattering.hh"
0085 #include "G4LivermoreComptonModel.hh"
0086 
0087 #include "G4GammaConversion.hh"
0088 #include "G4BetheHeitler5DModel.hh"
0089 
0090 #include "G4RayleighScattering.hh" 
0091 #include "G4LivermoreRayleighModel.hh"
0092 
0093 // e-
0094 #include "G4eMultipleScattering.hh"
0095 
0096 #include "G4eIonisation.hh"
0097 #include "G4LivermoreIonisationModel.hh"
0098 
0099 #include "G4eBremsstrahlung.hh"
0100 #include "G4UniversalFluctuation.hh"
0101 
0102 // e+
0103 #include "G4eIonisation.hh" 
0104 #include "G4eBremsstrahlung.hh" 
0105 #include "G4eplusAnnihilation.hh"
0106 
0107 // alpha and GenericIon and deuterons, triton, He3:
0108 //muon:
0109 #include "G4MuIonisation.hh"
0110 #include "G4MuBremsstrahlung.hh"
0111 #include "G4MuPairProduction.hh"
0112 #include "G4MuMultipleScattering.hh"
0113 #include "G4MuonMinusCapture.hh"
0114 
0115 //OTHERS:
0116 #include "G4hIonisation.hh" 
0117 #include "G4hMultipleScattering.hh"
0118 #include "G4hBremsstrahlung.hh"
0119 #include "G4ionIonisation.hh"
0120 #include "G4IonParametrisedLossModel.hh"
0121 #include "G4NuclearStopping.hh"
0122 
0123 //em process options to allow msc step-limitation to be switched off
0124 #include "G4EmParameters.hh"
0125 #include "G4VAtomDeexcitation.hh"
0126 #include "G4UAtomicDeexcitation.hh"
0127 #include "G4LossTableManager.hh"
0128 
0129 #include "G4Scintillation.hh"
0130 #include "G4OpAbsorption.hh"
0131 #include "G4OpBoundaryProcess.hh"
0132 #include "G4OpticalParameters.hh"
0133 
0134 // Elastic processes:
0135 #include "G4HadronElasticProcess.hh"
0136 #include "G4ChipsElasticModel.hh"
0137 #include "G4ElasticHadrNucleusHE.hh"
0138 
0139 // Inelastic processes:
0140 #include "G4HadronInelasticProcess.hh"
0141 
0142 // High energy FTFP model and Bertini cascade
0143 #include "G4FTFModel.hh"
0144 #include "G4LundStringFragmentation.hh"
0145 #include "G4ExcitedStringDecay.hh"
0146 #include "G4PreCompoundModel.hh"
0147 #include "G4GeneratorPrecompoundInterface.hh"
0148 #include "G4TheoFSGenerator.hh"
0149 #include "G4CascadeInterface.hh"
0150 
0151 // Cross sections
0152 #include "G4VCrossSectionDataSet.hh"
0153 #include "G4CrossSectionDataSetRegistry.hh"
0154 
0155 #include "G4CrossSectionElastic.hh"
0156 #include "G4CrossSectionInelastic.hh"
0157 #include "G4BGGPionElasticXS.hh"
0158 #include "G4BGGPionInelasticXS.hh"
0159 #include "G4AntiNuclElastic.hh"
0160 
0161 #include "G4CrossSectionInelastic.hh"
0162 #include "G4BGGNucleonInelasticXS.hh"
0163 #include "G4BGGNucleonElasticXS.hh"
0164 #include "G4NeutronInelasticXS.hh"
0165 #include "G4NeutronElasticXS.hh"
0166 #include "G4ComponentAntiNuclNuclearXS.hh"
0167 #include "G4ComponentGGNuclNuclXsc.hh"
0168 #include "G4ComponentGGHadronNucleusXsc.hh"
0169 
0170 #include "G4HadronElastic.hh"
0171 #include "G4NeutronCaptureProcess.hh"
0172 
0173 // Neutron high-precision models: <20 MeV
0174 #include "G4ParticleHPElastic.hh"
0175 #include "G4ParticleHPElasticData.hh"
0176 #include "G4ParticleHPCapture.hh"
0177 #include "G4ParticleHPCaptureData.hh"
0178 #include "G4ParticleHPInelastic.hh"
0179 #include "G4ParticleHPInelasticData.hh"
0180 
0181 // Stopping processes
0182 #include "G4HadronStoppingProcess.hh"
0183 #include "G4HadronicAbsorptionBertini.hh"
0184 #include "G4HadronicAbsorptionFritiof.hh"
0185 
0186 #include "G4HadronicParameters.hh"
0187 
0188 #include "G4Decay.hh"
0189 #include "G4RadioactiveDecay.hh"
0190 #include "G4PhysicsListHelper.hh"
0191 #include "G4NuclideTable.hh"
0192 #include "G4NuclearLevelData.hh"
0193 
0194 // Constructor /////////////////////////////////////////////////////////////
0195 DMXPhysicsList::DMXPhysicsList() : G4VUserPhysicsList() 
0196 {
0197 
0198   defaultCutValue     = 1.0*micrometer; //
0199   cutForGamma         = defaultCutValue;
0200   cutForElectron      = 1.0*nanometer;
0201   cutForPositron      = defaultCutValue;
0202 
0203   VerboseLevel = 1;
0204   OpVerbLevel = 0;
0205 
0206   //set a finer grid of the physic tables in order to improve precision
0207   //former LowEnergy models have 200 bins up to 100 GeV
0208   G4EmParameters* param = G4EmParameters::Instance();
0209   param->SetMaxEnergy(100*GeV);
0210   param->SetNumberOfBinsPerDecade(20);
0211   param->SetMscStepLimitType(fMinimal);
0212   param->SetFluo(true);
0213   param->SetPixe(true);
0214   param->SetAuger(true);
0215 
0216   G4EmParameters::Instance()->AddPhysics("World","G4RadioactiveDecay");
0217   G4DeexPrecoParameters* deex = G4NuclearLevelData::GetInstance()->GetParameters();
0218   deex->SetStoreICLevelData(true);
0219   deex->SetMaxLifeTime(G4NuclideTable::GetInstance()->GetThresholdOfHalfLife()
0220                        /std::log(2.));
0221   SetVerboseLevel(VerboseLevel);
0222 }
0223 
0224 // Destructor //////////////////////////////////////////////////////////////
0225 DMXPhysicsList::~DMXPhysicsList() 
0226 {;}
0227 
0228 // Construct Particles /////////////////////////////////////////////////////
0229 void DMXPhysicsList::ConstructParticle() 
0230 {
0231   // In this method, static member functions should be called
0232   // for all particles which you want to use.
0233   // This ensures that objects of these particle types will be
0234   // created in the program. 
0235 
0236   ConstructMyBosons();
0237   ConstructMyLeptons();
0238   ConstructMyHadrons();
0239   ConstructMyShortLiveds();
0240 }
0241 
0242 // construct Bosons://///////////////////////////////////////////////////
0243 void DMXPhysicsList::ConstructMyBosons()
0244 {
0245   // pseudo-particles
0246   G4Geantino::GeantinoDefinition();
0247   G4ChargedGeantino::ChargedGeantinoDefinition();
0248   
0249   // gamma
0250   G4Gamma::GammaDefinition();
0251 
0252   //OpticalPhotons
0253   G4OpticalPhoton::OpticalPhotonDefinition();
0254 
0255 }
0256 
0257 // construct Leptons://///////////////////////////////////////////////////
0258 void DMXPhysicsList::ConstructMyLeptons()
0259 {
0260   // leptons
0261   G4Electron::ElectronDefinition();
0262   G4Positron::PositronDefinition();
0263   G4MuonPlus::MuonPlusDefinition();
0264   G4MuonMinus::MuonMinusDefinition();
0265 
0266   G4NeutrinoE::NeutrinoEDefinition();
0267   G4AntiNeutrinoE::AntiNeutrinoEDefinition();
0268   G4NeutrinoMu::NeutrinoMuDefinition();
0269   G4AntiNeutrinoMu::AntiNeutrinoMuDefinition();
0270 }
0271 
0272 // construct Hadrons://///////////////////////////////////////////////////
0273 void DMXPhysicsList::ConstructMyHadrons()
0274 {
0275  //  mesons
0276   G4MesonConstructor mConstructor;
0277   mConstructor.ConstructParticle();
0278 
0279  //  baryons
0280   G4BaryonConstructor bConstructor;
0281   bConstructor.ConstructParticle();
0282 
0283  //  ions
0284   G4IonConstructor iConstructor;
0285   iConstructor.ConstructParticle();
0286 }
0287 
0288 // construct Shortliveds://///////////////////////////////////////////////////
0289 void DMXPhysicsList::ConstructMyShortLiveds()
0290 {
0291   G4ShortLivedConstructor slConstructor;
0292   slConstructor.ConstructParticle();
0293 }
0294 
0295 // Construct Processes //////////////////////////////////////////////////////
0296 void DMXPhysicsList::ConstructProcess() 
0297 {
0298   AddTransportation();
0299 
0300   ConstructEM();
0301 
0302   ConstructOp();
0303 
0304   ConstructHad();
0305 
0306   ConstructGeneral();
0307 }
0308 
0309 // Transportation ///////////////////////////////////////////////////////////
0310 void DMXPhysicsList::AddTransportation() {
0311   
0312   G4VUserPhysicsList::AddTransportation();
0313   
0314   auto particleIterator=GetParticleIterator();
0315   particleIterator->reset();
0316   while( (*particleIterator)() ){
0317     G4ParticleDefinition* particle = particleIterator->value();
0318     G4ProcessManager* pmanager = particle->GetProcessManager();
0319     G4String particleName = particle->GetParticleName();
0320     // time cuts for ONLY neutrons:
0321     if(particleName == "neutron") 
0322       pmanager->AddDiscreteProcess(new DMXMaxTimeCuts());
0323     // Energy cuts to kill charged (embedded in method) particles:
0324     pmanager->AddDiscreteProcess(new DMXMinEkineCuts());
0325 
0326     // Step limit applied to all particles:
0327     pmanager->AddDiscreteProcess(new G4StepLimiter);
0328   }           
0329 }
0330 
0331 // Electromagnetic Processes ////////////////////////////////////////////////
0332 // all charged particles
0333 void DMXPhysicsList::ConstructEM() {
0334   
0335   G4LossTableManager* man = G4LossTableManager::Instance();
0336   man->SetAtomDeexcitation(new G4UAtomicDeexcitation());
0337 
0338   G4EmParameters* em_params = G4EmParameters::Instance();
0339   
0340   auto particleIterator=GetParticleIterator();
0341   particleIterator->reset();
0342   while( (*particleIterator)() ){
0343     G4ParticleDefinition* particle = particleIterator->value();
0344     G4ProcessManager* pmanager = particle->GetProcessManager();
0345     G4String particleName = particle->GetParticleName();
0346     G4String particleType = particle->GetParticleType();
0347     G4double charge = particle->GetPDGCharge();
0348     
0349     if (particleName == "gamma") 
0350       {
0351     //gamma
0352     G4RayleighScattering* theRayleigh = new G4RayleighScattering();
0353     pmanager->AddDiscreteProcess(theRayleigh);
0354 
0355     G4PhotoElectricEffect* thePhotoElectricEffect = new G4PhotoElectricEffect();
0356     thePhotoElectricEffect->SetEmModel(new G4LivermorePhotoElectricModel());
0357     pmanager->AddDiscreteProcess(thePhotoElectricEffect);
0358     
0359     G4ComptonScattering* theComptonScattering = new G4ComptonScattering();
0360     theComptonScattering->SetEmModel(new G4LivermoreComptonModel());
0361     pmanager->AddDiscreteProcess(theComptonScattering);
0362     
0363     G4GammaConversion* theGammaConversion = new G4GammaConversion();
0364     theGammaConversion->SetEmModel(new G4BetheHeitler5DModel());
0365     pmanager->AddDiscreteProcess(theGammaConversion);
0366 
0367       } 
0368     else if (particleName == "e-") 
0369       {
0370     //electron
0371     // process ordering: AddProcess(name, at rest, along step, post step)
0372     // Multiple scattering
0373     G4eMultipleScattering* msc = new G4eMultipleScattering();
0374     em_params->SetMscStepLimitType(fUseDistanceToBoundary);
0375     pmanager->AddProcess(msc,-1, 1, -1);
0376 
0377     // Ionisation
0378     G4eIonisation* eIonisation = new G4eIonisation();
0379         G4VEmModel* theIoniLiv = new G4LivermoreIonisationModel();
0380         theIoniLiv->SetHighEnergyLimit(0.1*MeV); 
0381         eIonisation->AddEmModel(0, theIoniLiv, new G4UniversalFluctuation() );
0382     em_params->SetStepFunction(0.2, 100*um); //improved precision in tracking  
0383     pmanager->AddProcess(eIonisation,-1, 2, 1);
0384     
0385     // Bremsstrahlung
0386     G4eBremsstrahlung* eBremsstrahlung = new G4eBremsstrahlung();
0387     pmanager->AddProcess(eBremsstrahlung, -1,-3, 2);
0388       } 
0389     else if (particleName == "e+") 
0390       {
0391     //positron  
0392     G4eMultipleScattering* msc = new G4eMultipleScattering();
0393     msc->SetStepLimitType(fUseDistanceToBoundary);
0394     pmanager->AddProcess(msc,-1, 1, -1);
0395     
0396     // Ionisation
0397     G4eIonisation* eIonisation = new G4eIonisation();
0398     // eIonisation->SetStepFunction(0.2, 100*um); //     
0399     pmanager->AddProcess(eIonisation,                 -1, 2, 1);
0400 
0401     //Bremsstrahlung (use default, no low-energy available)
0402     pmanager->AddProcess(new G4eBremsstrahlung(), -1,-1, 2);
0403 
0404     //Annihilation
0405     pmanager->AddProcess(new G4eplusAnnihilation(),0,-1, 3);      
0406       } 
0407     else if( particleName == "mu+" || 
0408          particleName == "mu-"    ) 
0409       {
0410     //muon  
0411     pmanager->AddProcess(new G4MuMultipleScattering,    -1, 1,-1);
0412     pmanager->AddProcess(new G4MuIonisation(),          -1, 2, 1);
0413     pmanager->AddProcess(new G4MuBremsstrahlung(),      -1,-1, 2);
0414     pmanager->AddProcess(new G4MuPairProduction(),      -1,-1, 3);
0415     if( particleName == "mu-" )
0416       pmanager->AddProcess(new G4MuonMinusCapture(), 0,-1,-1);
0417       } 
0418     else if (particleName == "proton" || 
0419          particleName == "pi+" || 
0420          particleName == "pi-")
0421       {
0422     //multiple scattering
0423     pmanager->AddProcess(new G4hMultipleScattering, -1, 1, -1);
0424       
0425     //ionisation
0426     G4hIonisation* hIonisation = new G4hIonisation();
0427     em_params->SetStepFunctionMuHad(0.2, 50*um);
0428     pmanager->AddProcess(hIonisation,               -1, 2, 1);      
0429     
0430     //bremmstrahlung
0431     pmanager->AddProcess(new G4hBremsstrahlung,     -1,-3, 2);
0432       }
0433     else if(particleName == "alpha"      ||
0434          particleName == "deuteron"   ||
0435          particleName == "triton"     ||
0436          particleName == "He3")
0437       {
0438     //multiple scattering
0439     pmanager->AddProcess(new G4hMultipleScattering,-1,1,-1);
0440     
0441     //ionisation
0442     G4ionIonisation* ionIoni = new G4ionIonisation();
0443     em_params->SetStepFunctionLightIons(0.1, 1*CLHEP::um);
0444     pmanager->AddProcess(ionIoni,                   -1, 2, 1);
0445     pmanager->AddProcess(new G4NuclearStopping(),   -1, 3,-1);  
0446       }
0447     else if (particleName == "GenericIon")
0448       {
0449     // OBJECT may be dynamically created as either a GenericIon or nucleus
0450     // G4Nucleus exists and therefore has particle type nucleus
0451     // genericIon:
0452     
0453     //multiple scattering
0454     pmanager->AddProcess(new G4hMultipleScattering,-1,1,-1);
0455 
0456     //ionisation
0457     G4ionIonisation* ionIoni = new G4ionIonisation();
0458     em_params->SetStepFunctionIons(0.1, 1*CLHEP::um);
0459     pmanager->AddProcess(ionIoni,                   -1, 2, 1);
0460     pmanager->AddProcess(new G4NuclearStopping(),   -1, 3,-1);  
0461       } 
0462 
0463     else if ((!particle->IsShortLived()) &&
0464          (charge != 0.0) && 
0465          (particle->GetParticleName() != "chargedgeantino")) 
0466       {
0467     //all others charged particles except geantino
0468         G4hMultipleScattering* aMultipleScattering = new G4hMultipleScattering();
0469         G4hIonisation* ahadronIon = new G4hIonisation();
0470     
0471     //multiple scattering
0472     pmanager->AddProcess(aMultipleScattering,-1,1,-1);
0473 
0474     //ionisation
0475     pmanager->AddProcess(ahadronIon,       -1,2,1);      
0476       }
0477   }
0478 }
0479 
0480 // Optical Processes ////////////////////////////////////////////////////////
0481 void DMXPhysicsList::ConstructOp() 
0482 {
0483   G4OpticalParameters* opParams = G4OpticalParameters::Instance();
0484   G4Scintillation* theScintProcessDef = new G4Scintillation("Scintillation");
0485   opParams->SetScintTrackSecondariesFirst(true);
0486   opParams->SetScintByParticleType(true);
0487 
0488   // optical processes
0489   G4OpAbsorption* theAbsorptionProcess = new G4OpAbsorption();
0490   G4OpBoundaryProcess* theBoundaryProcess = new G4OpBoundaryProcess();
0491 
0492   auto particleIterator=GetParticleIterator();
0493   particleIterator->reset();
0494   while( (*particleIterator)() )
0495     {
0496       G4ParticleDefinition* particle = particleIterator->value();
0497       G4ProcessManager* pmanager = particle->GetProcessManager();
0498       G4String particleName = particle->GetParticleName();
0499       if (theScintProcessDef->IsApplicable(*particle)) {
0500         pmanager->AddProcess(theScintProcessDef);
0501         pmanager->SetProcessOrderingToLast(theScintProcessDef,idxAtRest);
0502         pmanager->SetProcessOrderingToLast(theScintProcessDef,idxPostStep);
0503       }
0504       
0505       if (particleName == "opticalphoton") {
0506     pmanager->AddDiscreteProcess(theAbsorptionProcess);
0507     pmanager->AddDiscreteProcess(theBoundaryProcess);
0508       }
0509     }
0510 }
0511 
0512 // Hadronic processes ////////////////////////////////////////////////////////
0513 
0514 void DMXPhysicsList::ConstructHad() 
0515 {
0516   //Elastic models
0517   G4HadronElastic* elastic_lhep0 = new G4HadronElastic();
0518   G4ChipsElasticModel* elastic_chip = new G4ChipsElasticModel();
0519   G4ElasticHadrNucleusHE* elastic_he = new G4ElasticHadrNucleusHE(); 
0520   
0521   // Inelastic scattering
0522   const G4double theFTFMin0 =    0.0*GeV;
0523   const G4double theFTFMin1 =    3.0*GeV;
0524   const G4double theFTFMax = G4HadronicParameters::Instance()->GetMaxEnergy();
0525   const G4double theBERTMin0 =   0.0*GeV;
0526   const G4double theBERTMin1 =  19.0*MeV;
0527   const G4double theBERTMax =    6.0*GeV;
0528   const G4double theHPMin =      0.0*GeV;
0529   const G4double theHPMax =     20.0*MeV;
0530 
0531   G4FTFModel * theStringModel = new G4FTFModel;
0532   G4ExcitedStringDecay * theStringDecay = new G4ExcitedStringDecay( new G4LundStringFragmentation );
0533   theStringModel->SetFragmentationModel( theStringDecay );
0534   G4PreCompoundModel * thePreEquilib = new G4PreCompoundModel( new G4ExcitationHandler );
0535   G4GeneratorPrecompoundInterface * theCascade = new G4GeneratorPrecompoundInterface( thePreEquilib );
0536 
0537   G4TheoFSGenerator * theFTFModel0 = new G4TheoFSGenerator( "FTFP" );
0538   theFTFModel0->SetHighEnergyGenerator( theStringModel );
0539   theFTFModel0->SetTransport( theCascade );
0540   theFTFModel0->SetMinEnergy( theFTFMin0 );
0541   theFTFModel0->SetMaxEnergy( theFTFMax );
0542 
0543   G4TheoFSGenerator * theFTFModel1 = new G4TheoFSGenerator( "FTFP" );
0544   theFTFModel1->SetHighEnergyGenerator( theStringModel );
0545   theFTFModel1->SetTransport( theCascade );
0546   theFTFModel1->SetMinEnergy( theFTFMin1 );
0547   theFTFModel1->SetMaxEnergy( theFTFMax );
0548 
0549   G4CascadeInterface * theBERTModel0 = new G4CascadeInterface;
0550   theBERTModel0->SetMinEnergy( theBERTMin0 );
0551   theBERTModel0->SetMaxEnergy( theBERTMax );
0552 
0553   G4CascadeInterface * theBERTModel1 = new G4CascadeInterface;
0554   theBERTModel1->SetMinEnergy( theBERTMin1 );
0555   theBERTModel1->SetMaxEnergy( theBERTMax );
0556 
0557   G4VCrossSectionDataSet * theAntiNucleonData = new G4CrossSectionInelastic( new G4ComponentAntiNuclNuclearXS );
0558   G4ComponentGGNuclNuclXsc * ggNuclNuclXsec = new G4ComponentGGNuclNuclXsc();
0559   G4VCrossSectionDataSet * theGGNuclNuclData = new G4CrossSectionInelastic(ggNuclNuclXsec);
0560   G4VCrossSectionDataSet * theGGNNEl = new G4CrossSectionElastic(ggNuclNuclXsec);
0561   G4ComponentGGHadronNucleusXsc * ggHNXsec = new G4ComponentGGHadronNucleusXsc();
0562   G4VCrossSectionDataSet * theGGHNEl = new G4CrossSectionElastic(ggHNXsec);
0563   G4VCrossSectionDataSet * theGGHNInel = new G4CrossSectionInelastic(ggHNXsec);
0564 
0565   auto particleIterator=GetParticleIterator();
0566   particleIterator->reset();
0567   while ((*particleIterator)())
0568     {
0569       G4ParticleDefinition* particle = particleIterator->value();
0570       G4ProcessManager* pmanager = particle->GetProcessManager();
0571       G4String particleName = particle->GetParticleName();
0572 
0573       if (particleName == "pi+") 
0574     {
0575       // Elastic scattering
0576           G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
0577           theElasticProcess->AddDataSet( new G4BGGPionElasticXS( particle ) );
0578           theElasticProcess->RegisterMe( elastic_he );
0579       pmanager->AddDiscreteProcess( theElasticProcess );
0580       //Inelastic scattering
0581       G4HadronInelasticProcess* theInelasticProcess = 
0582         new G4HadronInelasticProcess( "inelastic", G4PionPlus::Definition() );
0583       theInelasticProcess->AddDataSet( new G4BGGPionElasticXS( particle ) );
0584       theInelasticProcess->RegisterMe( theFTFModel1 );
0585           theInelasticProcess->RegisterMe( theBERTModel0 );
0586       pmanager->AddDiscreteProcess( theInelasticProcess );
0587     } 
0588 
0589       else if (particleName == "pi-") 
0590     {
0591       // Elastic scattering
0592           G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
0593           theElasticProcess->AddDataSet( new G4BGGPionElasticXS( particle ) );
0594           theElasticProcess->RegisterMe( elastic_he );
0595       pmanager->AddDiscreteProcess( theElasticProcess );
0596       //Inelastic scattering
0597       G4HadronInelasticProcess* theInelasticProcess = 
0598         new G4HadronInelasticProcess( "inelastic", G4PionMinus::Definition() );
0599       theInelasticProcess->AddDataSet( new G4BGGPionInelasticXS( particle ) );
0600       theInelasticProcess->RegisterMe( theFTFModel1 );
0601           theInelasticProcess->RegisterMe( theBERTModel0 );
0602       pmanager->AddDiscreteProcess( theInelasticProcess );    
0603       //Absorption
0604       pmanager->AddRestProcess(new G4HadronicAbsorptionBertini(G4PionMinus::Definition()), ordDefault);
0605     }
0606       else if (particleName == "kaon+") 
0607     {
0608       // Elastic scattering
0609           G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
0610       theElasticProcess->AddDataSet( theGGHNEl );
0611           theElasticProcess->RegisterMe( elastic_lhep0 );
0612       pmanager->AddDiscreteProcess( theElasticProcess );
0613           // Inelastic scattering   
0614       G4HadronInelasticProcess* theInelasticProcess = 
0615         new G4HadronInelasticProcess( "inelastic", G4KaonPlus::Definition() );
0616       theInelasticProcess->AddDataSet( theGGHNInel );
0617       theInelasticProcess->RegisterMe( theFTFModel1 );
0618           theInelasticProcess->RegisterMe( theBERTModel0 );
0619       pmanager->AddDiscreteProcess( theInelasticProcess );
0620     }      
0621       else if (particleName == "kaon0S") 
0622     {
0623       // Elastic scattering
0624           G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
0625       theElasticProcess->AddDataSet( theGGHNEl );
0626           theElasticProcess->RegisterMe( elastic_lhep0 );
0627       pmanager->AddDiscreteProcess( theElasticProcess );
0628           // Inelastic scattering    
0629       G4HadronInelasticProcess* theInelasticProcess = 
0630         new G4HadronInelasticProcess( "inelastic", G4KaonZeroShort::Definition() );
0631       theInelasticProcess->AddDataSet( theGGHNInel );
0632       theInelasticProcess->RegisterMe( theFTFModel1 );
0633           theInelasticProcess->RegisterMe( theBERTModel0 );
0634       pmanager->AddDiscreteProcess( theInelasticProcess );    
0635     }
0636 
0637       else if (particleName == "kaon0L") 
0638     {
0639       // Elastic scattering
0640           G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
0641       theElasticProcess->AddDataSet( theGGHNEl );
0642           theElasticProcess->RegisterMe( elastic_lhep0 );
0643       pmanager->AddDiscreteProcess( theElasticProcess );
0644       // Inelastic scattering
0645       G4HadronInelasticProcess* theInelasticProcess = 
0646         new G4HadronInelasticProcess( "inelastic", G4KaonZeroLong::Definition() );
0647       theInelasticProcess->AddDataSet( theGGHNInel );
0648       theInelasticProcess->RegisterMe( theFTFModel1 );
0649           theInelasticProcess->RegisterMe( theBERTModel0 ); 
0650       pmanager->AddDiscreteProcess( theInelasticProcess );    
0651     }
0652 
0653       else if (particleName == "kaon-") 
0654     {
0655       // Elastic scattering
0656           G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
0657       theElasticProcess->AddDataSet( theGGHNEl );
0658           theElasticProcess->RegisterMe( elastic_lhep0 );
0659       pmanager->AddDiscreteProcess( theElasticProcess );
0660           // Inelastic scattering
0661       G4HadronInelasticProcess* theInelasticProcess = 
0662         new G4HadronInelasticProcess( "inelastic", G4KaonMinus::Definition() ); 
0663           theInelasticProcess->AddDataSet( theGGHNInel );
0664       theInelasticProcess->RegisterMe( theFTFModel1 );
0665           theInelasticProcess->RegisterMe( theBERTModel0 );
0666       pmanager->AddDiscreteProcess( theInelasticProcess );
0667       pmanager->AddRestProcess(new G4HadronicAbsorptionBertini(G4KaonMinus::Definition()), ordDefault);
0668     }
0669 
0670       else if (particleName == "proton") 
0671     {
0672       // Elastic scattering
0673           G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
0674           theElasticProcess->AddDataSet( new G4BGGNucleonElasticXS( G4Proton::Proton() ) );
0675           theElasticProcess->RegisterMe( elastic_chip );
0676       pmanager->AddDiscreteProcess( theElasticProcess );
0677       // Inelastic scattering
0678       G4HadronInelasticProcess* theInelasticProcess = 
0679         new G4HadronInelasticProcess( "inelastic", G4Proton::Definition() );
0680       theInelasticProcess->AddDataSet( new G4BGGNucleonInelasticXS( G4Proton::Proton() ) );
0681       theInelasticProcess->RegisterMe( theFTFModel1 );
0682           theInelasticProcess->RegisterMe( theBERTModel0 );
0683       pmanager->AddDiscreteProcess( theInelasticProcess );
0684     }
0685       else if (particleName == "anti_proton") 
0686     {
0687       // Elastic scattering
0688           const G4double elastic_elimitAntiNuc = 100.0*MeV;
0689           G4AntiNuclElastic* elastic_anuc = new G4AntiNuclElastic();
0690           elastic_anuc->SetMinEnergy( elastic_elimitAntiNuc );
0691           G4CrossSectionElastic* elastic_anucxs = new G4CrossSectionElastic( elastic_anuc->GetComponentCrossSection() );
0692           G4HadronElastic* elastic_lhep2 = new G4HadronElastic();
0693           elastic_lhep2->SetMaxEnergy( elastic_elimitAntiNuc );
0694           G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
0695           theElasticProcess->AddDataSet( elastic_anucxs );
0696           theElasticProcess->RegisterMe( elastic_lhep2 );
0697           theElasticProcess->RegisterMe( elastic_anuc );
0698       pmanager->AddDiscreteProcess( theElasticProcess );
0699       // Inelastic scattering
0700       G4HadronInelasticProcess* theInelasticProcess = 
0701         new G4HadronInelasticProcess( "inelastic", G4AntiProton::Definition() );
0702       theInelasticProcess->AddDataSet( theAntiNucleonData );
0703       theInelasticProcess->RegisterMe( theFTFModel0 );
0704       pmanager->AddDiscreteProcess( theInelasticProcess );
0705       // Absorption
0706       pmanager->AddRestProcess(new G4HadronicAbsorptionFritiof(G4AntiProton::Definition()), ordDefault);
0707     }
0708       else if (particleName == "neutron") {
0709     // elastic scattering
0710     G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
0711         theElasticProcess->AddDataSet(new G4NeutronElasticXS());
0712         G4HadronElastic* elastic_neutronChipsModel = new G4ChipsElasticModel();
0713     elastic_neutronChipsModel->SetMinEnergy( 19.0*MeV );
0714         theElasticProcess->RegisterMe( elastic_neutronChipsModel );
0715     G4ParticleHPElastic * theElasticNeutronHP = new G4ParticleHPElastic;
0716         theElasticNeutronHP->SetMinEnergy( theHPMin );
0717         theElasticNeutronHP->SetMaxEnergy( theHPMax );
0718     theElasticProcess->RegisterMe( theElasticNeutronHP );
0719     theElasticProcess->AddDataSet( new G4ParticleHPElasticData );
0720     pmanager->AddDiscreteProcess( theElasticProcess );
0721     // inelastic scattering     
0722     G4HadronInelasticProcess* theInelasticProcess =
0723       new G4HadronInelasticProcess( "inelastic", G4Neutron::Definition() );
0724     theInelasticProcess->AddDataSet( new G4NeutronInelasticXS() );
0725     theInelasticProcess->RegisterMe( theFTFModel1 );
0726         theInelasticProcess->RegisterMe( theBERTModel1 );
0727     G4ParticleHPInelastic * theNeutronInelasticHPModel = new G4ParticleHPInelastic;
0728         theNeutronInelasticHPModel->SetMinEnergy( theHPMin );
0729         theNeutronInelasticHPModel->SetMaxEnergy( theHPMax );
0730     theInelasticProcess->RegisterMe( theNeutronInelasticHPModel );
0731     theInelasticProcess->AddDataSet( new G4ParticleHPInelasticData );
0732     pmanager->AddDiscreteProcess(theInelasticProcess);
0733     // capture
0734     G4NeutronCaptureProcess* theCaptureProcess =
0735       new G4NeutronCaptureProcess;
0736     G4ParticleHPCapture * theLENeutronCaptureModel = new G4ParticleHPCapture;
0737     theLENeutronCaptureModel->SetMinEnergy(theHPMin);
0738     theLENeutronCaptureModel->SetMaxEnergy(theHPMax);
0739     theCaptureProcess->RegisterMe(theLENeutronCaptureModel);
0740     theCaptureProcess->AddDataSet( new G4ParticleHPCaptureData);
0741     pmanager->AddDiscreteProcess(theCaptureProcess);
0742       }
0743       else if (particleName == "anti_neutron") 
0744     {
0745       // Elastic scattering
0746           G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
0747       theElasticProcess->AddDataSet( theGGHNEl );
0748           theElasticProcess->RegisterMe( elastic_lhep0 );
0749       pmanager->AddDiscreteProcess( theElasticProcess );
0750           // Inelastic scattering (include annihilation on-fly)
0751       G4HadronInelasticProcess* theInelasticProcess = 
0752         new G4HadronInelasticProcess( "inelastic", G4AntiNeutron::Definition() );
0753       theInelasticProcess->AddDataSet( theAntiNucleonData );
0754       theInelasticProcess->RegisterMe( theFTFModel0 );
0755       pmanager->AddDiscreteProcess( theInelasticProcess );    
0756     }
0757       else if (particleName == "deuteron") 
0758     {
0759       // Elastic scattering
0760           G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
0761       theElasticProcess->AddDataSet( theGGNNEl );
0762           theElasticProcess->RegisterMe( elastic_lhep0 );
0763       pmanager->AddDiscreteProcess( theElasticProcess );
0764           // Inelastic scattering
0765       G4HadronInelasticProcess* theInelasticProcess = 
0766         new G4HadronInelasticProcess( "inelastic", G4Deuteron::Definition() );
0767       theInelasticProcess->AddDataSet( theGGNuclNuclData );
0768       theInelasticProcess->RegisterMe( theFTFModel1 );
0769           theInelasticProcess->RegisterMe( theBERTModel0 );
0770       pmanager->AddDiscreteProcess( theInelasticProcess );
0771     }
0772       else if (particleName == "triton") 
0773     {
0774       // Elastic scattering
0775           G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
0776       theElasticProcess->AddDataSet( theGGNNEl );
0777           theElasticProcess->RegisterMe( elastic_lhep0 );
0778       pmanager->AddDiscreteProcess( theElasticProcess );
0779           // Inelastic scattering
0780       G4HadronInelasticProcess* theInelasticProcess = 
0781         new G4HadronInelasticProcess( "inelastic", G4Triton::Definition() );
0782       theInelasticProcess->AddDataSet( theGGNuclNuclData );
0783       theInelasticProcess->RegisterMe( theFTFModel1 );
0784           theInelasticProcess->RegisterMe( theBERTModel0 );
0785       pmanager->AddDiscreteProcess( theInelasticProcess );
0786     }
0787       else if (particleName == "alpha") 
0788     {
0789       // Elastic scattering
0790           G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
0791       theElasticProcess->AddDataSet( theGGNNEl );
0792           theElasticProcess->RegisterMe( elastic_lhep0 );
0793       pmanager->AddDiscreteProcess( theElasticProcess );
0794           // Inelastic scattering
0795       G4HadronInelasticProcess* theInelasticProcess = 
0796         new G4HadronInelasticProcess( "inelastic", G4Alpha::Definition() );  
0797           theInelasticProcess->AddDataSet( theGGNuclNuclData );
0798       theInelasticProcess->RegisterMe( theFTFModel1 );
0799           theInelasticProcess->RegisterMe( theBERTModel0 );
0800       pmanager->AddDiscreteProcess( theInelasticProcess );
0801     }
0802     }
0803 }
0804 
0805 // Decays ///////////////////////////////////////////////////////////////////
0806 void DMXPhysicsList::ConstructGeneral() {
0807 
0808   // Add Decay Process
0809   G4Decay* theDecayProcess = new G4Decay();
0810   auto particleIterator=GetParticleIterator();
0811   particleIterator->reset();
0812   while( (*particleIterator)() )
0813     {
0814       G4ParticleDefinition* particle = particleIterator->value();
0815       G4ProcessManager* pmanager = particle->GetProcessManager();
0816       
0817       if (theDecayProcess->IsApplicable(*particle) && !particle->IsShortLived()) 
0818     { 
0819       pmanager ->AddProcess(theDecayProcess);
0820       // set ordering for PostStepDoIt and AtRestDoIt
0821       pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep);
0822       pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest);
0823     }
0824     }
0825 
0826   // Declare radioactive decay to the GenericIon in the IonTable.
0827   G4LossTableManager* man = G4LossTableManager::Instance();
0828   G4VAtomDeexcitation* ad = man->AtomDeexcitation();
0829   if(!ad) {
0830     G4EmParameters::Instance()->SetAugerCascade(true);
0831     ad = new G4UAtomicDeexcitation();
0832     man->SetAtomDeexcitation(ad);
0833     ad->InitialiseAtomicDeexcitation();
0834   }
0835 
0836   G4PhysicsListHelper::GetPhysicsListHelper()->
0837     RegisterProcess(new G4RadioactiveDecay(), G4GenericIon::GenericIon());
0838 }
0839 
0840 // Cuts /////////////////////////////////////////////////////////////////////
0841 void DMXPhysicsList::SetCuts() 
0842 {
0843   
0844   if (verboseLevel >1)
0845     G4cout << "DMXPhysicsList::SetCuts:";
0846   
0847   if (verboseLevel>0){
0848     G4cout << "DMXPhysicsList::SetCuts:";
0849     G4cout << "CutLength : " 
0850        << G4BestUnit(defaultCutValue,"Length") << G4endl;
0851   }
0852 
0853   //special for low energy physics
0854   G4double lowlimit=250*eV;  
0855   G4ProductionCutsTable::GetProductionCutsTable()->SetEnergyRange(lowlimit,100.*GeV);
0856 
0857   // set cut values for gamma at first and for e- second and next for e+,
0858   // because some processes for e+/e- need cut values for gamma 
0859   SetCutValue(cutForGamma, "gamma");
0860   SetCutValue(cutForElectron, "e-");
0861   SetCutValue(cutForPositron, "e+");
0862   
0863   if (verboseLevel>0) DumpCutValuesTable();
0864 }
0865