Geant4 Cross Reference |
1 // 2 // ******************************************************************** 3 // * License and Disclaimer * 4 // * * 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. * 10 // * * 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitation of liability. * 17 // * * 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************************************** 25 // 26 /// \file electromagnetic/TestEm18/src/PhysListEmPenelope.cc 27 /// \brief Implementation of the PhysListEmPenelope class 28 // 29 // 30 // 31 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 32 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 33 34 #include "PhysListEmPenelope.hh" 35 36 #include "G4BuilderType.hh" 37 #include "G4ParticleDefinition.hh" 38 #include "G4PhysicsListHelper.hh" 39 #include "G4ProcessManager.hh" 40 41 // gamma 42 43 #include "G4ComptonScattering.hh" 44 #include "G4GammaConversion.hh" 45 #include "G4PenelopeComptonModel.hh" 46 #include "G4PenelopeGammaConversionModel.hh" 47 #include "G4PenelopePhotoElectricModel.hh" 48 #include "G4PenelopeRayleighModel.hh" 49 #include "G4PhotoElectricEffect.hh" 50 #include "G4RayleighScattering.hh" 51 52 // e- 53 54 #include "G4PenelopeBremsstrahlungModel.hh" 55 #include "G4PenelopeIonisationModel.hh" 56 #include "G4UniversalFluctuation.hh" 57 #include "G4eBremsstrahlung.hh" 58 #include "G4eIonisation.hh" 59 60 // e+ 61 62 #include "G4PenelopeAnnihilationModel.hh" 63 #include "G4eplusAnnihilation.hh" 64 65 // mu 66 67 #include "G4MuBremsstrahlung.hh" 68 #include "G4MuIonisation.hh" 69 #include "G4MuPairProduction.hh" 70 71 // hadrons, ions 72 73 #include "G4hIonisation.hh" 74 #include "G4ionIonisation.hh" 75 76 // deexcitation 77 78 #include "G4LossTableManager.hh" 79 #include "G4SystemOfUnits.hh" 80 #include "G4UAtomicDeexcitation.hh" 81 82 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 83 84 PhysListEmPenelope::PhysListEmPenelope(const G4String& name) : G4VPhysicsConstructor(name) 85 { 86 G4EmParameters* param = G4EmParameters::Instance(); 87 param->SetDefaults(); 88 param->SetMinEnergy(10 * eV); 89 param->SetMaxEnergy(10 * TeV); 90 param->SetNumberOfBinsPerDecade(10); 91 param->SetBuildCSDARange(true); 92 param->SetMaxEnergyForCSDARange(10 * TeV); 93 SetPhysicsType(bElectromagnetic); 94 95 param->SetVerbose(0); 96 param->Dump(); 97 } 98 99 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 100 101 void PhysListEmPenelope::ConstructProcess() 102 { 103 G4PhysicsListHelper* list = G4PhysicsListHelper::GetPhysicsListHelper(); 104 105 // Add standard EM Processes 106 107 auto particleIterator = GetParticleIterator(); 108 particleIterator->reset(); 109 while ((*particleIterator)()) { 110 G4ParticleDefinition* particle = particleIterator->value(); 111 G4String particleName = particle->GetParticleName(); 112 113 // Applicability range for Penelope models 114 // for higher energies, the Standard models are used 115 G4double highEnergyLimit = 1 * GeV; 116 117 if (particleName == "gamma") { 118 // gamma 119 120 G4PhotoElectricEffect* phot = new G4PhotoElectricEffect(); 121 G4PenelopePhotoElectricModel* photModel = new G4PenelopePhotoElectricModel(); 122 photModel->SetHighEnergyLimit(highEnergyLimit); 123 phot->AddEmModel(0, photModel); 124 list->RegisterProcess(phot, particle); 125 126 G4ComptonScattering* compt = new G4ComptonScattering(); 127 G4PenelopeComptonModel* comptModel = new G4PenelopeComptonModel(); 128 comptModel->SetHighEnergyLimit(highEnergyLimit); 129 compt->AddEmModel(0, comptModel); 130 list->RegisterProcess(compt, particle); 131 132 G4GammaConversion* conv = new G4GammaConversion(); 133 G4PenelopeGammaConversionModel* convModel = new G4PenelopeGammaConversionModel(); 134 convModel->SetHighEnergyLimit(highEnergyLimit); 135 conv->AddEmModel(0, convModel); 136 list->RegisterProcess(conv, particle); 137 138 G4RayleighScattering* rayl = new G4RayleighScattering(); 139 G4PenelopeRayleighModel* raylModel = new G4PenelopeRayleighModel(); 140 raylModel->SetHighEnergyLimit(highEnergyLimit); 141 rayl->AddEmModel(0, raylModel); 142 list->RegisterProcess(rayl, particle); 143 } 144 else if (particleName == "e-") { 145 // electron 146 147 G4eIonisation* eIoni = new G4eIonisation(); 148 G4PenelopeIonisationModel* eIoniModel = new G4PenelopeIonisationModel(); 149 eIoniModel->SetHighEnergyLimit(highEnergyLimit); 150 eIoni->AddEmModel(0, eIoniModel, new G4UniversalFluctuation()); 151 list->RegisterProcess(eIoni, particle); 152 153 G4eBremsstrahlung* eBrem = new G4eBremsstrahlung(); 154 G4PenelopeBremsstrahlungModel* eBremModel = new G4PenelopeBremsstrahlungModel(); 155 eBremModel->SetHighEnergyLimit(highEnergyLimit); 156 eBrem->AddEmModel(0, eBremModel); 157 list->RegisterProcess(eBrem, particle); 158 } 159 else if (particleName == "e+") { 160 // positron 161 G4eIonisation* eIoni = new G4eIonisation(); 162 G4PenelopeIonisationModel* eIoniModel = new G4PenelopeIonisationModel(); 163 eIoniModel->SetHighEnergyLimit(highEnergyLimit); 164 eIoni->AddEmModel(0, eIoniModel, new G4UniversalFluctuation()); 165 list->RegisterProcess(eIoni, particle); 166 167 G4eBremsstrahlung* eBrem = new G4eBremsstrahlung(); 168 G4PenelopeBremsstrahlungModel* eBremModel = new G4PenelopeBremsstrahlungModel(); 169 eBremModel->SetHighEnergyLimit(highEnergyLimit); 170 eBrem->AddEmModel(0, eBremModel); 171 list->RegisterProcess(eBrem, particle); 172 173 G4eplusAnnihilation* eAnni = new G4eplusAnnihilation(); 174 G4PenelopeAnnihilationModel* eAnniModel = new G4PenelopeAnnihilationModel(); 175 eAnniModel->SetHighEnergyLimit(highEnergyLimit); 176 eAnni->AddEmModel(0, eAnniModel); 177 list->RegisterProcess(eAnni, particle); 178 } 179 else if (particleName == "mu+" || particleName == "mu-") { 180 // muon 181 list->RegisterProcess(new G4MuIonisation, particle); 182 list->RegisterProcess(new G4MuBremsstrahlung, particle); 183 list->RegisterProcess(new G4MuPairProduction, particle); 184 } 185 else if (particleName == "alpha" || particleName == "GenericIon") { 186 list->RegisterProcess(new G4ionIonisation, particle); 187 } 188 else if ((!particle->IsShortLived()) && (particle->GetPDGCharge() != 0.0) 189 && (particle->GetParticleName() != "chargedgeantino")) 190 { 191 // all others charged particles except geantino 192 list->RegisterProcess(new G4hIonisation, particle); 193 } 194 } 195 // Deexcitation 196 // 197 G4VAtomDeexcitation* deex = new G4UAtomicDeexcitation(); 198 G4LossTableManager::Instance()->SetAtomDeexcitation(deex); 199 } 200 201 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 202