Version: SMASH-1.8
processstring.h
Go to the documentation of this file.
1 /*
2  *
3  * Copyright (c) 2017-
4  * SMASH Team
5  *
6  * GNU General Public License (GPLv3 or later)
7  *
8  */
9 
10 #ifndef SRC_INCLUDE_PROCESSSTRING_H_
11 #define SRC_INCLUDE_PROCESSSTRING_H_
12 
13 #include <memory>
14 #include <string>
15 #include <utility>
16 #include <vector>
17 
18 #include "Pythia8/Pythia.h"
19 
20 #include "constants.h"
21 #include "logging.h"
22 #include "particledata.h"
23 
24 namespace smash {
25 static constexpr int LPythia = LogArea::Pythia::id;
26 
27 // \todo Sangwook: make file (processstring) and class (StringProcess) name
28 // consistent
29 
48  private:
49  // The following 4 variables are in the center of mass frame
51  double PPosA_;
53  double PPosB_;
55  double PNegA_;
57  double PNegB_;
59  double massA_;
61  double massB_;
63  double sqrtsAB_;
65  std::array<PdgCode, 2> PDGcodes_;
67  std::array<FourVector, 2> plab_;
69  std::array<FourVector, 2> pcom_;
78  std::array<ThreeVector, 3> evecBasisAB_;
82  std::array<int, 2> NpartString_;
105  double sigma_qperp_;
136  double soft_t_form_;
156 
161 
166  ParticleList final_state_;
167 
169  std::unique_ptr<Pythia8::Pythia> pythia_parton_;
170 
172  std::unique_ptr<Pythia8::Pythia> pythia_hadron_;
173 
175  Pythia8::SigmaTotal pythia_sigmatot_;
176 
181  Pythia8::StringFlav pythia_stringflav_;
182 
187  Pythia8::Event event_intermediate_;
188 
189  public:
190  // clang-format off
191 
234  StringProcess(double string_tension, double time_formation,
235  double gluon_beta, double gluon_pmin,
236  double quark_alpha, double quark_beta,
237  double strange_supp, double diquark_supp,
238  double sigma_perp, double stringz_a_leading,
239  double stringz_b_leading, double stringz_a,
240  double stringz_b, double string_sigma_T,
241  double factor_t_form,
242  bool mass_dependent_formation_times,
243  double prob_proton_to_d_uu,
244  bool separate_fragment_baryon, double popcorn_rate);
245 
266  void common_setup_pythia(Pythia8::Pythia *pythia_in, double strange_supp,
267  double diquark_supp, double popcorn_rate,
268  double stringz_a, double stringz_b,
269  double string_sigma_T);
270 
278  const int seed_new =
280 
281  pythia_hadron_->rndm.init(seed_new);
282  logg[LPythia].debug("pythia_hadron_ : rndm is initialized with seed ", seed_new);
283  }
284 
285  // clang-format on
286 
291  Pythia8::Pythia *get_ptr_pythia_parton() { return pythia_parton_.get(); }
292 
302  std::array<double, 3> cross_sections_diffractive(int pdg_a, int pdg_b,
303  double sqrt_s) {
304  // This threshold magic is following Pythia. Todo(ryu): take care of this.
305  double sqrts_threshold = 2. * (1. + 1.0e-6);
306  /* In the case of mesons, the corresponding vector meson masses
307  * are used to evaluate the energy threshold. */
308  const int pdg_a_mod =
309  (std::abs(pdg_a) > 1000) ? pdg_a : 10 * (std::abs(pdg_a) / 10) + 3;
310  const int pdg_b_mod =
311  (std::abs(pdg_b) > 1000) ? pdg_b : 10 * (std::abs(pdg_b) / 10) + 3;
312  sqrts_threshold += pythia_hadron_->particleData.m0(pdg_a_mod) +
313  pythia_hadron_->particleData.m0(pdg_b_mod);
314  /* Constant cross-section for sub-processes below threshold equal to
315  * cross-section at the threshold. */
316  if (sqrt_s < sqrts_threshold) {
317  sqrt_s = sqrts_threshold;
318  }
319  pythia_sigmatot_.calc(pdg_a, pdg_b, sqrt_s);
320  return {pythia_sigmatot_.sigmaAX(), pythia_sigmatot_.sigmaXB(),
321  pythia_sigmatot_.sigmaXX()};
322  }
323 
338  void set_pmin_gluon_lightcone(double p_light_cone_min) {
339  pmin_gluon_lightcone_ = p_light_cone_min;
340  }
348  void set_pow_fgluon(double betapow) { pow_fgluon_beta_ = betapow; }
357  void set_pow_fquark(double alphapow, double betapow) {
358  pow_fquark_alpha_ = alphapow;
359  pow_fquark_beta_ = betapow;
360  }
365  void set_sigma_qperp_(double sigma_qperp) { sigma_qperp_ = sigma_qperp; }
370  void set_tension_string(double kappa_string) {
371  kappa_tension_string_ = kappa_string;
372  }
373 
374  // clang-format off
375 
383  void init(const ParticleList &incoming, double tcoll);
392  static void make_orthonormal_basis(ThreeVector &evec_polar,
393  std::array<ThreeVector, 3> &evec_basis);
409  const std::array<std::array<int, 2>, 2> &quarks,
410  const std::array<FourVector, 2> &pstr_com,
411  std::array<double, 2> &m_str,
412  std::array<ThreeVector, 2> &evec_str);
425  const std::array<std::array<int, 2>, 2> &quarks,
426  const std::array<FourVector, 2> &pstr_com,
427  const std::array<double, 2> &m_str,
428  const std::array<ThreeVector, 2> &evec_str,
429  bool flip_string_ends, bool separate_fragment_baryon);
430 
439  bool next_SDiff(bool is_AB_to_AX);
450  bool next_DDiff();
464  bool next_NDiffSoft();
470  bool next_NDiffHard();
481  bool next_BBbarAnn();
482 
498  static void find_excess_constituent(PdgCode &pdg_actual, PdgCode &pdg_mapped,
499  std::array<int, 5> &excess_quark,
500  std::array<int, 5> &excess_antiq);
522  void replace_constituent(Pythia8::Particle &particle,
523  std::array<int, 5> &excess_constituent);
524 
538  void find_total_number_constituent(Pythia8::Event &event_intermediate,
539  std::array<int, 5> &nquark_total,
540  std::array<int, 5> &nantiq_total);
541 
567  bool splitting_gluon_qqbar(Pythia8::Event &event_intermediate,
568  std::array<int, 5> &nquark_total, std::array<int, 5> &nantiq_total,
569  bool sign_constituent,
570  std::array<std::array<int, 5>, 2> &excess_constituent);
571 
599  void rearrange_excess(std::array<int, 5> &nquark_total,
600  std::array<std::array<int, 5>, 2> &excess_quark,
601  std::array<std::array<int, 5>, 2> &excess_antiq);
602 
640  bool restore_constituent(Pythia8::Event &event_intermediate,
641  std::array<std::array<int, 5>, 2> &excess_quark,
642  std::array<std::array<int, 5>, 2> &excess_antiq);
643 
668  void compose_string_parton(bool find_forward_string,
669  Pythia8::Event &event_intermediate,
670  Pythia8::Event &event_hadronize);
702  void compose_string_junction(bool &find_forward_string,
703  Pythia8::Event &event_intermediate,
704  Pythia8::Event &event_hadronize);
705 
727  void find_junction_leg(bool sign_color, std::vector<int> &col,
728  Pythia8::Event &event_intermediate,
729  Pythia8::Event &event_hadronize);
730 
747  int get_index_forward(bool find_forward, int np_end,
748  Pythia8::Event &event) {
749  int iforward = 1;
750  for (int ip = 2; ip < event.size() - np_end; ip++) {
751  const double y_quark_current = event[ip].y();
752  const double y_quark_forward = event[iforward].y();
753  if ((find_forward && y_quark_current > y_quark_forward) ||
754  (!find_forward && y_quark_current < y_quark_forward)) {
755  iforward = ip;
756  }
757  }
758  return iforward;
759  }
760 
766  ParticleList get_final_state() { return final_state_; }
767 
780  int append_final_state(ParticleList &intermediate_particles,
781  const FourVector &uString,
782  const ThreeVector &evecLong);
783 
792  static bool append_intermediate_list(int pdgid, FourVector momentum,
793  ParticleList &intermediate_particles) {
794  const std::string s = std::to_string(pdgid);
795  PdgCode pythia_code(s);
796  ParticleData new_particle(ParticleType::find(pythia_code));
797  new_particle.set_4momentum(momentum);
798  intermediate_particles.push_back(new_particle);
799  return true;
800  }
801 
806  static void convert_KaonLS(int &pythia_id) {
807  if (pythia_id == 310 || pythia_id == 130) {
808  pythia_id = (random::uniform_int(0, 1) == 0) ? 311 : -311;
809  }
810  }
811 
819  static void quarks_from_diquark(int diquark, int &q1, int &q2, int &deg_spin);
820 
827  static int diquark_from_quarks(int q1, int q2);
828 
836  static void make_string_ends(const PdgCode &pdgcode_in, int &idq1, int &idq2, double xi);
837 
844  static Pythia8::Vec4 set_Vec4(double energy, const ThreeVector &mom) {
845  return Pythia8::Vec4(mom.x1(), mom.x2(), mom.x3(), energy);
846  }
847 
859  static FourVector reorient(Pythia8::Particle &particle,
860  std::array<ThreeVector, 3> &evec_basis) {
861  ThreeVector three_momentum = evec_basis[0] * particle.pz() +
862  evec_basis[1] * particle.px() +
863  evec_basis[2] * particle.py();
864  return FourVector(particle.e(), three_momentum);
865  }
866 
884  int fragment_string(int idq1, int idq2, double mString,
885  ThreeVector &evecLong, bool flip_string_ends,
886  bool separate_fragment_baryon,
887  ParticleList &intermediate_particles);
888 
925  int fragment_off_hadron(bool from_forward,
926  bool separate_fragment_baryon,
927  std::array<ThreeVector, 3> &evec_basis,
928  double &ppos_string, double &pneg_string,
929  double &QTrx_string_pos, double &QTrx_string_neg,
930  double &QTry_string_pos, double &QTry_string_neg,
931  Pythia8::FlavContainer &flav_string_pos,
932  Pythia8::FlavContainer &flav_string_neg,
933  std::vector<int> &pdgid_frag,
934  std::vector<FourVector> &momentum_frag);
935 
946  int get_hadrontype_from_quark(int idq1, int idq2);
947 
956  int get_resonance_from_quark(int idq1, int idq2, double mass);
957 
976  bool separate_fragment_hadron,
977  double ppos_string, double pneg_string,
978  double mTrn_had_forward, double mTrn_had_backward,
979  double &ppos_had_forward, double &ppos_had_backward,
980  double &pneg_had_forward, double &pneg_had_backward);
981 
991  static double sample_zLund(double a, double b, double mTrn);
992 
1012  Pythia8::Event &event_fragments, std::array<ThreeVector, 3> &evec_basis,
1013  double ppos_string, double pneg_string,
1014  double QTrx_string, double QTry_string,
1015  double QTrx_add_pos, double QTry_add_pos,
1016  double QTrx_add_neg, double QTry_add_neg);
1017 
1029  Pythia8::Event &event, std::array<ThreeVector, 3> &evec_basis,
1030  double factor_yrapid, double diff_yrapid) {
1031  Pythia8::Vec4 pvec_string_now = Pythia8::Vec4(0., 0., 0., 0.);
1032  // loop over all particles in the record
1033  for (int ipyth = 1; ipyth < event.size(); ipyth++) {
1034  if (!event[ipyth].isFinal()) {
1035  continue;
1036  }
1037 
1038  FourVector p_frag = FourVector(
1039  event[ipyth].e(), event[ipyth].px(),
1040  event[ipyth].py(), event[ipyth].pz());
1041  const double E_frag = p_frag.x0();
1042  const double pl_frag = p_frag.threevec() * evec_basis[0];
1043  const double ppos_frag = (E_frag + pl_frag) * M_SQRT1_2;
1044  const double pneg_frag = (E_frag - pl_frag) * M_SQRT1_2;
1045  const double mTrn_frag = std::sqrt(2. * ppos_frag * pneg_frag);
1046  // evaluate the old momentum rapidity
1047  const double y_frag = 0.5 * std::log(ppos_frag / pneg_frag);
1048 
1049  // obtain the new momentum rapidity
1050  const double y_new_frag = factor_yrapid * y_frag + diff_yrapid;
1051  // compute the new four momentum
1052  const double E_new_frag = mTrn_frag * std::cosh(y_new_frag);
1053  const double pl_new_frag = mTrn_frag * std::sinh(y_new_frag);
1054  ThreeVector mom_new_frag =
1055  p_frag.threevec() + (pl_new_frag - pl_frag) * evec_basis[0];
1056  Pythia8::Vec4 pvec_new_frag = set_Vec4(E_new_frag, mom_new_frag);
1057  event[ipyth].p(pvec_new_frag);
1058  pvec_string_now += pvec_new_frag;
1059  }
1060  event[0].p(pvec_string_now);
1061  event[0].m(pvec_string_now.mCalc());
1062  }
1063 
1079  static void assign_all_scaling_factors(int baryon_string,
1080  ParticleList& outgoing_particles,
1081  const ThreeVector &evecLong,
1082  double suppression_factor);
1083 
1098  static std::pair<int, int> find_leading(int nq1, int nq2, ParticleList& list);
1099 
1112  static void assign_scaling_factor(int nquark, ParticleData& data,
1113  double suppression_factor);
1114 
1132  static int pdg_map_for_pythia(PdgCode &pdg);
1133 
1134 
1139  double getPPosA() { return PPosA_;}
1140 
1145  double getPNegA() { return PNegA_;}
1146 
1151  double getPPosB() { return PPosB_;}
1152 
1157  double getPnegB() { return PNegB_;}
1158 
1163  double get_massA() { return massA_;}
1164 
1169  double get_massB() { return massB_;}
1170 
1175  double get_sqrts() { return sqrtsAB_;}
1176 
1181  std::array<PdgCode, 2> get_PDGs() { return PDGcodes_;}
1182 
1187  std::array<FourVector, 2> get_plab() { return plab_;}
1188 
1193  std::array<FourVector, 2> get_pcom() { return pcom_;}
1194 
1200 
1206 
1211  double get_tcoll() { return time_collision_;}
1212 
1213  // clang-format on
1214 };
1215 
1216 } // namespace smash
1217 
1218 #endif // SRC_INCLUDE_PROCESSSTRING_H_
smash::StringProcess::get_tcoll
double get_tcoll()
Definition: processstring.h:1211
smash
Definition: action.h:24
smash::StringProcess::NpartString_
std::array< int, 2 > NpartString_
number of particles fragmented from strings
Definition: processstring.h:82
smash::StringProcess::set_pow_fquark
void set_pow_fquark(double alphapow, double betapow)
lightcone momentum fraction of quark is sampled according to probability distribution in non-diffrac...
Definition: processstring.h:357
smash::StringProcess::NpartFinal_
int NpartFinal_
total number of final state particles
Definition: processstring.h:80
smash::StringProcess::get_ucom
FourVector get_ucom()
Definition: processstring.h:1199
particledata.h
smash::StringProcess::get_final_state
ParticleList get_final_state()
a function to get the final state particle list which is called after the collision
Definition: processstring.h:766
smash::StringProcess::get_massA
double get_massA()
Definition: processstring.h:1163
smash::StringProcess::evecBasisAB_
std::array< ThreeVector, 3 > evecBasisAB_
Orthonormal basis vectors in the center of mass frame, where the 0th one is parallel to momentum of i...
Definition: processstring.h:78
smash::StringProcess::init_pythia_hadron_rndm
void init_pythia_hadron_rndm()
Set PYTHIA random seeds to be desired values.
Definition: processstring.h:277
smash::StringProcess::stringz_b_produce_
double stringz_b_produce_
parameter (StringZ:bLund) for the fragmentation function of other (produced) hadrons in soft non-diff...
Definition: processstring.h:125
smash::StringProcess::init
void init(const ParticleList &incoming, double tcoll)
initialization feed intial particles, time of collision and gamma factor of the center of mass.
Definition: processstring.cc:200
smash::StringProcess::append_final_state
int append_final_state(ParticleList &intermediate_particles, const FourVector &uString, const ThreeVector &evecLong)
compute the formation time and fill the arrays with final-state particles as described in Andersson:1...
Definition: processstring.cc:144
smash::StringProcess::stringz_b_leading_
double stringz_b_leading_
parameter (StringZ:bLund) for the fragmentation function of leading baryon in soft non-diffractive st...
Definition: processstring.h:115
smash::ParticleData
Definition: particledata.h:52
smash::StringProcess::compose_string_parton
void compose_string_parton(bool find_forward_string, Pythia8::Event &event_intermediate, Pythia8::Event &event_hadronize)
Identify a set of partons, which are connected to form a color-neutral string, from a given PYTHIA ev...
Definition: processstring.cc:1252
smash::StringProcess::get_resonance_from_quark
int get_resonance_from_quark(int idq1, int idq2, double mass)
Definition: processstring.cc:2638
smash::StringProcess::cross_sections_diffractive
std::array< double, 3 > cross_sections_diffractive(int pdg_a, int pdg_b, double sqrt_s)
Interface to pythia_sigmatot_ to compute cross-sections of A+B-> different final states Schuler:1993w...
Definition: processstring.h:302
smash::ThreeVector::x3
double x3() const
Definition: threevector.h:173
smash::StringProcess::next_NDiffHard
bool next_NDiffHard()
Hard Non-diffractive process is based on PYTHIA 8 with partonic showers and interactions.
Definition: processstring.cc:513
smash::StringProcess::massB_
double massB_
mass of incoming particle B [GeV]
Definition: processstring.h:61
smash::StringProcess::make_string_ends
static void make_string_ends(const PdgCode &pdgcode_in, int &idq1, int &idq2, double xi)
make a random selection to determine partonic contents at the string ends.
Definition: processstring.cc:1700
smash::StringProcess::plab_
std::array< FourVector, 2 > plab_
momenta of incoming particles in the lab frame [GeV]
Definition: processstring.h:67
smash::StringProcess::set_Vec4
static Pythia8::Vec4 set_Vec4(double energy, const ThreeVector &mom)
Easy setter of Pythia Vec4 from SMASH.
Definition: processstring.h:844
smash::StringProcess::pythia_sigmatot_
Pythia8::SigmaTotal pythia_sigmatot_
An object to compute cross-sections.
Definition: processstring.h:175
smash::StringProcess::find_leading
static std::pair< int, int > find_leading(int nq1, int nq2, ParticleList &list)
Find the leading string fragments.
Definition: processstring.cc:3175
smash::StringProcess::pcom_
std::array< FourVector, 2 > pcom_
momenta of incoming particles in the center of mass frame [GeV]
Definition: processstring.h:69
smash::StringProcess::set_sigma_qperp_
void set_sigma_qperp_(double sigma_qperp)
set the average amount of transverse momentum transfer sigma_qperp_.
Definition: processstring.h:365
smash::StringProcess::sqrtsAB_
double sqrtsAB_
sqrt of Mandelstam variable s of collision [GeV]
Definition: processstring.h:63
smash::StringProcess::prob_proton_to_d_uu_
double prob_proton_to_d_uu_
Probability of splitting a nucleon into the quark flavour it has only once and a diquark it has twice...
Definition: processstring.h:155
smash::StringProcess::separate_fragment_baryon_
bool separate_fragment_baryon_
Whether to use a separate fragmentation function for leading baryons.
Definition: processstring.h:158
smash::StringProcess::find_total_number_constituent
void find_total_number_constituent(Pythia8::Event &event_intermediate, std::array< int, 5 > &nquark_total, std::array< int, 5 > &nantiq_total)
Compute how many quarks and antiquarks we have in the system, and update the correspoing arrays with ...
Definition: processstring.cc:875
smash::StringProcess::get_massB
double get_massB()
Definition: processstring.h:1169
smash::StringProcess::rearrange_excess
void rearrange_excess(std::array< int, 5 > &nquark_total, std::array< std::array< int, 5 >, 2 > &excess_quark, std::array< std::array< int, 5 >, 2 > &excess_antiq)
Take total number of quarks and check if the system has enough constituents that need to be converted...
Definition: processstring.cc:1037
smash::StringProcess::replace_constituent
void replace_constituent(Pythia8::Particle &particle, std::array< int, 5 > &excess_constituent)
Convert a partonic PYTHIA particle into the desired species and update the excess of constituents.
Definition: processstring.cc:790
smash::StringProcess::time_formation_const_
double time_formation_const_
constant proper time in the case of constant formation time [fm]
Definition: processstring.h:134
smash::StringProcess::pythia_hadron_
std::unique_ptr< Pythia8::Pythia > pythia_hadron_
PYTHIA object used in fragmentation.
Definition: processstring.h:172
smash::StringProcess::pythia_parton_initialized_
bool pythia_parton_initialized_
Remembers if Pythia is initialized or not.
Definition: processstring.h:160
smash::random::uniform_int
T uniform_int(T min, T max)
Definition: random.h:100
smash::StringProcess::getPNegA
double getPNegA()
Definition: processstring.h:1145
smash::StringProcess::get_hadrontype_from_quark
int get_hadrontype_from_quark(int idq1, int idq2)
Determines hadron type from valence quark constituents.
Definition: processstring.cc:2550
smash::ParticleData::set_4momentum
void set_4momentum(const FourVector &momentum_vector)
Set the particle's 4-momentum directly.
Definition: particledata.h:148
smash::StringProcess::set_mass_and_direction_2strings
bool set_mass_and_direction_2strings(const std::array< std::array< int, 2 >, 2 > &quarks, const std::array< FourVector, 2 > &pstr_com, std::array< double, 2 > &m_str, std::array< ThreeVector, 2 > &evec_str)
Determine string masses and directions in which strings are stretched.
Definition: processstring.cc:311
smash::logg
std::array< einhard::Logger<>, std::tuple_size< LogArea::AreaTuple >::value > logg
An array that stores all pre-configured Logger objects.
Definition: logging.cc:39
smash::StringProcess::PNegA_
double PNegA_
backward lightcone momentum p^{-} of incoming particle A in CM-frame [GeV]
Definition: processstring.h:55
smash::StringProcess::pmin_gluon_lightcone_
double pmin_gluon_lightcone_
the minimum lightcone momentum scale carried by a gluon [GeV]
Definition: processstring.h:84
smash::StringProcess::soft_t_form_
double soft_t_form_
factor to be multiplied to formation times in soft strings
Definition: processstring.h:136
smash::StringProcess::getPPosB
double getPPosB()
Definition: processstring.h:1151
smash::StringProcess::PDGcodes_
std::array< PdgCode, 2 > PDGcodes_
PdgCodes of incoming particles.
Definition: processstring.h:65
smash::ParticleType::find
static const ParticleType & find(PdgCode pdgcode)
Returns the ParticleType object for the given pdgcode.
Definition: particletype.cc:105
smash::StringProcess::PPosB_
double PPosB_
forward lightcone momentum p^{+} of incoming particle B in CM-frame [GeV]
Definition: processstring.h:53
smash::StringProcess::compute_incoming_lightcone_momenta
void compute_incoming_lightcone_momenta()
compute the lightcone momenta of incoming particles where the longitudinal direction is set to be sam...
Definition: processstring.cc:1664
smash::StringProcess::shift_rapidity_event
void shift_rapidity_event(Pythia8::Event &event, std::array< ThreeVector, 3 > &evec_basis, double factor_yrapid, double diff_yrapid)
Shift the momentum rapidity of all particles in a given event record.
Definition: processstring.h:1028
smash::StringProcess
String excitation processes used in SMASH.
Definition: processstring.h:47
smash::StringProcess::quarks_from_diquark
static void quarks_from_diquark(int diquark, int &q1, int &q2, int &deg_spin)
find two quarks from a diquark.
Definition: processstring.cc:1671
smash::ThreeVector
Definition: threevector.h:31
smash::StringProcess::fragment_string
int fragment_string(int idq1, int idq2, double mString, ThreeVector &evecLong, bool flip_string_ends, bool separate_fragment_baryon, ParticleList &intermediate_particles)
perform string fragmentation to determine species and momenta of hadrons by implementing PYTHIA 8....
Definition: processstring.cc:1748
smash::StringProcess::get_plab
std::array< FourVector, 2 > get_plab()
Definition: processstring.h:1187
smash::StringProcess::diquark_from_quarks
static int diquark_from_quarks(int q1, int q2)
Construct diquark from two quarks.
Definition: processstring.cc:1687
smash::StringProcess::next_SDiff
bool next_SDiff(bool is_AB_to_AX)
Single-diffractive process is based on single pomeron exchange described in Ingelman:1984ns .
Definition: processstring.cc:229
smash::StringProcess::assign_scaling_factor
static void assign_scaling_factor(int nquark, ParticleData &data, double suppression_factor)
Assign a cross section scaling factor to the given particle.
Definition: processstring.cc:3160
smash::StringProcess::next_DDiff
bool next_DDiff()
Double-diffractive process ( A + B -> X + X ) is similar to the single-diffractive process,...
Definition: processstring.cc:381
smash::FourVector::x0
double x0() const
Definition: fourvector.h:303
smash::StringProcess::make_orthonormal_basis
static void make_orthonormal_basis(ThreeVector &evec_polar, std::array< ThreeVector, 3 > &evec_basis)
compute three orthonormal basis vectors from unit vector in the longitudinal direction
Definition: processstring.cc:1615
smash::StringProcess::set_pmin_gluon_lightcone
void set_pmin_gluon_lightcone(double p_light_cone_min)
set the minimum lightcone momentum scale carried by gluon.
Definition: processstring.h:338
smash::StringProcess::event_intermediate_
Pythia8::Event event_intermediate_
event record for intermediate partonic state in the hard string routine
Definition: processstring.h:187
smash::StringProcess::pythia_parton_
std::unique_ptr< Pythia8::Pythia > pythia_parton_
PYTHIA object used in hard string routine.
Definition: processstring.h:169
smash::StringProcess::set_pow_fgluon
void set_pow_fgluon(double betapow)
lightcone momentum fraction of gluon is sampled according to probability distribution P(x) = 1/x * (1...
Definition: processstring.h:348
smash::LPythia
static constexpr int LPythia
Definition: processstring.h:25
smash::StringProcess::find_junction_leg
void find_junction_leg(bool sign_color, std::vector< int > &col, Pythia8::Event &event_intermediate, Pythia8::Event &event_hadronize)
Identify partons, which are associated with junction legs, from a given PYTHIA event record.
Definition: processstring.cc:1442
smash::StringProcess::mass_dependent_formation_times_
bool mass_dependent_formation_times_
Whether the formation time should depend on the mass of the fragment according to Andersson:1983ia e...
Definition: processstring.h:150
smash::StringProcess::get_index_forward
int get_index_forward(bool find_forward, int np_end, Pythia8::Event &event)
Obtain index of the most forward or backward particle in a given PYTHIA event record.
Definition: processstring.h:747
smash::StringProcess::pdg_map_for_pythia
static int pdg_map_for_pythia(PdgCode &pdg)
Take pdg code and map onto particle specie which can be handled by PYTHIA.
Definition: processstring.cc:3239
smash::ThreeVector::x1
double x1() const
Definition: threevector.h:165
smash::StringProcess::get_vcom
ThreeVector get_vcom()
Definition: processstring.h:1205
smash::PdgCode
Definition: pdgcode.h:108
smash::StringProcess::PNegB_
double PNegB_
backward lightcone momentum p^{-} of incoming particle B in CM-frame [GeV]
Definition: processstring.h:57
smash::StringProcess::convert_KaonLS
static void convert_KaonLS(int &pythia_id)
convert Kaon-L or Kaon-S into K0 or Anti-K0
Definition: processstring.h:806
smash::StringProcess::assign_all_scaling_factors
static void assign_all_scaling_factors(int baryon_string, ParticleList &outgoing_particles, const ThreeVector &evecLong, double suppression_factor)
Assign a cross section scaling factor to all outgoing particles.
Definition: processstring.cc:3192
smash::StringProcess::get_ptr_pythia_parton
Pythia8::Pythia * get_ptr_pythia_parton()
Function to get the PYTHIA object for hard string routine.
Definition: processstring.h:291
smash::StringProcess::pow_fquark_beta_
double pow_fquark_beta_
parameter for the quark distribution function
Definition: processstring.h:99
smash::StringProcess::stringz_a_produce_
double stringz_a_produce_
parameter (StringZ:aLund) for the fragmentation function of other (produced) hadrons in soft non-diff...
Definition: processstring.h:120
smash::StringProcess::find_excess_constituent
static void find_excess_constituent(PdgCode &pdg_actual, PdgCode &pdg_mapped, std::array< int, 5 > &excess_quark, std::array< int, 5 > &excess_antiq)
Compare the valence quark contents of the actual and mapped hadrons and evaluate how many more consti...
Definition: processstring.cc:756
smash::StringProcess::sigma_qperp_
double sigma_qperp_
Transverse momentum spread of the excited strings.
Definition: processstring.h:105
smash::StringProcess::common_setup_pythia
void common_setup_pythia(Pythia8::Pythia *pythia_in, double strange_supp, double diquark_supp, double popcorn_rate, double stringz_a, double stringz_b, double string_sigma_T)
Common setup of PYTHIA objects for soft and hard string routines.
Definition: processstring.cc:80
smash::StringProcess::get_pcom
std::array< FourVector, 2 > get_pcom()
Definition: processstring.h:1193
smash::StringProcess::getPPosA
double getPPosA()
Definition: processstring.h:1139
smash::StringProcess::splitting_gluon_qqbar
bool splitting_gluon_qqbar(Pythia8::Event &event_intermediate, std::array< int, 5 > &nquark_total, std::array< int, 5 > &nantiq_total, bool sign_constituent, std::array< std::array< int, 5 >, 2 > &excess_constituent)
Take total number of quarks and check if the system has enough constituents that need to be converted...
Definition: processstring.cc:904
smash::StringProcess::ucomAB_
FourVector ucomAB_
velocity four vector of the center of mass in the lab frame
Definition: processstring.h:71
smash::StringProcess::reorient
static FourVector reorient(Pythia8::Particle &particle, std::array< ThreeVector, 3 > &evec_basis)
compute the four-momentum properly oriented in the lab frame.
Definition: processstring.h:859
constants.h
smash::StringProcess::sample_zLund
static double sample_zLund(double a, double b, double mTrn)
Sample lightcone momentum fraction according to the LUND fragmentation function.
Definition: processstring.cc:2829
smash::StringProcess::pythia_stringflav_
Pythia8::StringFlav pythia_stringflav_
An object for the flavor selection in string fragmentation in the case of separate fragmentation func...
Definition: processstring.h:181
logging.h
smash::StringProcess::additional_xsec_supp_
double additional_xsec_supp_
additional cross-section suppression factor to take coherence effect into account.
Definition: processstring.h:132
smash::StringProcess::make_final_state_2strings
bool make_final_state_2strings(const std::array< std::array< int, 2 >, 2 > &quarks, const std::array< FourVector, 2 > &pstr_com, const std::array< double, 2 > &m_str, const std::array< ThreeVector, 2 > &evec_str, bool flip_string_ends, bool separate_fragment_baryon)
Prepare kinematics of two strings, fragment them and append to final_state.
Definition: processstring.cc:348
smash::StringProcess::remake_kinematics_fragments
bool remake_kinematics_fragments(Pythia8::Event &event_fragments, std::array< ThreeVector, 3 > &evec_basis, double ppos_string, double pneg_string, double QTrx_string, double QTry_string, double QTrx_add_pos, double QTry_add_pos, double QTrx_add_neg, double QTry_add_neg)
Definition: processstring.cc:2861
smash::StringProcess::append_intermediate_list
static bool append_intermediate_list(int pdgid, FourVector momentum, ParticleList &intermediate_particles)
append new particle from PYTHIA to a specific particle list
Definition: processstring.h:792
smash::StringProcess::compose_string_junction
void compose_string_junction(bool &find_forward_string, Pythia8::Event &event_intermediate, Pythia8::Event &event_hadronize)
Identify a set of partons and junction(s), which are connected to form a color-neutral string,...
Definition: processstring.cc:1339
smash::StringProcess::stringz_a_leading_
double stringz_a_leading_
parameter (StringZ:aLund) for the fragmentation function of leading baryon in soft non-diffractive st...
Definition: processstring.h:110
smash::FourVector
Definition: fourvector.h:33
smash::StringProcess::time_collision_
double time_collision_
time of collision in the computational frame [fm]
Definition: processstring.h:138
smash::StringProcess::vcomAB_
ThreeVector vcomAB_
velocity three vector of the center of mass in the lab frame
Definition: processstring.h:73
smash::StringProcess::next_BBbarAnn
bool next_BBbarAnn()
Baryon-antibaryon annihilation process Based on what UrQMD Bass:1998ca , Bleicher:1999xi does,...
Definition: processstring.cc:1504
smash::StringProcess::PPosA_
double PPosA_
forward lightcone momentum p^{+} of incoming particle A in CM-frame [GeV]
Definition: processstring.h:51
smash::StringProcess::pow_fgluon_beta_
double pow_fgluon_beta_
parameter for the gluon distribution function
Definition: processstring.h:89
smash::StringProcess::make_lightcone_final_two
bool make_lightcone_final_two(bool separate_fragment_hadron, double ppos_string, double pneg_string, double mTrn_had_forward, double mTrn_had_backward, double &ppos_had_forward, double &ppos_had_backward, double &pneg_had_forward, double &pneg_had_backward)
Determines lightcone momenta of two final hadrons fragmented from a string in the same way as StringF...
Definition: processstring.cc:2764
smash::StringProcess::StringProcess
StringProcess(double string_tension, double time_formation, double gluon_beta, double gluon_pmin, double quark_alpha, double quark_beta, double strange_supp, double diquark_supp, double sigma_perp, double stringz_a_leading, double stringz_b_leading, double stringz_a, double stringz_b, double string_sigma_T, double factor_t_form, bool mass_dependent_formation_times, double prob_proton_to_d_uu, bool separate_fragment_baryon, double popcorn_rate)
Constructor, initializes pythia.
Definition: processstring.cc:21
smash::StringProcess::pow_fquark_alpha_
double pow_fquark_alpha_
parameter for the quark distribution function
Definition: processstring.h:94
smash::StringProcess::kappa_tension_string_
double kappa_tension_string_
string tension [GeV/fm]
Definition: processstring.h:127
smash::StringProcess::getPnegB
double getPnegB()
Definition: processstring.h:1157
smash::StringProcess::restore_constituent
bool restore_constituent(Pythia8::Event &event_intermediate, std::array< std::array< int, 5 >, 2 > &excess_quark, std::array< std::array< int, 5 >, 2 > &excess_antiq)
Take the intermediate partonic state from PYTHIA event with mapped hadrons and convert constituents i...
Definition: processstring.cc:1094
smash::StringProcess::fragment_off_hadron
int fragment_off_hadron(bool from_forward, bool separate_fragment_baryon, std::array< ThreeVector, 3 > &evec_basis, double &ppos_string, double &pneg_string, double &QTrx_string_pos, double &QTrx_string_neg, double &QTry_string_pos, double &QTry_string_neg, Pythia8::FlavContainer &flav_string_pos, Pythia8::FlavContainer &flav_string_neg, std::vector< int > &pdgid_frag, std::vector< FourVector > &momentum_frag)
Fragment one hadron from a given string configuration if the string mass is above threshold (given by...
Definition: processstring.cc:2219
smash::maximum_rndm_seed_in_pythia
constexpr int maximum_rndm_seed_in_pythia
The maximum value of the random seed used in PYTHIA.
Definition: constants.h:116
smash::FourVector::threevec
ThreeVector threevec() const
Definition: fourvector.h:319
smash::StringProcess::final_state_
ParticleList final_state_
final state array which must be accessed after the collision
Definition: processstring.h:166
smash::StringProcess::next_NDiffSoft
bool next_NDiffSoft()
Soft Non-diffractive process is modelled in accordance with dual-topological approach Capella:1978ig ...
Definition: processstring.cc:436
smash::StringProcess::get_sqrts
double get_sqrts()
Definition: processstring.h:1175
smash::StringProcess::get_PDGs
std::array< PdgCode, 2 > get_PDGs()
Definition: processstring.h:1181
smash::StringProcess::set_tension_string
void set_tension_string(double kappa_string)
set the string tension, which is used in append_final_state.
Definition: processstring.h:370
smash::ThreeVector::x2
double x2() const
Definition: threevector.h:169
smash::StringProcess::massA_
double massA_
mass of incoming particle A [GeV]
Definition: processstring.h:59