#include "HadronPhysics.h" #include "G4MesonConstructor.hh" #include "G4LeptonConstructor.hh" #include "G4BaryonConstructor.hh" #include "G4ShortLivedConstructor.hh" #include "G4IonConstructor.hh" // processes #include "G4ProcessManager.hh" #include "G4HadronElasticProcess.hh" #include "G4PionPlusInelasticProcess.hh" #include "G4PionMinusInelasticProcess.hh" #include "G4KaonPlusInelasticProcess.hh" #include "G4KaonMinusInelasticProcess.hh" #include "G4KaonZeroLInelasticProcess.hh" #include "G4KaonZeroSInelasticProcess.hh" #include "G4ProtonInelasticProcess.hh" #include "G4NeutronInelasticProcess.hh" #include "G4HadronFissionProcess.hh" #include "G4HadronCaptureProcess.hh" #include "G4AntiProtonInelasticProcess.hh" #include "G4AntiNeutronInelasticProcess.hh" #include "G4LambdaInelasticProcess.hh" #include "G4AntiLambdaInelasticProcess.hh" #include "G4SigmaPlusInelasticProcess.hh" #include "G4SigmaMinusInelasticProcess.hh" #include "G4AntiSigmaPlusInelasticProcess.hh" #include "G4AntiSigmaMinusInelasticProcess.hh" #include "G4PionMinusAbsorptionAtRest.hh" #include "G4KaonMinusAbsorption.hh" #include "G4AntiProtonAnnihilationAtRest.hh" #include "G4AntiNeutronAnnihilationAtRest.hh" // Lepton nuclear interations #include "G4PhotoNuclearProcess.hh" #include "G4GammaNuclearReaction.hh" #include "G4ElectroNuclearReaction.hh" #include "G4ElectronNuclearProcess.hh" // cross sections #include "G4PiNuclearCrossSection.hh" #include "G4ProtonInelasticCrossSection.hh" #include "G4NeutronInelasticCrossSection.hh" // models #include "G4LElastic.hh" #include "G4CascadeInterface.hh" #include "G4LEPionPlusInelastic.hh" #include "G4LEPionMinusInelastic.hh" #include "G4LEKaonPlusInelastic.hh" #include "G4LEKaonMinusInelastic.hh" #include "G4LEKaonZeroLInelastic.hh" #include "G4LEKaonZeroSInelastic.hh" #include "G4LEProtonInelastic.hh" #include "G4LENeutronInelastic.hh" #include "G4LFission.hh" #include "G4LCapture.hh" #include "G4LEAntiProtonInelastic.hh" #include "G4LEAntiNeutronInelastic.hh" #include "G4LELambdaInelastic.hh" #include "G4LEAntiLambdaInelastic.hh" #include "G4LESigmaPlusInelastic.hh" #include "G4LESigmaMinusInelastic.hh" #include "G4LEAntiSigmaPlusInelastic.hh" #include "G4LEAntiSigmaMinusInelastic.hh" HadronPhysics::HadronPhysics(gemc_opts Opt):G4VPhysicsConstructor("Hadron Physics") { gemcOpt = Opt; } HadronPhysics::~HadronPhysics(){;} void HadronPhysics::ConstructParticle() { // Construct all mesons G4MesonConstructor pMesonConstructor; pMesonConstructor.ConstructParticle(); // Construct all leptons G4LeptonConstructor pLeptonConstructor; pLeptonConstructor.ConstructParticle(); // Construct all baryons G4BaryonConstructor pBaryonConstructor; pBaryonConstructor.ConstructParticle(); // Construct resonaces and quarks G4ShortLivedConstructor pShortLivedConstructor; pShortLivedConstructor.ConstructParticle(); // Construct light ions (d, t, 3He, alpha, and generic ion) G4IonConstructor ionConstruct; ionConstruct.ConstructParticle(); } void HadronPhysics::ConstructProcess() { string hd_msg = gemcOpt.args["LOG_MSG"].args + " Hadron Physics List: <<< "; double PHY_VERBOSITY = gemcOpt.args["PHY_VERBOSITY"].arg ; cout << hd_msg << " Building Hadron Physics Processes " << endl; // Hadronic Elastic Process and Model (the same for all hadrons) G4HadronElasticProcess* elasticProcess = new G4HadronElasticProcess(); G4LElastic* elasticModel = new G4LElastic(); elasticProcess->RegisterMe(elasticModel); // %%%%%%%%%%%%%%%%%%%%%%%%% // Hadronic inelastic models // %%%%%%%%%%%%%%%%%%%%%%%%% // Bertini cascade model: use for p,n,pi+,pi- between 0 and 9.9 GeV G4CascadeInterface* bertiniModel = new G4CascadeInterface(); bertiniModel->SetMaxEnergy(9.9*GeV); // Low energy parameterized models : use between 9.5 and 25 GeV G4double LEPUpperLimit = 25*GeV; G4double LEPpnpiLimit = 9.5*GeV; G4LEKaonZeroLInelastic* LEPk0LModel = new G4LEKaonZeroLInelastic(); LEPk0LModel->SetMaxEnergy(LEPUpperLimit); G4LEKaonZeroSInelastic* LEPk0SModel = new G4LEKaonZeroSInelastic(); LEPk0SModel->SetMaxEnergy(LEPUpperLimit); G4ProcessManager * pmanager = 0; // Add Hadron Elastic Processes to all hadrons theParticleIterator->reset(); while( (*theParticleIterator)() ) { G4ParticleDefinition* particle = theParticleIterator->value(); pmanager = particle->GetProcessManager(); string pname = particle->GetParticleName(); if (elasticProcess->IsApplicable(*particle)) { if(PHY_VERBOSITY > 2) cout << hd_msg << " Adding Hadron Elastic Process for " << pname << endl; pmanager->AddDiscreteProcess(elasticProcess); } } // Adding inelastic nuclear scattering to leptons G4PhotoNuclearProcess* thePhotoNuclearProcess = new G4PhotoNuclearProcess; G4GammaNuclearReaction* theGammaReaction = new G4GammaNuclearReaction; theGammaReaction->SetMaxEnergy(10*GeV); thePhotoNuclearProcess->RegisterMe(theGammaReaction); G4ElectronNuclearProcess* theElectronNuclearProcess = new G4ElectronNuclearProcess; G4ElectroNuclearReaction* theElectroReaction = new G4ElectroNuclearReaction; theElectroReaction->SetMaxEnergy(10*GeV); theElectronNuclearProcess->RegisterMe(theElectroReaction); theParticleIterator->reset(); while( (*theParticleIterator)() ) { G4ParticleDefinition* particle = theParticleIterator->value(); pmanager = particle->GetProcessManager(); string pname = particle->GetParticleName(); if (pname == "gamma") { if(PHY_VERBOSITY > 2) cout << hd_msg << " Gamma: inelastic nuclear scattering" << endl; pmanager->AddDiscreteProcess(thePhotoNuclearProcess); } else if (pname == "e-") { if(PHY_VERBOSITY > 2) cout << hd_msg << " Electron: inelastic nuclear scattering" << endl; pmanager->AddDiscreteProcess(thePhotoNuclearProcess); } else if (pname == "e+") { if(PHY_VERBOSITY > 2) cout << hd_msg << " Positron: inelastic nuclear scattering" << endl; pmanager->AddDiscreteProcess(thePhotoNuclearProcess); } } // pi+ physics // pmanager = G4PionPlus::PionPlus()->GetProcessManager(); G4PionPlusInelasticProcess* pipinelProc = new G4PionPlusInelasticProcess(); G4PiNuclearCrossSection* pion_XC = new G4PiNuclearCrossSection(); pipinelProc->AddDataSet(pion_XC); pipinelProc->RegisterMe(bertiniModel); G4LEPionPlusInelastic* LEPpipModel = new G4LEPionPlusInelastic(); LEPpipModel->SetMinEnergy(LEPpnpiLimit); LEPpipModel->SetMaxEnergy(LEPUpperLimit); pipinelProc->RegisterMe(LEPpipModel); pmanager->AddDiscreteProcess(pipinelProc); // pi- physics // pmanager = G4PionMinus::PionMinus()->GetProcessManager(); G4PionMinusInelasticProcess* piminelProc = new G4PionMinusInelasticProcess(); piminelProc->AddDataSet(pion_XC); piminelProc->RegisterMe(bertiniModel); G4LEPionMinusInelastic* LEPpimModel = new G4LEPionMinusInelastic(); LEPpimModel->SetMinEnergy(LEPpnpiLimit); LEPpimModel->SetMaxEnergy(LEPUpperLimit); piminelProc->RegisterMe(LEPpimModel); pmanager->AddDiscreteProcess(piminelProc); // pi- absorption at rest G4PionMinusAbsorptionAtRest* pimAbsorb = new G4PionMinusAbsorptionAtRest(); pmanager->AddRestProcess(pimAbsorb); // K+ physics // pmanager = G4KaonPlus::KaonPlus()->GetProcessManager(); G4KaonPlusInelasticProcess* kpinelProc = new G4KaonPlusInelasticProcess(); G4LEKaonPlusInelastic* LEPkpModel = new G4LEKaonPlusInelastic(); LEPkpModel->SetMaxEnergy(LEPUpperLimit); kpinelProc->RegisterMe(LEPkpModel); pmanager->AddDiscreteProcess(kpinelProc); // K- physics // pmanager = G4KaonMinus::KaonMinus()->GetProcessManager(); G4KaonMinusInelasticProcess* kminelProc = new G4KaonMinusInelasticProcess(); G4LEKaonMinusInelastic* LEPkmModel = new G4LEKaonMinusInelastic(); LEPkmModel->SetMaxEnergy(LEPUpperLimit); kminelProc->RegisterMe(LEPkmModel); pmanager->AddDiscreteProcess(kminelProc); // K- absorption at rest G4KaonMinusAbsorption* kmAbsorb = new G4KaonMinusAbsorption(); pmanager->AddRestProcess(kmAbsorb); // K0L physics // pmanager = G4KaonZeroLong::KaonZeroLong()->GetProcessManager(); G4KaonZeroLInelasticProcess* k0LinelProc = new G4KaonZeroLInelasticProcess(); k0LinelProc->RegisterMe(LEPk0LModel); pmanager->AddDiscreteProcess(k0LinelProc); // K0S physics // pmanager = G4KaonZeroShort::KaonZeroShort()->GetProcessManager(); G4KaonZeroSInelasticProcess* k0SinelProc = new G4KaonZeroSInelasticProcess(); k0SinelProc->RegisterMe(LEPk0SModel); pmanager->AddDiscreteProcess(k0SinelProc); // Proton // pmanager = G4Proton::Proton()->GetProcessManager(); G4ProtonInelasticProcess* pinelProc = new G4ProtonInelasticProcess(); G4ProtonInelasticCrossSection* proton_XC = new G4ProtonInelasticCrossSection(); pinelProc->AddDataSet(proton_XC); pinelProc->RegisterMe(bertiniModel); G4LEProtonInelastic* LEPpModel = new G4LEProtonInelastic(); LEPpModel->SetMinEnergy(LEPpnpiLimit); LEPpModel->SetMaxEnergy(LEPUpperLimit); pinelProc->RegisterMe(LEPpModel); pmanager->AddDiscreteProcess(pinelProc); // Anti-Proton // pmanager = G4AntiProton::AntiProton()->GetProcessManager(); G4AntiProtonInelasticProcess* apinelProc = new G4AntiProtonInelasticProcess(); G4LEAntiProtonInelastic* LEPapModel = new G4LEAntiProtonInelastic(); apinelProc->RegisterMe(LEPapModel); pmanager->AddDiscreteProcess(apinelProc); // anti-proton annihilation at rest G4AntiProtonAnnihilationAtRest* apAnnihil = new G4AntiProtonAnnihilationAtRest(); pmanager->AddRestProcess(apAnnihil); // Neutron // pmanager = G4Neutron::Neutron()->GetProcessManager(); // hadron elastic pmanager->AddDiscreteProcess(elasticProcess); // hadron inelastic G4NeutronInelasticProcess* ninelProc = new G4NeutronInelasticProcess(); G4NeutronInelasticCrossSection* neutron_XC = new G4NeutronInelasticCrossSection(); ninelProc->AddDataSet(neutron_XC); ninelProc->RegisterMe(bertiniModel); G4LENeutronInelastic* LEPnModel = new G4LENeutronInelastic(); LEPnModel->SetMinEnergy(LEPpnpiLimit); LEPnModel->SetMaxEnergy(LEPUpperLimit); ninelProc->RegisterMe(LEPnModel); pmanager->AddDiscreteProcess(ninelProc); // neutron-induced fission G4HadronFissionProcess* neutronFission = new G4HadronFissionProcess(); G4LFission* neutronFissionModel = new G4LFission(); neutronFissionModel->SetMinEnergy(0.); neutronFissionModel->SetMaxEnergy(20*TeV); neutronFission->RegisterMe(neutronFissionModel); pmanager->AddDiscreteProcess(neutronFission); // neutron capture G4HadronCaptureProcess* neutronCapture = new G4HadronCaptureProcess(); G4LCapture* neutronCaptureModel = new G4LCapture(); neutronCaptureModel->SetMinEnergy(0.); neutronCaptureModel->SetMaxEnergy(20*TeV); neutronCapture->RegisterMe(neutronCaptureModel); pmanager->AddDiscreteProcess(neutronCapture); // Anti-Neutron // pmanager = G4AntiNeutron::AntiNeutron()->GetProcessManager(); G4AntiNeutronInelasticProcess* aninelProc = new G4AntiNeutronInelasticProcess(); G4LEAntiNeutronInelastic* LEPanModel = new G4LEAntiNeutronInelastic(); aninelProc->RegisterMe(LEPanModel); pmanager->AddDiscreteProcess(aninelProc); // anti-neutron annihilation at rest G4AntiNeutronAnnihilationAtRest* anAnnihil = new G4AntiNeutronAnnihilationAtRest(); pmanager->AddRestProcess(anAnnihil); // Lambda // pmanager = G4Lambda::Lambda()->GetProcessManager(); G4LambdaInelasticProcess* linelProc = new G4LambdaInelasticProcess(); G4LELambdaInelastic* LEPlModel = new G4LELambdaInelastic(); linelProc->RegisterMe(LEPlModel); pmanager->AddDiscreteProcess(linelProc); // Anti-Lambda // pmanager = G4AntiLambda::AntiLambda()->GetProcessManager(); G4AntiLambdaInelasticProcess* alinelProc = new G4AntiLambdaInelasticProcess(); G4LEAntiLambdaInelastic* LEPalModel = new G4LEAntiLambdaInelastic(); alinelProc->RegisterMe(LEPalModel); pmanager->AddDiscreteProcess(alinelProc); // Sigma- // pmanager = G4SigmaMinus::SigmaMinus()->GetProcessManager(); G4SigmaMinusInelasticProcess* sminelProc = new G4SigmaMinusInelasticProcess(); G4LESigmaMinusInelastic* LEPsmModel = new G4LESigmaMinusInelastic(); sminelProc->RegisterMe(LEPsmModel); pmanager->AddDiscreteProcess(sminelProc); // Anti-Sigma- // pmanager = G4AntiSigmaMinus::AntiSigmaMinus()->GetProcessManager(); G4AntiSigmaMinusInelasticProcess* asminelProc = new G4AntiSigmaMinusInelasticProcess(); G4LEAntiSigmaMinusInelastic* LEPasmModel = new G4LEAntiSigmaMinusInelastic(); asminelProc->RegisterMe(LEPasmModel); pmanager->AddDiscreteProcess(asminelProc); // Sigma+ // pmanager = G4SigmaPlus::SigmaPlus()->GetProcessManager(); G4SigmaPlusInelasticProcess* spinelProc = new G4SigmaPlusInelasticProcess(); G4LESigmaPlusInelastic* LEPspModel = new G4LESigmaPlusInelastic(); spinelProc->RegisterMe(LEPspModel); pmanager->AddDiscreteProcess(spinelProc); // Anti-Sigma+ // pmanager = G4AntiSigmaPlus::AntiSigmaPlus()->GetProcessManager(); G4AntiSigmaPlusInelasticProcess* aspinelProc = new G4AntiSigmaPlusInelasticProcess(); G4LEAntiSigmaPlusInelastic* LEPaspModel = new G4LEAntiSigmaPlusInelastic(); aspinelProc->RegisterMe(LEPaspModel); pmanager->AddDiscreteProcess(aspinelProc); }