PYTHIA  8.312
HepMC2.h
1 // HepMC2.h is a part of the PYTHIA event generator.
2 // Copyright (C) 2024 Torbjorn Sjostrand.
3 // PYTHIA is licenced under the GNU GPL v2 or later, see COPYING for details.
4 // Please respect the MCnet Guidelines, see GUIDELINES for details.
5 
6 // Author: Mikhail Kirsanov, Mikhail.Kirsanov@cern.ch
7 // Exception classes provided by James Monk, with minor changes.
8 // Header file and function definitions for the Pythia8ToHepMC class,
9 // which converts a PYTHIA event record to the standard HepMC format.
10 //
11 // Common wrapper for HepMC2/3 added by Leif Lönnblad.
12 
13 #ifndef Pythia8_HepMC2_H
14 #define Pythia8_HepMC2_H
15 
16 #ifdef Pythia8_HepMC3_H
17 #error Cannot include HepMC2.h if HepMC3.h has already been included.
18 #endif
19 
20 #include <exception>
21 #include <sstream>
22 #include <vector>
23 #include "HepMC/IO_BaseClass.h"
24 #include "HepMC/IO_GenEvent.h"
25 #include "HepMC/GenEvent.h"
26 #include "HepMC/Units.h"
27 #include "Pythia8/Pythia.h"
28 #include "Pythia8/HIInfo.h"
29 
30 namespace HepMC {
31 
32 //==========================================================================
33 
34 // Base exception for all exceptions that Pythia8ToHepMC might throw.
35 
36 class Pythia8ToHepMCException : public std::exception {
37 
38 public:
39  virtual const char* what() const throw() { return "Pythia8ToHepMCException";}
40 
41 };
42 
43 //--------------------------------------------------------------------------
44 
45 // Exception thrown when an undecayed parton is written into the record.
46 
48 
49 public:
50 
51  // Constructor and destructor.
53  iSave = i;
54  idSave = pdg_idIn;
55  std::stringstream ss;
56  ss << "Bad end vertex at " << i << " for flavour " << pdg_idIn;
57  msg = ss.str();
58  }
59  virtual ~PartonEndVertexException() throw() {}
60 
61  // Throw exception.
62  virtual const char* what() const throw() {return msg.c_str();}
63 
64  // Return info on location and flavour of bad parton.
65  int index() const {return iSave;}
66  int pdg_id() const {return idSave;}
67 
68 protected:
69 
70  std::string msg;
71  int iSave, idSave;
72 };
73 
74 //==========================================================================
75 
76 // The Pythia8ToHepMC class.
77 
78 class Pythia8ToHepMC : public IO_BaseClass {
79 
80 public:
81 
82  // Constructor and destructor.
83  Pythia8ToHepMC() : m_internal_event_number(0), m_print_inconsistency(true),
84  m_free_parton_exception(true), m_convert_gluon_to_0(false),
85  m_store_pdf(true), m_store_proc(true), m_store_xsec(true),
86  m_store_weights(true) {;}
87  virtual ~Pythia8ToHepMC() {;}
88 
89  // The recommended method to convert Pythia events into HepMC ones.
90  bool fill_next_event( Pythia8::Pythia& pythia, GenEvent& evt,
91  int ievnum = -1, bool append = false, GenParticle* rootParticle = 0,
92  int iBarcode = -1 ) {return fill_next_event( pythia.event, &evt, ievnum,
93  &pythia.info, &pythia.settings, append, rootParticle, iBarcode);}
94  bool fill_next_event( Pythia8::Pythia& pythia, GenEvent* evt,
95  int ievnum = -1, bool append = false, GenParticle* rootParticle = 0,
96  int iBarcode = -1 ) {return fill_next_event( pythia.event, evt, ievnum,
97  &pythia.info, &pythia.settings, append, rootParticle, iBarcode);}
98 
99  // Alternative method to convert Pythia events into HepMC ones.
100  bool fill_next_event( Pythia8::Event& pyev, GenEvent& evt,
101  int ievnum = -1, const Pythia8::Info* pyinfo = 0,
102  Pythia8::Settings* pyset = 0, bool append = false,
103  GenParticle* rootParticle = 0, int iBarcode = -1) {
104  return fill_next_event(pyev, &evt, ievnum, pyinfo, pyset,
105  append, rootParticle, iBarcode); }
106 
107  bool fill_next_event( Pythia8::Event& pyev, GenEvent* evt,
108  int ievnum = -1, const Pythia8::Info* pyinfo = 0,
109  Pythia8::Settings* pyset = 0, bool append = false,
110  GenParticle* rootParticle = 0, int iBarcode = -1);
111 
112  // Read out values for some switches.
113  bool print_inconsistency() const {return m_print_inconsistency;}
114  bool free_parton_exception() const {return m_free_parton_exception;}
115  bool free_parton_warnings() const {return m_free_parton_exception;}
116  bool convert_gluon_to_0() const {return m_convert_gluon_to_0;}
117  bool store_pdf() const {return m_store_pdf;}
118  bool store_proc() const {return m_store_proc;}
119  bool store_xsec() const {return m_store_xsec;}
120  bool store_weights() const {return m_store_weights;}
121 
122  // Set values for some switches.
123  void set_print_inconsistency(bool b = true) {m_print_inconsistency = b;}
124  void set_free_parton_exception(bool b = true) {m_free_parton_exception = b;}
125  void set_free_parton_warnings(bool b = true) {m_free_parton_exception = b;}
126  void set_convert_gluon_to_0(bool b = false) {m_convert_gluon_to_0 = b;}
127  void set_store_pdf(bool b = true) {m_store_pdf = b;}
128  void set_store_proc(bool b = true) {m_store_proc = b;}
129  void set_store_xsec(bool b = true) {m_store_xsec = b;}
130  void set_store_weights(bool b = true) {m_store_weights = b;}
131 
132 private:
133 
134  // Following methods are not implemented for this class.
135  virtual bool fill_next_event( GenEvent* ) { return 0; }
136  virtual void write_event( const GenEvent* ) {;}
137 
138  // Use of copy constructor is not allowed.
139  Pythia8ToHepMC( const Pythia8ToHepMC& ) : IO_BaseClass() {;}
140 
141  // Data members.
142  int m_internal_event_number;
143  bool m_print_inconsistency, m_free_parton_exception, m_convert_gluon_to_0,
144  m_store_pdf, m_store_proc, m_store_xsec, m_store_weights;
145 
146 };
147 
148 //==========================================================================
149 
150 // Main method for conversion from PYTHIA event to HepMC event.
151 // Read one event from Pythia8 and fill a new GenEvent, alternatively
152 // append to an existing GenEvent, and return T/F = success/failure.
153 
155  GenEvent* evt, int ievnum, const Pythia8::Info* pyinfo,
156  Pythia8::Settings* pyset, bool append, GenParticle* rootParticle,
157  int iBarcode) {
158 
159  // 1. Error if no event passed.
160  if (!evt) {
161  std::cout << " Pythia8ToHepMC::fill_next_event error: passed null event."
162  << std::endl;
163  return 0;
164  }
165 
166  // Update event number counter.
167  if (!append) {
168  if (ievnum >= 0) {
169  evt->set_event_number(ievnum);
170  m_internal_event_number = ievnum;
171  } else {
172  evt->set_event_number(m_internal_event_number);
173  ++m_internal_event_number;
174  }
175  }
176 
177  // Conversion factors from Pythia units GeV and mm to HepMC ones.
178  double momFac = HepMC::Units::conversion_factor(HepMC::Units::GEV,
179  evt->momentum_unit());
180  double lenFac = HepMC::Units::conversion_factor(HepMC::Units::MM,
181  evt->length_unit());
182 
183  // Set up for alternative to append to an existing event.
184  int iStart = 1;
185  int newBarcode = 0;
186  if (append) {
187  if (!rootParticle) {
188  std::cout << " Pythia8ToHepMC::fill_next_event error: passed null "
189  << "root particle in append mode." << std::endl;
190  return 0;
191  }
192  iStart = 2;
193  newBarcode = (iBarcode > -1) ? iBarcode : evt->particles_size();
194  // New vertex associated with appended particles.
195  GenVertex* prod_vtx0 = new GenVertex();
196  prod_vtx0->add_particle_in( rootParticle );
197  evt->add_vertex( prod_vtx0 );
198  }
199 
200  // 1a. If there is a HIInfo object fill info from that.
201  if ( pyinfo && pyinfo->hiInfo ) {
202  HepMC::HeavyIon ion;
203  ion.set_Ncoll_hard(pyinfo->hiInfo->nCollNDTot());
204  ion.set_Ncoll(pyinfo->hiInfo->nAbsProj() +
205  pyinfo->hiInfo->nDiffProj() +
206  pyinfo->hiInfo->nAbsTarg() +
207  pyinfo->hiInfo->nDiffTarg() -
208  pyinfo->hiInfo->nCollND() -
209  pyinfo->hiInfo->nCollDD());
210  ion.set_Npart_proj(pyinfo->hiInfo->nAbsProj() +
211  pyinfo->hiInfo->nDiffProj());
212  ion.set_Npart_targ(pyinfo->hiInfo->nAbsTarg() +
213  pyinfo->hiInfo->nDiffTarg());
214  ion.set_impact_parameter(pyinfo->hiInfo->b());
215  evt->set_heavy_ion(ion);
216  }
217 
218  // 2. Create a particle instance for each entry and fill a map, and
219  // a vector which maps from the particle index to the GenParticle address.
220  std::vector<GenParticle*> hepevt_particles( pyev.size() );
221  for (int i = iStart; i < pyev.size(); ++i) {
222 
223  // Fill the particle.
224  hepevt_particles[i] = new GenParticle(
225  FourVector( momFac * pyev[i].px(), momFac * pyev[i].py(),
226  momFac * pyev[i].pz(), momFac * pyev[i].e() ),
227  pyev[i].id(), pyev[i].statusHepMC() );
228  if (iBarcode != 0) ++newBarcode;
229  hepevt_particles[i]->suggest_barcode( (append) ? newBarcode : i );
230  hepevt_particles[i]->set_generated_mass( momFac * pyev[i].m() );
231 
232  // Colour flow uses index 1 and 2.
233  int colType = pyev[i].colType();
234  if (colType == 1 || colType == 2)
235  hepevt_particles[i]->set_flow(1, pyev[i].col());
236  if (colType == -1 || colType == 2)
237  hepevt_particles[i]->set_flow(2, pyev[i].acol());
238  }
239 
240  // Here we assume that the first two particles in the list
241  // are the incoming beam particles.
242  if (!append)
243  evt->set_beam_particles( hepevt_particles[1], hepevt_particles[2] );
244 
245  // 3. Loop over particles AGAIN, this time creating vertices.
246  // We build the production vertex for each entry in hepevt.
247  // The HEPEVT pointers are bi-directional, so gives decay vertices as well.
248  for (int i = iStart; i < pyev.size(); ++i) {
249  GenParticle* p = hepevt_particles[i];
250 
251  // 3a. Search to see if a production vertex already exists.
252  std::vector<int> mothers = pyev[i].motherList();
253  unsigned int imother = 0;
254  int mother = -1; // note that in Pythia8 there is a particle number 0!
255  if ( !mothers.empty() ) mother = mothers[imother];
256  GenVertex* prod_vtx = p->production_vertex();
257  while ( !prod_vtx && mother > 0 ) {
258  prod_vtx = (append && mother == 1) ? rootParticle->end_vertex()
259  : hepevt_particles[mother]->end_vertex();
260  if (prod_vtx) prod_vtx->add_particle_out( p );
261  mother = ( ++imother < mothers.size() ) ? mothers[imother] : -1;
262  }
263 
264  // 3b. If no suitable production vertex exists - and the particle has
265  // at least one mother or position information to store - make one.
266  FourVector prod_pos( lenFac * pyev[i].xProd(), lenFac * pyev[i].yProd(),
267  lenFac * pyev[i].zProd(), lenFac * pyev[i].tProd() );
268  if ( !prod_vtx && ( mothers.size() > 0 || prod_pos != FourVector() ) ) {
269  prod_vtx = new GenVertex();
270  prod_vtx->add_particle_out( p );
271  evt->add_vertex( prod_vtx );
272  }
273 
274  // 3c. If prod_vtx doesn't already have position specified, fill it.
275  if ( prod_vtx && prod_vtx->position() == FourVector() )
276  prod_vtx->set_position( prod_pos );
277 
278  // 3d. loop over mothers to make sure their end_vertices are consistent.
279  imother = 0;
280  mother = -1;
281  if ( !mothers.empty() ) mother = mothers[imother];
282  while ( prod_vtx && mother > 0 ) {
283 
284  // If end vertex of the mother isn't specified, do it now.
285  GenParticle* ppp = (append && mother == 1) ? rootParticle
286  : hepevt_particles[mother];
287  if ( !ppp->end_vertex() ) {
288  prod_vtx->add_particle_in( ppp );
289  } else if (ppp->end_vertex() != prod_vtx ) {
290  // Problem scenario: the mother already has a decay vertex which
291  // differs from the daughter's production vertex.
292  // Can happen with Vincia showers since antenna emissions have two
293  // parents. In that case, let vertex structure be defined by first
294  // parent, ignoring second parent.
295  if ( (pyev[i].statusAbs() == 43 || pyev[i].statusAbs() == 51
296  || pyev[i].statusAbs() == 53) && mother >= 1) break;
297  // Otherwise this means there is internal inconsistency in the
298  // HEPEVT event record. Print an error.
299  // Note: we could provide a fix by joining the two vertices with a
300  // dummy particle if the problem arises often.
301  if ( m_print_inconsistency ) std::cout
302  << " Pythia8ToHepMC::fill_next_event: inconsistent mother/daugher "
303  << "information in Pythia8 event " << std::endl
304  << "i = " << i << " mother = " << mother
305  << "\n This warning can be turned off with the "
306  << "Pythia8ToHepMC::print_inconsistency switch." << std::endl;
307  }
308 
309  // End of vertex-setting loops.
310  mother = ( ++imother < mothers.size() ) ? mothers[imother] : -1;
311  }
312  }
313 
314  // If hadronization switched on then no final coloured particles.
315  bool doHadr = (pyset == 0) ? m_free_parton_exception
316  : pyset->flag("HadronLevel:all") && pyset->flag("HadronLevel:Hadronize");
317 
318  // 4. Check for particles which come from nowhere, i.e. are without
319  // mothers or daughters. These need to be attached to a vertex, or else
320  // they will never become part of the event.
321  for (int i = iStart; i < pyev.size(); ++i) {
322  if ( !hepevt_particles[i]->end_vertex() &&
323  !hepevt_particles[i]->production_vertex() ) {
324  std::cout << " Pythia8ToHepMC::fill_next_event error: "
325  << "hanging particle " << i << std::endl;
326  GenVertex* prod_vtx = new GenVertex();
327  prod_vtx->add_particle_out( hepevt_particles[i] );
328  evt->add_vertex( prod_vtx );
329  }
330 
331  // Also check for free partons (= gluons and quarks; not diquarks?).
332  if ( doHadr && m_free_parton_exception ) {
333  int pdg_tmp = hepevt_particles[i]->pdg_id();
334  if ( (abs(pdg_tmp) <= 6 || pdg_tmp == 21)
335  && !hepevt_particles[i]->end_vertex() )
336  throw PartonEndVertexException(i, pdg_tmp);
337  }
338  }
339 
340  // Done if only appending to already existing event.
341  if (append) return true;
342 
343  // 5. Store PDF, weight, cross section and other event information.
344  // Flavours of incoming partons.
345  if (m_store_pdf && pyinfo != 0) {
346  int id1pdf = pyinfo->id1pdf();
347  int id2pdf = pyinfo->id2pdf();
348  if ( m_convert_gluon_to_0 ) {
349  if (id1pdf == 21) id1pdf = 0;
350  if (id2pdf == 21) id2pdf = 0;
351  }
352 
353  // Store PDF information.
354  evt->set_pdf_info( PdfInfo( id1pdf, id2pdf, pyinfo->x1pdf(),
355  pyinfo->x2pdf(), pyinfo->QFac(), pyinfo->pdf1(), pyinfo->pdf2() ) );
356  }
357 
358  // Store process code, scale, alpha_em, alpha_s.
359  if (m_store_proc && pyinfo != 0) {
360  evt->set_signal_process_id( pyinfo->code() );
361  evt->set_event_scale( pyinfo->QRen() );
362  if (evt->alphaQED() <= 0) evt->set_alphaQED( pyinfo->alphaEM() );
363  if (evt->alphaQCD() <= 0) evt->set_alphaQCD( pyinfo->alphaS() );
364  }
365 
366  // Store cross-section information in pb and event weight. The latter is
367  // usually dimensionless, but in units of pb for Les Houches strategies +-4.
368  if (m_store_xsec && pyinfo != 0) {
369  HepMC::GenCrossSection xsec;
370  xsec.set_cross_section( pyinfo->sigmaGen() * 1e9,
371  pyinfo->sigmaErr() * 1e9);
372  evt->set_cross_section(xsec);
373  // If multiweights with possibly different xsec, overwrite central value
374  std::vector<double> xsecVec = pyinfo->weightContainerPtr->getTotalXsec();
375  if (xsecVec.size() > 0) {
376  xsec.set_cross_section(xsecVec[0]*1e9);
377  evt->set_cross_section(xsec);
378  }
379  }
380 
381  if (m_store_weights && pyinfo != 0) {
382  for (int iweight = 0; iweight < pyinfo->numberOfWeights();
383  ++iweight) {
384  std::string name = pyinfo->weightNameByIndex(iweight);
385  double value = pyinfo->weightValueByIndex(iweight);
386  evt->weights()[name] = value;
387  }
388  }
389 
390  // Done for new event.
391  return true;
392 
393 }
394 
395 //==========================================================================
396 
397 } // end namespace HepMC
398 
399 namespace Pythia8 {
400 
401 //==========================================================================
402 // This a wrapper around HepMC::Pythia8ToHepMC in the Pythia8
403 // namespace that simplify the most common use cases. It stores the
404 // current GenEvent and output stream internally to avoid cluttering
405 // of user code. This class is also defined in HepMC3.h with the same
406 // signatures, and the user can therefore switch between HepMC version
407 // 2 and 3, by simply changing the include file.
409 
410 public:
411 
412  // We can either have standard ascii output or none at all.
413  enum OutputType { none, ascii2 };
414 
415  // Typedef for the version 2 specific classes used.
416  typedef HepMC::GenEvent GenEvent;
417  typedef shared_ptr<GenEvent> EventPtr;
418  typedef HepMC::IO_GenEvent Writer;
419  typedef shared_ptr<Writer> WriterPtr;
420 
421  // The empty constructor does not creat an aoutput stream.
423 
424  // Construct an object with an internal output stream.
425  Pythia8ToHepMC(string filename, OutputType ft = ascii2) {
426  setNewFile(filename, ft);
427  }
428 
429  // Open a new external output stream.
430  bool setNewFile(string filename, OutputType ft = ascii2) {
431  switch ( ft ) {
432  case ascii2:
433  writerPtr = make_shared<HepMC::IO_GenEvent>(filename);
434  break;
435  case none:
436  writerPtr = nullptr;
437  break;
438  }
439  return writerPtr != nullptr;
440  }
441 
442  // Create a new GenEvent object and fill it with information from
443  // the given Pythia object.
444  bool fillNextEvent(Pythia & pythia) {
445  geneve = make_shared<GenEvent>();
446  return fill_next_event(pythia, *geneve);
447  }
448 
449  // Write out the current GenEvent to the internal stream.
450  void writeEvent() {
451  writerPtr->write_event(&*geneve);
452  }
453 
454  // Create a new GenEvent object and fill it with information from
455  // the given Pythia object and write it out directly to the
456  // internal stream.
457  bool writeNextEvent(Pythia & pythia) {
458  if ( !fillNextEvent(pythia) ) return false;
459  writeEvent();
460  return true;
461  }
462 
463  // Get a reference to the current GenEvent.
464  GenEvent & event() {
465  return *geneve;
466  }
467 
468  // Get a pointer to the current GenEvent.
469  EventPtr getEventPtr() {
470  return geneve;
471  }
472 
473  // Get a reference to the internal stream.
474  Writer & output() {
475  return *writerPtr;
476  }
477 
478  // Get a pointer to the internal stream.
479  WriterPtr outputPtr() {
480  return writerPtr;
481  }
482 
483  // Set cross section information in the current GenEvent.
484  void setXSec(double xsec, double xsecerr) {
485  HepMC::GenCrossSection xs;
486  xs.set_cross_section(xsec, xsecerr);
487  geneve->set_cross_section(xs);
488  }
489 
490  // Update all weights in the current GenEvent.
491  void setWeights(const vector<double> & wv) {
492  if ( wv.size() != weightNames.size() )
493  geneve->weights() = wv;
494  else
495  for ( int i= 0, N = wv.size(); i < N; ++i )
496  geneve->weights()[weightNames[i]] = wv[i];
497  }
498 
499  // Set all weight names in the current run.
500  void setWeightNames(const vector<string> &wnv) {
501  weightNames = wnv;
502  }
503 
504  // Update the PDF information in the current GenEvent
505  void setPdfInfo(int id1, int id2, double x1, double x2,
506  double scale, double xf1, double xf2,
507  int pdf1 = 0, int pdf2 = 0) {
508  HepMC::PdfInfo pdf(id1, id2, x1, x2, scale, xf1, xf2, pdf1, pdf2);
509  geneve->set_pdf_info(pdf);
510  }
511 
512 private:
513 
514  // The current GenEvent.
515  EventPtr geneve = nullptr;
516 
517  // The output stream.
518  WriterPtr writerPtr = nullptr;
519 
520  // The weight names in the current run.
521  vector<string> weightNames;
522 
523 };
524 
525 }
526 
527 #endif // end Pythia8_HepMC2_H
Pythia8ToHepMC()
Constructor and destructor.
Definition: HepMC2.h:83
void setWeightNames(const vector< string > &wnv)
Set all weight names in the current run.
Definition: HepMC2.h:500
bool flag(string keyIn)
Give back current value, with check that key exists.
Definition: Settings.cc:1478
GenEvent & event()
Get a reference to the current GenEvent.
Definition: HepMC2.h:464
int nCollND() const
Definition: HIInfo.h:94
Definition: Info.h:45
The Event class holds all info on the generated event.
Definition: Event.h:453
Exception thrown when an undecayed parton is written into the record.
Definition: HepMC2.h:47
void setPdfInfo(int id1, int id2, double x1, double x2, double scale, double xf1, double xf2, int pdf1=0, int pdf2=0)
Update the PDF information in the current GenEvent.
Definition: HepMC2.h:505
Pythia8ToHepMC(string filename, OutputType ft=ascii2)
Construct an object with an internal output stream.
Definition: HepMC2.h:425
vector< double > getTotalXsec()
Return cross section estimate for total run.
Definition: Weights.cc:1025
Settings settings
Settings: databases of flags/modes/parms/words to control run.
Definition: Pythia.h:332
double b() const
The impact-parameter distance in the current event.
Definition: HIInfo.h:43
bool writeNextEvent(Pythia &pythia)
Definition: HepMC2.h:457
The Pythia8ToHepMC class.
Definition: HepMC2.h:78
double sigmaErr(int i=0) const
The uncertainty on the estimated cross-section in units of mb.
Definition: Info.h:314
Event event
The event record for the complete event history.
Definition: Pythia.h:323
Base exception for all exceptions that Pythia8ToHepMC might throw.
Definition: HepMC2.h:36
PartonEndVertexException(int i, int pdg_idIn)
Constructor and destructor.
Definition: HepMC2.h:52
int id1pdf(int i=0) const
Hard process flavours, x values, parton densities, couplings, Q2 scales.
Definition: Info.h:169
void setXSec(double xsec, double xsecerr)
Set cross section information in the current GenEvent.
Definition: HepMC2.h:484
void writeEvent()
Write out the current GenEvent to the internal stream.
Definition: HepMC2.h:450
virtual const char * what() const
Throw exception.
Definition: HepMC2.h:62
EventPtr getEventPtr()
Get a pointer to the current GenEvent.
Definition: HepMC2.h:469
Pythia8ToHepMC()
The empty constructor does not creat an aoutput stream.
Definition: HepMC2.h:422
HIInfo * hiInfo
Definition: Info.h:113
Definition: HepMC2.h:408
int nCollDD() const
Definition: HIInfo.h:109
bool fill_next_event(Pythia8::Event &pyev, GenEvent &evt, int ievnum=-1, const Pythia8::Info *pyinfo=0, Pythia8::Settings *pyset=0, bool append=false, GenParticle *rootParticle=0, int iBarcode=-1)
Alternative method to convert Pythia events into HepMC ones.
Definition: HepMC2.h:100
double sigmaGen(int i=0) const
The estimated cross-section in units of mb.
Definition: Info.h:306
bool setNewFile(string filename, OutputType ft=ascii2)
Open a new external output stream.
Definition: HepMC2.h:430
int size() const
Event record size.
Definition: Event.h:504
int nAbsTarg() const
Definition: HIInfo.h:139
WriterPtr outputPtr()
Get a pointer to the internal stream.
Definition: HepMC2.h:479
bool print_inconsistency() const
Read out values for some switches.
Definition: HepMC2.h:113
int nAbsProj() const
Definition: HIInfo.h:123
void setWeights(const vector< double > &wv)
Update all weights in the current GenEvent.
Definition: HepMC2.h:491
HepMC::GenEvent GenEvent
Typedef for the version 2 specific classes used.
Definition: HepMC2.h:416
void set_print_inconsistency(bool b=true)
Set values for some switches.
Definition: HepMC2.h:123
int nCollNDTot() const
The total number of non-diffractive sub collisions in the current event.
Definition: HIInfo.h:97
OutputType
We can either have standard ascii output or none at all.
Definition: HepMC2.h:413
Header for classes to set beam momentum and interaction vertex spread.
Definition: Analysis.h:20
The Pythia class contains the top-level routines to generate an event.
Definition: Pythia.h:71
bool fillNextEvent(Pythia &pythia)
Definition: HepMC2.h:444
const Info & info
Public information and statistic on the generation.
Definition: Pythia.h:326
Writer & output()
Get a reference to the internal stream.
Definition: HepMC2.h:474
int nDiffTarg() const
Definition: HIInfo.h:143
Definition: HepMC2.h:30
int index() const
Return info on location and flavour of bad parton.
Definition: HepMC2.h:65
int nDiffProj() const
Definition: HIInfo.h:127
bool fill_next_event(Pythia8::Pythia &pythia, GenEvent &evt, int ievnum=-1, bool append=false, GenParticle *rootParticle=0, int iBarcode=-1)
The recommended method to convert Pythia events into HepMC ones.
Definition: HepMC2.h:90
Definition: Settings.h:195