- Timestamp:
- Sep 30, 2010, 2:47:17 PM (14 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/examples/advanced/microbeam/src/MicrobeamPhysicsList.cc
r1230 r1337 24 24 // ******************************************************************** 25 25 // 26 // ------------------------------------------------------------------- 27 // $Id: MicrobeamPhysicsList.cc,v 1.8 2009/04/30 10:23:57 sincerti Exp $ 28 // ------------------------------------------------------------------- 29 30 #include "G4ParticleDefinition.hh" 26 // $Id: MicrobeamPhysicsList.cc,v 1.10 2010/06/10 09:54:05 vnivanch Exp $ 27 // GEANT4 tag $Name: geant4-09-04-beta-01 $ 28 // 29 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 30 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 31 32 #include "MicrobeamPhysicsList.hh" 33 #include "MicrobeamPhysicsListMessenger.hh" 34 35 #include "G4StepLimiter.hh" 36 37 #include "G4EmStandardPhysics.hh" 38 #include "G4EmStandardPhysics_option1.hh" 39 #include "G4EmStandardPhysics_option2.hh" 40 #include "G4EmStandardPhysics_option3.hh" 41 #include "G4EmLivermorePhysics.hh" 42 #include "G4EmPenelopePhysics.hh" 43 44 #include "G4DecayPhysics.hh" 45 46 #include "G4HadronElasticPhysics.hh" 47 #include "G4HadronDElasticPhysics.hh" 48 #include "G4HadronHElasticPhysics.hh" 49 #include "G4HadronQElasticPhysics.hh" 50 #include "G4HadronInelasticQBBC.hh" 51 #include "G4IonBinaryCascadePhysics.hh" 52 53 #include "G4LossTableManager.hh" 54 #include "G4EmConfigurator.hh" 55 #include "G4UnitsTable.hh" 56 31 57 #include "G4ProcessManager.hh" 32 #include "G4ParticleTypes.hh" 33 #include "G4StepLimiter.hh" 34 #include "G4BaryonConstructor.hh" 35 #include "G4IonConstructor.hh" 36 #include "G4MesonConstructor.hh" 37 38 #include "MicrobeamPhysicsList.hh" 39 40 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 41 42 MicrobeamPhysicsList::MicrobeamPhysicsList(): G4VUserPhysicsList() 43 { 58 #include "G4Decay.hh" 59 60 #include "G4IonFluctuations.hh" 61 #include "G4IonParametrisedLossModel.hh" 62 63 #include "G4BraggIonGasModel.hh" 64 #include "G4BetheBlochIonGasModel.hh" 65 66 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 67 68 MicrobeamPhysicsList::MicrobeamPhysicsList() : G4VModularPhysicsList() 69 { 70 G4LossTableManager::Instance(); 44 71 defaultCutValue = 0.01*micrometer; 45 72 cutForGamma = defaultCutValue; 46 73 cutForElectron = defaultCutValue; 47 74 cutForPositron = defaultCutValue; 75 76 helIsRegisted = false; 77 bicIsRegisted = false; 78 biciIsRegisted = false; 79 80 pMessenger = new MicrobeamPhysicsListMessenger(this); 81 48 82 SetVerboseLevel(1); 49 } 50 51 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 83 84 // EM physics 85 emName = G4String("emstandard_opt0"); 86 87 emPhysicsList = new G4EmStandardPhysics(1); 88 89 // Deacy physics and all particles 90 decPhysicsList = new G4DecayPhysics(); 91 } 92 93 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 52 94 53 95 MicrobeamPhysicsList::~MicrobeamPhysicsList() 54 {} 55 56 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 96 { 97 delete pMessenger; 98 delete emPhysicsList; 99 delete decPhysicsList; 100 for(size_t i=0; i<hadronPhys.size(); i++) {delete hadronPhys[i];} 101 } 102 103 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 57 104 58 105 void MicrobeamPhysicsList::ConstructParticle() 59 106 { 60 ConstructBosons(); 61 ConstructLeptons(); 62 ConstructBaryons(); 63 } 64 65 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 66 67 void MicrobeamPhysicsList::ConstructBosons() 68 { 69 // gamma 70 G4Gamma::GammaDefinition(); 71 72 // optical photon 73 G4OpticalPhoton::OpticalPhotonDefinition(); 74 } 75 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 76 77 void MicrobeamPhysicsList::ConstructLeptons() 78 { 79 // leptons 80 G4Electron::ElectronDefinition(); 81 G4Positron::PositronDefinition(); 82 } 83 84 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 85 86 void MicrobeamPhysicsList::ConstructBaryons() 87 { 88 // baryons 89 G4BaryonConstructor bConstructor; 90 bConstructor.ConstructParticle(); 91 92 G4IonConstructor iConstructor; 93 iConstructor.ConstructParticle(); 94 95 G4MesonConstructor mConstructor; 96 mConstructor.ConstructParticle(); 97 } 98 99 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 107 decPhysicsList->ConstructParticle(); 108 } 109 110 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 100 111 101 112 void MicrobeamPhysicsList::ConstructProcess() 102 113 { 114 // transportation 115 // 103 116 AddTransportation(); 104 ConstructEM(); 105 ConstructHad(); 106 ConstructGeneral(); 107 } 108 109 // *** Processes and models 110 111 // gamma 112 113 #include "G4PhotoElectricEffect.hh" 114 #include "G4LivermorePhotoElectricModel.hh" 115 116 #include "G4ComptonScattering.hh" 117 #include "G4LivermoreComptonModel.hh" 118 119 #include "G4GammaConversion.hh" 120 #include "G4LivermoreGammaConversionModel.hh" 121 122 #include "G4RayleighScattering.hh" 123 #include "G4LivermoreRayleighModel.hh" 124 125 // e- 126 127 #include "G4eMultipleScattering.hh" 128 #include "G4UniversalFluctuation.hh" 129 130 #include "G4eIonisation.hh" 131 #include "G4LivermoreIonisationModel.hh" 132 133 #include "G4eBremsstrahlung.hh" 134 #include "G4LivermoreBremsstrahlungModel.hh" 135 136 // e+ 137 138 #include "G4eplusAnnihilation.hh" 139 140 // mu 141 142 #include "G4MuIonisation.hh" 143 #include "G4MuBremsstrahlung.hh" 144 #include "G4MuPairProduction.hh" 145 146 // hadrons 147 148 #include "G4hMultipleScattering.hh" 149 #include "G4MscStepLimitType.hh" 150 151 #include "G4hBremsstrahlung.hh" 152 #include "G4hPairProduction.hh" 153 154 #include "G4hIonisation.hh" 155 #include "G4ionIonisation.hh" 156 #include "G4IonParametrisedLossModel.hh" 157 158 // 159 160 #include "G4LossTableManager.hh" 161 #include "G4EmProcessOptions.hh" 162 163 164 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 165 166 void MicrobeamPhysicsList::ConstructEM() 167 { 117 118 // electromagnetic physics list 119 // 120 emPhysicsList->ConstructProcess(); 121 122 // decay physics list 123 // 124 decPhysicsList->ConstructProcess(); 125 126 // hadronic physics lists 127 for(size_t i=0; i<hadronPhys.size(); i++) { 128 hadronPhys[i]->ConstructProcess(); 129 } 130 131 // step limitation (as a full process) 132 // 133 AddStepMax(); 134 135 } 136 137 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 138 139 void MicrobeamPhysicsList::AddPhysicsList(const G4String& name) 140 { 141 if (verboseLevel>1) { 142 G4cout << "PhysicsList::AddPhysicsList: <" << name << ">" << G4endl; 143 } 144 145 if (name == emName) return; 146 147 if (name == "emstandard_opt0") { 148 149 emName = name; 150 delete emPhysicsList; 151 emPhysicsList = new G4EmStandardPhysics(1); 152 153 } else if (name == "emstandard_opt1") { 154 155 emName = name; 156 delete emPhysicsList; 157 emPhysicsList = new G4EmStandardPhysics_option1(); 158 159 } else if (name == "emstandard_opt2") { 160 161 emName = name; 162 delete emPhysicsList; 163 emPhysicsList = new G4EmStandardPhysics_option2(); 164 165 } else if (name == "emstandard_opt3") { 166 167 emName = name; 168 delete emPhysicsList; 169 emPhysicsList = new G4EmStandardPhysics_option3(); 170 171 } else if (name == "emlivermore") { 172 emName = name; 173 delete emPhysicsList; 174 emPhysicsList = new G4EmLivermorePhysics(); 175 176 } else if (name == "empenelope") { 177 emName = name; 178 delete emPhysicsList; 179 emPhysicsList = new G4EmPenelopePhysics(); 180 181 } else if (name == "ionGasModels") { 182 183 AddPhysicsList("emstandard_opt3"); 184 emName = name; 185 AddIonGasModels(); 186 187 } else if (name == "elastic" && !helIsRegisted) { 188 hadronPhys.push_back( new G4HadronElasticPhysics()); 189 helIsRegisted = true; 190 191 } else if (name == "DElastic" && !helIsRegisted) { 192 hadronPhys.push_back( new G4HadronDElasticPhysics()); 193 helIsRegisted = true; 194 195 } else if (name == "HElastic" && !helIsRegisted) { 196 hadronPhys.push_back( new G4HadronHElasticPhysics()); 197 helIsRegisted = true; 198 199 } else if (name == "QElastic" && !helIsRegisted) { 200 hadronPhys.push_back( new G4HadronQElasticPhysics()); 201 helIsRegisted = true; 202 203 } else if (name == "binary" && !bicIsRegisted) { 204 hadronPhys.push_back(new G4HadronInelasticQBBC()); 205 bicIsRegisted = true; 206 207 } else if (name == "binary_ion" && !biciIsRegisted) { 208 hadronPhys.push_back(new G4IonBinaryCascadePhysics()); 209 biciIsRegisted = true; 210 211 } else { 212 213 G4cout << "PhysicsList::AddPhysicsList: <" << name << ">" 214 << " is not defined" 215 << G4endl; 216 } 217 } 218 219 220 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 221 222 void MicrobeamPhysicsList::AddStepMax() 223 { 224 // Step limitation seen as a process 225 stepMaxProcess = new G4StepLimiter(); 226 168 227 theParticleIterator->reset(); 169 170 while( (*theParticleIterator)() ){ 171 228 while ((*theParticleIterator)()){ 172 229 G4ParticleDefinition* particle = theParticleIterator->value(); 173 230 G4ProcessManager* pmanager = particle->GetProcessManager(); 174 G4String particleName = particle->GetParticleName(); 175 176 if (particleName == "gamma") { 177 178 179 G4PhotoElectricEffect* thePhotoElectricEffect = new G4PhotoElectricEffect(); 180 G4LivermorePhotoElectricModel* theLivermorePhotoElectricModel = new G4LivermorePhotoElectricModel(); 181 thePhotoElectricEffect->AddEmModel(0, theLivermorePhotoElectricModel); 182 pmanager->AddDiscreteProcess(thePhotoElectricEffect); 183 184 G4ComptonScattering* theComptonScattering = new G4ComptonScattering(); 185 G4LivermoreComptonModel* theLivermoreComptonModel = new G4LivermoreComptonModel(); 186 theComptonScattering->AddEmModel(0, theLivermoreComptonModel); 187 pmanager->AddDiscreteProcess(theComptonScattering); 188 189 G4GammaConversion* theGammaConversion = new G4GammaConversion(); 190 G4LivermoreGammaConversionModel* theLivermoreGammaConversionModel = new G4LivermoreGammaConversionModel(); 191 theGammaConversion->AddEmModel(0, theLivermoreGammaConversionModel); 192 pmanager->AddDiscreteProcess(theGammaConversion); 193 194 G4RayleighScattering* theRayleigh = new G4RayleighScattering(); 195 G4LivermoreRayleighModel* theRayleighModel = new G4LivermoreRayleighModel(); 196 theRayleigh->AddEmModel(0, theRayleighModel); 197 pmanager->AddDiscreteProcess(theRayleigh); 198 199 pmanager->AddProcess(new G4StepLimiter(), -1, -1, 5); 200 201 } else if (particleName == "e-") { 202 203 G4eMultipleScattering* msc = new G4eMultipleScattering(); 204 msc->SetStepLimitType(fUseDistanceToBoundary); 205 pmanager->AddProcess(msc, -1, 1, 1); 206 207 // Ionisation 208 G4eIonisation* eIoni = new G4eIonisation(); 209 eIoni->AddEmModel(0, new G4LivermoreIonisationModel(), new G4UniversalFluctuation() ); 210 eIoni->SetStepFunction(0.2, 100*um); // 211 pmanager->AddProcess(eIoni, -1, 2, 2); 212 213 // Bremsstrahlung 214 G4eBremsstrahlung* eBrem = new G4eBremsstrahlung(); 215 eBrem->AddEmModel(0, new G4LivermoreBremsstrahlungModel()); 216 pmanager->AddProcess(eBrem, -1,-3, 3); 217 218 pmanager->AddProcess(new G4StepLimiter(), -1, -1, 4); 219 220 } else if (particleName == "e+") { 221 222 // Identical to G4EmStandardPhysics_option3 223 224 G4eMultipleScattering* msc = new G4eMultipleScattering(); 225 msc->SetStepLimitType(fUseDistanceToBoundary); 226 pmanager->AddProcess(msc, -1, 1, 1); 227 228 G4eIonisation* eIoni = new G4eIonisation(); 229 eIoni->SetStepFunction(0.2, 100*um); 230 pmanager->AddProcess(eIoni, -1, 2, 2); 231 232 pmanager->AddProcess(new G4eBremsstrahlung, -1,-3, 3); 233 234 pmanager->AddProcess(new G4eplusAnnihilation,0,-1, 4); 235 236 pmanager->AddProcess(new G4StepLimiter(), -1, -1, 5); 237 238 } else if (particleName == "GenericIon") { 239 240 pmanager->AddProcess(new G4hMultipleScattering, -1, 1, 1); 241 242 G4ionIonisation* ionIoni = new G4ionIonisation(); 243 ionIoni->SetEmModel(new G4IonParametrisedLossModel()); 244 ionIoni->SetStepFunction(0.1, 20*um); 245 pmanager->AddProcess(ionIoni, -1, 2, 2); 246 247 pmanager->AddProcess(new G4StepLimiter(), -1, -1, 3); 248 249 } else if (particleName == "alpha" || 250 particleName == "He3" ) { 251 252 // Identical to G4EmStandardPhysics_option3 253 254 pmanager->AddProcess(new G4hMultipleScattering, -1, 1, 1); 255 256 G4ionIonisation* ionIoni = new G4ionIonisation(); 257 ionIoni->SetStepFunction(0.1, 20*um); 258 pmanager->AddProcess(ionIoni, -1, 2, 2); 259 260 pmanager->AddProcess(new G4StepLimiter(), -1, -1, 3); 261 } 262 263 // 264 265 } 266 } 267 268 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 269 270 #include "G4HadronElasticProcess.hh" 271 #include "G4LElastic.hh" 272 273 #include "G4AlphaInelasticProcess.hh" 274 #include "G4BinaryLightIonReaction.hh" 275 #include "G4TripathiCrossSection.hh" 276 #include "G4IonsShenCrossSection.hh" 277 #include "G4LEAlphaInelastic.hh" 278 279 void MicrobeamPhysicsList::ConstructHad() 280 { 281 282 G4HadronElasticProcess * theElasticProcess = new G4HadronElasticProcess; 283 theElasticProcess->RegisterMe( new G4LElastic() ); 284 285 theParticleIterator->reset(); 286 while( (*theParticleIterator)() ) 287 { 288 G4ParticleDefinition* particle = theParticleIterator->value(); 289 G4ProcessManager* pManager = particle->GetProcessManager(); 290 291 if (particle->GetParticleName() == "alpha") 292 { 293 294 // INELASTIC SCATTERING 295 // Binary Cascade 296 G4BinaryLightIonReaction* theBC = new G4BinaryLightIonReaction(); 297 theBC -> SetMinEnergy(80.*MeV); 298 theBC -> SetMaxEnergy(40.*GeV); 299 300 // TRIPATHI CROSS SECTION 301 // Implementation of formulas in analogy to NASA technical paper 3621 by 302 // Tripathi, et al. Cross-sections for ion ion scattering 303 G4TripathiCrossSection* TripathiCrossSection = new G4TripathiCrossSection; 304 305 // IONS SHEN CROSS SECTION 306 // Implementation of formulas 307 // Shen et al. Nuc. Phys. A 491 130 (1989) 308 // Total Reaction Cross Section for Heavy-Ion Collisions 309 G4IonsShenCrossSection* aShen = new G4IonsShenCrossSection; 310 311 // Final state production model for Alpha inelastic scattering below 20 GeV 312 G4LEAlphaInelastic* theAIModel = new G4LEAlphaInelastic; 313 theAIModel -> SetMaxEnergy(100.*MeV); 314 315 G4AlphaInelasticProcess * theIPalpha = new G4AlphaInelasticProcess; 316 theIPalpha->AddDataSet(TripathiCrossSection); 317 theIPalpha->AddDataSet(aShen); 318 319 // Register the Alpha Inelastic and Binary Cascade Model 320 theIPalpha->RegisterMe(theAIModel); 321 theIPalpha->RegisterMe(theBC); 322 323 // Activate the alpha inelastic scattering using the alpha inelastic and binary cascade model 324 pManager -> AddDiscreteProcess(theIPalpha); 325 326 // Activate the Hadron Elastic Process 327 pManager -> AddDiscreteProcess(theElasticProcess); 328 329 } 330 } 331 332 } 333 334 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 335 336 void MicrobeamPhysicsList::ConstructGeneral() 337 { } 338 339 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 231 232 if (stepMaxProcess->IsApplicable(*particle) && pmanager) 233 { 234 pmanager ->AddDiscreteProcess(stepMaxProcess); 235 } 236 } 237 } 238 239 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 340 240 341 241 void MicrobeamPhysicsList::SetCuts() 342 242 { 243 343 244 if (verboseLevel >0){ 344 G4cout << " MicrobeamPhysicsList::SetCuts:";245 G4cout << "PhysicsList::SetCuts:"; 345 246 G4cout << "CutLength : " << G4BestUnit(defaultCutValue,"Length") << G4endl; 346 } 347 247 } 248 348 249 // set cut values for gamma at first and for e- second and next for e+, 349 // because some processes for e+/e- need cut values for gamma 250 // because some processes for e+/e- need cut values for gamma 350 251 SetCutValue(cutForGamma, "gamma"); 351 252 SetCutValue(cutForElectron, "e-"); 352 253 SetCutValue(cutForPositron, "e+"); 353 254 354 255 if (verboseLevel>0) DumpCutValuesTable(); 355 356 } 357 256 } 257 258 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 259 260 void MicrobeamPhysicsList::SetCutForGamma(G4double cut) 261 { 262 cutForGamma = cut; 263 SetParticleCuts(cutForGamma, G4Gamma::Gamma()); 264 } 265 266 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 267 268 void MicrobeamPhysicsList::SetCutForElectron(G4double cut) 269 { 270 cutForElectron = cut; 271 SetParticleCuts(cutForElectron, G4Electron::Electron()); 272 } 273 274 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 275 276 void MicrobeamPhysicsList::SetCutForPositron(G4double cut) 277 { 278 cutForPositron = cut; 279 SetParticleCuts(cutForPositron, G4Positron::Positron()); 280 } 281 282 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 283 284 void MicrobeamPhysicsList::AddIonGasModels() 285 { 286 G4EmConfigurator* em_config = G4LossTableManager::Instance()->EmConfigurator(); 287 theParticleIterator->reset(); 288 while ((*theParticleIterator)()) 289 { 290 G4ParticleDefinition* particle = theParticleIterator->value(); 291 G4String partname = particle->GetParticleName(); 292 if(partname == "alpha" || partname == "He3" || partname == "GenericIon") { 293 G4BraggIonGasModel* mod1 = new G4BraggIonGasModel(); 294 G4BetheBlochIonGasModel* mod2 = new G4BetheBlochIonGasModel(); 295 G4double eth = 2.*MeV*particle->GetPDGMass()/proton_mass_c2; 296 em_config->SetExtraEmModel(partname,"ionIoni",mod1,"",0.0,eth, 297 new G4IonFluctuations()); 298 em_config->SetExtraEmModel(partname,"ionIoni",mod2,"",eth,100*TeV, 299 new G4UniversalFluctuation()); 300 301 } 302 } 303 } 304 305 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 306
Note: See TracChangeset
for help on using the changeset viewer.