Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer 3 // * License and Disclaimer * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. 9 // * include a list of copyright holders. * 10 // * 10 // * * 11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. * 17 // * 17 // * * 18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************* 24 // ******************************************************************** 25 // 25 // 26 // INCL++ intra-nuclear cascade model 26 // INCL++ intra-nuclear cascade model 27 // Alain Boudard, CEA-Saclay, France 27 // Alain Boudard, CEA-Saclay, France 28 // Joseph Cugnon, University of Liege, Belgium 28 // Joseph Cugnon, University of Liege, Belgium 29 // Jean-Christophe David, CEA-Saclay, France 29 // Jean-Christophe David, CEA-Saclay, France 30 // Pekka Kaitaniemi, CEA-Saclay, France, and H 30 // Pekka Kaitaniemi, CEA-Saclay, France, and Helsinki Institute of Physics, Finland 31 // Sylvie Leray, CEA-Saclay, France 31 // Sylvie Leray, CEA-Saclay, France 32 // Davide Mancusi, CEA-Saclay, France 32 // Davide Mancusi, CEA-Saclay, France 33 // 33 // 34 #define INCLXX_IN_GEANT4_MODE 1 34 #define INCLXX_IN_GEANT4_MODE 1 35 35 36 #include "globals.hh" 36 #include "globals.hh" 37 37 38 #include "G4INCLNNToNSKpiChannel.hh" 38 #include "G4INCLNNToNSKpiChannel.hh" 39 #include "G4INCLKinematicsUtils.hh" 39 #include "G4INCLKinematicsUtils.hh" 40 #include "G4INCLBinaryCollisionAvatar.hh" 40 #include "G4INCLBinaryCollisionAvatar.hh" 41 #include "G4INCLRandom.hh" 41 #include "G4INCLRandom.hh" 42 #include "G4INCLGlobals.hh" 42 #include "G4INCLGlobals.hh" 43 #include "G4INCLLogger.hh" 43 #include "G4INCLLogger.hh" 44 #include <algorithm> 44 #include <algorithm> 45 #include "G4INCLPhaseSpaceGenerator.hh" 45 #include "G4INCLPhaseSpaceGenerator.hh" 46 46 47 namespace G4INCL { 47 namespace G4INCL { 48 48 49 const G4double NNToNSKpiChannel::angularSlop 49 const G4double NNToNSKpiChannel::angularSlope = 2.; // What is the exact effect? Sould be check 50 50 51 NNToNSKpiChannel::NNToNSKpiChannel(Particle 51 NNToNSKpiChannel::NNToNSKpiChannel(Particle *p1, Particle *p2) 52 : particle1(p1), particle2(p2) 52 : particle1(p1), particle2(p2) 53 {} 53 {} 54 54 55 NNToNSKpiChannel::~NNToNSKpiChannel(){} 55 NNToNSKpiChannel::~NNToNSKpiChannel(){} 56 56 57 void NNToNSKpiChannel::fillFinalState(FinalS 57 void NNToNSKpiChannel::fillFinalState(FinalState *fs) { 58 58 59 // pp (36) pn (36) 59 // pp (36) pn (36) 60 // 60 // 61 // pp -> p pi+ S- K+ (9) 61 // pp -> p pi+ S- K+ (9) 62 // pp -> p pi+ S0 K0 (9) 62 // pp -> p pi+ S0 K0 (9) 63 // pp -> p pi0 S+ K0 (4) 63 // pp -> p pi0 S+ K0 (4) 64 // pp -> n pi+ S+ K0 (2) 64 // pp -> n pi+ S+ K0 (2) 65 // pp -> p pi0 S0 K+ (4) 65 // pp -> p pi0 S0 K+ (4) 66 // pp -> n pi+ S0 K+ (2) 66 // pp -> n pi+ S0 K+ (2) 67 // pp -> p pi- S+ K+ (2) 67 // pp -> p pi- S+ K+ (2) 68 // pp -> n pi0 S+ K+ (4) 68 // pp -> n pi0 S+ K+ (4) 69 69 70 // pn -> p pi0 S- K+ (4) 70 // pn -> p pi0 S- K+ (4) 71 // pn -> n pi+ S- K+ (2) 71 // pn -> n pi+ S- K+ (2) 72 // pn -> p pi0 S0 K0 (2) 72 // pn -> p pi0 S0 K0 (2) 73 // pn -> n pi+ S0 K0 (1) 73 // pn -> n pi+ S0 K0 (1) 74 // pn -> p pi+ S- K0 (9) 74 // pn -> p pi+ S- K0 (9) 75 75 76 const G4double sqrtS = KinematicsUtils::to 76 const G4double sqrtS = KinematicsUtils::totalEnergyInCM(particle1, particle2); 77 77 78 const G4int iso = ParticleTable::getIsospi 78 const G4int iso = ParticleTable::getIsospin(particle1->getType()) + ParticleTable::getIsospin(particle2->getType()); 79 79 80 ParticleType KaonType; 80 ParticleType KaonType; 81 ParticleType PionType; 81 ParticleType PionType; 82 82 83 G4double rdm = Random::shoot(); 83 G4double rdm = Random::shoot(); 84 84 85 if(iso == 2){ 85 if(iso == 2){ 86 if(rdm * 36. < 9.){ 86 if(rdm * 36. < 9.){ 87 KaonType = KPlus; 87 KaonType = KPlus; 88 PionType = PiPlus; 88 PionType = PiPlus; 89 particle2->setType(SigmaMinus); 89 particle2->setType(SigmaMinus); 90 } 90 } 91 else if(rdm * 36. < 18.){ 91 else if(rdm * 36. < 18.){ 92 KaonType = KZero; 92 KaonType = KZero; 93 PionType = PiPlus; 93 PionType = PiPlus; 94 particle2->setType(SigmaZero); 94 particle2->setType(SigmaZero); 95 } 95 } 96 else if(rdm * 36. < 22.){ 96 else if(rdm * 36. < 22.){ 97 KaonType = KZero; 97 KaonType = KZero; 98 PionType = PiZero; 98 PionType = PiZero; 99 particle2->setType(SigmaPlus); 99 particle2->setType(SigmaPlus); 100 } 100 } 101 else if(rdm * 36. < 24.){ 101 else if(rdm * 36. < 24.){ 102 KaonType = KZero; 102 KaonType = KZero; 103 PionType = PiPlus; 103 PionType = PiPlus; 104 particle1->setType(Neutron); 104 particle1->setType(Neutron); 105 particle2->setType(SigmaPlus); 105 particle2->setType(SigmaPlus); 106 } 106 } 107 else if(rdm * 36. < 28.){ 107 else if(rdm * 36. < 28.){ 108 KaonType = KPlus; 108 KaonType = KPlus; 109 PionType = PiZero; 109 PionType = PiZero; 110 particle2->setType(SigmaZero); 110 particle2->setType(SigmaZero); 111 } 111 } 112 else if(rdm * 36. < 30.){ 112 else if(rdm * 36. < 30.){ 113 KaonType = KPlus; 113 KaonType = KPlus; 114 PionType = PiPlus; 114 PionType = PiPlus; 115 particle1->setType(Neutron); 115 particle1->setType(Neutron); 116 particle2->setType(SigmaZero); 116 particle2->setType(SigmaZero); 117 } 117 } 118 else if(rdm * 36. < 32.){ 118 else if(rdm * 36. < 32.){ 119 KaonType = KPlus; 119 KaonType = KPlus; 120 PionType = PiMinus; 120 PionType = PiMinus; 121 particle2->setType(SigmaPlus); 121 particle2->setType(SigmaPlus); 122 } 122 } 123 else{ 123 else{ 124 KaonType = KPlus; 124 KaonType = KPlus; 125 PionType = PiZero; 125 PionType = PiZero; 126 particle1->setType(Neutron); 126 particle1->setType(Neutron); 127 particle2->setType(SigmaPlus); 127 particle2->setType(SigmaPlus); 128 } 128 } 129 129 130 } 130 } 131 else if(iso == -2){ 131 else if(iso == -2){ 132 if(rdm * 36. < 9.){ 132 if(rdm * 36. < 9.){ 133 KaonType = KZero; 133 KaonType = KZero; 134 PionType = PiMinus; 134 PionType = PiMinus; 135 particle2->setType(SigmaPlus); 135 particle2->setType(SigmaPlus); 136 } 136 } 137 else if(rdm * 36. < 18.){ 137 else if(rdm * 36. < 18.){ 138 KaonType = KPlus; 138 KaonType = KPlus; 139 PionType = PiMinus; 139 PionType = PiMinus; 140 particle2->setType(SigmaZero); 140 particle2->setType(SigmaZero); 141 } 141 } 142 else if(rdm * 36. < 22.){ 142 else if(rdm * 36. < 22.){ 143 KaonType = KPlus; 143 KaonType = KPlus; 144 PionType = PiZero; 144 PionType = PiZero; 145 particle2->setType(SigmaMinus); 145 particle2->setType(SigmaMinus); 146 } 146 } 147 else if(rdm * 36. < 24.){ 147 else if(rdm * 36. < 24.){ 148 KaonType = KPlus; 148 KaonType = KPlus; 149 PionType = PiMinus; 149 PionType = PiMinus; 150 particle1->setType(Proton); 150 particle1->setType(Proton); 151 particle2->setType(SigmaMinus); 151 particle2->setType(SigmaMinus); 152 } 152 } 153 else if(rdm * 36. < 28.){ 153 else if(rdm * 36. < 28.){ 154 KaonType = KZero; 154 KaonType = KZero; 155 PionType = PiZero; 155 PionType = PiZero; 156 particle2->setType(SigmaZero); 156 particle2->setType(SigmaZero); 157 } 157 } 158 else if(rdm * 36. < 30.){ 158 else if(rdm * 36. < 30.){ 159 KaonType = KZero; 159 KaonType = KZero; 160 PionType = PiMinus; 160 PionType = PiMinus; 161 particle1->setType(Proton); 161 particle1->setType(Proton); 162 particle2->setType(SigmaZero); 162 particle2->setType(SigmaZero); 163 } 163 } 164 else if(rdm * 36. < 32.){ 164 else if(rdm * 36. < 32.){ 165 KaonType = KZero; 165 KaonType = KZero; 166 PionType = PiPlus; 166 PionType = PiPlus; 167 particle2->setType(SigmaMinus); 167 particle2->setType(SigmaMinus); 168 } 168 } 169 else{ 169 else{ 170 KaonType = KZero; 170 KaonType = KZero; 171 PionType = PiZero; 171 PionType = PiZero; 172 particle1->setType(Proton); 172 particle1->setType(Proton); 173 particle2->setType(SigmaMinus); 173 particle2->setType(SigmaMinus); 174 } 174 } 175 175 176 } 176 } 177 else if(rdm*36. < 4.){ 177 else if(rdm*36. < 4.){ 178 KaonType = KPlus; 178 KaonType = KPlus; 179 PionType = PiZero; 179 PionType = PiZero; 180 particle1->setType(Proton); 180 particle1->setType(Proton); 181 particle2->setType(SigmaMinus); 181 particle2->setType(SigmaMinus); 182 } 182 } 183 else if(rdm*36. < 6.){ 183 else if(rdm*36. < 6.){ 184 KaonType = KZero; 184 KaonType = KZero; 185 PionType = PiZero; 185 PionType = PiZero; 186 particle1->setType(Neutron); 186 particle1->setType(Neutron); 187 particle2->setType(SigmaPlus); 187 particle2->setType(SigmaPlus); 188 } 188 } 189 else if(rdm*36. < 8.){ 189 else if(rdm*36. < 8.){ 190 KaonType = KPlus; 190 KaonType = KPlus; 191 PionType = PiPlus; 191 PionType = PiPlus; 192 particle1->setType(Neutron); 192 particle1->setType(Neutron); 193 particle2->setType(SigmaMinus); 193 particle2->setType(SigmaMinus); 194 } 194 } 195 else if(rdm*36. < 9.){ 195 else if(rdm*36. < 9.){ 196 KaonType = KZero; 196 KaonType = KZero; 197 PionType = PiMinus; 197 PionType = PiMinus; 198 particle1->setType(Proton); 198 particle1->setType(Proton); 199 particle2->setType(SigmaPlus); 199 particle2->setType(SigmaPlus); 200 } 200 } 201 else if(rdm*36. < 18.){ 201 else if(rdm*36. < 18.){ 202 KaonType = KZero; 202 KaonType = KZero; 203 PionType = PiZero; 203 PionType = PiZero; 204 particle1->setType(Proton); 204 particle1->setType(Proton); 205 particle2->setType(SigmaZero); 205 particle2->setType(SigmaZero); 206 } 206 } 207 else if(rdm*36. < 27.){ 207 else if(rdm*36. < 27.){ 208 KaonType = KPlus; 208 KaonType = KPlus; 209 PionType = PiZero; 209 PionType = PiZero; 210 particle1->setType(Neutron); 210 particle1->setType(Neutron); 211 particle2->setType(SigmaZero); 211 particle2->setType(SigmaZero); 212 } 212 } 213 else if(rdm*36. < 28.){ 213 else if(rdm*36. < 28.){ 214 KaonType = KZero; 214 KaonType = KZero; 215 PionType = PiPlus; 215 PionType = PiPlus; 216 particle1->setType(Neutron); 216 particle1->setType(Neutron); 217 particle2->setType(SigmaZero); 217 particle2->setType(SigmaZero); 218 } 218 } 219 else if(rdm*36. < 30.){ 219 else if(rdm*36. < 30.){ 220 KaonType = KPlus; 220 KaonType = KPlus; 221 PionType = PiMinus; 221 PionType = PiMinus; 222 particle1->setType(Proton); 222 particle1->setType(Proton); 223 particle2->setType(SigmaZero); 223 particle2->setType(SigmaZero); 224 } 224 } 225 else if(rdm*36. < 32.){ 225 else if(rdm*36. < 32.){ 226 KaonType = KZero; 226 KaonType = KZero; 227 PionType = PiPlus; 227 PionType = PiPlus; 228 particle1->setType(Proton); 228 particle1->setType(Proton); 229 particle2->setType(SigmaMinus); 229 particle2->setType(SigmaMinus); 230 } 230 } 231 else{ 231 else{ 232 KaonType = KPlus; 232 KaonType = KPlus; 233 PionType = PiMinus; 233 PionType = PiMinus; 234 particle1->setType(Neutron); 234 particle1->setType(Neutron); 235 particle2->setType(SigmaPlus); 235 particle2->setType(SigmaPlus); 236 } 236 } 237 237 238 ParticleList list; 238 ParticleList list; 239 list.push_back(particle1); 239 list.push_back(particle1); 240 list.push_back(particle2); 240 list.push_back(particle2); 241 const ThreeVector &rcol1 = particle1->getP 241 const ThreeVector &rcol1 = particle1->getPosition(); 242 const ThreeVector &rcol2 = particle2->getP 242 const ThreeVector &rcol2 = particle2->getPosition(); 243 const ThreeVector zero; 243 const ThreeVector zero; 244 Particle *pion = new Particle(PionType,zer 244 Particle *pion = new Particle(PionType,zero,rcol1); 245 Particle *kaon = new Particle(KaonType,zer 245 Particle *kaon = new Particle(KaonType,zero,rcol2); 246 list.push_back(kaon); 246 list.push_back(kaon); 247 list.push_back(pion); 247 list.push_back(pion); 248 248 249 if(Random::shoot()<0.5) PhaseSpaceGenerato << 249 PhaseSpaceGenerator::generateBiased(sqrtS, list, 0, angularSlope); 250 else PhaseSpaceGenerator::generateBiased(s << 251 250 252 INCL_DEBUG("NNToNSKpi " << (kaon->getMomen 251 INCL_DEBUG("NNToNSKpi " << (kaon->getMomentum().theta()) * 180. / G4INCL::Math::pi << '\n'); 253 252 254 fs->addModifiedParticle(particle1); 253 fs->addModifiedParticle(particle1); 255 fs->addModifiedParticle(particle2); 254 fs->addModifiedParticle(particle2); 256 fs->addCreatedParticle(kaon); 255 fs->addCreatedParticle(kaon); 257 fs->addCreatedParticle(pion); 256 fs->addCreatedParticle(pion); 258 257 259 } 258 } 260 } 259 } 261 260