source: trunk/examples/advanced/radioprotection/src/RemSimHadronicBinary.cc@ 1330

Last change on this file since 1330 was 1230, checked in by garnier, 16 years ago

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