![]() |
Project
|
#include <AODToHepMC.h>
Public Types | |
The containers we subscribe to | |
using | Headers = o2::aod::McCollisions |
using | Header = o2::aod::McCollision |
using | Tracks = o2::aod::McParticles |
using | Track = typename Tracks::iterator |
using | XSections = o2::aod::HepMCXSections |
using | PdfInfos = o2::aod::HepMCPdfInfos |
using | HeavyIons = o2::aod::HepMCHeavyIons |
using | XSection = typename XSections::iterator |
using | PdfInfo = typename PdfInfos::iterator |
using | HeavyIon = typename HeavyIons::iterator |
Types from HepMC3 | |
using | FourVector = HepMC3::FourVector |
using | ParticlePtr = HepMC3::GenParticlePtr |
using | VertexPtr = HepMC3::GenVertexPtr |
using | Event = HepMC3::GenEvent |
using | HeavyIonPtr = HepMC3::GenHeavyIonPtr |
using | CrossSectionPtr = HepMC3::GenCrossSectionPtr |
using | PdfInfoPtr = HepMC3::GenPdfInfoPtr |
using | ParticleMap = std::unordered_map< long, ParticlePtr > |
using | ParticleVector = std::vector< ParticlePtr > |
using | ParticleList = std::list< ParticlePtr > |
using | VertexVector = std::vector< VertexPtr > |
using | VertexList = std::list< VertexPtr > |
using | WriterAscii = HepMC3::WriterAscii |
using | WriterAsciiPtr = std::shared_ptr< WriterAscii > |
Public Member Functions | |
Interface member functions | |
virtual void | init () |
virtual void | startEvent () |
virtual void | process (Header const &collision, Tracks const &tracks) |
virtual void | process (Header const &collision, XSections const &xsections, PdfInfos const &pdfs, HeavyIons const &heavyions) |
virtual void | endEvent () |
virtual bool | postRun () |
Public Attributes | ||
struct { | ||
std::string dump {""} | ||
bool onlyGen {false} | ||
bool useTree {false} | ||
int precision {8} | ||
bool recenter {false} | ||
} | configs | |
HepMC3 objects | ||
Event | mEvent | |
CrossSectionPtr | mCrossSec = nullptr | |
HeavyIonPtr | mIon = nullptr | |
PdfInfoPtr | mPdf = nullptr | |
Containers etc. | ||
ParticleMap | mParticles | |
VertexList | mVertices | |
Cache of particles. | ||
ParticleVector | mBeams | |
Cache of vertices. | ||
ParticleList | mOrphans | |
Cache of beam particles. | ||
Options and such | ||
Cache of particles w/o mothers | ||
WriterAsciiPtr | mWriter = nullptr | |
int | mEventNo = 0 | |
int | mLastBC = -1 | |
bool | mOnlyGen = false | |
bool | mUseTree = true | |
std::ofstream * | mOutput = nullptr | |
int | mPrecision = 16 | |
bool | mRecenter = false | |
Protected Member Functions | |
void | enableDump (const std::string &dump) |
Actual worker member functions | |
virtual void | makeEvent (Header const &collision, Tracks const &tracks) |
virtual void | makeHeader (Header const &header) |
virtual void | makeXSection (XSections const &xsection) |
virtual void | makePdfInfo (PdfInfos const &pdf) |
virtual void | makeHeavyIon (HeavyIons const &heavyion, Header const &header) |
virtual void | makeBeams (Header const &header, const VertexPtr ip) |
virtual void | makeParticles (Tracks const &tracks) |
virtual ParticlePtr | getParticle (Track const &ref) const |
virtual bool | isIgnored (Track const &track) const |
virtual ParticlePtr | makeParticleRecursive (Track const &track, Tracks const &tracks, bool force=false) |
virtual ParticlePtr | makeParticle (const Track &track, Int_t &motherStatus, bool force) const |
virtual VertexPtr | makeVertex (const Track &track) |
virtual void | fleshOut (Tracks const &tracks) |
virtual void | fleshOutParticle (Track const &track, Tracks const &tracks) |
virtual void | recenter (Header const &collision) |
Convert AOD tables of MC information into a HepMC event structure.
The event structure is kept in memory.
The conventions used here are
o2::aod::McParticles
table (aliased as AODToHepMC::Tracks), and correspond to a track used during simulation transport. These are of type o2::aod::McParticle
aliased to AODToHepMC::Track.Track::producedByGenerator()
HepMC3::GenParticle
.HepMC3::GenVertex
.last - first + 1if both last and first are zero or greater.
An event header is the information stored in a row of o2::aod::McCollisions
table (aliased to AODToHepMC::Headers). Each event has one such header (aliased as AODToHepMC::Header).
The header stores information on the event such as event number, weight, primary interaction point (IP), and so on.
In addition, auxiliary information (impact parameter,
o2::aod::HepMCXSections
(aliased to AODToHepMC::XSections) stores the event cross section and weight.o2::aod::HepMCPdfInfos
(aliased to AODToHepMC::PdfInfos) stores the parton distribution function parameters used in the event.o2::aod::HepMCHeavyIons
(aliased to AODToHepMC::PdfHeavyIons) stores the heavy-ion auxiliary information, such as the event impact parameter A HepMC event has a simple header which contains the event number, number of particles, and number of vertexes in the event.
Other event information is stored in specific structures.
HepMC3::GenCrossSection
(aliased to AODToHepMC::CrossSectionPtr) object. The information to fill into that is read from AODToHepMC::XSections tableHepMC3::GenPdfInfo
(aliased to AODToHepMC::PdfInfoPtr) object. The information to fill into that is read from AODToHepMC::PdfInfos table.HepMC3::GenHeavyIon
(aliased to AODToHepMC::HeavyIonPtr) object. The information to fill into that is read from AODToHepMC::HeavyIons table.The conversion is done as follows:
--hepmc-only-generated
), i.e., tracks that were posted by the event generator, and this track is not such a track (AODToHepMC::isIgnored), nor forced (argument force
) then return nothing. Note that mothers of particles that we keep are always made (force
is set to true).--hepmc-recenter
--hepmc-dumb
filename) we may write the events to a file (AODToHepMC::mOutput). The event is still kept in memory.The utility o2-aod-mc-to-hepmc
will read in AODs and write out a HepMC event file (plain ASCII).
Definition at line 253 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::CrossSectionPtr = HepMC3::GenCrossSectionPtr |
Alias (smart-)pointer to HepMC cross section object
Definition at line 317 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::Event = HepMC3::GenEvent |
Alias HepMC eventt structure
Definition at line 313 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::FourVector = HepMC3::FourVector |
Alias HepMC four-vector
Definition at line 307 of file AODToHepMC.h.
Alias of MC collisions table row
Definition at line 282 of file AODToHepMC.h.
Alias of MC collisions table
Definition at line 280 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::HeavyIon = typename HeavyIons::iterator |
Alias row of auxiliary MC table of heavy-ion parameters
Definition at line 300 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::HeavyIonPtr = HepMC3::GenHeavyIonPtr |
Alias (smart-)pointer to HepMC heavy-ion object
Definition at line 315 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::HeavyIons = o2::aod::HepMCHeavyIons |
Alias auxiliary MC table of heavy-ion parameters
Definition at line 293 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::ParticleList = std::list<ParticlePtr> |
A container of pointers to particles
Definition at line 326 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::ParticleMap = std::unordered_map<long, ParticlePtr> |
Type used to map tracks to HepMC particles
Definition at line 322 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::ParticlePtr = HepMC3::GenParticlePtr |
Alias (smart-)pointer to HepMC particle
Definition at line 309 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::ParticleVector = std::vector<ParticlePtr> |
A container of pointers to particles
Definition at line 324 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::PdfInfo = typename PdfInfos::iterator |
Alias row of auxiliary MC table of parton distribution function parameters
Definition at line 298 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::PdfInfoPtr = HepMC3::GenPdfInfoPtr |
Alias (smart-)pointer to HepMC parton distribution function object
Definition at line 320 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::PdfInfos = o2::aod::HepMCPdfInfos |
Alias auxiliary MC table of parton distribution function parameters
Definition at line 291 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::Track = typename Tracks::iterator |
Alias MC particles (tracks) table row
Definition at line 286 of file AODToHepMC.h.
Alias MC particles (tracks) table
Definition at line 284 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::VertexList = std::list<VertexPtr> |
A container of pointers to vertexes
Definition at line 330 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::VertexPtr = HepMC3::GenVertexPtr |
Alias (smart-)pointer to HepMC vertex
Definition at line 311 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::VertexVector = std::vector<VertexPtr> |
A container of pointers to vertexes
Definition at line 328 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::WriterAscii = HepMC3::WriterAscii |
Alias of HepMC writer class
Definition at line 332 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::WriterAsciiPtr = std::shared_ptr<WriterAscii> |
The of pointer to HepMC writer class
Definition at line 334 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::XSection = typename XSections::iterator |
Alias row of auxiliary MC table of cross-sections
Definition at line 295 of file AODToHepMC.h.
using o2::eventgen::AODToHepMC::XSections = o2::aod::HepMCXSections |
Alias auxiliary MC table of cross-sections
Definition at line 288 of file AODToHepMC.h.
|
protected |
Open dump output, or close if an empty string was given.
dump |
Definition at line 561 of file AODToHepMC.cxx.
|
virtual |
Call after process an. Thisf finalises the event and optionally outputs to dump.
Definition at line 71 of file AODToHepMC.cxx.
Visit all tracks, but this time look for daughters and attach mothers as incoming particles if not already done.
Definition at line 436 of file AODToHepMC.cxx.
|
protectedvirtual |
Flesh out a single particle
track | The track for which we should flesh out the corresponding particle. |
Definition at line 443 of file AODToHepMC.cxx.
|
protectedvirtual |
Get particle corresponding to track no from particle cache
ref | Track reference |
Definition at line 292 of file AODToHepMC.cxx.
|
virtual |
Initialize the converter. Sets internal parameters based on the configurables.
Definition at line 21 of file AODToHepMC.cxx.
|
protectedvirtual |
Check if we are ignoring this track
track | Track to check |
Definition at line 302 of file AODToHepMC.cxx.
|
inlineprotectedvirtual |
This is supposed to make the beam particles from the information available. However, it seems like we really don't have enough information, so for now we will do nothing here. Perhaps the user will be forced to provide that information - either via the analyser configurables or somehow from somewhere else.
Definition at line 475 of file AODToHepMC.h.
|
protectedvirtual |
Generate the final event, including fleshing out the vertexes, and so on
collision | Header information |
tracks | Particle tracks |
Definition at line 82 of file AODToHepMC.cxx.
Set the various fields in the header of the HepMC3::GenEvent object
header | Header object |
Definition at line 130 of file AODToHepMC.cxx.
|
protectedvirtual |
Make heavy-ion collision information. If no header given, then fill in other reasonable values
Definition at line 204 of file AODToHepMC.cxx.
|
protectedvirtual |
Generate a HepMC particle from a track. Note that the job here is simply to make the object. The more complicated job of adding the track to the tree is done above in makeParticleRecursive
track | MC track |
mst | Mother status code - updated on return |
force | Force generation of particle even if onlyGen |
Definition at line 398 of file AODToHepMC.cxx.
|
protectedvirtual |
Truely make a particle, and its mother particles if any. We add vertexes as needed here too.
Note that this traverses the tree from the bottom up. That is, we check if a particle has any mothers, and if it does, we create those.
However, this can be a bit problematic, since the Kinematic tree (and thus the McParticles table) only allows for at most 2 mothers. In the case a vertex has 3 or more incoming particles, then some of the intermediate particles will be lost.
We remedy that problem by traversing the tree again, but this time from the bottom up - that is, we look for daughters of all particles, and if they are not registered with the out-going vertex, then we reattch the parent to the incoming vertex of the daughters. For this to work, we need to map from track index to HepMC3::GenParticle.
track | Current track |
tracks | MC tracks |
force | If true, do make the particle even if it would otherwise be skipped. |
Definition at line 313 of file AODToHepMC.cxx.
Make all particles. We loop through the MC particles, and for each of them create that particle and any possible mother particles (recursively). This allows us to traverse the data rather straight-forwardly.
tracks | The MC tracks |
Definition at line 282 of file AODToHepMC.cxx.
Make parton-distribition function information. If no entry in the table, then make dummy information
Definition at line 175 of file AODToHepMC.cxx.
|
protectedvirtual |
Generate vertex from production vertex of track
track | MC track |
Definition at line 591 of file AODToHepMC.cxx.
Make cross-section information. If no entry in the table, then make dummy information
Definition at line 151 of file AODToHepMC.cxx.
|
inlinevirtual |
End of run - closes output file if enabled. This is called via specialisation of o2::framework::OutputManager<AODToHepMC>.
Definition at line 426 of file AODToHepMC.h.
Process the collision header and tracks
collision | Header information |
tracks | Particle tracks |
Definition at line 47 of file AODToHepMC.cxx.
|
virtual |
Process collision header and HepMC auxiliary information
collision | Header information |
xsections | Cross-section table (possible no rows) |
pdfs | Parton-distribution function table (possible no rows) |
heavyions | Heavy ion collision table (possible no rows) |
Definition at line 60 of file AODToHepMC.cxx.
Recenter event to (0,0,0,0). This will use the vertex information from the event header to translate all vertexes in the event.
header | Event header |
|
virtual |
Call before starting to process an event. This clears the current event and internal data structures.
Definition at line 37 of file AODToHepMC.cxx.
struct { ... } o2::eventgen::AODToHepMC::configs |
Group of configurables which will be added to an program that uses this class. Note that it is really the specialisation of framework::OptionManager<AODToHepMC> that propagates the options to the program.
std::string o2::eventgen::AODToHepMC::dump {""} |
Option for dumping HepMC event structures to disk. Takes one argument - the name of the file to write to.
Definition at line 263 of file AODToHepMC.h.
ParticleVector o2::eventgen::AODToHepMC::mBeams |
CrossSectionPtr o2::eventgen::AODToHepMC::mCrossSec = nullptr |
Pointer to cross section-ion information
Definition at line 343 of file AODToHepMC.h.
Event o2::eventgen::AODToHepMC::mEvent |
The result of processing
Definition at line 341 of file AODToHepMC.h.
int o2::eventgen::AODToHepMC::mEventNo = 0 |
Current sequential event number
Definition at line 369 of file AODToHepMC.h.
HeavyIonPtr o2::eventgen::AODToHepMC::mIon = nullptr |
Pointer to heavy-ion information
Definition at line 345 of file AODToHepMC.h.
int o2::eventgen::AODToHepMC::mLastBC = -1 |
The last bunch crossing identifier
Definition at line 371 of file AODToHepMC.h.
bool o2::eventgen::AODToHepMC::mOnlyGen = false |
If true, only store particles from the generator
Definition at line 373 of file AODToHepMC.h.
ParticleList o2::eventgen::AODToHepMC::mOrphans |
std::ofstream* o2::eventgen::AODToHepMC::mOutput = nullptr |
Output stream if enabled
Definition at line 377 of file AODToHepMC.h.
ParticleMap o2::eventgen::AODToHepMC::mParticles |
Maps tracks to particles
Definition at line 354 of file AODToHepMC.h.
PdfInfoPtr o2::eventgen::AODToHepMC::mPdf = nullptr |
Pointer to parton distribution function information
Definition at line 347 of file AODToHepMC.h.
int o2::eventgen::AODToHepMC::mPrecision = 16 |
Precision used on the output stream
Definition at line 379 of file AODToHepMC.h.
bool o2::eventgen::AODToHepMC::mRecenter = false |
If true, recenter IP to (0,0,0,0)
Definition at line 381 of file AODToHepMC.h.
bool o2::eventgen::AODToHepMC::mUseTree = true |
If true, use HepMC tree parser
Definition at line 375 of file AODToHepMC.h.
VertexList o2::eventgen::AODToHepMC::mVertices |
WriterAsciiPtr o2::eventgen::AODToHepMC::mWriter = nullptr |
Output writer, if enabled
Definition at line 367 of file AODToHepMC.h.
bool o2::eventgen::AODToHepMC::onlyGen {false} |
Option for only storing particles from the event generator. Note, if a particle is stored down, then its mothers will also be stored.
Definition at line 267 of file AODToHepMC.h.
int o2::eventgen::AODToHepMC::precision {8} |
Floating point precision used when writing to disk
Definition at line 271 of file AODToHepMC.h.
void o2::eventgen::AODToHepMC::recenter {false} |
Recenter event at IP=(0,0,0,0).
Definition at line 273 of file AODToHepMC.h.
bool o2::eventgen::AODToHepMC::useTree {false} |
Use HepMC's tree parsing for building event structure
Definition at line 269 of file AODToHepMC.h.