| Geant4 Cross Reference |
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 // GEANT4 Class file
29 //
30 // For information related to this code contact:
31 //
32 // File name: G4XnpElasticLowE
33 //
34 // Author:
35 //
36 // Creation date: 15 April 1999
37 //
38 // Modifications:
39 //
40 // Neutron-Proton elastic cross section
41 // Linear interpolation from data table
42 //
43 // -------------------------------------------------------------------
44
45 #include "globals.hh"
46 #include "G4ios.hh"
47 #include "G4Exp.hh"
48 #include "G4Log.hh"
49 #include "G4SystemOfUnits.hh"
50 #include "G4XnpElasticLowE.hh"
51 #include "G4KineticTrack.hh"
52 #include "G4ParticleDefinition.hh"
53 #include "G4PhysicsVector.hh"
54 #include "G4PhysicsLogVector.hh"
55 #include "G4Proton.hh"
56 #include "G4Neutron.hh"
57
58 const G4double G4XnpElasticLowE::_lowLimit = 0.;
59 const G4double G4XnpElasticLowE::_highLimit = 3.*GeV;
60
61 // Low energy limit of the cross-section table (in GeV)
62 // Units are assigned when filling the PhysicsVector
63 const G4double G4XnpElasticLowE::_eMinTable = 1.8964808;
64 const G4double G4XnpElasticLowE::_eStepLog = 0.01;
65
66 // Cross-sections in mb
67 // Units are assigned when filling the PhysicsVector
68 const G4int G4XnpElasticLowE::_tableSize = 101;
69 const G4double G4XnpElasticLowE::_sigmaTable[101] =
70 {
71 1500.0, //GF was 0.0
72 248.20, 93.38, 55.26, 44.50, 41.33, 38.48, 37.20, 35.98,
73 35.02, 34.47, 32.48, 30.76, 29.46, 28.53, 27.84, 27.20,
74 26.53, 25.95, 25.59, 25.46, 25.00, 24.49, 24.08, 23.86,
75 23.17, 22.70, 21.88, 21.48, 20.22, 19.75, 18.97, 18.39,
76 17.98, 17.63, 17.21, 16.72, 16.68, 16.58, 16.42, 16.22,
77 15.98, 15.71, 15.42, 15.14, 14.87, 14.65, 14.44, 14.26,
78 14.10, 13.95, 13.80, 13.64, 13.47, 13.29, 13.09, 12.89,
79 12.68, 12.47, 12.27, 12.06, 11.84, 11.76, 11.69, 11.60,
80 11.50, 11.41, 11.29, 11.17, 11.06, 10.93, 10.81, 10.68,
81 10.56, 10.44, 10.33, 10.21, 10.12, 10.03, 9.96, 9.89,
82 9.83, 9.80, 9.77, 9.75, 9.74, 9.74, 9.74, 9.76,
83 9.73, 9.70, 9.68, 9.65, 9.63, 9.60, 9.57, 9.55,
84 9.52, 9.49, 9.46, 9.43
85 };
86
87
88 G4XnpElasticLowE::G4XnpElasticLowE()
89 {
90 // Cross-sections are available in the range (_eMin,_eMax)
91
92 _eMin = _eMinTable * GeV;
93 _eMin = G4Exp(G4Log(_eMinTable)-_eStepLog)*GeV;
94 _eMax = G4Exp(G4Log(_eMinTable) + _tableSize * _eStepLog) * GeV;
95
96 // Protections: validity limits must be compatible with available data
97
98 if (_eMin < _lowLimit)
99 throw G4HadronicException(__FILE__, __LINE__, "G4XnpElasticLowE::G4XnpElasticLowE - Low energy limit not valid");
100
101 if (_highLimit > _eMax)
102 throw G4HadronicException(__FILE__, __LINE__, "G4XnpElasticLowE::G4XnpElasticLowE - High energy limit not valid");
103
104 _sigma = new G4PhysicsLogVector(_eMin,_eMax,_tableSize);
105 G4int i;
106 for (i=0; i<_tableSize; i++)
107 {
108 G4double value = _sigmaTable[i] * millibarn;
109 _sigma->PutValue(i,value);
110 }
111 }
112
113
114 G4XnpElasticLowE::~G4XnpElasticLowE()
115 {
116 delete _sigma;
117 _sigma = 0;
118 }
119
120
121 G4bool G4XnpElasticLowE::operator==(const G4XnpElasticLowE &right) const
122 {
123 return (this == (G4XnpElasticLowE *) &right);
124 }
125
126
127 G4bool G4XnpElasticLowE::operator!=(const G4XnpElasticLowE &right) const
128 {
129 return (this != (G4XnpElasticLowE *) &right);
130 }
131
132
133 G4double G4XnpElasticLowE::CrossSection(const G4KineticTrack& trk1, const G4KineticTrack& trk2) const
134 {
135 G4double sigma = 0.;
136 G4double sqrtS = (trk1.Get4Momentum() + trk2.Get4Momentum()).mag();
137 G4bool dummy = false;
138
139 const G4ParticleDefinition* proton = G4Proton::ProtonDefinition();
140 const G4ParticleDefinition* neutron = G4Neutron::NeutronDefinition();
141
142 const G4ParticleDefinition* def1 = trk1.GetDefinition();
143 const G4ParticleDefinition* def2 = trk2.GetDefinition();
144 if ( (def1 == proton && def2 == neutron) ||
145 (def1 == neutron && def2 == proton) )
146 {
147 if (sqrtS >= _eMin && sqrtS <= _eMax)
148 {
149 sigma = _sigma->GetValue(sqrtS,dummy);
150 } else if ( sqrtS < _eMin )
151 {
152 sigma = _sigma->GetValue(_eMin,dummy);
153 }
154 }
155
156 return sigma;
157 }
158
159 void G4XnpElasticLowE::Print() const
160 {
161 // Dump the cross-section table
162 G4cout << Name() << "Cross-section table: " << G4endl;
163 G4bool dummy = false;
164 G4int i;
165
166 for (i=0; i<_tableSize; i++)
167 {
168 G4double e = _sigma->GetLowEdgeEnergy(i) / GeV;
169 G4double sigma = _sigma->GetValue(e,dummy) / millibarn;
170 G4cout << i << ") e = " << e << " GeV ---- Cross section = " << sigma << " mb " << G4endl;
171 }
172
173 G4VCrossSectionSource::Print();
174 }
175
176
177 G4String G4XnpElasticLowE::Name() const
178 {
179 G4String name("npElasticLowE");
180 return name;
181 }
182
183
184 G4bool G4XnpElasticLowE::IsValid(G4double e) const
185 {
186 G4bool answer = InLimits(e,_lowLimit,_highLimit);
187
188 return answer;
189 }
190
191
192