-
Notifications
You must be signed in to change notification settings - Fork 496
Expand file tree
/
Copy pathGeneratorFactory.cxx
More file actions
286 lines (272 loc) · 11.8 KB
/
GeneratorFactory.cxx
File metadata and controls
286 lines (272 loc) · 11.8 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
// 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 S. Wenzel - Mai 2018
#include <Generators/GeneratorFactory.h>
#include "FairPrimaryGenerator.h"
#include "FairGenerator.h"
#include "FairBoxGenerator.h"
#include <FairLogger.h>
#include <SimConfig/SimConfig.h>
#include <Generators/GeneratorFromFile.h>
#ifdef GENERATORS_WITH_PYTHIA6
#include <Generators/GeneratorPythia6.h>
#include <Generators/GeneratorPythia6Param.h>
#endif
#ifdef GENERATORS_WITH_PYTHIA8
#include <Generators/GeneratorPythia8.h>
#include <Generators/GeneratorPythia8Param.h>
#endif
#include <Generators/GeneratorTGenerator.h>
#include <Generators/GeneratorExternalParam.h>
#ifdef GENERATORS_WITH_HEPMC3
#include <Generators/GeneratorHepMC.h>
#include <Generators/GeneratorHepMCParam.h>
#endif
#include <Generators/BoxGunParam.h>
#include <Generators/PDG.h>
#include <Generators/TriggerParticle.h>
#include <Generators/TriggerExternalParam.h>
#include <Generators/TriggerParticleParam.h>
#include "CommonUtils/ConfigurationMacroHelper.h"
#include "TRandom.h"
namespace o2
{
namespace eventgen
{
// reusable helper class
// main purpose is to init a FairPrimGen given some (Sim)Config
void GeneratorFactory::setPrimaryGenerator(o2::conf::SimConfig const& conf, FairPrimaryGenerator* primGen)
{
if (!primGen) {
LOG(warning) << "No primary generator instance; Cannot setup";
return;
}
auto makeBoxGen = [](int pdgid, int mult, double etamin, double etamax, double pmin, double pmax, double phimin, double phimax, bool debug = false) {
auto gen = new FairBoxGenerator(pdgid, mult);
gen->SetEtaRange(etamin, etamax);
gen->SetPRange(pmin, pmax);
gen->SetPhiRange(phimin, phimax);
gen->SetDebug(debug);
return gen;
};
#ifdef GENERATORS_WITH_PYTHIA8
auto makePythia8Gen = [](std::string& config) {
auto gen = new o2::eventgen::GeneratorPythia8();
if (!config.empty()) {
LOG(info) << "Reading \'Pythia8\' base configuration: " << config << std::endl;
gen->readFile(config);
}
auto seed = (gRandom->GetSeed() % 900000000);
LOG(info) << "Using random seed from gRandom % 900000000: " << seed;
gen->readString("Random:setSeed on");
gen->readString("Random:seed " + std::to_string(seed));
return gen;
};
#endif
/** generators **/
o2::PDG::addParticlesToPdgDataBase();
auto genconfig = conf.getGenerator();
if (genconfig.compare("boxgen") == 0) {
// a simple "box" generator configurable via BoxGunparam
auto& boxparam = BoxGunParam::Instance();
LOG(info) << "Init generic box generator with following parameters";
LOG(info) << boxparam;
auto boxGen = makeBoxGen(boxparam.pdg, boxparam.number, boxparam.eta[0], boxparam.eta[1], boxparam.prange[0], boxparam.prange[1], boxparam.phirange[0], boxparam.phirange[1], boxparam.debug);
primGen->AddGenerator(boxGen);
} else if (genconfig.compare("fwmugen") == 0) {
// a simple "box" generator for forward muons
LOG(info) << "Init box forward muons generator";
auto boxGen = makeBoxGen(13, 1, -4, -2.5, 50., 50., 0., 360);
primGen->AddGenerator(boxGen);
} else if (genconfig.compare("hmpidgun") == 0) {
// a simple "box" generator for forward muons
LOG(info) << "Init hmpid gun generator";
auto boxGen = makeBoxGen(-211, 100, -0.5, -0.5, 2, 5, -5, 60);
primGen->AddGenerator(boxGen);
} else if (genconfig.compare("fwpigen") == 0) {
// a simple "box" generator for forward pions
LOG(info) << "Init box forward pions generator";
auto boxGen = makeBoxGen(-211, 10, -4, -2.5, 7, 7, 0, 360);
primGen->AddGenerator(boxGen);
} else if (genconfig.compare("fwrootino") == 0) {
// a simple "box" generator for forward rootinos
LOG(info) << "Init box forward rootinos generator";
auto boxGen = makeBoxGen(0, 1, -4, -2.5, 1, 5, 0, 360);
primGen->AddGenerator(boxGen);
} else if (genconfig.compare("zdcgen") == 0) {
// a simple "box" generator for forward neutrons
LOG(info) << "Init box forward/backward zdc generator";
auto boxGenC = makeBoxGen(2112 /*neutrons*/, 1, -8, -9999, 500, 1000, 0., 360.);
auto boxGenA = makeBoxGen(2112 /*neutrons*/, 1, 8, 9999, 500, 1000, 0., 360.);
primGen->AddGenerator(boxGenC);
primGen->AddGenerator(boxGenA);
} else if (genconfig.compare("emcgenele") == 0) {
// box generator with one electron per event
LOG(info) << "Init box generator for electrons in EMCAL";
// using phi range of emcal
auto elecgen = makeBoxGen(11, 1, -0.67, 0.67, 15, 15, 80, 187);
primGen->AddGenerator(elecgen);
} else if (genconfig.compare("emcgenphoton") == 0) {
LOG(info) << "Init box generator for photons in EMCAL";
auto photongen = makeBoxGen(22, 1, -0.67, 0.67, 15, 15, 80, 187);
primGen->AddGenerator(photongen);
} else if (genconfig.compare("fddgen") == 0) {
LOG(info) << "Init box FDD generator";
auto boxGenFDC = makeBoxGen(13, 1000, -7, -4.8, 10, 500, 0, 360.);
auto boxGenFDA = makeBoxGen(13, 1000, 4.9, 6.3, 10, 500, 0., 360);
primGen->AddGenerator(boxGenFDA);
primGen->AddGenerator(boxGenFDC);
} else if (genconfig.compare("extkin") == 0) {
// external kinematics
// needs precense of a kinematics file "Kinematics.root"
// TODO: make this configurable and check for presence
auto extGen = new o2::eventgen::GeneratorFromFile(conf.getExtKinematicsFileName().c_str());
extGen->SetStartEvent(conf.getStartEvent());
primGen->AddGenerator(extGen);
LOG(info) << "using external kinematics";
} else if (genconfig.compare("extkinO2") == 0) {
// external kinematics from previous O2 output
auto extGen = new o2::eventgen::GeneratorFromO2Kine(conf.getExtKinematicsFileName().c_str());
extGen->SetStartEvent(conf.getStartEvent());
primGen->AddGenerator(extGen);
LOG(info) << "using external O2 kinematics";
#ifdef GENERATORS_WITH_HEPMC3
} else if (genconfig.compare("hepmc") == 0) {
// external HepMC file
auto& param = GeneratorHepMCParam::Instance();
LOG(info) << "Init \'GeneratorHepMC\' with following parameters";
LOG(info) << param;
auto hepmcGen = new o2::eventgen::GeneratorHepMC();
hepmcGen->setFileName(param.fileName);
hepmcGen->setVersion(param.version);
primGen->AddGenerator(hepmcGen);
#endif
#ifdef GENERATORS_WITH_PYTHIA6
} else if (genconfig.compare("pythia6") == 0) {
// pythia6 pp
// configures pythia6 according to param
auto& param = GeneratorPythia6Param::Instance();
LOG(info) << "Init \'Pythia6\' generator with following parameters";
LOG(info) << param;
auto py6Gen = new o2::eventgen::GeneratorPythia6();
py6Gen->setConfig(param.config);
py6Gen->setFrame(param.frame);
py6Gen->setBeam(param.beam);
py6Gen->setTarget(param.target);
py6Gen->setWin(param.win);
primGen->AddGenerator(py6Gen);
#endif
#ifdef GENERATORS_WITH_PYTHIA8
} else if (genconfig.compare("alldets") == 0) {
// a simple generator for test purposes - making sure to generate hits
// in all detectors
// I compose it of:
// 1) pythia8
auto py8config = std::string(std::getenv("O2_ROOT")) + "/share/Generators/egconfig/pythia8_inel.cfg";
auto py8 = makePythia8Gen(py8config);
primGen->AddGenerator(py8);
// 2) forward muons
auto muon = makeBoxGen(13, 100, -2.5, -4.0, 100, 100, 0., 360);
primGen->AddGenerator(muon);
} else if (genconfig.compare("pythia8") == 0) {
auto py8config = std::string();
auto py8 = makePythia8Gen(py8config);
primGen->AddGenerator(py8);
} else if (genconfig.compare("pythia8pp") == 0) {
auto py8config = std::string(std::getenv("O2_ROOT")) + "/share/Generators/egconfig/pythia8_inel.cfg";
auto py8 = makePythia8Gen(py8config);
primGen->AddGenerator(py8);
} else if (genconfig.compare("pythia8hf") == 0) {
// pythia8 pp (HF production)
// configures pythia for HF production in pp collisions at 14 TeV
auto py8config = std::string(std::getenv("O2_ROOT")) + "/share/Generators/egconfig/pythia8_hf.cfg";
auto py8 = makePythia8Gen(py8config);
primGen->AddGenerator(py8);
} else if (genconfig.compare("pythia8hi") == 0) {
// pythia8 heavy-ion
// exploits pythia8 heavy-ion machinery (available from v8.230)
// configures pythia for min.bias Pb-Pb collisions at 5.52 TeV
auto py8config = std::string(std::getenv("O2_ROOT")) + "/share/Generators/egconfig/pythia8_hi.cfg";
auto py8 = makePythia8Gen(py8config);
primGen->AddGenerator(py8);
#endif
} else if (genconfig.compare("external") == 0 || genconfig.compare("extgen") == 0) {
// external generator via configuration macro
auto& params = GeneratorExternalParam::Instance();
LOG(info) << "Setting up external generator with following parameters";
LOG(info) << params;
auto extgen_filename = params.fileName;
auto extgen_func = params.funcName;
auto extgen = o2::conf::GetFromMacro<FairGenerator*>(extgen_filename, extgen_func, "FairGenerator*", "extgen");
if (!extgen) {
LOG(fatal) << "Failed to retrieve \'extgen\': problem with configuration ";
}
primGen->AddGenerator(extgen);
} else if (genconfig.compare("toftest") == 0) { // 1 muon per sector and per module
LOG(info) << "Init tof test generator -> 1 muon per sector and per module";
for (int i = 0; i < 18; i++) {
for (int j = 0; j < 5; j++) {
auto boxGen = new FairBoxGenerator(13, 1); /*protons*/
boxGen->SetEtaRange(-0.8 + 0.32 * j + 0.15, -0.8 + 0.32 * j + 0.17);
boxGen->SetPRange(9, 10);
boxGen->SetPhiRange(10 + 20. * i - 1, 10 + 20. * i + 1);
boxGen->SetDebug(kTRUE);
primGen->AddGenerator(boxGen);
}
}
} else {
LOG(fatal) << "Invalid generator";
}
/** triggers **/
Trigger trigger = nullptr;
DeepTrigger deeptrigger = nullptr;
auto trgconfig = conf.getTrigger();
if (trgconfig.empty()) {
return;
} else if (trgconfig.compare("particle") == 0) {
trigger = TriggerParticle(TriggerParticleParam::Instance());
} else if (trgconfig.compare("external") == 0) {
// external trigger via configuration macro
auto& params = TriggerExternalParam::Instance();
LOG(info) << "Setting up external trigger with following parameters";
LOG(info) << params;
auto external_trigger_filename = params.fileName;
auto external_trigger_func = params.funcName;
trigger = o2::conf::GetFromMacro<o2::eventgen::Trigger>(external_trigger_filename, external_trigger_func, "o2::eventgen::Trigger", "trigger");
if (!trigger) {
LOG(info) << "Trying to retrieve a \'o2::eventgen::DeepTrigger\' type" << std::endl;
deeptrigger = o2::conf::GetFromMacro<o2::eventgen::DeepTrigger>(external_trigger_filename, external_trigger_func, "o2::eventgen::DeepTrigger", "deeptrigger");
}
if (!trigger && !deeptrigger) {
LOG(fatal) << "Failed to retrieve \'external trigger\': problem with configuration ";
}
} else {
LOG(fatal) << "Invalid trigger";
}
/** add trigger to generators **/
auto generators = primGen->GetListOfGenerators();
for (int igen = 0; igen < generators->GetEntries(); ++igen) {
auto generator = dynamic_cast<o2::eventgen::Generator*>(generators->At(igen));
if (!generator) {
LOG(fatal) << "request to add a trigger to an unsupported generator";
return;
}
generator->setTriggerMode(o2::eventgen::Generator::kTriggerOR);
if (trigger) {
generator->addTrigger(trigger);
}
if (deeptrigger) {
generator->addDeepTrigger(deeptrigger);
}
}
}
} // end namespace eventgen
} // end namespace o2