![]() |
Project
|
#include <GeneratorPythia8.h>
Inherits o2::eventgen::Generator.
Public Member Functions | |
GeneratorPythia8 () | |
GeneratorPythia8 (Pythia8GenConfig const &) | |
GeneratorPythia8 (const Char_t *name, const Char_t *title="ALICEo2 Pythia8 Generator") | |
~GeneratorPythia8 () override=default | |
methods to override | |
Bool_t | Init () override |
Bool_t | generateEvent () override |
Bool_t | importParticles () override |
setters | |
void | setConfig (std::string val) |
void | setHooksFileName (std::string val) |
void | setHooksFuncName (std::string val) |
void | setUserHooks (Pythia8::UserHooks *hooks) |
Configuration methods | |
bool | readString (std::string val) |
bool | readFile (std::string val) |
![]() | |
Generator () | |
Generator (const Char_t *name, const Char_t *title="ALICEo2 Generator") | |
~Generator () override=default | |
Bool_t | Init () override |
Bool_t | ReadEvent (FairPrimaryGenerator *primGen) final |
Bool_t | triggerEvent () |
void | setMomentumUnit (double val) |
void | setEnergyUnit (double val) |
void | setPositionUnit (double val) |
void | setTimeUnit (double val) |
void | setBoost (Double_t val) |
void | setTriggerMode (ETriggerMode_t val) |
void | addTrigger (Trigger trigger) |
void | addDeepTrigger (DeepTrigger trigger) |
const std::vector< TParticle > & | getParticles () const |
void | clearParticles () |
virtual void | notifyEmbedding (const o2::dataformats::MCEventHeader *eventHeader) |
void | setTriggerOkHook (std::function< void(std::vector< TParticle > const &p, int eventCount)> f) |
void | setTriggerFalseHook (std::function< void(std::vector< TParticle > const &p, int eventCount)> f) |
Protected Types | |
Some function definitions | |
typedef UserFilterFcn | Select |
typedef std::vector< int >(* | GetRelatives) (const Pythia8::Particle &) |
typedef void(* | SetRelatives) (Pythia8::Particle &, int, int) |
typedef std::pair< int, int >(* | FirstLastRelative) (const Pythia8::Particle &) |
Protected Member Functions | |
void | seedGenerator () |
performs seeding of the random state of Pythia (called from Init) | |
void | initUserFilterCallback () |
ClassDefOverride (GeneratorPythia8, 1) | |
Methods that can be overridded | |
void | updateHeader (o2::dataformats::MCEventHeader *eventHeader) override |
Internal methods | |
Bool_t | importParticles (Pythia8::Event &event) |
void | pruneEvent (Pythia8::Event &event, Select select) |
void | investigateRelatives (Pythia8::Event &event, const std::vector< int > &old2New, size_t index, std::vector< bool > &done, GetRelatives getter, SetRelatives setter, FirstLastRelative firstLast, const std::string &what, const std::string &ind="") |
![]() | |
Generator (const Generator &) | |
Generator & | operator= (const Generator &) |
Bool_t | addTracks (FairPrimaryGenerator *primGen) |
Bool_t | boostEvent () |
void | addSubGenerator (int subGeneratorId, std::string const &subGeneratorDescription) |
void | notifySubGenerator (int subGeneratorId) |
Protected Attributes | |
Pythia8::Pythia | mPythia |
UserFilterFcn | mUserFilterFcn = [](Pythia8::Particle const&) -> bool { return true; } |
bool | mApplyPruning = false |
bool | mIsInitialized = false |
long | mInitialRNGSeed = -1 |
Pythia8GenConfig | mGenConfig |
int | mThisPythia8InstanceID = 0 |
Configurations | |
std::string | mConfig |
std::string | mHooksFileName |
std::string | mHooksFuncName |
![]() | |
void * | mInterface = nullptr |
std::string | mInterfaceName |
ETriggerMode_t | mTriggerMode = kTriggerOFF |
std::vector< Trigger > | mTriggers |
std::vector< DeepTrigger > | mDeepTriggers |
std::function< void(std::vector< TParticle > const &p, int eventCount)> | mTriggerOkHook = [](std::vector<TParticle> const& p, int eventCount) {} |
std::function< void(std::vector< TParticle > const &p, int eventCount)> | mTriggerFalseHook = [](std::vector<TParticle> const& p, int eventCount) {} |
int | mReadEventCounter = 0 |
double | mMomentumUnit = 1. |
double | mEnergyUnit = 1. |
double | mPositionUnit = 0.1 |
double | mTimeUnit = 3.3356410e-12 |
std::vector< TParticle > | mParticles |
Double_t | mBoost |
int | mThisInstanceID = 0 |
Static Protected Attributes | |
static std::atomic< int > | Pythia8InstanceCounter |
static constexpr long | MAX_SEED = 900000000 |
![]() | |
static std::atomic< int > | InstanceCounter {0} |
Utilities | |
typedef std::function< bool(const Pythia8::Particle &)> | UserFilterFcn |
void | getNcoll (int &nColl) |
void | getNpart (int &nPart) |
void | getNpart (int &nProtonProj, int &nNeutronProj, int &nProtonTarg, int &nNeutronTarg) |
void | getNremn (int &nProtonProj, int &nNeutronProj, int &nProtonTarg, int &nNeutronTarg) |
void | getNfreeSpec (int &nFreenProj, int &nFreepProj, int &nFreenTarg, int &nFreepTarg) |
bool | setInitialSeed (long seed) |
GeneratorPythia8 (const GeneratorPythia8 &) | |
GeneratorPythia8 & | operator= (const GeneratorPythia8 &) |
void | selectFromAncestor (int ancestor, Pythia8::Event &inputEvent, Pythia8::Event &outputEvent) |
void | getNcoll (const Pythia8::Info &info, int &nColl) |
void | getNpart (const Pythia8::Info &info, int &nPart) |
void | getNpart (const Pythia8::Info &info, int &nProtonProj, int &nNeutronProj, int &nProtonTarg, int &nNeutronTarg) |
void | getNremn (const Pythia8::Event &event, int &nProtonProj, int &nNeutronProj, int &nProtonTarg, int &nNeutronTarg) |
void | getNfreeSpec (const Pythia8::Info &info, int &nFreenProj, int &nFreepProj, int &nFreenTarg, int &nFreepTarg) |
Additional Inherited Members | |
![]() | |
enum | ETriggerMode_t { kTriggerOFF , kTriggerOR , kTriggerAND } |
![]() | |
static void | setTotalNEvents (unsigned int &n) |
static unsigned int | getTotalNEvents () |
Interface to the Pythia8 event generator. This generator is configured by configuration files (e.g., Generators/share/egconfig/pythia8_inel.cfg
for the type of events to generate.
The above file, for example, contains
* ### Beams, proton-proton collisions at sqrt(s)=14TeV * Beams:idA 2212 # proton * Beams:idB 2212 # proton * Beams:eCM 14000. # GeV * * ### Processes, min-bias inelastic events * SoftQCD:inelastic on # all inelastic processes * * ### Decays. Only decay until 1cm/c - corresponding to (physical) primaries * ParticleDecays:limitTau0 on * ParticleDecays:tau0Max 10. *
The configuration file to read is initially set in GeneratorFactory
, but an additional configuration file can be specified with the configuration key GeneratorPythia8::config
.
If the configuration key GeneratorPythia8::includePartonEvent
is set to false (default), then the event is pruned. That is, all particles that are not
are removed from the event before passing on to the simulation. The event structure is kept, so that we have a well-formed event structure. This reduces the event size by roughly 30%.
If the configuration key GeneratorPythia8::includePartonEvent
is true, then the full event is kept, including intermediate partons such as gluons, pomerons, diffractive hadrons, and so on.
In the future, the way to prune events may become more flexible. For example, we could have the configuration keys
The configuration key GeneratorPythia8::hooksFileName
allows the definition of a Pythia8 user hook. See for example Generators/share/egconfig/pythia8_userhooks_charm.C
. The file specified is interpreted via ROOT (i.e., a ROOT script), and the function name set via the configuration key GeneratorPythia8::hooksFuncName
(default pythia8_user_hooks
) is executed. That function must return a pointer to a Pythia8::UserHooks
object (see the Pythia8 manual for more on this).
Definition at line 86 of file GeneratorPythia8.h.
|
protected |
Get range of relatives (mothers or daughters) of a particle
Definition at line 191 of file GeneratorPythia8.h.
|
protected |
Get relatives (mothers or daughters) of a particle
Definition at line 187 of file GeneratorPythia8.h.
|
protected |
Select particles when pruning event
Definition at line 185 of file GeneratorPythia8.h.
|
protected |
Set relatives (mothers or daughters) of a particle
Definition at line 189 of file GeneratorPythia8.h.
typedef std::function<bool(const Pythia8::Particle&)> o2::eventgen::GeneratorPythia8::UserFilterFcn |
Definition at line 167 of file GeneratorPythia8.h.
o2::eventgen::GeneratorPythia8::GeneratorPythia8 | ( | ) |
default constructor
Definition at line 52 of file GeneratorPythia8.cxx.
o2::eventgen::GeneratorPythia8::GeneratorPythia8 | ( | Pythia8GenConfig const & | config | ) |
o2::eventgen::GeneratorPythia8::GeneratorPythia8 | ( | const Char_t * | name, |
const Char_t * | title = "ALICEo2 Pythia8 Generator" |
||
) |
|
overridedefault |
destructor
|
protected |
copy constructor
|
protected |
|
overridevirtual |
Generate a single event
generate event
[NOTE] The issue with large particle production vertex when running Pythia8 heavy-ion model (Angantyr) is solved in Pythia 8.3 series. For discussions about this issue, please refer to this JIRA ticket https://alice.its.cern.ch/jira/browse/O2-1382. The code remains within preprocessor directives, both for reference and in case future use demands to roll back to Pythia 8.2 series.
As we have inhibited all hadron decays before init, the event generation stops after hadronisation. We then pick all particles from here and force their production vertex to be (0,0,0,0). Afterwards we process the decays.
force production vertices to (0,0,0,0)
proceed with decays
success
Implements o2::eventgen::Generator.
Definition at line 239 of file GeneratorPythia8.cxx.
Get number of binary collisions. Note that this method deviates from how the Pythia authors count number of binary collisions
compute number of collisions from sub-collision information
Definition at line 825 of file GeneratorPythia8.cxx.
Get number of binary collisions. Note that this method deviates from how the Pythia authors count number of binary collisions
Definition at line 137 of file GeneratorPythia8.h.
|
protected |
Get number of spectators, split by nucleon type and origin. Note that this method deviates from how the Pythia authors count number of spectators
compute number of free spectator nucleons for ZDC response
Definition at line 1011 of file GeneratorPythia8.cxx.
|
inline |
Get number of spectators, split by nucleon type and origin. Note that this method deviates from how the Pythia authors count number of spectators
Definition at line 162 of file GeneratorPythia8.h.
Get number of participants. Note that this method deviates from how the Pythia authors count number of participants
compute number of participants as the sum of all participants nucleons
Definition at line 870 of file GeneratorPythia8.cxx.
|
protected |
Get number of participants, split by nucleon type and origin. Note that this method deviates from how the Pythia authors count number of participants
compute number of participants from sub-collision information
Definition at line 899 of file GeneratorPythia8.cxx.
Get number of participants. Note that this method deviates from how the Pythia authors count number of participants
Definition at line 143 of file GeneratorPythia8.h.
|
inline |
Get number of participants, split by nucleon type and origin. Note that this method deviates from how the Pythia authors count number of participants
Definition at line 150 of file GeneratorPythia8.h.
|
protected |
Get number of nuclei remnants, split by nucleon type and origin.
compute number of spectators from the nuclear remnant of the beams
Definition at line 960 of file GeneratorPythia8.cxx.
|
inline |
Get number of nuclei remnants, split by nucleon type and origin.
Definition at line 155 of file GeneratorPythia8.h.
|
inlineoverridevirtual |
Import particles from Pythia onto the simulation event stack
Implements o2::eventgen::Generator.
Definition at line 106 of file GeneratorPythia8.h.
|
protected |
Import particles from Pythia onto the simulation stack
import particles
success
Definition at line 616 of file GeneratorPythia8.cxx.
|
override |
Initialize the generator if needed
init
init base class
Seed the Pythia random number state. The user may override this seeding by providing separate Random:setSeed configurations in the configuration file.
read configuration
user hooks via configuration macro
[NOTE] The issue with large particle production vertex when running Pythia8 heavy-ion model (Angantyr) is solved in Pythia 8.3 series. For discussions about this issue, please refer to this JIRA ticket https://alice.its.cern.ch/jira/browse/O2-1382. The code remains within preprocessor directives, both for reference and in case future use demands to roll back to Pythia 8.2 series.
inhibit hadron decays
initialise
success
Definition at line 140 of file GeneratorPythia8.cxx.
|
protected |
Definition at line 593 of file GeneratorPythia8.cxx.
|
protected |
Investigate relatives (mothers or daughters) for particles to keep when pruning an event. This checks the current particle, identified by index, if any of its relatives (either mothers or daughters) are to be kept. If a relative is to be kept, then that is added to an internal list. If a relative is not to be kept, then that relatives relatives are queried (recursive call). The result of the recursive call is a list of relatives to the current particle which are to be kept. These are then also added to the internal list. The relatives that are found to be kept are then set to be relatives of the current particle. Note that this member function modifies the relatives of the passed particle, and thus modifies the passed event structure. Calculations are cached.
Definition at line 280 of file GeneratorPythia8.cxx.
|
protected |
operator=
Prune an event. Only particles for which the function select returns true are kept in the event record. The structure of the event is preserved, meaning that particles will point back to their ultimate (selected) mothers and, if select preserves the beam particles, the ultimate collision interaction.
Definition at line 399 of file GeneratorPythia8.cxx.
|
inline |
Read a Pythia8 configuration file
Definition at line 130 of file GeneratorPythia8.h.
|
inline |
Read a Pythia8 configuration string
Definition at line 128 of file GeneratorPythia8.h.
|
protected |
performs seeding of the random state of Pythia (called from Init)
Function is seeding the Pythia random numbers. In case a completely different logic is required by users, we could make this function virtual or allow to set/execute a user-given lambda function instead.
Note that this function is executed before the Pythia8 user config file is read. So the config file should either not contain seeding information ... or can be used to override seeding logic.
Definition at line 107 of file GeneratorPythia8.cxx.
|
protected |
Select from ancestor. Fills the output event with all particles related to an ancestor of the input event
Starting from ancestor, select all daughters (and their daughters recursively), and store them in the output event.
select from ancestor fills the output event with all particles related to an ancestor of the input event
Definition at line 783 of file GeneratorPythia8.cxx.
|
inline |
Set Pythia8 configuration file to read
Definition at line 112 of file GeneratorPythia8.h.
|
inline |
Set the ROOT script file name that defines a Pythia8::UserHooks object
Definition at line 116 of file GeneratorPythia8.h.
|
inline |
Function in ROOT script that returns a pointer to a Pythia8::UserHooks object
Definition at line 119 of file GeneratorPythia8.h.
bool o2::eventgen::GeneratorPythia8::setInitialSeed | ( | long | seed | ) |
A function allowing to set the initial value used in seeding Pythia. The function needs to be called before GeneratorPythia8::Init is invoked. The function value will be true upon success, or false if either Init has already been called or if the see is smaller than 0. For values of seed >= 0, a truncation to the range [0:90000000] will automatically take place via a modulus operation.
Definition at line 90 of file GeneratorPythia8.cxx.
void o2::eventgen::GeneratorPythia8::setUserHooks | ( | Pythia8::UserHooks * | hooks | ) |
Set the user hooks (defined in a Pythia8::UserHooks object) for the event generator.
Definition at line 227 of file GeneratorPythia8.cxx.
|
overrideprotectedvirtual |
Update the event header. This propagates all sorts of information from Pythia8 to the simulation event header, including parton distribution function parameters, event weight, cross-section information, heavy-ion collision information, and so on.
update header
set impact parameter
set event plane angle
set Ncoll, Npart and Nremn
Reimplemented from o2::eventgen::Generator.
Definition at line 682 of file GeneratorPythia8.cxx.
|
protected |
Definition at line 283 of file GeneratorPythia8.h.
|
staticconstexprprotected |
Definition at line 293 of file GeneratorPythia8.h.
|
protected |
configuration file to read
Definition at line 272 of file GeneratorPythia8.h.
|
protected |
Definition at line 288 of file GeneratorPythia8.h.
|
protected |
ROOT script defining a Pythia8::UserHooks object
Definition at line 274 of file GeneratorPythia8.h.
|
protected |
Function in mHooksFileName
to execute to return pointer to Pythia8::UserHooks object
Definition at line 277 of file GeneratorPythia8.h.
|
protected |
Definition at line 285 of file GeneratorPythia8.h.
|
protected |
Definition at line 284 of file GeneratorPythia8.h.
|
protected |
Pythia8
Definition at line 267 of file GeneratorPythia8.h.
|
protected |
Definition at line 291 of file GeneratorPythia8.h.
|
protected |
Definition at line 280 of file GeneratorPythia8.h.
|
staticprotected |
Definition at line 290 of file GeneratorPythia8.h.