Geant4 11.4.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4INCL::Particle Class Reference

#include <G4INCLParticle.hh>

Inheritance diagram for G4INCL::Particle:

Public Member Functions

 Particle ()
 Particle (ParticleType t, G4double energy, ThreeVector const &momentum, ThreeVector const &position)
 Particle (ParticleType t, ThreeVector const &momentum, ThreeVector const &position)
virtual ~Particle ()
 Particle (const Particle &rhs)
 Copy constructor.
Particleoperator= (const Particle &rhs)
 Assignment operator.
G4INCL::ParticleType getType () const
virtual G4INCL::ParticleSpecies getSpecies () const
 Get the particle species.
void setType (ParticleType t)
G4bool isNucleon () const
ParticipantType getParticipantType () const
void setParticipantType (ParticipantType const p)
G4bool isParticipant () const
G4bool isTargetSpectator () const
G4bool isProjectileSpectator () const
virtual void makeParticipant ()
virtual void makeTargetSpectator ()
virtual void makeProjectileSpectator ()
G4bool isPion () const
 Is this a pion?
G4bool isEta () const
 Is this an eta?
G4bool isOmega () const
 Is this an omega?
G4bool isEtaPrime () const
 Is this an etaprime?
G4bool isPhoton () const
 Is this a photon?
G4bool isResonance () const
 Is it a resonance?
G4bool isDelta () const
 Is it a Delta?
G4bool isSigma () const
 Is this a Sigma?
G4bool isKaon () const
 Is this a Kaon?
G4bool isAntiKaon () const
 Is this an antiKaon?
G4bool isLambda () const
 Is this a Lambda?
G4bool isNucleonorLambda () const
 Is this a Nucleon or a Lambda?
G4bool isHyperon () const
 Is this an Hyperon?
G4bool isMeson () const
 Is this a Meson?
G4bool isBaryon () const
 Is this a Baryon?
G4bool isStrange () const
 Is this a Strange?
G4bool isXi () const
 Is this a Xi?
G4bool isAntiNucleon () const
 Is this an antinucleon?
G4bool isAntiSigma () const
 Is this an antiSigma?
G4bool isAntiXi () const
 Is this an antiXi?
G4bool isAntiLambda () const
 Is this an antiLambda?
G4bool isAntiHyperon () const
 Is this an antiHyperon?
G4bool isAntiBaryon () const
 Is this an antiBaryon?
G4bool isAntiNucleonorAntiLambda () const
 Is this an antiNucleon or an antiLambda?
G4int getA () const
 Returns the baryon number.
G4int getZ () const
 Returns the charge number.
G4int getS () const
 Returns the strangeness number.
G4int getSrcPair () const
 Returns the strangeness number.
G4double getBeta () const
ThreeVector boostVector () const
void boost (const ThreeVector &aBoostVector)
void lorentzContract (const ThreeVector &aBoostVector, const ThreeVector &refPos)
 Lorentz-contract the particle position around some center.
G4double getMass () const
 Get the cached particle mass.
G4double getINCLMass () const
 Get the INCL particle mass.
virtual G4double getTableMass () const
 Get the tabulated particle mass.
G4double getRealMass () const
 Get the real particle mass.
void setRealMass ()
 Set the mass of the Particle to its real mass.
void setTableMass ()
 Set the mass of the Particle to its table mass.
void setINCLMass ()
 Set the mass of the Particle to its table mass.
G4double getEmissionQValueCorrection (const G4int AParent, const G4int ZParent) const
 Computes correction on the emission Q-value.
G4double getTransferQValueCorrection (const G4int AFrom, const G4int ZFrom, const G4int ATo, const G4int ZTo) const
 Computes correction on the transfer Q-value.
G4double getEmissionQValueCorrection (const G4int AParent, const G4int ZParent, const G4int SParent) const
 Computes correction on the emission Q-value for hypernuclei.
G4double getTransferQValueCorrection (const G4int AFrom, const G4int ZFrom, const G4int SFrom, const G4int ATo, const G4int ZTo, const G4int STo) const
 Computes correction on the transfer Q-value for hypernuclei.
G4double getInvariantMass () const
 Get the the particle invariant mass.
G4double getKineticEnergy () const
 Get the particle kinetic energy.
G4double getPotentialEnergy () const
 Get the particle potential energy.
void setPotentialEnergy (G4double v)
 Set the particle potential energy.
G4double getEnergy () const
void setMass (G4double mass)
void setEnergy (G4double energy)
const G4INCL::ThreeVectorgetMomentum () const
virtual G4INCL::ThreeVector getAngularMomentum () const
virtual void setMomentum (const G4INCL::ThreeVector &momentum)
const G4INCL::ThreeVectorgetPosition () const
virtual void setPosition (const G4INCL::ThreeVector &position)
G4double getHelicity ()
void setHelicity (G4double h)
void propagate (G4double step)
G4int getNumberOfCollisions () const
 Return the number of collisions undergone by the particle.
void setNumberOfCollisions (G4int n)
 Set the number of collisions undergone by the particle.
void incrementNumberOfCollisions ()
 Increment the number of collisions undergone by the particle.
G4int getNumberOfDecays () const
 Return the number of decays undergone by the particle.
void setNumberOfDecays (G4int n)
 Set the number of decays undergone by the particle.
void incrementNumberOfDecays ()
 Increment the number of decays undergone by the particle.
void setNumberOfSrcPair (int n)
 Set the number of srcpairs.
void setOutOfWell ()
 Mark the particle as out of its potential well.
G4bool isOutOfWell () const
 Check if the particle is out of its potential well.
void setSrcPartner ()
 Set and reset src partner.
void resetSrcPartner ()
G4bool isSrcPartner () const
 Check if the particle is a src partner.
void setEmissionTime (G4double t)
G4double getEmissionTime ()
ThreeVector getTransversePosition () const
 Transverse component of the position w.r.t. the momentum.
ThreeVector getLongitudinalPosition () const
 Longitudinal component of the position w.r.t. the momentum.
const ThreeVectoradjustMomentumFromEnergy ()
 Rescale the momentum to match the total energy.
G4double adjustEnergyFromMomentum ()
 Recompute the energy to match the momentum.
G4bool isCluster () const
void setFrozenMomentum (const ThreeVector &momentum)
 Set the frozen particle momentum.
void setFrozenEnergy (const G4double energy)
 Set the frozen particle momentum.
ThreeVector getFrozenMomentum () const
 Get the frozen particle momentum.
G4double getFrozenEnergy () const
 Get the frozen particle momentum.
ThreeVector getPropagationVelocity () const
 Get the propagation velocity of the particle.
void freezePropagation ()
 Freeze particle propagation.
void thawPropagation ()
 Unfreeze particle propagation.
virtual void rotatePositionAndMomentum (const G4double angle, const ThreeVector &axis)
 Rotate the particle position and momentum.
virtual void rotatePosition (const G4double angle, const ThreeVector &axis)
 Rotate the particle position.
virtual void rotateMomentum (const G4double angle, const ThreeVector &axis)
 Rotate the particle momentum.
std::string print () const
std::string dump () const
long getID () const
ParticleList const * getParticles () const
G4double getReflectionMomentum () const
 Return the reflection momentum.
void setUncorrelatedMomentum (const G4double p)
 Set the uncorrelated momentum.
void rpCorrelate ()
 Make the particle follow a strict r-p correlation.
void rpDecorrelate ()
 Make the particle not follow a strict r-p correlation.
G4double getCosRPAngle () const
 Get the cosine of the angle between position and momentum.
G4double getParticleBias () const
 Get the particle bias.
void setParticleBias (G4double ParticleBias)
 Set the particle bias.
std::vector< G4intgetBiasCollisionVector () const
 Get the vector list of biased vertices on the particle path.
void setBiasCollisionVector (std::vector< G4int > BiasCollisionVector)
 Set the vector list of biased vertices on the particle path.
G4int getNumberOfKaon () const
 Number of Kaon inside de nucleus.
void setNumberOfKaon (const G4int NK)
G4int getParentResonancePDGCode () const
void setParentResonancePDGCode (const G4int parentPDGCode)
G4int getParentResonanceID () const
void setParentResonanceID (const G4int parentID)

Static Public Member Functions

static G4double getTotalBias ()
 General bias vector function.
static void setINCLBiasVector (std::vector< G4double > NewVector)
static void FillINCLBiasVector (G4double newBias)
static G4double getBiasFromVector (std::vector< G4int > VectorBias)
static std::vector< G4intMergeVectorBias (Particle const *const p1, Particle const *const p2)
static std::vector< G4intMergeVectorBias (std::vector< G4int > p1, Particle const *const p2)

Static Public Attributes

static std::vector< G4doubleINCLBiasVector
 Time ordered vector of all bias applied.
static G4ThreadLocal G4int nextBiasedCollisionID = 0

Protected Member Functions

void swap (Particle &rhs)
 Helper method for the assignment operator.

Protected Attributes

G4int theZ
G4int theA
G4int theS
ParticipantType theParticipantType
G4INCL::ParticleType theType
G4double theEnergy
G4doublethePropagationEnergy
G4double theFrozenEnergy
G4INCL::ThreeVector theMomentum
G4INCL::ThreeVectorthePropagationMomentum
G4INCL::ThreeVector theFrozenMomentum
G4INCL::ThreeVector thePosition
G4int nCollisions
G4int nDecays
G4int nSrcPair
G4double thePotentialEnergy
long ID
G4bool rpCorrelated
G4double uncorrelatedMomentum
G4double theParticleBias
G4int theNKaon
 The number of Kaons inside the nucleus (update during the cascade).
G4int theParentResonancePDGCode
G4int theParentResonanceID

Detailed Description

Definition at line 75 of file G4INCLParticle.hh.

Constructor & Destructor Documentation

◆ Particle() [1/4]

G4INCL::Particle::Particle ( )

Definition at line 59 of file G4INCLParticle.cc.

60 : theZ(0), theA(0), theS(0),
63 theEnergy(0.0),
66 theMomentum(ThreeVector(0.,0.,0.)),
69 thePosition(ThreeVector(0.,0.,0.)),
70 nCollisions(0),
71 nDecays(0),
72 nSrcPair(0),
74 rpCorrelated(false),
77 theNKaon(0),
78#ifdef INCLXX_IN_GEANT4_MODE
81#endif
82 theHelicity(0.0),
83 emissionTime(0.0),
84 outOfWell(false),
85 theSrcPartner(false),
86 theMass(0.)
87 {
88 ID = nextID;
89 nextID++;
90 }
G4INCL::ThreeVector * thePropagationMomentum
G4INCL::ThreeVector theMomentum
ParticipantType theParticipantType
G4INCL::ParticleType theType
G4int theNKaon
The number of Kaons inside the nucleus (update during the cascade).
G4double uncorrelatedMomentum
G4double * thePropagationEnergy
G4INCL::ThreeVector thePosition
G4INCL::ThreeVector theFrozenMomentum

Referenced by G4INCL::Cluster::addParticle(), G4INCL::Cluster::Cluster(), G4INCL::Cluster::Cluster(), G4INCL::Cluster::Cluster(), MergeVectorBias(), MergeVectorBias(), operator=(), Particle(), G4INCL::Cluster::removeParticle(), G4INCL::ProjectileRemnant::reset(), G4INCL::ProjectileRemnant::storeComponents(), and swap().

◆ Particle() [2/4]

G4INCL::Particle::Particle ( ParticleType t,
G4double energy,
ThreeVector const & momentum,
ThreeVector const & position )

Definition at line 92 of file G4INCLParticle.cc.

94 : theEnergy(energy),
97 theMomentum(momentum),
100 thePosition(position),
101 nCollisions(0), nDecays(0), nSrcPair(0),
103 rpCorrelated(false),
105 theParticleBias(1.),
106 theNKaon(0),
107#ifdef INCLXX_IN_GEANT4_MODE
110#endif
111 theHelicity(0.0),
112 emissionTime(0.0), outOfWell(false), theSrcPartner(false)
113 {
115 ID = nextID;
116 nextID++;
117 if(theEnergy <= 0.0) {
118 INCL_WARN("Particle with energy " << theEnergy << " created." << '\n');
119 }
120 setType(t);
122 }
#define INCL_WARN(x)
void setMass(G4double mass)
G4double getInvariantMass() const
Get the the particle invariant mass.
void setType(ParticleType t)

◆ Particle() [3/4]

G4INCL::Particle::Particle ( ParticleType t,
ThreeVector const & momentum,
ThreeVector const & position )

Definition at line 124 of file G4INCLParticle.cc.

127 theMomentum(momentum),
130 thePosition(position),
131 nCollisions(0), nDecays(0),
132 nSrcPair(0),
134 rpCorrelated(false),
136 theParticleBias(1.),
137 theNKaon(0),
138#ifdef INCLXX_IN_GEANT4_MODE
141#endif
142 theHelicity(0.0),
143 emissionTime(0.0), outOfWell(false), theSrcPartner(false)
144 {
146 ID = nextID;
147 nextID++;
148 setType(t);
149 if( isResonance() ) {
150 INCL_ERROR("Cannot create resonance without specifying its momentum four-vector." << '\n');
151 }
152 G4double energy = std::sqrt(theMomentum.mag2() + theMass*theMass);
155 }
#define INCL_ERROR(x)
double G4double
Definition G4Types.hh:83
G4bool isResonance() const
Is it a resonance?
G4double energy(const ThreeVector &p, const G4double m)

◆ ~Particle()

virtual G4INCL::Particle::~Particle ( )
inlinevirtual

Definition at line 80 of file G4INCLParticle.hh.

80{}

◆ Particle() [4/4]

G4INCL::Particle::Particle ( const Particle & rhs)
inline

Copy constructor.

Does not copy the particle ID.

Definition at line 86 of file G4INCLParticle.hh.

86 :
87 theZ(rhs.theZ),
88 theA(rhs.theA),
89 theS(rhs.theS),
90 theParticipantType(rhs.theParticipantType),
91 theType(rhs.theType),
92 theEnergy(rhs.theEnergy),
93 theFrozenEnergy(rhs.theFrozenEnergy),
94 theMomentum(rhs.theMomentum),
95 theFrozenMomentum(rhs.theFrozenMomentum),
96 thePosition(rhs.thePosition),
97 nCollisions(rhs.nCollisions),
98 nDecays(rhs.nDecays),
99 nSrcPair(rhs.nSrcPair),
100 thePotentialEnergy(rhs.thePotentialEnergy),
101 rpCorrelated(rhs.rpCorrelated),
102 uncorrelatedMomentum(rhs.uncorrelatedMomentum),
103 theParticleBias(rhs.theParticleBias),
104 theNKaon(rhs.theNKaon),
105#ifdef INCLXX_IN_GEANT4_MODE
106 theParentResonancePDGCode(rhs.theParentResonancePDGCode),
107 theParentResonanceID(rhs.theParentResonanceID),
108#endif
109 theHelicity(rhs.theHelicity),
110 emissionTime(rhs.emissionTime),
111 outOfWell(rhs.outOfWell),
112 theSrcPartner(rhs.theSrcPartner),
113 theMass(rhs.theMass)
114 {
115 if(rhs.thePropagationEnergy == &(rhs.theFrozenEnergy))
117 else
119 if(rhs.thePropagationMomentum == &(rhs.theFrozenMomentum))
121 else
123 // ID intentionally not copied
124 ID = nextID++;
125
126 theBiasCollisionVector = rhs.theBiasCollisionVector;
127 }

Member Function Documentation

◆ adjustEnergyFromMomentum()

◆ adjustMomentumFromEnergy()

◆ boost()

void G4INCL::Particle::boost ( const ThreeVector & aBoostVector)
inline

Boost the particle using a boost vector.

Example (go to the particle rest frame): particle->boost(particle->boostVector());

Definition at line 517 of file G4INCLParticle.hh.

517 {
518 const G4double beta2 = aBoostVector.mag2();
519 const G4double gamma = 1.0 / std::sqrt(1.0 - beta2);
520 const G4double bp = theMomentum.dot(aBoostVector);
521 const G4double alpha = (gamma*gamma)/(1.0 + gamma);
522
523 theMomentum = theMomentum + aBoostVector * (alpha * bp - gamma * theEnergy);
524 theEnergy = gamma * (theEnergy - bp);
525 }

Referenced by G4INCL::Cluster::boost(), G4INCL::Nucleus::decayOutgoingPionResonances(), G4INCL::Nucleus::decayOutgoingSigmaZero(), and G4INCL::PhaseSpaceKopylov::generate().

◆ boostVector()

ThreeVector G4INCL::Particle::boostVector ( ) const
inline

Returns a three vector we can give to the boost() -method.

In order to go to the particle rest frame you need to multiply the boost vector by -1.0.

Definition at line 507 of file G4INCLParticle.hh.

507 {
508 return theMomentum / theEnergy;
509 }

Referenced by G4INCL::PhaseSpaceKopylov::generate(), G4INCL::StandardPropagationModel::shootComposite(), and G4INCL::StandardPropagationModel::shootParticle().

◆ dump()

std::string G4INCL::Particle::dump ( ) const
inline

Definition at line 1114 of file G4INCLParticle.hh.

1114 {
1115 std::stringstream ss;
1116 ss << "(particle " << ID << " ";
1118 ss << nSrcPair << " ";
1119 ss << '\n'
1120 << thePosition.dump()
1121 << '\n'
1122 << theMomentum.dump()
1123 << '\n'
1124 << theEnergy << ")" << '\n';
1125 return ss.str();
1126 };
std::string getName(const ParticleType t)
Get the native INCL name of the particle.

◆ FillINCLBiasVector()

void G4INCL::Particle::FillINCLBiasVector ( G4double newBias)
static

Definition at line 217 of file G4INCLParticle.cc.

217 {
218// assert(G4int(Particle::INCLBiasVector.size())==nextBiasedCollisionID);
219 //assert(G4int(Particle::INCLBiasVector.Size())==nextBiasedCollisionID);
220// assert(std::fabs(newBias - 1.) > 1E-6);
221 Particle::INCLBiasVector.push_back(newBias);
222 //Particle::INCLBiasVector.Push_back(newBias);
224 }
static std::vector< G4double > INCLBiasVector
Time ordered vector of all bias applied.
static G4ThreadLocal G4int nextBiasedCollisionID

Referenced by G4INCL::InteractionAvatar::postInteraction().

◆ freezePropagation()

void G4INCL::Particle::freezePropagation ( )
inline

Freeze particle propagation.

Make the particle use theFrozenMomentum and theFrozenEnergy for propagation. The normal state can be restored by calling the thawPropagation() method.

Definition at line 1052 of file G4INCLParticle.hh.

◆ getA()

◆ getAngularMomentum()

virtual G4INCL::ThreeVector G4INCL::Particle::getAngularMomentum ( ) const
inlinevirtual

Get the angular momentum w.r.t. the origin

Reimplemented in G4INCL::Cluster.

Definition at line 934 of file G4INCLParticle.hh.

935 {
936 return thePosition.vector(theMomentum);
937 };

Referenced by G4INCL::Cluster::getAngularMomentum(), G4INCL::StandardPropagationModel::shootCompositeAtrest(), and G4INCL::StandardPropagationModel::shootParticle().

◆ getBeta()

G4double G4INCL::Particle::getBeta ( ) const
inline

Definition at line 496 of file G4INCLParticle.hh.

496 {
497 const G4double P = theMomentum.mag();
498 return P/theEnergy;
499 }

◆ getBiasCollisionVector()

std::vector< G4int > G4INCL::Particle::getBiasCollisionVector ( ) const
inline

Get the vector list of biased vertices on the particle path.

Definition at line 1185 of file G4INCLParticle.hh.

1185{ return theBiasCollisionVector; }

Referenced by G4INCL::ClusterDecay::decay(), G4INCL::Nucleus::decayMe(), G4INCL::Nucleus::decayOutgoingClusters(), G4INCL::Nucleus::decayOutgoingPionResonances(), G4INCL::Nucleus::decayOutgoingSigmaZero(), MergeVectorBias(), and MergeVectorBias().

◆ getBiasFromVector()

G4double G4INCL::Particle::getBiasFromVector ( std::vector< G4int > VectorBias)
static

Definition at line 226 of file G4INCLParticle.cc.

226 {
227 if(VectorBias.empty()) return 1.;
228
229 G4double ParticleBias = 1.;
230
231 for(G4int i=0; i<G4int(VectorBias.size()); i++){
232 ParticleBias *= Particle::INCLBiasVector[G4int(VectorBias[i])];
233 }
234
235 return ParticleBias;
236 }
int G4int
Definition G4Types.hh:85

Referenced by G4INCL::BinaryCollisionAvatar::getChannel(), G4INCL::ParticleList::getParticleListBias(), and setBiasCollisionVector().

◆ getCosRPAngle()

G4double G4INCL::Particle::getCosRPAngle ( ) const
inline

Get the cosine of the angle between position and momentum.

Definition at line 1161 of file G4INCLParticle.hh.

1161 {
1162 const G4double norm = thePosition.mag2()*thePropagationMomentum->mag2();
1163 if(norm>0.)
1164 return thePosition.dot(*thePropagationMomentum) / std::sqrt(norm);
1165 else
1166 return 1.;
1167 }

◆ getEmissionQValueCorrection() [1/2]

G4double G4INCL::Particle::getEmissionQValueCorrection ( const G4int AParent,
const G4int ZParent ) const
inline

Computes correction on the emission Q-value.

Computes the correction that must be applied to INCL particles in order to obtain the correct Q-value for particle emission from a given nucleus. For absorption, the correction is obviously equal to minus the value returned by this function.

Parameters
AParentthe mass number of the emitting nucleus
ZParentthe charge number of the emitting nucleus
Returns
the correction

Definition at line 736 of file G4INCLParticle.hh.

736 {
737 const G4int SParent = 0;
738 const G4int ADaughter = AParent - theA;
739 const G4int ZDaughter = ZParent - theZ;
740 const G4int SDaughter = 0;
741
742 // Note the minus sign here
743 G4double theQValue;
744 if(isCluster())
745 theQValue = -ParticleTable::getTableQValue(theA, theZ, theS, ADaughter, ZDaughter, SDaughter);
746 else {
747 const G4double massTableParent = ParticleTable::getTableMass(AParent,ZParent,SParent);
748 const G4double massTableDaughter = ParticleTable::getTableMass(ADaughter,ZDaughter,SDaughter);
749 const G4double massTableParticle = getTableMass();
750 theQValue = massTableParent - massTableDaughter - massTableParticle;
751 }
752
753 const G4double massINCLParent = ParticleTable::getINCLMass(AParent,ZParent,SParent);
754 const G4double massINCLDaughter = ParticleTable::getINCLMass(ADaughter,ZDaughter,SDaughter);
755 const G4double massINCLParticle = getINCLMass();
756
757 // The rhs corresponds to the INCL Q-value
758 return theQValue - (massINCLParent-massINCLDaughter-massINCLParticle);
759 }
G4bool isCluster() const
G4double getINCLMass() const
Get the INCL particle mass.
virtual G4double getTableMass() const
Get the tabulated particle mass.
G4ThreadLocal NuclearMassFn getTableMass
Static pointer to the mass function for nuclei.
G4double getTableQValue(const G4int A1, const G4int Z1, const G4int S1, const G4int A2, const G4int Z2, const G4int S2)
Get Q-value (in MeV/c^2).
G4double getINCLMass(const G4int A, const G4int Z, const G4int S)
Get INCL nuclear mass (in MeV/c^2).

Referenced by G4INCL::Nucleus::emitInsideAntilambda(), G4INCL::Nucleus::emitInsideKaon(), G4INCL::Nucleus::emitInsideLambda(), G4INCL::Nucleus::emitInsidePions(), G4INCL::Nucleus::emitInsideStrangeParticles(), and G4INCL::SurfaceAvatar::getTransmissionProbability().

◆ getEmissionQValueCorrection() [2/2]

G4double G4INCL::Particle::getEmissionQValueCorrection ( const G4int AParent,
const G4int ZParent,
const G4int SParent ) const
inline

Computes correction on the emission Q-value for hypernuclei.

Computes the correction that must be applied to INCL particles in order to obtain the correct Q-value for particle emission from a given nucleus. For absorption, the correction is obviously equal to minus the value returned by this function.

Parameters
AParentthe mass number of the emitting nucleus
ZParentthe charge number of the emitting nucleus
SParentthe strangess number of the emitting nucleus
Returns
the correction

Definition at line 812 of file G4INCLParticle.hh.

812 {
813 const G4int ADaughter = AParent - theA;
814 const G4int ZDaughter = ZParent - theZ;
815 const G4int SDaughter = SParent - theS;
816
817 // Note the minus sign here
818 G4double theQValue;
819 if(isCluster())
820 theQValue = -ParticleTable::getTableQValue(theA, theZ, theS, ADaughter, ZDaughter, SDaughter);
821 else {
822 const G4double massTableParent = ParticleTable::getTableMass(AParent,ZParent,SParent);
823 const G4double massTableDaughter = ParticleTable::getTableMass(ADaughter,ZDaughter,SDaughter);
824 const G4double massTableParticle = getTableMass();
825 theQValue = massTableParent - massTableDaughter - massTableParticle;
826 }
827
828 const G4double massINCLParent = ParticleTable::getINCLMass(AParent,ZParent,SParent);
829 const G4double massINCLDaughter = ParticleTable::getINCLMass(ADaughter,ZDaughter,SDaughter);
830 const G4double massINCLParticle = getINCLMass();
831
832 // The rhs corresponds to the INCL Q-value
833 return theQValue - (massINCLParent-massINCLDaughter-massINCLParticle);
834 }

◆ getEmissionTime()

G4double G4INCL::Particle::getEmissionTime ( )
inline

Definition at line 1009 of file G4INCLParticle.hh.

1009{ return emissionTime; };

◆ getEnergy()

◆ getFrozenEnergy()

G4double G4INCL::Particle::getFrozenEnergy ( ) const
inline

Get the frozen particle momentum.

Definition at line 1041 of file G4INCLParticle.hh.

1041{ return theFrozenEnergy; }

◆ getFrozenMomentum()

ThreeVector G4INCL::Particle::getFrozenMomentum ( ) const
inline

Get the frozen particle momentum.

Definition at line 1038 of file G4INCLParticle.hh.

1038{ return theFrozenMomentum; }

◆ getHelicity()

G4double G4INCL::Particle::getHelicity ( )
inline

Definition at line 960 of file G4INCLParticle.hh.

960{ return theHelicity; };

◆ getID()

◆ getINCLMass()

G4double G4INCL::Particle::getINCLMass ( ) const
inline

Get the INCL particle mass.

Definition at line 549 of file G4INCLParticle.hh.

549 {
550 switch(theType) {
551 case Proton:
552 case Neutron:
553 case PiPlus:
554 case PiMinus:
555 case PiZero:
556 case Lambda:
557 case SigmaPlus:
558 case SigmaZero:
559 case SigmaMinus:
560 case antiProton:
561 case XiZero:
562 case XiMinus:
563 case antiNeutron:
564 case antiLambda:
565 case antiSigmaPlus:
566 case antiSigmaZero:
567 case antiSigmaMinus:
568 case antiXiZero:
569 case antiXiMinus:
570 case KPlus:
571 case KZero:
572 case KZeroBar:
573 case KShort:
574 case KLong:
575 case KMinus:
576 case Eta:
577 case Omega:
578 case EtaPrime:
579 case Photon:
581 break;
582
583 case DeltaPlusPlus:
584 case DeltaPlus:
585 case DeltaZero:
586 case DeltaMinus:
587 return theMass;
588 break;
589
590 case Composite:
592 break;
593 case antiComposite:
595 break;
596
597 default:
598 INCL_ERROR("Particle::getINCLMass: Unknown particle type." << '\n');
599 return 0.0;
600 break;
601 }
602 }

Referenced by getEmissionQValueCorrection(), getEmissionQValueCorrection(), G4INCL::StandardPropagationModel::getTime(), and setINCLMass().

◆ getInvariantMass()

G4double G4INCL::Particle::getInvariantMass ( ) const
inline

Get the the particle invariant mass.

Uses the relativistic invariant

\‍[ m = \sqrt{E^2 - {\vec p}^2}\‍]

Definition at line 882 of file G4INCLParticle.hh.

882 {
883 const G4double mass = std::pow(theEnergy, 2) - theMomentum.dot(theMomentum);
884 if(mass < 0.0) {
885 INCL_ERROR("E*E - p*p is negative." << '\n');
886 return 0.0;
887 } else {
888 return std::sqrt(mass);
889 }
890 };

Referenced by Particle().

◆ getKineticEnergy()

◆ getLongitudinalPosition()

ThreeVector G4INCL::Particle::getLongitudinalPosition ( ) const
inline

Longitudinal component of the position w.r.t. the momentum.

Definition at line 1017 of file G4INCLParticle.hh.

Referenced by getTransversePosition().

◆ getMass()

G4double G4INCL::Particle::getMass ( ) const
inline

Get the cached particle mass.

Definition at line 546 of file G4INCLParticle.hh.

546{ return theMass; }

Referenced by G4INCL::DeltaDecayChannel::computeDecayTime(), G4INCL::PionResonanceDecayChannel::computeDecayTime(), G4INCL::SigmaZeroDecayChannel::computeDecayTime(), G4INCL::Nucleus::emitInsideAntilambda(), G4INCL::Nucleus::emitInsideKaon(), G4INCL::Nucleus::emitInsideLambda(), G4INCL::Nucleus::emitInsidePions(), G4INCL::Nucleus::emitInsideStrangeParticles(), G4INCL::InteractionAvatar::enforceEnergyConservation(), G4INCL::CrossSectionsMultiPionsAndResonances::etaNToPiN(), G4INCL::EtaNElasticChannel::fillFinalState(), G4INCL::EtaNToPiNChannel::fillFinalState(), G4INCL::EtaOrOmegaNToLKChannel::fillFinalState(), G4INCL::EtaOrOmegaNToSKChannel::fillFinalState(), G4INCL::NSToNLChannel::fillFinalState(), G4INCL::NSToNSChannel::fillFinalState(), G4INCL::OmegaNElasticChannel::fillFinalState(), G4INCL::OmegaNToPiNChannel::fillFinalState(), G4INCL::Cluster::freezeInternalMotion(), G4INCL::PhaseSpaceKopylov::generate(), G4INCL::NuclearPotential::INuclearPotential::getFermiMomentum(), G4INCL::KinematicsUtils::getLocalEnergy(), G4INCL::SurfaceAvatar::getTransmissionProbability(), G4INCL::CrossSections::interactionDistanceKbarN(), G4INCL::CrossSections::interactionDistanceKN(), G4INCL::CrossSections::interactionDistanceNbarN(), G4INCL::CrossSections::interactionDistancenbarN(), G4INCL::CrossSections::interactionDistanceNN(), G4INCL::CrossSections::interactionDistancePiN(), G4INCL::CrossSections::interactionDistanceYN(), G4INCL::Cluster::internalBoostToCM(), G4INCL::KinematicsUtils::momentumInCM(), G4INCL::KinematicsUtils::momentumInLab(), G4INCL::CrossSectionsINCL46::NDeltaToNN(), G4INCL::CrossSectionsMultiPions::NDeltaToNN(), G4INCL::CrossSectionsStrangeness::NNToNLK2pi(), G4INCL::CrossSectionsStrangeness::NNToNLKpi(), G4INCL::CrossSectionsStrangeness::NNToNSK2pi(), G4INCL::CrossSectionsStrangeness::NNToNSKpi(), G4INCL::CrossSectionsStrangeness::omegaNToLK(), G4INCL::CrossSectionsMultiPionsAndResonances::omegaNToPiN(), G4INCL::CrossSectionsStrangeness::omegaNToSK(), G4INCL::CrossSectionsMultiPionsAndResonances::piMinuspToEtaN(), G4INCL::CrossSectionsMultiPionsAndResonances::piMinuspToOmegaN(), G4INCL::Cluster::print(), G4INCL::ProjectileRemnant::ProjectileRemnant(), G4INCL::StandardPropagationModel::shootCompositeAtrest(), and G4INCL::StandardPropagationModel::shootParticle().

◆ getMomentum()

const G4INCL::ThreeVector & G4INCL::Particle::getMomentum ( ) const
inline

◆ getNumberOfCollisions()

G4int G4INCL::Particle::getNumberOfCollisions ( ) const
inline

Return the number of collisions undergone by the particle.

Definition at line 968 of file G4INCLParticle.hh.

968{ return nCollisions; }

Referenced by G4INCL::Cluster::addParticle().

◆ getNumberOfDecays()

G4int G4INCL::Particle::getNumberOfDecays ( ) const
inline

Return the number of decays undergone by the particle.

Definition at line 977 of file G4INCLParticle.hh.

977{ return nDecays; }

◆ getNumberOfKaon()

G4int G4INCL::Particle::getNumberOfKaon ( ) const
inline

Number of Kaon inside de nucleus.

Put in the Particle class in order to calculate the "correct" mass of composit particle.

Definition at line 1200 of file G4INCLParticle.hh.

1200{ return theNKaon; };

◆ getParentResonanceID()

G4int G4INCL::Particle::getParentResonanceID ( ) const
inline

Definition at line 1206 of file G4INCLParticle.hh.

1206{ return theParentResonanceID; };

◆ getParentResonancePDGCode()

G4int G4INCL::Particle::getParentResonancePDGCode ( ) const
inline

Definition at line 1204 of file G4INCLParticle.hh.

1204{ return theParentResonancePDGCode; };

◆ getParticipantType()

ParticipantType G4INCL::Particle::getParticipantType ( ) const
inline

Definition at line 378 of file G4INCLParticle.hh.

378 {
379 return theParticipantType;
380 }

Referenced by G4INCL::StandardPropagationModel::generateBinaryCollisionAvatar().

◆ getParticleBias()

G4double G4INCL::Particle::getParticleBias ( ) const
inline

Get the particle bias.

Definition at line 1179 of file G4INCLParticle.hh.

1179{ return theParticleBias; };

◆ getParticles()

ParticleList const * G4INCL::Particle::getParticles ( ) const
inline

Return a NULL pointer

Definition at line 1133 of file G4INCLParticle.hh.

1133 {
1134 INCL_WARN("Particle::getParticles() method was called on a Particle object" << '\n');
1135 return 0;
1136 }

◆ getPosition()

◆ getPotentialEnergy()

◆ getPropagationVelocity()

ThreeVector G4INCL::Particle::getPropagationVelocity ( ) const
inline

Get the propagation velocity of the particle.

Definition at line 1044 of file G4INCLParticle.hh.

1044{ return (*thePropagationMomentum)/(*thePropagationEnergy); }

Referenced by G4INCL::CoulombNone::bringToSurface(), G4INCL::CoulombNone::bringToSurfaceAbar(), G4INCL::StandardPropagationModel::getReflectionTime(), and G4INCL::StandardPropagationModel::getTime().

◆ getRealMass()

G4double G4INCL::Particle::getRealMass ( ) const
inline

Get the real particle mass.

Definition at line 661 of file G4INCLParticle.hh.

661 {
662 switch(theType) {
663 case Proton:
664 case Neutron:
665 case PiPlus:
666 case PiMinus:
667 case PiZero:
668 case Lambda:
669 case SigmaPlus:
670 case SigmaZero:
671 case SigmaMinus:
672 case antiProton:
673 case XiZero:
674 case XiMinus:
675 case antiNeutron:
676 case antiLambda:
677 case antiSigmaPlus:
678 case antiSigmaZero:
679 case antiSigmaMinus:
680 case antiXiZero:
681 case antiXiMinus:
682 case KPlus:
683 case KZero:
684 case KZeroBar:
685 case KShort:
686 case KLong:
687 case KMinus:
688 case Eta:
689 case Omega:
690 case EtaPrime:
691 case Photon:
693 break;
694
695 case DeltaPlusPlus:
696 case DeltaPlus:
697 case DeltaZero:
698 case DeltaMinus:
699 return theMass;
700 break;
701
702 case Composite:
704 break;
705 case antiComposite:
707 break;
708
709 default:
710 INCL_ERROR("Particle::getRealMass: Unknown particle type." << '\n');
711 return 0.0;
712 break;
713 }
714 }
G4double getRealMass(const G4INCL::ParticleType t)
Get particle mass (in MeV/c^2).

Referenced by G4INCL::Cluster::getTableMass(), and setRealMass().

◆ getReflectionMomentum()

G4double G4INCL::Particle::getReflectionMomentum ( ) const
inline

Return the reflection momentum.

The reflection momentum is used by calls to getSurfaceRadius to compute the radius of the sphere where the nucleon moves. It is necessary to introduce fuzzy r-p correlations.

Definition at line 1144 of file G4INCLParticle.hh.

1144 {
1145 if(rpCorrelated)
1146 return theMomentum.mag();
1147 else
1148 return uncorrelatedMomentum;
1149 }

Referenced by G4INCL::KinematicsUtils::getLocalEnergy(), and G4INCL::Nucleus::getSurfaceRadius().

◆ getS()

◆ getSpecies()

virtual G4INCL::ParticleSpecies G4INCL::Particle::getSpecies ( ) const
inlinevirtual

Get the particle species.

Reimplemented in G4INCL::Cluster.

Definition at line 196 of file G4INCLParticle.hh.

196 {
197 return ParticleSpecies(theType);
198 };

Referenced by G4INCL::CoulombDistortion::maxImpactParameter(), and G4INCL::StandardPropagationModel::shootParticle().

◆ getSrcPair()

G4int G4INCL::Particle::getSrcPair ( ) const
inline

Returns the strangeness number.

Definition at line 494 of file G4INCLParticle.hh.

494{ return nSrcPair; }

◆ getTableMass()

virtual G4double G4INCL::Particle::getTableMass ( ) const
inlinevirtual

Get the tabulated particle mass.

Reimplemented in G4INCL::Cluster.

Definition at line 605 of file G4INCLParticle.hh.

605 {
606 switch(theType) {
607 case Proton:
608 case Neutron:
609 case PiPlus:
610 case PiMinus:
611 case PiZero:
612 case Lambda:
613 case SigmaPlus:
614 case SigmaZero:
615 case SigmaMinus:
616 case antiProton:
617 case XiZero:
618 case XiMinus:
619 case antiNeutron:
620 case antiLambda:
621 case antiSigmaPlus:
622 case antiSigmaZero:
623 case antiSigmaMinus:
624 case antiXiZero:
625 case antiXiMinus:
626 case KPlus:
627 case KZero:
628 case KZeroBar:
629 case KShort:
630 case KLong:
631 case KMinus:
632 case Eta:
633 case Omega:
634 case EtaPrime:
635 case Photon:
637 break;
638
639 case DeltaPlusPlus:
640 case DeltaPlus:
641 case DeltaZero:
642 case DeltaMinus:
643 return theMass;
644 break;
645
646 case Composite:
648 break;
649 case antiComposite:
651 break;
652
653 default:
654 INCL_ERROR("Particle::getTableMass: Unknown particle type." << '\n');
655 return 0.0;
656 break;
657 }
658 }
G4ThreadLocal ParticleMassFn getTableParticleMass
Static pointer to the mass function for particles.

Referenced by G4INCL::Nucleus::decayOutgoingPionResonances(), G4INCL::Nucleus::decayOutgoingSigmaZero(), getEmissionQValueCorrection(), getEmissionQValueCorrection(), getTransferQValueCorrection(), getTransferQValueCorrection(), and setTableMass().

◆ getTotalBias()

G4double G4INCL::Particle::getTotalBias ( )
static

◆ getTransferQValueCorrection() [1/2]

G4double G4INCL::Particle::getTransferQValueCorrection ( const G4int AFrom,
const G4int ZFrom,
const G4int ATo,
const G4int ZTo ) const
inline

Computes correction on the transfer Q-value.

Computes the correction that must be applied to INCL particles in order to obtain the correct Q-value for particle transfer from a given nucleus to another.

Assumes that the receving nucleus is INCL's target nucleus, with the INCL separation energy.

Parameters
AFromthe mass number of the donating nucleus
ZFromthe charge number of the donating nucleus
ATothe mass number of the receiving nucleus
ZTothe charge number of the receiving nucleus
Returns
the correction

Definition at line 776 of file G4INCLParticle.hh.

776 {
777 const G4int SFrom = 0;
778 const G4int STo = 0;
779 const G4int AFromDaughter = AFrom - theA;
780 const G4int ZFromDaughter = ZFrom - theZ;
781 const G4int SFromDaughter = 0;
782 const G4int AToDaughter = ATo + theA;
783 const G4int ZToDaughter = ZTo + theZ;
784 const G4int SToDaughter = 0;
785 const G4double theQValue = ParticleTable::getTableQValue(AToDaughter,ZToDaughter,SToDaughter,AFromDaughter,ZFromDaughter,SFromDaughter,AFrom,ZFrom,SFrom);
786
787 const G4double massINCLTo = ParticleTable::getINCLMass(ATo,ZTo,STo);
788 const G4double massINCLToDaughter = ParticleTable::getINCLMass(AToDaughter,ZToDaughter,SToDaughter);
789 /* Note that here we have to use the table mass in the INCL Q-value. We
790 * cannot use theMass, because at this stage the particle is probably
791 * still off-shell; and we cannot use getINCLMass(), because it leads to
792 * violations of global energy conservation.
793 */
794 const G4double massINCLParticle = getTableMass();
795
796 // The rhs corresponds to the INCL Q-value for particle absorption
797 return theQValue - (massINCLToDaughter-massINCLTo-massINCLParticle);
798 }

◆ getTransferQValueCorrection() [2/2]

G4double G4INCL::Particle::getTransferQValueCorrection ( const G4int AFrom,
const G4int ZFrom,
const G4int SFrom,
const G4int ATo,
const G4int ZTo,
const G4int STo ) const
inline

Computes correction on the transfer Q-value for hypernuclei.

Computes the correction that must be applied to INCL particles in order to obtain the correct Q-value for particle transfer from a given nucleus to another.

Assumes that the receving nucleus is INCL's target nucleus, with the INCL separation energy.

Parameters
AFromthe mass number of the donating nucleus
ZFromthe charge number of the donating nucleus
SFromthe strangess number of the donating nucleus
ATothe mass number of the receiving nucleus
ZTothe charge number of the receiving nucleus
STothe strangess number of the receiving nucleus
Returns
the correction

Definition at line 853 of file G4INCLParticle.hh.

853 {
854 const G4int AFromDaughter = AFrom - theA;
855 const G4int ZFromDaughter = ZFrom - theZ;
856 const G4int SFromDaughter = SFrom - theS;
857 const G4int AToDaughter = ATo + theA;
858 const G4int ZToDaughter = ZTo + theZ;
859 const G4int SToDaughter = STo + theS;
860 const G4double theQValue = ParticleTable::getTableQValue(AToDaughter,ZToDaughter,SFromDaughter,AFromDaughter,ZFromDaughter,SToDaughter,AFrom,ZFrom,SFrom);
861
862 const G4double massINCLTo = ParticleTable::getINCLMass(ATo,ZTo,STo);
863 const G4double massINCLToDaughter = ParticleTable::getINCLMass(AToDaughter,ZToDaughter,SToDaughter);
864 /* Note that here we have to use the table mass in the INCL Q-value. We
865 * cannot use theMass, because at this stage the particle is probably
866 * still off-shell; and we cannot use getINCLMass(), because it leads to
867 * violations of global energy conservation.
868 */
869 const G4double massINCLParticle = getTableMass();
870
871 // The rhs corresponds to the INCL Q-value for particle absorption
872 return theQValue - (massINCLToDaughter-massINCLTo-massINCLParticle);
873 }

◆ getTransversePosition()

ThreeVector G4INCL::Particle::getTransversePosition ( ) const
inline

Transverse component of the position w.r.t. the momentum.

Definition at line 1012 of file G4INCLParticle.hh.

1012 {
1014 }
ThreeVector getLongitudinalPosition() const
Longitudinal component of the position w.r.t. the momentum.

Referenced by G4INCL::StandardPropagationModel::shootComposite(), G4INCL::StandardPropagationModel::shootCompositeAtrest(), and G4INCL::StandardPropagationModel::shootParticle().

◆ getType()

G4INCL::ParticleType G4INCL::Particle::getType ( ) const
inline

Get the particle type.

See also
G4INCL::ParticleType

Definition at line 191 of file G4INCLParticle.hh.

191 {
192 return theType;
193 };

Referenced by G4INCL::PionResonanceDecayChannel::computeDecayTime(), G4INCL::NuclearPotential::INuclearPotential::computeKaonPotentialEnergy(), G4INCL::NuclearPotential::INuclearPotential::computePionPotentialEnergy(), G4INCL::NuclearPotential::INuclearPotential::computePionResonancePotentialEnergy(), G4INCL::NuclearPotential::NuclearPotentialConstant::computePotentialEnergy(), G4INCL::NuclearPotential::NuclearPotentialIsospin::computePotentialEnergy(), G4INCL::CrossSectionsINCL46::elasticNNLegacy(), G4INCL::NDeltaToDeltaLKChannel::fillFinalState(), G4INCL::NDeltaToDeltaSKChannel::fillFinalState(), G4INCL::NDeltaToNSKChannel::fillFinalState(), G4INCL::NKbToL2piChannel::fillFinalState(), G4INCL::NKbToLpiChannel::fillFinalState(), G4INCL::NKbToNKb2piChannel::fillFinalState(), G4INCL::NKbToNKbChannel::fillFinalState(), G4INCL::NKbToNKbpiChannel::fillFinalState(), G4INCL::NKbToS2piChannel::fillFinalState(), G4INCL::NKbToSpiChannel::fillFinalState(), G4INCL::NKToNK2piChannel::fillFinalState(), G4INCL::NKToNKChannel::fillFinalState(), G4INCL::NKToNKpiChannel::fillFinalState(), G4INCL::NpiToMissingStrangenessChannel::fillFinalState(), G4INCL::NSToNLChannel::fillFinalState(), G4INCL::NSToNSChannel::fillFinalState(), G4INCL::PauliStandard::getBlockingProbability(), G4INCL::NuclearPotential::INuclearPotential::getFermiEnergy(), G4INCL::NuclearPotential::INuclearPotential::getFermiMomentum(), G4INCL::KinematicsUtils::getLocalEnergy(), G4INCL::NuclearPotential::INuclearPotential::getSeparationEnergy(), G4INCL::Nucleus::getSurfaceRadius(), G4INCL::StandardPropagationModel::getTime(), G4INCL::NuclearDensity::getTransmissionRadius(), G4INCL::Nucleus::insertParticle(), G4INCL::ParticleConfig::isPair(), G4INCL::CrossSectionsStrangeness::NDeltaToDeltaLK(), G4INCL::CrossSectionsStrangeness::NDeltaToDeltaSK(), G4INCL::CrossSectionsStrangeness::NDeltaToNLK(), G4INCL::CrossSectionsINCL46::NDeltaToNN(), G4INCL::CrossSectionsMultiPions::NDeltaToNN(), G4INCL::CrossSectionsStrangeness::NDeltaToNNKKb(), G4INCL::CrossSectionsStrangeness::NDeltaToNSK(), G4INCL::CrossSectionsStrangeness::NKbToL2pi(), G4INCL::CrossSectionsStrangeness::NKbToLpi(), G4INCL::CrossSectionsStrangeness::NKbToNKb(), G4INCL::CrossSectionsStrangeness::NKbToNKb2pi(), G4INCL::CrossSectionsStrangeness::NKbToNKbpi(), G4INCL::CrossSectionsStrangeness::NKbToS2pi(), G4INCL::CrossSectionsStrangeness::NKbToSpi(), G4INCL::CrossSectionsStrangeness::NKToNK(), G4INCL::CrossSectionsStrangeness::NKToNK2pi(), G4INCL::CrossSectionsStrangeness::NKToNKpi(), G4INCL::CrossSectionsAntiparticles::NNbarCEX(), G4INCL::CrossSectionsAntiparticles::NNbarElastic(), G4INCL::CrossSectionsAntiparticles::NNbarToAnnihilation(), G4INCL::CrossSectionsAntiparticles::NNbarToLLbar(), G4INCL::CrossSectionsAntiparticles::NNbarToNNbar2pi(), G4INCL::CrossSectionsAntiparticles::NNbarToNNbar3pi(), G4INCL::CrossSectionsAntiparticles::NNbarToNNbarpi(), G4INCL::CrossSectionsMultiPions::NNElastic(), G4INCL::CrossSectionsMultiPions::NNOnePi(), G4INCL::CrossSectionsMultiPions::NNOnePiOrDelta(), G4INCL::CrossSectionsMultiPions::NNThreePi(), G4INCL::CrossSectionsStrangeness::NNToMissingStrangeness(), G4INCL::CrossSectionsINCL46::NNToNDelta(), G4INCL::CrossSectionsMultiPions::NNToNDelta(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNDeltaEta(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNDeltaOmega(), G4INCL::CrossSectionsStrangeness::NNToNLK(), G4INCL::CrossSectionsStrangeness::NNToNLK2pi(), G4INCL::CrossSectionsStrangeness::NNToNLKpi(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNEta(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNEtaExclu(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNEtaFourPi(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNEtaOnePi(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNEtaOnePiOrDelta(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNEtaThreePi(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNEtaTwoPi(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNEtaxPi(), G4INCL::CrossSectionsStrangeness::NNToNNKKb(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNOmega(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNOmegaExclu(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNOmegaFourPi(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNOmegaOnePi(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNOmegaOnePiOrDelta(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNOmegaThreePi(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNOmegaTwoPi(), G4INCL::CrossSectionsMultiPionsAndResonances::NNToNNOmegaxPi(), G4INCL::CrossSectionsStrangeness::NNToNSK(), G4INCL::CrossSectionsStrangeness::NNToNSK2pi(), G4INCL::CrossSectionsStrangeness::NNToNSKpi(), G4INCL::CrossSectionsMultiPions::NNTot(), G4INCL::CrossSectionsMultiPions::NNTwoPi(), G4INCL::CrossSectionsStrangeness::NpiToLK(), G4INCL::CrossSectionsStrangeness::NpiToLK2pi(), G4INCL::CrossSectionsStrangeness::NpiToLKpi(), G4INCL::CrossSectionsStrangeness::NpiToNKKb(), G4INCL::CrossSectionsStrangeness::NpiToSK(), G4INCL::CrossSectionsStrangeness::NpiToSK2pi(), G4INCL::CrossSectionsStrangeness::NpiToSKpi(), G4INCL::CrossSectionsStrangeness::NSToNL(), G4INCL::CrossSectionsStrangeness::NSToNS(), G4INCL::CrossSectionsStrangeness::p_pimToSzKz(), G4INCL::CrossSectionsINCL46::piNToDelta(), G4INCL::CrossSectionsMultiPions::piNToDelta(), G4INCL::CrossSectionsMultiPionsAndResonances::piNToEtaN(), G4INCL::CrossSectionsMultiPionsAndResonances::piNToOmegaN(), G4INCL::CrossSectionsMultiPions::piNTot(), and G4INCL::CrossSectionsAntiparticles::total().

◆ getZ()

◆ incrementNumberOfCollisions()

void G4INCL::Particle::incrementNumberOfCollisions ( )
inline

Increment the number of collisions undergone by the particle.

Definition at line 974 of file G4INCLParticle.hh.

974{ nCollisions++; }

◆ incrementNumberOfDecays()

void G4INCL::Particle::incrementNumberOfDecays ( )
inline

Increment the number of decays undergone by the particle.

Definition at line 983 of file G4INCLParticle.hh.

983{ nDecays++; }

◆ isAntiBaryon()

G4bool G4INCL::Particle::isAntiBaryon ( ) const
inline

Is this an antiBaryon?

Definition at line 479 of file G4INCLParticle.hh.

479{ return (isAntiNucleon() || isAntiHyperon()); }
G4bool isAntiHyperon() const
Is this an antiHyperon?
G4bool isAntiNucleon() const
Is this an antinucleon?

◆ isAntiHyperon()

G4bool G4INCL::Particle::isAntiHyperon ( ) const
inline

Is this an antiHyperon?

Definition at line 476 of file G4INCLParticle.hh.

476{ return (isAntiLambda() || isAntiSigma() || isAntiXi()); }
G4bool isAntiSigma() const
Is this an antiSigma?
G4bool isAntiXi() const
Is this an antiXi?
G4bool isAntiLambda() const
Is this an antiLambda?

Referenced by isAntiBaryon().

◆ isAntiKaon()

◆ isAntiLambda()

G4bool G4INCL::Particle::isAntiLambda ( ) const
inline

Is this an antiLambda?

Definition at line 473 of file G4INCLParticle.hh.

473{ return (theType == antiLambda); }

Referenced by isAntiHyperon(), and isAntiNucleonorAntiLambda().

◆ isAntiNucleon()

◆ isAntiNucleonorAntiLambda()

G4bool G4INCL::Particle::isAntiNucleonorAntiLambda ( ) const
inline

Is this an antiNucleon or an antiLambda?

Definition at line 482 of file G4INCLParticle.hh.

482{ return (isAntiNucleon() || isAntiLambda()); }

◆ isAntiSigma()

G4bool G4INCL::Particle::isAntiSigma ( ) const
inline

Is this an antiSigma?

Definition at line 467 of file G4INCLParticle.hh.

Referenced by isAntiHyperon().

◆ isAntiXi()

G4bool G4INCL::Particle::isAntiXi ( ) const
inline

Is this an antiXi?

Definition at line 470 of file G4INCLParticle.hh.

470{ return (theType == antiXiZero || theType == antiXiMinus); }

Referenced by isAntiHyperon().

◆ isBaryon()

G4bool G4INCL::Particle::isBaryon ( ) const
inline

Is this a Baryon?

Definition at line 455 of file G4INCLParticle.hh.

455{ return (isNucleon() || isResonance() || isHyperon()); }
G4bool isHyperon() const
Is this an Hyperon?
G4bool isNucleon() const

Referenced by G4INCL::CDPP::processOneParticle().

◆ isCluster()

G4bool G4INCL::Particle::isCluster ( ) const
inline

◆ isDelta()

◆ isEta()

◆ isEtaPrime()

G4bool G4INCL::Particle::isEtaPrime ( ) const
inline

◆ isHyperon()

G4bool G4INCL::Particle::isHyperon ( ) const
inline

Is this an Hyperon?

Definition at line 449 of file G4INCLParticle.hh.

449{ return (isLambda() || isSigma() ); } //|| isXi()
G4bool isLambda() const
Is this a Lambda?
G4bool isSigma() const
Is this a Sigma?

Referenced by G4INCL::CrossSectionsAntiparticles::elastic(), G4INCL::CrossSectionsStrangeness::elastic(), isBaryon(), isStrange(), and G4INCL::CrossSectionsStrangeness::NYelastic().

◆ isKaon()

◆ isLambda()

◆ isMeson()

G4bool G4INCL::Particle::isMeson ( ) const
inline

Is this a Meson?

Definition at line 452 of file G4INCLParticle.hh.

452{ return (isPion() || isKaon() || isAntiKaon() || isEta() || isEtaPrime() || isOmega()); }
G4bool isEtaPrime() const
Is this an etaprime?
G4bool isOmega() const
Is this an omega?
G4bool isEta() const
Is this an eta?
G4bool isPion() const
Is this a pion?
G4bool isAntiKaon() const
Is this an antiKaon?
G4bool isKaon() const
Is this a Kaon?

Referenced by G4INCL::StandardPropagationModel::generateBinaryCollisionAvatar(), G4INCL::InteractionAvatar::preInteractionLocalEnergy(), G4INCL::CDPP::processOneParticle(), and G4INCL::StandardPropagationModel::shootParticle().

◆ isNucleon()

G4bool G4INCL::Particle::isNucleon ( ) const
inline

Is this a nucleon?

Definition at line 371 of file G4INCLParticle.hh.

371 {
373 return true;
374 else
375 return false;
376 };

Referenced by G4INCL::NuclearPotential::NuclearPotentialEnergyIsospin::computePotentialEnergy(), G4INCL::NuclearPotential::NuclearPotentialEnergyIsospinSmooth::computePotentialEnergy(), G4INCL::CrossSectionsAntiparticles::elastic(), G4INCL::CrossSectionsINCL46::elastic(), G4INCL::CrossSectionsMultiPions::elastic(), G4INCL::CrossSectionsMultiPionsAndResonances::elastic(), G4INCL::CrossSectionsStrangeness::elastic(), G4INCL::CrossSectionsTruncatedMultiPions::elastic(), G4INCL::CrossSectionsINCL46::elasticNNLegacy(), G4INCL::StandardPropagationModel::generateBinaryCollisionAvatar(), G4INCL::Nucleus::getSurfaceRadius(), G4INCL::Nucleus::insertParticle(), isBaryon(), isNucleonorLambda(), G4INCL::CrossSectionsMultiPions::NNElastic(), G4INCL::CrossSectionsMultiPions::NNTot(), G4INCL::CrossSectionsStrangeness::omegaNToLK(), G4INCL::CrossSectionsStrangeness::omegaNToSK(), G4INCL::CrossSectionsMultiPions::piMinuspIne(), G4INCL::CrossSectionsMultiPions::piMinuspOnePi(), G4INCL::CrossSectionsMultiPions::piMinuspTwoPi(), G4INCL::CrossSectionsMultiPions::piNIne(), G4INCL::CrossSectionsMultiPions::piNOnePi(), G4INCL::CrossSectionsINCL46::piNToDelta(), G4INCL::CrossSectionsMultiPions::piNToxPiN(), G4INCL::CrossSectionsMultiPions::piNTwoPi(), G4INCL::CrossSectionsMultiPions::piPluspIne(), G4INCL::CrossSectionsMultiPions::piPluspOnePi(), G4INCL::CrossSectionsMultiPions::piPluspTwoPi(), G4INCL::CrossSectionsAntiparticles::total(), G4INCL::CrossSectionsINCL46::total(), G4INCL::CrossSectionsMultiPions::total(), G4INCL::CrossSectionsMultiPionsAndResonances::total(), and G4INCL::CrossSectionsStrangeness::total().

◆ isNucleonorLambda()

G4bool G4INCL::Particle::isNucleonorLambda ( ) const
inline

Is this a Nucleon or a Lambda?

Definition at line 446 of file G4INCLParticle.hh.

446{ return (isNucleon() || isLambda()); }

◆ isOmega()

◆ isOutOfWell()

G4bool G4INCL::Particle::isOutOfWell ( ) const
inline

◆ isParticipant()

G4bool G4INCL::Particle::isParticipant ( ) const
inline

◆ isPhoton()

G4bool G4INCL::Particle::isPhoton ( ) const
inline

◆ isPion()

◆ isProjectileSpectator()

G4bool G4INCL::Particle::isProjectileSpectator ( ) const
inline

Definition at line 394 of file G4INCLParticle.hh.

394 {
396 }

◆ isResonance()

G4bool G4INCL::Particle::isResonance ( ) const
inline

◆ isSigma()

◆ isSrcPartner()

G4bool G4INCL::Particle::isSrcPartner ( ) const
inline

Check if the particle is a src partner.

Definition at line 1006 of file G4INCLParticle.hh.

1006{ return theSrcPartner; }

◆ isStrange()

G4bool G4INCL::Particle::isStrange ( ) const
inline

Is this a Strange?

Definition at line 458 of file G4INCLParticle.hh.

458{ return (isKaon() || isAntiKaon() || isHyperon()); }

◆ isTargetSpectator()

G4bool G4INCL::Particle::isTargetSpectator ( ) const
inline

Definition at line 390 of file G4INCLParticle.hh.

390 {
392 }

Referenced by G4INCL::Nucleus::insertParticle().

◆ isXi()

G4bool G4INCL::Particle::isXi ( ) const
inline

Is this a Xi?

Definition at line 461 of file G4INCLParticle.hh.

461{ return (theType == XiZero || theType == XiMinus); }

◆ lorentzContract()

void G4INCL::Particle::lorentzContract ( const ThreeVector & aBoostVector,
const ThreeVector & refPos )
inline

Lorentz-contract the particle position around some center.

Apply Lorentz contraction to the position component along the direction of the boost vector.

Parameters
aBoostVectorthe boost vector (velocity) [c]
refPosthe reference position

Definition at line 535 of file G4INCLParticle.hh.

535 {
536 const G4double beta2 = aBoostVector.mag2();
537 const G4double gamma = 1.0 / std::sqrt(1.0 - beta2);
538 const ThreeVector theRelativePosition = thePosition - refPos;
539 const ThreeVector transversePosition = theRelativePosition - aBoostVector * (theRelativePosition.dot(aBoostVector) / aBoostVector.mag2());
540 const ThreeVector longitudinalPosition = theRelativePosition - transversePosition;
541
542 thePosition = refPos + transversePosition + longitudinalPosition / gamma;
543 }

◆ makeParticipant()

virtual void G4INCL::Particle::makeParticipant ( )
inlinevirtual

Reimplemented in G4INCL::Cluster.

Definition at line 398 of file G4INCLParticle.hh.

398 {
400 }

Referenced by G4INCL::Store::loadParticles(), and G4INCL::Cluster::makeParticipant().

◆ makeProjectileSpectator()

virtual void G4INCL::Particle::makeProjectileSpectator ( )
inlinevirtual

◆ makeTargetSpectator()

virtual void G4INCL::Particle::makeTargetSpectator ( )
inlinevirtual

Reimplemented in G4INCL::Cluster.

Definition at line 402 of file G4INCLParticle.hh.

402 {
404 }

Referenced by G4INCL::Cluster::makeTargetSpectator().

◆ MergeVectorBias() [1/2]

std::vector< G4int > G4INCL::Particle::MergeVectorBias ( Particle const *const p1,
Particle const *const p2 )
static

Definition at line 238 of file G4INCLParticle.cc.

238 {
239 std::vector<G4int> MergedVectorBias;
240 std::vector<G4int> VectorBias1 = p1->getBiasCollisionVector();
241 std::vector<G4int> VectorBias2 = p2->getBiasCollisionVector();
242 G4int i = 0;
243 G4int j = 0;
244 if(VectorBias1.size()==0 && VectorBias2.size()==0) return MergedVectorBias;
245 else if(VectorBias1.size()==0) return VectorBias2;
246 else if(VectorBias2.size()==0) return VectorBias1;
247
248 while(i < G4int(VectorBias1.size()) || j < G4int(VectorBias2.size())){
249 if(VectorBias1[i]==VectorBias2[j]){
250 MergedVectorBias.push_back(VectorBias1[i]);
251 i++;
252 j++;
253 if(i == G4int(VectorBias1.size())){
254 for(;j<G4int(VectorBias2.size());j++) MergedVectorBias.push_back(VectorBias2[j]);
255 }
256 else if(j == G4int(VectorBias2.size())){
257 for(;i<G4int(VectorBias1.size());i++) MergedVectorBias.push_back(VectorBias1[i]);
258 }
259 } else if(VectorBias1[i]<VectorBias2[j]){
260 MergedVectorBias.push_back(VectorBias1[i]);
261 i++;
262 if(i == G4int(VectorBias1.size())){
263 for(;j<G4int(VectorBias2.size());j++) MergedVectorBias.push_back(VectorBias2[j]);
264 }
265 }
266 else {
267 MergedVectorBias.push_back(VectorBias2[j]);
268 j++;
269 if(j == G4int(VectorBias2.size())){
270 for(;i<G4int(VectorBias1.size());i++) MergedVectorBias.push_back(VectorBias1[i]);
271 }
272 }
273 }
274 return MergedVectorBias;
275 }

Referenced by G4INCL::BinaryCollisionAvatar::getChannel(), G4INCL::ParticleList::getParticleListBias(), and G4INCL::ParticleList::getParticleListBiasVector().

◆ MergeVectorBias() [2/2]

std::vector< G4int > G4INCL::Particle::MergeVectorBias ( std::vector< G4int > p1,
Particle const *const p2 )
static

Definition at line 277 of file G4INCLParticle.cc.

277 {
278 std::vector<G4int> MergedVectorBias;
279 std::vector<G4int> VectorBias = p2->getBiasCollisionVector();
280 G4int i = 0;
281 G4int j = 0;
282 if(p1.size()==0 && VectorBias.size()==0) return MergedVectorBias;
283 else if(p1.size()==0) return VectorBias;
284 else if(VectorBias.size()==0) return p1;
285
286 while(i < G4int(p1.size()) || j < G4int(VectorBias.size())){
287 if(p1[i]==VectorBias[j]){
288 MergedVectorBias.push_back(p1[i]);
289 i++;
290 j++;
291 if(i == G4int(p1.size())){
292 for(;j<G4int(VectorBias.size());j++) MergedVectorBias.push_back(VectorBias[j]);
293 }
294 else if(j == G4int(VectorBias.size())){
295 for(;i<G4int(p1.size());i++) MergedVectorBias.push_back(p1[i]);
296 }
297 } else if(p1[i]<VectorBias[j]){
298 MergedVectorBias.push_back(p1[i]);
299 i++;
300 if(i == G4int(p1.size())){
301 for(;j<G4int(VectorBias.size());j++) MergedVectorBias.push_back(VectorBias[j]);
302 }
303 }
304 else {
305 MergedVectorBias.push_back(VectorBias[j]);
306 j++;
307 if(j == G4int(VectorBias.size())){
308 for(;i<G4int(p1.size());i++) MergedVectorBias.push_back(p1[i]);
309 }
310 }
311 }
312 return MergedVectorBias;
313 }

◆ operator=()

Particle & G4INCL::Particle::operator= ( const Particle & rhs)
inline

Assignment operator.

Does not copy the particle ID.

Definition at line 181 of file G4INCLParticle.hh.

181 {
182 Particle temporaryParticle(rhs);
183 swap(temporaryParticle);
184 return *this;
185 }
void swap(Particle &rhs)
Helper method for the assignment operator.

Referenced by G4INCL::Cluster::operator=().

◆ print()

std::string G4INCL::Particle::print ( ) const
inline

Definition at line 1097 of file G4INCLParticle.hh.

1097 {
1098 std::stringstream ss;
1099 ss << "Particle (ID = " << ID << ") type = ";
1101 ss << ", SRC pair = " << nSrcPair;
1102 ss << ", Potential energy = " << thePotentialEnergy;
1103 ss << '\n'
1104 << " energy = " << theEnergy << '\n'
1105 << " momentum = "
1106 << theMomentum.print()
1107 << '\n'
1108 << " position = "
1109 << thePosition.print()
1110 << '\n';
1111 return ss.str();
1112 };

Referenced by adjustMomentumFromEnergy(), G4INCL::Nucleus::emitInsideAnnihilationProducts(), G4INCL::Nucleus::emitInsideAntilambda(), G4INCL::Nucleus::emitInsideKaon(), G4INCL::Nucleus::emitInsideLambda(), G4INCL::Nucleus::emitInsidePions(), G4INCL::Nucleus::emitInsideStrangeParticles(), G4INCL::KinematicsUtils::getLocalEnergy(), G4INCL::StandardPropagationModel::getReflectionTime(), G4INCL::CrossSectionsMultiPions::piNTot(), G4INCL::ProjectileRemnant::removeParticle(), and G4INCL::Nucleus::restoreSrcPartner().

◆ propagate()

void G4INCL::Particle::propagate ( G4double step)
inline

Definition at line 963 of file G4INCLParticle.hh.

963 {
964 thePosition += ((*thePropagationMomentum)*(step/(*thePropagationEnergy)));
965 };

Referenced by G4INCL::StandardPropagationModel::generateBinaryCollisionAvatar().

◆ resetSrcPartner()

void G4INCL::Particle::resetSrcPartner ( )
inline

Definition at line 1003 of file G4INCLParticle.hh.

1003{ theSrcPartner = false; nSrcPair=0; }

◆ rotateMomentum()

virtual void G4INCL::Particle::rotateMomentum ( const G4double angle,
const ThreeVector & axis )
inlinevirtual

Rotate the particle momentum.

Parameters
anglethe rotation angle
axisa unit vector representing the rotation axis

Reimplemented in G4INCL::Cluster.

Definition at line 1092 of file G4INCLParticle.hh.

1092 {
1093 theMomentum.rotate(angle, axis);
1094 theFrozenMomentum.rotate(angle, axis);
1095 }
const axis_t axis_to_type< N >::axis
Definition pugixml.cc:9668

Referenced by G4INCL::Cluster::rotateMomentum(), and rotatePositionAndMomentum().

◆ rotatePosition()

virtual void G4INCL::Particle::rotatePosition ( const G4double angle,
const ThreeVector & axis )
inlinevirtual

Rotate the particle position.

Parameters
anglethe rotation angle
axisa unit vector representing the rotation axis

Reimplemented in G4INCL::Cluster.

Definition at line 1083 of file G4INCLParticle.hh.

1083 {
1084 thePosition.rotate(angle, axis);
1085 }

Referenced by G4INCL::Cluster::rotatePosition(), and rotatePositionAndMomentum().

◆ rotatePositionAndMomentum()

virtual void G4INCL::Particle::rotatePositionAndMomentum ( const G4double angle,
const ThreeVector & axis )
inlinevirtual

Rotate the particle position and momentum.

Parameters
anglethe rotation angle
axisa unit vector representing the rotation axis

Definition at line 1073 of file G4INCLParticle.hh.

1073 {
1074 rotatePosition(angle, axis);
1075 rotateMomentum(angle, axis);
1076 }
virtual void rotateMomentum(const G4double angle, const ThreeVector &axis)
Rotate the particle momentum.
virtual void rotatePosition(const G4double angle, const ThreeVector &axis)
Rotate the particle position.

◆ rpCorrelate()

void G4INCL::Particle::rpCorrelate ( )
inline

Make the particle follow a strict r-p correlation.

Definition at line 1155 of file G4INCLParticle.hh.

1155{ rpCorrelated = true; }

Referenced by G4INCL::InteractionAvatar::bringParticleInside(), and G4INCL::SurfaceAvatar::postInteraction().

◆ rpDecorrelate()

void G4INCL::Particle::rpDecorrelate ( )
inline

Make the particle not follow a strict r-p correlation.

Definition at line 1158 of file G4INCLParticle.hh.

1158{ rpCorrelated = false; }

◆ setBiasCollisionVector()

void G4INCL::Particle::setBiasCollisionVector ( std::vector< G4int > BiasCollisionVector)
inline

Set the vector list of biased vertices on the particle path.

Definition at line 1188 of file G4INCLParticle.hh.

1188 {
1189 this->theBiasCollisionVector = BiasCollisionVector;
1190 this->setParticleBias(Particle::getBiasFromVector(std::move(BiasCollisionVector)));
1191 }
void setParticleBias(G4double ParticleBias)
Set the particle bias.
static G4double getBiasFromVector(std::vector< G4int > VectorBias)

Referenced by G4INCL::Nucleus::decayOutgoingPionResonances(), G4INCL::Nucleus::decayOutgoingSigmaZero(), and G4INCL::SurfaceAvatar::postInteraction().

◆ setEmissionTime()

◆ setEnergy()

◆ setFrozenEnergy()

void G4INCL::Particle::setFrozenEnergy ( const G4double energy)
inline

Set the frozen particle momentum.

Definition at line 1035 of file G4INCLParticle.hh.

◆ setFrozenMomentum()

void G4INCL::Particle::setFrozenMomentum ( const ThreeVector & momentum)
inline

Set the frozen particle momentum.

Definition at line 1032 of file G4INCLParticle.hh.

1032{ theFrozenMomentum = momentum; }

◆ setHelicity()

void G4INCL::Particle::setHelicity ( G4double h)
inline

Definition at line 961 of file G4INCLParticle.hh.

961{ theHelicity = h; };

◆ setINCLBiasVector()

void G4INCL::Particle::setINCLBiasVector ( std::vector< G4double > NewVector)
static

Definition at line 321 of file G4INCLParticle.cc.

321 {
322 Particle::INCLBiasVector = std::move(NewVector);
323 }

◆ setINCLMass()

void G4INCL::Particle::setINCLMass ( )
inline

Set the mass of the Particle to its table mass.

Definition at line 723 of file G4INCLParticle.hh.

723{ setMass(getINCLMass()); }

Referenced by G4INCL::Cluster::Cluster(), G4INCL::Cluster::Cluster(), setType(), and G4INCL::StandardPropagationModel::shootParticle().

◆ setMass()

void G4INCL::Particle::setMass ( G4double mass)
inline

◆ setMomentum()

◆ setNumberOfCollisions()

void G4INCL::Particle::setNumberOfCollisions ( G4int n)
inline

Set the number of collisions undergone by the particle.

Definition at line 971 of file G4INCLParticle.hh.

◆ setNumberOfDecays()

void G4INCL::Particle::setNumberOfDecays ( G4int n)
inline

Set the number of decays undergone by the particle.

Definition at line 980 of file G4INCLParticle.hh.

980{ nDecays = n; }

◆ setNumberOfKaon()

void G4INCL::Particle::setNumberOfKaon ( const G4int NK)
inline

Definition at line 1201 of file G4INCLParticle.hh.

1201{ theNKaon = NK; }

◆ setNumberOfSrcPair()

void G4INCL::Particle::setNumberOfSrcPair ( int n)
inline

Set the number of srcpairs.

Definition at line 986 of file G4INCLParticle.hh.

986{ nSrcPair = n; }

◆ setOutOfWell()

void G4INCL::Particle::setOutOfWell ( )
inline

Mark the particle as out of its potential well.

This flag is used to control pions created outside their potential well in delta decay. The pion potential checks it and returns zero if it is true (necessary in order to correctly enforce energy conservation). The Nucleus::applyFinalState() method uses it to determine whether new avatars should be generated for the particle.

Definition at line 996 of file G4INCLParticle.hh.

996{ outOfWell = true; }

◆ setParentResonanceID()

void G4INCL::Particle::setParentResonanceID ( const G4int parentID)
inline

Definition at line 1207 of file G4INCLParticle.hh.

1207{ theParentResonanceID = parentID; };

◆ setParentResonancePDGCode()

void G4INCL::Particle::setParentResonancePDGCode ( const G4int parentPDGCode)
inline

Definition at line 1205 of file G4INCLParticle.hh.

1205{ theParentResonancePDGCode = parentPDGCode; };

◆ setParticipantType()

void G4INCL::Particle::setParticipantType ( ParticipantType const p)
inline

Definition at line 382 of file G4INCLParticle.hh.

382 {
384 }

◆ setParticleBias()

void G4INCL::Particle::setParticleBias ( G4double ParticleBias)
inline

Set the particle bias.

Definition at line 1182 of file G4INCLParticle.hh.

1182{ this->theParticleBias = ParticleBias; }

Referenced by setBiasCollisionVector().

◆ setPosition()

◆ setPotentialEnergy()

◆ setRealMass()

void G4INCL::Particle::setRealMass ( )
inline

Set the mass of the Particle to its real mass.

Definition at line 717 of file G4INCLParticle.hh.

717{ setMass(getRealMass()); }
G4double getRealMass() const
Get the real particle mass.

Referenced by G4INCL::ClusterDecay::decay().

◆ setSrcPartner()

void G4INCL::Particle::setSrcPartner ( )
inline

Set and reset src partner.

Definition at line 1002 of file G4INCLParticle.hh.

1002{ theSrcPartner = true; }

◆ setTableMass()

◆ setType()

void G4INCL::Particle::setType ( ParticleType t)
inline

Definition at line 200 of file G4INCLParticle.hh.

200 {
201 theType = t;
202 switch(theType)
203 {
204 case DeltaPlusPlus:
205 theA = 1;
206 theZ = 2;
207 theS = 0;
208 break;
209 case Proton:
210 case DeltaPlus:
211 theA = 1;
212 theZ = 1;
213 theS = 0;
214 break;
215 case Neutron:
216 case DeltaZero:
217 theA = 1;
218 theZ = 0;
219 theS = 0;
220 break;
221 case DeltaMinus:
222 theA = 1;
223 theZ = -1;
224 theS = 0;
225 break;
226 case PiPlus:
227 theA = 0;
228 theZ = 1;
229 theS = 0;
230 break;
231 case PiZero:
232 case Eta:
233 case Omega:
234 case EtaPrime:
235 case Photon:
236 theA = 0;
237 theZ = 0;
238 theS = 0;
239 break;
240 case PiMinus:
241 theA = 0;
242 theZ = -1;
243 theS = 0;
244 break;
245 case Lambda:
246 theA = 1;
247 theZ = 0;
248 theS = -1;
249 break;
250 case SigmaPlus:
251 theA = 1;
252 theZ = 1;
253 theS = -1;
254 break;
255 case SigmaZero:
256 theA = 1;
257 theZ = 0;
258 theS = -1;
259 break;
260 case SigmaMinus:
261 theA = 1;
262 theZ = -1;
263 theS = -1;
264 break;
265 case antiProton:
266 theA = -1;
267 theZ = -1;
268 theS = 0;
269 break;
270 case XiMinus:
271 theA = 1;
272 theZ = -1;
273 theS = -2;
274 break;
275 case XiZero:
276 theA = 1;
277 theZ = 0;
278 theS = -2;
279 break;
280 case antiNeutron:
281 theA = -1;
282 theZ = 0;
283 theS = 0;
284 break;
285 case antiLambda:
286 theA = -1;
287 theZ = 0;
288 theS = 1;
289 break;
290 case antiSigmaMinus:
291 theA = -1;
292 theZ = 1;
293 theS = 1;
294 break;
295 case antiSigmaPlus:
296 theA = -1;
297 theZ = -1;
298 theS = 1;
299 break;
300 case antiSigmaZero:
301 theA = -1;
302 theZ = 0;
303 theS = 1;
304 break;
305 case antiXiMinus:
306 theA = -1;
307 theZ = 1;
308 theS = 2;
309 break;
310 case antiXiZero:
311 theA = -1;
312 theZ = 0;
313 theS = 2;
314 break;
315 case KPlus:
316 theA = 0;
317 theZ = 1;
318 theS = 1;
319 break;
320 case KZero:
321 theA = 0;
322 theZ = 0;
323 theS = 1;
324 break;
325 case KZeroBar:
326 theA = 0;
327 theZ = 0;
328 theS = -1;
329 break;
330 case KShort:
331 theA = 0;
332 theZ = 0;
333// theS should not be defined
334 break;
335 case KLong:
336 theA = 0;
337 theZ = 0;
338// theS should not be defined
339 break;
340 case KMinus:
341 theA = 0;
342 theZ = -1;
343 theS = -1;
344 break;
345 case Composite:
346 // INCL_ERROR("Trying to set particle type to Composite! Construct a Cluster object instead" << '\n');
347 theA = 0;
348 theZ = 0;
349 theS = 0;
350 break;
351 case antiComposite:
352 theA = 0;
353 theZ = 0;
354 theS = 0;
355 break;
356 case UnknownParticle:
357 theA = 0;
358 theZ = 0;
359 theS = 0;
360 INCL_ERROR("Trying to set particle type to Unknown!" << '\n');
361 break;
362 }
363
364 if( !isResonance() && t!=Composite && t!=antiComposite )
365 setINCLMass();
366 }
void setINCLMass()
Set the mass of the Particle to its table mass.

Referenced by G4INCL::Cluster::Cluster(), G4INCL::Cluster::Cluster(), G4INCL::ClusterDecay::decay(), G4INCL::EtaNToPiNChannel::fillFinalState(), G4INCL::EtaNToPiPiNChannel::fillFinalState(), G4INCL::EtaOrOmegaNToLKChannel::fillFinalState(), G4INCL::EtaOrOmegaNToSKChannel::fillFinalState(), G4INCL::NDeltaToDeltaLKChannel::fillFinalState(), G4INCL::NKbToL2piChannel::fillFinalState(), G4INCL::NKbToLpiChannel::fillFinalState(), G4INCL::NKbToNKb2piChannel::fillFinalState(), G4INCL::NKbToNKbChannel::fillFinalState(), G4INCL::NKbToNKbpiChannel::fillFinalState(), G4INCL::NKbToS2piChannel::fillFinalState(), G4INCL::NKbToSpiChannel::fillFinalState(), G4INCL::NKToNK2piChannel::fillFinalState(), G4INCL::NKToNKChannel::fillFinalState(), G4INCL::NKToNKpiChannel::fillFinalState(), G4INCL::NNbarToNNbar2piChannel::fillFinalState(), G4INCL::NNToMissingStrangenessChannel::fillFinalState(), G4INCL::NpiToMissingStrangenessChannel::fillFinalState(), G4INCL::NSToNLChannel::fillFinalState(), G4INCL::NSToNSChannel::fillFinalState(), G4INCL::OmegaNToPiNChannel::fillFinalState(), G4INCL::OmegaNToPiPiNChannel::fillFinalState(), G4INCL::PiNToMultiPionsChannel::fillFinalState(), G4INCL::StrangeAbsorbtionChannel::fillFinalState(), Particle(), and Particle().

◆ setUncorrelatedMomentum()

void G4INCL::Particle::setUncorrelatedMomentum ( const G4double p)
inline

Set the uncorrelated momentum.

Definition at line 1152 of file G4INCLParticle.hh.

1152{ uncorrelatedMomentum = p; }

◆ swap()

void G4INCL::Particle::swap ( Particle & rhs)
inlineprotected

Helper method for the assignment operator.

Definition at line 131 of file G4INCLParticle.hh.

131 {
132 std::swap(theZ, rhs.theZ);
133 std::swap(theA, rhs.theA);
134 std::swap(theS, rhs.theS);
135 std::swap(theParticipantType, rhs.theParticipantType);
136 std::swap(theType, rhs.theType);
137 if(rhs.thePropagationEnergy == &(rhs.theFrozenEnergy))
139 else
141 std::swap(theEnergy, rhs.theEnergy);
142 std::swap(theFrozenEnergy, rhs.theFrozenEnergy);
143 if(rhs.thePropagationMomentum == &(rhs.theFrozenMomentum))
145 else
147 std::swap(theMomentum, rhs.theMomentum);
148 std::swap(theFrozenMomentum, rhs.theFrozenMomentum);
149 std::swap(thePosition, rhs.thePosition);
150 std::swap(nCollisions, rhs.nCollisions);
151 std::swap(nDecays, rhs.nDecays);
152 std::swap(nSrcPair, rhs.nSrcPair),
153 std::swap(thePotentialEnergy, rhs.thePotentialEnergy);
154 // ID intentionally not swapped
155
156#ifdef INCLXX_IN_GEANT4_MODE
157 std::swap(theParentResonancePDGCode, rhs.theParentResonancePDGCode);
158 std::swap(theParentResonanceID, rhs.theParentResonanceID);
159#endif
160
161 std::swap(theHelicity, rhs.theHelicity);
162 std::swap(emissionTime, rhs.emissionTime);
163 std::swap(outOfWell, rhs.outOfWell);
164 std::swap(theSrcPartner, rhs.theSrcPartner);
165
166 std::swap(theMass, rhs.theMass);
167 std::swap(rpCorrelated, rhs.rpCorrelated);
168 std::swap(uncorrelatedMomentum, rhs.uncorrelatedMomentum);
169
170 std::swap(theParticleBias, rhs.theParticleBias);
171 std::swap(theBiasCollisionVector, rhs.theBiasCollisionVector);
172
173 }

Referenced by operator=(), and G4INCL::Cluster::swap().

◆ thawPropagation()

void G4INCL::Particle::thawPropagation ( )
inline

Unfreeze particle propagation.

Make the particle use theMomentum and theEnergy for propagation. Call this method to restore the normal propagation if the freezePropagation() method has been called.

Definition at line 1063 of file G4INCLParticle.hh.

Member Data Documentation

◆ ID

long G4INCL::Particle::ID
protected

◆ INCLBiasVector

std::vector< G4double > G4INCL::Particle::INCLBiasVector
static

Time ordered vector of all bias applied.

/!\ Caution /!\ methods Assotiated to G4VectorCache<T> are: Push_back(…), operator[], Begin(), End(), Clear(), Size() and Pop_back()

Definition at line 1225 of file G4INCLParticle.hh.

Referenced by FillINCLBiasVector(), getBiasFromVector(), getTotalBias(), G4INCL::INCL::processEvent(), and setINCLBiasVector().

◆ nCollisions

◆ nDecays

G4int G4INCL::Particle::nDecays
protected

◆ nextBiasedCollisionID

G4ThreadLocal G4int G4INCL::Particle::nextBiasedCollisionID = 0
static

◆ nSrcPair

G4int G4INCL::Particle::nSrcPair
protected

◆ rpCorrelated

G4bool G4INCL::Particle::rpCorrelated
protected

◆ theA

◆ theEnergy

◆ theFrozenEnergy

G4double G4INCL::Particle::theFrozenEnergy
protected

◆ theFrozenMomentum

G4INCL::ThreeVector G4INCL::Particle::theFrozenMomentum
protected

◆ theMomentum

◆ theNKaon

G4int G4INCL::Particle::theNKaon
protected

The number of Kaons inside the nucleus (update during the cascade).

Definition at line 1255 of file G4INCLParticle.hh.

Referenced by G4INCL::Nucleus::emitInsideKaon(), getNumberOfKaon(), Particle(), Particle(), Particle(), Particle(), and setNumberOfKaon().

◆ theParentResonanceID

G4int G4INCL::Particle::theParentResonanceID
protected

◆ theParentResonancePDGCode

G4int G4INCL::Particle::theParentResonancePDGCode
protected

◆ theParticipantType

◆ theParticleBias

G4double G4INCL::Particle::theParticleBias
protected

◆ thePosition

◆ thePotentialEnergy

◆ thePropagationEnergy

G4double* G4INCL::Particle::thePropagationEnergy
protected

◆ thePropagationMomentum

◆ theS

◆ theType

◆ theZ

◆ uncorrelatedMomentum

G4double G4INCL::Particle::uncorrelatedMomentum
protected

The documentation for this class was generated from the following files: