CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
HydjetHadronizer.cc
Go to the documentation of this file.
1 
7 #include <iostream>
8 #include <cmath>
9 
16 
22 
23 #include "HepMC/IO_HEPEVT.h"
24 #include "HepMC/PythiaWrapper6_4.h"
25 #include "HepMC/GenEvent.h"
26 #include "HepMC/HeavyIon.h"
27 #include "HepMC/SimpleVector.h"
28 
33 
34 using namespace edm;
35 using namespace std;
36 using namespace gen;
37 
38 namespace {
39  int convertStatus(int st) {
40  if (st <= 0)
41  return 0;
42  if (st <= 10)
43  return 1;
44  if (st <= 20)
45  return 2;
46  if (st <= 30)
47  return 3;
48  else
49  return st;
50  }
51 } // namespace
52 
53 const std::vector<std::string> HydjetHadronizer::theSharedResources = {edm::SharedResourceNames::kPythia6,
55 
56 //_____________________________________________________________________
57 HydjetHadronizer::HydjetHadronizer(const ParameterSet& pset, edm::ConsumesCollector&& iC)
58  : BaseHadronizer(pset),
59  evt(nullptr),
60  pset_(pset),
61  abeamtarget_(pset.getParameter<double>("aBeamTarget")),
62  angularspecselector_(pset.getParameter<int>("angularSpectrumSelector")),
63  bfixed_(pset.getParameter<double>("bFixed")),
64  bmax_(pset.getParameter<double>("bMax")),
65  bmin_(pset.getParameter<double>("bMin")),
66  cflag_(pset.getParameter<int>("cFlag")),
67  embedding_(pset.getParameter<bool>("embeddingMode")),
68  comenergy(pset.getParameter<double>("comEnergy")),
69  doradiativeenloss_(pset.getParameter<bool>("doRadiativeEnLoss")),
70  docollisionalenloss_(pset.getParameter<bool>("doCollisionalEnLoss")),
71  fracsoftmult_(pset.getParameter<double>("fracSoftMultiplicity")),
72  hadfreeztemp_(pset.getParameter<double>("hadronFreezoutTemperature")),
73  hymode_(pset.getParameter<string>("hydjetMode")),
74  maxEventsToPrint_(pset.getUntrackedParameter<int>("maxEventsToPrint", 1)),
75  maxlongy_(pset.getParameter<double>("maxLongitudinalRapidity")),
76  maxtrany_(pset.getParameter<double>("maxTransverseRapidity")),
77  nsub_(0),
78  nhard_(0),
79  nmultiplicity_(pset.getParameter<int>("nMultiplicity")),
80  nsoft_(0),
81  nquarkflavor_(pset.getParameter<int>("qgpNumQuarkFlavor")),
82  pythiaPylistVerbosity_(pset.getUntrackedParameter<int>("pythiaPylistVerbosity", 0)),
83  qgpt0_(pset.getParameter<double>("qgpInitialTemperature")),
84  qgptau0_(pset.getParameter<double>("qgpProperTimeFormation")),
85  phi0_(0.),
86  sinphi0_(0.),
87  cosphi0_(1.),
88  rotate_(pset.getParameter<bool>("rotateEventPlane")),
89  shadowingswitch_(pset.getParameter<int>("shadowingSwitch")),
90  signn_(pset.getParameter<double>("sigmaInelNN")),
91  fVertex_(nullptr),
92  pythia6Service_(new Pythia6Service(pset)) {
93  // Default constructor
94 
95  if (pset.exists("signalVtx"))
96  signalVtx_ = pset.getUntrackedParameter<std::vector<double> >("signalVtx");
97 
98  if (signalVtx_.size() == 4) {
99  if (!fVertex_)
100  fVertex_ = new HepMC::FourVector();
101  LogDebug("EventSignalVertex") << "Setting event signal vertex "
102  << " x = " << signalVtx_.at(0) << " y = " << signalVtx_.at(1)
103  << " z= " << signalVtx_.at(2) << " t = " << signalVtx_.at(3) << endl;
104  fVertex_->set(signalVtx_.at(0), signalVtx_.at(1), signalVtx_.at(2), signalVtx_.at(3));
105  }
106 
107  // PYLIST Verbosity Level
108  // Valid PYLIST arguments are: 1, 2, 3, 5, 7, 11, 12, 13
109  pythiaPylistVerbosity_ = pset.getUntrackedParameter<int>("pythiaPylistVerbosity", 0);
110  LogDebug("PYLISTverbosity") << "Pythia PYLIST verbosity level = " << pythiaPylistVerbosity_;
111 
112  //Max number of events printed on verbosity level
113  maxEventsToPrint_ = pset.getUntrackedParameter<int>("maxEventsToPrint", 0);
114  LogDebug("Events2Print") << "Number of events to be printed = " << maxEventsToPrint_;
115 
116  if (embedding_) {
117  cflag_ = 0;
118  src_ = iC.consumes<CrossingFrame<edm::HepMCProduct> >(
119  pset.getUntrackedParameter<edm::InputTag>("backgroundLabel", edm::InputTag("mix", "generatorSmeared")));
120  }
121 
122  int cm = 1, va, vb, vc;
123  HYJVER(cm, va, vb, vc);
124  HepMC::HEPEVT_Wrapper::set_max_number_entries(4000);
125 }
126 
127 //_____________________________________________________________________
129  // destructor
130  call_pystat(1);
131  delete pythia6Service_;
132 }
133 
134 //_____________________________________________________________________
136 
137 //_____________________________________________________________________
139  // heavy ion record in the final CMSSW Event
140  double npart = hyfpar.npart;
141  int nproj = static_cast<int>(npart / 2);
142  int ntarg = static_cast<int>(npart - nproj);
143 
144  HepMC::HeavyIon* hi = new HepMC::HeavyIon(nsub_, // Ncoll_hard/N of SubEvents
145  nproj, // Npart_proj
146  ntarg, // Npart_targ
147  static_cast<int>(hyfpar.nbcol), // Ncoll
148  0, // spectator_neutrons
149  0, // spectator_protons
150  0, // N_Nwounded_collisions
151  0, // Nwounded_N_collisions
152  0, // Nwounded_Nwounded_collisions
153  hyfpar.bgen * nuclear_radius(), // impact_parameter in [fm]
154  phi0_, // event_plane_angle
155  0, //hypsi3.psi3, // eccentricity
156  hyjpar.sigin // sigma_inel_NN
157  );
158 
159  evt->set_heavy_ion(*hi);
160  delete hi;
161 }
162 
163 //___________________________________________________________________
165  // Build particle object corresponding to index in hyjets (soft+hard)
166  double x0 = hyjets.phj[0][index];
167  double y0 = hyjets.phj[1][index];
168 
169  double x = x0 * cosphi0_ - y0 * sinphi0_;
170  double y = y0 * cosphi0_ + x0 * sinphi0_;
171 
172  HepMC::GenParticle* p = new HepMC::GenParticle(HepMC::FourVector(x, // px
173  y, // py
174  hyjets.phj[2][index], // pz
175  hyjets.phj[3][index]), // E
176  hyjets.khj[1][index], // id
177  convertStatus(hyjets.khj[0][index] // status
178  ));
179 
180  p->suggest_barcode(barcode);
181  return p;
182 }
183 
184 //___________________________________________________________________
185 HepMC::GenVertex* HydjetHadronizer::build_hyjet_vertex(int i, int id) {
186  // build verteces for the hyjets stored events
187  double x0 = hyjets.vhj[0][i];
188  double y0 = hyjets.vhj[1][i];
189  double x = x0 * cosphi0_ - y0 * sinphi0_;
190  double y = y0 * cosphi0_ + x0 * sinphi0_;
191  double z = hyjets.vhj[2][i];
192  double t = hyjets.vhj[4][i];
193 
194  HepMC::GenVertex* vertex = new HepMC::GenVertex(HepMC::FourVector(x, y, z, t), id);
195  return vertex;
196 }
197 
198 //___________________________________________________________________
199 
202 
203  // generate single event
204  if (embedding_) {
205  const edm::Event& e = getEDMEvent();
206  HepMC::GenVertex* genvtx = nullptr;
207  const HepMC::GenEvent* inev = nullptr;
209  e.getByToken(src_, cf);
211  if (mix.size() < 1) {
212  throw cms::Exception("MatchVtx") << "Mixing has " << mix.size() << " sub-events, should have been at least 1"
213  << endl;
214  }
215  const HepMCProduct& bkg = mix.getObject(0);
216  if (!(bkg.isVtxGenApplied())) {
217  throw cms::Exception("MatchVtx") << "Input background does not have smeared vertex!" << endl;
218  } else {
219  inev = bkg.GetEvent();
220  }
221 
222  genvtx = inev->signal_process_vertex();
223 
224  if (!genvtx)
225  throw cms::Exception("MatchVtx") << "Input background does not have signal process vertex!" << endl;
226 
227  double aX, aY, aZ, aT;
228 
229  aX = genvtx->position().x();
230  aY = genvtx->position().y();
231  aZ = genvtx->position().z();
232  aT = genvtx->position().t();
233 
234  if (!fVertex_) {
235  fVertex_ = new HepMC::FourVector();
236  }
237  LogInfo("MatchVtx") << " setting vertex "
238  << " aX " << aX << " aY " << aY << " aZ " << aZ << " aT " << aT << endl;
239  fVertex_->set(aX, aY, aZ, aT);
240 
241  const HepMC::HeavyIon* hi = inev->heavy_ion();
242 
243  if (hi) {
244  bfixed_ = (hi->impact_parameter()) / nuclear_radius();
245  phi0_ = hi->event_plane_angle();
246  sinphi0_ = sin(phi0_);
247  cosphi0_ = cos(phi0_);
248  } else {
249  LogWarning("EventEmbedding") << "Background event does not have heavy ion record!";
250  }
251 
252  } else if (rotate_)
253  rotateEvtPlane();
254 
255  nsoft_ = 0;
256  nhard_ = 0;
257 
258  edm::LogInfo("HYDJETmode") << "##### HYDJET nhsel = " << hyjpar.nhsel;
259  edm::LogInfo("HYDJETfpart") << "##### HYDJET fpart = " << hyflow.fpart;
260  edm::LogInfo("HYDJETtf") << "##### HYDJET hadron freez-out temp, Tf = " << hyflow.Tf;
261  edm::LogInfo("HYDJETinTemp") << "##### HYDJET: QGP init temperature, T0 =" << pyqpar.T0u;
262  edm::LogInfo("HYDJETinTau") << "##### HYDJET: QGP formation time,tau0 =" << pyqpar.tau0u;
263 
264  int ntry = 0;
265  while (nsoft_ == 0 && nhard_ == 0) {
266  if (ntry > 100) {
267  edm::LogError("HydjetEmptyEvent") << "##### HYDJET: No Particles generated, Number of tries =" << ntry;
268 
269  // Throw an exception. Use the EventCorruption exception since it maps onto SkipEvent
270  // which is what we want to do here.
271 
272  std::ostringstream sstr;
273  sstr << "HydjetHadronizerProducer: No particles generated after " << ntry << " tries.\n";
274  edm::Exception except(edm::errors::EventCorruption, sstr.str());
275  throw except;
276  } else {
277  HYEVNT(bfixed_);
278  nsoft_ = hyfpar.nhyd;
279  nsub_ = hyjpar.njet;
280  nhard_ = hyfpar.npyt;
281  ++ntry;
282  }
283  }
284 
285  if (hyjpar.nhsel < 3)
286  nsub_++;
287 
288  // event information
290 
291  if (nhard_ > 0 || nsoft_ > 0)
292  get_particles(evt);
293 
294  evt->set_signal_process_id(pypars.msti[0]); // type of the process
295  evt->set_event_scale(pypars.pari[16]); // Q^2
296  add_heavy_ion_rec(evt);
297 
298  if (fVertex_) {
299  // generate new vertex & apply the shift
300 
301  // Copy the HepMC::GenEvent
302  std::unique_ptr<edm::HepMCProduct> HepMCEvt(new edm::HepMCProduct(evt));
303 
304  HepMCEvt->applyVtxGen(fVertex_);
305  evt = new HepMC::GenEvent((*HepMCEvt->GetEvent()));
306  }
307 
308  HepMC::HEPEVT_Wrapper::check_hepevt_consistency();
309  LogDebug("HEPEVT_info") << "Ev numb: " << HepMC::HEPEVT_Wrapper::event_number()
310  << " Entries number: " << HepMC::HEPEVT_Wrapper::number_entries() << " Max. entries "
311  << HepMC::HEPEVT_Wrapper::max_number_entries() << std::endl;
312 
313  event().reset(evt);
314  return true;
315 }
316 
317 //_____________________________________________________________________
319  // Hard particles. The first nhard_ lines from hyjets array.
320  // Pythia/Pyquen sub-events (sub-collisions) for a given event
321  // Return T/F if success/failure
322  // Create particles from lujet entries, assign them into vertices and
323  // put the vertices in the GenEvent, for each SubEvent
324  // The SubEvent information is kept by storing indeces of main vertices
325  // of subevents as a vector in GenHIEvent.
326 
327  LogDebug("SubEvent") << " Number of sub events " << nsub_;
328  LogDebug("Hydjet") << " Number of hard events " << hyjpar.njet;
329  LogDebug("Hydjet") << " Number of hard particles " << nhard_;
330  LogDebug("Hydjet") << " Number of soft particles " << nsoft_;
331  LogDebug("Hydjet") << " nhard_ + nsoft_ = " << nhard_ + nsoft_ << " hyjets.nhj = " << hyjets.nhj << endl;
332 
333  int ihy = 0;
334  int isub = -1;
335  int isub_l = -1;
336  int stab = 0;
337 
338  vector<HepMC::GenParticle*> primary_particle(hyjets.nhj);
339  vector<HepMC::GenParticle*> particle(hyjets.nhj);
340 
341  HepMC::GenVertex* sub_vertices = new HepMC::GenVertex(HepMC::FourVector(0, 0, 0, 0), 0); // just initialization
342 
343  // contain the last index in for each subevent
344  vector<int> index(nsub_);
345 
346  while (ihy < hyjets.nhj) {
347  isub = std::floor((hyjets.khj[2][ihy] / 50000));
348  int hjoffset = isub * 50000;
349 
350  if (isub != isub_l) {
351  sub_vertices = new HepMC::GenVertex(HepMC::FourVector(0, 0, 0, 0), isub);
352  evt->add_vertex(sub_vertices);
353  if (!evt->signal_process_vertex())
354  evt->set_signal_process_vertex(sub_vertices);
355  index[isub] = ihy - 1;
356  isub_l = isub;
357  }
358 
359  if (convertStatus(hyjets.khj[0][ihy]) == 1)
360  stab++;
361  LogDebug("Hydjet_array") << ihy << " MULTin ev.:" << hyjets.nhj << " SubEv.#" << isub << " Part #" << ihy + 1
362  << ", PDG: " << hyjets.khj[1][ihy] << " (st. " << convertStatus(hyjets.khj[0][ihy])
363  << ") mother=" << hyjets.khj[2][ihy] - (isub * 50000) + index[isub] + 1 << " ("
364  << hyjets.khj[2][ihy] << "), childs ("
365  << hyjets.khj[3][ihy] - (isub * 50000) + index[isub] + 1 << "-"
366  << hyjets.khj[4][ihy] - (isub * 50000) + index[isub] + 1 << "), vtx ("
367  << hyjets.vhj[0][ihy] << "," << hyjets.vhj[1][ihy] << "," << hyjets.vhj[2][ihy] << ") "
368  << std::endl;
369 
370  if (hyjets.khj[2][ihy] == 0) {
371  primary_particle[ihy] = build_hyjet(ihy, ihy + 1);
372  sub_vertices->add_particle_out(primary_particle[ihy]);
373  LogDebug("Hydjet_array") << " ---> " << ihy + 1 << std::endl;
374  } else {
375  particle[ihy] = build_hyjet(ihy, ihy + 1);
376  int mid = hyjets.khj[2][ihy] - hjoffset + index[isub];
377  int mid_t = mid;
378  while ((mid < ihy) && (hyjets.khj[1][mid] < 100) && (hyjets.khj[3][mid + 1] - hjoffset + index[isub] == ihy))
379  mid++;
380  if (hyjets.khj[1][mid] < 100)
381  mid = mid_t;
382 
383  HepMC::GenParticle* mother = primary_particle.at(mid);
384  HepMC::GenVertex* prods = build_hyjet_vertex(ihy, isub);
385 
386  if (!mother) {
387  mother = particle[mid];
388  primary_particle[mid] = mother;
389  }
390 
391  HepMC::GenVertex* prod_vertex = mother->end_vertex();
392  if (!prod_vertex) {
393  prod_vertex = prods;
394  prod_vertex->add_particle_in(mother);
395  LogDebug("Hydjet_array") << " <--- " << mid + 1 << std::endl;
396  evt->add_vertex(prod_vertex);
397  prods = nullptr;
398  }
399 
400  prod_vertex->add_particle_out(particle[ihy]);
401  LogDebug("Hydjet_array") << " ---" << mid + 1 << "---> " << ihy + 1 << std::endl;
402 
403  if (prods)
404  delete prods;
405  }
406  ihy++;
407  }
408  LogDebug("Hydjet_array") << " MULTin ev.:" << hyjets.nhj << ", last index: " << ihy - 1
409  << ", Sub events: " << isub + 1 << ", stable particles: " << stab << std::endl;
410 
411  return true;
412 }
413 
414 //______________________________________________________________
415 bool HydjetHadronizer::call_hyinit(double energy, double a, int ifb, double bmin, double bmax, double bfix, int nh) {
416  // initialize hydjet
417 
418  pydatr.mrpy[2] = 1;
419  HYINIT(energy, a, ifb, bmin, bmax, bfix, nh);
420  return true;
421 }
422 
423 //______________________________________________________________
425  // set hydjet options
426 
427  // hydjet running mode mode
428  // kHydroOnly --- nhsel=0 jet production off (pure HYDRO event), nhsel=0
429  // kHydroJets --- nhsle=1 jet production on, jet quenching off (HYDRO+njet*PYTHIA events)
430  // kHydroQJet --- nhsel=2 jet production & jet quenching on (HYDRO+njet*PYQUEN events)
431  // kJetsOnly --- nhsel=3 jet production on, jet quenching off, HYDRO off (njet*PYTHIA events)
432  // kQJetsOnly --- nhsel=4 jet production & jet quenching on, HYDRO off (njet*PYQUEN events)
433 
434  if (hymode_ == "kHydroOnly")
435  hyjpar.nhsel = 0;
436  else if (hymode_ == "kHydroJets")
437  hyjpar.nhsel = 1;
438  else if (hymode_ == "kHydroQJets")
439  hyjpar.nhsel = 2;
440  else if (hymode_ == "kJetsOnly")
441  hyjpar.nhsel = 3;
442  else if (hymode_ == "kQJetsOnly")
443  hyjpar.nhsel = 4;
444  else
445  hyjpar.nhsel = 2;
446 
447  // fraction of soft hydro induced multiplicity
448  hyflow.fpart = fracsoftmult_;
449 
450  // hadron freez-out temperature
451  hyflow.Tf = hadfreeztemp_;
452 
453  // maximum longitudinal collective rapidity
454  hyflow.ylfl = maxlongy_;
455 
456  // maximum transverse collective rapidity
457  hyflow.ytfl = maxtrany_;
458 
459  // shadowing on=1, off=0
460  hyjpar.ishad = shadowingswitch_;
461 
462  // set inelastic nucleon-nucleon cross section
463  hyjpar.sigin = signn_;
464 
465  // angular emitted gluon spectrum selection
466  pyqpar.ianglu = angularspecselector_;
467 
468  // number of active quark flavors in qgp
469  pyqpar.nfu = nquarkflavor_;
470 
471  // initial temperature of QGP
472  pyqpar.T0u = qgpt0_;
473 
474  // proper time of QGP formation
475  pyqpar.tau0u = qgptau0_;
476 
477  // type of medium induced partonic energy loss
479  edm::LogInfo("HydjetEnLoss") << "##### Radiative AND Collisional partonic energy loss ON ####";
480  pyqpar.ienglu = 0;
481  } else if (doradiativeenloss_) {
482  edm::LogInfo("HydjetenLoss") << "##### Only RADIATIVE partonic energy loss ON ####";
483  pyqpar.ienglu = 1;
484  } else if (docollisionalenloss_) {
485  edm::LogInfo("HydjetEnLoss") << "##### Only COLLISIONAL partonic energy loss ON ####";
486  pyqpar.ienglu = 2;
487  } else {
488  edm::LogInfo("HydjetEnLoss") << "##### Radiative AND Collisional partonic energy loss ON ####";
489  pyqpar.ienglu = 0;
490  }
491  return true;
492 }
493 
494 //_____________________________________________________________________
495 
499 
500  return true;
501 }
502 
503 //_____________________________________________________________________
504 
507  // pythia6Service_->setGeneralParams();
508 
509  // the input impact parameter (bxx_) is in [fm]; transform in [fm/RA] for hydjet usage
510  const float ra = nuclear_radius();
511  LogInfo("RAScaling") << "Nuclear radius(RA) = " << ra;
512  bmin_ /= ra;
513  bmax_ /= ra;
514  bfixed_ /= ra;
515 
516  // hydjet running options
518  // initialize hydjet
519  LogInfo("HYDJETinAction") << "##### Calling HYINIT(" << comenergy << "," << abeamtarget_ << "," << cflag_ << ","
520  << bmin_ << "," << bmax_ << "," << bfixed_ << "," << nmultiplicity_ << ") ####";
522  return true;
523 }
524 
525 bool HydjetHadronizer::declareStableParticles(const std::vector<int>& _pdg) {
526  std::vector<int> pdg = _pdg;
527  for (size_t i = 0; i < pdg.size(); i++) {
528  int pyCode = pycomp_(pdg[i]);
529  std::ostringstream pyCard;
530  pyCard << "MDCY(" << pyCode << ",1)=0";
531  std::cout << pyCard.str() << std::endl;
532  call_pygive(pyCard.str());
533  }
534  return true;
535 }
536 
537 //________________________________________________________________
539  const double pi = 3.14159265358979;
540  phi0_ = 2. * pi * gen::pyr_(nullptr) - pi;
541  sinphi0_ = sin(phi0_);
542  cosphi0_ = cos(phi0_);
543 }
544 
545 //________________________________________________________________
546 bool HydjetHadronizer::hadronize() { return false; }
547 
548 bool HydjetHadronizer::decay() { return true; }
549 
550 bool HydjetHadronizer::residualDecay() { return true; }
551 
553 
555 
556 const char* HydjetHadronizer::classname() const { return "gen::HydjetHadronizer"; }
#define HYEVNT
Definition: HydjetWrapper.h:20
double bfixed_
fixed impact param (fm); valid only if cflag_=0
T getUntrackedParameter(std::string const &, T const &) const
void doSetRandomEngine(CLHEP::HepRandomEngine *v) override
HepMC::GenEvent * evt
bool docollisionalenloss_
DEFAULT = true.
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:539
double pyr_(int *idummy)
bool hydjet_init(const edm::ParameterSet &pset)
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
bool call_pygive(const std::string &line)
bool exists(std::string const &parameterName) const
checks if a parameter exists
double npart
Definition: HydjetWrapper.h:46
void add_heavy_ion_rec(HepMC::GenEvent *evt)
#define HYJVER
Definition: HydjetWrapper.h:25
Log< level::Error, false > LogError
double comenergy
collision energy
#define hyfpar
Definition: HydjetWrapper.h:51
double v[5][pyjets_maxn]
const Double_t pi
static const std::string kPythia6
#define HYINIT
Definition: HydjetWrapper.h:13
double p[5][pyjets_maxn]
bool embedding_
Switch for embedding mode.
std::string hymode_
Hydjet running mode.
int nsoft_
multiplicity of HYDRO-induced particles in event
if(conf_.getParameter< bool >("UseStripCablingDB"))
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
const char * classname() const
unsigned int maxEventsToPrint_
Events to print if verbosity.
uint32_t nh
std::unique_ptr< HepMC::GenEvent > & event()
HepMC::GenParticle * build_hyjet(int index, int barcode)
Pythia6Service * pythia6Service_
#define hyjpar
Definition: HydjetWrapper.h:93
int nsub_
number of sub-events
std::vector< double > signalVtx_
Pset double vector to set event signal vertex.
Log< level::Info, false > LogInfo
static const std::string kFortranInstance
int nhard_
multiplicity of PYTHIA(+PYQUEN)-induced particles in event
int pycomp_(int &)
#define pypars
HepMC::GenVertex * build_hyjet_vertex(int i, int id)
const HepMC::GenEvent * GetEvent() const
Definition: HepMCProduct.h:37
T const * product() const
Definition: Handle.h:70
bool isVtxGenApplied() const
Definition: HepMCProduct.h:39
#define hyjets
Definition: HydjetWrapper.h:80
double nuclear_radius() const
bool call_hyinit(double energy, double a, int ifb, double bmin, double bmax, double bfix, int nh)
unsigned int pythiaPylistVerbosity_
pythia verbosity; def=1
void setRandomEngine(CLHEP::HepRandomEngine *v)
double abeamtarget_
beam/target atomic mass number
#define pyqpar
#define hyflow
Definition: HydjetWrapper.h:40
edm::EDGetTokenT< CrossingFrame< edm::HepMCProduct > > src_
double a
Definition: hdecay.h:119
edm::ParameterSet pset_
unsigned int shadowingswitch_
bool doradiativeenloss_
DEFAULT = true.
tuple cout
Definition: gather_cfg.py:144
int angularspecselector_
angular emitted gluon spectrum selection
Log< level::Warning, false > LogWarning
bool declareStableParticles(const std::vector< int > &)
bool rotate_
Switch to rotate event plane.
double phi0_
Event plane angle.
HepMC::FourVector * fVertex_
Event signal vertex.
bool get_particles(HepMC::GenEvent *evt)
edm::Event & getEDMEvent() const
#define LogDebug(id)