| 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 | // **************************************
|
|---|
| 27 | // * *
|
|---|
| 28 | // * RemSimHadronicBertini.cc *
|
|---|
| 29 | // * *
|
|---|
| 30 | // **************************************
|
|---|
| 31 | //
|
|---|
| 32 | // $Id: RemSimHadronicBertini.cc,v 1.5 2006/11/15 18:39:30 guatelli Exp $
|
|---|
| 33 | // GEANT4 tag $Name: $
|
|---|
| 34 | //
|
|---|
| 35 | // Author : Susanna Guatelli, guatelli@ge.infn.it
|
|---|
| 36 |
|
|---|
| 37 | // Code review: MGP, 7 November 2006 (to be completed)
|
|---|
| 38 | //
|
|---|
| 39 | #include "RemSimHadronicBertini.hh"
|
|---|
| 40 | #include "G4BinaryLightIonReaction.hh"
|
|---|
| 41 | #include "G4TripathiCrossSection.hh"
|
|---|
| 42 | #include "G4IonsShenCrossSection.hh"
|
|---|
| 43 | #include "G4ParticleDefinition.hh"
|
|---|
| 44 | #include "G4ProcessManager.hh"
|
|---|
| 45 | #include "G4LElastic.hh"
|
|---|
| 46 | #include "G4CascadeInterface.hh"
|
|---|
| 47 | //#include "G4PreCompoundModel.hh"
|
|---|
| 48 | //#include "G4QGSModel.hh"
|
|---|
| 49 | #include "G4LEProtonInelastic.hh"
|
|---|
| 50 | #include "G4LENeutronInelastic.hh"
|
|---|
| 51 | #include "G4LEPionPlusInelastic.hh"
|
|---|
| 52 | #include "G4LEPionMinusInelastic.hh"
|
|---|
| 53 | #include "G4LEAlphaInelastic.hh"
|
|---|
| 54 | #include "G4LFission.hh"
|
|---|
| 55 | #include "G4LCapture.hh"
|
|---|
| 56 | //#include "G4TheoFSGenerator.hh"
|
|---|
| 57 | //#include "G4GeneratorPrecompoundInterface.hh"
|
|---|
| 58 | //#include "G4PreCompoundModel.hh"
|
|---|
| 59 | //#include "G4QGSMFragmentation.hh"
|
|---|
| 60 | //#include "G4ExcitedStringDecay.hh"
|
|---|
| 61 | #include "G4HadronElasticProcess.hh"
|
|---|
| 62 | #include "G4HadronFissionProcess.hh"
|
|---|
| 63 | #include "G4HadronCaptureProcess.hh"
|
|---|
| 64 | #include "G4ProtonInelasticProcess.hh"
|
|---|
| 65 | #include "G4NeutronInelasticProcess.hh"
|
|---|
| 66 | #include "G4PionPlusInelasticProcess.hh"
|
|---|
| 67 | #include "G4PionMinusInelasticProcess.hh"
|
|---|
| 68 | #include "G4AlphaInelasticProcess.hh"
|
|---|
| 69 | #include "G4ExcitationHandler.hh"
|
|---|
| 70 |
|
|---|
| 71 | //
|
|---|
| 72 | // BERTINI PHYSICS LIST
|
|---|
| 73 | //
|
|---|
| 74 |
|
|---|
| 75 | RemSimHadronicBertini::RemSimHadronicBertini(const G4String& name):
|
|---|
| 76 | G4VPhysicsConstructor(name)
|
|---|
| 77 | {}
|
|---|
| 78 |
|
|---|
| 79 | RemSimHadronicBertini::~RemSimHadronicBertini()
|
|---|
| 80 | {}
|
|---|
| 81 |
|
|---|
| 82 | void RemSimHadronicBertini::ConstructProcess()
|
|---|
| 83 | {
|
|---|
| 84 | // Physics for proton, neutron, pion+ and pion-
|
|---|
| 85 |
|
|---|
| 86 | // Elastic scattering: LElastic model
|
|---|
| 87 | G4LElastic* elasticModel = new G4LElastic();
|
|---|
| 88 | G4HadronElasticProcess* elasticScattering = new G4HadronElasticProcess();
|
|---|
| 89 | elasticScattering->RegisterMe(elasticModel);
|
|---|
| 90 |
|
|---|
| 91 | // Inelastic scattering: Bertini model up to 3.2 GeV
|
|---|
| 92 | G4CascadeInterface* bertiniModel = new G4CascadeInterface();
|
|---|
| 93 | // Energy limit of the Bertini model
|
|---|
| 94 | G4double bertiniHighEnergyLimit = 3.2 * GeV;
|
|---|
| 95 | bertiniModel->SetMaxEnergy(bertiniHighEnergyLimit);
|
|---|
| 96 |
|
|---|
| 97 | // Inelastic scattering: LEP model between 2.8 * GeV and 25. * GeV
|
|---|
| 98 | G4LEProtonInelastic* LEPProtonModel = new G4LEProtonInelastic();
|
|---|
| 99 | G4LENeutronInelastic* LEPNeutronModel = new G4LENeutronInelastic();
|
|---|
| 100 | G4LEPionPlusInelastic* LEPPionPlusModel = new G4LEPionPlusInelastic();
|
|---|
| 101 | G4LEPionMinusInelastic* LEPPionMinusModel = new G4LEPionMinusInelastic();
|
|---|
| 102 | // Set the energy limits
|
|---|
| 103 | G4double LEPLowEnergyLimit = 2.8 * GeV;
|
|---|
| 104 | G4double LEPHighEnergyLimit = 25. * GeV;
|
|---|
| 105 | LEPProtonModel->SetMinEnergy(LEPLowEnergyLimit);
|
|---|
| 106 | LEPProtonModel->SetMaxEnergy(LEPHighEnergyLimit);
|
|---|
| 107 | LEPNeutronModel->SetMinEnergy(LEPLowEnergyLimit);
|
|---|
| 108 | LEPNeutronModel->SetMaxEnergy(LEPHighEnergyLimit);
|
|---|
| 109 | LEPPionPlusModel->SetMinEnergy(LEPLowEnergyLimit);
|
|---|
| 110 | LEPPionPlusModel->SetMaxEnergy(LEPHighEnergyLimit);
|
|---|
| 111 | LEPPionMinusModel->SetMinEnergy(LEPLowEnergyLimit);
|
|---|
| 112 | LEPPionMinusModel->SetMaxEnergy(LEPHighEnergyLimit);
|
|---|
| 113 |
|
|---|
| 114 | // Inelastic scattering: QGSP model between 20 GeV and 100 TeV
|
|---|
| 115 | QGSPModel = new G4TheoFSGenerator();
|
|---|
| 116 | // Set the energy limits
|
|---|
| 117 | G4double QGSPLowEnergyLimit = 20.* GeV;
|
|---|
| 118 | G4double QGSPHighEnergyLimit = 100.* GeV;
|
|---|
| 119 | QGSPModel->SetMinEnergy(QGSPLowEnergyLimit);
|
|---|
| 120 | QGSPModel->SetMaxEnergy(QGSPHighEnergyLimit);
|
|---|
| 121 |
|
|---|
| 122 | theCascade = new G4GeneratorPrecompoundInterface();
|
|---|
| 123 |
|
|---|
| 124 | thePreEquilib = new G4PreCompoundModel(&theHandler);
|
|---|
| 125 | theCascade->SetDeExcitation(thePreEquilib);
|
|---|
| 126 | QGSPModel->SetTransport(theCascade);
|
|---|
| 127 |
|
|---|
| 128 | // Set the fragmentation
|
|---|
| 129 | theFragmentation = new G4QGSMFragmentation();
|
|---|
| 130 | theStringDecay = new G4ExcitedStringDecay(theFragmentation);
|
|---|
| 131 | theStringModel.SetFragmentationModel(theStringDecay);
|
|---|
| 132 | QGSPModel->SetHighEnergyGenerator(&theStringModel);
|
|---|
| 133 |
|
|---|
| 134 | // ---------------------------------------------------------------------------------------------
|
|---|
| 135 | // Proton processes
|
|---|
| 136 | G4ParticleDefinition* proton = G4Proton::ProtonDefinition();
|
|---|
| 137 | G4ProcessManager* protonProcessManager = proton->GetProcessManager();
|
|---|
| 138 |
|
|---|
| 139 | // Proton inelastic scattering
|
|---|
| 140 | G4ProtonInelasticProcess* protonInelasticProcess = new G4ProtonInelasticProcess();
|
|---|
| 141 | // Activate the cross-sections for proton nuclear scattering up to 20 GeV
|
|---|
| 142 | //G4ProtonInelasticCrossSection protonCrossSection;
|
|---|
| 143 | protonInelasticProcess->AddDataSet(&protonCrossSection);
|
|---|
| 144 | // Set the models
|
|---|
| 145 | protonInelasticProcess->RegisterMe(bertiniModel);
|
|---|
| 146 | protonInelasticProcess->RegisterMe(LEPProtonModel);
|
|---|
| 147 |
|
|---|
| 148 | protonInelasticProcess->RegisterMe(QGSPModel);
|
|---|
| 149 | // Activate the inelastic scattering
|
|---|
| 150 | protonProcessManager->AddDiscreteProcess(protonInelasticProcess);
|
|---|
| 151 | // Activate the elastic scattering
|
|---|
| 152 | protonProcessManager->AddDiscreteProcess(elasticScattering);
|
|---|
| 153 |
|
|---|
| 154 | //------------------------------------------------------
|
|---|
| 155 | // Pion Plus processes
|
|---|
| 156 | G4ParticleDefinition* piPlus = G4PionPlus::PionPlusDefinition();
|
|---|
| 157 | G4ProcessManager* pionPlusProcessManager = piPlus->GetProcessManager();
|
|---|
| 158 |
|
|---|
| 159 | // Define the inelastic scattering for pion plus
|
|---|
| 160 | G4PionPlusInelasticProcess* pionPlusInelasticProcess = new G4PionPlusInelasticProcess();
|
|---|
| 161 | // Set the cross section
|
|---|
| 162 | // G4PiNuclearCrossSection pionCrossSection;
|
|---|
| 163 | pionPlusInelasticProcess->AddDataSet(&pionCrossSection);
|
|---|
| 164 | // Register the models
|
|---|
| 165 | pionPlusInelasticProcess->RegisterMe(bertiniModel);
|
|---|
| 166 | pionPlusInelasticProcess->RegisterMe(LEPPionPlusModel);
|
|---|
| 167 | pionPlusInelasticProcess->RegisterMe(QGSPModel);
|
|---|
| 168 | // Activate the inelastic scattering
|
|---|
| 169 | pionPlusProcessManager->AddDiscreteProcess(pionPlusInelasticProcess);
|
|---|
| 170 | // Activate the elastic process
|
|---|
| 171 | pionPlusProcessManager->AddDiscreteProcess(elasticScattering);
|
|---|
| 172 |
|
|---|
| 173 | //------------------------------------------------------------
|
|---|
| 174 | // Pion Minus processes
|
|---|
| 175 | G4ParticleDefinition* piMinus = G4PionMinus::PionMinusDefinition();
|
|---|
| 176 | G4ProcessManager* pionMinusProcessManager = piMinus->GetProcessManager();
|
|---|
| 177 |
|
|---|
| 178 | // Define the inelastic processes for pion minus
|
|---|
| 179 | G4PionMinusInelasticProcess* pionMinusInelasticProcess = new G4PionMinusInelasticProcess();
|
|---|
| 180 | // Set the cross section
|
|---|
| 181 | pionMinusInelasticProcess->AddDataSet(&pionCrossSection);
|
|---|
| 182 | // Register the models
|
|---|
| 183 | pionMinusInelasticProcess->RegisterMe(bertiniModel);
|
|---|
| 184 | pionMinusInelasticProcess->RegisterMe(LEPPionMinusModel);
|
|---|
| 185 | pionMinusInelasticProcess->RegisterMe(QGSPModel);
|
|---|
| 186 | // Activate the inelastic scattering
|
|---|
| 187 | pionMinusProcessManager->AddDiscreteProcess(pionMinusInelasticProcess);
|
|---|
| 188 | // Activate the elastic scattering
|
|---|
| 189 | pionMinusProcessManager->AddDiscreteProcess(elasticScattering);
|
|---|
| 190 |
|
|---|
| 191 | //-----------------------------------------------------
|
|---|
| 192 | // Neutron processes
|
|---|
| 193 | G4ParticleDefinition* neutron = G4Neutron::NeutronDefinition();
|
|---|
| 194 | G4ProcessManager* neutronProcessManager = neutron->GetProcessManager();
|
|---|
| 195 |
|
|---|
| 196 | // Inelastic process
|
|---|
| 197 | G4NeutronInelasticProcess* neutronInelasticProcess = new G4NeutronInelasticProcess();
|
|---|
| 198 | // Set the cross section
|
|---|
| 199 | //G4NeutronInelasticCrossSection neutronCrossSection;
|
|---|
| 200 | neutronInelasticProcess->AddDataSet(&neutronCrossSection);
|
|---|
| 201 | // Set the models
|
|---|
| 202 | neutronInelasticProcess->RegisterMe(bertiniModel);
|
|---|
| 203 | neutronInelasticProcess->RegisterMe(LEPNeutronModel);
|
|---|
| 204 | neutronInelasticProcess->RegisterMe(QGSPModel);
|
|---|
| 205 | // Activate the neutron inelastic scattering
|
|---|
| 206 | neutronProcessManager->AddDiscreteProcess(neutronInelasticProcess);
|
|---|
| 207 | // Activate the neutron elastic scattering
|
|---|
| 208 | neutronProcessManager->AddDiscreteProcess(elasticScattering);
|
|---|
| 209 |
|
|---|
| 210 | // Neutron capture process
|
|---|
| 211 | G4HadronCaptureProcess* neutronCapture = new G4HadronCaptureProcess();
|
|---|
| 212 | // Final state production model for capture of neutral hadrons in nuclei
|
|---|
| 213 | G4LCapture* captureModel = new G4LCapture();
|
|---|
| 214 | // Set the energy range for the capture model
|
|---|
| 215 | G4double neutronLowEnergyLimit = 0. * MeV;
|
|---|
| 216 | G4double neutronHighEnergyLimit = 100. * TeV;
|
|---|
| 217 | captureModel->SetMinEnergy(neutronLowEnergyLimit);
|
|---|
| 218 | captureModel->SetMaxEnergy(neutronHighEnergyLimit);
|
|---|
| 219 | // Activate the neutron capture model
|
|---|
| 220 | neutronCapture->RegisterMe(captureModel);
|
|---|
| 221 | // Activate the neutron capture process
|
|---|
| 222 | neutronProcessManager->AddDiscreteProcess(neutronCapture);
|
|---|
| 223 |
|
|---|
| 224 | // Process for induced fission
|
|---|
| 225 | G4HadronFissionProcess* fission = new G4HadronFissionProcess();
|
|---|
| 226 | //Final state production model for induced fission
|
|---|
| 227 | G4LFission* fissionModel = new G4LFission();
|
|---|
| 228 | // Set the energy range for the fission model
|
|---|
| 229 | fissionModel->SetMinEnergy(neutronLowEnergyLimit);
|
|---|
| 230 | fissionModel->SetMaxEnergy(neutronHighEnergyLimit);
|
|---|
| 231 | // Register the fission model
|
|---|
| 232 | fission->RegisterMe(fissionModel);
|
|---|
| 233 | // Activate the fission process
|
|---|
| 234 | neutronProcessManager->AddDiscreteProcess(fission);
|
|---|
| 235 |
|
|---|
| 236 | //--------------------------------------------------------
|
|---|
| 237 | // Physics for alpha particles
|
|---|
| 238 |
|
|---|
| 239 | G4ParticleDefinition* alpha = G4Alpha::AlphaDefinition();
|
|---|
| 240 | G4ProcessManager* alphaProcessManager = alpha->GetProcessManager();
|
|---|
| 241 |
|
|---|
| 242 | // Cross section data sets
|
|---|
| 243 |
|
|---|
| 244 | // TRIPATHI CROSS SECTION
|
|---|
| 245 | // Implementation of formulas taken from NASA technical paper 3621 by
|
|---|
| 246 | // Tripathi, et al. Cross-sections for ion ion scattering
|
|---|
| 247 | G4TripathiCrossSection* tripathi = new G4TripathiCrossSection();
|
|---|
| 248 |
|
|---|
| 249 | // IONS SHEN CROSS SECTION
|
|---|
| 250 | // Implementation of formulas
|
|---|
| 251 | // Shen et al. Nuc. Phys. A 491 130 (1989)
|
|---|
| 252 | // Total Reaction Cross Section for Heavy-Ion Collisions
|
|---|
| 253 | G4IonsShenCrossSection* shen = new G4IonsShenCrossSection();
|
|---|
| 254 |
|
|---|
| 255 | G4LEAlphaInelastic* LEPAlphaModel = new G4LEAlphaInelastic();
|
|---|
| 256 | // Energy limit of the LEP model for alpha particles
|
|---|
| 257 | G4double LEPAlphaHighLimit = 100 * MeV;
|
|---|
| 258 | LEPAlphaModel->SetMaxEnergy(LEPAlphaHighLimit);
|
|---|
| 259 |
|
|---|
| 260 | G4BinaryLightIonReaction* binaryIonModel = new G4BinaryLightIonReaction();
|
|---|
| 261 | // Energy limit of the binary ion model
|
|---|
| 262 | G4double binaryIonLowLimit = 80. * MeV;
|
|---|
| 263 | G4double binaryIonHighLimit = 40. * GeV;
|
|---|
| 264 | binaryIonModel->SetMinEnergy(binaryIonLowLimit);
|
|---|
| 265 | binaryIonModel->SetMaxEnergy(binaryIonHighLimit);
|
|---|
| 266 |
|
|---|
| 267 | // Define the alpha inelastic scattering
|
|---|
| 268 | G4AlphaInelasticProcess* alphaInelasticProcess = new G4AlphaInelasticProcess();
|
|---|
| 269 | // Activate the Tripathi and Shen Cross Section
|
|---|
| 270 | alphaInelasticProcess->AddDataSet(tripathi);
|
|---|
| 271 | alphaInelasticProcess->AddDataSet(shen);
|
|---|
| 272 | // Set the models
|
|---|
| 273 | alphaInelasticProcess->RegisterMe(LEPAlphaModel);
|
|---|
| 274 | alphaInelasticProcess->RegisterMe(binaryIonModel);
|
|---|
| 275 | // Activate the inelastic scattering
|
|---|
| 276 | alphaProcessManager->AddDiscreteProcess(alphaInelasticProcess);
|
|---|
| 277 | // Activate the elastic scattering
|
|---|
| 278 | alphaProcessManager->AddDiscreteProcess(elasticScattering);
|
|---|
| 279 | }
|
|---|
| 280 |
|
|---|
| 281 |
|
|---|
| 282 |
|
|---|
