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 // INCL++ intra-nuclear cascade model 27 // Alain Boudard, CEA-Saclay, France 28 // Joseph Cugnon, University of Liege, Belgium 29 // Jean-Christophe David, CEA-Saclay, France 30 // Pekka Kaitaniemi, CEA-Saclay, France, and H 31 // Sylvie Leray, CEA-Saclay, France 32 // Davide Mancusi, CEA-Saclay, France 33 // 34 #define INCLXX_IN_GEANT4_MODE 1 35 36 #include "globals.hh" 37 38 #include "G4INCLNNToMultiPionsChannel.hh" 39 #include "G4INCLKinematicsUtils.hh" 40 #include "G4INCLBinaryCollisionAvatar.hh" 41 #include "G4INCLRandom.hh" 42 #include "G4INCLGlobals.hh" 43 #include "G4INCLLogger.hh" 44 #include <algorithm> 45 #include "G4INCLPhaseSpaceGenerator.hh" 46 47 namespace G4INCL { 48 49 const G4double NNToMultiPionsChannel::angula 50 51 NNToMultiPionsChannel::NNToMultiPionsChannel 52 : npion(npi), 53 iso1(0), 54 iso2(0), 55 particle1(p1), 56 particle2(p2) 57 { 58 std::fill(isosp, isosp+4, 0); 59 } 60 61 NNToMultiPionsChannel::~NNToMultiPionsChanne 62 63 } 64 65 void NNToMultiPionsChannel::fillFinalState(F 66 // assert(npion > 0 && npion < 5); 67 68 iso1=ParticleTable::getIsospin(particle1 69 iso2=ParticleTable::getIsospin(particle2 70 71 ParticleList list; 72 list.push_back(particle1); 73 list.push_back(particle2); 74 fs->addModifiedParticle(particle1); 75 fs->addModifiedParticle(particle2); 76 77 isospinRepartition(); 78 79 const ParticleType tn1=ParticleTable::ge 80 particle1->setType(tn1); 81 const ParticleType tn2=ParticleTable::ge 82 particle2->setType(tn2); 83 const ThreeVector &rcolnucleon1 = partic 84 const ThreeVector &rcolnucleon2 = partic 85 const ThreeVector rcol = (rcolnucleon1+r 86 const ThreeVector zero; 87 for(G4int i=0; i<npion; ++i) { 88 const ParticleType pionType=ParticleTa 89 Particle *pion = new Particle(pionType 90 list.push_back(pion); 91 fs->addCreatedParticle(pion); 92 } 93 94 const G4double sqrtS = KinematicsUtils:: 95 G4int biasIndex = ((Random::shoot()<0.5) 96 PhaseSpaceGenerator::generateBiased(sqrt 97 98 } 99 100 void NNToMultiPionsChannel::isospinReparti 101 const G4double rjcd=Random::shoot(); 102 G4double p; 103 const G4int itot=iso1+iso2; 104 105 if (npion == 1) { 106 p=3.*rjcd; 107 if (p < 1.) pn_ppPim(); 108 else if (p < 2.) pn_pnPi0(); 109 else pn_nnPip(); 110 } 111 else if (npion == 2) { 112 if (itot == 2) { 113 p=20.*rjcd; 114 if (p >= 14.) pp_nnPipPip 115 else if (p >= 11.) pp_pnPipPi0 116 else if (p >= 7.) pp_ppPi0Pi0 117 else pp_ppPipPim( 118 } 119 else if (itot == -2) { 120 p=20.*rjcd; 121 if (p >= 14.) nn_ppPimPim 122 else if (p >= 11.) nn_pnPimPi0 123 else if (p >= 7.) nn_nnPi0Pi0 124 else nn_nnPipPim( 125 } 126 else { 127 G4double pp=Random::shoot(); 128 if (pp > 0.5) { 129 p=3.*rjcd; 130 if (p < 2.) { 131 pn_pnPipPim(); 132 } 133 else { 134 pn_pnPi0Pi0(); 135 } 136 } 137 else { 138 p=60.*rjcd; 139 if (p >= 51.) pn_nnPi 140 else if (p >= 33.) pn_pnPi 141 else if (p >= 9.) pn_pnPi 142 else pn_ppPim 143 } 144 } 145 } 146 else if (npion == 3) { 147 p=60.*rjcd; 148 if (itot == 2) { 149 if (p >= 42.) pp_nnPipPip 150 else if (p >= 39.) pp_pnPipPi0 151 else if (p >= 33.) pp_pnPipPip 152 else if (p >= 22.) pp_ppPi0Pi0 153 else pp_ppPipPimP 154 } 155 else if (itot == -2) { 156 if (p >= 42.) nn_ppPimPim 157 else if (p >= 39.) nn_pnPimPi0 158 else if (p >= 33.) nn_pnPipPim 159 else if (p >= 22.) nn_nnPi0Pi0 160 else nn_nnPipPimP 161 } 162 else { 163 if (p >= 57.) pn_nnPipPi0 164 else if (p >= 51.) pn_nnPipPip 165 else if (p >= 37.) pn_pnPi0Pi0 166 else if (p >= 9.) pn_pnPi0Pip 167 else if (p >= 6.) pn_ppPimPi0 168 else pn_ppPimPimP 169 170 } 171 } 172 else if (npion == 4) { 173 p=60.*rjcd; 174 if (itot == 2) { 175 if (p >= 48.) pp_nnPipPip 176 else if (p >= 42.) pp_nnPipPip 177 else if (p >= 36.) pp_pnPipPip 178 else if (p >= 33.) pp_pnPipPi0 179 else if (p >= 19.) pp_ppPipPip 180 else if (p >= 4.) pp_ppPipPi0 181 else pp_ppPi0Pi0P 182 } 183 else if (itot == -2) { 184 if (p >= 48.) nn_ppPipPim 185 else if (p >= 42.) nn_ppPi0Pi0 186 else if (p >= 36.) nn_pnPipPi0 187 else if (p >= 33.) nn_pnPi0Pi0 188 else if (p >= 19.) nn_nnPipPip 189 else if (p >= 4.) nn_nnPipPi0 190 else nn_nnPi0Pi0P 191 } 192 else { 193 G4double pp=Random::shoot(); 194 if (pp > 0.5) { 195 p=9.*rjcd; 196 if (p < 1.) pn_pnPi0P 197 else if (p < 5.) pn_pnPipP 198 else pn_pnPipPi 199 } 200 else { 201 if (p < 3.) pn_ppPi0 202 else if (p < 9.) pn_ppPip 203 else if (p < 15.) pn_pnPi0 204 else if (p < 35.) pn_pnPip 205 else if (p < 51.) pn_pnPip 206 else if (p < 54.) pn_nnPip 207 else pn_nnPipP 208 } 209 } 210 } 211 212 std::shuffle(isosp,isosp+npion,Random: 213 inter2Part(0.5); 214 } 215 216 217 void NNToMultiPionsChannel::pn_ppPim() { 218 isosp[0]=-2; 219 iso1=1; 220 iso2=1; 221 } 222 void NNToMultiPionsChannel::pn_pnPi0() { 223 isosp[0]=0; 224 } 225 void NNToMultiPionsChannel::pn_nnPip() { 226 isosp[0]=2; 227 iso1=-1; 228 iso2=-1; 229 } 230 void NNToMultiPionsChannel::pp_nnPipPip() 231 isosp[0]=2; 232 isosp[1]=2; 233 iso1=-1; 234 iso2=-1; 235 } 236 void NNToMultiPionsChannel::nn_ppPimPim() 237 isosp[0]=-2; 238 isosp[1]=-2; 239 iso1=1; 240 iso2=1; 241 } 242 void NNToMultiPionsChannel::pn_pnPipPim() 243 isosp[0]=2; 244 isosp[1]=-2; 245 } 246 void NNToMultiPionsChannel::pn_pnPi0Pi0() 247 isosp[0]=0; 248 isosp[1]=0; 249 } 250 void NNToMultiPionsChannel::pp_ppPipPim() 251 isosp[0]=2; 252 isosp[1]=-2; 253 } 254 void NNToMultiPionsChannel::nn_nnPipPim() 255 isosp[0]=2; 256 isosp[1]=-2; 257 } 258 void NNToMultiPionsChannel::pp_ppPi0Pi0() 259 isosp[0]=0; 260 isosp[1]=0; 261 } 262 void NNToMultiPionsChannel::nn_nnPi0Pi0() 263 isosp[0]=0; 264 isosp[1]=0; 265 } 266 void NNToMultiPionsChannel::pp_pnPipPi0() 267 isosp[0]=2; 268 isosp[1]=0; 269 iso1=1; 270 iso2=-1; 271 } 272 void NNToMultiPionsChannel::pn_ppPimPi0() 273 isosp[0]=-2; 274 isosp[1]=0; 275 iso1=1; 276 iso2=1; 277 } 278 void NNToMultiPionsChannel::pn_nnPipPi0() 279 isosp[0]=2; 280 isosp[1]=0; 281 iso1=-1; 282 iso2=-1; 283 } 284 void NNToMultiPionsChannel::nn_pnPimPi0() 285 isosp[0]=-2; 286 isosp[1]=0; 287 iso1=1; 288 iso2=-1; 289 } 290 void NNToMultiPionsChannel::pp_pnPipPi0Pi0 291 isosp[0]=2; 292 isosp[1]=0; 293 isosp[2]=0; 294 iso1=1; 295 iso2=-1; 296 } 297 void NNToMultiPionsChannel::nn_pnPimPi0Pi0 298 isosp[0]=-2; 299 isosp[1]=0; 300 isosp[2]=0; 301 iso1=1; 302 iso2=-1; 303 } 304 void NNToMultiPionsChannel::pn_nnPipPi0Pi0 305 isosp[0]=2; 306 isosp[1]=0; 307 isosp[2]=0; 308 iso1=-1; 309 iso2=-1; 310 } 311 void NNToMultiPionsChannel::pp_ppPipPimPi0 312 isosp[0]=2; 313 isosp[1]=-2; 314 isosp[2]=0; 315 } 316 void NNToMultiPionsChannel::nn_nnPipPimPi0 317 isosp[0]=2; 318 isosp[1]=-2; 319 isosp[2]=0; 320 } 321 void NNToMultiPionsChannel::pp_ppPi0Pi0Pi0 322 isosp[0]=0; 323 isosp[1]=0; 324 isosp[2]=0; 325 } 326 void NNToMultiPionsChannel::nn_nnPi0Pi0Pi0 327 isosp[0]=0; 328 isosp[1]=0; 329 isosp[2]=0; 330 } 331 void NNToMultiPionsChannel::pp_pnPipPipPim 332 isosp[0]=2; 333 isosp[1]=2; 334 isosp[2]=-2; 335 iso1=1; 336 iso2=-1; 337 } 338 void NNToMultiPionsChannel::pp_nnPipPipPi0 339 isosp[0]=2; 340 isosp[1]=2; 341 isosp[2]=0; 342 iso1=-1; 343 iso2=-1; 344 } 345 void NNToMultiPionsChannel::pn_ppPimPi0Pi0 346 isosp[0]=-2; 347 isosp[1]=0; 348 isosp[2]=0; 349 iso1=1; 350 iso2=1; 351 } 352 void NNToMultiPionsChannel::pn_ppPimPimPip 353 isosp[0]=-2; 354 isosp[1]=-2; 355 isosp[2]=2; 356 iso1=1; 357 iso2=1; 358 } 359 void NNToMultiPionsChannel::pn_pnPi0PipPim 360 isosp[0]=0; 361 isosp[1]=2; 362 isosp[2]=-2; 363 } 364 void NNToMultiPionsChannel::pn_pnPi0Pi0Pi0 365 isosp[0]=0; 366 isosp[1]=0; 367 isosp[2]=0; 368 } 369 void NNToMultiPionsChannel::pn_nnPipPipPim 370 isosp[0]=2; 371 isosp[1]=2; 372 isosp[2]=-2; 373 iso1=-1; 374 iso2=-1; 375 } 376 void NNToMultiPionsChannel::nn_pnPipPimPim 377 isosp[0]=2; 378 isosp[1]=-2; 379 isosp[2]=-2; 380 iso1=1; 381 iso2=-1; 382 } 383 void NNToMultiPionsChannel::nn_ppPimPimPi0 384 isosp[0]=-2; 385 isosp[1]=-2; 386 isosp[2]=0; 387 iso1=1; 388 iso2=1; 389 } 390 void NNToMultiPionsChannel::pp_nnPipPipPi0 391 isosp[0]=2; 392 isosp[1]=2; 393 isosp[2]=0; 394 isosp[3]=0; 395 iso1=-1; 396 iso2=-1; 397 } 398 void NNToMultiPionsChannel::pp_nnPipPipPip 399 isosp[0]=2; 400 isosp[1]=2; 401 isosp[2]=2; 402 isosp[3]=-2; 403 iso1=-1; 404 iso2=-1; 405 } 406 void NNToMultiPionsChannel::nn_ppPi0Pi0Pim 407 isosp[0]=0; 408 isosp[1]=0; 409 isosp[2]=-2; 410 isosp[3]=-2; 411 iso1=1; 412 iso2=1; 413 } 414 void NNToMultiPionsChannel::nn_ppPipPimPim 415 isosp[0]=2; 416 isosp[1]=-2; 417 isosp[2]=-2; 418 isosp[3]=-2; 419 iso1=1; 420 iso2=1; 421 } 422 void NNToMultiPionsChannel::pp_ppPi0Pi0Pi0 423 isosp[0]=0; 424 isosp[1]=0; 425 isosp[2]=0; 426 isosp[3]=0; 427 } 428 void NNToMultiPionsChannel::nn_nnPi0Pi0Pi0 429 isosp[0]=0; 430 isosp[1]=0; 431 isosp[2]=0; 432 isosp[3]=0; 433 } 434 void NNToMultiPionsChannel::pn_pnPi0Pi0Pi0 435 isosp[0]=0; 436 isosp[1]=0; 437 isosp[2]=0; 438 isosp[3]=0; 439 } 440 void NNToMultiPionsChannel::pp_ppPipPi0Pi0 441 isosp[0]=2; 442 isosp[1]=0; 443 isosp[2]=0; 444 isosp[3]=-2; 445 } 446 void NNToMultiPionsChannel::nn_nnPipPi0Pi0 447 isosp[0]=2; 448 isosp[1]=0; 449 isosp[2]=0; 450 isosp[3]=-2; 451 } 452 void NNToMultiPionsChannel::pn_pnPipPi0Pi0 453 isosp[0]=2; 454 isosp[1]=0; 455 isosp[2]=0; 456 isosp[3]=-2; 457 } 458 void NNToMultiPionsChannel::pp_ppPipPipPim 459 isosp[0]=2; 460 isosp[1]=2; 461 isosp[2]=-2; 462 isosp[3]=-2; 463 } 464 void NNToMultiPionsChannel::nn_nnPipPipPim 465 isosp[0]=2; 466 isosp[1]=2; 467 isosp[2]=-2; 468 isosp[3]=-2; 469 } 470 void NNToMultiPionsChannel::pn_pnPipPipPim 471 isosp[0]=2; 472 isosp[1]=2; 473 isosp[2]=-2; 474 isosp[3]=-2; 475 } 476 void NNToMultiPionsChannel::pp_pnPipPi0Pi0 477 isosp[0]=2; 478 isosp[1]=0; 479 isosp[2]=0; 480 isosp[3]=0; 481 iso1=1; 482 iso2=-1; 483 } 484 void NNToMultiPionsChannel::pn_nnPipPi0Pi0 485 isosp[0]=2; 486 isosp[1]=0; 487 isosp[2]=0; 488 isosp[3]=0; 489 iso1=-1; 490 iso2=-1; 491 } 492 void NNToMultiPionsChannel::pp_nnPipPi0Pi0 493 isosp[0]=2; 494 isosp[1]=0; 495 isosp[2]=0; 496 isosp[3]=0; 497 iso1=-1; 498 iso2=-1; 499 } 500 void NNToMultiPionsChannel::pp_pnPipPipPi0 501 isosp[0]=2; 502 isosp[1]=2; 503 isosp[2]=0; 504 isosp[3]=-2; 505 iso1=1; 506 iso2=-1; 507 } 508 void NNToMultiPionsChannel::pn_nnPipPipPi0 509 isosp[0]=2; 510 isosp[1]=2; 511 isosp[2]=0; 512 isosp[3]=-2; 513 iso1=-1; 514 iso2=-1; 515 } 516 void NNToMultiPionsChannel::pp_nnPipPipPi0 517 isosp[0]=2; 518 isosp[1]=2; 519 isosp[2]=0; 520 isosp[3]=-2; 521 iso1=-1; 522 iso2=-1; 523 } 524 void NNToMultiPionsChannel::nn_pnPi0Pi0Pi0 525 isosp[0]=0; 526 isosp[1]=0; 527 isosp[2]=0; 528 isosp[3]=-2; 529 iso1=1; 530 iso2=-1; 531 } 532 void NNToMultiPionsChannel::pn_ppPi0Pi0Pi0 533 isosp[0]=0; 534 isosp[1]=0; 535 isosp[2]=0; 536 isosp[3]=-2; 537 iso1=1; 538 iso2=1; 539 } 540 void NNToMultiPionsChannel::nn_pnPipPi0Pim 541 isosp[0]=2; 542 isosp[1]=0; 543 isosp[2]=-2; 544 isosp[3]=-2; 545 iso1=1; 546 iso2=-1; 547 } 548 void NNToMultiPionsChannel::pn_ppPipPi0Pim 549 isosp[0]=2; 550 isosp[1]=0; 551 isosp[2]=-2; 552 isosp[3]=-2; 553 iso1=1; 554 iso2=1; 555 } 556 557 void NNToMultiPionsChannel::inter2Part(con 558 559 if (Random::shoot() < p) std::swap(iso 560 561 } 562 563 564 } 565