Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4NucleonNuclearCrossSection.cc
Go to the documentation of this file.
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 // author: Vladimir.Grichine@cern.ch
27 //
28 // Implements data from: Barashenkov V.S., Nucleon-Nucleus Cross Section,
29 // Preprint JINR P2-89-770, p. 12, Dubna 1989 (scanned version from KEK)
30 // Based on G. Folger version of G4PiNuclearCrossSection class
31 //
32 // Modified:
33 // 05.03.07 V.Ivanchenko - add IfZAApplicable, remove "debug"
34 // 11.09.09 V.Ivanchenko - fixed bug in interpolation
35 //
36 
38 #include "G4SystemOfUnits.hh"
39 #include "G4DynamicParticle.hh"
40 #include "G4Neutron.hh"
41 #include "G4Proton.hh"
42 #include "G4Pow.hh"
43 
44 // factory
45 #include "G4CrossSectionFactory.hh"
46 //
48 
49 
50 // Group 1: He, Be, C for 44 energies
51 
52 const G4double G4NucleonNuclearCrossSection::e1[44] =
53 {
54  0.014, 0.015, 0.017, 0.02, 0.022, 0.025, 0.027, 0.03, 0.035, 0.04,
55  0.045, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1, 0.12, 0.14, 0.15,
56  0.16, 0.18, 0.20, 0.25, 0.30, 0.35, 0.4, 0.5, 0.6, 0.7,
57  0.8, 0.9, 1.0, 1.5, 2.0, 3.0, 5.0, 7.0, 10.0, 20.0,
58  50.0, 100.0, 500.0, 1000.0
59 };
60 
61 const G4double G4NucleonNuclearCrossSection::he_m_t[44] =
62 {
63  1090, 1020, 915, 800, 710, 640, 600, 560, 500, 440, 390, 360, 295, 256, 220, 192,
64  168, 136, 120, 116, 114, 110, 107, 104, 106, 108, 110, 120, 126, 135, 140, 144, 146,
65  148, 152, 150, 146, 142, 138, 132, 129, 126, 127, 128
66 };
67 const G4double G4NucleonNuclearCrossSection::he_m_in[44] =
68 {
69  0, 5, 10, 20, 35, 55, 70, 80, 90, 105, 115, 115, 100, 90, 86, 84, 84, 82, 80, 80, 80, 80,
70  79, 78, 80, 84, 88, 94, 100, 105, 108, 108, 108, 112, 114, 114, 112, 110, 108, 106, 104,
71  101, 102, 102
72 };
73 const G4double G4NucleonNuclearCrossSection::he_p_in[44] =
74 {
75  0, 2, 3, 13, 30, 50, 65, 77, 90, 105, 115, 115, 100, 90, 86, 84, 84, 82, 80, 80, 80, 80,
76  79, 78, 80, 84, 88, 94, 100, 105, 108, 108, 108, 112, 114, 114, 112, 110, 108, 106, 104,
77  101, 102, 102
78 };
79 
80 const G4double G4NucleonNuclearCrossSection::be_m_t[44] =
81 {
82  1490, 1460, 1400, 1350, 1270, 1200, 1160, 1100, 1000, 910, 810, 740, 625, 575, 455, 406,
83  365, 310, 275, 262, 255, 240, 235, 225, 225, 230, 238, 252, 270, 282, 288, 290, 294, 303,
84  303, 300, 292, 284, 277, 267, 263, 264, 268, 268
85 };
86 const G4double G4NucleonNuclearCrossSection::be_m_in[44] =
87 {
88  650, 640, 617, 595, 555, 520, 495, 470, 430, 385, 350, 320, 270, 250, 210, 190, 185, 178,
89  175, 175, 175, 175, 175, 170, 170, 172, 176, 184, 194, 200, 209, 213, 214, 216, 216, 212,
90  210, 210, 210, 210, 210, 210, 210, 210
91 };
92 const G4double G4NucleonNuclearCrossSection::be_p_in[44] =
93 {
94  490, 540, 580, 545, 525, 495, 470, 450, 420, 370, 340, 310, 262, 242, 205, 185, 180, 175,
95  172, 175, 175, 175, 175, 170, 170, 172, 176, 184, 194, 200, 209, 213, 214, 216, 216, 212,
96  210, 210, 210, 210, 210, 210, 210, 210
97 };
98 
99 const G4double G4NucleonNuclearCrossSection::c_m_t[44] =
100 {
101  1240, 1370, 1450, 1455, 1445, 1385, 1345, 1290, 1210, 1110, 1020, 940, 800, 700, 604, 530,
102  475, 396, 350, 336, 320, 303, 294, 280, 280, 286, 296, 314, 330, 344, 356, 360, 364, 384,
103  388, 384, 364, 352, 344, 330, 324, 324, 332, 332
104 };
105 const G4double G4NucleonNuclearCrossSection::c_m_in[44] =
106 {
107  590, 570, 542, 510, 500, 460, 445, 430, 395, 380, 350, 330, 295, 270, 255, 240, 228, 222,
108  216, 216, 210, 210, 210, 208, 210, 214, 216, 228, 240, 248, 254, 257, 260, 262, 260, 256,
109  252, 252, 250, 250, 248, 248, 248, 248
110 };
111 const G4double G4NucleonNuclearCrossSection::c_p_in[44] =
112 {
113  310, 330, 400, 440, 450, 435, 430, 420, 385, 370, 340, 320, 288, 263, 249, 234, 222, 216,
114  210, 211, 205, 208, 210, 208, 210, 214, 216, 228, 240, 248, 254, 257, 260, 262, 260, 256,
115  252, 252, 250, 250, 248, 248, 248, 248
116 };
117 
118 // Group 2: N, O, Na for 44 energies (e1=e2)
119 
120 const G4double G4NucleonNuclearCrossSection::e2[44] =
121 {
122  0.014, 0.015, 0.017, .02, 0.022, 0.025, 0.027, 0.03, 0.035, .04, 0.045, 0.05, .06, 0.07,
123  .08, 0.09, .1, .12, .14, .15, .16, .18, .20, .25, .30, .35, .4 , 0.5, 0.6, 0.7, 0.8,
124  0.9, 1, 1.5, 2, 3, 5, 7, 10,
125  20, 50, 100, 500, 1000
126 };
127 
128 const G4double G4NucleonNuclearCrossSection::n_m_t[44] =
129 {
130  1420,1480, 1537, 1550, 1525, 1500, 1480, 1425, 1340, 1260, 1175, 1090, 930, 805, 690, 612,
131  552, 462, 402, 384, 372, 350, 345, 326, 324, 328, 336, 356, 372, 388, 400, 408, 415, 430,
132  435, 432, 415, 402, 390, 375, 367, 370, 382, 385
133 };
134 const G4double G4NucleonNuclearCrossSection::n_m_in[44] =
135 {
136  680, 665, 625, 580, 562, 525, 510, 485, 450, 435, 410, 387, 340, 310, 290, 280, 276, 274,
137  260, 258, 254, 247, 245, 240, 240, 244, 250, 260, 268, 275, 280, 285, 290, 295, 300, 294,
138  292, 290, 285, 285, 282, 282, 282, 282
139 };
140 const G4double G4NucleonNuclearCrossSection::n_p_in[44] =
141 {
142  420, 440, 470, 490, 497, 500, 480, 462, 440, 425, 400, 377, 333, 303, 284, 274, 270, 268,
143  254, 252, 247, 245, 245, 240, 240, 244, 250, 260, 268, 275, 280, 285, 290, 295, 300, 294,
144  292, 290, 285, 285, 282, 282, 282, 282
145 };
146 
147 const G4double G4NucleonNuclearCrossSection::o_m_t[44] =
148 {
149  1520, 1570, 1630, 1660, 1647, 1623, 1595, 1555, 1475, 1395, 1290, 1207, 1035, 925, 816,
150  720, 645, 540, 462, 438, 415, 392, 378, 362, 361, 381, 390, 403, 417, 440, 460, 470,
151  479, 498, 504, 498, 477, 457, 443, 427, 420, 425, 429, 430
152 };
153 const G4double G4NucleonNuclearCrossSection::o_m_in[44] =
154 {
155  750, 740, 700, 650, 620, 575, 555, 530, 505, 462, 435, 420, 375, 345, 320, 310, 300, 293,
156  288, 282, 282, 280, 276, 270, 271, 275, 280, 290, 295, 304, 310, 315, 318, 332, 335, 330,
157  323, 320, 317, 315, 315, 315, 315, 315
158 };
159 const G4double G4NucleonNuclearCrossSection::o_p_in[44] =
160 {
161  460, 485, 510, 535, 537, 532, 520, 500, 460, 432, 405, 390, 350, 320, 310, 304, 293, 287,
162  283, 279, 279, 278, 276, 270, 271, 275, 280, 290, 295, 304, 310, 315, 318, 332, 335, 330,
163  323, 320, 317, 315, 315, 315, 315, 315
164 };
165 
166 const G4double G4NucleonNuclearCrossSection::na_m_t[44] =
167 {
168  1570, 1620, 1695, 1730, 1750, 1760, 1755, 1740, 1710, 1643, 1560, 1480, 1343, 1220, 1073,
169  953, 860, 720, 618, 582, 546, 522, 504, 484, 492, 500, 512, 538, 560, 586, 608, 622, 632,
170  660, 668, 664, 640, 616, 596, 568, 568, 568, 568, 568
171 };
172 const G4double G4NucleonNuclearCrossSection::na_m_in[44] =
173 {
174  960, 930, 890, 822, 790, 750, 725, 686, 620, 600, 575, 540, 497, 450, 414, 390, 380, 372,
175  354, 360, 355, 354, 350, 350, 350, 356, 364, 384, 392, 400, 408, 410, 420, 408, 412, 420,
176  411, 409, 407, 403, 400, 400, 400, 400
177 };
178 const G4double G4NucleonNuclearCrossSection::na_p_in[44] =
179 {
180  600, 617, 660, 675, 680, 680, 670, 650, 575, 550, 525, 490, 450, 420, 385, 367, 360, 350,
181  350, 350, 345, 347, 350, 350, 350, 356, 364, 384, 392, 400, 408, 410, 420, 408, 412, 420,
182  411, 409, 407, 403, 400, 400, 400, 400
183 };
184 
185 // Al, Si, Ca for 45 energies
186 
187 const G4double G4NucleonNuclearCrossSection::e3[45] =
188 {
189  0.014, 0.015, 0.016, 0.017, .02, 0.022, 0.025, 0.027, 0.03, 0.035, .04, 0.045, 0.05, .06, 0.07,
190  .08, 0.09, .1, .12, .14, .15, .16, .18, .20, .25, .30, .35, 0.4, 0.5, 0.6,
191  0.7, 0.8, 0.9, 1, 1.5, 2, 3, 5, 7, 10, 20, 50, 100, 500, 1000
192 };
193 
194 const G4double G4NucleonNuclearCrossSection::al_m_t[45] =
195 {
196  1735, 1750, 1760, 1795, 1830, 1855, 1885, 1895, 1900, 1870, 1835, 1785, 1710, 1522, 1350,
197  1212, 1080, 972, 816, 720, 678, 642, 600, 567, 558, 560, 578, 592, 616, 644,
198  672, 688, 708, 720, 736, 754, 736, 706, 680, 672, 646, 632, 632, 632, 632
199 };
200 const G4double G4NucleonNuclearCrossSection::al_m_in[45] =
201 {
202  1000, 990, 975, 950, 905, 875, 825, 800, 762, 690, 652, 610, 570, 495, 480,
203  456, 444, 432, 420, 420, 420, 420, 410, 410, 400, 402, 404, 408, 424, 438,
204  448, 450, 454, 456, 472, 480, 466, 456, 452, 448, 444, 440, 440, 440, 440
205 };
206 const G4double G4NucleonNuclearCrossSection::al_p_in[45] =
207 {
208  650, 682, 690, 715, 750, 762, 750, 740, 720, 655, 617, 575, 540, 470, 455,
209  // 532, 420, 408, 400, 403, 403, 408, 406, 404, 400, 402, 404, 408, 424, 438,
210  432, 420, 408, 400, 403, 403, 408, 406, 404, 400, 402, 404, 408, 424, 438,
211  448, 450, 454, 456, 472, 480, 466, 456, 452, 448, 444, 440, 440, 440, 440
212 };
213 
214 const G4double G4NucleonNuclearCrossSection::si_m_t[45] =
215 {
216  1810, 1833, 1850, 1872, 1920, 1950, 1995, 2020, 2035, 2000, 1930, 1850, 1760, 1570, 1400,
217  1255, 1110, 1008, 846, 742, 696, 671, 623, 588, 584, 584, 602, 618, 645, 679,
218  708, 727, 746, 757, 769, 782, 771, 734, 710, 698, 672, 654, 650, 650, 650
219 };
220 const G4double G4NucleonNuclearCrossSection::si_m_in[45] =
221 {
222  1060, 1035, 1015, 990, 935, 900, 860, 830, 790, 725, 665, 630, 600, 520, 504,
223  486, 470, 456, 444, 432, 432, 432, 418, 418, 415, 412, 416, 422, 440, 460,
224  472, 476, 479, 480, 492, 496, 488, 472, 472, 464, 460, 452, 448, 448, 448
225 };
226 const G4double G4NucleonNuclearCrossSection::si_p_in[45] =
227 {
228  670, 700, 725, 750, 780, 780, 770, 757, 735, 690, 635, 585, 570, 490, 475,
229  460, 446, 431, 423, 425, 425, 425, 425, 422, 422, 412, 416, 422, 440, 460,
230  472, 476, 479, 480, 492, 496, 488, 472, 472, 464, 460, 452, 448, 448, 448
231 };
232 
233 const G4double G4NucleonNuclearCrossSection::ca_m_t[45] =
234 {
235  2180, 2130, 2095, 2075, 2115, 2150, 2220, 2250, 2300, 2365, 2360, 2280, 2180, 2000,
236  1805, 1650, 1500, 1340, 1140, 990, 940, 890, 825, 790, 770, 773, 787, 800, 830, 870,
237  905, 930, 950, 965, 990, 1002, 990, 965, 945, 925, 892, 860, 860, 860, 860
238 };
239 const G4double G4NucleonNuclearCrossSection::ca_m_in[45] =
240 {
241  1240, 1225, 1200, 1180, 1125, 1090, 1045, 1020, 980, 925, 880, 825, 770, 680, 640,
242  620, 615, 600, 580, 565, 560, 560, 560, 550, 535, 530, 540, 550, 570, 595, 610, 615,
243  620, 622, 629, 630, 620, 612, 607, 592, 587, 580, 580, 580, 580
244 };
245 const G4double G4NucleonNuclearCrossSection::ca_p_in[45] =
246 {
247  770, 800, 823, 850, 900, 925, 935, 920, 895, 835, 800, 750, 715, 640, 605, 590, 588,
248  573, 555, 543, 540, 540, 540, 535, 530, 530, 540, 550, 570, 595, 610, 615,
249  620, 622, 629, 630, 620, 612, 607, 592, 587, 580, 580, 580, 580
250 };
251 
252 // Fe, Cu, Mo for 47 energies
253 
254 const G4double G4NucleonNuclearCrossSection::e4[47] =
255 {
256  0.014, 0.015, 0.017, .02, 0.022, 0.025, 0.027, 0.03, 0.033, 0.035, 0.037, .04, 0.045,
257  0.05, 0.055, .06, 0.07, .08, 0.09, .1, .12, .14, .15, .16, .18, .20, .25, .30, .35,
258  .4 , 0.5, 0.6, 0.7, 0.8, 0.9, 1, 1.5, 2, 3, 5, 7, 10,
259  20, 50, 100, 500, 1000
260 };
261 
262 const G4double G4NucleonNuclearCrossSection::fe_m_t[47] =
263 {
264  2580, 2490, 2370, 2282, 2275, 2285, 2320, 2370, 2432, 2445, 2460, 2485, 2530, 2540,
265  2517, 2480, 2290, 2110, 1940, 1790, 1510, 1290, 1220, 1150, 1070, 1030, 1013, 1020,
266  1030, 1043, 1075, 1110, 1133, 1163, 1185, 1225, 1252, 1260, 1260, 1233, 1207, 1185,
267  1140, 1110, 1110, 1110, 1110
268 };
269 const G4double G4NucleonNuclearCrossSection::fe_m_in[47] =
270 {
271  1440, 1433, 1390, 1325, 1280, 1260, 1215, 1180, 1140, 1110, 1080, 1040, 990, 955, 920,
272  885, 835, 800, 780, 765, 750, 725, 720, 720, 710, 700, 700, 700, 712, 705, 735, 750,
273  765, 775, 780, 795, 810, 813, 810, 784, 757, 743, 735, 720, 720, 720, 720
274 };
275 const G4double G4NucleonNuclearCrossSection::fe_p_in[47] =
276 {
277  900, 960, 1070, 1090, 1115, 1120, 1115, 1080, 1045, 1025, 1000, 960, 900, 885, 865, 790,
278  765, 740, 720, 700, 697, 697, 697, 697, 695, 690, 688, 690, 712, 705, 735, 750,
279  765, 775, 780, 795, 810, 813, 810, 784, 757, 743, 735, 720, 720, 720, 720
280 };
281 
282 const G4double G4NucleonNuclearCrossSection::cu_m_t[47] =
283 {
284  2920, 2800, 2615, 2480, 2455, 2430, 2440, 2460, 2500, 2530, 2560, 2615, 2690, 2720,
285  2700, 2645, 2500, 2320, 2140, 1970, 1670, 1460, 1380, 1285, 1200, 1160, 1140, 1147,
286  1163, 1170, 1200, 1237, 1265, 1285, 1305, 1328, 1375, 1390, 1395, 1370, 1335, 1315,
287  1270, 1230, 1230, 1230, 1230
288 };
289 const G4double G4NucleonNuclearCrossSection::cu_m_in[47] =
290 {
291  1540, 1535, 1500, 1445, 1407, 1380, 1330, 1300, 1285, 1270, 1240, 1190, 1090, 1010,
292  940, 920, 860, 835, 820, 810, 800, 780, 775, 770, 760, 760, 758, 765, 765, 770, 795,
293  810, 825, 830, 840, 848, 870, 870, 868, 840, 825, 810, 803, 795, 795, 795, 795
294 };
295 const G4double G4NucleonNuclearCrossSection::cu_p_in[47] =
296 {
297  935, 1000, 1060, 1190, 1220, 1250, 1240, 1210, 1150, 1130, 1115, 1050, 985, 950, 890,
298  870, 820, 800, 785, 780, 770, 750, 745, 740, 735, 735, 745, 760, 762, 770, 795,
299  810, 825, 830, 840, 848, 870, 870, 868, 840, 825, 810, 803, 795, 795, 795, 795
300 };
301 
302 const G4double G4NucleonNuclearCrossSection::mo_m_t[47] =
303 {
304  4150, 4040, 3800, 3490, 3300, 3060, 2960, 2845, 2785, 2820, 2850, 2980, 3170, 3230,
305  3270, 3280, 3225, 3075, 2895, 2710, 2355, 2060, 1925, 1800, 1630, 1560, 1540, 1550,
306  1570, 1590, 1650, 1685, 1715, 1740, 1760, 1780, 1850, 1880, 1858, 1815, 1790, 1782,
307  1720, 1690, 1690, 1690, 1690
308 };
309 const G4double G4NucleonNuclearCrossSection::mo_m_in[47] =
310 {
311  1790, 1775, 1740, 1680, 1640, 1580, 1550, 1510, 1460, 1440, 1418, 1380, 1330, 1280,
312  1240, 1200, 1155, 1140, 1110, 1110, 1080, 1065, 1050, 1050, 1025, 1020, 1015, 1020,
313  1022, 1026, 1060, 1085, 1100, 1110, 1120, 1127, 1150, 1160, 1140, 1100, 1085, 1080,
314  1070, 1070, 1070, 1070, 1070
315 };
316 const G4double G4NucleonNuclearCrossSection::mo_p_in[47] =
317 {
318  1025, 1080, 1190, 1380, 1440, 1495, 1475, 1420, 1350, 1310, 1300, 1290, 1250, 1200,
319  1170, 1130, 1095, 1060, 1040, 1022, 1020, 1016, 1016, 1016, 1016, 1012, 1005, 1005,
320  1005, 1010, 1060, 1085, 1100, 1110, 1120, 1127, 1150, 1160, 1140, 1100, 1085, 1080,
321  1070, 1070, 1070, 1070, 1070
322 };
323 
324 // Cd, Sn, W for 48 energies
325 
326 const G4double G4NucleonNuclearCrossSection::e5[48] =
327 {
328  0.014, 0.015, 0.017, 0.018, .02, 0.022, 0.025, 0.027, 0.03, 0.033, 0.035, .04,
329  0.045, 0.05, 0.055, .06, .065, 0.07, .08, 0.09, .1, .12, .14, .15, .16, .18,
330  .20, .25, .30, .35, .4 , 0.5, 0.6, 0.7, 0.8, 0.9, 1, 1.5, 2, 3, 5, 7, 10, 20,
331  50, 100, 500, 1000
332 };
333 
334 const G4double G4NucleonNuclearCrossSection::cd_m_t[48] =
335 {
336  4420, 4280, 4170, 4070, 3860, 3680, 3420, 3280, 3125, 3060, 3080, 3190, 3350, 3445,
337  3510, 3540, 3560, 3550, 3460, 3300, 3030, 2640, 2340, 2190, 2070, 1950, 1770, 1732,
338  1740, 1760, 1780, 1832, 1885, 1925, 1945, 1960, 1980, 2070, 2080, 2065, 2040, 2022,
339  1980, 1940, 1870, 1870, 1870, 1870
340 };
341 const G4double G4NucleonNuclearCrossSection::cd_m_in[48]=
342 {
343  1920, 1910, 1880, 1860, 1840, 1800, 1760, 1720, 1675, 1630, 1600, 1520, 1465, 1420,
344  1390, 1340, 1310, 1280, 1275, 1235, 1225, 1200, 1170, 1170, 1170, 1165, 1145, 1140,
345  1140, 1135, 1160, 1180, 1220, 1240, 1250, 1260, 1265, 1270, 1275, 1250, 1222, 1222,
346  1220, 1215, 1190, 1190, 1190, 1190
347 };
348 const G4double G4NucleonNuclearCrossSection::cd_p_in[48] =
349 {
350  1020, 1100, 1225, 1290, 1440, 1520, 1575, 1560, 1518, 1460, 1420, 1400, 1365, 1340,
351  1300, 1280, 1260, 1200, 1190, 1160, 1125, 1125, 1125, 1125, 1125, 1125, 1120, 1120,
352  1120, 1118, 1146, 1180, 1220, 1240, 1250, 1260, 1265, 1270, 1275, 1250, 1222, 1222,
353  1220, 1215, 1190, 1190, 1190, 1190
354 };
355 
356 const G4double G4NucleonNuclearCrossSection::sn_m_t[48] =
357 {
358  4420, 4400, 4260, 4150, 3980, 3770, 3530, 3370, 3245, 3180, 3170, 3260, 3400, 3500,
359  3560, 3610, 3650, 3680, 3580, 3390, 3190, 2760, 2430, 2295, 2175, 1990, 1880, 1810,
360  1820, 1840, 1865, 1940, 1985, 2020, 2040, 2060, 2080, 2160, 2185, 2180, 2110, 2105,
361  2080, 2050, 1980, 1980, 1980, 1980
362 };
363 const G4double G4NucleonNuclearCrossSection::sn_m_in[48] =
364 {
365  1945, 1940, 1905, 1890, 1860, 1830, 1780, 1755, 1717, 1680, 1645, 1570, 1500, 1455,
366  1410, 1370, 1340, 1320, 1290, 1285, 1260, 1240, 1235, 1212, 1200, 1200, 1200, 1190,
367  1190, 1200, 1210, 1240, 1270, 1285, 1300, 1300, 1310, 1320, 1320, 1290, 1240, 1240,
368  1240, 1240, 1240, 1240, 1240, 1240
369 };
370 const G4double G4NucleonNuclearCrossSection::sn_p_in[48] =
371 {
372  1020, 1080, 1270, 1335, 1465, 1505, 1610, 1610, 1550, 1535, 1500, 1440, 1407, 1370,
373  1340, 1300, 1285, 1260, 1230, 1215, 1200, 1180, 1170, 1170, 1165, 1165, 1170, 1165,
374  1165, 1183, 1195, 1240, 1270, 1285, 1300, 1300, 1310, 1320, 1320, 1290, 1240, 1240,
375  1240, 1240, 1240, 1240, 1240, 1240
376 };
377 
378 const G4double G4NucleonNuclearCrossSection::w_m_t[48] =
379 {
380  5320, 5430, 5480, 5450, 5330, 5190, 4960, 4790, 4550, 4340, 4200, 4070, 4000, 4030,
381  4125, 4220, 4270, 4390, 4440, 4360, 4200, 3800, 3380, 3200, 3040, 2790, 2660, 2575,
382  2575, 2600, 2640, 2690, 2755, 2790, 2812, 2837, 2850, 2950, 3000, 2970, 2940, 2910,
383  2880, 2820, 2730, 2730, 2730, 2730
384 };
385 const G4double G4NucleonNuclearCrossSection::w_m_in[48] =
386 {
387  2440, 2400, 2370, 2350, 2310, 2270, 2220, 2195, 2150, 2100, 2070, 2010, 1945, 1900,
388  1850, 1820, 1780, 1760, 1730, 1720, 1680, 1680, 1660, 1660, 1650, 1650, 1640, 1640,
389  1612, 1615, 1625, 1640, 1700, 1720, 1730, 1740, 1750, 1780, 1780, 1750, 1740, 1735,
390  1710, 1695, 1680, 1680, 1680, 1680
391 };
392 const G4double G4NucleonNuclearCrossSection::w_p_in[48] =
393 {
394  950, 1020, 1240, 1400, 1560, 1670, 1760, 1830, 1850, 1855, 1870, 1840, 1800, 1770,
395  1740, 1715, 1680, 1670, 1650, 1620, 1610, 1600, 1600, 1600, 1600, 1600, 1600, 1595,
396  1585, 1595, 1615, 1640, 1700, 1720, 1730, 1740, 1750, 1780, 1780, 1750, 1740, 1735,
397  1710, 1695, 1680, 1680, 1680, 1680
398 };
399 
400 // Pb, U for 46 energies
401 
402 const G4double G4NucleonNuclearCrossSection::e6[46] =
403 {
404  0.014, 0.015, 0.017, 0.019, 0.02, 0.022, 0.025, 0.027, 0.03, 0.035,
405  0.04, 0.045, 0.05, 0.055, 0.06, 0.07, 0.08, 0.09, 0.1, 0.12,
406  0.14, 0.15, 0.16, 0.18, 0.20, 0.25, 0.30, 0.35, 0.4 , 0.5,
407  0.6, 0.7, 0.8, 0.9, 1.0, 1.5, 2.0, 3.0, 5.0, 7.0,
408  10.0, 20.0, 50.0, 100.0, 500.0, 1000.0
409 };
410 
411 const G4double G4NucleonNuclearCrossSection::pb_m_t[46] =
412 {
413  5300, 5440, 5720, 5880, 5765, 5745, 5480, 5280, 4970, 4550, 4390, 4300, 4265, 4325,
414  4450, 4540, 4740, 4710, 4600, 4100, 3660, 3480, 3300, 3000, 2890, 2865, 2855, 2850,
415  2865, 2920, 2955, 3000, 3030, 3060, 3105, 3240, 3290, 3270, 3240, 3180, 3090, 3060,
416  2970, 2970, 2970, 2970
417 
418 };
419 const G4double G4NucleonNuclearCrossSection::pb_m_in[46] =
420 {
421  2580, 2550, 2505, 2462, 2460, 2435, 2380, 2355, 2280, 2180, 2170, 2130, 2080, 2035,
422  1980, 1940, 1900, 1870, 1840, 1800, 1800, 1800, 1780, 1760, 1760, 1740, 1730, 1725,
423  1740, 1785, 1815, 1835, 1860, 1890, 1895, 1920, 1920, 1890, 1850, 1835, 1830, 1830,
424  1830, 1830, 1830, 1830
425 };
426 const G4double G4NucleonNuclearCrossSection::pb_p_in[46] =
427 {
428  900, 1060, 1200, 1420, 1515, 1620, 1750, 1800, 1915, 2030, 1960, 1940, 1910, 1860,
429  1840, 1780, 1770, 1760, 1740, 1720, 1725, 1740, 1740, 1730, 1720, 1700, 1710, 1720,
430  1730, 1740, 1815, 1835, 1860, 1890, 1895, 1920, 1920, 1890, 1850, 1835, 1830, 1830,
431  1830, 1830, 1830, 1830
432 };
433 
434 const G4double G4NucleonNuclearCrossSection::u_m_t[46] =
435 {
436  5800, 5940, 6160, 6345, 6360, 6350, 6170, 6020, 5760, 5350, 4990, 4800, 4710, 4690,
437  4760, 5040, 5190, 5200, 5080, 4600, 4120, 3920, 3720, 3420, 3240, 3150, 3160, 3180,
438  3210, 3240, 3280, 3350, 3390, 3435, 3480, 3560, 3585, 3580, 3540, 3500, 3470, 3410,
439  3335, 3335, 3335, 3335
440 };
441 const G4double G4NucleonNuclearCrossSection::u_m_in[46] =
442 {
443  2820, 2770, 2700, 2660, 2645, 2620, 2580, 2550, 2515, 2450, 2390, 2320, 2260, 2225,
444  2200, 2140, 2080, 2060, 2040, 2000, 1980, 1965, 1960, 1930, 1920, 1890, 1905, 1920,
445  1945, 1970, 1985, 2010, 2040, 2070, 2080, 2090, 2095, 2080, 2063, 2060, 2050, 2040,
446  2005, 2005, 2005, 2005
447 };
448 const G4double G4NucleonNuclearCrossSection::u_p_in[46] =
449 {
450  800, 900, 1100, 1300, 1410, 1510, 1680, 1800, 2000, 2200, 2080, 2060, 2035, 2100,
451  2030, 2030, 2000, 1960, 1960, 1960, 1940, 1925, 1920, 1905, 1890, 1860, 1880, 1910,
452  1930, 1945, 1985, 2010, 2040, 2070, 2080, 2090, 2095, 2080, 2063, 2060, 2050, 2040,
453  2005, 2005, 2005, 2005
454 };
455 
456 using namespace std;
457 
459 
461  : G4VCrossSectionDataSet(Default_Name()),
462  fTotalXsc(0.0), fElasticXsc(0.0)
463 {
464  theNeutron = G4Neutron::Neutron();
465  theProton = G4Proton::Proton();
466 
467  // He, Be, C
468 
469  thePimData.push_back(new G4PiData(he_m_t, he_m_in, e1, 44));
470  thePipData.push_back(new G4PiData(he_m_t, he_p_in, e1, 44));
471 
472  thePimData.push_back(new G4PiData(be_m_t, be_m_in, e1, 44));
473  thePipData.push_back(new G4PiData(be_m_t, be_p_in, e1, 44));
474 
475  thePimData.push_back(new G4PiData(c_m_t, c_m_in, e1, 44));
476  thePipData.push_back(new G4PiData(c_m_t, c_p_in, e1, 44));
477 
478  // N, O, Na
479 
480  thePimData.push_back(new G4PiData(n_m_t, n_m_in, e2, 44));
481  thePipData.push_back(new G4PiData(n_m_t, n_p_in, e2, 44));
482 
483  thePimData.push_back(new G4PiData(o_m_t, o_m_in, e2, 44));
484  thePipData.push_back(new G4PiData(o_m_t, o_p_in, e2, 44));
485 
486  thePimData.push_back(new G4PiData(na_m_t, na_m_in, e2, 44));
487  thePipData.push_back(new G4PiData(na_m_t, na_p_in, e2, 44));
488 
489  // Al, Si, Ca
490 
491  thePimData.push_back(new G4PiData(al_m_t, al_m_in, e3, 45));
492  thePipData.push_back(new G4PiData(al_m_t, al_p_in, e3, 45));
493 
494  thePimData.push_back(new G4PiData(si_m_t, si_m_in, e3, 45));
495  thePipData.push_back(new G4PiData(si_m_t, si_p_in, e3, 45));
496 
497  thePimData.push_back(new G4PiData(ca_m_t, ca_m_in, e3, 45));
498  thePipData.push_back(new G4PiData(ca_m_t, ca_p_in, e3, 45));
499 
500  // Fe, Cu, Mo
501 
502  thePimData.push_back(new G4PiData(fe_m_t, fe_m_in, e4, 47));
503  thePipData.push_back(new G4PiData(fe_m_t, fe_p_in, e4, 47));
504 
505  thePimData.push_back(new G4PiData(cu_m_t, cu_m_in, e4, 47));
506  thePipData.push_back(new G4PiData(cu_m_t, cu_p_in, e4, 47));
507 
508  thePimData.push_back(new G4PiData(mo_m_t, mo_m_in, e4, 47));
509  thePipData.push_back(new G4PiData(mo_m_t, mo_p_in, e4, 47));
510 
511  // Cd, Sn, W
512 
513  thePimData.push_back(new G4PiData(cd_m_t, cd_m_in, e5, 48));
514  thePipData.push_back(new G4PiData(cd_m_t, cd_p_in, e5, 48));
515 
516  thePimData.push_back(new G4PiData(sn_m_t, sn_m_in, e5, 48));
517  thePipData.push_back(new G4PiData(sn_m_t, sn_p_in, e5, 48));
518 
519  thePimData.push_back(new G4PiData(w_m_t, w_m_in, e5, 48));
520  thePipData.push_back(new G4PiData(w_m_t, w_p_in, e5, 48));
521 
522  // Pb, U
523 
524  thePimData.push_back(new G4PiData(pb_m_t, pb_m_in, e6, 46));
525  thePipData.push_back(new G4PiData(pb_m_t, pb_p_in, e6, 46));
526 
527  thePimData.push_back(new G4PiData(u_m_t, u_m_in, e6, 46));
528  thePipData.push_back(new G4PiData(u_m_t, u_p_in, e6, 46));
529 
530  theZ.push_back(2); // He
531  theZ.push_back(4); // Be
532  theZ.push_back(6); // C
533  theZ.push_back(7); // N
534  theZ.push_back(8); // O
535  theZ.push_back(11); // Na
536  theZ.push_back(13); // Al
537  theZ.push_back(14); // Si
538  theZ.push_back(20); // Ca
539  theZ.push_back(26); // Fe
540  theZ.push_back(29); // Cu
541  theZ.push_back(42); // Mo
542  theZ.push_back(48); // Cd
543  theZ.push_back(50); // Sn
544  theZ.push_back(74); // W
545  theZ.push_back(82); // Pb
546  theZ.push_back(92); // U
547 
548 }
549 
551 //
552 
554 {
555  std::for_each(thePimData.begin(), thePimData.end(), G4PiData::Delete());
556  std::for_each(thePipData.begin(), thePipData.end(), G4PiData::Delete());
557 }
558 
560 //
561 
562 G4bool
564  G4int Z,
565  const G4Material*)
566 {
567  G4bool result = false;
568  if(aParticle->GetDefinition() == theNeutron ) result = true;
569  if(aParticle->GetDefinition() == theProton) result = true;
570  if(Z < 2) result = false;
571  if(aParticle->GetKineticEnergy() > 999.9*GeV) result = false;
572  return result;
573 }
574 
576 //
577 
578 G4double
580  G4int Z, const G4Material*)
581 {
582  G4double kineticEnergy = aParticle->GetKineticEnergy();
583 
584  G4double result = 0;
585  // G4cout<<"Z = "<<Z<<G4endl;
586 
587  size_t it = 0;
588  size_t itmax = theZ.size() - 1;
589  for(; it <= itmax; ++it) { if(Z <= theZ[it]) { break; } }
590  if( it > itmax ) { it = itmax; }
591  G4int Z1, Z2;
592  G4double x1, x2, xt1, xt2;
593 
594  std::vector<G4PiData *> * theData = &thePimData;
595  if(aParticle->GetDefinition() == theProton) { theData = &thePipData; }
596 
597  if( theZ[it] == Z )
598  {
599  result = (*theData)[it]->ReactionXSection(kineticEnergy);
600  fTotalXsc = (*theData)[it]->TotalXSection(kineticEnergy);
601  }
602  else
603  {
604  if(0 == it) { it = 1; }
605  x1 = (*theData)[it-1]->ReactionXSection(kineticEnergy);
606  xt1 = (*theData)[it-1]->TotalXSection(kineticEnergy);
607  Z1 = theZ[it-1];
608  x2 = (*theData)[it]->ReactionXSection(kineticEnergy);
609  xt2 = (*theData)[it]->TotalXSection(kineticEnergy);
610  Z2 = theZ[it];
611 
612  result = Interpolate(Z1, Z2, Z, x1, x2);
613  fTotalXsc = Interpolate(Z1, Z2, Z, xt1, xt2);
614  }
615 
616  fElasticXsc = fTotalXsc - result;
617  if( fElasticXsc < 0.) { fElasticXsc = 0.; }
618 
619  return result;
620 }
621 
623 //
624 
625 G4double G4NucleonNuclearCrossSection::
626 Interpolate(G4int Z1, G4int Z2, G4int Z, G4double x1, G4double x2)
627 {
628 // Nucleon numbers obtained from G4NistManager G4 8.0
629 
630  static const G4double alpha = 2./3.;
631 
632  static const G4double A[92] =
633  {
634  1.0001, 4.0000, 6.9241, 9.0000, 10.801, 12.011, 14.004, 16.004, 19.000, 20.188,
635  23.000, 24.320, 27.000, 28.109, 31.000, 32.094, 35.484, 39.985, 39.135, 40.116,
636  45.000, 47.918, 50.998, 52.055, 55.000, 55.910, 59.000, 58.760, 63.617, 65.468,
637  69.798, 72.691, 75.000, 79.042, 79.986, 83.887, 85.557, 87.710, 89.000, 91.318,
638  93.000, 96.025, 98.000, 101.16, 103.00, 106.51, 107.96, 112.51, 114.91, 118.81,
639  121.86, 127.70, 127.00, 131.39, 133.00, 137.42, 139.00, 140.21, 141.00, 144.32,
640  145.00, 150.45, 152.04, 157.33, 159.00, 162.57, 165.00, 167.32, 169.00, 173.10,
641  175.03, 178.54, 181.00, 183.89, 186.25, 190.27, 192.25, 195.11, 197.00, 200.63,
642  204.41, 207.24, 209.00, 209.00, 210.00, 222.00, 223.00, 226.00, 227.00, 232.00,
643  231.00, 237.98
644  };
645  static G4ThreadLocal G4bool NeedInit = true;
646 
647  static G4ThreadLocal G4double A75[92];
648 
649  if ( NeedInit )
650  {
651  for (G4int i=0; i<92; ++i)
652  {
653  A75[i] = G4Pow::GetInstance()->powA(A[i], alpha); // interpolate by square ~ A^(2/3)
654  }
655  NeedInit=false;
656  }
657 
658  // for tabulated data, cross section scales with A^(2/3)
659  G4double r1 = x1 / A75[Z1-1] * A75[Z-1];
660  G4double r2 = x2 / A75[Z2-1] * A75[Z-1];
661  G4double result = 0.5*(r1+r2);
662 
663  // More precise average
664  if(Z1 != Z2) {
665  G4double alp1 = (A[Z-1] - A[Z1-1]);
666  G4double alp2 = (A[Z2-1] - A[Z-1]);
667  result = (r1*alp2 + r2*alp1)/(alp1 + alp2);
668  }
669  // G4cout << "x1/2, z1/2 z" <<x1<<" "<<x2<<" "<<Z1<<" "<<Z2<<" "<<Z<<G4endl;
670  // G4cout << "res1/2 " << r1 <<" " << r2 <<" " << result<< G4endl;
671  return result;
672 }
673 
674 void
676 {
677  outFile << "G4NucleonNuclearCrossSection is a variant of the Barashenkov\n"
678  << "cross section parameterization to be used of protons and\n"
679  << "nucleons on targets heavier than hydrogen. It is intended for\n"
680  << "use as a cross section component and is currently used by\n"
681  << "G4BGGNucleonInelasticXS. It is valid for incident energies up\n"
682  << "to 1 TeV.\n";
683 }
684 
G4double G4ParticleHPJENDLHEData::G4double result
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
G4double GetKineticEnergy() const
virtual G4bool IsElementApplicable(const G4DynamicParticle *aParticle, G4int Z, const G4Material *mat=0)
virtual G4double GetElementCrossSection(const G4DynamicParticle *aParticle, G4int Z, const G4Material *mat=0)
G4ParticleDefinition * GetDefinition() const
virtual void CrossSectionDescription(std::ostream &) const
#define G4ThreadLocal
Definition: tls.hh:89
int G4int
Definition: G4Types.hh:78
double A(double temperature)
bool G4bool
Definition: G4Types.hh:79
static G4Proton * Proton()
Definition: G4Proton.cc:93
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
#define G4_DECLARE_XS_FACTORY(cross_section)
static constexpr double GeV
Definition: G4SIunits.hh:217
double G4double
Definition: G4Types.hh:76
static const G4double alpha