HepMC3 event record library
|
For the full list of the installation options including the description of the flags to build the HepMC3 from the sources see at HepMC3 page at CERN GitLab.
A quick minimalist build that requires only C++11 compiler and a recent version of CMake (3.15+) is described below. To build HepMC3 is using CMake on the supported platforms, C++11 compiler and a recent CMake is needed (3.15). The commands executed in the unpacked source tarball of HepMC3
will configure the HepMC3 sources, compile them and install the library into "desired_installation_path".
The following is a list of main differences that should be taken into account when transitioning from HepMC2 to HepMC3.
Following changes in header files have been applied:
The structure of GenCrossSection class has been changed to handle multiple values of cross-sections. The cross-section values and errors (uncertainties) can be accessed only by public functions, the corresponding data members are private. By default, the number of cross-section values in every event is equal to the number of event weights. Accordingly, each cross-section value can be accessed using the corresponding event weight name (std::string) or event weight index (int).
Following header files are no longer available:
An example of interface to Pythia6 Fortran blocks is given in the examples. Please note that the provided interface Pythia6ToHepMC3.cc and Pythia6ToHepMC3.inc is an interface for HepMC3 from Pythia6 and not an interface to Pythia6 from HepMC, as it was in the case of the HepMC2.
Multiple file formats are supported. The implementation of reading and writing is separated in HepMC3. All the reading operations are performed in "reader" objects inherited from HepMC::Reader and the writing operations in the "writer" objects inherited from HePMC::Writer. Therefore it is to use the desired headers explicitly, as needed.
The IO_GenEvent.h header is not available anymore. to write and/or read HepMC2 files the following includes
should be used instead of
Please note that HepMC2 format is outdated and is not able to contain a lot of information stored into event record by the modern Monte Carlo event generators. It is recommended to use HepMC3 native event record in plain text or in ROOT TTree format. The corresponding readers and writers are
and
Implementation of custom Reader and Writer objects is possible as well.
Please, note the difference in the behaviour of default Readers with respect to HepMC2. when reading files with multiple headers. The ASCII files with multiple headers ( e.g. obtained with cat 1.hepmc 2.hepmc > 12.hepmc) will be processed by the readers only till the first occurrence of END_EVENT_LISTING.
In addition to the standard readers, starting for the version 3.2.5 HepMC3 provides as set of templates readers/writers to handle the zip-,lzma-,bz2-compressed files (ReaderGZ and WriterGZ) and to perform multithread reading (ReaderMT).
Particles and vertices are managed using shared pointers, so they should not be created through the call to 'new'.
The iterator-bases classes and access functions from HepMC2, e.g.
were removed. The C++11 iterations should be used instead, e.g. instead of
one should use
or alternatively
Particles and vertices in HepMC3 are stored in topological order. This means that when creating vertices, incoming particles must have id lower than any of the outgoing particles.
This forces the tree structure to be constructed top-to-bottom and disallows creating loops.
Deleting a particle using GenEvent::remove_particle() will also remove its end_vertex if this is the only particle that is on this vertex particles_in() list.
Deleting a vertex will delete all of its outgoing particles. (and subsequently, all of their decays).
The "barcode" integer in HepMC2 was an uncomfortable object, simultaneously declared in the code documentation to be a meaningless unique identifier for vertex and particle objects, and set to specific ranges by experiments' production systems to encode information about a particle's origins. It proved impossible to satisfactorily reconcile these twin uses, and experiments' demands for particle provenance information have exceeded the capacity of an int (or even a long int).
Hence, barcodes are no longer available. Use attributes to provide additional information that was previously encoded using barcodes (see Attributes).
The unique identifier of particles and vertices is now called id() to separate its role from barcodes. Id is set automatically and cannot be changed. Id is not permanently attached to particle/vertex. When a particle or vertex is removed from the event, id's of other particles or vertices may change.
The Flow class has been removed, since it was unused by any widespread event generator, and to our knowledge the only active use-case is an abuse of it to provide more ints in which to encode provenance information. As this is now done via attributes, there is no case for Flow's continued existence. No backward compatibility Flow class is provided since this usage is extremely localised in one piece of user code and migration to the newer scheme should be simple.
The default units are set to GEV and MM. They can be provided as constructor parameters or changed later using HepMC::GenEvent::set_units
A lot of HepMC2 functions has been declared obsolete and are marked as deprecated. Warnings displayed at compilation time hint to what functions or classes should be used instead.
For the user convenience and backward compatibility the following standard attributes are supported for the
GenEvent
GenVertex
GenParticle
The presence of cycles in the event structure is indicated with an attribute
Note that attributes belong to the event, therefore these can be set only for particles and vertices that belong to a GenEvent object.
The most recent versions of HepMC3 has multiple implementations of the interfaces to HEPEVT Fortran common block. These are
define
directive. The block can hold float/double precision momenta. This implementation is not compiled into any library. All functions and variables are static, so only one instance of the interface can exists.A new class has been provided to store run-level information, such as weight names, names and description of tools used to generate the event, global attributes such as LHE run information or any other run information provided by user. See HepMC::GenRunInfo class description for details.
Attributes can be attached to GenEvent, GenParticle or GenVertex and they can have any format defined by the user (see Writing custom attributes). An attribute is accessed through a shared pointer and identified by its name.
Example of reading an attribute from the event:
Example of adding an attribute to the event:
Adding and getting attributes of a vertex or particle uses the same principles.
Note: An event (or particle or vertex) can have more than one attribute of the same type distinguished by different names. This might be useful in some applications, however, we strongly encourage to use just one instance named by its class name, as in these examples.
Any class that derives from HepMC::Attribute class can be used as an attribute that can be attached to the event, vertex or particle.
User has to provide two abstract methods from HepMC::Attribute used to parse the class content from/to string.
Example:
For other examples see attributes provided in the HepMC3 package.
The main HepMC3 library contains the classes for the I/O of multiple event formats.
Optionally the I/O capabilities can be implemented as plugin Reader/Writer classes compiled separately into dynamically loadable libraries and used via RearedPlugin and WriterPlugin classes. Please note that all required libraries/dlls should be loadable. See examples for details.
In some cases the fine tuning of the Reader/Writer classes behavior can be done using a map of string "options"
The options for ReaderAsciiHepMC2
The option for WriterAscii and WriterAsciiHepMC2
The relations between vertices and particles in GenEventData are encoded via members links1 and links2, wich are std::vector<int> containing object ids. Direct manipulations with links1 and links2 can be useful. For instance, when the events are saved in ROOT format, one can extract the information from links1 and links2 without reading whole event. In case links1[i] is particle, links2[i] is end vertex. In case links1[i] is vertex, links2[i] is outgoing particle. An example of usage is given below.
HepMC3 comes with an optional Search library for finding particles related to other particles or vertices. It provides a set of functions to perform simple search operations e.g.
and interfaces for a more advanced usage. For the latter two main interfaces are defined: Relatives, for finding a particular type of relative, and Feature, for generating filters based on Features extracted from particles. In addition, operator on Filters are also defined.
The Relatives interface is defined within search/include/HepMC3/Relatives.h. Classes that obey this interface must provide a set of operator functions that take either a GenParticlePtr, ConstGenParticlePtr, GenVertexPtr or ConstGenVertexPtr and return a vector of either GenParticlePtr or ConstGenParticlePtr. Note that the versions of the operator functions that receive a consted input parameter also return a vector<ConstGenParticlePtr>, while the versions that take non-consted input return non-const output. This ensures consistency with the rule that consted objects may only return pointers to const objects.
The Relatives base class is abstract, and has a concrete implementation in the templated RelativesInterface class. The RelativesInterface uses type erasure so that any class that obeys the defined Relatives interface (i.e that has the necessary four operator functions) can be wrapped in the RelativesInterface without needing to inherit from Relatives directly.
For example, if class foo implements the four necessary functions then the following will work
The purpose of Relatives is to be able to wrap any viable class in a common interface for finding relatives from a particle or vertex. Examples are provided in the form of the _parents and _children classes. These do not inherit from Raltives, but they do implement the necessary functions. The _parents and _children class are not intended to be used directly, but they are aliased by wrapping in the RelativesInterface:
Note as well that the _parents and _children classes use some utility aliases to help declare the appropriately consted return type. For example
has a return type GenParticles_type that is a vector of GenParticlePtr that is consted if GenObject_type is const, but is not consted if GenObject_type is not consted. Note as well the use of enable_if so that a single implementation can be used for both the const and non-const version of the functions. For the simple case of _parents the four required funcs could have been implemented directly without such templating, but for more complicated relatives it avoids duplicated code.
In addition to the RelativesInterface wrapper, Relatives.h also contains a Recursive class that can wrap the underlying relation in recursion. For example, recursion applied to the parents relationship provides all of the ancestors, i.e. parents repeatedly applied to the output of parents. The only additional requirement to use the Recursive wrapper is that the underlying class must implement a vertex(GenParticlePtr) method that returns the appropriate vertex to follow from a given GenParticle. As long as a class has such a method, it is possible to make a recursive version of it
The Relatives class contains static implementations of the Parents, Children, Ancestors and Descendants relatives, which can be accessed and used as follows
A Filter is any object that has an operator that takes as input a ConstGenParticlePtr and returns a bool that reflects whether the input particle passes the filter requirements or not. Filter is defined in Filter.h as
Filter.h also contains some logical operators that allow filters to be combined to create new filters, for example
Filter.h additionally contains a dummy filter that always accepts every possible particle. This may be needed in functions that require a default filter. The dummy filter is accessed as
It is possible to define a Filter by hand. However, there are some utility classes to define Filters based on features that can be obtained from GenParticles
The Feature interface is defined in Feature.h. The interface is templated on a Feature_type that is any type that can be extracted from a GenParticle. This is very flexible, and the only criteria is that the Feature must have the set of comparison operators. While the templated Feature is general enough to be used with any type of Feature, there are specialisations for both integral and floating point features. The specialisations will cover the vast majority of Features that are likely to be useful, although note that Attributes may be a source of more complicated Features.
To create a Feature, one need only wrap a lambda expression in the Feature interface. For example, to create a Feature based on particle status or pT:
The more general form for any type of Feature would be
Having created a Feature, it can be used to create Filters for particle selection. Applying operators to Features creates the Filter, which is a functor that evaluates on a particle. For example
It is also possible to make a new Feature from the absolute value of a previous Feature, e.g.
Some standard features are contained within the non-templated Selector class
Selector is a simplified interface that contains some predefined Features that can be used to search. Selector defines comparisons operators for both integral and floating point types, as well as the following selection Features:
So, for example, a filter can be defined as follows
As with Feature, it is possible to take tbe absolute value of a Selector. However, note that while Featue is templated, Selector is abstract and so it is not possible for abs() to return a Selector object directly, only a pointer
Note that the ATTRIBUTE selection is different from the others and does not have the full set of comparison operators. This is a current limitation of the Attributes, which are not guaranteed to offer all comparisons. ATTRIBUTE takes a string, which is the name of the attribute, and permits the equality operator and the method exists, which checks if the attribute is even present
The function applyFilter is used to apply the Filter to a set of particles. See for example examples/BasicExamples/basic_tree.cc
HepMC3 includes Python bindings codes suitable for compilation of python modules for Python3.
The binding codes were generated automatically using the binder utility version 1.4.0 created by Sergey Lyskov (Johns Hopkins University) et al.. See
The binding codes use the pybind11 library version 2.6.0 by Wenzel Jakob, EPFL's School of Computer and Communication Sciences. See
The Python bindings together with the HepMC3 itself can be installed from PyPy and multiple other repositories. Please see HepMC3 page at CERN GitLab for details.
To turn on the compilation of bindings use -DHEPMC3_ENABLE_PYTHON = ON. By default the python modules will be generated for python3 if these are found in the system. In case the test suite is enabled, tests of python bindings with all the enabled versions will run as well.
Despite not recommended, it should be possible to compile the python bindings using the installed version of HepMC3. To do this, copy the python directory outside of source tree, uncomment #project(pyHepMC3 CXX) in python/CMakeLists.txt and run CMake inside python directory with -DUSE_INSTALLED_HEPMC3=ON option.
In general, the syntax used in the python bindings is exactly the same as in the C++ code. However, some C++ classes and routines are don't have their Python equivalents, namsly:
deduce_reader
was binded and uses the internal python compression libraries, i.e. it has no dependence on zlib, zstd etc. The only requirement is that the corresponding module is available.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 perfected 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 structure 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), MCatNLO, 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)<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 errors and warnings since HepMC 3.3.0 have the following categories:
The categories are numbered approximatelly according to their importance and if the current warning/error level is set below the warning level of the category, the warnigs/errors from the category will not be printed.
The warning level can be obtained/set with int Setup::warnings_level()/void Setup::set_warnings_level(const int)
functions. The functions bool Setup::print_warnings()/void Setup::set_print_errors(const bool flag)
obtain/set the global flag for printing the warnings but don't change the warning levels. The default warning level is 750.
The error level can be obtained/set with int Setup::errors_level()/void Setup::set_errors_level(const int)
functions. The functions bool Setup::print_errors()/void Setup::set_print_errors(const bool flag)
obtain/set the global flag for printing the errors but don't change the warning levels. The default error level is 1000.