Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4INCLNNToMissingStrangenessChannel.cc
이 파일의 문서화 페이지로 가기
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 // 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 Helsinki Institute of Physics, Finland
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 
39 #include "G4INCLKinematicsUtils.hh"
41 #include "G4INCLRandom.hh"
42 #include "G4INCLGlobals.hh"
43 #include "G4INCLLogger.hh"
44 #include <algorithm>
46 
47 namespace G4INCL {
48 
50 
52  : particle1(p1), particle2(p2)
53  {}
54 
56 
58 
60  const G4double pLab = 0.001*KinematicsUtils::momentumInLab(particle1, particle2); // GeV
61 // assert(sqrtS > 3100); // ! > 3047.34. Not supposed to be under 3,626 GeV.
62 
64 // assert(iso == -2 || iso == 0 || iso == 2);
65  G4int iso_system = 0;
66  G4int available_iso = 0;
67  G4int nbr_pions = 0;
68  G4int min_pions = 0;
69  G4int max_pions = 0;
70 
71  G4double rdm = Random::shoot();
72 
73  G4int nbr_particle = 3;
74 
75  if(rdm < 0.35){
76  // N-K-Lambda chosen
78  available_iso = 2;
79  min_pions = 3;
81  }
82  else if((iso == 0 && rdm < 0.55) || rdm < 0.5){
83  // N-N-K-Kb chosen
84  nbr_particle++;
85  available_iso = 4;
86  min_pions = 1;
88  }
89  else{
90  // N-K-Sigma chosen
91  available_iso = 4;
92  min_pions = 3;
94  }
95  /* Gaussian noise + mean value nbr pions fonction energy (choice)*/
96  G4double intermediaire = min_pions + Random::gauss(2) + std::sqrt(pLab-5);
97  nbr_pions = std::min(max_pions,std::max(min_pions,G4int(intermediaire )));
98 
99  available_iso += nbr_pions*2;
100  nbr_particle += nbr_pions;
101 
102  ParticleList list;
103  ParticleType PionType = PiZero;
104  const ThreeVector &rcol1 = particle1->getPosition();
105  const ThreeVector zero;
106  Particle *kaon = new Particle(KPlus,zero,rcol1);
107 
108  for(Int_t i=0; i<nbr_pions; i++){
109  Particle *pion = new Particle(PionType,zero,rcol1);
110  if(available_iso-std::abs(iso-iso_system) >= 4){ // pp(pn) pip:0.40(0.3) piz:0.35(0.4) pim:0.25(0.3)
111  rdm = Random::shoot();
112  if(((iso == 0) && rdm < 0.3) || ((iso == -2) && rdm < 0.40) || rdm < 0.25){
113  pion->setType(PiMinus);
114  iso_system -= 2;
115  available_iso -= 2;
116  }
117  else if(((iso == 0) && rdm < 0.7) || ((iso == -2) && rdm < 0.75) || rdm < 0.60){
118  pion->setType(PiZero);
119  available_iso -= 2;
120  }
121  else{
122  pion->setType(PiPlus);
123  iso_system += 2;
124  available_iso -= 2;
125  }
126  }
127  else if(available_iso-std::abs(iso-iso_system) == 2){
128  rdm = Random::shoot();
129  // pn pp too high (nn too low) -> PiMinus(PiPlus) pp too low (nn too high) -> PiPlus(PiMinus)
130  if(((iso == 0) && (rdm*0.7 < 0.3)) || ((rdm*0.60 < 0.25) && (Math::sign(iso-iso_system)*2-iso != 0)) || ((rdm*0.75 < 0.40) && (Math::sign(iso-iso_system)*2+iso != 0) && (iso != 0))){
131  pion->setType(ParticleTable::getPionType(Math::sign(iso-iso_system)*2));
132  iso_system += Math::sign(iso-iso_system)*2;
133  available_iso -= 2;
134  }
135  else{
136  PionType = PiZero;
137  available_iso -= 2;
138  }
139  }
140  else if(available_iso-std::abs(iso-iso_system) == 0){
141  pion->setType(ParticleTable::getPionType(Math::sign(iso-iso_system)*2));
142  iso_system += Math::sign(iso-iso_system)*2;
143  available_iso -= 2;
144  }
145  else INCL_ERROR("Pion Generation Problem in NNToMissingStrangenessChannel" << '\n');
146  list.push_back(pion);
147  }
148 
149  if(particle2->isLambda()){ // N-K-Lambda
150 // assert(available_iso == 2);
151  if(std::abs(iso-iso_system) == 2){
152  particle1->setType(ParticleTable::getNucleonType((iso-iso_system)/2));
153  kaon->setType(ParticleTable::getKaonType((iso-iso_system)/2));
154  }
155  else if(std::abs(iso-iso_system) == 0){
156  rdm = G4int(Random::shoot()*2.)*2-1;
159  }
160  else INCL_ERROR("Isospin non-conservation in NNToMissingStrangenessChannel" << '\n');
161  }
162  else if(min_pions == 1){ // N-N-K-Kb chosen
163 // assert(available_iso == 4);
164  Particle *antikaon = new Particle(KMinus,zero,rcol1);
165  if(std::abs(iso-iso_system) == 4){
166  particle1->setType(ParticleTable::getNucleonType((iso-iso_system)/4));
167  particle2->setType(ParticleTable::getNucleonType((iso-iso_system)/4));
168  kaon->setType(ParticleTable::getKaonType((iso-iso_system)/4));
169  antikaon->setType(ParticleTable::getAntiKaonType((iso-iso_system)/4));
170  }
171  else if(std::abs(iso-iso_system) == 2){ // choice: kaon type free, nucleon type fixed
172  rdm = G4int(Random::shoot()*2.)*2-1;
173  particle1->setType(ParticleTable::getNucleonType((iso-iso_system)/2));
174  particle2->setType(ParticleTable::getNucleonType((iso-iso_system)/2));
176  antikaon->setType(ParticleTable::getAntiKaonType(G4int(-rdm)));
177  }
178  else if(std::abs(iso-iso_system) == 0){ // particle1 3/4 proton, 1/4 neutron; particle2 1/4 proton, 3/4 neutron
179  rdm = G4int(Random::shoot()*2.)*2-1;
180  G4double rdm2 = G4int(Random::shoot()*2.)*2-1;
184  antikaon->setType(ParticleTable::getAntiKaonType(-G4int(rdm2)));
185  }
186  else INCL_ERROR("Isospin non-conservation in NNToMissingStrangenessChannel" << '\n');
187  list.push_back(antikaon);
188  nbr_pions += 1; // need for addCreatedParticle loop
189  }
190  else{// N-K-Sigma
191 // assert(available_iso == 4);
192  if(std::abs(iso-iso_system) == 4){
193  particle1->setType(ParticleTable::getNucleonType((iso-iso_system)/4));
194  particle2->setType(ParticleTable::getSigmaType((iso-iso_system)/2));
195  kaon->setType(ParticleTable::getKaonType((iso-iso_system)/4));
196  }
197  else if(std::abs(iso-iso_system) == 2){ // choice: sigma quasi-free, kaon type semi-free, nucleon type fixed
198  rdm = Random::shoot();
199  // pp(pn) sp:0.42(0.23) sz:0.51(0.54) sm:0.07(0.23)
200  if(((iso == 0) && (rdm*0.77 < 0.23)) || ((rdm*0.58 < 0.07) && (Math::sign(iso-iso_system)*2-iso != 0)) || ((rdm*0.93 < 0.42) && (Math::sign(iso-iso_system)*2+iso != 0) && (iso != 0))){
202  rdm = G4int(Random::shoot()*2.)*2-1;
205  }
206  else{
208  particle1->setType(ParticleTable::getNucleonType((iso-iso_system)/2));
209  kaon->setType(ParticleTable::getKaonType((iso-iso_system)/2));
210  }
211  }
212  else if(std::abs(iso-iso_system) == 0){ // choice: sigma free, kaontype semi-free, nucleon type fixed
213  if(((iso == 0) && rdm < 0.23) || ((iso == 2) && rdm < 0.42) || rdm < 0.07){
216  kaon->setType(KZero);
217  }
218  else if(((iso == 0) && rdm < 0.77) || ((iso == 2) && rdm < 0.93) || rdm < 0.58){
220  rdm = G4int(Random::shoot()*2.)*2-1;
223  }
224  else{
227  kaon->setType(KPlus);
228  }
229  }
230  else INCL_ERROR("Isospin non-conservation in NNToMissingStrangenessChannel" << '\n');
231  }
232 
233  list.push_back(particle1);
234  list.push_back(particle2);
235  list.push_back(kaon);
236 
237  PhaseSpaceGenerator::generateBiased(sqrtS, list, list.size()-3, angularSlope);
238 
239  // INCL_DEBUG("NN Missing " << (kaon->getMomentum().theta()) * 180. / G4INCL::Math::pi << '\n');
240 
243  fs->addCreatedParticle(kaon);
244  for(Int_t i=0; i<nbr_pions; i++) fs->addCreatedParticle(list[i]);
245 
246  }
247 }
void generateBiased(const G4double sqrtS, ParticleList &particles, const size_t index, const G4double slope)
Generate a biased event in the CM system.
G4bool pion(G4int ityp)
T max(const T t1, const T t2)
brief Return the largest of the two arguments
void addModifiedParticle(Particle *p)
const G4INCL::ThreeVector & getPosition() const
ParticleType getNucleonType(const G4int isosp)
Get the type of nucleon.
G4double totalEnergyInCM(Particle const *const p1, Particle const *const p2)
ParticleType getAntiKaonType(const G4int isosp)
Get the type of antikaon.
G4bool isLambda() const
Is this a Lambda?
double G4double
Definition: G4Types.hh:76
G4int Int_t
G4double shoot()
Definition: G4INCLRandom.cc:93
#define INCL_ERROR(x)
G4double getINCLMass(const G4int A, const G4int Z)
Get INCL nuclear mass (in MeV/c^2)
ParticleType getKaonType(const G4int isosp)
Get the type of kaon.
void setType(ParticleType t)
G4int sign(const T t)
int G4int
Definition: G4Types.hh:78
G4int getIsospin(const ParticleType t)
Get the isospin of a particle.
G4INCL::ParticleType getType() const
ParticleType getSigmaType(const G4int isosp)
Get the type of sigma.
G4double momentumInLab(Particle const *const p1, Particle const *const p2)
gives the momentum in the lab frame of two particles.
void addCreatedParticle(Particle *p)
G4double gauss(G4double sigma=1.)
ParticleType getPionType(const G4int isosp)
Get the type of pion.
T min(const T t1, const T t2)
brief Return the smallest of the two arguments