LHAinit and LHAevnt classes are base classes,
containing reading and printout functions, plus a pure virtual function
each. Derived classes have to provide these two virtual functions to do
the actual work. The existing derived classes are for reading information
from a Les Houches Event File (LHEF) or from the respective Fortran
commonblock.
Normally, pointers to objects of the derived classes should be handed
in with the
pythia.init( LHAinit*, LHAevnt*) method. However,
with the LHEF format a filename can replace the two pointers, see below.
LHAinit class stores information equivalent to the
/HEPRUP/ commonblock, as required to initialize the event
generation chain. The main difference is that the vector container
now allows a flexible number of subprocesses to be defined. For the
rest, names have been modified, since the 6-character-limit does not
apply, and variables have been regrouped for clarity, but nothing
fundamental is changed.
The pure virtual function set() has to be implemented in the
derived class, to set relevant information when called. It should
return false if it fails to set the info.
Inside set(), such information can be set by the following
methods:
method beamA( identity, energy, pdfGroup, pdfSet) IDBMUP(1), EBMUP(1), PDFGUP(1), PDFSUP(1)), and similarly
a beamB method exists. The parton distribution information
defaults to zero, meaning that internal sets are used.
method strategy( strategyIn) argument strategyIn :
chosen strategy (cf. IDWTUP; see [Sjo06]
section 9.9.1 for extensive comments).
argumentoption 1 : events come with non-negative weight, given in units
of pb, with an average that converges towards the cross section of the
process. PYTHIA is in charge of the event mixing, i.e. for each new
try decides which process should be generated, and then decides whether
is should be kept based on a comparison with xMax.
Accepted events therefore have unit weight.
argumentoption -1 : as option 1, except that cross sections can now be
negative and events after unweighting have weight +-1. You can use
pythia.info.weight()
to find the weight of the current event. A correct event mixing requires
that a process that can take both signs should be split in two, one limited
to positive or zero and the other to negative or zero values, with
xMax chosen appropriately for the two.
argumentoption 2 : events come with non-negative weight, in unspecified
units, but such that xMax can be used to unweight the events
to unit weight. Again PYTHIA is in charge of the event mixing.
The total cross section of a process is stored in
xSec.
argumentoption -2 : as option 2, except that cross sections can now be
negative and events after unweighting have weight +-1. As for option -1
processes with indeterminate sign should be split in two.
argumentoption 3 : events come with unit weight, and are thus accepted
as is. The total cross section of the process is stored in
xSec.
argumentoption -3 : as option 3, except that events now come with weight
+-1. Unlike options -1 and -2 processes with indeterminate sign need not be
split in two, unless you intend to mix with internal PYTHIA processes
(see below).
argumentoption 4 : events come with non-negative weight, given in units
of pb, with an average that converges towards the cross section of the
process, like for option 1. No attempt is made to unweight the events,
however, but all are generated in full, and retain their original weight.
For consistency with normal PYTHIA units, the weight stored in
pythia.info.weight() has been converted to mb, however.
argumentoption -4 : as option 4, except that events now can come
either with positive or negative weights.
method process( idProcess, xSec, xErr, xMax) LPRUP, XSECUP, XERRUP,
XMAXUP).
Each new call will append one more entry to the list of processes.
The choice of strategy determines which quantities are mandatory:
xSec for strategies +-2 and +-3,
xErr never, and
xMax for strategies +-1 and +-2.
xErr values, but calculates a statistical cross section
error based on the spread of event-to-event weights. This should work
fine for strategy options +-1, but not for the others. Specifically,
for options +-2 and +-3 the weight spread may well vanish, and anyway
is likely to be an underestimate of the true error. If the author of the
LHA input information does provide error information you may use that -
this information is displayed at initialization. If not, then a relative
error decreasing like 1/sqrt(n_acc), where n_acc
is the number of accepted events, should offer a reasonable estimate.
Information is handed back by the following methods:
method idBeamA(), eBeamA(), pdfGroupBeamA(), pdfSetBeamA() method strategy() method size() method idProcess(i), xSec(i), xErr(i), xMax(i) i in the range 0 <= i <
size().
The information can also be printed using the list()
method, e.g. LHAinitObject->list().
This is automatically done by the pythia.init call.
LHAevnt class stores information equivalent to the
/HEPEUP/ commonblock, as required to hand in the next
parton-level configuration for complete event generation. The main
difference is that the vector container now allows a flexible number
of partons to be defined. For the rest, names have been modified,
since the 6-character-limit does not apply, and variables have been
regrouped for clarity, but nothing fundamental is changed.
The LHA standard is based on Fortran arrays beginning with
index 1, and mother information is defined accordingly. In order to
be compatible with this convention, the zeroth line of the C++ particle
array is kept empty, so that index 1 also here corresponds to the first
particle. One small incompatibility is that the size()
method returns the full size of the particle array, including the
empty zeroth line, and thus is one larger than the true number of
particles (NUP).
The pure virtual function set(idProcess) has to be
implemented in the derived class, to set relevant information when
called. For strategy options +-1 and +-2 the input
idProcess value specifies which process that should be
generated by set(...), while idProcess is
irrelevant for strategies +-3 and +-4. The set(...)
function should return false if it fails to set the info,
i.e. normally that the supply of events in a file is exhausted.
If so, no event is generated, and pythia.next() returns
false. You can then interrogate
pythia.info.atEndOfFile()
to confirm that indeed the failure is caused in the set(...)
function, and decide to break out of the event generation loop.
Inside a normal set(...) call, information can be set by
the following methods:
method process( idProcess, weight, scale, alphaQED, alphaQCD) IDPRUP, XWTGUP, SCALUP, AQEDUP, AQCDUP). This method
also resets the size of the particle list, and adds the empty zeroth
line, so it has to be called before the particle method below.
method particle( id, status, mother1, mother2, colourTag1, colourTag2, p_x, p_y, p_z, e, m, tau, spin) IDUP, ISTUP,
MOTHUP(1,..), MOTHUP(2,..), ICOLUP(1,..), ICOLUP(2,..), PUP(J,..),
VTIMUP, SPINUP) .
Information is handed back by the following methods:
method idProcess(), weight(), scale(), alphaQED(), alphaQCD() pythia.info.weight() as a rule
is not the same as the above weight(): the method here gives
the value before unweighting while the one in info gives
the one after unweighting and thus normally is 1 or -1. Only with strategy
options +-4 would the value in info be the same as here,
except for a conversion from pb to mb.
method size() method id(i), status(i), mother1(i), mother2(i), col1(i), col2(i),px(i), py(i), pz(i), e(i), m(i), tau(i), spin(i) i in the range
0 <= i < size(). (But again note that
i = 0 is an empty line, so the true range begins at 1.)
In the LHEF description [Alw06] an extension to
include information on the parton densities of the colliding partons
is suggested. This optional further information can be set by
method pdf( id1, id2, x1, x2, scalePDF, xpdf1, xpdf2) method pdfIsSet(), id1(), id2(), x1(), x2(), scalePDF(), xpdf1(), xpdf2() pdf(...) must be called after the
process(...) call of the event for this to work.)
The information can also be printed using the list()
method, e.g. LHAevntObject->list().
In cases where the LHAevntObject is not available to the
user, the pythia.LHAevntList() method can be used, which
is a wrapper for the above.
The LHA expects the decay of resonances to be included as part of the
hard process, i.e. if unstable particles are produced in a process then
their decays are also described. This includes Z^0, W^+-, H^0
and other short-lived particles in models beyond the Standard Model.
Should this not be the case then PYTHIA will perform the decays of all
resonances it knows how to do, in the same way as for internal processes.
Note that you will be on slippery ground if you then restrict the decay of
these resonances to specific allowed channels since, if this is not what
was intended, you will obtain the wrong cross section and potentially the
wrong mix of different event types. (Since the original intention is
unknown, the cross section will not be corrected for the fraction of
open channels, i.e. the procedure used for internal processes is not
applied in this case.)
pythia.init(fileName)
initialization option exists, where the LHEF name is provided as single
input. Internally this name is then used to create instances of two derived
classes, LHAinitLHEF and LHAevntLHEF. Both of
them are allowed to read from the same LHEF, first the former and then
the latter.
An example how to generate events from a LHEF is found in
main12.cc. Note the use of
pythia.info.atEndOfFile() to find out when the whole
LHEF has been processed.
To allow the sequential use of several event files the
init(...) method has an optional second argument:
pythia.init(fileName, bool skipInit = false).
If called with this argument true then there will
be no initialization, except that the existing LHAinitLHEF
and LHAevntLHEF classes will be deleted and replaced
by ones pointing to the new file. It is assumed (but never checked) that
the initialization information is identical, and that the new file simply
contains further events of exactly the same kind as the previous one.
An example of this possibility, and the option to mix with internal
processes, is found in main13.cc.
LHAFortran.h file, that is not included in any of the
other library files. Instead it should be included in the
user-supplied main program, together with the implementation of two
methods below that call the Fortran program to do its part of the job.
The LHAinitFortran class derives from LHAinit.
It reads initialization information from the LHA standard
Fortran commonblock, assuming this commonblock behaves like an
extern "C" struct named heprup_. (Note the final
underscore, to match how the gcc compiler internally names Fortran
files.)
Initialization is with
LHAinitFortran lhaInit();
i.e. does not require any arguments.
The user has to supply an implementation of the fillHepRup()
method, that is to do the actual calling of the external Fortran routine(s)
that fills the HEPRUP commonblock. The translation of this
information to the C++ structure is provided by the existing
set() code.
The LHAevntFortran class derives from LHAevnt.
It reads information on the next event, stored in the LHA
standard Fortran commonblock, assuming this commonblock behaves like
an extern "C" struct named hepeup_.
Initialization is with
LHAevntFortran lhaEvnt();
i.e. does not require any arguments.
The user has to supply an implementation of the fillHepEup()
method, that is to do the actual calling of the external Fortran routine(s)
that fills the HEPEUP commonblock. The translation of this
information to the C++ structure is provided by the existing
set() code.
See further
here for
information how PYTHIA 6.4 can be linked to make use of this facility.