HepMC3 event record library
|
Namespaces | |
LHEF | |
Les Houches event file classes. | |
This module contains helper classes and Reader and Writer classes for handling Les Houches event files - LHEF.
The Les Houches accord on an event file format (LHEF) to be used for passing events from a matrix element generator program (MEG) to an event generator program (EG) implementing parton showers, underlying event models, and hadronisation models etc., was not originally included in the HepMC event record format. But as the demand for more information to be included in HepMC, it was decided to allow HepMC to include also the original information from a MEG in the LHEF format (see the run attribute HepMC3::HEPRUPAttribute and event attribute HepMC3::HEPEUPAttribute). A separate /standard/ implementation in C++ of the LHEF format had already been maintained by Leif Lönnblad, and it was decided to include this (header only - LHEF.h) as a part of HepMC3. This will both be used in above mentioned HepMC3::Attribute classes and as a kind of definition of the LHEF format, which so far has not been extremely well documented. From now on these pages will serve as the defining information about the LHEF format.
The original Les Houches accord for communicating between MEGs and EGs was agreed upon in 2001 [arXiv:hep-ph/0109068] and consisted of two simple FORTRAN common blocks. In fact this structure survived in the LHEF format, which was introduced in 2006 [arXiv:hep-ph/0609017], and is still there after the updated versions 2 in 2009 [arXiv:1003.1643], and 3 in 2013 [arXiv:1405.1067], and in the current proposal developed at the Les Houches workshop on TeV Colliders 2015 .
As the methods for combining MEGs and EGs has advanced since the first accord, from the tree-level merging methods and NLO matching at the turn of the millennium, to the multi-jet (N)NLO matching and merging methods being perfeted to day, the LHEF format has developed and a lot of optional information can be passed beyond the original common block structures. In the following all features included will be described, also those that were added a bit prematurely and later became deprecated.
The LHEF format is based on XML, but has some oddities that goes beyond pure XML. As the name indicates, XML is extensible, and anyone writing a LHEF file can add whatever information she or he wants, however the following basic stucture must be observed.
This looks like fairly normal XML tags, and indeed they are. The only addition to the structure is that the init
and event
(and their respective end tags) are required to be alone on a line, and the content of these blocks are required to start with a number of lines on a specific format that follows exactly the structure of the fortran common block in original Les Houches Accord. This means that the first line in the init
block must start with a line containing the numbers
and the following NPRUP
lines should be numbers on the form
where the different variable names are defined as follows:
IDBMUP
: the PDG numbers for the incoming beams; EBMUP
: the energy (in GeV) of the incoming beams; PDFGUP
and PDFSUP
: the LHAPDF group and set numbers for the corresponding parton densities used; IDWTUP
: the weight strategy used; NPRUP
The number of different processes included in the file; and for each process IPR
:
XSECUP
, XERRUP
: the Monte Carlo integrated cross section and error estimate; XMAXUP
: the overestimated cross section such that the sum over the individual weights in each <event>
times this file times XMAXUP
equals XSECUP
times the number of events; LPRUP
: a unique integer identifying the corresponding process. In the LHEF::Reader and LHEF::Writer classes this information is available as the public heprup
member of class LHEF::HEPRUP with public members mimicking the Fortran common block variables.
Similarly, every event
block must start with a line containing then numbers
and the following NUP
lines should be numbers describing each particle on the form
where the different variable names are defined as follows:
NUP
: the number of particles in the event; IDPRUP
: the integer identifying the process; XWGTUP
: the (default) weight of the event SCALUP
: the scale of the hard process in GeV; AQEDUP
: the value of αEM used; AQCDUP
: the value of αS used; and for each particle I
:
IDUP
: the PDG code for the particle; ISTUP
: the status code for the particle; MOTHUP
: the line numbers of two particles considered to be the mothers of this particle; ICOLUP
: indices of the colour and anti-colour lines connected to this particle; PUP
: the x, y, z and energy component of the 4-momentum and invariant masss (in units of GeV); VTIMUP
; the proper lifetime of this particle; SPINUP
: the spin. In the LHEF::Reader and LHEF::Writer classes this information is available as the public hepeup
member of class LHEF::HEPEUP with public members mimicking the Fortran common block variables.
Over the years several additional XML-tags has been formalised to specify information on top of what is given in the original Les Houches accord common block. These are listed below. In most cases the tag name corresponds to a class with a corresponding name available as suitably named public members in the LHEF::HEPRUP and LHEF::HEPEUP class respectively.
Note that a tag may contain attributes in the following ways:
where the second case is a tag with only attributes an no contents. In the following an attribute may be described as required (R) or optional with a default value (D).
The <init>
block contains information about the full run (similar to the information contained in HepMC3::GenRunInfo). The following tags are defined.
<generator>
(optional, multiple, see LHEF::HEPRUP::generators)version
can be given with a string containg a version string. The content of the tag can include any generator specific information. <xsecinfo>
(required, multiple, see LHEF::HEPRUP::xsecinfos)weightname
: in case of multiple weights for each event several xsecinfo
tags can be given, in which case This should give the corresponding weight name. If missing This is assumed to describe the default weight. neve
(R): the number of events in the file totxsec
(R): the total cross section (in units of pb) of all processes in the file maxweight
(D=1): the maximum weight of any event in the file (in an arbitrary unit) minweight
: if the file contains negative weights, the minweight may contain the most negative of the negative weights in the file. If not given, minweight is assumed to be -maxweight
. meanweight
(D=1): The average weight of the events in the file (same unit as maxweight
). negweights
(D=no): If yes, the file may contain events with negative weights. varweights
(D=no): If yes, the file may contain varying weights (if no, all events are weighted with maxweight
(or minweight
)). <cutsinfo>
(optional, multiple, see LHEF::HEPRUP::cuts)cut
and ptype
tags can be supplied. <ptype>
(optional, multiple, see LHEF::HEPRUP::ptypes)cutinfo
tag to define a group of particles to which a cut can be applied. Its contents should be a white-space separated list of PDG particle codes. The only attribute is name
(R): the string used to represent this ptype
in a cut
. <cut>
(optional, multiple, see LHEF::HEPRUP::cuts)
This tag has information of a particular cut used. The information needed is which particle(s) are affected, what variable is used and the maximum and/or minimum value of that parameter. The contents should be the minimum and maximum allowed values of this variable (in that order). If only one number is given, there is no maximum. If the minumum is larger or equal to the maximum, there is no minimum. The attributes are:
p1
(D=0): The particles for which this cut applies. This can either be a number corresponding to a given particle PDG code or 0 for any particle or the name of a previously defined ptype
tag. p2
(D=0): If cut is between pairs of particles, this attribute should be non-zero. Allowed values are the same as for p1
. type
(R) This defines the variable which is cut. The following values are standardised (the lab frame is assumed in all cases): <procinfo>
(optional, multiple, see LHEF::HEPRUP::procinfo)iproc
(D=0): The process number for which the information is given. 0 means all processes. loops
: The number of loops used in calculating this process. qcdorder
: The number of QCD vertices used in calculating this process. eworder
: The number of electro-weak vertices used in calculating this process. rscheme
: The renormalization scheme used, if applicable. fscheme
: The factorization scheme used, if applicable. scheme
: Information about the scheme used to calculate the matrix elements. If absent, a pure tree-level calculation is assumed. Other possible values could be CSdipole (NLO calculation with Catani-Seymour subtraction), MC, POWHEG and NLOexclusive (NLO calculation according to the exclusive cross section withing the given cuts). <mergeinfo>
(DEPRECATED, multiple, see LHEF::HEPRUP::mergeinfo)iproc
(D=0): The process number for which the information is given. "0" means all processes. mergingscale
(R): The value of the merging scale in GeV. maxmult
(D=no): If yes the corresponding process is reweighted as if it is the maximum multiplicity process (ie. the Sudakov for the last step down to the merging scale is not included. <eventfile>
(optional, multiple, see LHEF::HEPRUP::eventfiles)init
block, eg. when the event files become extremely long or when one wants to write the init
block after having generated all events in order to collect the correct information from the run. The names of these files are then specified in eventfile
tags with attributes: name
(R): the file name, interpreted as an absolute path if starting with "/", otherwise as a relative path to where the file with the init
block is located; neve
: the number of events in the file ntries
: the number of attempts the generator needed to produce the neve
events (for statistics purposes. <weightinfo>
(optional, multiple, see LHEF::HEPRUP::weightinfo)weightinfo
tag. The default weight (as given by the LHEF:HEPEUP::XWGTUP variable) is always treated as index 0 and given the name "Default", while the additional weights are indexed by the order of the weightinfo
tags. The attributes are: name
(R): the name of the weight (in the best of all worlds this would be standardised); muf
: he factor multiplying the nominal factorisation scale of the event for the given weight; mur
: he factor multiplying the nominal renormalisation scale of the event for the given weight; pdf
(D=0): the LHAPDF code used for the given weight (where 0 corresponds to the default PDF of the run); pdf2
(D=pdf): the LHAPDF code used for the second beam for the given weight; <weightgroup>
(optional, multiple, see LHEF::HEPRUP::weightgroup)weightinfo
tags together. The only well defined attribute is name
giving a string which will be combined with the name
attribute of the included weightinfo
tags to give the HepMC weight names. Also other attributes can be included to include information about the weight variation used and is available in LHEF::WeightGroup::attributes. <initrwgt>
(optional, multiple, see LHEF::HEPRUP::weightinfo)weightinfo
tag. It accepts the same attributes as weightinfo
, except that the name attribute is named id
. Note that some generators puts these tags in the header
block. The current implementation of LHEF::Reader cannot handle this. Note that this is handled by the same classes (LHEF::WeightInfo and LHEF::WeightGroup) in LHEF::Reader and LHEF::Writer. After the init
block any number of events can be given. In addition events can be given in separate files declared with eventfile
tags in the init
block.
The main tag here is simply called event
and can (for statistics purposes) have an attribute ntries
specifying how many attempts the generator needed to produce the event. Also other attributes can be given (and will be stored in the LHEF::HEPEUP::attributes member variable).
The event
tags may be grouped together in a eventgroup
tag. This is useful mainly for NLO generators that produces a (number) "real" event(s) accompanied with a number of "counter" events, where these events should be treated together for statistics purposes. For this reason the eventgroup
tag can be provided with the optional tags nreal
and ncounter
to indicate the number of event
tags included of each type.
As indicated above the block must start with required information corresponding to the original Les Houches Accord Fortran common block. Here is a list of additional tags that may be provided.
<weights>
(optional, see LHEF::HEPEUP::weights)weightinfo
tags in the init
block. If for some obscure reason a weight is absent, it should be anyway be included as a zero weight. Note that the weight for the nominal event is still given by XWGTUP
. This tag has no attributes. <rwgt>
(optional, see LHEF::HEPEUP::weights)initrwgt
this should contain the weights for this event. The difference wrt. the weights
tag is that each weight needs to be given in an wgt
tag. <wgt>
(optional, see LHEF::HEPEUP::weights)initrwgt
tag. The only attribute is the id
of the weight defined in the corresponing weightinfo
tag. <scales>
(optional, see LHEF::HEPEUP::scales)muf
(D=SCLAUP): the factorisation scale (in GeV); mur
(D=SCLAUP): the renormalisation scale (in GeV); mups
(D=SCLAUP): the suggested parton shower starting scale (in GeV). <scale>
(optional inside a <scales>
tag)
It is also possible to specify an individual scale for any particle in the event in a <scale>
tag. This is typically used to define a starting scale for a parton shower in resonance decays. This tag can also be used to specify the scale of a set of particle types. The following attributes can be given:
pos
(optional): the index of the particle for which this scale applies, optionally followed by a space-separated list of indices of recoilers involved when the particle was created. etype
(optional): a space-separated list of PDG codes for which this scale applies. The short-hand notation "EW" can be used to specify all leptons and electro-weak bosons, and "QCD" can be used to specify the guon and all quarks. The contents of the tag gives the scale in GeV. Last update 27 Oct 2020