00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027 #include "G4PiNuclearCrossSection.hh"
00028 #include "G4SystemOfUnits.hh"
00029 #include "G4DynamicParticle.hh"
00030 #include "G4HadronicException.hh"
00031 #include "G4HadTmpUtil.hh"
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045 const G4double G4PiNuclearCrossSection::e1[38] = {
00046 .02, .04, .06, .08, .1, .12, .13, .14, .15, .16, .17, .18, .19, .20,
00047 .22, .24, .26, .28, .30, .35, .40, .45, 0.5, 0.55, 0.6, 0.7, 0.8, 0.9,
00048 1, 2, 3, 5, 10, 20, 50, 100, 500, 100000};
00049
00050 const G4double G4PiNuclearCrossSection::he_t[38] = {
00051 40, 70, 108, 152, 208, 276, 300, 320, 329, 333, 332, 328, 322, 310, 288,
00052 260, 240, 216, 196, 144, 125, 112,108.5, 109, 110.5, 117, 123,128.5, 135,
00053 110, 96, 87, 85, 83.5, 83.5, 83.5, 83.5, 83.5};
00054
00055 const G4double G4PiNuclearCrossSection::he_in[38] = {
00056 18, 38, 62, 98, 136, 176, 190, 200, 209, 212, 212, 208, 204, 196,
00057 176, 164, 150, 134, 124,97.5, 90, 85, 82.5, 83.5, 86.5, 93, 97.5, 100,
00058 102, 83, 77, 75, 74, 72.5, 72.5, 72.5, 72.5, 72.5};
00059
00060 const G4double G4PiNuclearCrossSection::be_m_t[38] = {
00061 150, 210, 294, 396, 520, 600, 623, 635, 642, 640, 630, 615, 600, 576, 540,
00062 504, 470, 435, 400, 340, 294, 258, 236, 230, 233, 244, 257, 270, 276, 250,
00063 230, 215, 205, 194, 188, 186, 186, 186};
00064
00065 const G4double G4PiNuclearCrossSection::be_m_in[38] = {
00066 90, 126, 177, 240, 320, 380, 400, 410, 414, 410, 400, 387, 371, 360, 333,
00067 312, 285, 260, 237, 216, 198, 187, 182, 180, 182, 187, 193, 203, 207, 179,
00068 172, 165, 159, 155, 144, 144, 144, 144};
00069
00070 const G4double G4PiNuclearCrossSection::be_p_t[24] = {
00071 96, 150, 222, 320, 430, 514, 545, 565, 574, 574, 564, 552, 535, 522, 490,
00072 462, 432, 398, 367, 314, 276, 248, 232, 230};
00073
00074 const G4double G4PiNuclearCrossSection::be_p_in[24] = {
00075 60, 95, 142, 194, 262, 319, 345, 361, 364, 364, 354, 350, 330, 319, 298,
00076 280, 258, 237, 216, 200, 189, 183, 182, 180};
00077
00078 const G4double G4PiNuclearCrossSection::e2[39] = {
00079 .02, .04, .06, .08, .10, .11, .12, .13, .14, .15, .16, .17, .18, .20, .22,
00080 .24, .26, .28, .30, .35, .40, .45, .50, .55, .575, .60, .70, .80, .90, 1,
00081 2, 3, 5, 10, 20, 50, 100, 500, 100000};
00082
00083 const G4double G4PiNuclearCrossSection::c_m_t[39] = {
00084 204, 260, 366, 517, 630, 673, 694, 704, 710, 711, 706, 694, 676, 648, 616,
00085 584, 548, 518, 489, 426, 376, 342, 323, 310, 312, 313, 319, 333, 342, 348,
00086 310, 290, 268, 250, 245, 237, 234, 234, 234};
00087
00088 const G4double G4PiNuclearCrossSection::c_m_in[39] = {
00089 128, 160, 224, 315, 388, 416, 430, 438, 444, 445, 440, 432, 416, 400, 380,
00090 354, 320, 304, 288, 264, 246, 240, 233, 232, 233, 234, 238, 246, 252, 256,
00091 220, 210, 198, 187, 183, 176, 174, 174, 174};
00092
00093 const G4double G4PiNuclearCrossSection::c_p_t[24] = {
00094 140, 192, 294, 428, 594, 642, 662, 687, 685, 688, 684, 672, 656, 630, 598,
00095 567, 533, 504, 474, 416, 369, 336, 319, 310};
00096
00097 const G4double G4PiNuclearCrossSection::c_p_in[24] = {
00098 94, 132, 184, 260, 370, 398, 408, 420, 426, 428, 424, 416, 400, 386, 366,
00099 340, 308, 294, 280, 257, 241, 236, 231, 232};
00100
00101 const G4double G4PiNuclearCrossSection::n_m_t[39] = {
00102 246, 308, 424, 590, 729, 776, 800, 821, 822, 817, 800, 778, 768, 728, 690,
00103 654, 615, 584, 556, 480, 430, 393, 373, 367, 368, 370, 375, 388, 390, 397,
00104 364, 337, 310, 291, 275, 268, 268, 268, 268};
00105
00106 const G4double G4PiNuclearCrossSection::n_m_in[39] = {
00107 155, 188, 256, 360, 456, 492, 512, 526, 526, 520, 504, 491, 475, 450, 425,
00108 396, 376, 360, 340, 300, 282, 270, 265, 265, 266, 268, 273, 280, 288, 288,
00109 256, 237, 226, 218, 208, 202, 202, 202, 202};
00110
00111 const G4double G4PiNuclearCrossSection::n_p_t[27] = {
00112 150, 212, 328, 500, 680, 735, 762, 781, 782, 779, 770, 748, 740, 706, 672,
00113 633, 600, 569, 541, 467, 419, 385, 368, 364, 366, 368, 375};
00114
00115 const G4double G4PiNuclearCrossSection::n_p_in[27] = {
00116 90, 140, 208, 300, 426, 467, 490, 504, 504, 500, 484, 474, 460, 437, 413,
00117 381, 365, 350, 330, 292, 276, 267, 263, 264, 265, 267, 273};
00118
00119 const G4double G4PiNuclearCrossSection::e3[31] = {
00120 .02, .04, .06, .08, .10, .12, .14, .16, .18, .20, .22, .25, .30, .35, .40,
00121 .45, .50, .60, .70, .80, .90, 1, 2, 3, 5, 10, 20, 50, 100, 500,
00122 100000};
00123
00124 const G4double G4PiNuclearCrossSection::o_m_t[31] = {
00125 280, 360, 500, 685, 812, 861, 870, 865, 835, 800, 755, 700, 600, 537, 493,
00126 468, 441, 436, 443, 449, 460, 463, 432, 385, 350, 325, 312, 307, 303, 303,
00127 303};
00128
00129 const G4double G4PiNuclearCrossSection::o_m_in[31] = {
00130 190, 207, 300, 420, 500, 540, 550, 542, 520, 490, 460, 423, 360, 339, 321,
00131 314, 312, 314, 319, 324, 328, 330, 300, 275, 250, 240, 229, 225, 222, 222,
00132 222};
00133
00134 const G4double G4PiNuclearCrossSection::o_p_t[20] = {
00135 170, 240, 390, 570, 740, 818, 830, 822, 800, 765, 725, 675, 585, 525, 483,
00136 458, 444, 447, 453, 449};
00137
00138 const G4double G4PiNuclearCrossSection::o_p_in[20] = {
00139 100, 145, 240, 340, 470, 518, 530, 522, 505, 477, 448, 412, 350, 330, 316,
00140 310, 308, 311, 317, 324};
00141
00142 const G4double G4PiNuclearCrossSection::na_m_t[31] = {
00143 450, 545, 705, 910, 1020, 1075, 1087, 1080, 1042, 987, 943, 885, 790, 700,
00144 650, 610, 585, 575, 585, 595, 600, 610, 556, 524, 494, 458, 445, 429,
00145 427, 427, 427};
00146
00147 const G4double G4PiNuclearCrossSection::na_m_in[31] = {
00148 275, 315, 413, 545, 620, 660, 670, 662, 630, 593, 570, 520, 465, 420, 410,
00149 395, 390, 400, 410, 418, 420, 422, 372, 348, 330, 320, 310, 294, 292, 292,
00150 292};
00151
00152 const G4double G4PiNuclearCrossSection::na_p_t[22] = {
00153 210, 320, 530, 795, 960, 1035, 1050, 1040, 1007, 957, 918, 865, 773, 685,
00154 636, 598, 575, 565, 578, 590, 598, 610};
00155
00156 const G4double G4PiNuclearCrossSection::na_p_in[22] = {
00157 115, 210, 340, 495, 585, 630, 645, 637, 605, 572, 550, 505, 455, 410, 401,
00158 388, 383, 393, 405, 414, 418, 422};
00159
00160 const G4double G4PiNuclearCrossSection::e3_1[31] = {
00161 0.02, 0.04, 0.06, 0.08, 0.10, 0.12, 0.14, 0.16, 0.18, 0.20,
00162 0.22, 0.25, 0.30, 0.35, 0.40, 0.45, 0.50, 0.60, 0.70, 0.80,
00163 0.90, 1.0, 2.0, 3.0, 5.0, 10.0, 20.0, 50.0, 100.0, 500.0, 100000.0};
00164
00165 const G4double G4PiNuclearCrossSection::al_m_t[31] = {
00166 532, 637, 832, 1057, 1207, 1230, 1210, 1174, 1133, 1095,
00167 1038, 970, 890, 807, 750, 710, 675, 665, 670, 673,
00168 678, 682, 618, 574, 546, 520, 507, 495, 488, 488, 488};
00169
00170 const G4double G4PiNuclearCrossSection::al_m_in[31] = {
00171 300, 360, 495, 665, 750, 765, 750, 730, 700, 660, 615, 570, 520, 490, 470,
00172 450, 448, 450, 450, 452, 456, 460, 408, 392, 376, 356, 347, 338, 332, 332,
00173 332};
00174
00175 const G4double G4PiNuclearCrossSection::al_p_t[21] = {
00176 225, 350, 616, 945, 1122, 1175, 1157, 1128, 1088, 1045,
00177 988, 935, 870, 787, 730, 690, 660, 652, 660, 668, 678};
00178
00179 const G4double G4PiNuclearCrossSection::al_p_in[21] = {
00180 120, 238, 390, 610, 712, 735, 720, 703, 655, 635, 590, 550, 505, 475, 455,
00181 438, 440, 445, 445, 450, 456};
00182
00183 const G4double G4PiNuclearCrossSection::ca_m_t[31] = {
00184 800, 980, 1240, 1460, 1570, 1600, 1580, 1535, 1475, 1425,
00185 1375, 1295, 1200, 1083, 1000, 948, 915, 895, 900, 908,
00186 915, 922, 856, 795, 740, 705, 682, 660, 660, 660, 660};
00187
00188 const G4double G4PiNuclearCrossSection::ca_m_in[31] = {
00189 470, 550, 620, 860, 955, 980, 960, 920, 860, 820, 780, 740, 665, 637, 615,
00190 600, 590, 590, 600, 608, 610, 615, 550, 525, 510, 488, 470, 450, 450, 450,
00191 450};
00192
00193 const G4double G4PiNuclearCrossSection::ca_p_t[23] = {
00194 275, 445, 790, 1195, 1440, 1485, 1475, 1435, 1385, 1335, 1295, 1245, 1160, 1050, 970,
00195 923, 895, 877, 887, 897, 904, 913, 855};
00196
00197 const G4double G4PiNuclearCrossSection::ca_p_in[23] = {
00198 160, 315, 500, 745, 870, 905, 900, 860, 810, 770, 740, 710, 640, 617, 595,
00199 585, 575, 575, 590, 600, 602, 608, 510};
00200
00201
00202 const G4double G4PiNuclearCrossSection::e4[32] = {
00203 0.02, 0.04, 0.06, 0.08, 0.10, 0.12, 0.14, 0.16, 0.18, 0.20, 0.22, 0.25, 0.30, 0.35, 0.40,
00204 0.45, 0.50, 0.55, 0.60, 0.70, 0.80, 0.90, 1, 2, 3, 5, 10, 20, 50, 100,
00205 500, 100000};
00206
00207 const G4double G4PiNuclearCrossSection::fe_m_t[32] = {
00208 1175, 1363, 1670, 1950, 2050, 2040, 1975, 1886, 1834, 1773, 1720, 1635,
00209 1474, 1380, 1269, 1225, 1182, 1162, 1159, 1162, 1178, 1190, 1197, 1102,
00210 1135, 975, 945, 925, 905, 905, 905, 905};
00211
00212 const G4double G4PiNuclearCrossSection::fe_m_in[32] = {
00213 625, 725, 910, 1180, 1275, 1250, 1200, 1150, 1100, 1040, 995, 925,
00214 825, 810, 780, 760, 745, 740, 740, 740, 750, 760, 765, 690,
00215 660, 635, 615, 600, 585, 585, 585, 585};
00216
00217 const G4double G4PiNuclearCrossSection::fe_p_t[25] = {
00218 330, 575, 1010, 1500, 1837, 1875, 1820, 1751, 1691, 1636, 1690, 1450,
00219 1396, 1305, 1219, 1190, 1148, 1138, 1134, 1144, 1163, 1175, 1183, 1198,
00220 1135};
00221
00222 const G4double G4PiNuclearCrossSection::fe_p_in[25] = {
00223 210, 410, 707, 1010, 1125, 1150, 1100, 1070, 1010, 960, 920, 776,
00224 780, 760, 750, 740, 720, 725, 725, 730, 740, 750, 755, 690,
00225 660};
00226
00227 const G4double G4PiNuclearCrossSection::cu_m_t[32] = {
00228 1400, 1600, 1875, 2088, 2200, 2220, 2175, 2125, 2075, 2012, 1950, 1855,
00229 1670, 1530, 1430, 1370, 1315, 1315, 1315, 1330, 1345, 1360, 1365, 1250,
00230 1185, 1128, 1070, 1035, 1010, 1010, 1010, 1010};
00231
00232 const G4double G4PiNuclearCrossSection::cu_m_in[32] = {
00233 725, 840, 1020, 1200, 1295, 1300, 1267, 1240, 1213, 1175, 1125, 1042,
00234 950, 900, 860, 840, 830, 832, 835, 840, 850, 860, 865, 785,
00235 735, 705, 680, 650, 630, 630, 630, 630};
00236
00237 const G4double G4PiNuclearCrossSection::cu_p_t[25] = {
00238 355, 605, 1120, 1630, 1940, 2010, 2010, 1980, 1925, 1895, 1830, 1730,
00239 1585, 1490, 1400, 1340, 1290, 1290, 1290, 1310, 1330, 1345, 1350, 1240,
00240 1185};
00241
00242 const G4double G4PiNuclearCrossSection::cu_p_in[25] = {
00243 230, 425, 780, 1025, 1155, 1190, 1190, 1180, 1125, 1100, 1050, 1000,
00244 900, 870, 835, 815, 810, 812, 815, 825, 840, 850, 855, 780,
00245 735};
00246
00247 const G4double G4PiNuclearCrossSection::e5[34] = {
00248 0.02, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.10, 0.12, 0.14, 0.16, 0.18, 0.20, 0.22, 0.25,
00249 0.30, 0.35, 0.40, 0.45, 0.50, 0.60, 0.70, 0.80, 0.90, 1, 2, 3, 5, 10, 20,
00250 50, 100, 500, 100000};
00251
00252 const G4double G4PiNuclearCrossSection::mo_m_t[34] = {
00253 2430, 2610, 2710, 2790, 2880, 2940, 2965, 2970, 2970, 2920, 2840, 2720,
00254 2570, 2500, 2365, 2200, 2050, 1926, 1825, 1768, 1749, 1750, 1778, 1789,
00255 1808, 1690, 1645, 1530, 1492, 1450, 1425, 1425, 1425, 1425};
00256
00257 const G4double G4PiNuclearCrossSection::mo_m_in[34] = {
00258 925, 1125, 1250, 1375, 1500, 1600, 1680, 1750, 1770, 1730, 1660, 1580,
00259 1500, 1450, 1330, 1250, 1190, 1140, 1100, 1075, 1075, 1070, 1088, 1095,
00260 1110, 1035, 1005, 940, 917, 880, 860, 860, 860, 860};
00261
00262 const G4double G4PiNuclearCrossSection::mo_p_t[27] = {
00263 410, 730, 1110, 1530, 1920, 2200, 2385, 2520, 2600, 2630, 2575, 2470,
00264 2320, 2285, 2185, 2053, 1945, 1852, 1776, 1719, 1710, 1716, 1746, 1759,
00265 1778, 1675, 1645};
00266
00267 const G4double G4PiNuclearCrossSection::mo_p_in[27] = {
00268 270, 540, 825, 975, 1140, 1285, 1400, 1480, 1555, 1580, 1525, 1470,
00269 1360, 1340, 1255, 1160, 1120, 1085, 1060, 1045, 1045, 1045, 1065, 1075,
00270 1090, 1025, 1005};
00271
00272 const G4double G4PiNuclearCrossSection::cd_m_t[34] = {
00273 3060, 3125, 3170, 3220, 3255, 3280, 3290, 3260, 3270, 3200, 3120, 3080,
00274 3090, 2920, 2810, 2640, 2362, 2230, 2115, 2050, 2020, 2025, 2040, 2070,
00275 2100, 1900, 1795, 1740, 1675, 1645, 1625, 1620, 1620, 1620};
00276
00277 const G4double G4PiNuclearCrossSection::cd_m_in[34]= {
00278 1025, 1275, 1440, 1625, 1740, 1800, 1880, 1920, 1980, 1920, 1850, 1810,
00279 1720, 1650, 1560, 1450, 1330, 1290, 1245, 1210, 1200, 1200, 1205, 1205,
00280 1230, 1130, 1085, 1060, 1000, 985, 975, 970, 970, 970};
00281
00282 const G4double G4PiNuclearCrossSection::cd_p_t[28] = {
00283 455, 780, 1170, 1700, 2120, 2400, 2600, 2720, 2820, 2840, 2800, 2760,
00284 2720, 2640, 2560, 2450, 2252, 2130, 2035, 1985, 1970, 1975, 2005, 2035,
00285 2070, 1880, 1795, 1740};
00286
00287 const G4double G4PiNuclearCrossSection::cd_p_in[28] = {
00288 310, 580, 880, 1060, 1270, 1400, 1530, 1610, 1660, 1680, 1640, 1600,
00289 1560, 1500, 1430, 1330, 1280, 1230, 1200, 1180, 1170, 1175, 1180, 1180,
00290 1210, 1120, 1085, 1060};
00291
00292 const G4double G4PiNuclearCrossSection::e6[35] = {
00293 0.02, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.10, 0.12, 0.14, 0.16, 0.18,
00294 0.20, 0.22, 0.25, 0.30, 0.35, 0.40, 0.45, 0.50, 0.55, 0.60, 0.70, 0.80,
00295 0.90, 1.0, 2.0, 3.0, 5.0, 10.0, 20.0, 50.0, 100.0, 500.0, 100000.0};
00296
00297 const G4double G4PiNuclearCrossSection::sn_m_t[35] = {
00298 3000, 3180, 3250, 3300, 3300, 3410, 3470, 3450, 3410, 3350, 3280, 3200,
00299 3120, 3050, 2900, 2630, 2500, 2325, 2190, 2100, 2060, 2055, 2055, 2055,
00300 2067, 2085, 2000, 1900, 1835, 1770, 1720, 1700, 1695, 1695, 1695};
00301
00302 const G4double G4PiNuclearCrossSection::sn_m_in[35] = {
00303 1050, 1350, 1520, 1650, 1800, 1980, 2070, 2120, 2090, 2050, 1980, 1920,
00304 1830, 1770, 1670, 1500, 1435, 1350, 1300, 1230, 1220, 1235, 1235, 1235,
00305 1237, 1240, 1160, 1120, 1090, 1065, 1040, 1020, 1015, 1015, 1015};
00306
00307 const G4double G4PiNuclearCrossSection::sn_p_t[29] = {
00308 465, 800, 1200, 1760, 2170, 2480, 2730, 2885, 2970, 2980, 2970, 2890,
00309 2840, 2790, 2620, 2450, 2335, 2205, 2080, 2020, 2010, 1990, 1990, 2015,
00310 2030, 2045, 1980, 1890, 1835};
00311
00312 const G4double G4PiNuclearCrossSection::sn_p_in[29] = {
00313 315, 590, 880, 1220, 1460, 1580, 1700, 1770, 1810, 1810, 1800, 1730,
00314 1680, 1630, 1530, 1400, 1335, 1270, 1210, 1180, 1190, 1190, 1190, 1205,
00315 1210, 1210, 1150, 1115, 1090};
00316
00317 const G4double G4PiNuclearCrossSection::w_m_t[35] = {
00318 5200, 5115, 5025, 4975, 4900, 4850, 4780, 4725, 4600, 4490, 4355, 4255,
00319 4125, 4040, 3830, 3580, 3330, 3110, 2955, 2860, 2852, 2845, 2885, 2900,
00320 2915, 2940, 2800, 2660, 2570, 2490, 2460, 2425, 2420, 2420, 2420};
00321
00322 const G4double G4PiNuclearCrossSection::w_m_in[35] = {
00323 1450, 1850, 2100, 2350, 2550, 2700, 2825, 2900, 2850, 2750, 2630, 2525,
00324 2400, 2300, 2200, 2070, 1880, 1770, 1715, 1680, 1680, 1680, 1685, 1690,
00325 1700, 1720, 1635, 1560, 1530, 1460, 1440, 1410, 1410, 1410, 1410};
00326
00327 const G4double G4PiNuclearCrossSection::w_p_t[30] = {
00328 480, 900, 1500, 2350, 3020, 3420, 3650, 3775, 3875, 3830, 3750, 3700,
00329 3630, 3550, 3550, 3290, 3070, 2890, 2840, 2730, 2725, 2720, 2770, 2805,
00330 2828, 2865, 2770, 2640, 2570, 2490};
00331
00332 const G4double G4PiNuclearCrossSection::w_p_in[30] = {
00333 325, 680, 990, 1500, 1850, 2150, 2250, 2300, 2350, 2330, 2280, 2230,
00334 2200, 2120, 2130, 1900, 1780, 1670, 1635, 1600, 1602, 1605, 1610, 1615,
00335 1630, 1660, 1620, 1550, 1530, 1460};
00336
00337 const G4double G4PiNuclearCrossSection::e7[35] = {
00338 0.02, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.10, 0.12, 0.14, 0.16, 0.18,
00339 0.20, 0.22, 0.25, 0.30, 0.35, 0.40, 0.45, 0.50, 0.55, 0.60, 0.70, 0.80,
00340 0.90, 1, 2, 3, 5, 10, 20, 50, 100, 500, 100000};
00341
00342 const G4double G4PiNuclearCrossSection::pb_m_t[35] = {
00343 5890, 5700, 5610, 5580, 5550, 5480, 5400, 5300, 5100, 4930, 4750, 4600,
00344 4400, 4280, 4170, 3915, 3650, 3470, 3260, 3150, 3120, 3070, 3085, 3100,
00345 3120, 3160, 3070, 2930, 2820, 2750, 2710, 2655, 2640, 2640, 2640};
00346
00347 const G4double G4PiNuclearCrossSection::pb_m_in[35] = {
00348 1575, 2025, 2300, 2575, 2850, 3000, 3115, 3180, 3080, 2940, 2800, 2670, 2550, 2450, 2370,
00349 2220, 2110, 2000, 1920, 1880, 1850, 1800, 1805, 1810, 1820, 1840, 1800, 1720, 1640, 1620,
00350 1570, 1530, 1530, 1530, 1530};
00351
00352 const G4double G4PiNuclearCrossSection::pb_p_t[30] = {
00353 515, 940, 1500, 2400, 3270, 3750, 4050, 4140, 4260, 4200, 4080, 3990, 3990, 3810, 3730,
00354 3520, 3370, 3186, 3110, 3010, 2990, 2985, 3005, 3020, 3040, 3080, 3020, 2905, 2790, 2750};
00355
00356 const G4double G4PiNuclearCrossSection::pb_p_in[30] = {
00357 348, 707, 1040, 1650, 2100, 2400, 2580, 2640, 2650, 2520, 2410, 2300, 2250, 2190, 2130,
00358 2000, 1930, 1870, 1830, 1790, 1770, 1765, 1775, 1780, 1790, 1800, 1775, 1710, 1620, 1620};
00359
00360 const G4double G4PiNuclearCrossSection::u_m_t[35] = {
00361 7080, 6830, 6650, 6530, 6400, 6280, 6100, 5840, 5660, 5520, 5330, 5160,
00362 4990, 4810, 4630, 4323, 4130, 3870, 3700, 3550, 3490, 3465, 3467, 3475,
00363 3495, 3515, 3440, 3360, 3150, 3040, 2985, 2955, 2940, 2940, 2940};
00364
00365 const G4double G4PiNuclearCrossSection::u_m_in[35] = {
00366 1740, 2220, 2500, 2820, 3080, 3300, 3420, 3500, 3420, 3330, 3200, 3060,
00367 2940, 2850, 2710, 2470, 2380, 2250, 2160, 2080, 2040, 2045, 2047, 2050,
00368 2055, 2060, 2010, 1980, 1830, 1780, 1735, 1710, 1700, 1700, 1700};
00369
00370 const G4double G4PiNuclearCrossSection::u_p_t[30] = {
00371 485, 960, 1580, 2700, 3550, 4050, 4320, 4420, 4620, 4660, 4580, 4470,
00372 4350, 4295, 4187, 3938, 3755, 3573, 3450, 3342, 3310, 3295, 3310, 3330,
00373 3375, 3405, 3350, 3338, 3135, 3040};
00374
00375 const G4double G4PiNuclearCrossSection::u_p_in[30] = {
00376 334, 720, 1020, 1560, 2100, 2300, 2550, 2700, 2880, 2880, 2760, 2660,
00377 2550, 2510, 2430, 2270, 2130, 2060, 2000, 1970, 1950, 1950, 1960, 1960,
00378 1970, 1980, 1950, 1978, 1830, 1780};
00379
00380
00381 G4PiNuclearCrossSection::G4PiNuclearCrossSection()
00382 : G4VCrossSectionDataSet("G4PiNuclearCrossSection"),
00383 fTotalXsc(0.0), fElasticXsc(0.0)
00384 {
00385 SetMinKinEnergy(0.0);
00386 SetMaxKinEnergy(99.9*TeV);
00387
00388 thePimData.push_back(new G4PiData(he_t, he_in, e1, 38));
00389 thePipData.push_back(new G4PiData(he_t, he_in, e1, 38));
00390 thePimData.push_back(new G4PiData(be_m_t, be_m_in, e1, 38));
00391 thePipData.push_back(new G4PiData(be_p_t, be_p_in, e1, 24));
00392 thePimData.push_back(new G4PiData(c_m_t, c_m_in, e2, 39));
00393 thePipData.push_back(new G4PiData(c_p_t, c_p_in, e2, 24));
00394 thePimData.push_back(new G4PiData(n_m_t, n_m_in, e2, 39));
00395 thePipData.push_back(new G4PiData(n_p_t, n_p_in, e2, 27));
00396 thePimData.push_back(new G4PiData(o_m_t, o_m_in, e3, 31));
00397 thePipData.push_back(new G4PiData(o_p_t, o_p_in, e3, 20));
00398 thePimData.push_back(new G4PiData(na_m_t, na_m_in, e3, 31));
00399 thePipData.push_back(new G4PiData(na_p_t, na_p_in, e3, 22));
00400 thePimData.push_back(new G4PiData(al_m_t, al_m_in, e3_1, 31));
00401 thePipData.push_back(new G4PiData(al_p_t, al_p_in, e3_1, 21));
00402 thePimData.push_back(new G4PiData(ca_m_t, ca_m_in, e3_1, 31));
00403 thePipData.push_back(new G4PiData(ca_p_t, ca_p_in, e3_1, 23));
00404 thePimData.push_back(new G4PiData(fe_m_t, fe_m_in, e4, 32));
00405 thePipData.push_back(new G4PiData(fe_p_t, fe_p_in, e4, 25));
00406 thePimData.push_back(new G4PiData(cu_m_t, cu_m_in, e4, 32));
00407 thePipData.push_back(new G4PiData(cu_p_t, cu_p_in, e4, 25));
00408 thePimData.push_back(new G4PiData(mo_m_t, mo_m_in, e5, 34));
00409 thePipData.push_back(new G4PiData(mo_p_t, mo_p_in, e5, 27));
00410 thePimData.push_back(new G4PiData(cd_m_t, cd_m_in, e5, 34));
00411 thePipData.push_back(new G4PiData(cd_p_t, cd_p_in, e5, 28));
00412 thePimData.push_back(new G4PiData(sn_m_t, sn_m_in, e6, 35));
00413 thePipData.push_back(new G4PiData(sn_p_t, sn_p_in, e6, 29));
00414 thePimData.push_back(new G4PiData(w_m_t, w_m_in, e6, 35));
00415 thePipData.push_back(new G4PiData(w_p_t, w_p_in, e6, 30));
00416 thePimData.push_back(new G4PiData(pb_m_t, pb_m_in, e7, 35));
00417 thePipData.push_back(new G4PiData(pb_p_t, pb_p_in, e7, 30));
00418 thePimData.push_back(new G4PiData(u_m_t, u_m_in, e7, 35));
00419 thePipData.push_back(new G4PiData(u_p_t, u_p_in, e7, 30));
00420
00421 theZ.push_back(2);
00422 theZ.push_back(4);
00423 theZ.push_back(6);
00424 theZ.push_back(7);
00425 theZ.push_back(8);
00426 theZ.push_back(11);
00427 theZ.push_back(13);
00428 theZ.push_back(20);
00429 theZ.push_back(26);
00430 theZ.push_back(29);
00431 theZ.push_back(42);
00432 theZ.push_back(48);
00433 theZ.push_back(50);
00434 theZ.push_back(74);
00435 theZ.push_back(82);
00436 theZ.push_back(92);
00437 }
00438
00439 G4PiNuclearCrossSection::
00440 ~G4PiNuclearCrossSection()
00441 {
00442 std::for_each(thePimData.begin(), thePimData.end(), G4PiData::Delete());
00443 std::for_each(thePipData.begin(), thePipData.end(), G4PiData::Delete());
00444 }
00445
00446 void
00447 G4PiNuclearCrossSection::CrossSectionDescription(std::ostream& outFile) const
00448 {
00449 outFile << "G4PiNuclearCrossSection calculates the pion inelastic cross\n"
00450 << "section for all nuclei heavier than hydrogen. It uses the\n"
00451 << "Barashenkov cross sections and is valid for all incident\n"
00452 << "energies.\n";
00453 }
00454
00455
00456 G4bool
00457 G4PiNuclearCrossSection::IsElementApplicable(const G4DynamicParticle*,
00458 G4int Z, const G4Material*)
00459 {
00460 return (1 < Z);
00461 }
00462
00463
00464 void G4PiNuclearCrossSection::BuildPhysicsTable(const G4ParticleDefinition& p)
00465 {
00466 if(&p == G4PionMinus::PionMinus() || &p == G4PionPlus::PionPlus()) { return; }
00467 throw G4HadronicException(__FILE__, __LINE__,"Is applicable only for pions");
00468 }
00469
00470 G4double
00471 G4PiNuclearCrossSection::GetElementCrossSection(const G4DynamicParticle* particle,
00472 G4int Z, const G4Material*)
00473 {
00474 G4double charge = particle->GetDefinition()->GetPDGCharge();
00475 G4double kineticEnergy = particle->GetKineticEnergy();
00476
00477
00478
00479 G4double result = 0;
00480
00481 size_t it = 0;
00482
00483 while(it < theZ.size() && Z > theZ[it]) it++;
00484
00485
00486
00487
00488 if(Z > theZ[it])
00489 {
00490 throw G4HadronicException(__FILE__, __LINE__,
00491 "Called G4PiNuclearCrossSection outside parametrization");
00492 }
00493 G4int Z1, Z2;
00494 G4double x1, x2, xt1, xt2;
00495 if( charge < 0 )
00496 {
00497 if( theZ[it] == Z )
00498 {
00499 result = thePimData[it]->ReactionXSection(kineticEnergy);
00500 fTotalXsc = thePimData[it]->TotalXSection(kineticEnergy);
00501
00502
00503
00504
00505 }
00506 else
00507 {
00508 x1 = thePimData[it-1]->ReactionXSection(kineticEnergy);
00509 xt1 = thePimData[it-1]->TotalXSection(kineticEnergy);
00510 Z1 = theZ[it-1];
00511 x2 = thePimData[it]->ReactionXSection(kineticEnergy);
00512 xt2 = thePimData[it]->TotalXSection(kineticEnergy);
00513 Z2 = theZ[it];
00514
00515 result = Interpolate(Z1, Z2, Z, x1, x2);
00516 fTotalXsc = Interpolate(Z1, Z2, Z, xt1, xt2);
00517
00518
00519
00520
00521
00522
00523
00524
00525
00526
00527 }
00528 }
00529 else
00530 {
00531 if(theZ[it]==Z)
00532 {
00533
00534 std::vector<G4PiData *> * theData = &thePimData;
00535 if(thePipData[it]->AppliesTo(kineticEnergy))
00536 {
00537 theData = &thePipData;
00538 }
00539 result = theData->operator[](it)->ReactionXSection(kineticEnergy);
00540 fTotalXsc = theData->operator[](it)->TotalXSection(kineticEnergy);
00541
00542
00543
00544
00545 }
00546 else
00547 {
00548 std::vector<G4PiData *> * theLData = &thePimData;
00549 if(thePipData[it-1]->AppliesTo(kineticEnergy))
00550 {
00551 theLData = &thePipData;
00552 }
00553 std::vector<G4PiData *> * theHData = &thePimData;
00554 if(thePipData[it]->AppliesTo(kineticEnergy))
00555 {
00556 theHData = &thePipData;
00557 }
00558 x1 = theLData->operator[](it-1)->ReactionXSection(kineticEnergy);
00559 xt1 = theLData->operator[](it-1)->TotalXSection(kineticEnergy);
00560 Z1 = theZ[it-1];
00561 x2 = theHData->operator[](it)->ReactionXSection(kineticEnergy);
00562 xt2 = theHData->operator[](it)->TotalXSection(kineticEnergy);
00563 Z2 = theZ[it];
00564
00565 result = Interpolate(Z1, Z2, Z, x1, x2);
00566 fTotalXsc = Interpolate(Z1, Z2, Z, xt1, xt2);
00567
00568
00569
00570
00571
00572
00573
00574
00575
00576
00577 }
00578 }
00579
00580
00581 fElasticXsc = fTotalXsc - result;
00582 if( fElasticXsc < 0.) fElasticXsc = 0.;
00583
00584 return result;
00585 }
00586
00587
00588 G4double G4PiNuclearCrossSection::
00589 Interpolate(G4int Z1, G4int Z2, G4int Z, G4double x1, G4double x2)
00590 {
00591
00592 static const G4double A[92] = {
00593 1.0001, 4.0000, 6.9241, 9.000, 10.801, 12.011, 14.004, 16.004, 19.000,
00594 20.188, 23.000, 24.320, 27.000, 28.109, 31.000, 32.094, 35.484, 39.985,
00595 39.135, 40.116, 45.000, 47.918, 50.998, 52.055, 55.000, 55.910, 59.000,
00596 58.760, 63.617, 65.468, 69.798, 72.691, 75.000, 79.042, 79.986, 83.887,
00597 85.557, 87.710, 89.000, 91.318, 93.000, 96.025, 98.000, 101.16, 103.00,
00598 106.51, 107.96, 112.51, 114.91, 118.81, 121.86, 127.70, 127.00, 131.39,
00599 133.00, 137.42, 139.00, 140.21, 141.00, 144.32, 145.00, 150.45, 152.04,
00600 157.33, 159.00, 162.57, 165.00, 167.32, 169.00, 173.10, 175.03, 178.54,
00601 181.00, 183.89, 186.25, 190.27, 192.25, 195.11, 197.00, 200.63, 204.41,
00602 207.24, 209.00, 209.00, 210.00, 222.00, 223.00, 226.00, 227.00, 232.00,
00603 231.00, 237.98};
00604
00605 static G4bool NeedInit=true;
00606 static G4double A75[92];
00607 if ( NeedInit )
00608 {
00609 for (G4int i=0; i<92; ++i)
00610 {
00611 A75[i]=std::pow(A[i],0.75);
00612 }
00613 NeedInit=false;
00614 }
00615
00616
00617 G4double r1 = x1 / A75[Z1-1] * A75[Z-1];
00618 G4double r2 = x2 / A75[Z2-1] * A75[Z-1];
00619 G4double result=0.5*(r1+r2);
00620
00621
00622 return result;
00623 }