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

Last change on this file since 1321 was 1230, checked in by garnier, 14 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.