1 <chapter name="HepMC Interface">
3 <h2>HepMC Interface</h2>
5 An interface to the HepMC <ref>Dob01</ref> standard event record
6 format has been provided by M. Kirsanov. To use it, the relevant
7 libraries need to be linked, as explained in the <code>README</code>
8 file. Only version 2 of HepMC is supported. (Version 1 requires
9 a different interface structure, which was only supported up until
13 The (simple) procedure to translate PYTHIA 8 events into HepMC ones
14 is illustrated in the <code>main41.cc</code>, <code>main42.cc</code>
15 <code>main61.cc</code> and <code>main62.cc</code>
16 main programs. At the core is a call to the
18 HepMC::I_Pythia8::fill_next_event( pythia, hepmcevt, ievnum = -1, convertGluonTo0 = false )
20 which takes a reference of the generator object and uses it, on the one
21 hand, to read out and covert the event record in <code>pythia.event</code>
22 and, on the other hand, to extract and store parton-density (PDF)
23 information for the hard subprocess from <code>pythia.info</code>.
24 The optional last argument, if <code>true</code>, allows you to store
25 gluons as "PDG" code 0 rather than the normal 21; this only applies to
26 the PDF information, not the event record.
29 The earlier version of this routine,
31 HepMC::I_Pythia8::fill_next_event( pythia.event, hepmcevt, ievnum = -1 )
33 is retained (for now) for backwards compatibility. It takes a PYTHIA event
34 as input and returns a HepMC one, but without storing the PDF information.
35 The latter could then instead be stored by a separate call
37 HepMC::I_Pythia8::pdf_put_info( hepmcevt, pythia, convertGluonTo0 = false )
42 While PYTHIA stores momenta in units of GeV, with <ei>c = 1</ei>,
43 the HepMC standard originally did not specify which units to use,
44 unfortunately. Later versions allow units to be specified either as
45 MeV or as GeV. The dividing line goes with the introduction of the
46 <code>HEPMC_HAS_UNITS</code> environment variable in HepMC 2.04.
47 When PYTHIA is linked to older versions, which do not have this
48 environment variable, the default behaviour of
49 <code>fill_next_event</code> is to convert momenta to MeV.
50 To store momenta in GeV instead, you must call
52 HepMC::I_Pythia8::set_convert_to_mev( false )
54 beforehand. When linked to newer HepMC versions, PYTHIA will adapt
55 to the unit specified for the HepMC event record. A default unit
56 is chosen when HepMC is built, but this choice can be overridden
57 in the constructor of a <code>HepMC::GenEvent</code>, to be either
58 GeV or MeV. The <code>main41.cc</code> code illustrates these
62 Also for length units there could exist ambiguities, but the mm
63 choice of PYTHIA seems to be shared by most other programs.
64 When linked to older HepMC versions there is never a conversion,
65 while for newer versions the PYTHIA interface will convert to the
66 units set for the HepMC event record, as for momenta.
69 By agreement with the LHC experimental community, support for the
70 older behaviour will be dropped for all versions released after
74 The status code is now based on the new standard for HepMC 2.05,
75 see the <aloc href="EventRecord">Event::statusHepMC(...)</aloc>
76 conversion routine for details. The earlier behaviour, where all
77 final particles had status 1 and all initial or intermediate ones
78 status 2, is available as a commented-out line in the
79 <code>I_Pythia8::fill_next_event(...)</code> method, and so is simple
83 In HepMC 2.05 it also becomes possible to store the generated cross
84 section and its error. The environment variable
85 <code>HEPMC_HAS_CROSS_SECTION</code> is used to check whether this
86 possibility exists and, if it does the current values from
87 <code>pythia.info.sigmaGen()</code> and
88 <code>pythia.info.sigmaErr()</code> are stored for each event,
89 multiplied by <ei>10^9</ei> to convert from mb to pb. Note that
90 PYTHIA improves its accuracy by Monte Carlo integration in the course
91 of the run, so the values associated with the last generated event
92 should be the most accurate ones. (If events also come with a dimensional
93 weight, like in some Les Houches strategies, conversion from mb to fb
94 for that weight must be set by hand, see the last method below.)
96 <h2>The public methods</h2>
98 Here comes a complete list of all public methods of the
99 <code>I_Pythia8</code> class in the <code>HepMC</code>
100 (<i>not</i> <code>Pythia8</code>!) namespace.
102 <method name="I_Pythia8::I_Pythia8()">
104 <methodmore name="virtual I_Pythia8::~I_Pythia8()">
105 the constructor and destructor take no arguments.
108 <method name="bool I_Pythia8::fill_next_event( Pythia8::Pythia& pythia,
109 GenEvent* evt, int ievnum = -1, bool convertGluonTo0 = false)">
110 convert a <code>Pythia</code> event to a <code>HepMC</code> one.
111 Will return true if succeeded.
112 <argument name="pythia">
113 the input <code>Pythia</code> generator object, from which both the
114 event and the parton density information can be obtained.
116 <argument name="evt">
117 the output <code>GenEvt</code> event, in its standard form.
119 <argument name="iev">
120 set the event number of the current event. If negative then the
121 internal event number is used, which is incremented by one for
124 <argument name="convertGluonTo0">
125 the normal gluon identity code 21 is used also when parton density
126 information is stored, unless this optional argument is set true to
127 have gluons represented by a 0. This choice does not affect the
128 normal event record, where a gluon is always 21.
132 <method name="bool I_Pythia8::fill_next_event( Pythia8::Event& pyev,
133 GenEvent* evt, int ievnum = -1 )">
134 convert a <code>Pythia</code> event to a <code>HepMC</code> one.
135 Will return true if succeeded. Do not store parton-density information.
136 <argument name="pyev">
137 the input <code>Pythia</code> event, in its standard form.
139 <argument name="evt">
140 the output <code>GenEvt</code> event, in its standard form.
142 <argument name="iev">
143 set the event number of the current event. If negative then the
144 internal event number is used, which is incremented by one for
149 <method name="bool I_Pythia8::put_pdf_info( GenEvent* evt,
150 Pythia8::Pythia& pythia, bool convertGluonTo0 = false )">
151 append parton-density information to an event already stored
152 by the previous method.
153 <argument name="evt">
154 the output <code>GenEvt</code> event record, in its standard form.
156 <argument name="pythia">
157 the input <code>Pythia</code> generator object, from which both the
158 event and the parton density information can be obtained.
160 <argument name="convertGluonTo0">
161 the normal gluon identity code 21 is used also when parton density
162 information is stored, unless this optional argument is set true to
163 have gluons represented by a 0.
168 The following methods can be used to set, and in some cases interrogate,
169 the status of some switches that can be used to modify the behaviour
170 of the conversion routine. The <code>set</code> methods have the
171 same default input values as the internal initialization ones, so
172 that a call without an argument (re)stores the default.
174 <method name="void I_Pythia8::set_trust_mothers_before_daughters(
177 <methodmore name="bool I_Pythia8::trust_mothers_before_daughters()">
178 if there is a conflict in the history information, then trust the
179 information on mothers above that on daughters. Currently this is
180 the only option implemented.
183 <method name="void I_Pythia8::set_trust_both_mothers_and_daughters(
186 <methodmore name="bool I_Pythia8::trust_both_mothers_and_daughters()">
187 currently dummy methods intended to resolve conflicts in the event
191 <method name="void I_Pythia8::set_print_inconsistency_errors(
194 <methodmore name="bool I_Pythia8::print_inconsistency_errors()">
195 print a warning line, on <code>cerr</code>, when inconsistent mother
196 and daughter information is encountered.
199 <method name="void I_Pythia8::set_crash_on_problem( bool b = false )">
200 if problems are encountered then the run is interrupted by an
201 <code>exit(1)</code> command. Default is not to crash.
204 <method name="void I_Pythia8::set_freepartonwarnings( bool b = true )">
205 interrupt the run by an <code>exit(1)</code> command if unhadronized
206 gluons or quarks are encountered in the event record, unless
207 hadronization is switched off. Default is to crash.
210 <method name="void I_Pythia8::set_convert_to_mev( bool b = false )">
211 convert the normal GeV energies, momenta and masses to MeV.
216 <!-- Copyright (C) 2012 Torbjorn Sjostrand -->