Geant4 Cross Reference |
1 // 1 2 // ******************************************* 3 // * License and Disclaimer 4 // * 5 // * The Geant4 software is copyright of th 6 // * the Geant4 Collaboration. It is provided 7 // * conditions of the Geant4 Software License 8 // * LICENSE and available at http://cern.ch/ 9 // * include a list of copyright holders. 10 // * 11 // * Neither the authors of this software syst 12 // * institutes,nor the agencies providing fin 13 // * work make any representation or warran 14 // * regarding this software system or assum 15 // * use. Please see the license in the file 16 // * for the full disclaimer and the limitatio 17 // * 18 // * This code implementation is the result 19 // * technical work of the GEANT4 collaboratio 20 // * By using, copying, modifying or distri 21 // * any work based on the software) you ag 22 // * use in resulting scientific publicati 23 // * acceptance of all terms of the Geant4 Sof 24 // ******************************************* 25 // 26 // 27 // Calculation of the total, elastic and inela 28 // based on Barashenkov parametrisations of pi 29 // 30 // 16.08.06 V.Ivanchenko - first implementatio 31 // J.P Wellisch class 32 // 22.01.07 V.Ivanchenko - add cross section i 33 // 05.03.07 V.Ivanchenko - fix weight for inte 34 // 13.03.07 V.Ivanchenko - cleanup at low ener 35 // 11.09.09 V.Ivanchenko - fixed bug in interp 36 // 37 38 #include "G4UPiNuclearCrossSection.hh" 39 #include "G4SystemOfUnits.hh" 40 #include "G4PhysicsFreeVector.hh" 41 #include "G4PionMinus.hh" 42 #include "G4PionPlus.hh" 43 #include "G4PhysicsTable.hh" 44 #include "G4NistManager.hh" 45 #include "G4HadronXSDataTable.hh" 46 47 G4int G4UPiNuclearCrossSection::theZ[NZ] = 48 {2,4,6,7,8,11,13,20,26,29,42,48,50,74,82,92}; 49 G4int G4UPiNuclearCrossSection::idxZ[93] = {0} 50 G4double G4UPiNuclearCrossSection::theA[NZ] 51 G4double G4UPiNuclearCrossSection::APower[93] 52 53 G4PhysicsTable* G4UPiNuclearCrossSection::piPl 54 G4PhysicsTable* G4UPiNuclearCrossSection::piPl 55 G4PhysicsTable* G4UPiNuclearCrossSection::piMi 56 G4PhysicsTable* G4UPiNuclearCrossSection::piMi 57 58 G4UPiNuclearCrossSection::G4UPiNuclearCrossSec 59 : G4VCrossSectionDataSet("G4UPiNuclearCrossSe 60 { 61 piPlus = G4PionPlus::PionPlus(); 62 piMinus = G4PionMinus::PionMinus(); 63 elow = 20.0*CLHEP::MeV; 64 65 if (idxZ[0] == 0) { LoadData(); } 66 } 67 68 G4bool 69 G4UPiNuclearCrossSection::IsElementApplicable( 70 G4int Z, const G4Material*) 71 { 72 return (1 < Z); 73 } 74 75 G4double G4UPiNuclearCrossSection::Interpolate 76 G4int Z, G4int A, G4double e, const G 77 { 78 G4double res = 0.0; 79 G4double ekin = std::max(e, elow); 80 G4int iz = std::min(Z, 92); 81 G4int idx = idxZ[iz]; 82 std::size_t jdx = (std::size_t)(std::max(eki 83 //G4cout << "Interpolate: Z= " << iz << " A= 84 // << " jdx= " << jdx << " Ekin= " << ekin < 85 if(idx < 0 || 2 == iz) { 86 res = ((*table)[std::abs(idx)])->Value(eki 87 //G4cout << "1: jdx= " << jdx << G4endl; 88 } else { 89 G4int iz2 = theZ[idx]; 90 G4double x2 = (((*table)[idx])->Value(ekin 91 //G4cout << "2: jdx= " << jdx << G4endl; 92 G4int iz1 = theZ[idx-1]; 93 G4double x1 = (((*table)[idx-1])->Value(ek 94 G4double w1 = ((G4double)A - theA[idx-1])/ 95 res = w1*x2 + (1.0 - w1)*x1; 96 } 97 //G4cout << " res(nb)= " << res/CLHEP::barn 98 return res; 99 } 100 101 void G4UPiNuclearCrossSection::AddDataSet(cons 102 const G4double* tot, 103 const G4double* in, 104 const G4double* e, 105 G4int n) 106 { 107 G4PhysicsFreeVector* pvin = new G4PhysicsFre 108 G4PhysicsFreeVector* pvel = new G4PhysicsFre 109 for (G4int i=0; i<n; ++i) { 110 pvin->PutValues(i, e[i]*CLHEP::GeV, in[i]* 111 pvel->PutValues(i, e[i]*CLHEP::GeV, std::m 112 } 113 if (spline) { 114 pvin->FillSecondDerivatives(); 115 pvel->FillSecondDerivatives(); 116 } 117 if (p == "pi+") { 118 piPlusInelastic->push_back(pvin); 119 piPlusElastic->push_back(pvel); 120 } else { 121 piMinusInelastic->push_back(pvin); 122 piMinusElastic->push_back(pvel); 123 } 124 } 125 126 void G4UPiNuclearCrossSection::DumpPhysicsTabl 127 { 128 if(&p == piPlus) { 129 G4cout << "### G4UPiNuclearCrossSection El 130 G4cout << *piPlusElastic << G4endl; 131 G4cout << "### G4UPiNuclearCrossSection In 132 G4cout << *piPlusInelastic << G4endl; 133 } else if(&p == piMinus) { 134 G4cout << "### G4UPiNuclearCrossSection El 135 G4cout << *piMinusElastic << G4endl; 136 G4cout << "### G4UPiNuclearCrossSection In 137 G4cout << *piMinusInelastic << G4endl; 138 } 139 } 140 141 void G4UPiNuclearCrossSection::BuildPhysicsTab 142 { 143 if(&p != piPlus && &p != piMinus) { 144 G4ExceptionDescription ed; 145 ed << "This cross section is applicable on 146 << p.GetParticleName() << G4endl; 147 G4Exception("G4UPiNuclearCrossSection::Bui 148 FatalException, ed); 149 } 150 } 151 152 void G4UPiNuclearCrossSection::LoadData() 153 { 154 idxZ[0] = 1; 155 idxZ[1] = idxZ[2] = 0; 156 G4NistManager* nist = G4NistManager::Instanc 157 G4Pow* g4pow = G4Pow::GetInstance(); 158 for(G4int i=0; i<NZ; ++i) { 159 theA[i] = nist->GetAtomicMassAmu(theZ[i]); 160 } 161 for(G4int i=1; i<93; ++i) { 162 APower[i] = g4pow->powA(nist->GetAtomicMas 163 } 164 G4int idx = 1; 165 for(G4int i=3; i<93; ++i) { 166 if(theZ[idx] == i) { 167 idxZ[i] = -idx; 168 ++idx; 169 } else { 170 idxZ[i] = idx; 171 } 172 } 173 174 piPlusElastic = new G4PhysicsTable(); 175 piPlusInelastic = new G4PhysicsTable(); 176 piMinusElastic = new G4PhysicsTable(); 177 piMinusInelastic = new G4PhysicsTable(); 178 auto ptr = G4HadronXSDataTable::Instance(); 179 ptr->AddTable(piPlusElastic); 180 ptr->AddTable(piPlusInelastic); 181 ptr->AddTable(piMinusElastic); 182 ptr->AddTable(piMinusInelastic); 183 184 static const G4double e1[38] = { 185 0.02, 0.04, 0.06, 0.08, 0.1, 0.12, 0.13, 0 186 0.17, 0.18, 0.19, 0.2, 0.22,0.24, 0.26, 0 187 0.4, 0.45, 0.5, 0.55, 0.6, 0.7, 0.8, 0 188 3, 5, 10, 20, 50, 100, 500, 10 189 static const G4double e2[39] = { 190 0.02, 0.04, 0.06, 0.08, 0.1, 0.11, 0.12, 191 0.16, 0.17, 0.18, 0.2, 0.22, 0.24, 0.26, 192 0.4, 0.45, 0.5, 0.55, 0.575,0.6, 0.7, 193 2, 3, 5, 10, 20, 50, 100, 194 static const G4double e3[31] = { 195 0.02, 0.04, 0.06, 0.08, 0.1, 0.12, 0.14, 196 0.22, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 197 0.9, 1, 2, 3, 5, 10, 20, 198 static const G4double e4[32] = { 199 0.02, 0.04, 0.06, 0.08, 0.1, 0.12, 0.14, 200 0.22, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0 201 0.8, 0.9, 1, 2, 3, 5, 10, 202 static const G4double e5[34] = { 203 0.02, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 204 0.16, 0.18, 0.2, 0.22, 0.25, 0.3, 0.35, 205 0.6, 0.7, 0.8, 0.9, 1, 2, 3, 206 static const G4double e6[35] = { 207 0.02, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 208 0.16, 0.18, 0.2, 0.22, 0.25, 0.3, 0.35, 209 0.55, 0.6, 0.7, 0.8, 0.9, 1, 2, 210 211 static const G4double he_t[38] = { 212 40, 70, 108, 152, 208, 276, 300, 320, 213 332, 328, 322, 310, 288, 260, 240, 216, 214 125, 112, 108.5, 109, 110.5, 117, 123, 128.5 215 96, 87, 85, 83.5, 83.5, 83.5, 83.5, 83.5 216 static const G4double he_in[38] = { 217 18, 38, 62, 98, 136, 176, 190, 200, 218 212, 208, 204, 196, 176, 164, 150, 134, 219 90, 85,82.5, 83.5, 86.5, 93, 97.5,100, 220 77, 75, 74, 72.5, 72.5, 72.5, 72.5, 72.5 221 static const G4double be_m_t[38] = { 222 150, 210, 294, 396, 520, 600, 623, 635, 223 630, 615, 600, 576, 540, 504, 470, 435, 224 294, 258, 236, 230, 233, 244, 257, 270, 225 230, 215, 205, 194, 188, 186, 186, 186}; 226 static const G4double be_m_in[38] = { 227 90, 126, 177, 240, 320, 380, 400, 410, 228 400, 387, 371, 360, 333, 312, 285, 260, 229 198, 187, 182, 180, 182, 187, 193, 203, 230 172, 165, 159, 155, 144, 144, 144, 144}; 231 static const G4double be_p_t[38] = { 232 96, 150, 222, 320, 430, 514, 545, 565, 233 564, 552, 535, 522, 490, 462, 432, 398, 234 276, 248, 232, 230, 233, 244, 257, 270, 235 230, 215, 205, 194, 188, 186, 186, 186}; 236 static const G4double be_p_in[38] = { 237 60, 95, 142, 194, 262, 319, 345, 361, 238 354, 350, 330, 319, 298, 280, 258, 237, 239 189, 183, 182, 180, 182, 187, 193, 203, 240 172, 165, 159, 155, 144, 144, 144, 144}; 241 242 static const G4double c_m_t[39] = { 243 204, 260, 366, 517, 630, 673, 694, 704, 244 706, 694, 676, 648, 616, 584, 548, 518, 245 376, 342, 323, 310, 312, 313, 319, 333, 246 310, 290, 268, 250, 245, 237, 234, 234, 247 static const G4double c_m_in[39] = { 248 128, 160, 224, 315, 388, 416, 430, 438, 249 440, 432, 416, 400, 380, 354, 320, 304, 250 246, 240, 233, 232, 233, 234, 238, 246, 251 220, 210, 198, 187, 183, 176, 174, 174, 252 static const G4double c_p_t[39] = { 253 140, 192, 294, 428, 594, 642, 662, 687, 254 684, 672, 656, 630, 598, 567, 533, 504, 255 369, 336, 319, 310, 312, 313, 319, 333, 256 310, 290, 268, 250, 245, 237, 234, 234, 257 static const G4double c_p_in[39] = { 258 94, 132, 184, 260, 370, 398, 408, 420, 259 424, 416, 400, 386, 366, 340, 308, 294, 260 241, 236, 231, 232, 233, 234, 238, 246, 261 220, 210, 198, 187, 183, 176, 174, 174, 262 static const G4double n_m_t[39] = { 263 246, 308, 424, 590, 729, 776, 800, 821, 264 800, 778, 768, 728, 690, 654, 615, 584, 265 430, 393, 373, 367, 368, 370, 375, 388, 266 364, 337, 310, 291, 275, 268, 268, 268, 267 static const G4double n_m_in[39] = { 268 155, 188, 256, 360, 456, 492, 512, 526, 269 504, 491, 475, 450, 425, 396, 376, 360, 270 282, 270, 265, 265, 266, 268, 273, 280, 271 256, 237, 226, 218, 208, 202, 202, 202, 272 static const G4double n_p_t[39] = { 273 150, 212, 328, 500, 680, 735, 762, 781, 274 770, 748, 740, 706, 672, 633, 600, 569, 275 419, 385, 368, 364, 366, 368, 375, 388, 276 364, 337, 310, 291, 275, 268, 268, 268, 277 static const G4double n_p_in[39] = { 278 90, 140, 208, 300, 426, 467, 490, 504, 279 484, 474, 460, 437, 413, 381, 365, 350, 280 276, 267, 263, 264, 265, 267, 273, 280, 281 256, 237, 226, 218, 208, 202, 202, 202, 282 283 static const G4double o_m_t[31] = { 284 280, 360, 500, 685, 812, 861, 870, 865, 285 755, 700, 600, 537, 493, 468, 441, 436, 286 460, 463, 432, 385, 350, 325, 312, 307, 287 static const G4double o_m_in[31] = { 288 190, 207, 300, 420, 500, 540, 550, 542, 289 460, 423, 360, 339, 321, 314, 312, 314, 290 328, 330, 300, 275, 250, 240, 229, 225, 291 static const G4double o_p_t[31] = { 292 170, 240, 390, 570, 740, 818, 830, 822, 293 725, 675, 585, 525, 483, 458, 444, 447, 294 460, 463, 432, 385, 350, 325, 312, 307, 295 static const G4double o_p_in[31] = { 296 100, 145, 240, 340, 470, 518, 530, 522, 297 448, 412, 350, 330, 316, 310, 308, 311, 298 328, 330, 300, 275, 250, 240, 229, 225, 299 static const G4double na_m_t[31] = { 300 450, 545, 705, 910, 1020, 1075, 1087, 1080 301 943, 885, 790, 700, 650, 610, 585, 575 302 600, 610, 556, 524, 494, 458, 445, 429 303 static const G4double na_m_in[31] = { 304 275, 315, 413, 545, 620, 660, 670, 662 305 570, 520, 465, 420, 410, 395, 390, 400 306 420, 422, 372, 348, 330, 320, 310, 294 307 static const G4double na_p_t[31] = { 308 210, 320, 530, 795, 960, 1035, 1050, 1040 309 918, 865, 773, 685, 636, 598, 575, 565 310 598, 610, 556, 524, 494, 458, 445, 429 311 static const G4double na_p_in[31] = { 312 115, 210, 340, 495, 585, 630, 645, 637 313 550, 505, 455, 410, 401, 388, 383, 393 314 418, 422, 372, 348, 330, 320, 310, 294 315 static const G4double al_m_t[31] = { 316 532, 637, 832, 1057, 1207, 1230, 1210, 1174 317 1038, 970, 890, 807, 750, 710, 675, 665 318 678, 682, 618, 574, 546, 520, 507, 495 319 static const G4double al_m_in[31] = { 320 300, 360, 495, 665, 750, 765, 750, 730 321 615, 570, 520, 490, 470, 450, 448, 450 322 456, 460, 408, 392, 376, 356, 347, 338 323 static const G4double al_p_t[31] = { 324 225, 350, 616, 945, 1122, 1175, 1157, 1128 325 988, 935, 870, 787, 730, 690, 660, 652 326 678, 682, 618, 574, 546, 520, 507, 495 327 static const G4double al_p_in[31] = { 328 120, 238, 390, 610, 712, 735, 720, 703 329 590, 550, 505, 475, 455, 438, 440, 445 330 456, 460, 408, 392, 376, 356, 347, 338 331 332 static const G4double ca_m_t[31] = { 333 800, 980, 1240, 1460, 1570, 1600, 1580, 1535 334 1375,1295, 1200, 1083, 1000, 948, 915, 895 335 915, 922, 856, 795, 740, 705, 682, 660 336 static const G4double ca_m_in[31] = { 337 470, 550, 620, 860, 955, 960, 920, 860 338 740, 665, 637, 615, 600, 590, 580, 580 339 610, 615, 550, 525, 510, 488, 470, 450 340 static const G4double ca_p_t[31] = { 341 275, 445, 790, 1195, 1440, 1485, 1475, 1435 342 1295,1245,1160, 1050, 970, 923, 895, 877 343 904, 913, 855, 795, 740, 705, 682, 660 344 static const G4double ca_p_in[31] = { 345 160, 315, 500, 745, 870, 905, 900, 860 346 740, 710, 640, 617, 595, 585, 575, 575 347 602, 608, 550, 525, 510, 488, 470, 450 348 // ca data may have typo 349 350 static const G4double fe_m_t[32] = { 351 1175, 1363, 1670, 1950, 2050, 2040, 1975, 188 352 1720, 1635, 1474, 1380, 1269, 1225, 1182, 116 353 1178, 1190, 1197, 1102, 1135, 975, 945, 92 354 905, 905}; 355 static const G4double fe_m_in[32] = { 356 625, 725, 910, 1180, 1275, 1250, 1200, 115 357 995, 925, 825, 810, 780, 760, 745, 74 358 750, 760, 765, 690, 660, 635, 615, 60 359 585, 585}; 360 static const G4double fe_p_t[32] = { 361 330, 575, 1010, 1500, 1837, 1875, 1820, 175 362 1690, 1450, 1396, 1305, 1219, 1190, 1148, 113 363 1163, 1175, 1183, 1198, 1135, 975, 945, 92 364 905, 905}; 365 static const G4double fe_p_in[32] = { 366 210, 410, 707, 1010, 1125, 1150, 1100, 107 367 920, 776, 780, 760, 750, 740, 720, 72 368 740, 750, 755, 690, 660, 635, 615, 60 369 585, 585}; 370 static const G4double cu_m_t[32] = { 371 1400, 1600, 1875, 2088, 2200, 2220, 2175, 212 372 1950, 1855, 1670, 1530, 1430, 1370, 1315, 131 373 1345, 1360, 1365, 1250, 1185, 1128, 1070, 103 374 1010, 1010}; 375 static const G4double cu_m_in[32] = { 376 725, 840, 1020, 1200, 1295, 1300, 1267, 124 377 1125, 1042, 950, 900, 860, 840, 830, 83 378 850, 860, 865, 785, 735, 705, 680, 65 379 630, 630}; 380 static const G4double cu_p_t[32] = { 381 355, 605, 1120, 1630, 1940, 2010, 2010, 198 382 1830, 1730, 1585, 1490, 1400, 1340, 1290, 129 383 1330, 1345, 1350, 1240, 1185, 1128, 1070, 103 384 1010, 1010}; 385 static const G4double cu_p_in[32] = { 386 230, 425, 780, 1025, 1155, 1190, 1190, 11 387 1050, 1000, 900, 870, 835, 815, 810, 8 388 840, 850, 855, 780, 735, 705, 680, 6 389 630, 630}; 390 391 static const G4double mo_m_t[34] = { 392 2430, 2610, 2710, 2790, 2880, 2940, 2965, 29 393 2840, 2720, 2570, 2500, 2365, 2200, 2050, 19 394 1749, 1750, 1778, 1789, 1808, 1690, 1645, 15 395 1425, 1425, 1425, 1425}; 396 static const G4double mo_m_in[34] = { 397 925, 1125, 1250, 1375, 1500, 1600, 1680, 17 398 1660, 1580, 1500, 1450, 1330, 1250, 1190, 11 399 1075, 1070, 1088, 1095, 1110, 1035, 1005, 9 400 860, 860, 860, 860}; 401 static const G4double mo_p_t[34] = { 402 410, 730, 1110, 1530, 1920, 2200, 2385, 25 403 2575, 2470, 2320, 2285, 2185, 2053, 1945, 18 404 1710, 1716, 1746, 1759, 1778, 1675, 1645, 15 405 1425, 1425, 1425, 1425}; 406 static const G4double mo_p_in[34] = { 407 270, 540, 825, 975, 1140, 1285, 1400, 14 408 1525, 1470, 1360, 1340, 1255, 1160, 1120, 10 409 1045, 1045, 1065, 1075, 1090, 1025, 1005, 9 410 860, 860, 860, 860}; 411 static const G4double cd_m_t[34] = { 412 3060, 3125, 3170, 3220, 3255, 3280, 3290, 32 413 3120, 3080, 3090, 2920, 2810, 2640, 2362, 22 414 2020, 2025, 2040, 2070, 2100, 1900, 1795, 17 415 1625, 1620, 1620, 1620}; 416 static const G4double cd_m_in[34]= { 417 1025, 1275, 1440, 1625, 1740, 1800, 1880, 19 418 1850, 1810, 1720, 1650, 1560, 1450, 1330, 12 419 1200, 1200, 1205, 1205, 1230, 1130, 1085, 10 420 975, 970, 970, 970}; 421 static const G4double cd_p_t[34] = { 422 455, 780, 1170, 1700, 2120, 2400, 2600, 27 423 2800, 2760, 2720, 2640, 2560, 2450, 2252, 21 424 1970, 1975, 2005, 2035, 2070, 1880, 1795, 17 425 1625, 1620, 1620, 1620}; 426 static const G4double cd_p_in[34] = { 427 310, 580, 880, 1060, 1270, 1400, 1530, 16 428 1640, 1600, 1560, 1500, 1430, 1330, 1280, 12 429 1170, 1175, 1180, 1180, 1210, 1120, 1085, 10 430 975, 970, 970, 970}; 431 432 static const G4double sn_m_t[35] = { 433 3000, 3180, 3250, 3300, 3300, 3410, 3470, 34 434 3280, 3200, 3120, 3050, 2900, 2630, 2500, 23 435 2060, 2055, 2055, 2055, 2067, 2085, 2000, 19 436 1720, 1700, 1695, 1695, 1695}; 437 static const G4double sn_m_in[35] = { 438 1050, 1350, 1520, 1650, 1800, 1980, 2070, 21 439 1980, 1920, 1830, 1770, 1670, 1500, 1435, 13 440 1220, 1235, 1235, 1235, 1237, 1240, 1160, 11 441 1040, 1020, 1015, 1015, 1015}; 442 static const G4double sn_p_t[35] = { 443 465, 800, 1200, 1760, 2170, 2480, 2730, 28 444 2970, 2890, 2840, 2790, 2620, 2450, 2335, 22 445 2010, 1990, 1990, 2015, 2030, 2045, 1980, 18 446 1720, 1700, 1695, 1695, 1695}; 447 static const G4double sn_p_in[35] = { 448 315, 590, 880, 1220, 1460, 1580, 1700, 17 449 1800, 1730, 1680, 1630, 1530, 1400, 1335, 12 450 1190, 1190, 1190, 1205, 1210, 1210, 1150, 11 451 1040, 1020, 1015, 1015, 1015}; 452 static const G4double w_m_t[35] = { 453 5200, 5115, 5025, 4975, 4900, 4850, 4780, 47 454 4355, 4255, 4125, 4040, 3830, 3580, 3330, 31 455 2852, 2845, 2885, 2900, 2915, 2940, 2800, 26 456 2460, 2425, 2420, 2420, 2420}; 457 static const G4double w_m_in[35] = { 458 1450, 1850, 2100, 2350, 2550, 2700, 2825, 29 459 2630, 2525, 2400, 2300, 2200, 2070, 1880, 17 460 1680, 1680, 1685, 1690, 1700, 1720, 1635, 15 461 1440, 1410, 1410, 1410, 1410}; 462 static const G4double w_p_t[35] = { 463 480, 900, 1500, 2350, 3020, 3420, 3650, 37 464 3750, 3700, 3630, 3550, 3550, 3290, 3070, 28 465 2725, 2720, 2770, 2805, 2828, 2865, 2770, 26 466 2460, 2425, 2420, 2420, 2420}; 467 static const G4double w_p_in[35] = { 468 325, 680, 990, 1500, 1850, 2150, 2250, 23 469 2280, 2230, 2200, 2120, 2130, 1900, 1780, 16 470 1602, 1605, 1610, 1615, 1630, 1660, 1620, 15 471 1440, 1410, 1410, 1410, 1410}; 472 473 static const G4double pb_m_t[35] = { 474 5890, 5700, 5610, 5580, 5550, 5480, 5400, 53 475 4750, 4600, 4400, 4280, 4170, 3915, 3650, 34 476 3120, 3070, 3085, 3100, 3120, 3160, 3070, 29 477 2710, 2655, 2640, 2640, 2640}; 478 static const G4double pb_m_in[35] = { 479 1575, 2025, 2300, 2575, 2850, 3000, 3115, 31 480 2800, 2670, 2550, 2450, 2370, 2220, 2110, 20 481 1850, 1800, 1805, 1810, 1820, 1840, 1800, 17 482 1570, 1530, 1530, 1530, 1530}; 483 static const G4double pb_p_t[35] = { 484 515, 940, 1500, 2400, 3270, 3750, 4050, 41 485 4080, 3990, 3990, 3810, 3730, 3520, 3370, 31 486 2990, 2985, 3005, 3020, 3040, 3080, 3020, 29 487 2710, 2655, 2640, 2640, 2640}; 488 static const G4double pb_p_in[35] = { 489 348, 707, 1040, 1650, 2100, 2400, 2580, 26 490 2410, 2300, 2250, 2190, 2130, 2000, 1930, 18 491 1770, 1765, 1775, 1780, 1790, 1800, 1775, 17 492 1570, 1530, 1530, 1530, 1530}; 493 static const G4double u_m_t[35] = { 494 7080, 6830, 6650, 6530, 6400, 6280, 6100, 58 495 5330, 5160, 4990, 4810, 4630, 4323, 4130, 38 496 3490, 3465, 3467, 3475, 3495, 3515, 3440, 33 497 2985, 2955, 2940, 2940, 2940}; 498 static const G4double u_m_in[35] = { 499 1740, 2220, 2500, 2820, 3080, 3300, 3420, 35 500 3200, 3060, 2940, 2850, 2710, 2470, 2380, 22 501 2040, 2045, 2047, 2050, 2055, 2060, 2010, 19 502 1735, 1710, 1700, 1700, 1700}; 503 static const G4double u_p_t[35] = { 504 485, 960, 1580, 2700, 3550, 4050, 4320, 44 505 4580, 4470, 4350, 4295, 4187, 3938, 3755, 35 506 3310, 3295, 3310, 3330, 3375, 3405, 3350, 33 507 2985, 2955, 2940, 2940, 2940}; 508 static const G4double u_p_in[35] = { 509 334, 720, 1020, 1560, 2100, 2300, 2550, 27 510 2760, 2660, 2550, 2510, 2430, 2270, 2130, 20 511 1950, 1950, 1960, 1960, 1970, 1980, 1950, 19 512 1735, 1710, 1700, 1700, 1700}; 513 514 AddDataSet("pi-",he_t, he_in, e1, 38); 515 AddDataSet("pi+",he_t, he_in, e1, 38); 516 AddDataSet("pi-",be_m_t, be_m_in, e1, 38); 517 AddDataSet("pi+",be_p_t, be_p_in, e1, 38); 518 AddDataSet("pi-",c_m_t, c_m_in, e2, 39); 519 AddDataSet("pi+",c_p_t, c_p_in, e2, 39); 520 AddDataSet("pi-",n_m_t, n_m_in, e2, 39); 521 AddDataSet("pi+",n_p_t, n_p_in, e2, 39); 522 AddDataSet("pi-",o_m_t, o_m_in, e3, 31); 523 AddDataSet("pi+",o_p_t, o_p_in, e3, 31); 524 AddDataSet("pi-",na_m_t, na_m_in, e3, 31); 525 AddDataSet("pi+",na_p_t, na_p_in, e3, 31); 526 AddDataSet("pi-",al_m_t, al_m_in, e3, 31); 527 AddDataSet("pi+",al_p_t, al_p_in, e3, 31); 528 AddDataSet("pi-",ca_m_t, ca_m_in, e3, 31); 529 AddDataSet("pi+",ca_p_t, ca_p_in, e3, 31); 530 AddDataSet("pi-",fe_m_t, fe_m_in, e4, 32); 531 AddDataSet("pi+",fe_p_t, fe_p_in, e4, 32); 532 AddDataSet("pi-",cu_m_t, cu_m_in, e4, 32); 533 AddDataSet("pi+",cu_p_t, cu_p_in, e4, 32); 534 AddDataSet("pi-",mo_m_t, mo_m_in, e5, 34); 535 AddDataSet("pi+",mo_p_t, mo_p_in, e5, 34); 536 AddDataSet("pi-",cd_m_t, cd_m_in, e5, 34); 537 AddDataSet("pi+",cd_p_t, cd_p_in, e5, 34); 538 AddDataSet("pi-",sn_m_t, sn_m_in, e6, 35); 539 AddDataSet("pi+",sn_p_t, sn_p_in, e6, 35); 540 AddDataSet("pi-",w_m_t, w_m_in, e6, 35); 541 AddDataSet("pi+",w_p_t, w_p_in, e6, 35); 542 AddDataSet("pi-",pb_m_t, pb_m_in, e6, 35); 543 AddDataSet("pi+",pb_p_t, pb_p_in, e6, 35); 544 AddDataSet("pi-",u_m_t, u_m_in, e6, 35); 545 AddDataSet("pi+",u_p_t, u_p_in, e6, 35); 546 } 547 548 void G4UPiNuclearCrossSection::CrossSectionDes 549 { 550 outFile << "G4UPiNuclearCrossSection calcula 551 << "inelastic cross sections for pio 552 << "heavier than hydrogen. It is ba 553 << "parameterization and is valid fo 554 } 555 556