Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4INCLParticle.hh
이 파일의 문서화 페이지로 가기
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // INCL++ intra-nuclear cascade model
27 // Alain Boudard, CEA-Saclay, France
28 // Joseph Cugnon, University of Liege, Belgium
29 // Jean-Christophe David, CEA-Saclay, France
30 // Pekka Kaitaniemi, CEA-Saclay, France, and Helsinki Institute of Physics, Finland
31 // Sylvie Leray, CEA-Saclay, France
32 // Davide Mancusi, CEA-Saclay, France
33 //
34 #define INCLXX_IN_GEANT4_MODE 1
35 
36 #include "globals.hh"
37 
38 /*
39  * G4INCLParticle.hh
40  *
41  * \date Jun 5, 2009
42  * \author Pekka Kaitaniemi
43  */
44 
45 #ifndef PARTICLE_HH_
46 #define PARTICLE_HH_
47 
48 #include "G4INCLThreeVector.hh"
49 #include "G4INCLParticleTable.hh"
50 #include "G4INCLParticleType.hh"
51 #include "G4INCLParticleSpecies.hh"
52 #include "G4INCLLogger.hh"
53 #include "G4INCLUnorderedVector.hh"
54 #include "G4INCLAllocationPool.hh"
55 #include <sstream>
56 #include <string>
57 
58 namespace G4INCL {
59 
60  class Particle;
61 
62  class ParticleList : public UnorderedVector<Particle*> {
63  public:
64  void rotatePositionAndMomentum(const G4double angle, const ThreeVector &axis) const;
65  void rotatePosition(const G4double angle, const ThreeVector &axis) const;
66  void rotateMomentum(const G4double angle, const ThreeVector &axis) const;
67  void boost(const ThreeVector &b) const;
69  std::vector<G4int> getParticleListBiasVector() const;
70  };
71 
72  typedef ParticleList::const_iterator ParticleIter;
73  typedef ParticleList::iterator ParticleMutableIter;
74 
75  class Particle {
76  public:
77  Particle();
79  Particle(ParticleType t, ThreeVector const &momentum, ThreeVector const &position);
80  virtual ~Particle() {}
81 
86  Particle(const Particle &rhs) :
87  theZ(rhs.theZ),
88  theA(rhs.theA),
89  theS(rhs.theS),
91  theType(rhs.theType),
92  theEnergy(rhs.theEnergy),
98  nDecays(rhs.nDecays),
105  outOfWell(rhs.outOfWell),
106  theMass(rhs.theMass)
107  {
108  if(rhs.thePropagationEnergy == &(rhs.theFrozenEnergy))
110  else
112  if(rhs.thePropagationMomentum == &(rhs.theFrozenMomentum))
114  else
116  // ID intentionally not copied
117  ID = nextID++;
118  }
119 
120  protected:
122  void swap(Particle &rhs) {
123  std::swap(theZ, rhs.theZ);
124  std::swap(theA, rhs.theA);
125  std::swap(theS, rhs.theS);
126  std::swap(theParticipantType, rhs.theParticipantType);
127  std::swap(theType, rhs.theType);
128  if(rhs.thePropagationEnergy == &(rhs.theFrozenEnergy))
130  else
132  std::swap(theEnergy, rhs.theEnergy);
133  std::swap(theFrozenEnergy, rhs.theFrozenEnergy);
134  if(rhs.thePropagationMomentum == &(rhs.theFrozenMomentum))
136  else
138  std::swap(theMomentum, rhs.theMomentum);
139  std::swap(theFrozenMomentum, rhs.theFrozenMomentum);
140  std::swap(thePosition, rhs.thePosition);
141  std::swap(nCollisions, rhs.nCollisions);
142  std::swap(nDecays, rhs.nDecays);
143  std::swap(thePotentialEnergy, rhs.thePotentialEnergy);
144  // ID intentionally not swapped
145 
146  std::swap(theHelicity, rhs.theHelicity);
147  std::swap(emissionTime, rhs.emissionTime);
148  std::swap(outOfWell, rhs.outOfWell);
149 
150  std::swap(theMass, rhs.theMass);
151  std::swap(rpCorrelated, rhs.rpCorrelated);
153  }
154 
155  public:
156 
161  Particle &operator=(const Particle &rhs) {
162  Particle temporaryParticle(rhs);
163  swap(temporaryParticle);
164  return *this;
165  }
166 
172  return theType;
173  };
174 
177  return ParticleSpecies(theType);
178  };
179 
181  theType = t;
182  switch(theType)
183  {
184  case DeltaPlusPlus:
185  theA = 1;
186  theZ = 2;
187  theS = 0;
188  break;
189  case Proton:
190  case DeltaPlus:
191  theA = 1;
192  theZ = 1;
193  theS = 0;
194  break;
195  case Neutron:
196  case DeltaZero:
197  theA = 1;
198  theZ = 0;
199  theS = 0;
200  break;
201  case DeltaMinus:
202  theA = 1;
203  theZ = -1;
204  theS = 0;
205  break;
206  case PiPlus:
207  theA = 0;
208  theZ = 1;
209  theS = 0;
210  break;
211  case PiZero:
212  case Eta:
213  case Omega:
214  case EtaPrime:
215  case Photon:
216  theA = 0;
217  theZ = 0;
218  theS = 0;
219  break;
220  case PiMinus:
221  theA = 0;
222  theZ = -1;
223  theS = 0;
224  break;
225  case Lambda:
226  theA = 1;
227  theZ = 0;
228  theS = -1;
229  break;
230  case SigmaPlus:
231  theA = 1;
232  theZ = 1;
233  theS = -1;
234  break;
235  case SigmaZero:
236  theA = 1;
237  theZ = 0;
238  theS = -1;
239  break;
240  case SigmaMinus:
241  theA = 1;
242  theZ = -1;
243  theS = -1;
244  break;
245  case KPlus:
246  theA = 0;
247  theZ = 1;
248  theS = 1;
249  break;
250  case KZero:
251  theA = 0;
252  theZ = 0;
253  theS = 1;
254  break;
255  case KZeroBar:
256  theA = 0;
257  theZ = 0;
258  theS = -1;
259  break;
260  case KShort:
261  theA = 0;
262  theZ = 0;
263  theS = -99;
264  break;
265  case KLong:
266  theA = 0;
267  theZ = 0;
268  theS = 99;
269  break;
270  case KMinus:
271  theA = 0;
272  theZ = -1;
273  theS = -1;
274  break;
275  case Composite:
276  // INCL_ERROR("Trying to set particle type to Composite! Construct a Cluster object instead" << '\n');
277  theA = 0;
278  theZ = 0;
279  theS = 0;
280  break;
281  case UnknownParticle:
282  theA = 0;
283  theZ = 0;
284  theS = 0;
285  INCL_ERROR("Trying to set particle type to Unknown!" << '\n');
286  break;
287  }
288 
289  if( !isResonance() && t!=Composite )
290  setINCLMass();
291  }
292 
296  G4bool isNucleon() const {
298  return true;
299  else
300  return false;
301  };
302 
304  return theParticipantType;
305  }
306 
309  }
310 
313  }
314 
317  }
318 
321  }
322 
323  virtual void makeParticipant() {
325  }
326 
327  virtual void makeTargetSpectator() {
329  }
330 
331  virtual void makeProjectileSpectator() {
333  }
334 
336  G4bool isPion() const { return (theType == PiPlus || theType == PiZero || theType == PiMinus); }
337 
339  G4bool isEta() const { return (theType == Eta); }
340 
342  G4bool isOmega() const { return (theType == Omega); }
343 
345  G4bool isEtaPrime() const { return (theType == EtaPrime); }
346 
348  G4bool isPhoton() const { return (theType == Photon); }
349 
351  inline G4bool isResonance() const { return isDelta(); }
352 
354  inline G4bool isDelta() const {
355  return (theType==DeltaPlusPlus || theType==DeltaPlus ||
357 
359  G4bool isSigma() const { return (theType == SigmaPlus || theType == SigmaZero || theType == SigmaMinus); }
360 
362  G4bool isKaon() const { return (theType == KPlus || theType == KZero); }
363 
365  G4bool isAntiKaon() const { return (theType == KZeroBar || theType == KMinus); }
366 
368  G4bool isLambda() const { return (theType == Lambda); }
369 
371  G4bool isNucleonorLambda() const { return (isNucleon() || isLambda()); }
372 
374  G4bool isHyperon() const { return (isLambda() || isSigma()); }
375 
377  G4bool isMeson() const { return (isPion() || isKaon() || isAntiKaon() || isEta() || isEtaPrime() || isOmega()); }
378 
380  G4bool isBaryon() const { return (isNucleon() || isResonance() || isHyperon()); }
381 
383  G4bool isStrange() const { return (isKaon() || isAntiKaon() || isHyperon()); }
384 
386  G4int getA() const { return theA; }
387 
389  G4int getZ() const { return theZ; }
390 
392  G4int getS() const { return theS; }
393 
394  G4double getBeta() const {
395  const G4double P = theMomentum.mag();
396  return P/theEnergy;
397  }
398 
406  return theMomentum / theEnergy;
407  }
408 
415  void boost(const ThreeVector &aBoostVector) {
416  const G4double beta2 = aBoostVector.mag2();
417  const G4double gamma = 1.0 / std::sqrt(1.0 - beta2);
418  const G4double bp = theMomentum.dot(aBoostVector);
419  const G4double alpha = (gamma*gamma)/(1.0 + gamma);
420 
421  theMomentum = theMomentum + aBoostVector * (alpha * bp - gamma * theEnergy);
422  theEnergy = gamma * (theEnergy - bp);
423  }
424 
433  void lorentzContract(const ThreeVector &aBoostVector, const ThreeVector &refPos) {
434  const G4double beta2 = aBoostVector.mag2();
435  const G4double gamma = 1.0 / std::sqrt(1.0 - beta2);
436  const ThreeVector theRelativePosition = thePosition - refPos;
437  const ThreeVector transversePosition = theRelativePosition - aBoostVector * (theRelativePosition.dot(aBoostVector) / aBoostVector.mag2());
438  const ThreeVector longitudinalPosition = theRelativePosition - transversePosition;
439 
440  thePosition = refPos + transversePosition + longitudinalPosition / gamma;
441  }
442 
444  inline G4double getMass() const { return theMass; }
445 
447  inline G4double getINCLMass() const {
448  switch(theType) {
449  case Proton:
450  case Neutron:
451  case PiPlus:
452  case PiMinus:
453  case PiZero:
454  case Lambda:
455  case SigmaPlus:
456  case SigmaZero:
457  case SigmaMinus:
458  case KPlus:
459  case KZero:
460  case KZeroBar:
461  case KShort:
462  case KLong:
463  case KMinus:
464  case Eta:
465  case Omega:
466  case EtaPrime:
467  case Photon:
469  break;
470 
471  case DeltaPlusPlus:
472  case DeltaPlus:
473  case DeltaZero:
474  case DeltaMinus:
475  return theMass;
476  break;
477 
478  case Composite:
480  break;
481 
482  default:
483  INCL_ERROR("Particle::getINCLMass: Unknown particle type." << '\n');
484  return 0.0;
485  break;
486  }
487  }
488 
490  inline virtual G4double getTableMass() const {
491  switch(theType) {
492  case Proton:
493  case Neutron:
494  case PiPlus:
495  case PiMinus:
496  case PiZero:
497  case Lambda:
498  case SigmaPlus:
499  case SigmaZero:
500  case SigmaMinus:
501  case KPlus:
502  case KZero:
503  case KZeroBar:
504  case KShort:
505  case KLong:
506  case KMinus:
507  case Eta:
508  case Omega:
509  case EtaPrime:
510  case Photon:
512  break;
513 
514  case DeltaPlusPlus:
515  case DeltaPlus:
516  case DeltaZero:
517  case DeltaMinus:
518  return theMass;
519  break;
520 
521  case Composite:
523  break;
524 
525  default:
526  INCL_ERROR("Particle::getTableMass: Unknown particle type." << '\n');
527  return 0.0;
528  break;
529  }
530  }
531 
533  inline G4double getRealMass() const {
534  switch(theType) {
535  case Proton:
536  case Neutron:
537  case PiPlus:
538  case PiMinus:
539  case PiZero:
540  case Lambda:
541  case SigmaPlus:
542  case SigmaZero:
543  case SigmaMinus:
544  case KPlus:
545  case KZero:
546  case KZeroBar:
547  case KShort:
548  case KLong:
549  case KMinus:
550  case Eta:
551  case Omega:
552  case EtaPrime:
553  case Photon:
555  break;
556 
557  case DeltaPlusPlus:
558  case DeltaPlus:
559  case DeltaZero:
560  case DeltaMinus:
561  return theMass;
562  break;
563 
564  case Composite:
566  break;
567 
568  default:
569  INCL_ERROR("Particle::getRealMass: Unknown particle type." << '\n');
570  return 0.0;
571  break;
572  }
573  }
574 
577 
580 
583 
595  G4double getEmissionQValueCorrection(const G4int AParent, const G4int ZParent) const {
596  const G4int ADaughter = AParent - theA;
597  const G4int ZDaughter = ZParent - theZ;
598 
599  // Note the minus sign here
600  G4double theQValue;
601  if(isCluster())
602  theQValue = -ParticleTable::getTableQValue(theA, theZ, ADaughter, ZDaughter);
603  else {
604  const G4double massTableParent = ParticleTable::getTableMass(AParent,ZParent);
605  const G4double massTableDaughter = ParticleTable::getTableMass(ADaughter,ZDaughter);
606  const G4double massTableParticle = getTableMass();
607  theQValue = massTableParent - massTableDaughter - massTableParticle;
608  }
609 
610  const G4double massINCLParent = ParticleTable::getINCLMass(AParent,ZParent);
611  const G4double massINCLDaughter = ParticleTable::getINCLMass(ADaughter,ZDaughter);
612  const G4double massINCLParticle = getINCLMass();
613 
614  // The rhs corresponds to the INCL Q-value
615  return theQValue - (massINCLParent-massINCLDaughter-massINCLParticle);
616  }
617 
633  G4double getTransferQValueCorrection(const G4int AFrom, const G4int ZFrom, const G4int ATo, const G4int ZTo) const {
634  const G4int AFromDaughter = AFrom - theA;
635  const G4int ZFromDaughter = ZFrom - theZ;
636  const G4int AToDaughter = ATo + theA;
637  const G4int ZToDaughter = ZTo + theZ;
638  const G4double theQValue = ParticleTable::getTableQValue(AToDaughter,ZToDaughter,AFromDaughter,ZFromDaughter,AFrom,ZFrom);
639 
640  const G4double massINCLTo = ParticleTable::getINCLMass(ATo,ZTo);
641  const G4double massINCLToDaughter = ParticleTable::getINCLMass(AToDaughter,ZToDaughter);
642  /* Note that here we have to use the table mass in the INCL Q-value. We
643  * cannot use theMass, because at this stage the particle is probably
644  * still off-shell; and we cannot use getINCLMass(), because it leads to
645  * violations of global energy conservation.
646  */
647  const G4double massINCLParticle = getTableMass();
648 
649  // The rhs corresponds to the INCL Q-value for particle absorption
650  return theQValue - (massINCLToDaughter-massINCLTo-massINCLParticle);
651  }
652 
659  const G4double mass = std::pow(theEnergy, 2) - theMomentum.dot(theMomentum);
660  if(mass < 0.0) {
661  INCL_ERROR("E*E - p*p is negative." << '\n');
662  return 0.0;
663  } else {
664  return std::sqrt(mass);
665  }
666  };
667 
669  inline G4double getKineticEnergy() const { return theEnergy - theMass; }
670 
672  inline G4double getPotentialEnergy() const { return thePotentialEnergy; }
673 
676 
681  {
682  return theEnergy;
683  };
684 
688  void setMass(G4double mass)
689  {
690  this->theMass = mass;
691  }
692 
696  void setEnergy(G4double energy)
697  {
698  this->theEnergy = energy;
699  };
700 
705  {
706  return theMomentum;
707  };
708 
711  {
713  };
714 
718  virtual void setMomentum(const G4INCL::ThreeVector &momentum)
719  {
720  this->theMomentum = momentum;
721  };
722 
727  {
728  return thePosition;
729  };
730 
731  virtual void setPosition(const G4INCL::ThreeVector &position)
732  {
733  this->thePosition = position;
734  };
735 
737  void setHelicity(G4double h) { theHelicity = h; };
738 
739  void propagate(G4double step) {
740  thePosition += ((*thePropagationMomentum)*(step/(*thePropagationEnergy)));
741  };
742 
745 
748 
751 
753  G4int getNumberOfDecays() const { return nDecays; }
754 
757 
760 
769  void setOutOfWell() { outOfWell = true; }
770 
772  G4bool isOutOfWell() const { return outOfWell; }
773 
776 
780  }
781 
785  }
786 
789 
792 
793  G4bool isCluster() const {
794  return (theType == Composite);
795  }
796 
798  void setFrozenMomentum(const ThreeVector &momentum) { theFrozenMomentum = momentum; }
799 
801  void setFrozenEnergy(const G4double energy) { theFrozenEnergy = energy; }
802 
805 
808 
810  ThreeVector getPropagationVelocity() const { return (*thePropagationMomentum)/(*thePropagationEnergy); }
811 
821  }
822 
832  }
833 
839  virtual void rotatePositionAndMomentum(const G4double angle, const ThreeVector &axis) {
840  rotatePosition(angle, axis);
841  rotateMomentum(angle, axis);
842  }
843 
849  virtual void rotatePosition(const G4double angle, const ThreeVector &axis) {
850  thePosition.rotate(angle, axis);
851  }
852 
858  virtual void rotateMomentum(const G4double angle, const ThreeVector &axis) {
859  theMomentum.rotate(angle, axis);
860  theFrozenMomentum.rotate(angle, axis);
861  }
862 
863  std::string print() const {
864  std::stringstream ss;
865  ss << "Particle (ID = " << ID << ") type = ";
867  ss << '\n'
868  << " energy = " << theEnergy << '\n'
869  << " momentum = "
870  << theMomentum.print()
871  << '\n'
872  << " position = "
873  << thePosition.print()
874  << '\n';
875  return ss.str();
876  };
877 
878  std::string dump() const {
879  std::stringstream ss;
880  ss << "(particle " << ID << " ";
882  ss << '\n'
883  << thePosition.dump()
884  << '\n'
885  << theMomentum.dump()
886  << '\n'
887  << theEnergy << ")" << '\n';
888  return ss.str();
889  };
890 
891  long getID() const { return ID; };
892 
896  ParticleList const *getParticles() const {
897  INCL_WARN("Particle::getParticles() method was called on a Particle object" << '\n');
898  return 0;
899  }
900 
908  if(rpCorrelated)
909  return theMomentum.mag();
910  else
911  return uncorrelatedMomentum;
912  }
913 
916 
918  void rpCorrelate() { rpCorrelated = true; }
919 
921  void rpDecorrelate() { rpCorrelated = false; }
922 
926  if(norm>0.)
927  return thePosition.dot(*thePropagationMomentum) / std::sqrt(norm);
928  else
929  return 1.;
930  }
931 
933  static G4double getTotalBias();
934  static void setINCLBiasVector(std::vector<G4double> NewVector);
935  static void FillINCLBiasVector(G4double newBias);
936  static G4double getBiasFromVector(std::vector<G4int> VectorBias);
937 
938  static std::vector<G4int> MergeVectorBias(Particle const * const p1, Particle const * const p2);
939  static std::vector<G4int> MergeVectorBias(std::vector<G4int> p1, Particle const * const p2);
940 
943 
945  void setParticleBias(G4double ParticleBias) { this->theParticleBias = ParticleBias; }
946 
948  std::vector<G4int> getBiasCollisionVector() const { return theBiasCollisionVector; }
949 
951  void setBiasCollisionVector(std::vector<G4int> BiasCollisionVector) {
952  this->theBiasCollisionVector = BiasCollisionVector;
953  this->setParticleBias(Particle::getBiasFromVector(BiasCollisionVector));
954  }
955 
956  public:
970 #ifdef INCLXX_IN_GEANT4_MODE
971  static std::vector<G4double> INCLBiasVector;
972  //static G4VectorCache<G4double> INCLBiasVector;
973 #else
974  static G4ThreadLocal std::vector<G4double> INCLBiasVector;
975  //static G4VectorCache<G4double> INCLBiasVector;
976 #endif
978 
979  protected:
993  long ID;
994 
997 
999 
1000  private:
1004 
1006  std::vector<G4int> theBiasCollisionVector;
1007 
1009  static G4ThreadLocal long nextID;
1010 
1012  };
1013 }
1014 
1015 #endif /* PARTICLE_HH_ */
std::vector< G4int > getBiasCollisionVector() const
Get the vector list of biased vertices on the particle path.
ParticleList::const_iterator ParticleIter
void rotatePositionAndMomentum(const G4double angle, const ThreeVector &axis) const
virtual G4double getTableMass() const
Get the tabulated particle mass.
G4bool isNucleonorLambda() const
Is this a Nucleon or a Lambda?
void propagate(G4double step)
const G4INCL::ThreeVector & getPosition() const
void thawPropagation()
Unfreeze particle propagation.
virtual G4INCL::ParticleSpecies getSpecies() const
Get the particle species.
std::string dump() const
G4double * thePropagationEnergy
G4double getCosRPAngle() const
Get the cosine of the angle between position and momentum.
void incrementNumberOfCollisions()
Increment the number of collisions undergone by the particle.
ThreeVector vector(const ThreeVector &v) const
G4double getEnergy() const
const ThreeVector & adjustMomentumFromEnergy()
Rescale the momentum to match the total energy.
ThreeVector getFrozenMomentum() const
Get the frozen particle momentum.
void freezePropagation()
Freeze particle propagation.
void rpDecorrelate()
Make the particle not follow a strict r-p correlation.
void setTableMass()
Set the mass of the Particle to its table mass.
static std::vector< G4double > INCLBiasVector
Time ordered vector of all bias applied.
#define INCL_DECLARE_ALLOCATION_POOL(T)
G4double getKineticEnergy() const
Get the particle kinetic energy.
static void FillINCLBiasVector(G4double newBias)
void setMass(G4double mass)
G4bool isOmega() const
Is this an omega?
const G4INCL::ThreeVector & getMomentum() const
G4int getS() const
Returns the strangeness number.
G4INCL::ThreeVector theMomentum
G4INCL::ThreeVector theFrozenMomentum
void setNumberOfCollisions(G4int n)
Set the number of collisions undergone by the particle.
G4double dot(const ThreeVector &v) const
const char * p
Definition: xmltok.h:285
ParticleList const * getParticles() const
long getID() const
ParticipantType theParticipantType
G4double getHelicity()
ThreeVector getPropagationVelocity() const
Get the propagation velocity of the particle.
void setINCLMass()
Set the mass of the Particle to its table mass.
G4double uncorrelatedMomentum
G4bool isHyperon() const
Is this an Hyperon?
G4double getMass() const
Get the cached particle mass.
G4bool isSigma() const
Is this a Sigma?
virtual void makeTargetSpectator()
void setParticipantType(ParticipantType const p)
void setBiasCollisionVector(std::vector< G4int > BiasCollisionVector)
Set the vector list of biased vertices on the particle path.
G4bool isResonance() const
Is it a resonance?
static G4double getBiasFromVector(std::vector< G4int > VectorBias)
#define G4ThreadLocal
Definition: tls.hh:69
void setHelicity(G4double h)
void boost(const ThreeVector &b) const
G4bool isOutOfWell() const
Check if the particle is out of its potential well.
G4bool isLambda() const
Is this a Lambda?
G4double getRealMass(const G4INCL::ParticleType t)
Get particle mass (in MeV/c^2)
G4double getTransferQValueCorrection(const G4int AFrom, const G4int ZFrom, const G4int ATo, const G4int ZTo) const
Computes correction on the transfer Q-value.
static G4double getTotalBias()
General bias vector function.
G4bool isBaryon() const
Is this a Baryon?
void rotatePosition(const G4double angle, const ThreeVector &axis) const
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
G4double getEmissionTime()
G4INCL::ParticleType theType
#define position
Definition: xmlparse.cc:622
G4double getPotentialEnergy() const
Get the particle potential energy.
G4bool isDelta() const
Is it a Delta?
G4ThreadLocal NuclearMassFn getTableMass
Static pointer to the mass function for nuclei.
std::string print() const
virtual void rotatePositionAndMomentum(const G4double angle, const ThreeVector &axis)
Rotate the particle position and momentum.
G4int getNumberOfDecays() const
Return the number of decays undergone by the particle.
virtual void setPosition(const G4INCL::ThreeVector &position)
void rotate(const G4double angle, const ThreeVector &axis)
Rotate the vector by a given angle around a given axis.
G4double getInvariantMass() const
Get the the particle invariant mass.
static std::vector< G4int > MergeVectorBias(Particle const *const p1, Particle const *const p2)
G4double getReflectionMomentum() const
Return the reflection momentum.
#define INCL_ERROR(x)
G4double getINCLMass(const G4int A, const G4int Z)
Get INCL nuclear mass (in MeV/c^2)
G4bool isPhoton() const
Is this a photon?
static const G4double alpha
virtual void rotateMomentum(const G4double angle, const ThreeVector &axis)
Rotate the particle momentum.
G4INCL::ThreeVector thePosition
G4double thePotentialEnergy
G4double getBeta() const
double energy
Definition: plottest35.C:25
void setPotentialEnergy(G4double v)
Set the particle potential energy.
std::string print() const
G4int getZ() const
Returns the charge number.
void setRealMass()
Set the mass of the Particle to its real mass.
void setFrozenEnergy(const G4double energy)
Set the frozen particle momentum.
G4double theFrozenEnergy
Singleton for recycling allocation of instances of a given class.
G4int getA() const
Returns the baryon number.
void rpCorrelate()
Make the particle follow a strict r-p correlation.
void setType(ParticleType t)
G4bool isParticipant() const
G4bool isTargetSpectator() const
std::string getName(const ParticleType t)
Get the native INCL name of the particle.
G4double getINCLMass() const
Get the INCL particle mass.
virtual void rotatePosition(const G4double angle, const ThreeVector &axis)
Rotate the particle position.
static double P[]
virtual void setMomentum(const G4INCL::ThreeVector &momentum)
void swap(Particle &rhs)
Helper method for the assignment operator.
static const G4double bp
void incrementNumberOfDecays()
Increment the number of decays undergone by the particle.
void setParticleBias(G4double ParticleBias)
Set the particle bias.
std::vector< G4int > theBiasCollisionVector
Time ordered vector of all biased vertices on the particle path.
#define INCL_WARN(x)
Particle(const Particle &rhs)
Copy constructor.
static void setINCLBiasVector(std::vector< G4double > NewVector)
int G4int
Definition: G4Types.hh:78
virtual void makeParticipant()
G4double getRealMass() const
Get the real particle mass.
G4double getTableQValue(const G4int A1, const G4int Z1, const G4int A2, const G4int Z2)
Get Q-value (in MeV/c^2)
ThreeVector getLongitudinalPosition() const
Longitudinal component of the position w.r.t. the momentum.
virtual void makeProjectileSpectator()
G4int getNumberOfCollisions() const
Return the number of collisions undergone by the particle.
Float_t norm
void lorentzContract(const ThreeVector &aBoostVector, const ThreeVector &refPos)
Lorentz-contract the particle position around some center.
G4bool isMeson() const
Is this a Meson?
void setUncorrelatedMomentum(const G4double p)
Set the uncorrelated momentum.
G4INCL::ParticleType getType() const
void boost(const ThreeVector &aBoostVector)
G4double mag2() const
G4bool isKaon() const
Is this a Kaon?
G4double getEmissionQValueCorrection(const G4int AParent, const G4int ZParent) const
Computes correction on the emission Q-value.
void setEnergy(G4double energy)
G4bool isEta() const
Is this an eta?
G4double mag() const
G4double adjustEnergyFromMomentum()
Recompute the energy to match the momentum.
void setNumberOfDecays(G4int n)
Set the number of decays undergone by the particle.
G4bool isCluster() const
G4bool isPion() const
Is this a pion?
static G4ThreadLocal long nextID
static G4ThreadLocal G4int nextBiasedCollisionID
G4INCL::ThreeVector * thePropagationMomentum
G4bool isAntiKaon() const
Is this an antiKaon?
std::string dump() const
Particle & operator=(const Particle &rhs)
Assignment operator.
Char_t n[5]
virtual G4INCL::ThreeVector getAngularMomentum() const
void setEmissionTime(G4double t)
G4bool isStrange() const
Is this an Strange?
void setFrozenMomentum(const ThreeVector &momentum)
Set the frozen particle momentum.
G4bool isProjectileSpectator() const
G4bool isNucleon() const
G4ThreadLocal ParticleMassFn getTableParticleMass
Static pointer to the mass function for particles.
G4double getParticleBias() const
Get the particle bias.
void setOutOfWell()
Mark the particle as out of its potential well.
ThreeVector getTransversePosition() const
Transverse component of the position w.r.t. the momentum.
G4double getFrozenEnergy() const
Get the frozen particle momentum.
G4bool isEtaPrime() const
Is this an etaprime?
ThreeVector boostVector() const
G4double getParticleListBias() const
ParticleList::iterator ParticleMutableIter
ParticipantType getParticipantType() const
std::vector< G4int > getParticleListBiasVector() const
void rotateMomentum(const G4double angle, const ThreeVector &axis) const
G4double theParticleBias