gammaray_telescope/src/GammaRayTelPhysicsList.cc

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 // Authors: Francesco Longo, franzlongo1969@gmail.com
0027 //
0028 // Code based on the hadrontherapy && radioprotection advanced example
0029
0030 #include "GammaRayTelPhysicsList.hh"
0031 #include "GammaRayTelPhysicsListMessenger.hh"
0032 #include "G4PhysListFactory.hh"
0033 #include "G4VPhysicsConstructor.hh"
0034
0035 // Physics lists (contained inside the Geant4 distribution)
0036 #include "G4EmLivermorePhysics.hh"
0037 #include "G4EmLivermorePolarizedPhysics.hh"
0038 #include "G4EmPenelopePhysics.hh"
0039 #include "G4EmStandardPhysics_option3.hh"
0040 #include "G4EmStandardPhysics_option4.hh" // to treat the new polarized process
0041
0042 #include "G4Decay.hh"
0043 #include "G4DecayPhysics.hh"
0044 #include "G4HadronDElasticPhysics.hh"
0045 #include "G4HadronElasticPhysics.hh"
0046 #include "G4HadronElasticPhysicsHP.hh"
0047 #include "G4HadronPhysicsQGSP_BIC_HP.hh"
0048 #include "G4IonBinaryCascadePhysics.hh"
0049 #include "G4IonFluctuations.hh"
0050 #include "G4IonParametrisedLossModel.hh"
0051 #include "G4LossTableManager.hh"
0052 #include "G4ProcessManager.hh"
0053 #include "G4RadioactiveDecayPhysics.hh"
0054 #include "G4StepLimiter.hh"
0055 #include "G4SystemOfUnits.hh"
0056 #include "G4UnitsTable.hh"
0057
0058 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0059
0060 GammaRayTelPhysicsList::GammaRayTelPhysicsList() {
0061     G4LossTableManager::Instance();
0062
0063     constexpr auto DEFAULT_CUT_VALUE{100 * micrometer};
0064     defaultCutValue = DEFAULT_CUT_VALUE;
0065
0066     SetVerboseLevel(1);
0067
0068     constexpr auto ENERGY_LOWER_BOUND{250 * eV};
0069     constexpr auto ENERGY_UPPER_BOUND{1 * GeV};
0070
0071     G4ProductionCutsTable::GetProductionCutsTable()->SetEnergyRange(ENERGY_LOWER_BOUND, ENERGY_UPPER_BOUND);
0072     SetDefaultCutValue (defaultCutValue);
0073     DumpCutValuesTable();
0074
0075     helIsRegisted = false;
0076     bicIsRegisted = false;
0077     biciIsRegisted = false;
0078     locIonIonInelasticIsRegistered = false;
0079     radioactiveDecayIsRegisted = false;
0080
0081     pMessenger = new GammaRayTelPhysicsListMessenger(this);
0082
0083     SetVerboseLevel(1);
0084
0085     // EM physics
0086     emPhysicsList = new G4EmStandardPhysics_option3(1);
0087     emName = G4String("emstandard_opt3");
0088
0089     // Decay physics and all particles
0090     decPhysicsList = new G4DecayPhysics();
0091 }
0092
0093 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0094
0095 GammaRayTelPhysicsList::~GammaRayTelPhysicsList() {
0096     delete pMessenger;
0097     delete emPhysicsList;
0098     delete decPhysicsList;
0099
0100     for (auto &hadronPhy : hadronPhys) {
0101         delete hadronPhy;
0102     }
0103 }
0104
0105 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0106
0107 void GammaRayTelPhysicsList::AddPackage(const G4String &name) {
0108     G4PhysListFactory factory;
0109     auto *phys = factory.GetReferencePhysList(name);
0110
0111     G4int i = 0;
0112     const G4VPhysicsConstructor *element = phys->GetPhysics(i);
0113     auto *tmp = const_cast<G4VPhysicsConstructor*>(element);
0114
0115     while (element != nullptr) {
0116         RegisterPhysics(tmp);
0117         element = phys->GetPhysics(++i);
0118         tmp = const_cast<G4VPhysicsConstructor*>(element);
0119     }
0120 }
0121
0122 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0123
0124 void GammaRayTelPhysicsList::ConstructParticle() {
0125     decPhysicsList->ConstructParticle();
0126 }
0127
0128 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0129
0130 void GammaRayTelPhysicsList::ConstructProcess() {
0131     // transportation
0132     //
0133     AddTransportation();
0134
0135     // electromagnetic physics list
0136     //
0137     emPhysicsList->ConstructProcess();
0138     emConfigurator.AddModels();
0139
0140     // decay physics list
0141     //
0142     decPhysicsList->ConstructProcess();
0143
0144     // hadronic physics lists
0145     for (auto &hadronPhy : hadronPhys) {
0146         hadronPhy->ConstructProcess();
0147     }
0148
0149     // step limitation (as a full process)
0150     // AddStepMax();
0151 }
0152
0153 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0154
0155 void GammaRayTelPhysicsList::AddPhysicsList(const G4String &name) {
0156     if (verboseLevel > 1) {
0157         G4cout << "PhysicsList::AddPhysicsList: <" << name << ">" << G4endl;
0158     }
0159     if (name == emName) {
0160         return;
0161     }
0162
0163     ////////////////////////////////
0164     //   ELECTROMAGNETIC MODELS   //
0165     ////////////////////////////////
0166
0167     if (name == "standard_opt3") {
0168         emName = name;
0169         delete emPhysicsList;
0170         emPhysicsList = new G4EmStandardPhysics_option3();
0171         G4cout << "THE FOLLOWING ELECTROMAGNETIC PHYSICS LIST HAS BEEN ACTIVATED: G4EmStandardPhysics_option3" << G4endl;
0172     } else if (name == "LowE_Livermore") {
0173         emName = name;
0174         delete emPhysicsList;
0175         emPhysicsList = new G4EmLivermorePhysics();
0176         G4cout << "THE FOLLOWING ELECTROMAGNETIC PHYSICS LIST HAS BEEN ACTIVATED: G4EmLivermorePhysics" << G4endl;
0177     } else if (name == "LowE_Penelope") {
0178         emName = name;
0179         delete emPhysicsList;
0180         emPhysicsList = new G4EmPenelopePhysics();
0181         G4cout << "THE FOLLOWING ELECTROMAGNETIC PHYSICS LIST HAS BEEN ACTIVATED: G4EmLivermorePhysics" << G4endl;
0182     } else if (name == "LowE_Polarized") {
0183         emName = name;
0184         delete emPhysicsList;
0185         emPhysicsList = new G4EmLivermorePolarizedPhysics();
0186         G4cout << "THE FOLLOWING ELECTROMAGNETIC PHYSICS LIST HAS BEEN ACTIVATED: G4EmLivermorePhysics" << G4endl;
0187     } else if (name == "standard_opt4") {
0188         emName = name;
0189         delete emPhysicsList;
0190         emPhysicsList = new G4EmStandardPhysics_option4();
0191         G4cout << "THE FOLLOWING ELECTROMAGNETIC PHYSICS LIST HAS BEEN ACTIVATED: G4EmStandardOption_4" << G4endl;
0192
0193         /////////////////////////
0194         //   HADRONIC MODELS   //
0195         /////////////////////////
0196
0197     } else if (name == "elastic" && !helIsRegisted) {
0198         G4cout << "THE FOLLOWING HADRONIC ELASTIC PHYSICS LIST HAS BEEN ACTIVATED: G4HadronElasticPhysics" << G4endl;
0199         hadronPhys.push_back(new G4HadronElasticPhysics());
0200         helIsRegisted = true;
0201     } else if (name == "DElastic" && !helIsRegisted) {
0202         hadronPhys.push_back(new G4HadronDElasticPhysics());
0203         helIsRegisted = true;
0204     } else if (name == "HPElastic" && !helIsRegisted) {
0205         hadronPhys.push_back(new G4HadronElasticPhysicsHP());
0206         helIsRegisted = true;
0207     } else if (name == "binary" && !bicIsRegisted) {
0208         hadronPhys.push_back(new G4HadronPhysicsQGSP_BIC_HP());
0209         bicIsRegisted = true;
0210         G4cout << "THE FOLLOWING HADRONIC INELASTIC PHYSICS LIST HAS BEEN ACTIVATED: HadronPhysicsQGSP_BIC_HP" << G4endl;
0211     } else if (name == "binary_ion" && !biciIsRegisted) {
0212         hadronPhys.push_back(new G4IonBinaryCascadePhysics());
0213         biciIsRegisted = true;
0214         G4cout << "THE FOLLOWING HADRONIC INELASTIC PHYSICS LIST HAS BEEN ACTIVATED: G4IonBinaryCascadePhysics" << G4endl;
0215     } else if (name == "radioactive_decay" && !radioactiveDecayIsRegisted) {
0216         hadronPhys.push_back(new G4RadioactiveDecayPhysics());
0217         radioactiveDecayIsRegisted = true;
0218         G4cout << "THE FOLLOWING HADRONIC INELASTIC PHYSICS LIST HAS BEEN ACTIVATED: G4RadioactiveDecayPhysics" << G4endl;
0219     } else {
0220         G4cout << "PhysicsList::AddPhysicsList: <" << name << "> is not defined" << G4endl;
0221     }
0222 }
0223
0224 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0225
0226 void GammaRayTelPhysicsList::SetCutForGamma(G4double cut) {
0227     SetParticleCuts(cut, G4Gamma::Gamma());
0228 }
0229
0230 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0231
0232 void GammaRayTelPhysicsList::SetCutForElectron(G4double cut) {
0233     SetParticleCuts(cut, G4Electron::Electron());
0234 }
0235
0236 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0237
0238 void GammaRayTelPhysicsList::SetCutForPositron(G4double cut) {
0239     SetParticleCuts(cut, G4Positron::Positron());
0240 }