The Les Houches Accord

The Les Houches Accord for user processes [Boo01] is the standard way to input parton-level information from a matrix-elements based generator into PYTHIA. The conventions for which information should be stored has been defined in a Fortran context, as two commonblocks. Here a C++ equivalent is defined, as two separate classes.

The 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 the respective Fortran commonblock or from a Les Houches Event File.

Normally, pointers to objects of the derived classes should be handed in with the pythia.init( LHAinit*, LHAevnt*) method. However, if you use the LHA runtime interface to PYTHIA 6.4 this is taken care of internally, so no pointers need to be handed in. Also, with the new Les Houches Event File format a filename can replace the two pointers, see below.

Initialization

The 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 name="beamA( identity, energy, pdfGroup, pdfSet)"
sets the properties of the first incoming beam (cf. the Fortran 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 name="strategy( choice)"
sets the event weighting and cross section strategy (cf. IDWTUP).

method name="process( idProcess, crossSection, crossSectionError, crossSectionMaximum)"
sets info on an allowed process (cf. LPRUP, XSECUP, XERRUP, XMAXUP). Each new call will append one more entry to the list of processes.

Information is handed back by the following methods:

method name="idBeamA(), eBeamA(), pdfGroupBeamA(), pdfSetBeamA()"
and similarly with A -> B, for the two beam particles.

method name="strategy()"
for the strategy choice.

method name="size()"
for the number of subprocesses.

method name="idProcess(i), xSec(i), xErr(i), xMax(i)"
for process 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, unless the runtime interface to PYTHIA 6 is being used, in which case that program is intended to print the information.

Event input

The 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 Les Houches 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() 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, e.g. if the supply of events in a file is exhausted.

Inside set(), such information can be set by the following methods:

method name="process( idProcess, weight, scale, alphaQED, alphaQCD)"
tells which kind of process occured, with what weight, at what scale, and which alpha_EM and alpha_strong were used (cf. 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 name="particle( id, status, mother1, mother2, colourTag1, colourTag2, p_x, p_y, p_z, e, m, tau, spin)"
gives the properties of the next particle handed in (cf. IDUP, ISTUP, MOTHUP(1,..), MOTHUP(2,..), ICOLUP(1,..), ICOLUP(2,..), PUP(J,..), VTIMUP, SPINUP) .

Information is handed back by the following methods:

method name="idProc(), weight(), scale(), alphaQED(), alphaQCD()".

method name="size()"
for the size of the particle array, which is one larger than the number of particles in the event, since the zeroth entry is kept empty (see above).

method name="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)"
for particle 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 Les Houches Event File proposal [Alw06] an extension to include information on the parton densities of the colliding partons was suggested. This optional further information can be set by

method name="pdf( id1, id2, x1, x2, scalePDF, xpdf1, xpdf2)"
which gives the flavours , the x and the Q scale (in GeV) at which the parton densities x*f_i(x, Q) have been evaluated.

This information is returned by the methods

method name="pdfIsSet(), id1(), id2(), x1(), x2(), scalePDF(), xpdf1(), xpdf2()"
where the first one tells whether this optional information has been set for the current event. (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.

A runtime Fortran interface

The LHAinitFortran class derives from LHAinit. It reads initialization information from the Les Houches 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 LHAevntFortran class derives from LHAevnt. It reads information on the next event, stored in the Les Houches 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.

See further here for information how PYTHIA 6.4 can be linked to make use of this facility. Several of the example main programs illustrate how it can be used.

An interface to Les Houches Event Files

The new Les Houches Event File (LHEF) standard [Alw06] specifies a format where a single file packs initialization and event information. This is likely to become the most frequently used procedure to process external parton-level events in Pythia. Therefore a special % 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 main14.cc.

Other examples

A special strategy = 10 (not present in the IDWTUP specification) has been added. It takes a given partonic input, no questions asked, and hadronizes it, i.e. does string fragmentation and decay. Thereby the normal process-level and parton-level machineries are bypassed, to the largest extent possible. (Some parts are used, e.g. first to translate the Les Houches event to the process record and later to the event record.) Such an option can therefore be used to feed in ready-made parton-level configurations, without needing to specify where these come from, i.e. there need be no beams or any such explicit information, but of course the user must have taken care of it beforehand.

An example how this can be used for toy studies is found in main16.cc.