Skip to content

Commit 305f55e

Browse files
glegrasglegras
committed
possibility to pass rejection factors by string
1 parent 10623d9 commit 305f55e

3 files changed

Lines changed: 156 additions & 42 deletions

File tree

‎MC/config/PWGDQ/EvtGen/DecayTablesEvtgen/BTOPSIJPSITODIELECTRON.DEC‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -371,7 +371,7 @@ Enddecay
371371

372372
Decay psi(2S)
373373
0.00795 e+ e- PHOTOS VLL;
374-
0.3478 J/psi pi+ pi- VVPIPI;
374+
0.0208 J/psi pi+ pi- VVPIPI; # multiplied by J/psi->ee BR to keep relative statistics, since J/psi->ee itself is also enhanced
375375
Enddecay
376376

377377
Decay J/psi

‎MC/config/PWGDQ/external/generator/generator_pythia8_HadronTriggered_withGap.C‎

Lines changed: 148 additions & 41 deletions
Original file line numberDiff line numberDiff line change
@@ -19,14 +19,14 @@ class GeneratorPythia8HadronTriggeredWithGap : public o2::eventgen::GeneratorPyt
1919
public:
2020

2121
/// constructor
22-
GeneratorPythia8HadronTriggeredWithGap(int inputTriggerRatio = 5, bool useOniaShower = false) {
22+
GeneratorPythia8HadronTriggeredWithGap(int inputTriggerRatio = 5) {
2323

2424
mGeneratedEvents = 0;
2525
mInverseTriggerRatio = inputTriggerRatio;
2626
// define minimum bias event generator
2727
auto seed = (gRandom->TRandom::GetSeed() % 900000000);
2828
// main physics option for the min bias pythia events: SoftQCD:Inelastic
29-
TString pathconfigMB = useOniaShower ? gSystem->ExpandPathName("${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_oniaAll_triggerGap.cfg") : gSystem->ExpandPathName("${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_triggerGap.cfg");
29+
TString pathconfigMB = gSystem->ExpandPathName("${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_triggerGap.cfg");
3030
pythiaMBgen.readFile(pathconfigMB.Data());
3131
pythiaMBgen.readString("Random:setSeed on");
3232
pythiaMBgen.readString("Random:seed " + std::to_string(seed));
@@ -39,7 +39,23 @@ public:
3939
/// Destructor
4040
~GeneratorPythia8HadronTriggeredWithGap() = default;
4141

42-
void addHadronPDGs(int pdg) { mHadronsPDGs.push_back(pdg); };
42+
void addHadronPDGs(int pdg) { mHadronsPDGs.push_back(pdg); mRejFactorPrompt.push_back(1.0); mRejFactorNonPrompt.push_back(1.0);}
43+
44+
void setRejFactorPrompt(int pdg, float rejFactor) {
45+
for (size_t i = 0; i < mHadronsPDGs.size(); i++) {
46+
if (pdg == mHadronsPDGs[i]) {
47+
mRejFactorPrompt[i] = rejFactor;
48+
}
49+
}
50+
}
51+
52+
void setRejFactorNonPrompt(int pdg, float rejFactor) {
53+
for (size_t i = 0; i < mHadronsPDGs.size(); i++) {
54+
if (pdg == mHadronsPDGs[i]) {
55+
mRejFactorNonPrompt[i] = rejFactor;
56+
}
57+
}
58+
}
4359

4460
void setRapidityRange(double valMin, double valMax)
4561
{
@@ -93,28 +109,72 @@ bool Init() override {
93109
addSubGenerator(1, "Hadron triggered");
94110
GeneratorPythia8::Init();
95111
pythiaMBgen.init();
112+
113+
for (size_t i = 0; i < mHadronsPDGs.size(); i++) {
114+
LOGF(info, "triggering hadron %d with rejection factor (prompt/non-prompt) %f/%f", mHadronsPDGs[i], mRejFactorPrompt[i], mRejFactorNonPrompt[i]);
115+
}
116+
96117
return true;
97118
}
98119

120+
121+
bool isOpenBhadron(int pdg) {
122+
// all open beauty hadrons, no upsilon
123+
return ((abs(pdg) >= 500 && abs(pdg) < 599) || (abs(pdg) >= 5000 && abs(pdg) < 5999)) && pdg != 553;
124+
}
125+
99126
// search for the presence of at least one of the required hadrons in a selected rapidity window
100127
bool findHadrons(Pythia8::Event& event) {
101-
128+
int ihad = 0;
102129
for (int ipa = 0; ipa < event.size(); ++ipa) {
103130

104131
auto daughterList = event[ipa].daughterList();
105132

106133
for (auto ida : daughterList) {
134+
ihad = 0;
107135
for (int pdg : mHadronsPDGs) { // check that at least one of the pdg code is found in the event
108136
if (event[ida].id() == pdg) {
109137
if ((event[ida].y() > mRapidityMin) && (event[ida].y() < mRapidityMax)) {
110-
cout << "============= Found jpsi y,pt " << event[ida].y() << ", " << event[ida].pT() << endl;
138+
cout << "============= Found jpsi y,pt,pdg " << event[ida].y() << ", " << event[ida].pT() << ", " << event[ida].pdg() << endl;
111139
std::vector<int> daughters = event[ida].daughterList();
112140
for (int d : daughters) {
113141
cout << "###### daughter " << d << ": code " << event[d].id() << ", pt " << event[d].pT() << endl;
114142
}
115-
return true;
143+
144+
// check whether particle is prompt or non-prompt, since rejection factor can depend on it
145+
bool isNonPrompt = false;
146+
if (isOpenBhadron(pdg)) {
147+
isNonPrompt = true;
148+
cout << "particle is non-prompt" << endl
149+
} else {
150+
// we check the mother
151+
int mother = event[ida].mother1();
152+
if (mother >= 0 && isOpenBhadron(event[mother].id())) {
153+
isNonPrompt = true;
154+
cout << "particle is non-prompt, mother pdg: " << event[mother].id() << endl;
155+
}
156+
if (mother >= 0 && !isOpenBhadron(event[mother].id())) {
157+
// we check the grand-mother
158+
int grandmother = event[mother].mother1();
159+
if (grandmother >= 0 && isOpenBhadron(event[grandmother].id())) {
160+
isNonPrompt = true;
161+
cout << "particle is non-prompt, mother pdg: " << event[mother].id() << ", grand-mother pdg: "<< event[grandmother].id() << endl;
162+
}
163+
if (grandmother >= 0 && !isOpenBhadron(event[grandmother].id())) {
164+
isNonPrompt = false;
165+
cout << "particle is prompt, mother pdg: " << event[mother].id() << ", grand-mother pdg: "<< event[grandmother].id() << endl;
166+
}
167+
}
168+
}
169+
170+
// rejection factor given in the ini file
171+
float randomNumber = gRandom->Rndm();
172+
if ((isPrompt && (randomNumber <= mRejFactorPrompt[ihad])) || (isNonPrompt && (randomNumber <= mRejFactorNonPrompt[ihad]))) {
173+
return true;
174+
}
116175
}
117176
}
177+
ihad++;
118178
}
119179
}
120180
}
@@ -133,6 +193,8 @@ private:
133193
Pythia8::Pythia pythiaMBgen; // minimum bias event
134194
TString mConfigMBdecays;
135195
std::vector<int> mHadronsPDGs;
196+
std::vector<float> mRejFactorPrompt; // rejection factors for possibility to trigger a given particle only a fraction of the time, 1 by default
197+
std::vector<float> mRejFactorNonPrompt;
136198
double mRapidityMin;
137199
double mRapidityMax;
138200
bool mVerbose;
@@ -251,52 +313,97 @@ GeneratorInclusiveJpsiPsi2SChiC_EvtGenMidY(int triggerGap, double rapidityMin =
251313
return gen;
252314
}
253315
FairGenerator *
254-
GeneratorInclusiveAllQuarkonia_EvtGenMidY(int triggerGap, double rapidityMin = -1.0, double rapidityMax = 1.0, bool verbose = false)
316+
GeneratorInclusiveAllQuarkonia_EvtGenMidY(int triggerGap, double rapidityMin = -1.0, double rapidityMax = 1.0, TString rejFactors = "", bool verbose = false)
255317
{
256-
auto gen = new o2::eventgen::GeneratorEvtGen<o2::eventgen::GeneratorPythia8HadronTriggeredWithGap>(triggerGap, true);
318+
319+
int particleList[16] = {443, // Jpsi
320+
100443, // psi(2S)
321+
10441, // chic0
322+
20443, // chic1
323+
445, // chic2
324+
553, // upsilon(1S)
325+
100553, // upsilon(2S)
326+
200553, // upsilon(3S)
327+
// we also add B hadrons to trigger correct rapidity range (e.g. B is within |y|<1 but non-prompt J/psi has |y|>1)
328+
511, // B0
329+
521, // B+
330+
531, // Bs
331+
541, // Bc
332+
5122, // Lambdab
333+
5132, // Xib+
334+
5232, // Xib0
335+
5332 // Omegab
336+
};
337+
338+
auto gen = new o2::eventgen::GeneratorEvtGen<o2::eventgen::GeneratorPythia8HadronTriggeredWithGap>();
257339
gen->setTriggerGap(triggerGap);
258340
gen->setRapidityRange(rapidityMin, rapidityMax);
259-
gen->addHadronPDGs(443); // Jpsi
260-
gen->addHadronPDGs(100443); // psi(2S)
261-
gen->addHadronPDGs(10441); // chic0
262-
gen->addHadronPDGs(20443); // chic1
263-
gen->addHadronPDGs(445); // chic2
264-
gen->addHadronPDGs(553); // upsilon(1S)
265-
gen->addHadronPDGs(100553); // upsilon(2S)
266-
gen->addHadronPDGs(200553); // upsilon(3S)
267-
// we also add B hadrons to trigger correct rapidity range (e.g. B is within |y|<1 but non-prompt J/psi has |y|>1)
268-
gen->addHadronPDGs(511); // B0
269-
gen->addHadronPDGs(521); // B+
270-
gen->addHadronPDGs(531); // Bs
271-
gen->addHadronPDGs(541); // Bc
272-
gen->addHadronPDGs(5122); // Lambdab
273-
gen->addHadronPDGs(5132); // Xib+
274-
gen->addHadronPDGs(5232); // Xib0
275-
gen->addHadronPDGs(5332); // Omegab
341+
// specify particles to be triggered
342+
for (int i = 0; i < 16; i++) {
343+
gen->addHadronPDGs(particleList[i]);
344+
}
276345
gen->setVerbose(verbose);
346+
347+
// possibility to enhance a particle compared to another (or completely reject one particle) using rejection factors configured from a string
348+
// the rejection factors can be kept in a comma separated list (e.g. "pdg1:rejFactor1,pdg2:rejFactor2")
349+
// also the keywords "prompt" and "non-prompt" can be used to modify all prompt and all non-prompt (e.g. "prompt:rejFactor1,non-prompt:rejFactor2")
350+
// or the keyword can be used for only one particle (e.g. "pdg1:rejFactor1:prompt,pdg1:rejFactor2:non-prompt")
351+
TObjArray* objArray = rejFactors.Tokenize(",");
352+
for (int i = 0; i < objArray->GetEntries(); i++) {
353+
TString rejStr = objArray->At(i);
354+
TObjArray* objArrayCurrent = rejStr.Tokenize(":");
355+
if (objArrayCurrent->GetEntries() != 2 && objArrayCurrent->GetEntries() != 3) {
356+
LOGF(fatal, "Problem when configuring string for particle rejection factors: %s, incorrect length", rejStr.Data());
357+
}
358+
if (!objArrayCurrent[1].IsFloat()) {
359+
LOGF(fatal, "Problem when configuring string for particle rejection factors: %s, is not float", rejStr.Data());
360+
}
361+
if (objArrayCurrent[0].CompareTo("prompt") == 0) {
362+
// Common switch for all prompt particles
363+
for (int ihad = 0; ihad < 16; i++) {
364+
gen->setRejFactorPrompt(particleList[i], objArrayCurrent[1].Atof());
365+
}
366+
continue;
367+
}
368+
if (objArrayCurrent[0].CompareTo("non-prompt") == 0) {
369+
// Common switch for all non-prompt particles
370+
for (int ihad = 0; ihad < 16; i++) {
371+
gen->setRejFactorNonPrompt(particleList[i], objArrayCurrent[1].Atof());
372+
}
373+
continue;
374+
}
375+
if (objArrayCurrent[0].IsDigit()) {
376+
// Setting the rejection factor for a specific particle
377+
if (objArrayCurrent->GetEntries() == 2) {
378+
gen->setRejFactorPrompt(objArrayCurrent[0].Atoi(), objArrayCurrent[1].Atof());
379+
gen->setRejFactorNonPrompt(objArrayCurrent[0].Atoi(), objArrayCurrent[1].Atof());
380+
continue;
381+
}
382+
else {
383+
if (objArrayCurrent[2].CompareTo("prompt") == 0) {
384+
gen->setRejFactorPrompt(objArrayCurrent[0].Atoi(), objArrayCurrent[1].Atof());
385+
continue;
386+
}
387+
if (objArrayCurrent[2].CompareTo("non-prompt") == 0) {
388+
gen->setRejFactorNonPrompt(objArrayCurrent[0].Atoi(), objArrayCurrent[1].Atof());
389+
continue;
390+
}
391+
}
392+
}
393+
LOGF(fatal, "Problem when configuring string for particle rejection factors: %s, incorrect template", rejStr.Data());
394+
}
395+
277396

278397
TString pathO2table = gSystem->ExpandPathName("${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/decayer/switchOffAllQuarkonia.cfg");
279398
gen->readFile(pathO2table.Data());
280399
gen->setConfigMBdecays(pathO2table);
281400
gen->PrintDebug(true);
282401

402+
// specify particles to be decayed with EvtGen
283403
gen->SetSizePdg(16);
284-
gen->AddPdg(443, 0);
285-
gen->AddPdg(100443, 1);
286-
gen->AddPdg(10441, 2);
287-
gen->AddPdg(20443, 3);
288-
gen->AddPdg(445, 4);
289-
gen->AddPdg(553, 5);
290-
gen->AddPdg(100553, 6);
291-
gen->AddPdg(200553, 7);
292-
gen->AddPdg(511, 8);
293-
gen->AddPdg(521, 9);
294-
gen->AddPdg(531, 10);
295-
gen->AddPdg(541, 11);
296-
gen->AddPdg(5122, 12);
297-
gen->AddPdg(5132, 13);
298-
gen->AddPdg(5232, 14);
299-
gen->AddPdg(5332, 15);
404+
for (int i = 0; i < 16; i++) {
405+
gen->AddPdg(particleList[i], i);
406+
}
300407

301408
gen->SetForceDecay(kEvtBPsiAndJpsiDiElectron);
302409

‎MC/config/PWGDQ/ini/Generator_TriggeredInclusiveQuarkonia_Pythia8_TriggerGap.ini‎

Lines changed: 7 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,7 @@
1+
### The external generator derives from GeneratorPythia8.
2+
[GeneratorExternal]
3+
fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/external/generator/generator_pythia8_HadronTriggered_withGap.C
4+
funcName=GeneratorInclusiveAllQuarkonia_EvtGenMidY(1,-1.0,1.0)
5+
6+
[GeneratorPythia8]
7+
config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_oniaAll_triggerGap.cfg

0 commit comments

Comments
 (0)