Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer 3 // * License and Disclaimer * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. 9 // * include a list of copyright holders. * 10 // * 10 // * * 11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. * 17 // * 17 // * * 18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************* 24 // ******************************************************************** 25 // 25 // 26 // Hadrontherapy advanced example for Geant4 << 26 // $Id: HadrontherapyPhysicsList.cc,v 1.0 27 // See more at: https://twiki.cern.ch/twiki/bi << 27 // ---------------------------------------------------------------------------- >> 28 // GEANT 4 - Hadrontherapy example >> 29 // ---------------------------------------------------------------------------- >> 30 // Code developed by: 28 // 31 // 29 // Using the builder concepts of Geant4 we ass << 32 // G.A.P. Cirrone(a)*, F. Di Rosa(a), S. Guatelli(b), G. Russo(a) 30 // Physics Lists that are particuilarly suited << 33 // 31 // << 34 // (a) Laboratori Nazionali del Sud 32 // 'HADRONTHERAPY_1' is more suited for proton << 35 // of the National Institute for Nuclear Physics, Catania, Italy 33 // 'HADRONTHERAPY_2' is suggested for better p << 36 // (b) National Institute for Nuclear Physics Section of Genova, genova, Italy 34 // << 37 // 35 // The Reference physics lists (already presen << 38 // * cirrone@lns.infn.it 36 // be used as well. In this case the more suit << 39 // ---------------------------------------------------------------------------- 37 // "QBBC", "QGSP_BIC", "Shielding", "QGSP_BERT << 38 // "QGSP_BIC_AllHP" and "QGSP_BIC_HP" << 39 // << 40 // NOTE: to activate the "_HP" physics you hav << 41 // variable pointing to the external dataset n << 42 // << 43 // All the lists can be activated inside any << 44 // /Physics/addPhysics << 45 // << 46 // Examples of usage are: << 47 // /Physics/addPhysics HADRONTHERAPY_1 or /Phy << 48 40 49 #include "G4SystemOfUnits.hh" << 41 #include "globals.hh" 50 #include "G4RunManager.hh" << 42 #include "G4ProcessManager.hh" 51 #include "G4Region.hh" 43 #include "G4Region.hh" 52 #include "G4RegionStore.hh" 44 #include "G4RegionStore.hh" 53 #include "HadrontherapyPhysicsList.hh" 45 #include "HadrontherapyPhysicsList.hh" 54 #include "HadrontherapyPhysicsListMessenger.hh 46 #include "HadrontherapyPhysicsListMessenger.hh" 55 #include "HadrontherapyStepMax.hh" << 47 #include "HadrontherapyParticles.hh" 56 #include "G4PhysListFactory.hh" << 48 #include "HadrontherapyPhotonStandard.hh" 57 #include "G4VPhysicsConstructor.hh" << 49 #include "HadrontherapyPhotonEPDL.hh" 58 #include "G4HadronPhysicsQGSP_BIC_HP.hh" << 50 #include "HadrontherapyPhotonPenelope.hh" 59 #include "G4HadronPhysicsQGSP_BIC.hh" << 51 #include "HadrontherapyElectronStandard.hh" 60 #include "G4EmStandardPhysics_option4.hh" << 52 #include "HadrontherapyElectronEEDL.hh" 61 #include "G4EmStandardPhysics.hh" << 53 #include "HadrontherapyElectronPenelope.hh" 62 #include "G4EmExtraPhysics.hh" << 54 #include "HadrontherapyPositronStandard.hh" 63 #include "G4StoppingPhysics.hh" << 55 #include "HadrontherapyPositronPenelope.hh" 64 #include "G4DecayPhysics.hh" << 56 #include "HadrontherapyIonLowE.hh" 65 #include "G4HadronElasticPhysics.hh" << 57 #include "HadrontherapyIonLowEZiegler1977.hh" 66 #include "G4HadronElasticPhysicsHP.hh" << 58 #include "HadrontherapyIonLowEZiegler1985.hh" 67 #include "G4RadioactiveDecayPhysics.hh" << 59 #include "HadrontherapyIonLowEZiegler2000.hh" 68 #include "G4IonBinaryCascadePhysics.hh" << 60 #include "HadrontherapyIonStandard.hh" 69 #include "G4DecayPhysics.hh" << 61 #include "HadrontherapyProtonPrecompound.hh" 70 #include "G4NeutronTrackingCut.hh" << 62 #include "HadrontherapyProtonPrecompoundFermi.hh" 71 #include "G4LossTableManager.hh" << 63 #include "HadrontherapyProtonPrecompoundGEM.hh" 72 #include "G4UnitsTable.hh" << 64 #include "HadrontherapyProtonPrecompoundGEMFermi.hh" 73 #include "G4ProcessManager.hh" << 65 #include "HadrontherapyProtonBertini.hh" 74 #include "G4IonFluctuations.hh" << 66 #include "HadrontherapyProtonBinary.hh" 75 #include "G4IonParametrisedLossModel.hh" << 67 #include "HadrontherapyMuonStandard.hh" 76 #include "G4EmParameters.hh" << 68 #include "HadrontherapyDecay.hh" 77 #include "G4ParallelWorldPhysics.hh" << 69 #include "HadrontherapyParticles.hh" 78 #include "G4EmLivermorePhysics.hh" << 70 #include "G4ParticleDefinition.hh" 79 #include "G4AutoDelete.hh" << 71 #include "G4ParticleTypes.hh" 80 #include "G4HadronPhysicsQGSP_BIC_AllHP.hh" << 72 #include "G4ParticleTable.hh" 81 #include "QGSP_BIC_HP.hh" << 73 HadrontherapyPhysicsList::HadrontherapyPhysicsList(): G4VModularPhysicsList(), 82 #include "QGSP_BIC.hh" << 74 electronIsRegistered(false), 83 #include "G4HadronPhysicsQGSP_BERT.hh" << 75 positronIsRegistered(false), 84 #include "G4HadronPhysicsQGSP_BERT_HP.hh" << 76 photonIsRegistered(false), 85 #include "G4ParallelWorldPhysics.hh" << 77 ionIsRegistered(false), 86 // Physics List << 78 protonHadronicIsRegistered(false), 87 #include "QBBC.hh" << 79 muonIsRegistered(false), 88 #include "QGSP_BIC.hh" << 80 decayIsRegistered(false) 89 #include "Shielding.hh" << 81 { 90 #include "QGSP_BERT.hh" << 82 // 91 #include "QGSP_BIC_AllHP.hh" << 83 // The threshold of production of secondaries is fixed to 10. mm 92 #include "QGSP_BIC_HP.hh" << 84 // for all the particles, in all the experimental set-up >> 85 // The phantom is defined as a Geant4 Region. Here the cut is fixed to 0.001 * mm. >> 86 // >> 87 >> 88 defaultCutValue = 10. * mm; 93 89 >> 90 // Messenger: it is possible to activate interactively physics processes and models >> 91 messenger = new HadrontherapyPhysicsListMessenger(this); 94 92 >> 93 SetVerboseLevel(1); 95 94 96 ////////////////////////////////////////////// << 95 // Register all the particles involved in the experimental set-up 97 HadrontherapyPhysicsList::HadrontherapyPhysics << 96 RegisterPhysics( new HadrontherapyParticles("particles") ); 98 { << 99 G4LossTableManager::Instance(); << 100 defaultCutValue = 1.*mm; << 101 cutForGamma = defaultCutValue; << 102 cutForElectron = defaultCutValue; << 103 cutForPositron = defaultCutValue; << 104 << 105 pMessenger = new HadrontherapyPhysicsListM << 106 SetVerboseLevel(1); << 107 decay_List = new G4DecayPhysics(); << 108 // Elecromagnetic physics << 109 // << 110 emPhysicsList = new G4EmStandardPhysics_op << 111 } 97 } 112 98 113 ////////////////////////////////////////////// << 114 HadrontherapyPhysicsList::~HadrontherapyPhysic 99 HadrontherapyPhysicsList::~HadrontherapyPhysicsList() >> 100 { >> 101 delete messenger; >> 102 } >> 103 >> 104 void HadrontherapyPhysicsList::AddPhysicsList(const G4String& name) 115 { 105 { 116 delete pMessenger; << 106 G4cout << "Adding PhysicsList component " << name << G4endl; 117 delete emPhysicsList; << 107 118 delete decay_List; << 108 //*************************// 119 //delete radioactiveDecay_List; << 109 // Electromagnetic physics // 120 hadronPhys.clear(); << 110 //*************************// 121 for(size_t i=0; i<hadronPhys.size(); i++) << 111 // >> 112 // The user can choose three alternative approaches: >> 113 // Standard, Low Energy based on the Livermore libraries and Low Energy Penelope >> 114 // >> 115 >> 116 // ******** PHOTONS ********// >> 117 >> 118 // Register standard processes for photons >> 119 if (name == "photon-standard") 122 { 120 { 123 delete hadronPhys[i]; << 121 if (photonIsRegistered) >> 122 { >> 123 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 124 << " cannot be registered ---- photon List already existing" >> 125 << G4endl; >> 126 } >> 127 else >> 128 { >> 129 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 130 << " is registered" << G4endl; >> 131 RegisterPhysics( new HadrontherapyPhotonStandard(name) ); >> 132 photonIsRegistered = true; >> 133 } 124 } 134 } 125 } << 126 135 127 ////////////////////////////////////////////// << 136 // Register LowE-EPDL processes for photons 128 void HadrontherapyPhysicsList::ConstructPartic << 137 if (name == "photon-epdl") 129 { << 138 { 130 decay_List -> ConstructParticle(); << 139 if (photonIsRegistered) 131 << 140 { 132 } << 141 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 142 << " cannot be registered ---- photon List already existing" >> 143 << G4endl; >> 144 } >> 145 else >> 146 { >> 147 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 148 << " is registered" << G4endl; >> 149 RegisterPhysics( new HadrontherapyPhotonEPDL(name) ); >> 150 photonIsRegistered = true; >> 151 } >> 152 } 133 153 134 ////////////////////////////////////////////// << 154 // Register processes a' la Penelope for photons 135 void HadrontherapyPhysicsList::ConstructProces << 155 if (name == "photon-penelope") 136 { << 137 // Transportation << 138 // << 139 AddTransportation(); << 140 << 141 decay_List -> ConstructProcess(); << 142 emPhysicsList -> ConstructProcess(); << 143 << 144 << 145 //em_config.AddModels(); << 146 << 147 // Hadronic physics << 148 // << 149 for(size_t i=0; i < hadronPhys.size(); i++ << 150 { 156 { 151 hadronPhys[i] -> ConstructProcess(); << 157 if (photonIsRegistered) >> 158 { >> 159 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 160 << " cannot be registered ---- photon List already existing" >> 161 << G4endl; >> 162 } >> 163 else >> 164 { >> 165 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 166 << " is registered" << G4endl; >> 167 RegisterPhysics( new HadrontherapyPhotonPenelope(name) ); >> 168 photonIsRegistered = true; >> 169 } 152 } 170 } 153 << 154 // step limitation (as a full process) << 155 // << 156 AddStepMax(); << 157 << 158 //Parallel world sensitivity << 159 // << 160 G4ParallelWorldPhysics* pWorld = new G4Par << 161 pWorld->ConstructProcess(); << 162 << 163 return; << 164 } << 165 171 166 ////////////////////////////////////////////// << 172 // ******** Electrons ********// 167 void HadrontherapyPhysicsList::AddPhysicsList( << 173 168 { << 174 // Register standard processes for electrons 169 if (verboseLevel>1) { << 175 if (name == "electron-standard") 170 G4cout << "PhysicsList::AddPhysicsList << 176 { >> 177 if (electronIsRegistered) >> 178 { >> 179 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 180 << " cannot be registered ---- electron List already existing" >> 181 << G4endl; >> 182 } >> 183 else >> 184 { >> 185 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 186 << " is registered" << G4endl; >> 187 RegisterPhysics( new HadrontherapyElectronStandard(name) ); >> 188 electronIsRegistered = true; >> 189 } 171 } 190 } 172 if (name == emName) return; << 191 173 << 192 // Register LowE-EEDL processes for electrons 174 /////////////////////////////////// << 193 if (name == "electron-eedl") 175 // ELECTROMAGNETIC MODELS << 194 { 176 /////////////////////////////////// << 195 if (electronIsRegistered) 177 if (name == "standard_opt4") { << 196 { 178 emName = name; << 197 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name 179 delete emPhysicsList; << 198 << " cannot be registered ---- electron List already existing" 180 hadronPhys.clear(); << 199 << G4endl; 181 emPhysicsList = new G4EmStandardPhysic << 200 } 182 G4RunManager::GetRunManager() -> Physi << 201 else 183 G4cout << "THE FOLLOWING ELECTROMAGNET << 202 { 184 << 203 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name 185 ////////////////////////////////////// << 204 << " is registered" << G4endl; 186 // ELECTROMAGNETIC + HADRONIC MODELS << 205 RegisterPhysics( new HadrontherapyElectronEEDL(name) ); 187 ////////////////////////////////////// << 206 electronIsRegistered = true; 188 << 207 } 189 } else if (name == "HADRONTHERAPY_1") { << 208 } 190 << 209 191 AddPhysicsList("standard_opt4"); << 210 // Register processes a' la Penelope for electrons 192 hadronPhys.push_back( new G4DecayPhysi << 211 if (name == "electron-penelope") 193 hadronPhys.push_back( new G4Radioactiv << 212 { 194 hadronPhys.push_back( new G4IonBinaryC << 213 if (electronIsRegistered) 195 hadronPhys.push_back( new G4EmExtraPhy << 214 { 196 hadronPhys.push_back( new G4HadronElas << 215 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name 197 hadronPhys.push_back( new G4StoppingPh << 216 << " cannot be registered ---- electron List already existing" 198 hadronPhys.push_back( new G4HadronPhys << 217 << G4endl; 199 hadronPhys.push_back( new G4NeutronTra << 218 } 200 << 219 else 201 G4cout << "HADRONTHERAPY_1 PHYSICS LIS << 220 { >> 221 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 222 << " is registered" << G4endl; >> 223 RegisterPhysics( new HadrontherapyElectronPenelope(name) ); >> 224 electronIsRegistered = true; >> 225 } 202 } 226 } 203 << 227 204 else if (name == "HADRONTHERAPY_2") { << 228 // ******** Positrons ********// 205 << 229 // Register standard processes for positrons 206 AddPhysicsList("standard_opt4"); << 230 if (name == "positron-standard") 207 hadronPhys.push_back( new G4DecayPhysi << 231 { 208 hadronPhys.push_back( new G4Radioactiv << 232 if (positronIsRegistered) 209 hadronPhys.push_back( new G4IonBinaryC << 233 { 210 hadronPhys.push_back( new G4EmExtraPhy << 234 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name 211 hadronPhys.push_back( new G4HadronElas << 235 << " cannot be registered ---- positron List already existing" 212 hadronPhys.push_back( new G4StoppingPh << 236 << G4endl; 213 hadronPhys.push_back( new G4HadronPhys << 237 } 214 hadronPhys.push_back( new G4NeutronTra << 238 else 215 << 239 { 216 G4cout << "HADRONTHERAPY_2 PHYSICS LIS << 240 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name 217 << 241 << " is registered" << G4endl; 218 } << 242 RegisterPhysics( new HadrontherapyPositronStandard(name) ); 219 << 243 positronIsRegistered = true; 220 else if (name == "QGSP_BIC"){ << 244 } 221 auto physicsList = new QGSP_BIC; << 222 G4RunManager::GetRunManager() -> SetUs << 223 G4RunManager::GetRunManager() -> Physi << 224 physicsList -> RegisterPhysics(new G4P << 225 } 245 } 226 << 246 // Register penelope processes for positrons 227 else if (name == "QGSP_BERT"){ << 247 if (name == "positron-penelope") 228 auto physicsList = new QGSP_BERT; << 248 { 229 G4RunManager::GetRunManager() -> SetUs << 249 if (positronIsRegistered) 230 G4RunManager::GetRunManager() -> Physi << 250 { 231 physicsList -> RegisterPhysics(new G4P << 251 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 252 << " cannot be registered ---- positron List already existing" << G4endl; >> 253 } >> 254 else >> 255 { >> 256 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 257 << " is registered" << G4endl; >> 258 RegisterPhysics( new HadrontherapyPositronPenelope(name) ); >> 259 positronIsRegistered = true; >> 260 } >> 261 } >> 262 >> 263 // ******** HADRONS AND IONS ********// >> 264 // Register Low Energy processes for protons and ions >> 265 // Stopping power parameterisation: ICRU49 (default model) >> 266 >> 267 if (name == "ion-LowE") >> 268 { >> 269 if (ionIsRegistered) >> 270 { >> 271 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 272 << " cannot be registered ---- proton List already existing" >> 273 << G4endl; >> 274 } >> 275 else >> 276 { >> 277 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 278 << " is registered" << G4endl; >> 279 RegisterPhysics( new HadrontherapyIonLowE(name) ); >> 280 ionIsRegistered = true; >> 281 } 232 } 282 } 233 << 283 234 else if (name == "QGSP_BIC_AllHP"){ << 284 // Register Low Energy processes for protons and ions 235 auto physicsList = new QGSP_BIC_AllHP; << 285 // Stopping power parameterisation: Ziegler 1977 236 G4RunManager::GetRunManager() -> SetUs << 286 if (name == "ion-LowE-ziegler1977") 237 G4RunManager::GetRunManager() -> Physi << 287 { 238 physicsList -> RegisterPhysics(new G4P << 288 if (ionIsRegistered) >> 289 { >> 290 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 291 << " cannot be registered ---- proton List already existing" >> 292 << G4endl; >> 293 } >> 294 else >> 295 { >> 296 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 297 << " is registered" << G4endl; >> 298 RegisterPhysics( new HadrontherapyIonLowEZiegler1977(name) ); >> 299 ionIsRegistered = true; >> 300 } 239 } 301 } 240 << 302 241 else if (name == "QGSP_BIC_HP"){ << 303 242 auto physicsList = new QGSP_BIC_HP; << 304 // Register Low Energy processes for protons and ions 243 G4RunManager::GetRunManager() -> SetUs << 305 // Stopping power parameterisation: Ziegler 1985 244 G4RunManager::GetRunManager() -> Physi << 306 if (name == "ion-LowE-ziegler1985") 245 physicsList -> RegisterPhysics(new G4P << 307 { >> 308 if (ionIsRegistered) >> 309 { >> 310 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 311 << " cannot be registered ---- proton List already existing" >> 312 << G4endl; >> 313 } >> 314 else >> 315 { >> 316 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 317 << " is registered" << G4endl; >> 318 RegisterPhysics( new HadrontherapyIonLowEZiegler1985(name) ); >> 319 ionIsRegistered = true; >> 320 } 246 } 321 } 247 << 322 248 else if (name == "Shielding"){ << 323 249 auto physicsList = new Shielding; << 324 // Register Low Energy processes for protons and ions 250 G4RunManager::GetRunManager() -> SetUs << 325 // Stopping power parameterisation: SRIM2000 251 G4RunManager::GetRunManager() -> Physi << 326 if (name == "ion-LowE-ziegler2000") 252 physicsList -> RegisterPhysics(new G4P << 327 { 253 } << 328 if (ionIsRegistered) 254 << 329 { 255 else if (name == "QBBC"){ << 330 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name 256 auto physicsList = new QBBC; << 331 << " cannot be registered ---- proton List already existing" 257 G4RunManager::GetRunManager() -> SetUs << 332 << G4endl; 258 G4RunManager::GetRunManager() -> Physi << 333 } 259 physicsList -> RegisterPhysics(new G4P << 334 else >> 335 { >> 336 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 337 << " is registered" << G4endl; >> 338 RegisterPhysics( new HadrontherapyIonLowEZiegler2000(name) ); >> 339 ionIsRegistered = true; >> 340 } 260 } 341 } 261 << 342 262 else { << 343 // Register Standard processes for protons and ions 263 G4cout << "PhysicsList::AddPhysicsList << 344 264 << " is not defined" << 345 if (name == "ion-standard") 265 << G4endl; << 346 { >> 347 if (ionIsRegistered) >> 348 { >> 349 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 350 << " cannot be registered ---- ion List already existing" >> 351 << G4endl; >> 352 } >> 353 else >> 354 { >> 355 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 356 << " is registered" << G4endl; >> 357 RegisterPhysics( new HadrontherapyIonStandard(name) ); >> 358 ionIsRegistered = true; >> 359 } 266 } 360 } 267 << 361 >> 362 // Register the Standard processes for muons >> 363 if (name == "muon-standard") >> 364 { >> 365 if (muonIsRegistered) >> 366 { >> 367 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 368 << " cannot be registered ---- decay List already existing" >> 369 << G4endl; >> 370 } >> 371 else >> 372 { >> 373 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 374 << " is registered" << G4endl; >> 375 RegisterPhysics( new HadrontherapyMuonStandard(name) ); >> 376 muonIsRegistered = true; >> 377 } >> 378 } >> 379 >> 380 // Register the decay process >> 381 if (name == "decay") >> 382 { >> 383 if (decayIsRegistered) >> 384 { >> 385 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 386 << " cannot be registered ---- decay List already existing" >> 387 << G4endl; >> 388 } >> 389 else >> 390 { >> 391 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 392 << " is registered" << G4endl; >> 393 RegisterPhysics( new HadrontherapyDecay(name) ); >> 394 decayIsRegistered = true; >> 395 } >> 396 } >> 397 >> 398 >> 399 //-------------------------------------------------------------------------------------------- >> 400 // Begin Hadronic Precompound models >> 401 //-------------------------------------------------------------------------------------------- >> 402 // >> 403 // Register the hadronic physics for protons, neutrons, ions >> 404 // >> 405 >> 406 //Precompound Default Evaporation >> 407 >> 408 if (name == "proton-precompound") >> 409 { >> 410 if (protonHadronicIsRegistered) >> 411 { >> 412 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 413 << " cannot be registered ---- decay List already existing" >> 414 << G4endl; >> 415 } >> 416 else >> 417 { >> 418 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 419 << " is registered" << G4endl; >> 420 RegisterPhysics( new HadrontherapyProtonPrecompound(name) ); >> 421 protonHadronicIsRegistered = true; >> 422 } >> 423 } >> 424 >> 425 >> 426 //Precompond Default Evaporation + Fermi Breck-up Model >> 427 >> 428 if (name == "proton-precompoundFermi") >> 429 { >> 430 if (protonHadronicIsRegistered) >> 431 { >> 432 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 433 << " cannot be registered ---- decay List already existing" >> 434 << G4endl; >> 435 } >> 436 else >> 437 { >> 438 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 439 << " is registered" << G4endl; >> 440 RegisterPhysics( new HadrontherapyProtonPrecompoundFermi(name) ); >> 441 protonHadronicIsRegistered = true; >> 442 } >> 443 } >> 444 >> 445 //Precompound GEM Evaporation >> 446 >> 447 if (name == "proton-precompoundGEM") >> 448 { >> 449 if (protonHadronicIsRegistered) >> 450 { >> 451 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 452 << " cannot be registered ---- decay List already existing" >> 453 << G4endl; >> 454 } >> 455 else >> 456 { >> 457 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 458 << " is registered" << G4endl; >> 459 RegisterPhysics( new HadrontherapyProtonPrecompoundGEM(name) ); >> 460 protonHadronicIsRegistered = true; >> 461 } >> 462 >> 463 } >> 464 >> 465 //Precompound GEM Evaporation + Fermi Breck-up Model >> 466 >> 467 if (name == "proton-precompoundGEMFermi") >> 468 { >> 469 if (protonHadronicIsRegistered) >> 470 { >> 471 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 472 << " cannot be registered ---- decay List already existing" >> 473 << G4endl; >> 474 } >> 475 else >> 476 { >> 477 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 478 << " is registered" << G4endl; >> 479 RegisterPhysics( new HadrontherapyProtonPrecompoundGEMFermi(name) ); >> 480 protonHadronicIsRegistered = true; >> 481 } >> 482 >> 483 } >> 484 >> 485 //------------------------------------------------------------------------------------------------- >> 486 // End Hadronic Precompound models >> 487 //------------------------------------------------------------------------------------------------- >> 488 >> 489 //-------------------------------------------------------------------------------------------- >> 490 //Begin Hadronic Binary models >> 491 //-------------------------------------------------------------------------------------------- >> 492 >> 493 // Binary cascade model with the default precompound >> 494 >> 495 >> 496 if (name == "proton-precompound-binary") >> 497 { >> 498 if (protonHadronicIsRegistered) >> 499 { >> 500 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 501 << " cannot be registered ---- decay List already existing" >> 502 << G4endl; >> 503 } >> 504 else >> 505 { >> 506 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 507 << " is registered" << G4endl; >> 508 RegisterPhysics( new HadrontherapyProtonBinary(name) ); >> 509 protonHadronicIsRegistered = true; >> 510 } >> 511 >> 512 } >> 513 >> 514 //-------------------------------------------------------------------------------------------- >> 515 // Begin Hadronic Bertini model >> 516 //-------------------------------------------------------------------------------------------- >> 517 >> 518 // Bertini model for protons, pions and neutrons >> 519 >> 520 if (name == "proton-bertini") >> 521 { >> 522 if (protonHadronicIsRegistered) >> 523 { >> 524 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 525 << " cannot be registered ---- decay List already existing" >> 526 << G4endl; >> 527 } >> 528 else >> 529 { >> 530 G4cout << "HadrontherapyPhysicsList::AddPhysicsList: " << name >> 531 << " is registered" << G4endl; >> 532 RegisterPhysics( new HadrontherapyProtonBertini(name) ); >> 533 protonHadronicIsRegistered = true; >> 534 } >> 535 >> 536 } >> 537 >> 538 //-------------------------------------------------------------------------------------------- >> 539 // End Hadronic models >> 540 //-------------------------------------------------------------------------------------------- >> 541 >> 542 if (electronIsRegistered && positronIsRegistered && photonIsRegistered && >> 543 ionIsRegistered) >> 544 { >> 545 G4cout << >> 546 "Electromagnetic physics is registered for electron, positron, photons, protons" >> 547 << G4endl; >> 548 } >> 549 >> 550 if (protonHadronicIsRegistered && muonIsRegistered && decayIsRegistered) >> 551 { >> 552 G4cout << " Hadronic physics is registered" << G4endl; >> 553 } 268 } 554 } 269 555 270 ////////////////////////////////////////////// << 556 void HadrontherapyPhysicsList::SetCuts() 271 void HadrontherapyPhysicsList::AddStepMax() << 557 { 272 { << 558 // Set the threshold of production equal to the defaultCutValue 273 // Step limitation seen as a process << 559 // in the experimental set-up 274 // This process must exist in all threads. << 560 G4VUserPhysicsList::SetCutsWithDefault(); 275 // << 276 HadrontherapyStepMax* stepMaxProcess = ne << 277 561 278 << 562 // Definition of a smaller threshold of production in the phantom region 279 auto particleIterator = GetParticleIterato << 563 // where high accuracy is required in the energy deposit calculation 280 particleIterator->reset(); << 564 281 while ((*particleIterator)()){ << 565 G4String regionName = "PhantomLog"; 282 G4ParticleDefinition* particle = parti << 566 G4Region* region = G4RegionStore::GetInstance()->GetRegion(regionName); 283 G4ProcessManager* pmanager = particle- << 567 G4ProductionCuts* cuts = new G4ProductionCuts ; 284 << 568 G4double regionCut = 0.001*mm; 285 if (stepMaxProcess->IsApplicable(*part << 569 cuts -> SetProductionCut(regionCut,G4ProductionCuts::GetIndex("gamma")); 286 { << 570 cuts -> SetProductionCut(regionCut,G4ProductionCuts::GetIndex("e-")); 287 pmanager ->AddDiscreteProcess(step << 571 cuts -> SetProductionCut(regionCut,G4ProductionCuts::GetIndex("e+")); 288 } << 572 cuts -> SetProductionCut(regionCut,G4ProductionCuts::GetIndex("proton")); 289 } << 573 cuts -> SetProductionCut(regionCut,G4ProductionCuts::GetIndex("genericIons")); >> 574 region -> SetProductionCuts(cuts); >> 575 >> 576 if (verboseLevel>0) DumpCutValuesTable(); 290 } 577 } >> 578 >> 579 291 580