// // ******************************************************************** // * License and Disclaimer * // * * // * The Geant4 software is copyright of the Copyright Holders of * // * the Geant4 Collaboration. It is provided under the terms and * // * conditions of the Geant4 Software License, included in the file * // * LICENSE and available at http://cern.ch/geant4/license . These * // * include a list of copyright holders. * // * * // * Neither the authors of this software system, nor their employing * // * institutes,nor the agencies providing financial support for this * // * work make any representation or warranty, express or implied, * // * regarding this software system or assume any liability for its * // * use. Please see the license in the file LICENSE and URL above * // * for the full disclaimer and the limitation of liability. * // * * // * This code implementation is the result of the scientific and * // * technical work of the GEANT4 collaboration. * // * By using, copying, modifying or distributing the software (or * // * any work based on the software) you agree to acknowledge its * // * use in resulting scientific publications, and indicate your * // * acceptance of all terms of the Geant4 Software license. * // ******************************************************************** // // Taken from /example/Extended/OpNovice // gpaterno: add the fIWantOpticalPhotons variable and Set method // gpaterno: add radioactivedecay (see exaples/extended/radioactivedecay/rdecay01) // gpaterno: add setCuts method for the messenger // gpaterno: add Atomic deexcitation processes // gpaterno: add Modular PhysicsLists for Hadronic Processes // gpaterno: substituted custom EM PhysicsLists with Modular ones // gpaterno: added commands to select PE or scattering (Penelope) only // //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... #include "PhysicsList.hh" #include "PhysicsListMessenger.hh" #include "G4UnitsTable.hh" #include "G4SystemOfUnits.hh" #include "globals.hh" #include "G4ProcessManager.hh" #include "G4Threading.hh" //Particles #include "G4ParticleDefinition.hh" #include "G4ParticleTypes.hh" #include "G4ParticleTable.hh" //EM Physics Lists #include "G4EmStandardPhysics.hh" #include "G4EmLivermorePhysics.hh" #include "G4EmPenelopePhysics.hh" #include "G4EmPenelopePhysicsMI.hh" #include "G4EmLowEPPhysics.hh" #include "G4EmStandardPhysics_option4.hh" //EM options #include "G4EmSaturation.hh" #include "G4LossTableManager.hh" #include "G4UAtomicDeexcitation.hh" //Hadronic and Extra Physics Lists #include "G4EmExtraPhysics.hh" #include "G4HadronPhysicsQGSP_BIC_HP.hh" #include "G4HadronElasticPhysicsHP.hh" #include "G4IonPhysics.hh" #include "G4StoppingPhysics.hh" //Optical processes #include "G4Cerenkov.hh" #include "G4Scintillation.hh" #include "G4OpAbsorption.hh" #include "G4OpRayleigh.hh" #include "G4OpMieHG.hh" #include "G4OpBoundaryProcess.hh" //Decays #include "G4Decay.hh" #include "G4DecayPhysics.hh" #include "G4PhysicsListHelper.hh" #include "G4Radioactivation.hh" #include "G4NuclideTable.hh" #include "G4NuclearLevelData.hh" #include "G4DeexPrecoParameters.hh" #include "G4PhysListUtil.hh" //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... G4ThreadLocal G4Cerenkov* PhysicsList::fCerenkovProcess = 0; G4ThreadLocal G4Scintillation* PhysicsList::fScintillationProcess = 0; G4ThreadLocal G4OpAbsorption* PhysicsList::fAbsorptionProcess = 0; G4ThreadLocal G4OpRayleigh* PhysicsList::fRayleighScatteringProcess = 0; G4ThreadLocal G4OpMieHG* PhysicsList::fMieHGScatteringProcess = 0; G4ThreadLocal G4OpBoundaryProcess* PhysicsList::fBoundaryProcess = 0; //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... PhysicsList::PhysicsList(): G4VUserPhysicsList(), fPMessenger(0) { G4cout << "### PhysicsList instantiated ###" << G4endl; G4LossTableManager::Instance(); //Hadronic and Extra Physics option fIWantHadExtraPhysics = true; //optical photon options fIWantOpticalPhotons = false; fOPVerboseLevel = 0; fMaxNumPhotonStep = 30; fIWantCerernkov = false; //set default cuts value defaultCutValue = 0.001*mm; //define the messenger fPMessenger = new PhysicsListMessenger(this); //set verbosity SetVerboseLevel(1); //particle list fParticleList = new G4DecayPhysics(verboseLevel); //EM Physics fEmPhysicsList = new G4EmStandardPhysics_option4(verboseLevel); //Hadronic and Extra Physics fEmExtraPhysicsList = new G4EmExtraPhysics(verboseLevel); fHadPhysicsList = new G4HadronPhysicsQGSP_BIC_HP(verboseLevel); fHadElPhysicsList = new G4HadronElasticPhysicsHP(verboseLevel); fIonPhysicsList = new G4IonPhysics(verboseLevel); fStoppingPhysicsList = new G4StoppingPhysics(verboseLevel); } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... PhysicsList::~PhysicsList() {} //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void PhysicsList::ConstructParticle() { fParticleList->ConstructParticle(); } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void PhysicsList::ConstructProcess() { //transportation AddTransportation(); //EM physics fEmPhysicsList->ConstructProcess(); //Hadronic and Extra Physics if (fIWantHadExtraPhysics) { fEmExtraPhysicsList->ConstructProcess(); fHadPhysicsList->ConstructProcess(); fHadElPhysicsList->ConstructProcess(); fIonPhysicsList->ConstructProcess(); fStoppingPhysicsList->ConstructProcess(); } //Atomic deexcitation G4VAtomDeexcitation* de = new G4UAtomicDeexcitation(); de->SetFluo(true); //To activate deexcitation processes and fluorescence de->SetAuger(true); //To activate Auger effect if deexcitation is activated de->SetAugerCascade(true); //To activate Auger Cascade if deexcitation is activated de->SetPIXE(true); //To activate Particle Induced X-Ray Emission (PIXE) G4LossTableManager::Instance()->SetAtomDeexcitation(de); //optical physics activation if (fIWantOpticalPhotons) { ConstructOp(); G4cout << "### Optical Physics activated ###" << G4endl; } //instantiate Physics List infrastructure G4PhysListUtil::InitialiseParameters(); //update G4NuclideTable time limit const G4double meanLife = 0.1*picosecond; G4NuclideTable::GetInstance()->SetMeanLifeThreshold(meanLife); G4NuclideTable::GetInstance()->SetLevelTolerance(1.0*eV); //define flags for nuclear gamma de-excitation model G4DeexPrecoParameters* deex = G4NuclearLevelData::GetInstance()->GetParameters(); deex->SetCorrelatedGamma(false); deex->SetStoreAllLevels(true); deex->SetInternalConversionFlag(true); deex->SetIsomerProduction(true); deex->SetMaxLifeTime(meanLife); //decay ConstructDecay(); } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void PhysicsList::ConstructDecay() { G4cout << G4endl << "### Decay processes activated ###" << G4endl; G4Decay* theDecayProcess = new G4Decay(); GetParticleIterator()->reset(); while ( (*GetParticleIterator())() ) { G4ParticleDefinition* particle = GetParticleIterator()->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); } } G4Radioactivation* radioactiveDecay = new G4Radioactivation("Radioactivation", 1.0e+60*CLHEP::year); radioactiveDecay->SetARM(true); //Atomic Rearangement G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper(); ph->RegisterProcess(radioactiveDecay, G4GenericIon::GenericIon()); } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void PhysicsList::ConstructOp() { if (fIWantCerernkov) { fCerenkovProcess = new G4Cerenkov("Cerenkov"); fCerenkovProcess->SetMaxNumPhotonsPerStep(fMaxNumPhotonStep); fCerenkovProcess->SetMaxBetaChangePerStep(10.0); fCerenkovProcess->SetTrackSecondariesFirst(true); fCerenkovProcess->SetVerboseLevel(fOPVerboseLevel); G4cout << "### Cerenkov Process activated ###" << G4endl; } fScintillationProcess = new G4Scintillation("Scintillation"); //fScintillationProcess->SetScintillationYieldFactor(1.); fScintillationProcess->SetTrackSecondariesFirst(true); fAbsorptionProcess = new G4OpAbsorption(); fRayleighScatteringProcess = new G4OpRayleigh(); fMieHGScatteringProcess = new G4OpMieHG(); fBoundaryProcess = new G4OpBoundaryProcess(); fScintillationProcess->SetVerboseLevel(fOPVerboseLevel); fAbsorptionProcess->SetVerboseLevel(fOPVerboseLevel); fRayleighScatteringProcess->SetVerboseLevel(fOPVerboseLevel); fMieHGScatteringProcess->SetVerboseLevel(fOPVerboseLevel); fBoundaryProcess->SetVerboseLevel(fOPVerboseLevel); /* //Use Birks Correction in the Scintillation process if(G4Threading::IsMasterThread()) { G4EmSaturation* emSaturation = G4LossTableManager::Instance()->EmSaturation(); fScintillationProcess->AddSaturation(emSaturation); } */ GetParticleIterator()->reset(); while ((*GetParticleIterator())()) { G4ParticleDefinition* particle = GetParticleIterator()->value(); G4ProcessManager* pmanager = particle->GetProcessManager(); G4String particleName = particle->GetParticleName(); if (fIWantCerernkov && fCerenkovProcess->IsApplicable(*particle)) { pmanager->AddProcess(fCerenkovProcess); pmanager->SetProcessOrdering(fCerenkovProcess,idxPostStep); } if (fScintillationProcess->IsApplicable(*particle)) { pmanager->AddProcess(fScintillationProcess); pmanager->SetProcessOrderingToLast(fScintillationProcess, idxAtRest); pmanager->SetProcessOrderingToLast(fScintillationProcess, idxPostStep); } if (particleName == "opticalphoton") { G4cout << "### AddDiscreteProcess to OpticalPhoton ###" << G4endl; pmanager->AddDiscreteProcess(fAbsorptionProcess); pmanager->AddDiscreteProcess(fRayleighScatteringProcess); pmanager->AddDiscreteProcess(fMieHGScatteringProcess); pmanager->AddDiscreteProcess(fBoundaryProcess); } } } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void PhysicsList::SelectPhysicsList(const G4String& name) { if (verboseLevel > 1) { G4cout << "### PhysicsList::SelectPhysicsList: <" << name << "> ###" << G4endl; } if (name == "standard") { delete fEmPhysicsList; fEmPhysicsList = new G4EmStandardPhysics(verboseLevel); G4cout << "### selected Standard PhysicsList ###" << G4endl; } else if (name == "livermore") { delete fEmPhysicsList; fEmPhysicsList = new G4EmLivermorePhysics(verboseLevel); G4cout << "### selected Livermore PhysicsList ###" << G4endl; } else if (name == "penelope") { delete fEmPhysicsList; fEmPhysicsList = new G4EmPenelopePhysics(verboseLevel); G4cout << "### selected Penelope PhysicsList ###" << G4endl; } else if (name == "penelopeMI") { delete fEmPhysicsList; fEmPhysicsList = new G4EmPenelopePhysicsMI(verboseLevel); G4cout << "### selected PenelopeMI PhysicsList ###" << G4endl; } else if (name == "penelopeMI_PE") { delete fEmPhysicsList; fEmPhysicsList = new G4EmPenelopePhysicsMI(verboseLevel, "G4EmPenelopeMI", false, true, false); G4cout << "### selected PenelopeMI PhysicsList with PE effect only ###" << G4endl; } else if (name == "penelopeMI_scatt") { delete fEmPhysicsList; fEmPhysicsList = new G4EmPenelopePhysicsMI(verboseLevel, "G4EmPenelopeMI", true, false, true); G4cout << "### selected PenelopeMI PhysicsList with scattering only ###" << G4endl; } else if (name == "LowEP") { delete fEmPhysicsList; fEmPhysicsList = new G4EmLowEPPhysics(1,name); G4cout << "### selected LowEP PhysicsList ###" << G4endl; } else { G4cout << "### PhysicsList::SelectPhysicsList: <" << name << ">"<< " is not defined ###" << G4endl; } } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void PhysicsList::SetVerbose(G4int verbose) { fOPVerboseLevel = verbose; } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void PhysicsList::SetNbOfPhotonsCerenkov(G4int MaxNumber) { fMaxNumPhotonStep = MaxNumber; } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void PhysicsList::SetCuts() { //set the default cut value for all particle types SetCutsWithDefault(); if (verboseLevel>0) DumpCutValuesTable(); } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void PhysicsList::SetOpPhotonPhysics(G4int val) { if (val == 1) fIWantOpticalPhotons = true; else if (val == 0) fIWantOpticalPhotons = false; } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void PhysicsList::SetDefaultCutsValue(G4double value) { defaultCutValue = value; } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... void PhysicsList::SetCerenkovProcess(G4int val) { if (val == 1) fIWantCerernkov = true; else if (val == 0) fIWantCerernkov = false; } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......