CMS 3D CMS Logo

List of all members | Public Member Functions | Private Member Functions | Private Attributes | Static Private Attributes
gen::Hydjet2Hadronizer Class Reference

#include <Hydjet2Hadronizer.h>

Inheritance diagram for gen::Hydjet2Hadronizer:
gen::BaseHadronizer

Public Member Functions

const char * classname () const
 
bool decay ()
 
bool declareSpecialSettings (const std::vector< std::string > &)
 
bool declareStableParticles (const std::vector< int > &)
 
void finalizeEvent ()
 
bool generatePartonsAndHadronize ()
 
bool hadronize ()
 
 Hydjet2Hadronizer (const edm::ParameterSet &, edm::ConsumesCollector &&)
 
bool initializeForExternalPartons ()
 
bool initializeForInternalPartons ()
 
bool readSettings (int)
 
bool residualDecay ()
 
void statistics ()
 
 ~Hydjet2Hadronizer () override
 
- Public Member Functions inherited from gen::BaseHadronizer
 BaseHadronizer (edm::ParameterSet const &ps)
 
void cleanLHE ()
 
void generateLHE (edm::LuminosityBlock const &lumi, CLHEP::HepRandomEngine *rengine, unsigned int ncpu)
 
edm::EventgetEDMEvent () const
 
std::unique_ptr< HepMC::GenEventgetGenEvent ()
 
std::unique_ptr< GenEventInfoProductgetGenEventInfo ()
 
virtual std::unique_ptr< GenLumiInfoHeadergetGenLumiInfoHeader () const
 
GenRunInfoProductgetGenRunInfo ()
 
std::unique_ptr< lhef::LHEEventgetLHEEvent ()
 
const std::shared_ptr< lhef::LHERunInfo > & getLHERunInfo () const
 
const std::string & gridpackPath () const
 
int randomIndex () const
 
const std::string & randomInitConfigDescription () const
 
void randomizeIndex (edm::LuminosityBlock const &lumi, CLHEP::HepRandomEngine *rengine)
 
void resetEvent (std::unique_ptr< HepMC::GenEvent > event)
 
void resetEventInfo (std::unique_ptr< GenEventInfoProduct > eventInfo)
 
virtual bool select (HepMC::GenEvent *) const
 
void setEDMEvent (edm::Event &event)
 
void setLHEEvent (std::unique_ptr< lhef::LHEEvent > event)
 
void setLHERunInfo (std::unique_ptr< lhef::LHERunInfo > runInfo)
 
void setRandomEngine (CLHEP::HepRandomEngine *v)
 
std::vector< std::string > const & sharedResources () const
 
virtual ~BaseHadronizer () noexcept(false)
 

Private Member Functions

void add_heavy_ion_rec (HepMC::GenEvent *evt)
 
HepMC::GenParticle * build_hyjet2 (int index, int barcode)
 
HepMC::GenVertex * build_hyjet2_vertex (int i, int id)
 
int convertStatus (int)
 
int convertStatusForComponents (int, int)
 
void doSetRandomEngine (CLHEP::HepRandomEngine *v) override
 
std::vector< std::string > const & doSharedResources () const override
 
bool get_particles (HepMC::GenEvent *evt)
 
double nuclear_radius () const
 
void rotateEvtPlane ()
 

Private Attributes

double cosphi0_
 
bool ev = false
 
HepMC::GenEventevt
 
InitialParamsHydjet_t fParams
 
HepMC::FourVector * fVertex_
 
Hydjet2 * hj2
 
unsigned int maxEventsToPrint_
 
int nhard_
 
int nsoft_
 
int nsub_
 
double phi0_
 
edm::ParameterSet pset
 
Pythia6Servicepythia6Service_
 
unsigned int pythiaPylistVerbosity_
 
bool rotate_
 
bool separateHydjetComponents_
 
double Sigin
 
double Sigjet
 
std::vector< double > signalVtx_
 
double sinphi0_
 
edm::EDGetTokenT< CrossingFrame< edm::HepMCProduct > > src_
 

Static Private Attributes

static const std::vector< std::string > theSharedResources = {edm::SharedResourceNames::kPythia6}
 

Additional Inherited Members

- Protected Member Functions inherited from gen::BaseHadronizer
std::unique_ptr< HepMC::GenEvent > & event ()
 
std::unique_ptr< GenEventInfoProduct > & eventInfo ()
 
lhef::LHEEventlheEvent ()
 
lhef::LHERunInfolheRunInfo ()
 
GenRunInfoProductrunInfo ()
 
- Protected Attributes inherited from gen::BaseHadronizer
std::string lheFile_
 
int randomIndex_
 

Detailed Description

Definition at line 45 of file Hydjet2Hadronizer.h.

Constructor & Destructor Documentation

◆ Hydjet2Hadronizer()

Hydjet2Hadronizer::Hydjet2Hadronizer ( const edm::ParameterSet pset,
edm::ConsumesCollector &&  iC 
)

Definition at line 74 of file Hydjet2Hadronizer.cc.

References edm::ParameterSet::exists(), DeadROC_duringRun::f1, DeadROC_duringRun::f2, fParams, fVertex_, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), HLT_2022v12_cff::InputTag, LogDebug, maxEventsToPrint_, pset, pythiaPylistVerbosity_, separateHydjetComponents_, signalVtx_, and src_.

76  rotate_(pset.getParameter<bool>("rotateEventPlane")),
77  evt(nullptr),
78  nsub_(0),
79  nhard_(0),
80  nsoft_(0),
81  phi0_(0.),
82  sinphi0_(0.),
83  cosphi0_(1.),
84  fVertex_(nullptr),
86 
87 {
88  fParams.doPrintInfo = false;
89  fParams.allowEmptyEvent = false;
90  fParams.fNevnt = 0; //not used in CMSSW
91  fParams.femb = pset.getParameter<int>("embeddingMode"); //
92  fParams.fSqrtS = pset.getParameter<double>("fSqrtS"); // C.m.s. energy per nucleon pair
93  fParams.fAw = pset.getParameter<double>("fAw"); // Atomic weigth of nuclei, fAw
94  fParams.fIfb = pset.getParameter<int>(
95  "fIfb"); // Flag of type of centrality generation, fBfix (=0 is fixed by fBfix, >0 distributed [fBfmin, fBmax])
96  fParams.fBmin = pset.getParameter<double>("fBmin"); // Minimum impact parameter in units of nuclear radius, fBmin
97  fParams.fBmax = pset.getParameter<double>("fBmax"); // Maximum impact parameter in units of nuclear radius, fBmax
98  fParams.fBfix = pset.getParameter<double>("fBfix"); // Fixed impact parameter in units of nuclear radius, fBfix
99  fParams.fT = pset.getParameter<double>("fT"); // Temperature at chemical freeze-out, fT [GeV]
100  fParams.fMuB = pset.getParameter<double>("fMuB"); // Chemical baryon potential per unit charge, fMuB [GeV]
101  fParams.fMuS = pset.getParameter<double>("fMuS"); // Chemical strangeness potential per unit charge, fMuS [GeV]
102  fParams.fMuC = pset.getParameter<double>(
103  "fMuC"); // Chemical charm potential per unit charge, fMuC [GeV] (used if charm production is turned on)
104  fParams.fMuI3 = pset.getParameter<double>("fMuI3"); // Chemical isospin potential per unit charge, fMuI3 [GeV]
105  fParams.fThFO = pset.getParameter<double>("fThFO"); // Temperature at thermal freeze-out, fThFO [GeV]
106  fParams.fMu_th_pip =
107  pset.getParameter<double>("fMu_th_pip"); // Chemical potential of pi+ at thermal freeze-out, fMu_th_pip [GeV]
108  fParams.fTau = pset.getParameter<double>(
109  "fTau"); // Proper time proper at thermal freeze-out for central collisions, fTau [fm/c]
110  fParams.fSigmaTau = pset.getParameter<double>(
111  "fSigmaTau"); // Duration of emission at thermal freeze-out for central collisions, fSigmaTau [fm/c]
112  fParams.fR = pset.getParameter<double>(
113  "fR"); // Maximal transverse radius at thermal freeze-out for central collisions, fR [fm]
114  fParams.fYlmax =
115  pset.getParameter<double>("fYlmax"); // Maximal longitudinal flow rapidity at thermal freeze-out, fYlmax
116  fParams.fUmax = pset.getParameter<double>(
117  "fUmax"); // Maximal transverse flow rapidity at thermal freeze-out for central collisions, fUmax
118  fParams.frhou2 = pset.getParameter<double>("fRhou2"); //parameter to swich ON/OFF = 0) rhou2
119  fParams.frhou3 = pset.getParameter<double>("fRhou3"); //parameter to swich ON/OFF(0) rhou3
120  fParams.frhou4 = pset.getParameter<double>("fRhou4"); //parameter to swich ON/OFF(0) rhou4
121  fParams.fDelta =
122  pset.getParameter<double>("fDelta"); // Momentum azimuthal anizotropy parameter at thermal freeze-out, fDelta
123  fParams.fEpsilon =
124  pset.getParameter<double>("fEpsilon"); // Spatial azimuthal anisotropy parameter at thermal freeze-out, fEpsilon
125  fParams.fv2 = pset.getParameter<double>("fKeps2"); //parameter to swich ON/OFF(0) epsilon2 fluctuations
126  fParams.fv3 = pset.getParameter<double>("fKeps3"); //parameter to swich ON/OFF(0) epsilon3 fluctuations
127  fParams.fIfDeltaEpsilon = pset.getParameter<double>(
128  "fIfDeltaEpsilon"); // Flag to specify fDelta and fEpsilon values, fIfDeltaEpsilon (=0 user's ones, >=1 calculated)
129  fParams.fDecay =
130  pset.getParameter<int>("fDecay"); // Flag to switch on/off hadron decays, fDecay (=0 decays off, >=1 decays on)
131  fParams.fWeakDecay = pset.getParameter<double>(
132  "fWeakDecay"); // Low decay width threshold fWeakDecay[GeV]: width<fWeakDecay decay off, width>=fDecayWidth decay on; can be used to switch off weak decays
133  fParams.fEtaType = pset.getParameter<double>(
134  "fEtaType"); // Flag to choose longitudinal flow rapidity distribution, fEtaType (=0 uniform, >0 Gaussian with the dispersion Ylmax)
135  fParams.fTMuType = pset.getParameter<double>(
136  "fTMuType"); // Flag to use calculated T_ch, mu_B and mu_S as a function of fSqrtS, fTMuType (=0 user's ones, >0 calculated)
137  fParams.fCorrS = pset.getParameter<double>(
138  "fCorrS"); // Strangeness supression factor gamma_s with fCorrS value (0<fCorrS <=1, if fCorrS <= 0 then it is calculated)
139  fParams.fCharmProd = pset.getParameter<int>(
140  "fCharmProd"); // Flag to include thermal charm production, fCharmProd (=0 no charm production, >=1 charm production)
141  fParams.fCorrC = pset.getParameter<double>(
142  "fCorrC"); // Charmness enhancement factor gamma_c with fCorrC value (fCorrC >0, if fCorrC<0 then it is calculated)
143  fParams.fNhsel = pset.getParameter<int>(
144  "fNhsel"); //Flag to include jet (J)/jet quenching (JQ) and hydro (H) state production, fNhsel (0 H on & J off, 1 H/J on & JQ off, 2 H/J/HQ on, 3 J on & H/JQ off, 4 H off & J/JQ on)
145  fParams.fPyhist = pset.getParameter<int>(
146  "fPyhist"); // Flag to suppress the output of particle history from PYTHIA, fPyhist (=1 only final state particles; =0 full particle history from PYTHIA)
147  fParams.fIshad = pset.getParameter<int>(
148  "fIshad"); // Flag to switch on/off nuclear shadowing, fIshad (0 shadowing off, 1 shadowing on)
149  fParams.fPtmin =
150  pset.getParameter<double>("fPtmin"); // Minimal pt of parton-parton scattering in PYTHIA event, fPtmin [GeV/c]
151  fParams.fT0 = pset.getParameter<double>(
152  "fT0"); // Initial QGP temperature for central Pb+Pb collisions in mid-rapidity, fT0 [GeV]
153  fParams.fTau0 = pset.getParameter<double>("fTau0"); // Proper QGP formation time in fm/c, fTau0 (0.01<fTau0<10)
154  fParams.fNf = pset.getParameter<int>("fNf"); // Number of active quark flavours in QGP, fNf (0, 1, 2 or 3)
155  fParams.fIenglu = pset.getParameter<int>(
156  "fIenglu"); // Flag to fix type of partonic energy loss, fIenglu (0 radiative and collisional loss, 1 radiative loss only, 2 collisional loss only)
157  fParams.fIanglu = pset.getParameter<int>(
158  "fIanglu"); // Flag to fix type of angular distribution of in-medium emitted gluons, fIanglu (0 small-angular, 1 wide-angular, 2 collinear).
159 
160  edm::FileInPath f1("externals/hydjet2/particles.data");
161  strcpy(fParams.partDat, (f1.fullPath()).c_str());
162 
163  edm::FileInPath f2("externals/hydjet2/tabledecay.txt");
164  strcpy(fParams.tabDecay, (f2.fullPath()).c_str());
165 
166  fParams.fPythiaTune = false;
167 
168  if (pset.exists("signalVtx"))
169  signalVtx_ = pset.getUntrackedParameter<std::vector<double>>("signalVtx");
170 
171  if (signalVtx_.size() == 4) {
172  if (!fVertex_)
173  fVertex_ = new HepMC::FourVector();
174  LogDebug("EventSignalVertex") << "Setting event signal vertex "
175  << " x = " << signalVtx_.at(0) << " y = " << signalVtx_.at(1)
176  << " z= " << signalVtx_.at(2) << " t = " << signalVtx_.at(3) << endl;
177  fVertex_->set(signalVtx_.at(0), signalVtx_.at(1), signalVtx_.at(2), signalVtx_.at(3));
178  }
179 
180  // PYLIST Verbosity Level
181  // Valid PYLIST arguments are: 1, 2, 3, 5, 7, 11, 12, 13
182  pythiaPylistVerbosity_ = pset.getUntrackedParameter<int>("pythiaPylistVerbosity", 0);
183  LogDebug("PYLISTverbosity") << "Pythia PYLIST verbosity level = " << pythiaPylistVerbosity_;
184  //Max number of events printed on verbosity level
185  maxEventsToPrint_ = pset.getUntrackedParameter<int>("maxEventsToPrint", 0);
186  LogDebug("Events2Print") << "Number of events to be printed = " << maxEventsToPrint_;
187  if (fParams.femb == 1) {
188  fParams.fIfb = 0;
190  pset.getUntrackedParameter<edm::InputTag>("backgroundLabel", edm::InputTag("mix", "generatorSmeared")));
191  }
192 
193  separateHydjetComponents_ = pset.getUntrackedParameter<bool>("separateHydjetComponents", false);
194 }
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
edm::ParameterSet pset
unsigned int maxEventsToPrint_
InitialParamsHydjet_t fParams
BaseHadronizer(edm::ParameterSet const &ps)
unsigned int pythiaPylistVerbosity_
bool exists(std::string const &parameterName) const
checks if a parameter exists
HepMC::FourVector * fVertex_
std::vector< double > signalVtx_
edm::EDGetTokenT< CrossingFrame< edm::HepMCProduct > > src_
T getUntrackedParameter(std::string const &, T const &) const
Pythia6Service * pythia6Service_
HepMC::GenEvent * evt
#define LogDebug(id)

◆ ~Hydjet2Hadronizer()

Hydjet2Hadronizer::~Hydjet2Hadronizer ( )
override

Definition at line 196 of file Hydjet2Hadronizer.cc.

References pythia6Service_.

196  {
197  call_pystat(1);
198  delete pythia6Service_;
199 }
Pythia6Service * pythia6Service_

Member Function Documentation

◆ add_heavy_ion_rec()

void Hydjet2Hadronizer::add_heavy_ion_rec ( HepMC::GenEvent evt)
private

Definition at line 504 of file Hydjet2Hadronizer.cc.

References evt, hj2, nsub_, nuclear_radius(), phi0_, and Sigin.

Referenced by generatePartonsAndHadronize().

504  {
505  // heavy ion record in the final CMSSW Event
506  int nproj = static_cast<int>((hj2->GetNpart()) / 2);
507  int ntarg = static_cast<int>((hj2->GetNpart()) - nproj);
508 
509  HepMC::HeavyIon *hi = new HepMC::HeavyIon(nsub_, // Ncoll_hard/N of SubEvents
510  nproj, // Npart_proj
511  ntarg, // Npart_targ
512  hj2->GetNbcol(), // Ncoll
513  0, // spectator_neutrons
514  0, // spectator_protons
515  0, // N_Nwounded_collisions
516  0, // Nwounded_N_collisions
517  0, // Nwounded_Nwounded_collisions
518  hj2->GetBgen() * nuclear_radius(), // impact_parameter in [fm]
519  phi0_, // event_plane_angle
520  hj2->GetPsiv3(), // eccentricity <<<---- psi for v3!!!
521  Sigin // sigma_inel_NN
522  );
523 
524  evt->set_heavy_ion(*hi);
525  delete hi;
526 }
double nuclear_radius() const
Definition: EPCuts.h:4
HepMC::GenEvent * evt

◆ build_hyjet2()

HepMC::GenParticle * Hydjet2Hadronizer::build_hyjet2 ( int  index,
int  barcode 
)
private

Definition at line 466 of file Hydjet2Hadronizer.cc.

References convertStatus(), convertStatusForComponents(), cosphi0_, GenParticle::GenParticle, hj2, gen::p, multPhiCorr_741_25nsDY_cfi::px, multPhiCorr_741_25nsDY_cfi::py, and sinphi0_.

Referenced by get_particles().

466  {
467  // Build particle object corresponding to index in hyjets (soft+hard)
468  double px0 = (hj2->GetPx()).at(index);
469  double py0 = (hj2->GetPy()).at(index);
470 
471  double px = px0 * cosphi0_ - py0 * sinphi0_;
472  double py = py0 * cosphi0_ + px0 * sinphi0_;
473 
475  new HepMC::GenParticle(HepMC::FourVector(px, // px
476  py, // py
477  (hj2->GetPz()).at(index), // pz
478  (hj2->GetE()).at(index)), // E
479  (hj2->GetPdg()).at(index), // id
481  (hj2->GetType()).at(index))) // status
482  );
483 
484  p->suggest_barcode(barcode);
485  return p;
486 }
double p[5][pyjets_maxn]
int convertStatusForComponents(int, int)

◆ build_hyjet2_vertex()

HepMC::GenVertex * Hydjet2Hadronizer::build_hyjet2_vertex ( int  i,
int  id 
)
private

Definition at line 489 of file Hydjet2Hadronizer.cc.

References cosphi0_, hj2, mps_fire::i, sinphi0_, submitPVValidationJobs::t, bphysicsOniaDQM_cfi::vertex, and x.

Referenced by get_particles().

489  {
490  // build verteces for the hyjets stored events
491  double x0 = (hj2->GetX()).at(i);
492  double y0 = (hj2->GetY()).at(i);
493  double x = x0 * cosphi0_ - y0 * sinphi0_;
494  double y = y0 * cosphi0_ + x0 * sinphi0_;
495  double z = (hj2->GetZ()).at(i);
496  double t = (hj2->GetT()).at(i);
497 
498  HepMC::GenVertex *vertex = new HepMC::GenVertex(HepMC::FourVector(x, y, z, t), id);
499 
500  return vertex;
501 }

◆ classname()

const char * Hydjet2Hadronizer::classname ( ) const

Definition at line 377 of file Hydjet2Hadronizer.cc.

377 { return "gen::Hydjet2Hadronizer"; }

◆ convertStatus()

int Hydjet2Hadronizer::convertStatus ( int  st)
private

Definition at line 57 of file Hydjet2Hadronizer.cc.

Referenced by build_hyjet2(), and get_particles().

57  {
59  if (st <= 0)
60  return 0;
61  if (st <= 10)
62  return 1;
63  if (st <= 20)
64  return 2;
65  if (st <= 30)
66  return 3;
67  }
68  return st;
69 }

◆ convertStatusForComponents()

int Hydjet2Hadronizer::convertStatusForComponents ( int  sta,
int  typ 
)
private

Definition at line 43 of file Hydjet2Hadronizer.cc.

Referenced by build_hyjet2().

43  {
44  if (sta == 1 && typ == 0)
45  return 6;
46  if (sta == 1 && typ == 1)
47  return 7;
48  if (sta == 2 && typ == 0)
49  return 16;
50  if (sta == 2 && typ == 1)
51  return 17;
52 
53  else
54  return sta;
55 }

◆ decay()

bool Hydjet2Hadronizer::decay ( )

Definition at line 373 of file Hydjet2Hadronizer.cc.

373 { return true; }

◆ declareSpecialSettings()

bool gen::Hydjet2Hadronizer::declareSpecialSettings ( const std::vector< std::string > &  )
inline

Definition at line 51 of file Hydjet2Hadronizer.h.

51 { return true; }

◆ declareStableParticles()

bool Hydjet2Hadronizer::declareStableParticles ( const std::vector< int > &  _pdg)

Definition at line 360 of file Hydjet2Hadronizer.cc.

References gen::call_pygive(), gather_cfg::cout, mps_fire::i, and gen::pycomp_().

360  {
361  std::vector<int> pdg = _pdg;
362  for (size_t i = 0; i < pdg.size(); i++) {
363  int pyCode = pycomp_(pdg[i]);
364  std::ostringstream pyCard;
365  pyCard << "MDCY(" << pyCode << ",1)=0";
366  std::cout << pyCard.str() << std::endl;
367  call_pygive(pyCard.str());
368  }
369  return true;
370 }
bool call_pygive(const std::string &line)
int pycomp_(int &)

◆ doSetRandomEngine()

void Hydjet2Hadronizer::doSetRandomEngine ( CLHEP::HepRandomEngine *  v)
overrideprivatevirtual

Reimplemented from gen::BaseHadronizer.

Definition at line 202 of file Hydjet2Hadronizer.cc.

References hjRandomEngine, pythia6Service_, gen::Pythia6Service::setRandomEngine(), and gen::v.

202  {
204  hjRandomEngine = v;
205 }
double v[5][pyjets_maxn]
Pythia6Service * pythia6Service_
CLHEP::HepRandomEngine * hjRandomEngine
Interface to the HYDJET++ (Hydjet2) generator (since core v. 2.4.3), produces HepMC events...
void setRandomEngine(CLHEP::HepRandomEngine *v)

◆ doSharedResources()

std::vector<std::string> const& gen::Hydjet2Hadronizer::doSharedResources ( ) const
inlineoverrideprivatevirtual

Reimplemented from gen::BaseHadronizer.

Definition at line 72 of file Hydjet2Hadronizer.h.

References theSharedResources.

72 { return theSharedResources; }
static const std::vector< std::string > theSharedResources

◆ finalizeEvent()

void Hydjet2Hadronizer::finalizeEvent ( )

Definition at line 375 of file Hydjet2Hadronizer.cc.

375 {}

◆ generatePartonsAndHadronize()

bool Hydjet2Hadronizer::generatePartonsAndHadronize ( )

Definition at line 236 of file Hydjet2Hadronizer.cc.

References add_heavy_ion_rec(), funct::cos(), cosphi0_, MillePedeFileConverter_cfg::e, ev, gen::BaseHadronizer::event(), edm::errors::EventCorruption, evt, Exception, fParams, fVertex_, get_particles(), gen::BaseHadronizer::getEDMEvent(), edm::HepMCProduct::GetEvent(), hj2, edm::HepMCProduct::isVtxGenApplied(), LogDebug, GeneratorMix_cff::mix, eostools::move(), nhard_, nsoft_, nsub_, nuclear_radius(), phi0_, edm::Handle< T >::product(), pypars, pythia6Service_, rotate_, rotateEvtPlane(), Sigin, Sigjet, funct::sin(), sinphi0_, and src_.

236  {
237  Pythia6Service::InstanceWrapper guard(pythia6Service_);
238 
239  // generate single event
240  if (fParams.femb == 1) {
241  const edm::Event &e = getEDMEvent();
242  HepMC::GenVertex *genvtx = nullptr;
243  const HepMC::GenEvent *inev = nullptr;
245  e.getByToken(src_, cf);
247  if (mix.size() < 1) {
248  throw cms::Exception("MatchVtx") << "Mixing has " << mix.size() << " sub-events, should have been at least 1"
249  << endl;
250  }
251  const HepMCProduct &bkg = mix.getObject(0);
252  if (!(bkg.isVtxGenApplied())) {
253  throw cms::Exception("MatchVtx") << "Input background does not have smeared vertex!" << endl;
254  } else {
255  inev = bkg.GetEvent();
256  }
257 
258  genvtx = inev->signal_process_vertex();
259 
260  if (!genvtx)
261  throw cms::Exception("MatchVtx") << "Input background does not have signal process vertex!" << endl;
262 
263  double aX, aY, aZ, aT;
264 
265  aX = genvtx->position().x();
266  aY = genvtx->position().y();
267  aZ = genvtx->position().z();
268  aT = genvtx->position().t();
269 
270  if (!fVertex_) {
271  fVertex_ = new HepMC::FourVector();
272  }
273  LogInfo("MatchVtx") << " setting vertex "
274  << " aX " << aX << " aY " << aY << " aZ " << aZ << " aT " << aT << endl;
275  fVertex_->set(aX, aY, aZ, aT);
276 
277  const HepMC::HeavyIon *hi = inev->heavy_ion();
278 
279  if (hi) {
280  fParams.fBfix = (hi->impact_parameter()) / nuclear_radius();
281  phi0_ = hi->event_plane_angle();
282  sinphi0_ = sin(phi0_);
283  cosphi0_ = cos(phi0_);
284  } else {
285  LogWarning("EventEmbedding") << "Background event does not have heavy ion record!";
286  }
287 
288  } else if (rotate_)
289  rotateEvtPlane();
290 
291  nsoft_ = 0;
292  nhard_ = 0;
293 
294  // generate one HYDJET event
295  int ntry = 0;
296 
297  while (nsoft_ == 0 && nhard_ == 0) {
298  if (ntry > 100) {
299  LogError("Hydjet2EmptyEvent") << "##### HYDJET2: No Particles generated, Number of tries =" << ntry;
300  // Throw an exception. Use the EventCorruption exception since it maps onto SkipEvent
301  // which is what we want to do here.
302  std::ostringstream sstr;
303  sstr << "Hydjet2HadronizerProducer: No particles generated after " << ntry << " tries.\n";
304  edm::Exception except(edm::errors::EventCorruption, sstr.str());
305  throw except;
306  } else {
307  hj2->GenerateEvent(fParams.fBfix);
308 
309  if (hj2->IsEmpty()) {
310  continue;
311  }
312 
313  nsoft_ = hj2->GetNhyd();
314  nsub_ = hj2->GetNjet();
315  nhard_ = hj2->GetNpyt();
316 
317  //100 trys
318  ++ntry;
319  }
320  }
321 
322  if (ev == 0) {
323  Sigin = hj2->GetSigin();
324  Sigjet = hj2->GetSigjet();
325  }
326  ev = true;
327 
328  if (fParams.fNhsel < 3)
329  nsub_++;
330 
331  // event information
332  std::unique_ptr<HepMC::GenEvent> evt = std::make_unique<HepMC::GenEvent>();
333  std::unique_ptr<edm::HepMCProduct> HepMCEvt = std::make_unique<edm::HepMCProduct>();
334 
335  if (nhard_ > 0 || nsoft_ > 0)
336  get_particles(evt.get());
337 
338  evt->set_signal_process_id(pypars.msti[0]); // type of the process
339  evt->set_event_scale(pypars.pari[16]); // Q^2
340  add_heavy_ion_rec(evt.get());
341 
342  if (fVertex_) {
343  // generate new vertex & apply the shift
344  // Copy the HepMC::GenEvent
345  HepMCEvt = std::make_unique<edm::HepMCProduct>(evt.get());
346  HepMCEvt->applyVtxGen(fVertex_);
347  evt = std::make_unique<HepMC::GenEvent>(*HepMCEvt->GetEvent());
348  }
349 
350  HepMC::HEPEVT_Wrapper::check_hepevt_consistency();
351  LogDebug("HEPEVT_info") << "Ev numb: " << HepMC::HEPEVT_Wrapper::event_number()
352  << " Entries number: " << HepMC::HEPEVT_Wrapper::number_entries() << " Max. entries "
353  << HepMC::HEPEVT_Wrapper::max_number_entries() << std::endl;
354 
355  event() = std::move(evt);
356  return kTRUE;
357 }
InitialParamsHydjet_t fParams
bool isVtxGenApplied() const
Definition: HepMCProduct.h:39
void add_heavy_ion_rec(HepMC::GenEvent *evt)
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
T const * product() const
Definition: Handle.h:70
HepMC::FourVector * fVertex_
Log< level::Error, false > LogError
edm::EDGetTokenT< CrossingFrame< edm::HepMCProduct > > src_
double nuclear_radius() const
Pythia6Service * pythia6Service_
edm::Event & getEDMEvent() const
bool get_particles(HepMC::GenEvent *evt)
Definition: EPCuts.h:4
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
std::unique_ptr< HepMC::GenEvent > & event()
const HepMC::GenEvent * GetEvent() const
Definition: HepMCProduct.h:37
Log< level::Info, false > LogInfo
#define pypars
Log< level::Warning, false > LogWarning
HepMC::GenEvent * evt
def move(src, dest)
Definition: eostools.py:511
#define LogDebug(id)

◆ get_particles()

bool Hydjet2Hadronizer::get_particles ( HepMC::GenEvent evt)
private

Definition at line 388 of file Hydjet2Hadronizer.cc.

References funct::abs(), build_hyjet2(), build_hyjet2_vertex(), convertStatus(), evt, GenParticle::GenParticle, hj2, LogDebug, nhard_, nsoft_, and nsub_.

Referenced by generatePartonsAndHadronize().

388  {
389  LogDebug("Hydjet2") << " Number of sub events " << nsub_;
390  LogDebug("Hydjet2") << " Number of hard events " << hj2->GetNjet();
391  LogDebug("Hydjet2") << " Number of hard particles " << nhard_;
392  LogDebug("Hydjet2") << " Number of soft particles " << nsoft_;
393  LogDebug("Hydjet2") << " nhard_ + nsoft_ = " << nhard_ + nsoft_ << " Ntot = " << hj2->GetNtot() << endl;
394 
395  int ihy = 0;
396  int isub_l = -1;
397  int stab = 0;
398 
399  vector<HepMC::GenParticle *> particle(hj2->GetNtot());
400  HepMC::GenVertex *sub_vertices = nullptr;
401 
402  while (ihy < hj2->GetNtot()) {
403  if ((hj2->GetiJet().at(ihy)) != isub_l) {
404  sub_vertices = new HepMC::GenVertex(HepMC::FourVector(0, 0, 0, 0), hj2->GetiJet().at(ihy));
405  evt->add_vertex(sub_vertices);
406  if (!evt->signal_process_vertex())
407  evt->set_signal_process_vertex(sub_vertices);
408  isub_l = hj2->GetiJet().at(ihy);
409  }
410 
411  if ((hj2->GetFinal().at(ihy)) == 1) //convertStatus(hj2->GetPythiaStatus().at(ihy)) == 1)
412  stab++;
413  LogDebug("Hydjet2_array") << ihy << " MULTin ev.:" << hj2->GetNtot() << " SubEv.#" << hj2->GetiJet().at(ihy)
414  << " Part #" << ihy + 1 << ", PDG: " << hj2->GetPdg().at(ihy) << " (st. "
415  << convertStatus(hj2->GetPythiaStatus().at(ihy))
416  << ") mother=" << hj2->GetMotherIndex().at(ihy) + 1 << ", childs ("
417  << hj2->GetFirstDaughterIndex().at(ihy) + 1 << "-"
418  << hj2->GetLastDaughterIndex().at(ihy) + 1 << "), vtx (" << hj2->GetX().at(ihy) << ","
419  << hj2->GetY().at(ihy) << "," << hj2->GetZ().at(ihy) << ") " << std::endl;
420 
421  if ((hj2->GetMotherIndex().at(ihy)) <= 0) {
422  particle.at(ihy) = build_hyjet2(ihy, ihy + 1);
423  if (!sub_vertices)
424  LogError("Hydjet2_array") << "##### HYDJET2: Vertex not initialized!";
425  else
426  sub_vertices->add_particle_out(particle.at(ihy));
427  LogDebug("Hydjet2_array") << " ---> " << ihy + 1 << std::endl;
428  } else {
429  particle.at(ihy) = build_hyjet2(ihy, ihy + 1);
430  int mid = hj2->GetMotherIndex().at(ihy);
431 
432  while (((mid + 1) < ihy) && (std::abs(hj2->GetPdg().at(mid)) < 100) &&
433  ((hj2->GetFirstDaughterIndex().at(mid + 1)) <= ihy)) {
434  mid++;
435  LogDebug("Hydjet2_array") << "======== MID changed to " << mid
436  << " ======== PDG(mid) = " << hj2->GetPdg().at(mid) << std::endl;
437  }
438 
439  if (std::abs(hj2->GetPdg().at(mid)) < 100) {
440  mid = hj2->GetMotherIndex().at(ihy);
441  LogDebug("Hydjet2_array") << "======== MID changed BACK to " << mid
442  << " ======== PDG(mid) = " << hj2->GetPdg().at(mid) << std::endl;
443  }
444 
445  HepMC::GenParticle *mother = particle.at(mid);
446  HepMC::GenVertex *prod_vertex = mother->end_vertex();
447 
448  if (!prod_vertex) {
449  prod_vertex = build_hyjet2_vertex(ihy, (hj2->GetiJet().at(ihy)));
450  prod_vertex->add_particle_in(mother);
451  LogDebug("Hydjet2_array") << " <--- " << mid + 1 << std::endl;
452  evt->add_vertex(prod_vertex);
453  }
454  prod_vertex->add_particle_out(particle.at(ihy));
455  LogDebug("Hydjet2_array") << " ---" << mid + 1 << "---> " << ihy + 1 << std::endl;
456  }
457  ihy++;
458  }
459 
460  LogDebug("Hydjet2_array") << " MULTin ev.:" << hj2->GetNtot() << ", last index: " << ihy - 1
461  << ", stable particles: " << stab << std::endl;
462  return kTRUE;
463 }
Log< level::Error, false > LogError
HepMC::GenVertex * build_hyjet2_vertex(int i, int id)
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
HepMC::GenParticle * build_hyjet2(int index, int barcode)
HepMC::GenEvent * evt
#define LogDebug(id)

◆ hadronize()

bool Hydjet2Hadronizer::hadronize ( )

Definition at line 372 of file Hydjet2Hadronizer.cc.

372 { return false; }

◆ initializeForExternalPartons()

bool gen::Hydjet2Hadronizer::initializeForExternalPartons ( )

◆ initializeForInternalPartons()

bool Hydjet2Hadronizer::initializeForInternalPartons ( )

Definition at line 220 of file Hydjet2Hadronizer.cc.

References fParams, hj2, nuclear_radius(), and pythia6Service_.

220  {
221  Pythia6Service::InstanceWrapper guard(pythia6Service_);
222 
223  // the input impact parameter (bxx_) is in [fm]; transform in [fm/RA] for hydjet usage
224  const double ra = nuclear_radius();
225  LogInfo("Hydjet2Hadronizer|RAScaling") << "Nuclear radius(RA) = " << ra;
226  fParams.fBmin /= ra;
227  fParams.fBmax /= ra;
228  fParams.fBfix /= ra;
229 
230  hj2 = new Hydjet2(fParams);
231 
232  return kTRUE;
233 }
InitialParamsHydjet_t fParams
double nuclear_radius() const
Pythia6Service * pythia6Service_
Log< level::Info, false > LogInfo

◆ nuclear_radius()

double gen::Hydjet2Hadronizer::nuclear_radius ( ) const
inlineprivate

Definition at line 108 of file Hydjet2Hadronizer.h.

References fParams, and funct::pow().

Referenced by add_heavy_ion_rec(), generatePartonsAndHadronize(), and initializeForInternalPartons().

108  {
109  // Return the nuclear radius derived from the
110  // beam/target atomic mass number.
111 
112  return 1.15 * pow((double)fParams.fAw, 1. / 3.);
113  }
InitialParamsHydjet_t fParams
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29

◆ readSettings()

bool Hydjet2Hadronizer::readSettings ( int  )

Definition at line 208 of file Hydjet2Hadronizer.cc.

References fParams, hjRandomEngine, pythia6Service_, and gen::Pythia6Service::setGeneralParams().

208  {
209  Pythia6Service::InstanceWrapper guard(pythia6Service_);
211 
212  fParams.fSeed = hjRandomEngine->CLHEP::HepRandomEngine::getSeed();
213  LogInfo("Hydjet2Hadronizer|GenSeed") << "Seed for random number generation: "
214  << hjRandomEngine->CLHEP::HepRandomEngine::getSeed();
215 
216  return kTRUE;
217 }
InitialParamsHydjet_t fParams
Pythia6Service * pythia6Service_
CLHEP::HepRandomEngine * hjRandomEngine
Interface to the HYDJET++ (Hydjet2) generator (since core v. 2.4.3), produces HepMC events...
Log< level::Info, false > LogInfo

◆ residualDecay()

bool Hydjet2Hadronizer::residualDecay ( )

Definition at line 374 of file Hydjet2Hadronizer.cc.

374 { return true; }

◆ rotateEvtPlane()

void Hydjet2Hadronizer::rotateEvtPlane ( )
private

Definition at line 380 of file Hydjet2Hadronizer.cc.

References funct::cos(), cosphi0_, phi0_, pi, gen::pyr_(), funct::sin(), and sinphi0_.

Referenced by generatePartonsAndHadronize().

380  {
381  const double pi = 3.14159265358979;
382  phi0_ = 2. * pi * gen::pyr_(nullptr) - pi;
383  sinphi0_ = sin(phi0_);
384  cosphi0_ = cos(phi0_);
385 }
double pyr_(int *idummy)
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
const Double_t pi
Cos< T >::type cos(const T &t)
Definition: Cos.h:22

◆ statistics()

void Hydjet2Hadronizer::statistics ( )

Definition at line 376 of file Hydjet2Hadronizer.cc.

376 {}

Member Data Documentation

◆ cosphi0_

double gen::Hydjet2Hadronizer::cosphi0_
private

◆ ev

bool gen::Hydjet2Hadronizer::ev = false
private

Definition at line 83 of file Hydjet2Hadronizer.h.

Referenced by generatePartonsAndHadronize().

◆ evt

HepMC::GenEvent* gen::Hydjet2Hadronizer::evt
private

◆ fParams

InitialParamsHydjet_t gen::Hydjet2Hadronizer::fParams
private

◆ fVertex_

HepMC::FourVector* gen::Hydjet2Hadronizer::fVertex_
private

Definition at line 100 of file Hydjet2Hadronizer.h.

Referenced by generatePartonsAndHadronize(), and Hydjet2Hadronizer().

◆ hj2

Hydjet2* gen::Hydjet2Hadronizer::hj2
private

◆ maxEventsToPrint_

unsigned int gen::Hydjet2Hadronizer::maxEventsToPrint_
private

Definition at line 95 of file Hydjet2Hadronizer.h.

Referenced by Hydjet2Hadronizer().

◆ nhard_

int gen::Hydjet2Hadronizer::nhard_
private

Definition at line 88 of file Hydjet2Hadronizer.h.

Referenced by generatePartonsAndHadronize(), and get_particles().

◆ nsoft_

int gen::Hydjet2Hadronizer::nsoft_
private

Definition at line 89 of file Hydjet2Hadronizer.h.

Referenced by generatePartonsAndHadronize(), and get_particles().

◆ nsub_

int gen::Hydjet2Hadronizer::nsub_
private

◆ phi0_

double gen::Hydjet2Hadronizer::phi0_
private

◆ pset

edm::ParameterSet gen::Hydjet2Hadronizer::pset
private

Definition at line 97 of file Hydjet2Hadronizer.h.

Referenced by Hydjet2Hadronizer().

◆ pythia6Service_

Pythia6Service* gen::Hydjet2Hadronizer::pythia6Service_
private

◆ pythiaPylistVerbosity_

unsigned int gen::Hydjet2Hadronizer::pythiaPylistVerbosity_
private

Definition at line 94 of file Hydjet2Hadronizer.h.

Referenced by Hydjet2Hadronizer().

◆ rotate_

bool gen::Hydjet2Hadronizer::rotate_
private

Definition at line 85 of file Hydjet2Hadronizer.h.

Referenced by generatePartonsAndHadronize().

◆ separateHydjetComponents_

bool gen::Hydjet2Hadronizer::separateHydjetComponents_
private

Definition at line 84 of file Hydjet2Hadronizer.h.

Referenced by Hydjet2Hadronizer().

◆ Sigin

double gen::Hydjet2Hadronizer::Sigin
private

Definition at line 98 of file Hydjet2Hadronizer.h.

Referenced by add_heavy_ion_rec(), and generatePartonsAndHadronize().

◆ Sigjet

double gen::Hydjet2Hadronizer::Sigjet
private

Definition at line 98 of file Hydjet2Hadronizer.h.

Referenced by generatePartonsAndHadronize().

◆ signalVtx_

std::vector<double> gen::Hydjet2Hadronizer::signalVtx_
private

Definition at line 102 of file Hydjet2Hadronizer.h.

Referenced by Hydjet2Hadronizer().

◆ sinphi0_

double gen::Hydjet2Hadronizer::sinphi0_
private

◆ src_

edm::EDGetTokenT<CrossingFrame<edm::HepMCProduct> > gen::Hydjet2Hadronizer::src_
private

Definition at line 105 of file Hydjet2Hadronizer.h.

Referenced by generatePartonsAndHadronize(), and Hydjet2Hadronizer().

◆ theSharedResources

const std::vector< std::string > Hydjet2Hadronizer::theSharedResources = {edm::SharedResourceNames::kPythia6}
staticprivate

Definition at line 73 of file Hydjet2Hadronizer.h.

Referenced by doSharedResources().