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