| 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 /* 38 /*
39 * G4INCLParticle.hh 39 * G4INCLParticle.hh
40 * 40 *
41 * \date Jun 5, 2009 41 * \date Jun 5, 2009
42 * \author Pekka Kaitaniemi 42 * \author Pekka Kaitaniemi
43 */ 43 */
44 44
45 #ifndef PARTICLE_HH_ 45 #ifndef PARTICLE_HH_
46 #define PARTICLE_HH_ 46 #define PARTICLE_HH_
47 47
48 #include "G4INCLThreeVector.hh" 48 #include "G4INCLThreeVector.hh"
49 #include "G4INCLParticleTable.hh" 49 #include "G4INCLParticleTable.hh"
50 #include "G4INCLParticleType.hh" 50 #include "G4INCLParticleType.hh"
51 #include "G4INCLParticleSpecies.hh" 51 #include "G4INCLParticleSpecies.hh"
52 #include "G4INCLLogger.hh" 52 #include "G4INCLLogger.hh"
53 #include "G4INCLUnorderedVector.hh" 53 #include "G4INCLUnorderedVector.hh"
54 #include "G4INCLAllocationPool.hh" 54 #include "G4INCLAllocationPool.hh"
55 #include <sstream> 55 #include <sstream>
56 #include <string> 56 #include <string>
57 57
58 namespace G4INCL { 58 namespace G4INCL {
59 59
60 class Particle; 60 class Particle;
61 61
62 class ParticleList : public UnorderedVector< 62 class ParticleList : public UnorderedVector<Particle*> {
63 public: 63 public:
64 void rotatePositionAndMomentum(const G4d 64 void rotatePositionAndMomentum(const G4double angle, const ThreeVector &axis) const;
65 void rotatePosition(const G4double angle 65 void rotatePosition(const G4double angle, const ThreeVector &axis) const;
66 void rotateMomentum(const G4double angle 66 void rotateMomentum(const G4double angle, const ThreeVector &axis) const;
67 void boost(const ThreeVector &b) const; 67 void boost(const ThreeVector &b) const;
68 G4double getParticleListBias() const; 68 G4double getParticleListBias() const;
69 std::vector<G4int> getParticleListBiasVe 69 std::vector<G4int> getParticleListBiasVector() const;
70 }; 70 };
71 71
72 typedef ParticleList::const_iterator Particl 72 typedef ParticleList::const_iterator ParticleIter;
73 typedef ParticleList::iterator Particl 73 typedef ParticleList::iterator ParticleMutableIter;
74 74
75 class Particle { 75 class Particle {
76 public: 76 public:
77 Particle(); 77 Particle();
78 Particle(ParticleType t, G4double energy, 78 Particle(ParticleType t, G4double energy, ThreeVector const &momentum, ThreeVector const &position);
79 Particle(ParticleType t, ThreeVector const 79 Particle(ParticleType t, ThreeVector const &momentum, ThreeVector const &position);
80 virtual ~Particle() {} 80 virtual ~Particle() {}
81 81
82 /** \brief Copy constructor 82 /** \brief Copy constructor
83 * 83 *
84 * Does not copy the particle ID. 84 * Does not copy the particle ID.
85 */ 85 */
86 Particle(const Particle &rhs) : 86 Particle(const Particle &rhs) :
87 theZ(rhs.theZ), 87 theZ(rhs.theZ),
88 theA(rhs.theA), 88 theA(rhs.theA),
89 theS(rhs.theS), 89 theS(rhs.theS),
90 theParticipantType(rhs.theParticipantTyp 90 theParticipantType(rhs.theParticipantType),
91 theType(rhs.theType), 91 theType(rhs.theType),
92 theEnergy(rhs.theEnergy), 92 theEnergy(rhs.theEnergy),
93 theFrozenEnergy(rhs.theFrozenEnergy), 93 theFrozenEnergy(rhs.theFrozenEnergy),
94 theMomentum(rhs.theMomentum), 94 theMomentum(rhs.theMomentum),
95 theFrozenMomentum(rhs.theFrozenMomentum) 95 theFrozenMomentum(rhs.theFrozenMomentum),
96 thePosition(rhs.thePosition), 96 thePosition(rhs.thePosition),
97 nCollisions(rhs.nCollisions), 97 nCollisions(rhs.nCollisions),
98 nDecays(rhs.nDecays), 98 nDecays(rhs.nDecays),
99 thePotentialEnergy(rhs.thePotentialEnerg 99 thePotentialEnergy(rhs.thePotentialEnergy),
100 rpCorrelated(rhs.rpCorrelated), 100 rpCorrelated(rhs.rpCorrelated),
101 uncorrelatedMomentum(rhs.uncorrelatedMom 101 uncorrelatedMomentum(rhs.uncorrelatedMomentum),
102 theParticleBias(rhs.theParticleBias), 102 theParticleBias(rhs.theParticleBias),
103 theNKaon(rhs.theNKaon), 103 theNKaon(rhs.theNKaon),
104 #ifdef INCLXX_IN_GEANT4_MODE <<
105 theParentResonancePDGCode(rhs.theParentR <<
106 theParentResonanceID(rhs.theParentResona <<
107 #endif <<
108 theHelicity(rhs.theHelicity), 104 theHelicity(rhs.theHelicity),
109 emissionTime(rhs.emissionTime), 105 emissionTime(rhs.emissionTime),
110 outOfWell(rhs.outOfWell), 106 outOfWell(rhs.outOfWell),
111 theMass(rhs.theMass) 107 theMass(rhs.theMass)
112 { 108 {
113 if(rhs.thePropagationEnergy == &(rhs.t 109 if(rhs.thePropagationEnergy == &(rhs.theFrozenEnergy))
114 thePropagationEnergy = &theFrozenEne 110 thePropagationEnergy = &theFrozenEnergy;
115 else 111 else
116 thePropagationEnergy = &theEnergy; 112 thePropagationEnergy = &theEnergy;
117 if(rhs.thePropagationMomentum == &(rhs 113 if(rhs.thePropagationMomentum == &(rhs.theFrozenMomentum))
118 thePropagationMomentum = &theFrozenM 114 thePropagationMomentum = &theFrozenMomentum;
119 else 115 else
120 thePropagationMomentum = &theMomentu 116 thePropagationMomentum = &theMomentum;
121 // ID intentionally not copied 117 // ID intentionally not copied
122 ID = nextID++; 118 ID = nextID++;
123 119
124 theBiasCollisionVector = rhs.theBiasCo 120 theBiasCollisionVector = rhs.theBiasCollisionVector;
125 } 121 }
126 122
127 protected: 123 protected:
128 /// \brief Helper method for the assignmen 124 /// \brief Helper method for the assignment operator
129 void swap(Particle &rhs) { 125 void swap(Particle &rhs) {
130 std::swap(theZ, rhs.theZ); 126 std::swap(theZ, rhs.theZ);
131 std::swap(theA, rhs.theA); 127 std::swap(theA, rhs.theA);
132 std::swap(theS, rhs.theS); 128 std::swap(theS, rhs.theS);
133 std::swap(theParticipantType, rhs.thePar 129 std::swap(theParticipantType, rhs.theParticipantType);
134 std::swap(theType, rhs.theType); 130 std::swap(theType, rhs.theType);
135 if(rhs.thePropagationEnergy == &(rhs.the 131 if(rhs.thePropagationEnergy == &(rhs.theFrozenEnergy))
136 thePropagationEnergy = &theFrozenEnerg 132 thePropagationEnergy = &theFrozenEnergy;
137 else 133 else
138 thePropagationEnergy = &theEnergy; 134 thePropagationEnergy = &theEnergy;
139 std::swap(theEnergy, rhs.theEnergy); 135 std::swap(theEnergy, rhs.theEnergy);
140 std::swap(theFrozenEnergy, rhs.theFrozen 136 std::swap(theFrozenEnergy, rhs.theFrozenEnergy);
141 if(rhs.thePropagationMomentum == &(rhs.t 137 if(rhs.thePropagationMomentum == &(rhs.theFrozenMomentum))
142 thePropagationMomentum = &theFrozenMom 138 thePropagationMomentum = &theFrozenMomentum;
143 else 139 else
144 thePropagationMomentum = &theMomentum; 140 thePropagationMomentum = &theMomentum;
145 std::swap(theMomentum, rhs.theMomentum); 141 std::swap(theMomentum, rhs.theMomentum);
146 std::swap(theFrozenMomentum, rhs.theFroz 142 std::swap(theFrozenMomentum, rhs.theFrozenMomentum);
147 std::swap(thePosition, rhs.thePosition); 143 std::swap(thePosition, rhs.thePosition);
148 std::swap(nCollisions, rhs.nCollisions); 144 std::swap(nCollisions, rhs.nCollisions);
149 std::swap(nDecays, rhs.nDecays); 145 std::swap(nDecays, rhs.nDecays);
150 std::swap(thePotentialEnergy, rhs.thePot 146 std::swap(thePotentialEnergy, rhs.thePotentialEnergy);
151 // ID intentionally not swapped 147 // ID intentionally not swapped
152 148
153 #ifdef INCLXX_IN_GEANT4_MODE <<
154 std::swap(theParentResonancePDGCode, rhs <<
155 std::swap(theParentResonanceID, rhs.theP <<
156 #endif <<
157 <<
158 std::swap(theHelicity, rhs.theHelicity); 149 std::swap(theHelicity, rhs.theHelicity);
159 std::swap(emissionTime, rhs.emissionTime 150 std::swap(emissionTime, rhs.emissionTime);
160 std::swap(outOfWell, rhs.outOfWell); 151 std::swap(outOfWell, rhs.outOfWell);
161 152
162 std::swap(theMass, rhs.theMass); 153 std::swap(theMass, rhs.theMass);
163 std::swap(rpCorrelated, rhs.rpCorrelated 154 std::swap(rpCorrelated, rhs.rpCorrelated);
164 std::swap(uncorrelatedMomentum, rhs.unco 155 std::swap(uncorrelatedMomentum, rhs.uncorrelatedMomentum);
165 156
166 std::swap(theParticleBias, rhs.thePartic 157 std::swap(theParticleBias, rhs.theParticleBias);
167 std::swap(theBiasCollisionVector, rhs.th 158 std::swap(theBiasCollisionVector, rhs.theBiasCollisionVector);
168 159
169 } 160 }
170 161
171 public: 162 public:
172 163
173 /** \brief Assignment operator 164 /** \brief Assignment operator
174 * 165 *
175 * Does not copy the particle ID. 166 * Does not copy the particle ID.
176 */ 167 */
177 Particle &operator=(const Particle &rhs) { 168 Particle &operator=(const Particle &rhs) {
178 Particle temporaryParticle(rhs); 169 Particle temporaryParticle(rhs);
179 swap(temporaryParticle); 170 swap(temporaryParticle);
180 return *this; 171 return *this;
181 } 172 }
182 173
183 /** 174 /**
184 * Get the particle type. 175 * Get the particle type.
185 * @see G4INCL::ParticleType 176 * @see G4INCL::ParticleType
186 */ 177 */
187 G4INCL::ParticleType getType() const { 178 G4INCL::ParticleType getType() const {
188 return theType; 179 return theType;
189 }; 180 };
190 181
191 /// \brief Get the particle species 182 /// \brief Get the particle species
192 virtual G4INCL::ParticleSpecies getSpecies 183 virtual G4INCL::ParticleSpecies getSpecies() const {
193 return ParticleSpecies(theType); 184 return ParticleSpecies(theType);
194 }; 185 };
195 186
196 void setType(ParticleType t) { 187 void setType(ParticleType t) {
197 theType = t; 188 theType = t;
198 switch(theType) 189 switch(theType)
199 { 190 {
200 case DeltaPlusPlus: 191 case DeltaPlusPlus:
201 theA = 1; 192 theA = 1;
202 theZ = 2; 193 theZ = 2;
203 theS = 0; 194 theS = 0;
204 break; 195 break;
205 case Proton: 196 case Proton:
206 case DeltaPlus: 197 case DeltaPlus:
207 theA = 1; 198 theA = 1;
208 theZ = 1; 199 theZ = 1;
209 theS = 0; 200 theS = 0;
210 break; 201 break;
211 case Neutron: 202 case Neutron:
212 case DeltaZero: 203 case DeltaZero:
213 theA = 1; 204 theA = 1;
214 theZ = 0; 205 theZ = 0;
215 theS = 0; 206 theS = 0;
216 break; 207 break;
217 case DeltaMinus: 208 case DeltaMinus:
218 theA = 1; 209 theA = 1;
219 theZ = -1; 210 theZ = -1;
220 theS = 0; 211 theS = 0;
221 break; 212 break;
222 case PiPlus: 213 case PiPlus:
223 theA = 0; 214 theA = 0;
224 theZ = 1; 215 theZ = 1;
225 theS = 0; 216 theS = 0;
226 break; 217 break;
227 case PiZero: 218 case PiZero:
228 case Eta: 219 case Eta:
229 case Omega: 220 case Omega:
230 case EtaPrime: 221 case EtaPrime:
231 case Photon: 222 case Photon:
232 theA = 0; 223 theA = 0;
233 theZ = 0; 224 theZ = 0;
234 theS = 0; 225 theS = 0;
235 break; 226 break;
236 case PiMinus: 227 case PiMinus:
237 theA = 0; 228 theA = 0;
238 theZ = -1; 229 theZ = -1;
239 theS = 0; 230 theS = 0;
240 break; 231 break;
241 case Lambda: 232 case Lambda:
242 theA = 1; 233 theA = 1;
243 theZ = 0; 234 theZ = 0;
244 theS = -1; 235 theS = -1;
245 break; 236 break;
246 case SigmaPlus: 237 case SigmaPlus:
247 theA = 1; 238 theA = 1;
248 theZ = 1; 239 theZ = 1;
249 theS = -1; 240 theS = -1;
250 break; 241 break;
251 case SigmaZero: 242 case SigmaZero:
252 theA = 1; 243 theA = 1;
253 theZ = 0; 244 theZ = 0;
254 theS = -1; 245 theS = -1;
255 break; 246 break;
256 case SigmaMinus: 247 case SigmaMinus:
257 theA = 1; 248 theA = 1;
258 theZ = -1; 249 theZ = -1;
259 theS = -1; 250 theS = -1;
260 break; <<
261 case antiProton: <<
262 theA = -1; <<
263 theZ = -1; <<
264 theS = 0; <<
265 break; <<
266 case XiMinus: <<
267 theA = 1; <<
268 theZ = -1; <<
269 theS = -2; <<
270 break; <<
271 case XiZero: <<
272 theA = 1; <<
273 theZ = 0; <<
274 theS = -2; <<
275 break; <<
276 case antiNeutron: <<
277 theA = -1; <<
278 theZ = 0; <<
279 theS = 0; <<
280 break; 251 break;
281 case antiLambda: <<
282 theA = -1; <<
283 theZ = 0; <<
284 theS = 1; <<
285 break; <<
286 case antiSigmaMinus: <<
287 theA = -1; <<
288 theZ = 1; <<
289 theS = 1; <<
290 break; <<
291 case antiSigmaPlus: <<
292 theA = -1; <<
293 theZ = -1; <<
294 theS = 1; <<
295 break; <<
296 case antiSigmaZero: <<
297 theA = -1; <<
298 theZ = 0; <<
299 theS = 1; <<
300 break; <<
301 case antiXiMinus: <<
302 theA = -1; <<
303 theZ = 1; <<
304 theS = 2; <<
305 break; <<
306 case antiXiZero: <<
307 theA = -1; <<
308 theZ = 0; <<
309 theS = 2; <<
310 break; <<
311 case KPlus: 252 case KPlus:
312 theA = 0; 253 theA = 0;
313 theZ = 1; 254 theZ = 1;
314 theS = 1; 255 theS = 1;
315 break; 256 break;
316 case KZero: 257 case KZero:
317 theA = 0; 258 theA = 0;
318 theZ = 0; 259 theZ = 0;
319 theS = 1; 260 theS = 1;
320 break; 261 break;
321 case KZeroBar: 262 case KZeroBar:
322 theA = 0; 263 theA = 0;
323 theZ = 0; 264 theZ = 0;
324 theS = -1; 265 theS = -1;
325 break; 266 break;
326 case KShort: 267 case KShort:
327 theA = 0; 268 theA = 0;
328 theZ = 0; 269 theZ = 0;
329 // theS should not be defined 270 // theS should not be defined
330 break; 271 break;
331 case KLong: 272 case KLong:
332 theA = 0; 273 theA = 0;
333 theZ = 0; 274 theZ = 0;
334 // theS should not be defined 275 // theS should not be defined
335 break; 276 break;
336 case KMinus: 277 case KMinus:
337 theA = 0; 278 theA = 0;
338 theZ = -1; 279 theZ = -1;
339 theS = -1; 280 theS = -1;
340 break; 281 break;
341 case Composite: 282 case Composite:
342 // INCL_ERROR("Trying to set particle 283 // INCL_ERROR("Trying to set particle type to Composite! Construct a Cluster object instead" << '\n');
343 theA = 0; 284 theA = 0;
344 theZ = 0; 285 theZ = 0;
345 theS = 0; 286 theS = 0;
346 break; << 287 break;
347 case UnknownParticle: 288 case UnknownParticle:
348 theA = 0; 289 theA = 0;
349 theZ = 0; 290 theZ = 0;
350 theS = 0; 291 theS = 0;
351 INCL_ERROR("Trying to set particle t 292 INCL_ERROR("Trying to set particle type to Unknown!" << '\n');
352 break; 293 break;
353 } 294 }
354 295
355 if( !isResonance() && t!=Composite ) 296 if( !isResonance() && t!=Composite )
356 setINCLMass(); 297 setINCLMass();
357 } 298 }
358 299
359 /** 300 /**
360 * Is this a nucleon? 301 * Is this a nucleon?
361 */ 302 */
362 G4bool isNucleon() const { 303 G4bool isNucleon() const {
363 if(theType == G4INCL::Proton || theType 304 if(theType == G4INCL::Proton || theType == G4INCL::Neutron)
364 return true; 305 return true;
365 else 306 else
366 return false; 307 return false;
367 }; 308 };
368 309
369 ParticipantType getParticipantType() const 310 ParticipantType getParticipantType() const {
370 return theParticipantType; 311 return theParticipantType;
371 } 312 }
372 313
373 void setParticipantType(ParticipantType co 314 void setParticipantType(ParticipantType const p) {
374 theParticipantType = p; 315 theParticipantType = p;
375 } 316 }
376 317
377 G4bool isParticipant() const { 318 G4bool isParticipant() const {
378 return (theParticipantType==Participant) 319 return (theParticipantType==Participant);
379 } 320 }
380 321
381 G4bool isTargetSpectator() const { 322 G4bool isTargetSpectator() const {
382 return (theParticipantType==TargetSpecta 323 return (theParticipantType==TargetSpectator);
383 } 324 }
384 325
385 G4bool isProjectileSpectator() const { 326 G4bool isProjectileSpectator() const {
386 return (theParticipantType==ProjectileSp 327 return (theParticipantType==ProjectileSpectator);
387 } 328 }
388 329
389 virtual void makeParticipant() { 330 virtual void makeParticipant() {
390 theParticipantType = Participant; 331 theParticipantType = Participant;
391 } 332 }
392 333
393 virtual void makeTargetSpectator() { 334 virtual void makeTargetSpectator() {
394 theParticipantType = TargetSpectator; 335 theParticipantType = TargetSpectator;
395 } 336 }
396 337
397 virtual void makeProjectileSpectator() { 338 virtual void makeProjectileSpectator() {
398 theParticipantType = ProjectileSpectator 339 theParticipantType = ProjectileSpectator;
399 } 340 }
400 341
401 /** \brief Is this a pion? */ 342 /** \brief Is this a pion? */
402 G4bool isPion() const { return (theType == 343 G4bool isPion() const { return (theType == PiPlus || theType == PiZero || theType == PiMinus); }
403 344
404 /** \brief Is this an eta? */ 345 /** \brief Is this an eta? */
405 G4bool isEta() const { return (theType == 346 G4bool isEta() const { return (theType == Eta); }
406 347
407 /** \brief Is this an omega? */ 348 /** \brief Is this an omega? */
408 G4bool isOmega() const { return (theType = 349 G4bool isOmega() const { return (theType == Omega); }
409 350
410 /** \brief Is this an etaprime? */ 351 /** \brief Is this an etaprime? */
411 G4bool isEtaPrime() const { return (theTyp 352 G4bool isEtaPrime() const { return (theType == EtaPrime); }
412 353
413 /** \brief Is this a photon? */ 354 /** \brief Is this a photon? */
414 G4bool isPhoton() const { return (theType 355 G4bool isPhoton() const { return (theType == Photon); }
415 356
416 /** \brief Is it a resonance? */ 357 /** \brief Is it a resonance? */
417 inline G4bool isResonance() const { return 358 inline G4bool isResonance() const { return isDelta(); }
418 359
419 /** \brief Is it a Delta? */ 360 /** \brief Is it a Delta? */
420 inline G4bool isDelta() const { 361 inline G4bool isDelta() const {
421 return (theType==DeltaPlusPlus || theTyp 362 return (theType==DeltaPlusPlus || theType==DeltaPlus ||
422 theType==DeltaZero || theType==Delta 363 theType==DeltaZero || theType==DeltaMinus); }
423 364
424 /** \brief Is this a Sigma? */ 365 /** \brief Is this a Sigma? */
425 G4bool isSigma() const { return (theType = << 366 G4bool isSigma() const { return (theType == SigmaPlus || theType == SigmaZero || theType == SigmaMinus); }
426 367
427 /** \brief Is this a Kaon? */ 368 /** \brief Is this a Kaon? */
428 G4bool isKaon() const { return (theType == << 369 G4bool isKaon() const { return (theType == KPlus || theType == KZero); }
429 370
430 /** \brief Is this an antiKaon? */ 371 /** \brief Is this an antiKaon? */
431 G4bool isAntiKaon() const { return (theTyp 372 G4bool isAntiKaon() const { return (theType == KZeroBar || theType == KMinus); }
432 373
433 /** \brief Is this a Lambda? */ 374 /** \brief Is this a Lambda? */
434 G4bool isLambda() const { return (theType 375 G4bool isLambda() const { return (theType == Lambda); }
435 376
436 /** \brief Is this a Nucleon or a Lambda? 377 /** \brief Is this a Nucleon or a Lambda? */
437 G4bool isNucleonorLambda() const { return 378 G4bool isNucleonorLambda() const { return (isNucleon() || isLambda()); }
438 379
439 /** \brief Is this an Hyperon? */ 380 /** \brief Is this an Hyperon? */
440 G4bool isHyperon() const { return (isLambd << 381 G4bool isHyperon() const { return (isLambda() || isSigma()); }
441 382
442 /** \brief Is this a Meson? */ 383 /** \brief Is this a Meson? */
443 G4bool isMeson() const { return (isPion() 384 G4bool isMeson() const { return (isPion() || isKaon() || isAntiKaon() || isEta() || isEtaPrime() || isOmega()); }
444 385
445 /** \brief Is this a Baryon? */ 386 /** \brief Is this a Baryon? */
446 G4bool isBaryon() const { return (isNucleo 387 G4bool isBaryon() const { return (isNucleon() || isResonance() || isHyperon()); }
447 388
448 /** \brief Is this a Strange? */ << 389 /** \brief Is this an Strange? */
449 G4bool isStrange() const { return (isKaon( 390 G4bool isStrange() const { return (isKaon() || isAntiKaon() || isHyperon()); }
450 <<
451 /** \brief Is this a Xi? */ <<
452 G4bool isXi() const { return (theType == X <<
453 <<
454 /** \brief Is this an antinucleon? */ <<
455 G4bool isAntiNucleon() const { return (the <<
456 <<
457 /** \brief Is this an antiSigma? */ <<
458 G4bool isAntiSigma() const { return (theTy <<
459 <<
460 /** \brief Is this an antiXi? */ <<
461 G4bool isAntiXi() const { return (theType <<
462 <<
463 /** \brief Is this an antiLambda? */ <<
464 G4bool isAntiLambda() const { return (theT <<
465 <<
466 /** \brief Is this an antiHyperon? */ <<
467 G4bool isAntiHyperon() const { return (isA <<
468 <<
469 /** \brief Is this an antiBaryon? */ <<
470 G4bool isAntiBaryon() const { return (isAn <<
471 <<
472 /** \brief Is this an antiNucleon or an an <<
473 G4bool isAntiNucleonorAntiLambda() const { <<
474 391
475 /** \brief Returns the baryon number. */ 392 /** \brief Returns the baryon number. */
476 G4int getA() const { return theA; } 393 G4int getA() const { return theA; }
477 394
478 /** \brief Returns the charge number. */ 395 /** \brief Returns the charge number. */
479 G4int getZ() const { return theZ; } 396 G4int getZ() const { return theZ; }
480 397
481 /** \brief Returns the strangeness number. 398 /** \brief Returns the strangeness number. */
482 G4int getS() const { return theS; } 399 G4int getS() const { return theS; }
483 400
484 G4double getBeta() const { 401 G4double getBeta() const {
485 const G4double P = theMomentum.mag(); 402 const G4double P = theMomentum.mag();
486 return P/theEnergy; 403 return P/theEnergy;
487 } 404 }
488 405
489 /** 406 /**
490 * Returns a three vector we can give to t 407 * Returns a three vector we can give to the boost() -method.
491 * 408 *
492 * In order to go to the particle rest fra 409 * In order to go to the particle rest frame you need to multiply
493 * the boost vector by -1.0. 410 * the boost vector by -1.0.
494 */ 411 */
495 ThreeVector boostVector() const { 412 ThreeVector boostVector() const {
496 return theMomentum / theEnergy; 413 return theMomentum / theEnergy;
497 } 414 }
498 415
499 /** 416 /**
500 * Boost the particle using a boost vector 417 * Boost the particle using a boost vector.
501 * 418 *
502 * Example (go to the particle rest frame) 419 * Example (go to the particle rest frame):
503 * particle->boost(particle->boostVector() 420 * particle->boost(particle->boostVector());
504 */ 421 */
505 void boost(const ThreeVector &aBoostVector 422 void boost(const ThreeVector &aBoostVector) {
506 const G4double beta2 = aBoostVector.mag2 423 const G4double beta2 = aBoostVector.mag2();
507 const G4double gamma = 1.0 / std::sqrt(1 424 const G4double gamma = 1.0 / std::sqrt(1.0 - beta2);
508 const G4double bp = theMomentum.dot(aBoo 425 const G4double bp = theMomentum.dot(aBoostVector);
509 const G4double alpha = (gamma*gamma)/(1. 426 const G4double alpha = (gamma*gamma)/(1.0 + gamma);
510 427
511 theMomentum = theMomentum + aBoostVector 428 theMomentum = theMomentum + aBoostVector * (alpha * bp - gamma * theEnergy);
512 theEnergy = gamma * (theEnergy - bp); 429 theEnergy = gamma * (theEnergy - bp);
513 } 430 }
514 431
515 /** \brief Lorentz-contract the particle p 432 /** \brief Lorentz-contract the particle position around some center
516 * 433 *
517 * Apply Lorentz contraction to the positi 434 * Apply Lorentz contraction to the position component along the
518 * direction of the boost vector. 435 * direction of the boost vector.
519 * 436 *
520 * \param aBoostVector the boost vector (v 437 * \param aBoostVector the boost vector (velocity) [c]
521 * \param refPos the reference position 438 * \param refPos the reference position
522 */ 439 */
523 void lorentzContract(const ThreeVector &aB 440 void lorentzContract(const ThreeVector &aBoostVector, const ThreeVector &refPos) {
524 const G4double beta2 = aBoostVector.mag2 441 const G4double beta2 = aBoostVector.mag2();
525 const G4double gamma = 1.0 / std::sqrt(1 442 const G4double gamma = 1.0 / std::sqrt(1.0 - beta2);
526 const ThreeVector theRelativePosition = 443 const ThreeVector theRelativePosition = thePosition - refPos;
527 const ThreeVector transversePosition = t 444 const ThreeVector transversePosition = theRelativePosition - aBoostVector * (theRelativePosition.dot(aBoostVector) / aBoostVector.mag2());
528 const ThreeVector longitudinalPosition = 445 const ThreeVector longitudinalPosition = theRelativePosition - transversePosition;
529 446
530 thePosition = refPos + transversePositio 447 thePosition = refPos + transversePosition + longitudinalPosition / gamma;
531 } 448 }
532 449
533 /** \brief Get the cached particle mass. * 450 /** \brief Get the cached particle mass. */
534 inline G4double getMass() const { return t 451 inline G4double getMass() const { return theMass; }
535 452
536 /** \brief Get the INCL particle mass. */ 453 /** \brief Get the INCL particle mass. */
537 inline G4double getINCLMass() const { 454 inline G4double getINCLMass() const {
538 switch(theType) { 455 switch(theType) {
539 case Proton: 456 case Proton:
540 case Neutron: 457 case Neutron:
541 case PiPlus: 458 case PiPlus:
542 case PiMinus: 459 case PiMinus:
543 case PiZero: 460 case PiZero:
544 case Lambda: 461 case Lambda:
545 case SigmaPlus: 462 case SigmaPlus:
546 case SigmaZero: 463 case SigmaZero:
547 case SigmaMinus: << 464 case SigmaMinus:
548 case antiProton: <<
549 case XiZero: <<
550 case XiMinus: <<
551 case antiNeutron: <<
552 case antiLambda: <<
553 case antiSigmaPlus: <<
554 case antiSigmaZero: <<
555 case antiSigmaMinus: <<
556 case antiXiZero: <<
557 case antiXiMinus: <<
558 case KPlus: 465 case KPlus:
559 case KZero: 466 case KZero:
560 case KZeroBar: 467 case KZeroBar:
561 case KShort: 468 case KShort:
562 case KLong: 469 case KLong:
563 case KMinus: 470 case KMinus:
564 case Eta: 471 case Eta:
565 case Omega: 472 case Omega:
566 case EtaPrime: 473 case EtaPrime:
567 case Photon: << 474 case Photon:
568 return ParticleTable::getINCLMass(th 475 return ParticleTable::getINCLMass(theType);
569 break; 476 break;
570 477
571 case DeltaPlusPlus: 478 case DeltaPlusPlus:
572 case DeltaPlus: 479 case DeltaPlus:
573 case DeltaZero: 480 case DeltaZero:
574 case DeltaMinus: 481 case DeltaMinus:
575 return theMass; 482 return theMass;
576 break; 483 break;
577 484
578 case Composite: 485 case Composite:
579 return ParticleTable::getINCLMass(th 486 return ParticleTable::getINCLMass(theA,theZ,theS);
580 break; 487 break;
581 488
582 default: 489 default:
583 INCL_ERROR("Particle::getINCLMass: U 490 INCL_ERROR("Particle::getINCLMass: Unknown particle type." << '\n');
584 return 0.0; 491 return 0.0;
585 break; 492 break;
586 } 493 }
587 } 494 }
588 495
589 /** \brief Get the tabulated particle mass 496 /** \brief Get the tabulated particle mass. */
590 inline virtual G4double getTableMass() con 497 inline virtual G4double getTableMass() const {
591 switch(theType) { 498 switch(theType) {
592 case Proton: 499 case Proton:
593 case Neutron: 500 case Neutron:
594 case PiPlus: 501 case PiPlus:
595 case PiMinus: 502 case PiMinus:
596 case PiZero: 503 case PiZero:
597 case Lambda: 504 case Lambda:
598 case SigmaPlus: 505 case SigmaPlus:
599 case SigmaZero: 506 case SigmaZero:
600 case SigmaMinus: << 507 case SigmaMinus:
601 case antiProton: <<
602 case XiZero: <<
603 case XiMinus: <<
604 case antiNeutron: <<
605 case antiLambda: <<
606 case antiSigmaPlus: <<
607 case antiSigmaZero: <<
608 case antiSigmaMinus: <<
609 case antiXiZero: <<
610 case antiXiMinus: <<
611 case KPlus: 508 case KPlus:
612 case KZero: 509 case KZero:
613 case KZeroBar: 510 case KZeroBar:
614 case KShort: 511 case KShort:
615 case KLong: 512 case KLong:
616 case KMinus: 513 case KMinus:
617 case Eta: 514 case Eta:
618 case Omega: 515 case Omega:
619 case EtaPrime: 516 case EtaPrime:
620 case Photon: << 517 case Photon:
621 return ParticleTable::getTablePartic 518 return ParticleTable::getTableParticleMass(theType);
622 break; 519 break;
623 520
624 case DeltaPlusPlus: 521 case DeltaPlusPlus:
625 case DeltaPlus: 522 case DeltaPlus:
626 case DeltaZero: 523 case DeltaZero:
627 case DeltaMinus: 524 case DeltaMinus:
628 return theMass; 525 return theMass;
629 break; 526 break;
630 527
631 case Composite: 528 case Composite:
632 return ParticleTable::getTableMass(t 529 return ParticleTable::getTableMass(theA,theZ,theS);
633 break; 530 break;
634 531
635 default: 532 default:
636 INCL_ERROR("Particle::getTableMass: 533 INCL_ERROR("Particle::getTableMass: Unknown particle type." << '\n');
637 return 0.0; 534 return 0.0;
638 break; 535 break;
639 } 536 }
640 } 537 }
641 538
642 /** \brief Get the real particle mass. */ 539 /** \brief Get the real particle mass. */
643 inline G4double getRealMass() const { 540 inline G4double getRealMass() const {
644 switch(theType) { 541 switch(theType) {
645 case Proton: 542 case Proton:
646 case Neutron: 543 case Neutron:
647 case PiPlus: 544 case PiPlus:
648 case PiMinus: 545 case PiMinus:
649 case PiZero: 546 case PiZero:
650 case Lambda: 547 case Lambda:
651 case SigmaPlus: 548 case SigmaPlus:
652 case SigmaZero: 549 case SigmaZero:
653 case SigmaMinus: << 550 case SigmaMinus:
654 case antiProton: <<
655 case XiZero: <<
656 case XiMinus: <<
657 case antiNeutron: <<
658 case antiLambda: <<
659 case antiSigmaPlus: <<
660 case antiSigmaZero: <<
661 case antiSigmaMinus: <<
662 case antiXiZero: <<
663 case antiXiMinus: <<
664 case KPlus: 551 case KPlus:
665 case KZero: 552 case KZero:
666 case KZeroBar: 553 case KZeroBar:
667 case KShort: 554 case KShort:
668 case KLong: 555 case KLong:
669 case KMinus: 556 case KMinus:
670 case Eta: 557 case Eta:
671 case Omega: 558 case Omega:
672 case EtaPrime: 559 case EtaPrime:
673 case Photon: << 560 case Photon:
674 return ParticleTable::getRealMass(th 561 return ParticleTable::getRealMass(theType);
675 break; 562 break;
676 563
677 case DeltaPlusPlus: 564 case DeltaPlusPlus:
678 case DeltaPlus: 565 case DeltaPlus:
679 case DeltaZero: 566 case DeltaZero:
680 case DeltaMinus: 567 case DeltaMinus:
681 return theMass; 568 return theMass;
682 break; 569 break;
683 570
684 case Composite: 571 case Composite:
685 return ParticleTable::getRealMass(th 572 return ParticleTable::getRealMass(theA,theZ,theS);
686 break; 573 break;
687 574
688 default: 575 default:
689 INCL_ERROR("Particle::getRealMass: U 576 INCL_ERROR("Particle::getRealMass: Unknown particle type." << '\n');
690 return 0.0; 577 return 0.0;
691 break; 578 break;
692 } 579 }
693 } 580 }
694 581
695 /// \brief Set the mass of the Particle to 582 /// \brief Set the mass of the Particle to its real mass
696 void setRealMass() { setMass(getRealMass() 583 void setRealMass() { setMass(getRealMass()); }
697 584
698 /// \brief Set the mass of the Particle to 585 /// \brief Set the mass of the Particle to its table mass
699 void setTableMass() { setMass(getTableMass 586 void setTableMass() { setMass(getTableMass()); }
700 587
701 /// \brief Set the mass of the Particle to 588 /// \brief Set the mass of the Particle to its table mass
702 void setINCLMass() { setMass(getINCLMass() 589 void setINCLMass() { setMass(getINCLMass()); }
703 590
704 /**\brief Computes correction on the emiss 591 /**\brief Computes correction on the emission Q-value
705 * 592 *
706 * Computes the correction that must be ap 593 * Computes the correction that must be applied to INCL particles in
707 * order to obtain the correct Q-value for 594 * order to obtain the correct Q-value for particle emission from a given
708 * nucleus. For absorption, the correction 595 * nucleus. For absorption, the correction is obviously equal to minus
709 * the value returned by this function. 596 * the value returned by this function.
710 * 597 *
711 * \param AParent the mass number of the e 598 * \param AParent the mass number of the emitting nucleus
712 * \param ZParent the charge number of the 599 * \param ZParent the charge number of the emitting nucleus
713 * \return the correction 600 * \return the correction
714 */ 601 */
715 G4double getEmissionQValueCorrection(const 602 G4double getEmissionQValueCorrection(const G4int AParent, const G4int ZParent) const {
716 const G4int SParent = 0; 603 const G4int SParent = 0;
717 const G4int ADaughter = AParent - theA; 604 const G4int ADaughter = AParent - theA;
718 const G4int ZDaughter = ZParent - theZ; 605 const G4int ZDaughter = ZParent - theZ;
719 const G4int SDaughter = 0; 606 const G4int SDaughter = 0;
720 607
721 // Note the minus sign here 608 // Note the minus sign here
722 G4double theQValue; 609 G4double theQValue;
723 if(isCluster()) 610 if(isCluster())
724 theQValue = -ParticleTable::getTableQV 611 theQValue = -ParticleTable::getTableQValue(theA, theZ, theS, ADaughter, ZDaughter, SDaughter);
725 else { 612 else {
726 const G4double massTableParent = Parti 613 const G4double massTableParent = ParticleTable::getTableMass(AParent,ZParent,SParent);
727 const G4double massTableDaughter = Par 614 const G4double massTableDaughter = ParticleTable::getTableMass(ADaughter,ZDaughter,SDaughter);
728 const G4double massTableParticle = get 615 const G4double massTableParticle = getTableMass();
729 theQValue = massTableParent - massTabl 616 theQValue = massTableParent - massTableDaughter - massTableParticle;
730 } 617 }
731 618
732 const G4double massINCLParent = Particle 619 const G4double massINCLParent = ParticleTable::getINCLMass(AParent,ZParent,SParent);
733 const G4double massINCLDaughter = Partic 620 const G4double massINCLDaughter = ParticleTable::getINCLMass(ADaughter,ZDaughter,SDaughter);
734 const G4double massINCLParticle = getINC 621 const G4double massINCLParticle = getINCLMass();
735 622
736 // The rhs corresponds to the INCL Q-val 623 // The rhs corresponds to the INCL Q-value
737 return theQValue - (massINCLParent-massI 624 return theQValue - (massINCLParent-massINCLDaughter-massINCLParticle);
738 } 625 }
739 626
740 G4double getEmissionPbarQvalueCorrection(c <<
741 G4int SParent = 0; <<
742 G4int SDaughter = 0; <<
743 G4int ADaughter = AParent - 1; <<
744 G4int ZDaughter; <<
745 G4bool isProton = Victim; <<
746 if(isProton){ //proton is annihilate <<
747 ZDaughter = ZParent - 1; <<
748 } <<
749 else { //neutron is annihilated <<
750 ZDaughter = ZParent; <<
751 } <<
752 <<
753 G4double theQValue; //same procedure as <<
754 <<
755 const G4double massTableParent = Particl <<
756 const G4double massTableDaughter = Parti <<
757 const G4double massTableParticle = getTa <<
758 theQValue = massTableParent - massTableD <<
759 <<
760 const G4double massINCLParent = Particle <<
761 const G4double massINCLDaughter = Partic <<
762 const G4double massINCLParticle = getINC <<
763 <<
764 return theQValue - (massINCLParent-massI <<
765 } <<
766 <<
767 /**\brief Computes correction on the trans 627 /**\brief Computes correction on the transfer Q-value
768 * 628 *
769 * Computes the correction that must be ap 629 * Computes the correction that must be applied to INCL particles in
770 * order to obtain the correct Q-value for 630 * order to obtain the correct Q-value for particle transfer from a given
771 * nucleus to another. 631 * nucleus to another.
772 * 632 *
773 * Assumes that the receving nucleus is IN 633 * Assumes that the receving nucleus is INCL's target nucleus, with the
774 * INCL separation energy. 634 * INCL separation energy.
775 * 635 *
776 * \param AFrom the mass number of the don 636 * \param AFrom the mass number of the donating nucleus
777 * \param ZFrom the charge number of the d 637 * \param ZFrom the charge number of the donating nucleus
778 * \param ATo the mass number of the recei 638 * \param ATo the mass number of the receiving nucleus
779 * \param ZTo the charge number of the rec 639 * \param ZTo the charge number of the receiving nucleus
780 * \return the correction 640 * \return the correction
781 */ 641 */
782 G4double getTransferQValueCorrection(const 642 G4double getTransferQValueCorrection(const G4int AFrom, const G4int ZFrom, const G4int ATo, const G4int ZTo) const {
783 const G4int SFrom = 0; 643 const G4int SFrom = 0;
784 const G4int STo = 0; 644 const G4int STo = 0;
785 const G4int AFromDaughter = AFrom - theA 645 const G4int AFromDaughter = AFrom - theA;
786 const G4int ZFromDaughter = ZFrom - theZ 646 const G4int ZFromDaughter = ZFrom - theZ;
787 const G4int SFromDaughter = 0; 647 const G4int SFromDaughter = 0;
788 const G4int AToDaughter = ATo + theA; 648 const G4int AToDaughter = ATo + theA;
789 const G4int ZToDaughter = ZTo + theZ; 649 const G4int ZToDaughter = ZTo + theZ;
790 const G4int SToDaughter = 0; 650 const G4int SToDaughter = 0;
791 const G4double theQValue = ParticleTable 651 const G4double theQValue = ParticleTable::getTableQValue(AToDaughter,ZToDaughter,SToDaughter,AFromDaughter,ZFromDaughter,SFromDaughter,AFrom,ZFrom,SFrom);
792 652
793 const G4double massINCLTo = ParticleTabl 653 const G4double massINCLTo = ParticleTable::getINCLMass(ATo,ZTo,STo);
794 const G4double massINCLToDaughter = Part 654 const G4double massINCLToDaughter = ParticleTable::getINCLMass(AToDaughter,ZToDaughter,SToDaughter);
795 /* Note that here we have to use the tab 655 /* Note that here we have to use the table mass in the INCL Q-value. We
796 * cannot use theMass, because at this s 656 * cannot use theMass, because at this stage the particle is probably
797 * still off-shell; and we cannot use ge 657 * still off-shell; and we cannot use getINCLMass(), because it leads to
798 * violations of global energy conservat 658 * violations of global energy conservation.
799 */ 659 */
800 const G4double massINCLParticle = getTab 660 const G4double massINCLParticle = getTableMass();
801 661
802 // The rhs corresponds to the INCL Q-val 662 // The rhs corresponds to the INCL Q-value for particle absorption
803 return theQValue - (massINCLToDaughter-m 663 return theQValue - (massINCLToDaughter-massINCLTo-massINCLParticle);
804 } 664 }
805 665
806 /**\brief Computes correction on the emiss 666 /**\brief Computes correction on the emission Q-value for hypernuclei
807 * 667 *
808 * Computes the correction that must be ap 668 * Computes the correction that must be applied to INCL particles in
809 * order to obtain the correct Q-value for 669 * order to obtain the correct Q-value for particle emission from a given
810 * nucleus. For absorption, the correction 670 * nucleus. For absorption, the correction is obviously equal to minus
811 * the value returned by this function. 671 * the value returned by this function.
812 * 672 *
813 * \param AParent the mass number of the e 673 * \param AParent the mass number of the emitting nucleus
814 * \param ZParent the charge number of the 674 * \param ZParent the charge number of the emitting nucleus
815 * \param SParent the strangess number of 675 * \param SParent the strangess number of the emitting nucleus
816 * \return the correction 676 * \return the correction
817 */ 677 */
818 G4double getEmissionQValueCorrection(const 678 G4double getEmissionQValueCorrection(const G4int AParent, const G4int ZParent, const G4int SParent) const {
819 const G4int ADaughter = AParent - theA; 679 const G4int ADaughter = AParent - theA;
820 const G4int ZDaughter = ZParent - theZ; 680 const G4int ZDaughter = ZParent - theZ;
821 const G4int SDaughter = SParent - theS; 681 const G4int SDaughter = SParent - theS;
822 682
823 // Note the minus sign here 683 // Note the minus sign here
824 G4double theQValue; 684 G4double theQValue;
825 if(isCluster()) 685 if(isCluster())
826 theQValue = -ParticleTable::getTableQV 686 theQValue = -ParticleTable::getTableQValue(theA, theZ, theS, ADaughter, ZDaughter, SDaughter);
827 else { 687 else {
828 const G4double massTableParent = Parti 688 const G4double massTableParent = ParticleTable::getTableMass(AParent,ZParent,SParent);
829 const G4double massTableDaughter = Par 689 const G4double massTableDaughter = ParticleTable::getTableMass(ADaughter,ZDaughter,SDaughter);
830 const G4double massTableParticle = get 690 const G4double massTableParticle = getTableMass();
831 theQValue = massTableParent - massTabl 691 theQValue = massTableParent - massTableDaughter - massTableParticle;
832 } 692 }
833 693
834 const G4double massINCLParent = Particle 694 const G4double massINCLParent = ParticleTable::getINCLMass(AParent,ZParent,SParent);
835 const G4double massINCLDaughter = Partic 695 const G4double massINCLDaughter = ParticleTable::getINCLMass(ADaughter,ZDaughter,SDaughter);
836 const G4double massINCLParticle = getINC 696 const G4double massINCLParticle = getINCLMass();
837 697
838 // The rhs corresponds to the INCL Q-val 698 // The rhs corresponds to the INCL Q-value
839 return theQValue - (massINCLParent-massI 699 return theQValue - (massINCLParent-massINCLDaughter-massINCLParticle);
840 } 700 }
841 701
842 /**\brief Computes correction on the trans 702 /**\brief Computes correction on the transfer Q-value for hypernuclei
843 * 703 *
844 * Computes the correction that must be ap 704 * Computes the correction that must be applied to INCL particles in
845 * order to obtain the correct Q-value for 705 * order to obtain the correct Q-value for particle transfer from a given
846 * nucleus to another. 706 * nucleus to another.
847 * 707 *
848 * Assumes that the receving nucleus is IN 708 * Assumes that the receving nucleus is INCL's target nucleus, with the
849 * INCL separation energy. 709 * INCL separation energy.
850 * 710 *
851 * \param AFrom the mass number of the don 711 * \param AFrom the mass number of the donating nucleus
852 * \param ZFrom the charge number of the d 712 * \param ZFrom the charge number of the donating nucleus
853 * \param SFrom the strangess number of th 713 * \param SFrom the strangess number of the donating nucleus
854 * \param ATo the mass number of the recei 714 * \param ATo the mass number of the receiving nucleus
855 * \param ZTo the charge number of the rec 715 * \param ZTo the charge number of the receiving nucleus
856 * \param STo the strangess number of the 716 * \param STo the strangess number of the receiving nucleus
857 * \return the correction 717 * \return the correction
858 */ 718 */
859 G4double getTransferQValueCorrection(const 719 G4double getTransferQValueCorrection(const G4int AFrom, const G4int ZFrom, const G4int SFrom, const G4int ATo, const G4int ZTo , const G4int STo) const {
860 const G4int AFromDaughter = AFrom - theA 720 const G4int AFromDaughter = AFrom - theA;
861 const G4int ZFromDaughter = ZFrom - theZ 721 const G4int ZFromDaughter = ZFrom - theZ;
862 const G4int SFromDaughter = SFrom - theS 722 const G4int SFromDaughter = SFrom - theS;
863 const G4int AToDaughter = ATo + theA; 723 const G4int AToDaughter = ATo + theA;
864 const G4int ZToDaughter = ZTo + theZ; 724 const G4int ZToDaughter = ZTo + theZ;
865 const G4int SToDaughter = STo + theS; 725 const G4int SToDaughter = STo + theS;
866 const G4double theQValue = ParticleTable 726 const G4double theQValue = ParticleTable::getTableQValue(AToDaughter,ZToDaughter,SFromDaughter,AFromDaughter,ZFromDaughter,SToDaughter,AFrom,ZFrom,SFrom);
867 727
868 const G4double massINCLTo = ParticleTabl 728 const G4double massINCLTo = ParticleTable::getINCLMass(ATo,ZTo,STo);
869 const G4double massINCLToDaughter = Part 729 const G4double massINCLToDaughter = ParticleTable::getINCLMass(AToDaughter,ZToDaughter,SToDaughter);
870 /* Note that here we have to use the tab 730 /* Note that here we have to use the table mass in the INCL Q-value. We
871 * cannot use theMass, because at this s 731 * cannot use theMass, because at this stage the particle is probably
872 * still off-shell; and we cannot use ge 732 * still off-shell; and we cannot use getINCLMass(), because it leads to
873 * violations of global energy conservat 733 * violations of global energy conservation.
874 */ 734 */
875 const G4double massINCLParticle = getTab 735 const G4double massINCLParticle = getTableMass();
876 736
877 // The rhs corresponds to the INCL Q-val 737 // The rhs corresponds to the INCL Q-value for particle absorption
878 return theQValue - (massINCLToDaughter-m 738 return theQValue - (massINCLToDaughter-massINCLTo-massINCLParticle);
879 } 739 }
880 740
881 741
882 742
883 /** \brief Get the the particle invariant 743 /** \brief Get the the particle invariant mass.
884 * 744 *
885 * Uses the relativistic invariant 745 * Uses the relativistic invariant
886 * \f[ m = \sqrt{E^2 - {\vec p}^2}\f] 746 * \f[ m = \sqrt{E^2 - {\vec p}^2}\f]
887 **/ 747 **/
888 G4double getInvariantMass() const { 748 G4double getInvariantMass() const {
889 const G4double mass = std::pow(theEnergy 749 const G4double mass = std::pow(theEnergy, 2) - theMomentum.dot(theMomentum);
890 if(mass < 0.0) { 750 if(mass < 0.0) {
891 INCL_ERROR("E*E - p*p is negative." << 751 INCL_ERROR("E*E - p*p is negative." << '\n');
892 return 0.0; 752 return 0.0;
893 } else { 753 } else {
894 return std::sqrt(mass); 754 return std::sqrt(mass);
895 } 755 }
896 }; 756 };
897 757
898 /// \brief Get the particle kinetic energy 758 /// \brief Get the particle kinetic energy.
899 inline G4double getKineticEnergy() const { 759 inline G4double getKineticEnergy() const { return theEnergy - theMass; }
900 760
901 /// \brief Get the particle potential ener 761 /// \brief Get the particle potential energy.
902 inline G4double getPotentialEnergy() const 762 inline G4double getPotentialEnergy() const { return thePotentialEnergy; }
903 763
904 /// \brief Set the particle potential ener 764 /// \brief Set the particle potential energy.
905 inline void setPotentialEnergy(G4double v) 765 inline void setPotentialEnergy(G4double v) { thePotentialEnergy = v; }
906 766
907 /** 767 /**
908 * Get the energy of the particle in MeV. 768 * Get the energy of the particle in MeV.
909 */ 769 */
910 G4double getEnergy() const 770 G4double getEnergy() const
911 { 771 {
912 return theEnergy; 772 return theEnergy;
913 }; 773 };
914 774
915 /** 775 /**
916 * Set the mass of the particle in MeV/c^2 776 * Set the mass of the particle in MeV/c^2.
917 */ 777 */
918 void setMass(G4double mass) 778 void setMass(G4double mass)
919 { 779 {
920 this->theMass = mass; 780 this->theMass = mass;
921 } 781 }
922 782
923 /** 783 /**
924 * Set the energy of the particle in MeV. 784 * Set the energy of the particle in MeV.
925 */ 785 */
926 void setEnergy(G4double energy) 786 void setEnergy(G4double energy)
927 { 787 {
928 this->theEnergy = energy; 788 this->theEnergy = energy;
929 }; 789 };
930 790
931 /** 791 /**
932 * Get the momentum vector. 792 * Get the momentum vector.
933 */ 793 */
934 const G4INCL::ThreeVector &getMomentum() c 794 const G4INCL::ThreeVector &getMomentum() const
935 { 795 {
936 return theMomentum; 796 return theMomentum;
937 }; 797 };
938 798
939 /** Get the angular momentum w.r.t. the or 799 /** Get the angular momentum w.r.t. the origin */
940 virtual G4INCL::ThreeVector getAngularMome 800 virtual G4INCL::ThreeVector getAngularMomentum() const
941 { 801 {
942 return thePosition.vector(theMomentum); 802 return thePosition.vector(theMomentum);
943 }; 803 };
944 804
945 /** 805 /**
946 * Set the momentum vector. 806 * Set the momentum vector.
947 */ 807 */
948 virtual void setMomentum(const G4INCL::Thr 808 virtual void setMomentum(const G4INCL::ThreeVector &momentum)
949 { 809 {
950 this->theMomentum = momentum; 810 this->theMomentum = momentum;
951 }; 811 };
952 812
953 /** 813 /**
954 * Set the position vector. 814 * Set the position vector.
955 */ 815 */
956 const G4INCL::ThreeVector &getPosition() c 816 const G4INCL::ThreeVector &getPosition() const
957 { 817 {
958 return thePosition; 818 return thePosition;
959 }; 819 };
960 820
961 virtual void setPosition(const G4INCL::Thr 821 virtual void setPosition(const G4INCL::ThreeVector &position)
962 { 822 {
963 this->thePosition = position; 823 this->thePosition = position;
964 }; 824 };
965 825
966 G4double getHelicity() { return theHelicit 826 G4double getHelicity() { return theHelicity; };
967 void setHelicity(G4double h) { theHelicity 827 void setHelicity(G4double h) { theHelicity = h; };
968 828
969 void propagate(G4double step) { 829 void propagate(G4double step) {
970 thePosition += ((*thePropagationMomentum 830 thePosition += ((*thePropagationMomentum)*(step/(*thePropagationEnergy)));
971 }; 831 };
972 832
973 /** \brief Return the number of collisions 833 /** \brief Return the number of collisions undergone by the particle. **/
974 G4int getNumberOfCollisions() const { retu 834 G4int getNumberOfCollisions() const { return nCollisions; }
975 835
976 /** \brief Set the number of collisions un 836 /** \brief Set the number of collisions undergone by the particle. **/
977 void setNumberOfCollisions(G4int n) { nCol 837 void setNumberOfCollisions(G4int n) { nCollisions = n; }
978 838
979 /** \brief Increment the number of collisi 839 /** \brief Increment the number of collisions undergone by the particle. **/
980 void incrementNumberOfCollisions() { nColl 840 void incrementNumberOfCollisions() { nCollisions++; }
981 841
982 /** \brief Return the number of decays und 842 /** \brief Return the number of decays undergone by the particle. **/
983 G4int getNumberOfDecays() const { return n 843 G4int getNumberOfDecays() const { return nDecays; }
984 844
985 /** \brief Set the number of decays underg 845 /** \brief Set the number of decays undergone by the particle. **/
986 void setNumberOfDecays(G4int n) { nDecays 846 void setNumberOfDecays(G4int n) { nDecays = n; }
987 847
988 /** \brief Increment the number of decays 848 /** \brief Increment the number of decays undergone by the particle. **/
989 void incrementNumberOfDecays() { nDecays++ 849 void incrementNumberOfDecays() { nDecays++; }
990 850
991 /** \brief Mark the particle as out of its 851 /** \brief Mark the particle as out of its potential well
992 * 852 *
993 * This flag is used to control pions crea 853 * This flag is used to control pions created outside their potential well
994 * in delta decay. The pion potential chec 854 * in delta decay. The pion potential checks it and returns zero if it is
995 * true (necessary in order to correctly e 855 * true (necessary in order to correctly enforce energy conservation). The
996 * Nucleus::applyFinalState() method uses 856 * Nucleus::applyFinalState() method uses it to determine whether new
997 * avatars should be generated for the par 857 * avatars should be generated for the particle.
998 */ 858 */
999 void setOutOfWell() { outOfWell = true; } 859 void setOutOfWell() { outOfWell = true; }
1000 860
1001 /// \brief Check if the particle is out o 861 /// \brief Check if the particle is out of its potential well
1002 G4bool isOutOfWell() const { return outOf 862 G4bool isOutOfWell() const { return outOfWell; }
1003 863
1004 void setEmissionTime(G4double t) { emissi 864 void setEmissionTime(G4double t) { emissionTime = t; }
1005 G4double getEmissionTime() { return emiss 865 G4double getEmissionTime() { return emissionTime; };
1006 866
1007 /** \brief Transverse component of the po 867 /** \brief Transverse component of the position w.r.t. the momentum. */
1008 ThreeVector getTransversePosition() const 868 ThreeVector getTransversePosition() const {
1009 return thePosition - getLongitudinalPos 869 return thePosition - getLongitudinalPosition();
1010 } 870 }
1011 871
1012 /** \brief Longitudinal component of the 872 /** \brief Longitudinal component of the position w.r.t. the momentum. */
1013 ThreeVector getLongitudinalPosition() con 873 ThreeVector getLongitudinalPosition() const {
1014 return *thePropagationMomentum * (thePo 874 return *thePropagationMomentum * (thePosition.dot(*thePropagationMomentum)/thePropagationMomentum->mag2());
1015 } 875 }
1016 876
1017 /** \brief Rescale the momentum to match 877 /** \brief Rescale the momentum to match the total energy. */
1018 const ThreeVector &adjustMomentumFromEner 878 const ThreeVector &adjustMomentumFromEnergy();
1019 879
1020 /** \brief Recompute the energy to match 880 /** \brief Recompute the energy to match the momentum. */
1021 G4double adjustEnergyFromMomentum(); 881 G4double adjustEnergyFromMomentum();
1022 882
1023 G4bool isCluster() const { 883 G4bool isCluster() const {
1024 return (theType == Composite); 884 return (theType == Composite);
1025 } 885 }
1026 886
1027 /// \brief Set the frozen particle moment 887 /// \brief Set the frozen particle momentum
1028 void setFrozenMomentum(const ThreeVector 888 void setFrozenMomentum(const ThreeVector &momentum) { theFrozenMomentum = momentum; }
1029 889
1030 /// \brief Set the frozen particle moment 890 /// \brief Set the frozen particle momentum
1031 void setFrozenEnergy(const G4double energ 891 void setFrozenEnergy(const G4double energy) { theFrozenEnergy = energy; }
1032 892
1033 /// \brief Get the frozen particle moment 893 /// \brief Get the frozen particle momentum
1034 ThreeVector getFrozenMomentum() const { r 894 ThreeVector getFrozenMomentum() const { return theFrozenMomentum; }
1035 895
1036 /// \brief Get the frozen particle moment 896 /// \brief Get the frozen particle momentum
1037 G4double getFrozenEnergy() const { return 897 G4double getFrozenEnergy() const { return theFrozenEnergy; }
1038 898
1039 /// \brief Get the propagation velocity o 899 /// \brief Get the propagation velocity of the particle
1040 ThreeVector getPropagationVelocity() cons 900 ThreeVector getPropagationVelocity() const { return (*thePropagationMomentum)/(*thePropagationEnergy); }
1041 901
1042 /** \brief Freeze particle propagation 902 /** \brief Freeze particle propagation
1043 * 903 *
1044 * Make the particle use theFrozenMomentu 904 * Make the particle use theFrozenMomentum and theFrozenEnergy for
1045 * propagation. The normal state can be r 905 * propagation. The normal state can be restored by calling the
1046 * thawPropagation() method. 906 * thawPropagation() method.
1047 */ 907 */
1048 void freezePropagation() { 908 void freezePropagation() {
1049 thePropagationMomentum = &theFrozenMome 909 thePropagationMomentum = &theFrozenMomentum;
1050 thePropagationEnergy = &theFrozenEnergy 910 thePropagationEnergy = &theFrozenEnergy;
1051 } 911 }
1052 912
1053 /** \brief Unfreeze particle propagation 913 /** \brief Unfreeze particle propagation
1054 * 914 *
1055 * Make the particle use theMomentum and 915 * Make the particle use theMomentum and theEnergy for propagation. Call
1056 * this method to restore the normal prop 916 * this method to restore the normal propagation if the
1057 * freezePropagation() method has been ca 917 * freezePropagation() method has been called.
1058 */ 918 */
1059 void thawPropagation() { 919 void thawPropagation() {
1060 thePropagationMomentum = &theMomentum; 920 thePropagationMomentum = &theMomentum;
1061 thePropagationEnergy = &theEnergy; 921 thePropagationEnergy = &theEnergy;
1062 } 922 }
1063 923
1064 /** \brief Rotate the particle position a 924 /** \brief Rotate the particle position and momentum
1065 * 925 *
1066 * \param angle the rotation angle 926 * \param angle the rotation angle
1067 * \param axis a unit vector representing 927 * \param axis a unit vector representing the rotation axis
1068 */ 928 */
1069 virtual void rotatePositionAndMomentum(co 929 virtual void rotatePositionAndMomentum(const G4double angle, const ThreeVector &axis) {
1070 rotatePosition(angle, axis); 930 rotatePosition(angle, axis);
1071 rotateMomentum(angle, axis); 931 rotateMomentum(angle, axis);
1072 } 932 }
1073 933
1074 /** \brief Rotate the particle position 934 /** \brief Rotate the particle position
1075 * 935 *
1076 * \param angle the rotation angle 936 * \param angle the rotation angle
1077 * \param axis a unit vector representing 937 * \param axis a unit vector representing the rotation axis
1078 */ 938 */
1079 virtual void rotatePosition(const G4doubl 939 virtual void rotatePosition(const G4double angle, const ThreeVector &axis) {
1080 thePosition.rotate(angle, axis); 940 thePosition.rotate(angle, axis);
1081 } 941 }
1082 942
1083 /** \brief Rotate the particle momentum 943 /** \brief Rotate the particle momentum
1084 * 944 *
1085 * \param angle the rotation angle 945 * \param angle the rotation angle
1086 * \param axis a unit vector representing 946 * \param axis a unit vector representing the rotation axis
1087 */ 947 */
1088 virtual void rotateMomentum(const G4doubl 948 virtual void rotateMomentum(const G4double angle, const ThreeVector &axis) {
1089 theMomentum.rotate(angle, axis); 949 theMomentum.rotate(angle, axis);
1090 theFrozenMomentum.rotate(angle, axis); 950 theFrozenMomentum.rotate(angle, axis);
1091 } 951 }
1092 952
1093 std::string print() const { 953 std::string print() const {
1094 std::stringstream ss; 954 std::stringstream ss;
1095 ss << "Particle (ID = " << ID << ") typ 955 ss << "Particle (ID = " << ID << ") type = ";
1096 ss << ParticleTable::getName(theType); 956 ss << ParticleTable::getName(theType);
1097 ss << '\n' 957 ss << '\n'
1098 << " energy = " << theEnergy << '\n 958 << " energy = " << theEnergy << '\n'
1099 << " momentum = " 959 << " momentum = "
1100 << theMomentum.print() 960 << theMomentum.print()
1101 << '\n' 961 << '\n'
1102 << " position = " 962 << " position = "
1103 << thePosition.print() 963 << thePosition.print()
1104 << '\n'; 964 << '\n';
1105 return ss.str(); 965 return ss.str();
1106 }; 966 };
1107 967
1108 std::string dump() const { 968 std::string dump() const {
1109 std::stringstream ss; 969 std::stringstream ss;
1110 ss << "(particle " << ID << " "; 970 ss << "(particle " << ID << " ";
1111 ss << ParticleTable::getName(theType); 971 ss << ParticleTable::getName(theType);
1112 ss << '\n' 972 ss << '\n'
1113 << thePosition.dump() 973 << thePosition.dump()
1114 << '\n' 974 << '\n'
1115 << theMomentum.dump() 975 << theMomentum.dump()
1116 << '\n' 976 << '\n'
1117 << theEnergy << ")" << '\n'; 977 << theEnergy << ")" << '\n';
1118 return ss.str(); 978 return ss.str();
1119 }; 979 };
1120 980
1121 long getID() const { return ID; }; 981 long getID() const { return ID; };
1122 982
1123 /** 983 /**
1124 * Return a NULL pointer 984 * Return a NULL pointer
1125 */ 985 */
1126 ParticleList const *getParticles() const 986 ParticleList const *getParticles() const {
1127 INCL_WARN("Particle::getParticles() met 987 INCL_WARN("Particle::getParticles() method was called on a Particle object" << '\n');
1128 return 0; 988 return 0;
1129 } 989 }
1130 990
1131 /** \brief Return the reflection momentum 991 /** \brief Return the reflection momentum
1132 * 992 *
1133 * The reflection momentum is used by cal 993 * The reflection momentum is used by calls to getSurfaceRadius to compute
1134 * the radius of the sphere where the nuc 994 * the radius of the sphere where the nucleon moves. It is necessary to
1135 * introduce fuzzy r-p correlations. 995 * introduce fuzzy r-p correlations.
1136 */ 996 */
1137 G4double getReflectionMomentum() const { 997 G4double getReflectionMomentum() const {
1138 if(rpCorrelated) 998 if(rpCorrelated)
1139 return theMomentum.mag(); 999 return theMomentum.mag();
1140 else 1000 else
1141 return uncorrelatedMomentum; 1001 return uncorrelatedMomentum;
1142 } 1002 }
1143 1003
1144 /// \brief Set the uncorrelated momentum 1004 /// \brief Set the uncorrelated momentum
1145 void setUncorrelatedMomentum(const G4doub 1005 void setUncorrelatedMomentum(const G4double p) { uncorrelatedMomentum = p; }
1146 1006
1147 /// \brief Make the particle follow a str 1007 /// \brief Make the particle follow a strict r-p correlation
1148 void rpCorrelate() { rpCorrelated = true; 1008 void rpCorrelate() { rpCorrelated = true; }
1149 1009
1150 /// \brief Make the particle not follow a 1010 /// \brief Make the particle not follow a strict r-p correlation
1151 void rpDecorrelate() { rpCorrelated = fal 1011 void rpDecorrelate() { rpCorrelated = false; }
1152 1012
1153 /// \brief Get the cosine of the angle be 1013 /// \brief Get the cosine of the angle between position and momentum
1154 G4double getCosRPAngle() const { 1014 G4double getCosRPAngle() const {
1155 const G4double norm = thePosition.mag2( 1015 const G4double norm = thePosition.mag2()*thePropagationMomentum->mag2();
1156 if(norm>0.) 1016 if(norm>0.)
1157 return thePosition.dot(*thePropagatio 1017 return thePosition.dot(*thePropagationMomentum) / std::sqrt(norm);
1158 else 1018 else
1159 return 1.; 1019 return 1.;
1160 } 1020 }
1161 1021
1162 /// \brief General bias vector function 1022 /// \brief General bias vector function
1163 static G4double getTotalBias(); 1023 static G4double getTotalBias();
1164 static void setINCLBiasVector(std::vector 1024 static void setINCLBiasVector(std::vector<G4double> NewVector);
1165 static void FillINCLBiasVector(G4double n 1025 static void FillINCLBiasVector(G4double newBias);
1166 static G4double getBiasFromVector(std::ve 1026 static G4double getBiasFromVector(std::vector<G4int> VectorBias);
1167 1027
1168 static std::vector<G4int> MergeVectorBias 1028 static std::vector<G4int> MergeVectorBias(Particle const * const p1, Particle const * const p2);
1169 static std::vector<G4int> MergeVectorBias 1029 static std::vector<G4int> MergeVectorBias(std::vector<G4int> p1, Particle const * const p2);
1170 1030
1171 /// \brief Get the particle bias. 1031 /// \brief Get the particle bias.
1172 G4double getParticleBias() const { return 1032 G4double getParticleBias() const { return theParticleBias; };
1173 1033
1174 /// \brief Set the particle bias. 1034 /// \brief Set the particle bias.
1175 void setParticleBias(G4double ParticleBia 1035 void setParticleBias(G4double ParticleBias) { this->theParticleBias = ParticleBias; }
1176 1036
1177 /// \brief Get the vector list of biased 1037 /// \brief Get the vector list of biased vertices on the particle path.
1178 std::vector<G4int> getBiasCollisionVector 1038 std::vector<G4int> getBiasCollisionVector() const { return theBiasCollisionVector; }
1179 1039
1180 /// \brief Set the vector list of biased 1040 /// \brief Set the vector list of biased vertices on the particle path.
1181 void setBiasCollisionVector(std::vector<G 1041 void setBiasCollisionVector(std::vector<G4int> BiasCollisionVector) {
1182 this->theBiasCollisionVector = BiasCollis 1042 this->theBiasCollisionVector = BiasCollisionVector;
1183 this->setParticleBias(Particle::getBiasFr << 1043 this->setParticleBias(Particle::getBiasFromVector(BiasCollisionVector));
1184 } 1044 }
1185 1045
1186 /** \brief Number of Kaon inside de nucle 1046 /** \brief Number of Kaon inside de nucleus
1187 * 1047 *
1188 * Put in the Particle class in order to 1048 * Put in the Particle class in order to calculate the
1189 * "correct" mass of composit particle. 1049 * "correct" mass of composit particle.
1190 * 1050 *
1191 */ 1051 */
1192 1052
1193 G4int getNumberOfKaon() const { return th 1053 G4int getNumberOfKaon() const { return theNKaon; };
1194 void setNumberOfKaon(const G4int NK) { th 1054 void setNumberOfKaon(const G4int NK) { theNKaon = NK; }
1195 <<
1196 #ifdef INCLXX_IN_GEANT4_MODE <<
1197 G4int getParentResonancePDGCode() const { <<
1198 void setParentResonancePDGCode(const G4in <<
1199 G4int getParentResonanceID() const { retu <<
1200 void setParentResonanceID(const G4int par <<
1201 #endif <<
1202 1055
1203 public: 1056 public:
1204 /** \brief Time ordered vector of all bia 1057 /** \brief Time ordered vector of all bias applied
1205 * 1058 *
1206 * /!\ Caution /!\ 1059 * /!\ Caution /!\
1207 * methods Assotiated to G4VectorCache<T> 1060 * methods Assotiated to G4VectorCache<T> are:
1208 * Push_back(…), 1061 * Push_back(…),
1209 * operator[], 1062 * operator[],
1210 * Begin(), 1063 * Begin(),
1211 * End(), 1064 * End(),
1212 * Clear(), 1065 * Clear(),
1213 * Size() and 1066 * Size() and
1214 * Pop_back() 1067 * Pop_back()
1215 * 1068 *
1216 */ 1069 */
1217 #ifdef INCLXX_IN_GEANT4_MODE 1070 #ifdef INCLXX_IN_GEANT4_MODE
1218 static std::vector<G4double> INCLBiasVe 1071 static std::vector<G4double> INCLBiasVector;
1219 //static G4VectorCache<G4double> INCLBi 1072 //static G4VectorCache<G4double> INCLBiasVector;
1220 #else 1073 #else
1221 static G4ThreadLocal std::vector<G4doub 1074 static G4ThreadLocal std::vector<G4double> INCLBiasVector;
1222 //static G4VectorCache<G4double> INCLBi 1075 //static G4VectorCache<G4double> INCLBiasVector;
1223 #endif 1076 #endif
1224 static G4ThreadLocal G4int nextBiasedColl 1077 static G4ThreadLocal G4int nextBiasedCollisionID;
1225 1078
1226 protected: 1079 protected:
1227 G4int theZ, theA, theS; 1080 G4int theZ, theA, theS;
1228 ParticipantType theParticipantType; 1081 ParticipantType theParticipantType;
1229 G4INCL::ParticleType theType; 1082 G4INCL::ParticleType theType;
1230 G4double theEnergy; 1083 G4double theEnergy;
1231 G4double *thePropagationEnergy; 1084 G4double *thePropagationEnergy;
1232 G4double theFrozenEnergy; 1085 G4double theFrozenEnergy;
1233 G4INCL::ThreeVector theMomentum; 1086 G4INCL::ThreeVector theMomentum;
1234 G4INCL::ThreeVector *thePropagationMoment 1087 G4INCL::ThreeVector *thePropagationMomentum;
1235 G4INCL::ThreeVector theFrozenMomentum; 1088 G4INCL::ThreeVector theFrozenMomentum;
1236 G4INCL::ThreeVector thePosition; 1089 G4INCL::ThreeVector thePosition;
1237 G4int nCollisions; 1090 G4int nCollisions;
1238 G4int nDecays; 1091 G4int nDecays;
1239 G4double thePotentialEnergy; 1092 G4double thePotentialEnergy;
1240 long ID; 1093 long ID;
1241 1094
1242 G4bool rpCorrelated; 1095 G4bool rpCorrelated;
1243 G4double uncorrelatedMomentum; 1096 G4double uncorrelatedMomentum;
1244 1097
1245 G4double theParticleBias; 1098 G4double theParticleBias;
1246 /// \brief The number of Kaons inside the 1099 /// \brief The number of Kaons inside the nucleus (update during the cascade)
1247 G4int theNKaon; 1100 G4int theNKaon;
1248 <<
1249 #ifdef INCLXX_IN_GEANT4_MODE <<
1250 G4int theParentResonancePDGCode; <<
1251 G4int theParentResonanceID; <<
1252 #endif <<
1253 1101
1254 private: 1102 private:
1255 G4double theHelicity; 1103 G4double theHelicity;
1256 G4double emissionTime; 1104 G4double emissionTime;
1257 G4bool outOfWell; 1105 G4bool outOfWell;
1258 1106
1259 /// \brief Time ordered vector of all bia 1107 /// \brief Time ordered vector of all biased vertices on the particle path
1260 std::vector<G4int> theBiasCollisionVector 1108 std::vector<G4int> theBiasCollisionVector;
1261 1109
1262 G4double theMass; 1110 G4double theMass;
1263 static G4ThreadLocal long nextID; 1111 static G4ThreadLocal long nextID;
1264 1112
1265 INCL_DECLARE_ALLOCATION_POOL(Particle) 1113 INCL_DECLARE_ALLOCATION_POOL(Particle)
1266 }; 1114 };
1267 } 1115 }
1268 1116
1269 #endif /* PARTICLE_HH_ */ 1117 #endif /* PARTICLE_HH_ */
1270 1118