STCyclotron/src/STCyclotronPhysicsList.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 // Author: F. Poignant, floriane.poignant@gmail.com
0027 //
0028
0029 #include "STCyclotronPhysicsList.hh"
0030 #include "G4Proton.hh"
0031 #include "G4Neutron.hh"
0032 #include "G4Electron.hh"
0033 #include "G4Positron.hh"
0034 #include "G4Gamma.hh"
0035 #include "G4Deuteron.hh"
0036
0037 //Physic Lists (contained inside the Geant 4 distribution)
0038 #include "G4EmStandardPhysics_option3.hh"
0039 #include "G4DecayPhysics.hh"
0040 #include "G4Decay.hh"
0041 #include "G4StepLimiter.hh"
0042 #include "G4LossTableManager.hh"
0043 #include "G4UnitsTable.hh"
0044 #include "G4SystemOfUnits.hh"
0045 #include "G4ProcessManager.hh"
0046 #include "G4Region.hh"
0047 #include "G4RegionStore.hh"
0048
0049 #include "G4EmExtraPhysics.hh"
0050 #include "G4EmParameters.hh"
0051 #include "G4NuclideTable.hh"
0052
0053 #include "G4HadronPhysicsQGSP_BERT.hh"
0054 #include "G4HadronPhysicsQGSP_BIC.hh"
0055 #include "G4HadronPhysicsQGSP_BIC_HP.hh"
0056 #include "G4HadronPhysicsQGSP_BIC_AllHP.hh"
0057 #include "G4RadioactiveDecayPhysics.hh"
0058 //#include "QGSP_BIC_HP.hh"
0059 #include "G4PhysListFactory.hh"
0060 #include "G4DeexPrecoParameters.hh"
0061 #include "G4NuclideTable.hh"
0062
0063 STCyclotronPhysicsList::STCyclotronPhysicsList(STCyclotronDetectorConstruction* det)
0064   : G4VModularPhysicsList()
0065 {
0066   // Mandatory for G4NuclideTable
0067   // Half-life threshold must be set small or many short-lived isomers
0068   // will not be assigned life times (default to 0)
0069   G4NuclideTable::GetInstance()->SetThresholdOfHalfLife(0.1*second);
0070   G4NuclideTable::GetInstance()->SetLevelTolerance(1.0*eV);
0071
0072   //---
0073   fDetector = det;
0074
0075   G4LossTableManager::Instance();
0076
0077   defaultCutValue = 0.1*mm;
0078   fCutForGamma     = defaultCutValue;
0079   fCutForElectron  = defaultCutValue;
0080   fCutForPositron  = defaultCutValue;
0081
0082   fThickness_foil = defaultCutValue;
0083   fThickness_target = defaultCutValue;
0084
0085   fCutTargetProton  = fThickness_target;
0086   fCutTargetElectron = fThickness_target;
0087   fCutTargetPositron = fThickness_target;
0088   fCutTargetGamma = fThickness_target;
0089
0090   fCutFoilProton = fThickness_foil;
0091   fCutFoilElectron = fThickness_foil;
0092   fCutFoilPositron = fThickness_foil;
0093   fCutFoilGamma = fThickness_foil;
0094
0095   //EM physics
0096   fEmPhysicsList = new G4EmStandardPhysics_option3(0);
0097   fEmName = G4String("emstandard_opt3");
0098
0099   //Decay physics and all particles
0100   fDecPhysicsList = new G4DecayPhysics(0);
0101   fRaddecayList = new G4RadioactiveDecayPhysics(0);
0102
0103   //Hadron physics
0104   fHadPhysicsList = new G4HadronPhysicsQGSP_BIC_AllHP(0);
0105   //fHadPhysicsList = new G4HadronPhysicsQGSP_BIC(0);
0106
0107
0108
0109
0110
0111 }
0112
0113 STCyclotronPhysicsList::~STCyclotronPhysicsList()
0114 {
0115
0116   delete fEmPhysicsList ;
0117   delete fDecPhysicsList;
0118   delete fRaddecayList;
0119   delete fHadPhysicsList;
0120
0121
0122 }
0123
0124 void STCyclotronPhysicsList::ConstructParticle()
0125 {
0126
0127   G4Proton::ProtonDefinition();
0128   G4Gamma::GammaDefinition();
0129   G4Electron::ElectronDefinition();
0130   G4Positron::PositronDefinition();
0131   G4Neutron::NeutronDefinition();
0132   G4Deuteron::DeuteronDefinition();
0133
0134   fDecPhysicsList->ConstructParticle();
0135
0136
0137 }
0138
0139 void STCyclotronPhysicsList::ConstructProcess()
0140 {
0141   // Define transportation process
0142
0143
0144   AddTransportation();
0145
0146   //electromagnetic physics list
0147   fEmPhysicsList->ConstructProcess();
0148   //em_config.AddModels();
0149
0150   //decay physics list
0151   fDecPhysicsList->ConstructProcess();
0152   fRaddecayList->ConstructProcess();
0153
0154   //hadronic physics lists
0155   fHadPhysicsList->ConstructProcess();
0156
0157   //Get the value of the fThickness of foil and target
0158   fThickness_foil = fDetector->GetFoilThickness()*mm;
0159   fThickness_target = fDetector->GetTargetThickness()*mm;
0160
0161   //Update the cuts with the 1/2 of the thickness of the foil/target
0162   //SetCuts();
0163
0164   SetCutTarget(0.01,fThickness_target/2.,fThickness_target/2.,fThickness_target/2.);
0165
0166
0167
0168 }
0169
0170 void STCyclotronPhysicsList::SetCuts()
0171 {
0172
0173   if (verboseLevel >0){
0174     G4cout << "PhysicsList::SetCuts:";
0175     G4cout << "CutLength : " << G4BestUnit(defaultCutValue,"Length") << G4endl;
0176   }
0177
0178   // set cut values for gamma at first and for e- second and next for e+,
0179   // because some processes for e+/e- need cut values for gamma
0180   SetCutValue(fCutForGamma, "gamma");
0181   SetCutValue(fCutForElectron, "e-");
0182   SetCutValue(fCutForPositron, "e+");
0183
0184   // Set cuts for detector
0185   SetCutFoil(0.1,fThickness_foil/2.,fThickness_foil/2.,fThickness_foil/2.);
0186   SetCutTarget(0.01,fThickness_target/2.,fThickness_target/2.,fThickness_target/2.);
0187   if (verboseLevel>0) DumpCutValuesTable();
0188 }
0189
0190 void STCyclotronPhysicsList::SetCutForGamma(G4double cut)
0191 {
0192   fCutForGamma = cut;
0193   SetParticleCuts(fCutForGamma, G4Gamma::Gamma());
0194 }
0195
0196 void STCyclotronPhysicsList::SetCutForElectron(G4double cut)
0197 {
0198   fCutForElectron = cut;
0199   SetParticleCuts(fCutForElectron, G4Electron::Electron());
0200 }
0201
0202 void STCyclotronPhysicsList::SetCutForPositron(G4double cut)
0203 {
0204   fCutForPositron = cut;
0205   SetParticleCuts(fCutForPositron, G4Positron::Positron());
0206 }
0207
0208 void STCyclotronPhysicsList::SetCutTarget(G4double cutProton, G4double cutElectron, G4double cutPositron, G4double cutGamma){
0209
0210   fCutTargetProton = cutProton*mm;
0211   fCutTargetElectron = cutElectron*mm;
0212   fCutTargetPositron = cutPositron*mm;
0213   fCutTargetGamma = cutGamma*mm;
0214
0215   G4String regionNameTarget = "Target";
0216   G4Region* regionTarget = G4RegionStore::GetInstance()->GetRegion(regionNameTarget);
0217
0218   G4ProductionCuts* cutsTarget = new G4ProductionCuts ;
0219   cutsTarget -> SetProductionCut(fCutTargetGamma,"gamma");
0220   cutsTarget -> SetProductionCut(fCutTargetElectron,"e-");
0221   cutsTarget -> SetProductionCut(fCutTargetPositron,"e+");
0222   cutsTarget -> SetProductionCut(fCutTargetProton,"proton");
0223
0224   regionTarget -> SetProductionCuts(cutsTarget);
0225
0226 }
0227
0228 void STCyclotronPhysicsList::SetCutFoil(G4double cutProton, G4double cutElectron, G4double cutPositron, G4double cutGamma){
0229
0230   fCutFoilProton = cutProton*mm;
0231   fCutFoilElectron = cutElectron*mm;
0232   fCutFoilPositron = cutPositron*mm;
0233   fCutFoilGamma = cutGamma*mm;
0234
0235   G4RegionStore::GetInstance()->GetRegion("Foil");
0236   G4String regionNameFoil = "Foil";
0237   G4Region* regionFoil = G4RegionStore::GetInstance()->GetRegion(regionNameFoil);
0238
0239   G4ProductionCuts* cutsFoil = new G4ProductionCuts ;
0240   cutsFoil -> SetProductionCut(fCutFoilGamma,"gamma");
0241   cutsFoil -> SetProductionCut(fCutFoilElectron,"e-");
0242   cutsFoil -> SetProductionCut(fCutFoilPositron,"e+");
0243   cutsFoil -> SetProductionCut(fCutFoilProton,"proton");
0244
0245   regionFoil -> SetProductionCuts(cutsFoil);
0246
0247 }
0248