TestEm1/src/PhysicsList.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 /// \file electromagnetic/TestEm1/src/PhysicsList.cc
0027 /// \brief Implementation of the PhysicsList class
0028 //
0029 //
0030 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0031 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0032
0033 #include "PhysicsList.hh"
0034
0035 #include "DetectorConstruction.hh"
0036 #include "PhysListEmStandard.hh"
0037 #include "PhysicsListMessenger.hh"
0038
0039 #include "G4EmLivermorePhysics.hh"
0040 #include "G4EmLowEPPhysics.hh"
0041 #include "G4EmParameters.hh"
0042 #include "G4EmPenelopePhysics.hh"
0043 #include "G4EmStandardPhysics.hh"
0044 #include "G4EmStandardPhysicsGS.hh"
0045 #include "G4EmStandardPhysicsSS.hh"
0046 #include "G4EmStandardPhysicsWVI.hh"
0047 #include "G4EmStandardPhysics_option1.hh"
0048 #include "G4EmStandardPhysics_option2.hh"
0049 #include "G4EmStandardPhysics_option3.hh"
0050 #include "G4EmStandardPhysics_option4.hh"
0051 #include "G4LossTableManager.hh"
0052 #include "G4SystemOfUnits.hh"
0053 #include "G4UnitsTable.hh"
0054
0055 // particles
0056
0057 #include "StepMax.hh"
0058
0059 #include "G4BaryonConstructor.hh"
0060 #include "G4BosonConstructor.hh"
0061 #include "G4Decay.hh"
0062 #include "G4Electron.hh"
0063 #include "G4Gamma.hh"
0064 #include "G4GenericIon.hh"
0065 #include "G4IonConstructor.hh"
0066 #include "G4LeptonConstructor.hh"
0067 #include "G4Material.hh"
0068 #include "G4MesonConstructor.hh"
0069 #include "G4NuclideTable.hh"
0070 #include "G4PhysicsListHelper.hh"
0071 #include "G4ProcessManager.hh"
0072 #include "G4Proton.hh"
0073 #include "G4RadioactiveDecay.hh"
0074 #include "G4ShortLivedConstructor.hh"
0075
0076 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0077
0078 PhysicsList::PhysicsList(DetectorConstruction* det) : fDet(det)
0079 {
0080   fMessenger = new PhysicsListMessenger(this);
0081   SetVerboseLevel(1);
0082
0083   // EM physics
0084   AddPhysicsList("emstandard_opt3");
0085
0086   // fix lower limit for cut
0087   G4ProductionCutsTable::GetProductionCutsTable()->SetEnergyRange(10 * eV, 1 * GeV);
0088   SetDefaultCutValue(1 * mm);
0089 }
0090
0091 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0092
0093 PhysicsList::~PhysicsList()
0094 {
0095   delete fMessenger;
0096 }
0097
0098 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0099
0100 void PhysicsList::ConstructParticle()
0101 {
0102   G4BosonConstructor pBosonConstructor;
0103   pBosonConstructor.ConstructParticle();
0104
0105   G4LeptonConstructor pLeptonConstructor;
0106   pLeptonConstructor.ConstructParticle();
0107
0108   G4MesonConstructor pMesonConstructor;
0109   pMesonConstructor.ConstructParticle();
0110
0111   G4BaryonConstructor pBaryonConstructor;
0112   pBaryonConstructor.ConstructParticle();
0113
0114   G4IonConstructor pIonConstructor;
0115   pIonConstructor.ConstructParticle();
0116
0117   G4ShortLivedConstructor pShortLivedConstructor;
0118   pShortLivedConstructor.ConstructParticle();
0119 }
0120
0121 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0122
0123 void PhysicsList::ConstructProcess()
0124 {
0125   // Transportation
0126   //
0127   AddTransportation();
0128
0129   // Electromagnetic physics list
0130   //
0131   fEmPhysicsList->ConstructProcess();
0132
0133   // Decay Process
0134   //
0135   AddDecay();
0136
0137   // Decay Process
0138   //
0139   AddRadioactiveDecay();
0140
0141   // step limitation (as a full process)
0142   //
0143   AddStepMax();
0144
0145   // example of Get process
0146   auto process = GetProcess("RadioactiveDecay");
0147   if (process != nullptr) {
0148     G4cout << "\n  GetProcess : " << process->GetProcessName() << G4endl;
0149   }
0150 }
0151
0152 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0153
0154 void PhysicsList::AddPhysicsList(const G4String& name)
0155 {
0156   if (verboseLevel > 0) {
0157     G4cout << "PhysicsList::AddPhysicsList: <" << name << ">" << G4endl;
0158   }
0159
0160   if (name == fEmName) return;
0161
0162   if (name == "local") {
0163     fEmName = name;
0164     delete fEmPhysicsList;
0165     fEmPhysicsList = new PhysListEmStandard(name);
0166   }
0167   else if (name == "emstandard_opt0") {
0168     fEmName = name;
0169     delete fEmPhysicsList;
0170     fEmPhysicsList = new G4EmStandardPhysics();
0171   }
0172   else if (name == "emstandard_opt1") {
0173     fEmName = name;
0174     delete fEmPhysicsList;
0175     fEmPhysicsList = new G4EmStandardPhysics_option1();
0176   }
0177   else if (name == "emstandard_opt2") {
0178     fEmName = name;
0179     delete fEmPhysicsList;
0180     fEmPhysicsList = new G4EmStandardPhysics_option2();
0181   }
0182   else if (name == "emstandard_opt3") {
0183     fEmName = name;
0184     delete fEmPhysicsList;
0185     fEmPhysicsList = new G4EmStandardPhysics_option3();
0186   }
0187   else if (name == "emstandard_opt4") {
0188     fEmName = name;
0189     delete fEmPhysicsList;
0190     fEmPhysicsList = new G4EmStandardPhysics_option4();
0191   }
0192   else if (name == "emstandardSS") {
0193     fEmName = name;
0194     delete fEmPhysicsList;
0195     fEmPhysicsList = new G4EmStandardPhysicsSS();
0196   }
0197   else if (name == "emstandardGS") {
0198     fEmName = name;
0199     delete fEmPhysicsList;
0200     fEmPhysicsList = new G4EmStandardPhysicsGS();
0201   }
0202   else if (name == "emstandardWVI") {
0203     fEmName = name;
0204     delete fEmPhysicsList;
0205     fEmPhysicsList = new G4EmStandardPhysicsWVI();
0206   }
0207   else if (name == "emlivermore") {
0208     fEmName = name;
0209     delete fEmPhysicsList;
0210     fEmPhysicsList = new G4EmLivermorePhysics();
0211   }
0212   else if (name == "empenelope") {
0213     fEmName = name;
0214     delete fEmPhysicsList;
0215     fEmPhysicsList = new G4EmPenelopePhysics();
0216   }
0217   else if (name == "emlowenergy") {
0218     fEmName = name;
0219     delete fEmPhysicsList;
0220     fEmPhysicsList = new G4EmLowEPPhysics();
0221   }
0222   else {
0223     G4cout << "PhysicsList::AddPhysicsList: <" << name << ">"
0224            << " is not defined" << G4endl;
0225   }
0226
0227   // Em options
0228   //
0229   G4EmParameters::Instance()->SetBuildCSDARange(true);
0230   G4EmParameters::Instance()->SetGeneralProcessActive(false);
0231 }
0232
0233 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0234
0235 void PhysicsList::AddDecay()
0236 {
0237   G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
0238
0239   // Decay Process
0240   //
0241   G4Decay* fDecayProcess = new G4Decay();
0242
0243   auto particleIterator = GetParticleIterator();
0244   particleIterator->reset();
0245   while ((*particleIterator)()) {
0246     G4ParticleDefinition* particle = particleIterator->value();
0247     if (fDecayProcess->IsApplicable(*particle) && !particle->IsShortLived())
0248       ph->RegisterProcess(fDecayProcess, particle);
0249   }
0250 }
0251
0252 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0253
0254 void PhysicsList::AddRadioactiveDecay()
0255 {
0256   G4RadioactiveDecay* radioactiveDecay = new G4RadioactiveDecay();
0257
0258   G4bool armFlag = false;
0259   radioactiveDecay->SetARM(armFlag);  // Atomic Rearangement
0260
0261   // atomic de-excitation module
0262   if (armFlag) {
0263     G4EmParameters::Instance()->SetAuger(true);
0264     G4EmParameters::Instance()->SetDeexcitationIgnoreCut(true);
0265   }
0266
0267   G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
0268   ph->RegisterProcess(radioactiveDecay, G4GenericIon::GenericIon());
0269
0270   // mandatory for G4NuclideTable
0271   //
0272   const G4double meanLife = 1 * picosecond, halfLife = meanLife * std::log(2);
0273   G4NuclideTable::GetInstance()->SetThresholdOfHalfLife(halfLife);
0274 }
0275
0276 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0277
0278 void PhysicsList::AddStepMax()
0279 {
0280   // Step limitation seen as a process
0281   StepMax* stepMaxProcess = new StepMax();
0282
0283   auto particleIterator = GetParticleIterator();
0284   particleIterator->reset();
0285   while ((*particleIterator)()) {
0286     G4ParticleDefinition* particle = particleIterator->value();
0287     G4ProcessManager* pmanager = particle->GetProcessManager();
0288
0289     if (stepMaxProcess->IsApplicable(*particle) && !particle->IsShortLived())
0290       pmanager->AddDiscreteProcess(stepMaxProcess);
0291   }
0292 }
0293
0294 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0295
0296 void PhysicsList::GetRange(G4double val)
0297 {
0298   G4LogicalVolume* lBox = fDet->GetWorld()->GetLogicalVolume();
0299   const G4MaterialCutsCouple* couple = lBox->GetMaterialCutsCouple();
0300   const G4Material* currMat = lBox->GetMaterial();
0301
0302   G4ParticleDefinition* part;
0303   G4double cut;
0304   part = G4Electron::Electron();
0305   cut = G4LossTableManager::Instance()->GetRange(part, val, couple);
0306   G4cout << "material : " << currMat->GetName() << G4endl;
0307   G4cout << "particle : " << part->GetParticleName() << G4endl;
0308   G4cout << "energy   : " << G4BestUnit(val, "Energy") << G4endl;
0309   G4cout << "range    : " << G4BestUnit(cut, "Length") << G4endl;
0310 }
0311
0312 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0313
0314 G4VProcess* PhysicsList::GetProcess(const G4String& processName) const
0315 {
0316   G4ParticleDefinition* particle = G4GenericIon::GenericIon();
0317   G4ProcessVector* procList = particle->GetProcessManager()->GetProcessList();
0318   G4int nbProc = particle->GetProcessManager()->GetProcessListLength();
0319   for (G4int k = 0; k < nbProc; k++) {
0320     G4VProcess* process = (*procList)[k];
0321     if (process->GetProcessName() == processName) return process;
0322   }
0323   return nullptr;
0324 }
0325
0326 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......