#include "myDetPhysicsList.hh" #include "G4ProcessManager.hh" #include "G4ParticleDefinition.hh" #include "G4ParticleWithCuts.hh" #include "G4ParticleTypes.hh" #include "G4ParticleTable.hh" //particle constructors #include "G4LeptonConstructor.hh" #include "G4IonConstructor.hh" #include "G4MesonConstructor.hh" #include "G4BaryonConstructor.hh" #include "G4ShortLivedConstructor.hh" //EM interactions: //gamma #include "G4LowEnergyCompton.hh" #include "G4LowEnergyGammaConversion.hh" #include "G4LowEnergyPhotoElectric.hh" #include "G4LowEnergyRayleigh.hh" //e+e- #include "G4LowEnergyIonisation.hh" #include "G4LowEnergyBremsstrahlung.hh" #include "G4eMultipleScattering.hh" #include "G4eIonisation.hh" #include "G4eBremsstrahlung.hh" #include "G4eplusAnnihilation.hh" //mu+,mu- #include "G4MuIonisation.hh" #include "G4MuBremsstrahlung.hh" #include "G4MuPairProduction.hh" //proton, neutron, other baryons #include "G4hIonisation.hh" #include "G4hBremsstrahlung.hh" #include "G4hPairProduction.hh" #include "G4hMultipleScattering.hh" //ions #include "G4ionIonisation.hh" //Hadronic physics: // Hadronics #include "G4HadronElasticProcess.hh" #include "G4HadronInelasticProcess.hh" #include "G4HadronCaptureProcess.hh" #include "G4HadronFissionProcess.hh" //Elastic model #include "G4LElastic.hh" //Inelastic Processes #include "G4ProtonInelasticProcess.hh" #include "G4NeutronInelasticProcess.hh" #include "G4DeuteronInelasticProcess.hh" #include "G4TritonInelasticProcess.hh" #include "G4AlphaInelasticProcess.hh" #include "G4IonInelasticProcess.hh" //Inelastic models and data #include "G4Evaporation.hh" #include "G4FermiBreakUp.hh" #include "G4StatMF.hh" #include "G4ExcitationHandler.hh" #include "G4LEDeuteronInelastic.hh" #include "G4LETritonInelastic.hh" #include "G4LEAlphaInelastic.hh" #include "G4WilsonAbrasionModel.hh" #include "G4WilsonAblationModel.hh" // Kinetic Model #include "G4BinaryCascade.hh" #include "G4BinaryLightIonReaction.hh" #include "G4TripathiCrossSection.hh" #include "G4IonsShenCrossSection.hh" //Neutron HP Models: Thermal to < 20 MeV #include "G4NeutronHPElasticData.hh" #include "G4NeutronHPJENDLHEData.hh" #include "G4NeutronHPorLElasticModel.hh" #include "G4NeutronHPCapture.hh" #include "G4NeutronHPCaptureData.hh" #include "G4NeutronHPorLCapture.hh" #include "G4NeutronHPFission.hh" #include "G4NeutronHPFissionData.hh" #include "G4NeutronHPorLFission.hh" #include "G4NeutronHPInelastic.hh" #include "G4NeutronHPInelasticData.hh" #include "G4NeutronHPorLEInelastic.hh" #include "G4NeutronHPorLEInelasticData.hh" #include "G4NeutronHPThermalScattering.hh" #include "G4NeutronHPThermalScatteringData.hh" //decay #include "G4Decay.hh" //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... myDetPhysicsList::myDetPhysicsList(): G4VUserPhysicsList() { SetVerboseLevel(1); } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... myDetPhysicsList::~myDetPhysicsList() {} //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void myDetPhysicsList::ConstructParticle() { // In this method, static member functions should be called // for all particles which you want to use. // This ensures that objects of these particle types will be // created in the program. //Bosons G4Gamma::GammaDefinition(); //Leptons G4LeptonConstructor lConstructor; lConstructor.ConstructParticle(); //Ions G4IonConstructor iConstructor; iConstructor.ConstructParticle(); //Mesons G4MesonConstructor mConstructor; mConstructor.ConstructParticle(); //Baryons G4BaryonConstructor bConstructor; bConstructor.ConstructParticle(); //short-lived G4ShortLivedConstructor shortLived; shortLived.ConstructParticle(); } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void myDetPhysicsList::ConstructProcess() { AddTransportation(); ConstructGeneral(); ConstructEM(); ConstructHadronic(); } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void myDetPhysicsList::AddTransportation() { G4VUserPhysicsList::AddTransportation(); } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void myDetPhysicsList::ConstructGeneral() { // Add Decay Process G4Decay* theDecayProcess = new G4Decay(); theParticleIterator->reset(); while( (*theParticleIterator)() ){ G4ParticleDefinition* particle = theParticleIterator->value(); G4ProcessManager* pmanager = particle->GetProcessManager(); if (theDecayProcess->IsApplicable(*particle)) { pmanager ->AddProcess(theDecayProcess); // set ordering for PostStepDoIt and AtRestDoIt pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep); pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest); } } } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void myDetPhysicsList::ConstructEM() { theParticleIterator->reset(); while( (*theParticleIterator)() ){ G4ParticleDefinition* particle = theParticleIterator->value(); G4ProcessManager* pmanager = particle->GetProcessManager(); G4String particleName = particle->GetParticleName(); if (particleName == "gamma") { // gamma pmanager->AddDiscreteProcess(new G4LowEnergyPhotoElectric()); pmanager->AddDiscreteProcess(new G4LowEnergyCompton()); pmanager->AddDiscreteProcess(new G4LowEnergyGammaConversion()); pmanager->AddDiscreteProcess(new G4LowEnergyRayleigh()); } else if (particleName == "e-") { //electron pmanager->AddProcess(new G4eMultipleScattering(),-1, 1, 1); pmanager->AddProcess(new G4LowEnergyIonisation(), -1, 2, 2); pmanager->AddProcess(new G4LowEnergyBremsstrahlung(), -1, -1, 3); } else if (particleName == "e+") { //positron pmanager->AddProcess(new G4eMultipleScattering(),-1, 1, 1); pmanager->AddProcess(new G4eIonisation(), -1, 2, 2); pmanager->AddProcess(new G4eBremsstrahlung(), -1, -1, 3); pmanager->AddProcess(new G4eplusAnnihilation(), 0,-1, 4); }else if( particleName == "mu+" || particleName == "mu-" ) { //muon pmanager->AddProcess(new G4hMultipleScattering, -1, 1, 1); pmanager->AddProcess(new G4MuIonisation, -1, 2, 2); pmanager->AddProcess(new G4MuBremsstrahlung, -1, 3, 3); pmanager->AddProcess(new G4MuPairProduction, -1, 4, 4); } else if( particleName == "proton" || particleName == "pi-" || particleName == "pi+" ) { //proton pmanager->AddProcess(new G4hMultipleScattering, -1, 1, 1); pmanager->AddProcess(new G4hIonisation, -1, 2, 2); pmanager->AddProcess(new G4hBremsstrahlung, -1, 3, 3); pmanager->AddProcess(new G4hPairProduction, -1, 4, 4); } else if( particleName == "alpha" || particleName == "He3" || particleName == "GenericIon" ) { //Ions pmanager->AddProcess(new G4hMultipleScattering, -1, 1, 1); pmanager->AddProcess(new G4ionIonisation, -1, 2, 2); } else if ((!particle->IsShortLived()) && (particle->GetPDGCharge() != 0.0)) { //all others charged particles pmanager->AddProcess(new G4hMultipleScattering, -1, 1, 1); pmanager->AddProcess(new G4hIonisation, -1, 2, 2); } } } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void myDetPhysicsList::ConstructHadronic() { // Construct processes: //---------------------- // ***************************// // *** elastic scattering *** // // ***************************// G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess(); G4LElastic* theElasticModel = new G4LElastic(); theElasticProcess->RegisterMe(theElasticModel); // ***************************// // *** inelastic processes *** // // ***************************// G4ProtonInelasticProcess* theProtonInelasticProcess; G4DeuteronInelasticProcess* theDeuteronInelasticProcess; G4TritonInelasticProcess* theTritonInelasticProcess; G4AlphaInelasticProcess* theAlphaInelasticProcess; G4HadronInelasticProcess* theHe3InelasticProcess; G4HadronInelasticProcess* theGenericIonInelasticProcess; // all models for treatment of thermal nucleus G4Evaporation* theEvaporation = new G4Evaporation; G4FermiBreakUp* theFermiBreakUp = new G4FermiBreakUp; G4StatMF* theMF = new G4StatMF; // evaporation logic G4ExcitationHandler* theHandler = new G4ExcitationHandler; theHandler->SetEvaporation(theEvaporation); theHandler->SetFermiModel(theFermiBreakUp); theHandler->SetMultiFragmentation(theMF); theHandler->SetMaxAandZForFermiBreakUp(12, 6); theHandler->SetMinEForMultiFrag(3.*MeV); // // deuteron G4LEDeuteronInelastic* theLEDeuteronInelasticModel = new G4LEDeuteronInelastic(); theLEDeuteronInelasticModel->SetMaxEnergy(100*MeV); // // triton G4LETritonInelastic* theLETritonInelasticModel = new G4LETritonInelastic(); theLETritonInelasticModel->SetMaxEnergy(100*MeV); // // alpha G4LEAlphaInelastic* theLEAlphaInelasticModel = new G4LEAlphaInelastic(); theLEAlphaInelasticModel->SetMaxEnergy(100*MeV); // pre-equilibrium stage (model for proton inelastic) G4PreCompoundModel* thePreEquilib = new G4PreCompoundModel(theHandler); thePreEquilib->SetMaxEnergy(70*MeV); // Cascade for generic ions and He3 G4BinaryLightIonReaction* theGenIonCascade = new G4BinaryLightIonReaction; theGenIonCascade->SetMinEnergy(0*MeV); theGenIonCascade->SetMaxEnergy(30*GeV); G4TripathiCrossSection* theTripathiCrossSection = new G4TripathiCrossSection; G4IonsShenCrossSection* theShenCrossSection = new G4IonsShenCrossSection; //attach process to particles: //---------------------------- while( (*theParticleIterator)() ){ G4ParticleDefinition* particle = theParticleIterator->value(); G4ProcessManager* pmanager = particle->GetProcessManager(); G4String particleName = particle->GetParticleName(); if( particle->GetParticleName() == "neutron" ) { // elastic scattering G4HadronElasticProcess* NeutronElasticProcess = new G4HadronElasticProcess(); G4NeutronHPElasticData* HPElasticData = new G4NeutronHPElasticData(); G4NeutronHPJENDLHEData* JENDLHEElasticData = new G4NeutronHPJENDLHEData(); G4NeutronHPorLElasticModel* NeutronElasticModel = new G4NeutronHPorLElasticModel(); NeutronElasticProcess->AddDataSet(NeutronElasticModel->GiveHPXSectionDataSet()); NeutronElasticProcess->AddDataSet(HPElasticData); NeutronElasticProcess->AddDataSet(JENDLHEElasticData); NeutronElasticModel->SetMinEnergy(4.0*eV); NeutronElasticProcess->RegisterMe(NeutronElasticModel); G4NeutronHPThermalScatteringData* HPThermalScatteringData = new G4NeutronHPThermalScatteringData(); NeutronElasticProcess->AddDataSet(HPThermalScatteringData); G4NeutronHPThermalScattering* NeutronThermalElasticModel = new G4NeutronHPThermalScattering(); NeutronThermalElasticModel->SetMaxEnergy(4.0*eV); NeutronElasticProcess->RegisterMe(NeutronThermalElasticModel); pmanager->AddDiscreteProcess(NeutronElasticProcess); //Inelastic scattering G4NeutronInelasticProcess* NeutronInelasticProcess = new G4NeutronInelasticProcess(); G4NeutronHPInelasticData* HPInelasticData = new G4NeutronHPInelasticData(); G4NeutronHPorLEInelastic* NeutronInelasticModel = new G4NeutronHPorLEInelastic(); NeutronInelasticProcess->AddDataSet(NeutronInelasticModel->GiveXSectionDataSet()); NeutronInelasticProcess->AddDataSet(HPInelasticData); NeutronInelasticProcess->RegisterMe(NeutronInelasticModel); pmanager->AddDiscreteProcess(NeutronInelasticProcess); //Capture G4HadronCaptureProcess* NeutronCapture = new G4HadronCaptureProcess(); G4NeutronHPCaptureData* HPCaptureData = new G4NeutronHPCaptureData(); G4NeutronHPorLCapture* NeutronCapModel = new G4NeutronHPorLCapture(); NeutronCapture->AddDataSet(NeutronCapModel->GiveXSectionDataSet()); NeutronCapture->AddDataSet(HPCaptureData); NeutronCapture->RegisterMe(NeutronCapModel); pmanager->AddDiscreteProcess(NeutronCapture); //Fission G4HadronFissionProcess* NeutronFission = new G4HadronFissionProcess(); G4NeutronHPorLFission* FissionModel = new G4NeutronHPorLFission(); G4NeutronHPFissionData* FissionData = new G4NeutronHPFissionData(); NeutronFission->AddDataSet(FissionModel->GiveXSectionDataSet()); NeutronFission->AddDataSet(FissionData); NeutronFission->RegisterMe(FissionModel); pmanager->AddDiscreteProcess(NeutronFission); } else if (particleName == "proton") { pmanager->AddDiscreteProcess(theElasticProcess); theProtonInelasticProcess = new G4ProtonInelasticProcess; theProtonInelasticProcess->RegisterMe(thePreEquilib); pmanager->AddDiscreteProcess(theProtonInelasticProcess); } else if (particleName == "deuteron") { pmanager->AddDiscreteProcess(theElasticProcess); theDeuteronInelasticProcess = new G4DeuteronInelasticProcess; theDeuteronInelasticProcess->AddDataSet(theTripathiCrossSection); theDeuteronInelasticProcess->AddDataSet(theShenCrossSection); theDeuteronInelasticProcess->RegisterMe(theLEDeuteronInelasticModel); pmanager->AddDiscreteProcess(theDeuteronInelasticProcess); } else if (particleName == "triton") { pmanager->AddDiscreteProcess(theElasticProcess); theTritonInelasticProcess = new G4TritonInelasticProcess; theTritonInelasticProcess->AddDataSet(theTripathiCrossSection); theTritonInelasticProcess->AddDataSet(theShenCrossSection); theTritonInelasticProcess->RegisterMe(theLETritonInelasticModel); pmanager->AddDiscreteProcess(theTritonInelasticProcess); } else if (particleName == "alpha") { pmanager->AddDiscreteProcess(theElasticProcess); theAlphaInelasticProcess = new G4AlphaInelasticProcess; theAlphaInelasticProcess->AddDataSet(theTripathiCrossSection); theAlphaInelasticProcess->AddDataSet(theShenCrossSection); theAlphaInelasticProcess->RegisterMe(theLEAlphaInelasticModel); pmanager->AddDiscreteProcess(theAlphaInelasticProcess); } else if (particleName == "He3") { // NOTE elastic scattering does not stick to He3! pmanager->AddDiscreteProcess(theElasticProcess); theHe3InelasticProcess = new G4HadronInelasticProcess ("He3Inelastic", G4He3::He3()); theHe3InelasticProcess->AddDataSet(theTripathiCrossSection); theHe3InelasticProcess->AddDataSet(theShenCrossSection); theHe3InelasticProcess->RegisterMe(theGenIonCascade); pmanager->AddDiscreteProcess(theHe3InelasticProcess); } else if (particleName == "GenericIon") { pmanager->AddDiscreteProcess(theElasticProcess); theGenericIonInelasticProcess = new G4HadronInelasticProcess ("IonInelastic", G4GenericIon::GenericIon()); G4WilsonAblationModel* AblModel = new G4WilsonAblationModel(); G4WilsonAbrasionModel* AbrModel = new G4WilsonAbrasionModel(AblModel); theGenericIonInelasticProcess->RegisterMe(AbrModel); theGenericIonInelasticProcess->AddDataSet(theTripathiCrossSection); theGenericIonInelasticProcess->AddDataSet(theShenCrossSection); theGenericIonInelasticProcess->RegisterMe(theGenIonCascade); pmanager->AddDiscreteProcess(theGenericIonInelasticProcess); } } } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void myDetPhysicsList::SetCuts() { //G4VUserPhysicsList::SetCutsWithDefault method sets //the default cut value for all particle types // SetCutsWithDefault(); if (verboseLevel>0) DumpCutValuesTable(); } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......