Version: SMASH-3.1
smash::ParticleData Class Reference

#include <particledata.h>

ParticleData contains the dynamic information of a certain particle.

Each particle has its momentum, position and other relevant physical data entry.

Definition at line 58 of file particledata.h.

Public Member Functions

 ParticleData (const ParticleType &particle_type, int unique_id=-1)
 Create a new particle with the given particle_type and optionally a specific unique_id. More...
 
int32_t id () const
 Get the id of the particle. More...
 
void set_id (int i)
 Set id of the particle. More...
 
PdgCode pdgcode () const
 Get the pdgcode of the particle. More...
 
bool is_hadron () const
 
bool is_baryon () const
 
bool is_nucleus () const
 
bool is_rho () const
 
bool is_proton () const
 
bool is_neutron () const
 
bool is_pion () const
 
double pole_mass () const
 Get the particle's pole mass ("on-shell"). More...
 
double effective_mass () const
 Get the particle's effective mass. More...
 
const ParticleTypetype () const
 Get the type of the particle. More...
 
uint32_t id_process () const
 Get the id of the last action. More...
 
HistoryData get_history () const
 Get history information. More...
 
void set_history (int ncoll, uint32_t pid, ProcessType pt, double time_of_or, const ParticleList &plist)
 Store history information. More...
 
const FourVectormomentum () const
 Get the particle's 4-momentum. More...
 
void set_4momentum (const FourVector &momentum_vector)
 Set the particle's 4-momentum directly. More...
 
void set_4momentum (double mass, const ThreeVector &mom)
 Set the momentum of the particle given its mass and momentum three-vector. More...
 
void set_4momentum (double mass, double px, double py, double pz)
 Set the momentum of the particle. More...
 
void set_3momentum (const ThreeVector &mom)
 Set the momentum of the particle without modifying the energy. More...
 
const FourVectorposition () const
 Get the particle's position in Minkowski space. More...
 
void set_4position (const FourVector &pos)
 Set the particle's 4-position directly. More...
 
void set_3position (const ThreeVector &pos)
 Set particle's 3-position. More...
 
ParticleData translated (const ThreeVector &delta) const
 Translate the particle position. More...
 
double formation_time () const
 Get the absolute formation time of the particle. More...
 
double begin_formation_time () const
 Get the absolute time, where the cross section scaling factor slowly starts increasing from the given scaling factor to 1. More...
 
void set_formation_time (const double &form_time)
 Set the absolute formation time. More...
 
void set_slow_formation_times (double begin_form_time, double form_time)
 Set the time, when the cross section scaling factor begins, and finishes to increase from the given cross section scaling factor to 1. More...
 
const double & initial_xsec_scaling_factor () const
 Get the initially assigned cross section scaling factor. More...
 
void set_cross_section_scaling_factor (const double &xsec_scal)
 Set the particle's initial cross_section_scaling_factor. More...
 
ThreeVector velocity () const
 Get the velocity 3-vector. More...
 
double inverse_gamma () const
 Get the inverse of the gamma factor from the current velocity of the particle. More...
 
void boost (const ThreeVector &v)
 Apply a full Lorentz boost of momentum and position. More...
 
void boost_momentum (const ThreeVector &v)
 Apply a Lorentz-boost to only the momentum. More...
 
void set_belongs_to (BelongsTo label)
 Setter for belongs_to label. More...
 
BelongsTo belongs_to () const
 Getter for belongs_to label. More...
 
bool operator== (const ParticleData &a) const
 Check whether two particles have the same id. More...
 
bool operator< (const ParticleData &a) const
 Check if this particle has a smaller id than another particle. More...
 
bool operator== (int id_a) const
 Check if the particle has a given id. More...
 
bool operator< (int id_a) const
 Check whether the particle's id is smaller than the given id. More...
 
 ParticleData (const ParticleType &ptype, int uid, int index)
 Construct a particle with the given type, id and index in Particles. More...
 
double xsec_scaling_factor (double delta_time=0.) const
 Return the cross section scaling factor at a given time. More...
 

Static Public Attributes

static double formation_power_ = 0.0
 Power with which the cross section scaling factor grows in time. More...
 

Private Member Functions

 ParticleData ()=default
 Default constructor. More...
 
void copy_to (ParticleData &dst) const
 Copies some information of the particle to the given particle dst. More...
 

Private Attributes

int32_t id_ = -1
 Each particle has a unique identifier. More...
 
unsigned index_ = std::numeric_limits<unsigned>::max()
 Internal index in the Particles list. More...
 
ParticleTypePtr type_
 A reference to the ParticleType object for this particle (this contains all the static information). More...
 
bool hole_ = false
 If true, the object is an entry in Particles::data_ and does not hold valid particle data. More...
 
FourVector momentum_
 momenta of the particle: x0, x1, x2, x3 as E, px, py, pz More...
 
FourVector position_
 position in space: x0, x1, x2, x3 as t, x, y, z More...
 
double formation_time_ = 0.0
 Formation time at which the particle is fully formed given as an absolute value in the computational frame. More...
 
double begin_formation_time_ = 0.0
 time when the cross section scaling factor starts to increase to 1 More...
 
double initial_xsec_scaling_factor_ = 1.0
 Initial cross section scaling factor. More...
 
HistoryData history_
 history information More...
 
BelongsTo belongs_to_ = BelongsTo::Nothing
 is it part of projectile or target nuclei? More...
 

Friends

class Particles
 

Constructor & Destructor Documentation

◆ ParticleData() [1/3]

smash::ParticleData::ParticleData ( const ParticleType particle_type,
int  unique_id = -1 
)
inlineexplicit

Create a new particle with the given particle_type and optionally a specific unique_id.

All other values are initialized to unphysical values.

Parameters
[in]particle_typeType of particle to be created
[in]unique_idid of particle to be created

Definition at line 69 of file particledata.h.

70  : id_(unique_id), type_(&particle_type) {}
int32_t id_
Each particle has a unique identifier.
Definition: particledata.h:420
ParticleTypePtr type_
A reference to the ParticleType object for this particle (this contains all the static information).
Definition: particledata.h:436

◆ ParticleData() [2/3]

smash::ParticleData::ParticleData ( const ParticleType ptype,
int  uid,
int  index 
)
inline

Construct a particle with the given type, id and index in Particles.

This constructor may only be called (directly or indirectly) from Particles. This constructor should be private, but can't be in order to support vector::emplace_back.

Parameters
[in]ptypeType of the particle to be constructed
[in]uidid of the particle to be constructed
[in]indexindex of the particle to be constructed

Definition at line 376 of file particledata.h.

377  : id_(uid), index_(index), type_(&ptype) {}
unsigned index_
Internal index in the Particles list.
Definition: particledata.h:430

◆ ParticleData() [3/3]

smash::ParticleData::ParticleData ( )
privatedefault

Default constructor.

Member Function Documentation

◆ id()

int32_t smash::ParticleData::id ( ) const
inline

Get the id of the particle.

Returns
particle id

Definition at line 76 of file particledata.h.

76 { return id_; }

◆ set_id()

void smash::ParticleData::set_id ( int  i)
inline

Set id of the particle.

Parameters
[in]iid to be assigned to the particle

Definition at line 81 of file particledata.h.

81 { id_ = i; }

◆ pdgcode()

PdgCode smash::ParticleData::pdgcode ( ) const
inline

Get the pdgcode of the particle.

Returns
pdgcode of the particle

Definition at line 87 of file particledata.h.

87 { return type_->pdgcode(); }
PdgCode pdgcode() const
Definition: particletype.h:157

◆ is_hadron()

bool smash::ParticleData::is_hadron ( ) const
inline

Returns
true if this is a baryon, antibaryon or meson.

Definition at line 91 of file particledata.h.

91 { return type_->is_hadron(); }
bool is_hadron() const
Definition: particletype.h:198

◆ is_baryon()

bool smash::ParticleData::is_baryon ( ) const
inline

Returns
whether this PDG code identifies a baryon.

Definition at line 94 of file particledata.h.

94 { return pdgcode().is_baryon(); }
PdgCode pdgcode() const
Get the pdgcode of the particle.
Definition: particledata.h:87
bool is_baryon() const
Definition: pdgcode.h:391

◆ is_nucleus()

bool smash::ParticleData::is_nucleus ( ) const
inline

Returns
true if this is a nucleus, false otherwise

Definition at line 97 of file particledata.h.

97 { return pdgcode().is_nucleus(); }
bool is_nucleus() const
Definition: pdgcode.h:364

◆ is_rho()

bool smash::ParticleData::is_rho ( ) const
inline

Returns
whether this is a rho meson (rho+/rho0/rho-)

Definition at line 100 of file particledata.h.

100 { return type_->is_rho(); }
bool is_rho() const
Definition: particletype.h:228

◆ is_proton()

bool smash::ParticleData::is_proton ( ) const
inline

Returns
whether this is a proton/anti-proton

Definition at line 103 of file particledata.h.

103 { return pdgcode().is_proton(); }
bool is_proton() const
Definition: pdgcode.h:403

◆ is_neutron()

bool smash::ParticleData::is_neutron ( ) const
inline

Returns
whether this is a neutron/anti-neutron

Definition at line 106 of file particledata.h.

106 { return pdgcode().is_neutron(); }
bool is_neutron() const
Definition: pdgcode.h:409

◆ is_pion()

bool smash::ParticleData::is_pion ( ) const
inline

Returns
whether this is a pion (pi+/pi0/pi-)

Definition at line 109 of file particledata.h.

109 { return pdgcode().is_pion(); }
bool is_pion() const
Definition: pdgcode.h:457

◆ pole_mass()

double smash::ParticleData::pole_mass ( ) const
inline

Get the particle's pole mass ("on-shell").

Returns
pole mass of the particle [GeV]

Definition at line 115 of file particledata.h.

115 { return type_->mass(); }
double mass() const
Definition: particletype.h:145

◆ effective_mass()

double smash::ParticleData::effective_mass ( ) const

Get the particle's effective mass.

Determined from the 4-momentum \(m=\sqrt{p_\mu p^\mu}\). Possibly "off-shell".

Returns
Effective mass [GeV]

Definition at line 24 of file particledata.cc.

24  {
25  const double m_pole = pole_mass();
26  if (m_pole < really_small) {
27  // prevent numerical problems with massless or very light particles
28  return m_pole;
29  } else {
30  return momentum().abs();
31  }
32 }
double abs() const
calculate the lorentz invariant absolute value
Definition: fourvector.h:464
const FourVector & momentum() const
Get the particle's 4-momentum.
Definition: particledata.h:158
double pole_mass() const
Get the particle's pole mass ("on-shell").
Definition: particledata.h:115
constexpr double really_small
Numerical error tolerance.
Definition: constants.h:37

◆ type()

const ParticleType& smash::ParticleData::type ( ) const
inline

Get the type of the particle.

Returns
ParticleType object associated to this particle.

Definition at line 128 of file particledata.h.

128 { return *type_; }

◆ id_process()

uint32_t smash::ParticleData::id_process ( ) const
inline

Get the id of the last action.

Returns
id of particle's latest collision

Definition at line 134 of file particledata.h.

134 { return history_.id_process; }
HistoryData history_
history information
Definition: particledata.h:468
int32_t id_process
id of the last action
Definition: particledata.h:34

◆ get_history()

HistoryData smash::ParticleData::get_history ( ) const
inline

Get history information.

Returns
particle history struct

Definition at line 139 of file particledata.h.

139 { return history_; }

◆ set_history()

void smash::ParticleData::set_history ( int  ncoll,
uint32_t  pid,
ProcessType  pt,
double  time_of_or,
const ParticleList &  plist 
)

Store history information.

The history contains the type of process and possibly the PdgCodes of the parent particles (plist). Note that history is not set for dileptons and photons.

Parameters
[in]ncollparticle's number of collisions
[in]pidid of the particle's latest process
[in]ptprocess type of the particle's latest process
[in]time_of_ortime of latest collision [fm]
[in]plistlist of parent particles

Definition at line 34 of file particledata.cc.

36  {
37  if (pt != ProcessType::Wall) {
39  history_.time_last_collision = time_last_coll;
40  }
41  history_.id_process = pid;
43  switch (pt) {
44  case ProcessType::Decay:
45  case ProcessType::Wall:
46  // only store one parent
47  history_.p1 = plist[0].pdgcode();
48  history_.p2 = 0x0;
49  break;
53  // Parent particles are not updated by the elastic scatterings,
54  // hypersurface crossings or failed string processes
55  break;
68  // store two parent particles
69  history_.p1 = plist[0].pdgcode();
70  history_.p2 = plist[1].pdgcode();
71  break;
78  case ProcessType::None:
79  // nullify parents
80  history_.p1 = 0x0;
81  history_.p2 = 0x0;
82  break;
83  }
84 }
@ FailedString
See here for a short description.
@ TwoToOne
See here for a short description.
@ MultiParticleThreeToTwo
See here for a short description.
@ StringSoftDoubleDiffractive
See here for a short description.
@ Bremsstrahlung
See here for a short description.
@ Thermalization
See here for a short description.
@ Freeforall
See here for a short description.
@ Decay
See here for a short description.
@ TwoToFive
See here for a short description.
@ None
See here for a short description.
@ StringSoftSingleDiffractiveXB
See here for a short description.
@ TwoToTwo
See here for a short description.
@ Wall
See here for a short description.
@ Elastic
See here for a short description.
@ TwoToFour
See here for a short description.
@ StringSoftAnnihilation
See here for a short description.
@ MultiParticleThreeMesonsToOne
See here for a short description.
@ StringSoftNonDiffractive
See here for a short description.
@ MultiParticleFourToTwo
See here for a short description.
@ StringSoftSingleDiffractiveAX
See here for a short description.
@ StringHard
See here for a short description.
@ HyperSurfaceCrossing
See here for a short description.
@ TwoToThree
See here for a short description.
@ MultiParticleFiveToTwo
See here for a short description.
double time_last_collision
Time of the last action (excluding walls), time of kinetic freeze_out for HBT analysis this time shou...
Definition: particledata.h:43
PdgCode p2
PdgCode of the second parent particles.
Definition: particledata.h:47
PdgCode p1
PdgCode of the first parent particles.
Definition: particledata.h:45
int32_t collisions_per_particle
Collision counter per particle, zero only for initially present particles.
Definition: particledata.h:32
ProcessType process_type
type of the last action
Definition: particledata.h:36

◆ momentum()

const FourVector& smash::ParticleData::momentum ( ) const
inline

Get the particle's 4-momentum.

Returns
particle's 4-momentum [GeV]

Definition at line 158 of file particledata.h.

158 { return momentum_; }
FourVector momentum_
momenta of the particle: x0, x1, x2, x3 as E, px, py, pz
Definition: particledata.h:453

◆ set_4momentum() [1/3]

void smash::ParticleData::set_4momentum ( const FourVector momentum_vector)
inline

Set the particle's 4-momentum directly.

Parameters
[in]momentum_vector4-vector \(p^\mu = (E,\mathbf{p})^T\)

Definition at line 164 of file particledata.h.

164  {
165  momentum_ = momentum_vector;
166  }

◆ set_4momentum() [2/3]

void smash::ParticleData::set_4momentum ( double  mass,
const ThreeVector mom 
)
inline

Set the momentum of the particle given its mass and momentum three-vector.

Parameters
[in]massthe mass of the particle (without E_kin contribution) [GeV]
[in]momthe three-momentum of the particle [GeV]

Definition at line 174 of file particledata.h.

174  {
175  momentum_ = FourVector(std::sqrt(mass * mass + mom * mom), mom);
176  }

◆ set_4momentum() [3/3]

void smash::ParticleData::set_4momentum ( double  mass,
double  px,
double  py,
double  pz 
)
inline

Set the momentum of the particle.

Parameters
[in]massthe mass of the particle (without E_kin contribution) [GeV]
[in]pxx-component of the momentum [GeV]
[in]pyy-component of the momentum [GeV]
[in]pzz-component of the momentum [GeV]

Definition at line 186 of file particledata.h.

186  {
187  momentum_ = FourVector(std::sqrt(mass * mass + px * px + py * py + pz * pz),
188  px, py, pz);
189  }

◆ set_3momentum()

void smash::ParticleData::set_3momentum ( const ThreeVector mom)
inline

Set the momentum of the particle without modifying the energy.

WARNING: Mass gets modified.

Parameters
[in]mommomentum 3-vector [GeV]

Definition at line 196 of file particledata.h.

196  {
197  momentum_ = FourVector(momentum_.x0(), mom);
198  }
double x0() const
Definition: fourvector.h:313

◆ position()

const FourVector& smash::ParticleData::position ( ) const
inline

Get the particle's position in Minkowski space.

Returns
particle's position 4-vector

Definition at line 204 of file particledata.h.

204 { return position_; }
FourVector position_
position in space: x0, x1, x2, x3 as t, x, y, z
Definition: particledata.h:455

◆ set_4position()

void smash::ParticleData::set_4position ( const FourVector pos)
inline

Set the particle's 4-position directly.

Parameters
[in]posposition 4-vector

Definition at line 209 of file particledata.h.

209 { position_ = pos; }

◆ set_3position()

void smash::ParticleData::set_3position ( const ThreeVector pos)
inline

Set particle's 3-position.

The time component is not changed

Parameters
[in]posposition 3-vector

Definition at line 216 of file particledata.h.

216  {
217  position_ = FourVector(position_.x0(), pos);
218  }

◆ translated()

ParticleData smash::ParticleData::translated ( const ThreeVector delta) const
inline

Translate the particle position.

Parameters
[in]delta3-vector by which the particle is translated [fm]

Definition at line 224 of file particledata.h.

224  {
225  ParticleData p = *this;
226  p.position_[1] += delta[0];
227  p.position_[2] += delta[1];
228  p.position_[3] += delta[2];
229  return p;
230  }
ParticleData()=default
Default constructor.
constexpr int p
Proton.

◆ formation_time()

double smash::ParticleData::formation_time ( ) const
inline

Get the absolute formation time of the particle.

Returns
particle's formation time

Definition at line 236 of file particledata.h.

236 { return formation_time_; }
double formation_time_
Formation time at which the particle is fully formed given as an absolute value in the computational ...
Definition: particledata.h:459

◆ begin_formation_time()

double smash::ParticleData::begin_formation_time ( ) const
inline

Get the absolute time, where the cross section scaling factor slowly starts increasing from the given scaling factor to 1.

Returns
time, when scaling factor starts increasing

Definition at line 242 of file particledata.h.

242 { return begin_formation_time_; }
double begin_formation_time_
time when the cross section scaling factor starts to increase to 1
Definition: particledata.h:461

◆ set_formation_time()

void smash::ParticleData::set_formation_time ( const double &  form_time)
inline

Set the absolute formation time.

The particle's cross section scaling factor will be a Heavyside fuction of time.

Parameters
[in]form_timeabsolute formation time

Definition at line 251 of file particledata.h.

251  {
252  formation_time_ = form_time;
253  // cross section scaling factor will be a step function in time
254  begin_formation_time_ = form_time;
255  }

◆ set_slow_formation_times()

void smash::ParticleData::set_slow_formation_times ( double  begin_form_time,
double  form_time 
)
inline

Set the time, when the cross section scaling factor begins, and finishes to increase from the given cross section scaling factor to 1.

The cross section will only grow slowly, if the option is used.

Parameters
[in]begin_form_timetime when the cross section starts to increase
[in]form_timetime when the cross section reaches 1

Definition at line 265 of file particledata.h.

265  {
266  begin_formation_time_ = begin_form_time;
267  formation_time_ = form_time;
268  }

◆ initial_xsec_scaling_factor()

const double& smash::ParticleData::initial_xsec_scaling_factor ( ) const
inline

Get the initially assigned cross section scaling factor.

Depending on the config, the cross section scaling factor might change with time, while this value will not be updated.

Returns
particle's initially assigned cross section scaling factor

Definition at line 278 of file particledata.h.

278  {
280  }
double initial_xsec_scaling_factor_
Initial cross section scaling factor.
Definition: particledata.h:466

◆ set_cross_section_scaling_factor()

void smash::ParticleData::set_cross_section_scaling_factor ( const double &  xsec_scal)
inline

Set the particle's initial cross_section_scaling_factor.

All cross sections of this particle are scaled down by this factor until the formation time is over.

If the particle formation power is set to be positive, this will only be the initial scaling factor, while the actual scaling factor grows with time.

Parameters
[in]xsec_scalcross section scaling factor

Definition at line 293 of file particledata.h.

293  {
294  initial_xsec_scaling_factor_ = xsec_scal;
295  }

◆ velocity()

ThreeVector smash::ParticleData::velocity ( ) const
inline

Get the velocity 3-vector.

Returns
3-velocity of the particle

Definition at line 301 of file particledata.h.

301 { return momentum_.velocity(); }
ThreeVector velocity() const
Get the velocity (3-vector divided by zero component).
Definition: fourvector.h:333

◆ inverse_gamma()

double smash::ParticleData::inverse_gamma ( ) const
inline

Get the inverse of the gamma factor from the current velocity of the particle.

\[\frac{1}{\gamma}=\sqrt{1-v^2}\]

This functions is more efficient than calculating the gamma factor from velocity, since the velocity function must execute three divisions (for every space component of the momentum vector).

Returns
inverse gamma factor

Definition at line 315 of file particledata.h.

315  {
316  return std::sqrt(1. - momentum_.sqr3() / (momentum_.x0() * momentum_.x0()));
317  }
double sqr3() const
calculate the square of the spatial three-vector
Definition: fourvector.h:474

◆ boost()

void smash::ParticleData::boost ( const ThreeVector v)
inline

Apply a full Lorentz boost of momentum and position.

Parameters
[in]vboost 3-velocity

Definition at line 323 of file particledata.h.

323  {
326  }
FourVector lorentz_boost(const ThreeVector &v) const
Returns the FourVector boosted with velocity v.
Definition: fourvector.cc:17
void set_4momentum(const FourVector &momentum_vector)
Set the particle's 4-momentum directly.
Definition: particledata.h:164
void set_4position(const FourVector &pos)
Set the particle's 4-position directly.
Definition: particledata.h:209

◆ boost_momentum()

void smash::ParticleData::boost_momentum ( const ThreeVector v)
inline

Apply a Lorentz-boost to only the momentum.

Parameters
[in]vboost 3-veloctity

Definition at line 332 of file particledata.h.

332  {
334  }

◆ set_belongs_to()

void smash::ParticleData::set_belongs_to ( BelongsTo  label)
inline

Setter for belongs_to label.

Definition at line 337 of file particledata.h.

337 { belongs_to_ = label; }
BelongsTo belongs_to_
is it part of projectile or target nuclei?
Definition: particledata.h:470

◆ belongs_to()

BelongsTo smash::ParticleData::belongs_to ( ) const
inline

Getter for belongs_to label.

Definition at line 339 of file particledata.h.

339 { return belongs_to_; }

◆ operator==() [1/2]

bool smash::ParticleData::operator== ( const ParticleData a) const
inline

Check whether two particles have the same id.

Parameters
[in]aparticle to compare to
Returns
whether the particles have the same id

Definition at line 346 of file particledata.h.

346 { return this->id_ == a.id_; }

◆ operator<() [1/2]

bool smash::ParticleData::operator< ( const ParticleData a) const
inline

Check if this particle has a smaller id than another particle.

Parameters
[in]aparticle to compare to
Returns
whether this particle has a smaller id than other particle

Definition at line 352 of file particledata.h.

352 { return this->id_ < a.id_; }

◆ operator==() [2/2]

bool smash::ParticleData::operator== ( int  id_a) const
inline

Check if the particle has a given id.

Parameters
[in]id_aid to compare to
Returns
whether the particle has the given id

Definition at line 359 of file particledata.h.

359 { return this->id_ == id_a; }

◆ operator<() [2/2]

bool smash::ParticleData::operator< ( int  id_a) const
inline

Check whether the particle's id is smaller than the given id.

Parameters
[in]id_anumber to compare particle's id to

Definition at line 364 of file particledata.h.

364 { return this->id_ < id_a; }

◆ xsec_scaling_factor()

double smash::ParticleData::xsec_scaling_factor ( double  delta_time = 0.) const

Return the cross section scaling factor at a given time.

Parameters
[in]delta_timescaling factor at current time plus this time will be returned.
Returns
the cross section scaling factor at a specified time.

Definition at line 86 of file particledata.cc.

86  {
87  double time_of_interest = position_.x0() + delta_time;
88  // cross section scaling factor at the time_of_interest
89  double scaling_factor;
90 
91  if (formation_power_ <= 0.) {
92  // use a step function to form particles
93  if (time_of_interest < formation_time_) {
94  // particles will not be fully formed at time of interest
95  scaling_factor = initial_xsec_scaling_factor_;
96  } else {
97  // particles are fully formed at time of interest
98  scaling_factor = 1.;
99  }
100  } else {
101  // use smooth function to scale cross section (unless particles are already
102  // fully formed at desired time or will start to form later)
103  if (formation_time_ <= time_of_interest) {
104  // particles are fully formed when colliding
105  scaling_factor = 1.;
106  } else if (begin_formation_time_ >= time_of_interest) {
107  // particles will start formimg later
108  scaling_factor = initial_xsec_scaling_factor_;
109  } else {
110  // particles are in the process of formation at the given time
111  scaling_factor =
114  std::pow((time_of_interest - begin_formation_time_) /
117  }
118  }
119  return scaling_factor;
120 }
static double formation_power_
Power with which the cross section scaling factor grows in time.
Definition: particledata.h:389

◆ copy_to()

void smash::ParticleData::copy_to ( ParticleData dst) const
inlineprivate

Copies some information of the particle to the given particle dst.

Specifically it avoids to copy id_, index_, and type_.

Parameters
[in]dstparticle values are copied to

Definition at line 402 of file particledata.h.

402  {
403  dst.history_ = history_;
404  dst.momentum_ = momentum_;
405  dst.position_ = position_;
406  dst.formation_time_ = formation_time_;
407  dst.initial_xsec_scaling_factor_ = initial_xsec_scaling_factor_;
408  dst.begin_formation_time_ = begin_formation_time_;
409  dst.belongs_to_ = belongs_to_;
410  }

Friends And Related Function Documentation

◆ Particles

friend class Particles
friend

Definition at line 392 of file particledata.h.

Member Data Documentation

◆ formation_power_

double smash::ParticleData::formation_power_ = 0.0
static

Power with which the cross section scaling factor grows in time.

Definition at line 389 of file particledata.h.

◆ id_

int32_t smash::ParticleData::id_ = -1
private

Each particle has a unique identifier.

This identifier is used for identifying the particle in the output files. It is specifically not used for searching for ParticleData objects in lists of particles, though it may be used to identify two ParticleData objects as referencing the same particle. This is why the comparison operators depend only on the id_ member.

Definition at line 420 of file particledata.h.

◆ index_

unsigned smash::ParticleData::index_ = std::numeric_limits<unsigned>::max()
private

Internal index in the Particles list.

This number is used to find the Experiment-wide original of this copy.

The value is read and written from the Particles class.

See also
Particles::data_

Definition at line 430 of file particledata.h.

◆ type_

ParticleTypePtr smash::ParticleData::type_
private

A reference to the ParticleType object for this particle (this contains all the static information).

Default-initialized with an invalid index.

Definition at line 436 of file particledata.h.

◆ hole_

bool smash::ParticleData::hole_ = false
private

If true, the object is an entry in Particles::data_ and does not hold valid particle data.

Specifically iterations over Particles must skip objects with hole_ == true. All other ParticleData instances should set this member to false.

See also
Particles::data_

Definition at line 450 of file particledata.h.

◆ momentum_

FourVector smash::ParticleData::momentum_
private

momenta of the particle: x0, x1, x2, x3 as E, px, py, pz

Definition at line 453 of file particledata.h.

◆ position_

FourVector smash::ParticleData::position_
private

position in space: x0, x1, x2, x3 as t, x, y, z

Definition at line 455 of file particledata.h.

◆ formation_time_

double smash::ParticleData::formation_time_ = 0.0
private

Formation time at which the particle is fully formed given as an absolute value in the computational frame.

Definition at line 459 of file particledata.h.

◆ begin_formation_time_

double smash::ParticleData::begin_formation_time_ = 0.0
private

time when the cross section scaling factor starts to increase to 1

Definition at line 461 of file particledata.h.

◆ initial_xsec_scaling_factor_

double smash::ParticleData::initial_xsec_scaling_factor_ = 1.0
private

Initial cross section scaling factor.

1 by default, since a particle is fully formed in this case.

Definition at line 466 of file particledata.h.

◆ history_

HistoryData smash::ParticleData::history_
private

history information

Definition at line 468 of file particledata.h.

◆ belongs_to_

BelongsTo smash::ParticleData::belongs_to_ = BelongsTo::Nothing
private

is it part of projectile or target nuclei?

Definition at line 470 of file particledata.h.


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