-
Notifications
You must be signed in to change notification settings - Fork 496
Expand file tree
/
Copy pathGeneratorPythia8.h
More file actions
308 lines (280 loc) · 12.9 KB
/
GeneratorPythia8.h
File metadata and controls
308 lines (280 loc) · 12.9 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
// All rights not expressly granted are reserved.
//
// This software is distributed under the terms of the GNU General Public
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
//
// In applying this license CERN does not waive the privileges and immunities
// granted to it by virtue of its status as an Intergovernmental Organization
// or submit itself to any jurisdiction.
/// \author R+Preghenella - January 2020
#ifndef ALICEO2_EVENTGEN_GENERATORPYTHIA8_H_
#define ALICEO2_EVENTGEN_GENERATORPYTHIA8_H_
#include "Generators/Generator.h"
#include "Pythia8/Pythia.h"
#include <functional>
#include "Generators/GeneratorPythia8Param.h"
namespace o2
{
namespace eventgen
{
/*****************************************************************/
/*****************************************************************/
/** 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
*
* @verbatim
* ### 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.
* @endverbatim
*
* 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
*
* - beam particles (HepMC status = 4),
* - decayed particles (HepMC status = 2), nor
* - final state partcles (HepMC status = 1)
*
* 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
*
* - GeneratorPythia8::onlyStatus a list of HepMC status codes to accept
* - GeneratorPythia8::onlyPDGs a list of PDG particle codes to accept
*
* 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).
*/
class GeneratorPythia8 : public Generator
{
public:
/** default constructor **/
GeneratorPythia8();
/** constructor **/
GeneratorPythia8(Pythia8GenConfig const&);
/** constructor **/
GeneratorPythia8(const Char_t* name, const Char_t* title = "ALICEo2 Pythia8 Generator");
/** destructor **/
~GeneratorPythia8() override = default;
/** @{
@name methods to override **/
/** Initialize the generator if needed **/
Bool_t Init() override;
/** Generate a single event */
Bool_t generateEvent() override;
/** Import particles from Pythia onto the simulation event stack */
Bool_t importParticles() override { return importParticles(mPythia.event); };
/** @} */
/** @{
* @name setters **/
/** Set Pythia8 configuration file to read */
void setConfig(std::string val) { mConfig = val; };
/**
* Set the ROOT script file name that defines a Pythia8::UserHooks
* object */
void setHooksFileName(std::string val) { mHooksFileName = val; };
/** Function in ROOT script that returns a pointer to a
* Pythia8::UserHooks object */
void setHooksFuncName(std::string val) { mHooksFuncName = val; };
/** Set the user hooks (defined in a Pythia8::UserHooks object) for
* the event generator. */
void setUserHooks(Pythia8::UserHooks* hooks);
/** @} */
/** @{
* @name Configuration methods **/
/** Read a Pythia8 configuration string */
bool readString(std::string val) { return mPythia.readString(val, true); };
/** Read a Pythia8 configuration file */
bool readFile(std::string val) { return mPythia.readFile(val, true); };
/** @} */
/** @{
* @name Utilities **/
/** Get number of binary collisions. Note that this method deviates
* from how the Pythia authors count number of binary collisions */
void getNcoll(int& nColl)
{
getNcoll(mPythia.info, nColl);
};
/** Get number of participants. Note that this method deviates
* from how the Pythia authors count number of participants */
void getNpart(int& nPart)
{
getNpart(mPythia.info, nPart);
};
/** Get number of participants, split by nucleon type and origin.
* Note that this method deviates from how the Pythia authors count
* number of participants */
void getNpart(int& nProtonProj, int& nNeutronProj, int& nProtonTarg, int& nNeutronTarg)
{
getNpart(mPythia.info, nProtonProj, nNeutronProj, nProtonTarg, nNeutronTarg);
};
/** Get number of nuclei remnants, split by nucleon type and origin. */
void getNremn(int& nProtonProj, int& nNeutronProj, int& nProtonTarg, int& nNeutronTarg)
{
getNremn(mPythia.event, nProtonProj, nNeutronProj, nProtonTarg, nNeutronTarg);
};
/** Get number of spectators, split by nucleon type and origin.
* Note that this method deviates from how the Pythia authors count
* number of spectators */
void getNfreeSpec(int& nFreenProj, int& nFreepProj, int& nFreenTarg, int& nFreepTarg)
{
getNfreeSpec(mPythia.info, nFreenProj, nFreepProj, nFreenTarg, nFreepTarg);
};
typedef std::function<bool(const Pythia8::Particle&)> UserFilterFcn;
/// 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.
bool setInitialSeed(long seed);
protected:
/** copy constructor **/
GeneratorPythia8(const GeneratorPythia8&);
/** operator= **/
GeneratorPythia8& operator=(const GeneratorPythia8&);
/** @{
* @name Some function definitions
*/
/** Select particles when pruning event */
typedef UserFilterFcn Select;
/** Get relatives (mothers or daughters) of a particle */
typedef std::vector<int> (*GetRelatives)(const Pythia8::Particle&);
/** Set relatives (mothers or daughters) of a particle */
typedef void (*SetRelatives)(Pythia8::Particle&, int, int);
/** Get range of relatives (mothers or daughters) of a particle */
typedef std::pair<int, int> (*FirstLastRelative)(const Pythia8::Particle&);
/** @} */
/** @{
* @name Methods that can be overridded **/
/** 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. */
void updateHeader(o2::dataformats::MCEventHeader* eventHeader) override;
/** @} */
/** @{
* @name Internal methods **/
/** Import particles from Pythia onto the simulation stack */
Bool_t importParticles(Pythia8::Event& event);
/** 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. */
void pruneEvent(Pythia8::Event& event, Select select);
/** 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. */
void investigateRelatives(Pythia8::Event& event, // Event
const std::vector<int>& old2New, // Map from old to new idx
size_t index, // Current particle
std::vector<bool>& done, // cache flag
GetRelatives getter, // get relatives
SetRelatives setter, // set relatives
FirstLastRelative firstLast, // get first and last relative
const std::string& what, // what are we looking for
const std::string& ind = ""); // logging indent
/** @{
* @name utilities **/
/** 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.
**/
void selectFromAncestor(int ancestor, Pythia8::Event& inputEvent, Pythia8::Event& outputEvent);
/** Get number of binary collisions. Note that this method deviates
* from how the Pythia authors count number of binary collisions */
void getNcoll(const Pythia8::Info& info, int& nColl);
/** Get number of participants. Note that this method deviates
* from how the Pythia authors count number of participants */
void getNpart(const Pythia8::Info& info, int& nPart);
/** Get number of participants, split by nucleon type and origin.
* Note that this method deviates from how the Pythia authors count
* number of participants */
void getNpart(const Pythia8::Info& info, int& nProtonProj, int& nNeutronProj, int& nProtonTarg, int& nNeutronTarg);
/** Get number of nuclei remnants, split by nucleon type and origin. */
void getNremn(const Pythia8::Event& event, int& nProtonProj, int& nNeutronProj, int& nProtonTarg, int& nNeutronTarg);
/** Get number of spectators, split by nucleon type and origin.
* Note that this method deviates from how the Pythia authors count
* number of spectators */
void getNfreeSpec(const Pythia8::Info& info, int& nFreenProj, int& nFreepProj, int& nFreenTarg, int& nFreepTarg);
/** @} */
/// performs seeding of the random state of Pythia (called from Init)
void seedGenerator();
/** Pythia8 **/
Pythia8::Pythia mPythia; //!
/** @{
* @name Configurations */
/** configuration file to read **/
std::string mConfig;
/** ROOT script defining a Pythia8::UserHooks object */
std::string mHooksFileName;
/** Function in `mHooksFileName` to execute to return pointer to
* Pythia8::UserHooks object */
std::string mHooksFuncName;
/** @} */
UserFilterFcn mUserFilterFcn = [](Pythia8::Particle const&) -> bool { return true; };
void initUserFilterCallback();
bool mApplyPruning = false;
bool mIsInitialized = false; // if Init function has been called
long mInitialRNGSeed = -1; // initial seed for Pythia random number state;
// will be transported to Pythia in the Init function through the Pythia::readString("Random:seed") mechanism.
// Value of -1 means unitialized; 0 will be time-dependent and values >1 <= MAX_SEED concrete reproducible seeding
Pythia8GenConfig mGenConfig; // configuration object
static std::atomic<int> Pythia8InstanceCounter;
int mThisPythia8InstanceID = 0;
constexpr static long MAX_SEED = 900000000;
ClassDefOverride(GeneratorPythia8, 1);
}; /** class GeneratorPythia8 **/
/*****************************************************************/
/*****************************************************************/
} // namespace eventgen
} // namespace o2
#endif /* ALICEO2_EVENTGEN_GENERATORPYTHIA8_H_ */
//
// EOF
//