Skip to content

Commit 72014f6

Browse files
lucamicheletti93lucamicheletti
authored and
chiarazampolli
committed
[PWGDQ] Adding Upsilon Injected Generator (#1787)
* Adding Upsilon Generator * fixing ini fil bug * Fix the test --------- Co-authored-by: Lucamicheletti93 <luca.mike93@gmail.com> (cherry picked from commit 06628bc)
1 parent c790ffa commit 72014f6

File tree

4 files changed

+466
-0
lines changed

4 files changed

+466
-0
lines changed

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

Lines changed: 280 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,280 @@
1+
//
2+
// generators for bottomonia considering at midrapidity and forward rapidity
3+
//
4+
5+
R__ADD_INCLUDE_PATH($O2DPG_ROOT/MC/config/PWGDQ/EvtGen)
6+
R__LOAD_LIBRARY(libpythia6)
7+
#include "GeneratorCocktail.C"
8+
#include "GeneratorEvtGen.C"
9+
10+
namespace o2
11+
{
12+
namespace eventgen
13+
{
14+
15+
/////////////////////////////////////////////////////////////////////////////
16+
class O2_GeneratorParamUpsilon1SFwdY : public GeneratorTGenerator
17+
{
18+
19+
public:
20+
O2_GeneratorParamUpsilon1SFwdY() : GeneratorTGenerator("ParamUpsilon1S")
21+
{
22+
paramUpsilon1S = new GeneratorParam(1, -1, PtUpsilon1Spp13TeV, YUpsilon1Spp13TeV, V2Upsilon1Spp13TeV, IpUpsilon1Spp13TeV);
23+
paramUpsilon1S->SetMomentumRange(0., 1.e6);
24+
paramUpsilon1S->SetPtRange(0, 999.);
25+
paramUpsilon1S->SetYRange(-4.2, -2.3);
26+
paramUpsilon1S->SetPhiRange(0., 360.);
27+
paramUpsilon1S->SetDecayer(new TPythia6Decayer());
28+
paramUpsilon1S->SetForceDecay(kNoDecay); // particle left undecayed
29+
// - - paramUpsilon1S->SetTrackingFlag(1); // (from AliGenParam) -> check this
30+
setTGenerator(paramUpsilon1S);
31+
};
32+
33+
~O2_GeneratorParamUpsilon1SFwdY()
34+
{
35+
delete paramUpsilon1S;
36+
};
37+
38+
Bool_t Init() override
39+
{
40+
GeneratorTGenerator::Init();
41+
paramUpsilon1S->Init();
42+
return true;
43+
}
44+
45+
void SetNSignalPerEvent(Int_t nsig) { paramUpsilon1S->SetNumberParticles(nsig); }
46+
47+
//-------------------------------------------------------------------------//
48+
static Double_t PtUpsilon1Spp13TeV(const Double_t* px, const Double_t* /*dummy*/)
49+
{
50+
// Upsilon1S pt shape from LHCb pp@13TeV arXiv:1804.09214
51+
Double_t x = *px;
52+
Float_t p0, p1, p2, p3;
53+
54+
p0 = 4.11195e+02;
55+
p1 = 1.03097e+01;
56+
p2 = 1.62309e+00;
57+
p3 = 4.84709e+00;
58+
59+
return (p0 * x / TMath::Power(1. + TMath::Power(x / p1, p2), p3));
60+
}
61+
62+
//-------------------------------------------------------------------------//
63+
static Double_t YUpsilon1Spp13TeV(const Double_t* py, const Double_t* /*dummy*/)
64+
{
65+
// Upsilon1S y shape from LHCb pp@13TeV arXiv:1804.09214
66+
Double_t x = *py;
67+
Float_t p0, p1;
68+
69+
p0 = 3.07931e+03;
70+
p1 = -3.53102e-02;
71+
72+
return (p0 * (1. + p1 * x * x));
73+
}
74+
75+
//-------------------------------------------------------------------------//
76+
static Double_t V2Upsilon1Spp13TeV(const Double_t* /*dummy*/, const Double_t* /*dummy*/)
77+
{
78+
// Upsilon(1S) v2
79+
return 0.;
80+
}
81+
82+
//-------------------------------------------------------------------------//
83+
static Int_t IpUpsilon1Spp13TeV(TRandom*)
84+
{
85+
return 553;
86+
}
87+
88+
private:
89+
GeneratorParam* paramUpsilon1S = nullptr;
90+
};
91+
92+
/////////////////////////////////////////////////////////////////////////////
93+
class O2_GeneratorParamUpsilon2SFwdY : public GeneratorTGenerator
94+
{
95+
96+
public:
97+
O2_GeneratorParamUpsilon2SFwdY() : GeneratorTGenerator("ParamUpsilon2S")
98+
{
99+
paramUpsilon2S = new GeneratorParam(1, -1, PtUpsilon2Spp13TeV, YUpsilon2Spp13TeV, V2Upsilon2Spp13TeV, IpUpsilon2Spp13TeV);
100+
paramUpsilon2S->SetMomentumRange(0., 1.e6);
101+
paramUpsilon2S->SetPtRange(0, 999.);
102+
paramUpsilon2S->SetYRange(-4.2, -2.3);
103+
paramUpsilon2S->SetPhiRange(0., 360.);
104+
paramUpsilon2S->SetDecayer(new TPythia6Decayer());
105+
paramUpsilon2S->SetForceDecay(kNoDecay); // particle left undecayed
106+
// - - paramUpsilon2S->SetTrackingFlag(1); // (from AliGenParam) -> check this
107+
setTGenerator(paramUpsilon2S);
108+
};
109+
110+
~O2_GeneratorParamUpsilon2SFwdY()
111+
{
112+
delete paramUpsilon2S;
113+
};
114+
115+
Bool_t Init() override
116+
{
117+
GeneratorTGenerator::Init();
118+
paramUpsilon2S->Init();
119+
return true;
120+
}
121+
122+
void SetNSignalPerEvent(Int_t nsig) { paramUpsilon2S->SetNumberParticles(nsig); }
123+
124+
//-------------------------------------------------------------------------//
125+
static Double_t PtUpsilon2Spp13TeV(const Double_t* px, const Double_t* /*dummy*/)
126+
{
127+
// Upsilon2S pt shape from LHCb pp@13TeV arXiv:1804.09214
128+
Double_t x = *px;
129+
Float_t p0, p1, p2, p3;
130+
131+
p0 = 8.15699e+01;
132+
p1 = 1.48060e+01;
133+
p2 = 1.50018e+00;
134+
p3 = 6.34208e+00;
135+
136+
return (p0 * x / TMath::Power(1. + TMath::Power(x / p1, p2), p3));
137+
}
138+
139+
//-------------------------------------------------------------------------//
140+
static Double_t YUpsilon2Spp13TeV(const Double_t* py, const Double_t* /*dummy*/)
141+
{
142+
// Upsilon2s y shape from LHCb pp@13TeV arXiv:1804.09214
143+
Double_t x = *py;
144+
Float_t p0, p1;
145+
146+
p0 = 7.50409e+02;
147+
p1 = -3.57039e-02;
148+
149+
return (p0 * (1. + p1 * x * x));
150+
}
151+
152+
//-------------------------------------------------------------------------//
153+
static Double_t V2Upsilon2Spp13TeV(const Double_t* /*dummy*/, const Double_t* /*dummy*/)
154+
{
155+
// Upsilon(2S) v2
156+
return 0.;
157+
}
158+
159+
//-------------------------------------------------------------------------//
160+
static Int_t IpUpsilon2Spp13TeV(TRandom*)
161+
{
162+
return 100553;
163+
}
164+
165+
private:
166+
GeneratorParam* paramUpsilon2S = nullptr;
167+
};
168+
169+
/////////////////////////////////////////////////////////////////////////////
170+
class O2_GeneratorParamUpsilon3SFwdY : public GeneratorTGenerator
171+
{
172+
173+
public:
174+
O2_GeneratorParamUpsilon3SFwdY() : GeneratorTGenerator("ParamUpsilon3S")
175+
{
176+
paramUpsilon3S = new GeneratorParam(1, -1, PtUpsilon3Spp13TeV, YUpsilon3Spp13TeV, V2Upsilon3Spp13TeV, IpUpsilon3Spp13TeV);
177+
paramUpsilon3S->SetMomentumRange(0., 1.e6);
178+
paramUpsilon3S->SetPtRange(0, 999.);
179+
paramUpsilon3S->SetYRange(-4.2, -2.3);
180+
paramUpsilon3S->SetPhiRange(0., 360.);
181+
paramUpsilon3S->SetDecayer(new TPythia6Decayer());
182+
paramUpsilon3S->SetForceDecay(kNoDecay); // particle left undecayed
183+
// - - paramUpsilon3S->SetTrackingFlag(1); // (from AliGenParam) -> check this
184+
setTGenerator(paramUpsilon3S);
185+
};
186+
187+
~O2_GeneratorParamUpsilon3SFwdY()
188+
{
189+
delete paramUpsilon3S;
190+
};
191+
192+
Bool_t Init() override
193+
{
194+
GeneratorTGenerator::Init();
195+
paramUpsilon3S->Init();
196+
return true;
197+
}
198+
199+
void SetNSignalPerEvent(Int_t nsig) { paramUpsilon3S->SetNumberParticles(nsig); }
200+
201+
//-------------------------------------------------------------------------//
202+
static Double_t PtUpsilon3Spp13TeV(const Double_t* px, const Double_t* /*dummy*/)
203+
{
204+
// Upsilon3S pt shape from LHCb pp@13TeV arXiv:1804.09214
205+
Double_t x = *px;
206+
Float_t p0, p1, p2, p3;
207+
208+
p0 = 3.51590e+01;
209+
p1 = 2.30813e+01;
210+
p2 = 1.40822e+00;
211+
p3 = 9.38026e+00;
212+
213+
return (p0 * x / TMath::Power(1. + TMath::Power(x / p1, p2), p3));
214+
}
215+
216+
//-------------------------------------------------------------------------//
217+
static Double_t YUpsilon3Spp13TeV(const Double_t* py, const Double_t* /*dummy*/)
218+
{
219+
// Upsilon3s y shape from LHCb pp@13TeV arXiv:1804.09214
220+
Double_t x = *py;
221+
Float_t p0, p1;
222+
223+
p0 = 3.69961e+02;
224+
p1 = -3.54650e-02;
225+
226+
return (p0 * (1. + p1 * x * x));
227+
}
228+
229+
//-------------------------------------------------------------------------//
230+
static Double_t V2Upsilon3Spp13TeV(const Double_t* /*dummy*/, const Double_t* /*dummy*/)
231+
{
232+
// Upsilon(3S) v2
233+
return 0.;
234+
}
235+
236+
//-------------------------------------------------------------------------//
237+
static Int_t IpUpsilon3Spp13TeV(TRandom*)
238+
{
239+
return 200553;
240+
}
241+
242+
private:
243+
GeneratorParam* paramUpsilon3S = nullptr;
244+
};
245+
246+
247+
} // namespace eventgen
248+
} // namespace o2
249+
250+
FairGenerator* GeneratorCocktailBottomoniaToMuonEvtGen_pp13TeV()
251+
{
252+
253+
auto genCocktailEvtGen = new o2::eventgen::GeneratorEvtGen<GeneratorCocktail>();
254+
255+
auto genUpsilon1S = new o2::eventgen::O2_GeneratorParamUpsilon1SFwdY;
256+
genUpsilon1S->SetNSignalPerEvent(1); // 1 Upsilon(1S) generated per event by GeneratorParam
257+
258+
auto genUpsilon2S = new o2::eventgen::O2_GeneratorParamUpsilon2SFwdY;
259+
genUpsilon2S->SetNSignalPerEvent(1); // 1 Upsilon(2S) generated per event by GeneratorParam
260+
261+
auto genUpsilon3S = new o2::eventgen::O2_GeneratorParamUpsilon3SFwdY;
262+
genUpsilon3S->SetNSignalPerEvent(1); // 1 Upsilon(3S) generated per event by GeneratorParam
263+
264+
genCocktailEvtGen->AddGenerator(genUpsilon1S, 1); // add Upsilon(1S) generator
265+
genCocktailEvtGen->AddGenerator(genUpsilon2S, 1); // add Upsilon(2S) generator
266+
genCocktailEvtGen->AddGenerator(genUpsilon3S, 1); // add Upsilon(3S) generator
267+
268+
TString pdgs = "553;100553;200553";
269+
std::string spdg;
270+
TObjArray* obj = pdgs.Tokenize(";");
271+
genCocktailEvtGen->SetSizePdg(obj->GetEntriesFast());
272+
for (int i = 0; i < obj->GetEntriesFast(); i++) {
273+
spdg = obj->At(i)->GetName();
274+
genCocktailEvtGen->AddPdg(std::stoi(spdg), i);
275+
printf("PDG %d \n", std::stoi(spdg));
276+
}
277+
genCocktailEvtGen->SetForceDecay(kEvtDiMuon);
278+
279+
return genCocktailEvtGen;
280+
}

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

Lines changed: 76 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,76 @@
1+
#include "FairGenerator.h"
2+
#include "Generators/GeneratorPythia8.h"
3+
#include "Pythia8/Pythia.h"
4+
#include "TRandom.h"
5+
#include "GeneratorBottomonia.C"
6+
#include <string>
7+
8+
using namespace o2::eventgen;
9+
using namespace Pythia8;
10+
11+
class GeneratorPythia8BottomoniaInjectedGapTriggeredDQ : public o2::eventgen::GeneratorPythia8 {
12+
public:
13+
14+
/// default constructor
15+
GeneratorPythia8BottomoniaInjectedGapTriggeredDQ() = default;
16+
17+
/// constructor
18+
GeneratorPythia8BottomoniaInjectedGapTriggeredDQ(int inputTriggerRatio = 5, int gentype = 0) {
19+
20+
mGeneratedEvents = 0;
21+
mGeneratorParam = 0x0;
22+
mInverseTriggerRatio = inputTriggerRatio;
23+
switch (gentype) {
24+
case 0: // generate bottomonia cocktail at forward rapidity
25+
mGeneratorParam = (Generator*)GeneratorCocktailBottomoniaToMuonEvtGen_pp13TeV();
26+
break;
27+
}
28+
mGeneratorParam->Init();
29+
}
30+
31+
/// Destructor
32+
~GeneratorPythia8BottomoniaInjectedGapTriggeredDQ() = default;
33+
34+
protected:
35+
36+
Bool_t importParticles() override
37+
{
38+
GeneratorPythia8::importParticles();
39+
bool genOk = false;
40+
if (mGeneratedEvents % mInverseTriggerRatio == 0) { // add injected prompt signals to the stack
41+
bool genOk = false;
42+
while (!genOk) {
43+
genOk = (mGeneratorParam->generateEvent() && mGeneratorParam->importParticles()) ? true : false ;
44+
}
45+
int originalSize = mParticles.size();
46+
for (int ipart=0; ipart < mGeneratorParam->getParticles().size(); ipart++) {
47+
TParticle part = TParticle(mGeneratorParam->getParticles().at(ipart));
48+
if(part.GetFirstMother() >= 0) part.SetFirstMother(part.GetFirstMother() + originalSize);
49+
if(part.GetFirstDaughter() >= 0) part.SetFirstDaughter(part.GetFirstDaughter() + originalSize);
50+
if(part.GetLastDaughter() >= 0) part.SetLastDaughter(part.GetLastDaughter() + originalSize);
51+
mParticles.push_back(part);
52+
// encodeParticleStatusAndTracking method already called in GeneratorEvtGen.C
53+
}
54+
mGeneratorParam->clearParticles();
55+
}
56+
57+
mGeneratedEvents++;
58+
return true;
59+
}
60+
61+
62+
private:
63+
Generator* mGeneratorParam;
64+
// Control gap-triggering
65+
unsigned long long mGeneratedEvents;
66+
int mInverseTriggerRatio;
67+
};
68+
69+
// Predefined generators:
70+
FairGenerator *GeneratorPythia8InjectedBottomoniaGapTriggered(int inputTriggerRatio, int gentype) {
71+
auto myGen = new GeneratorPythia8BottomoniaInjectedGapTriggeredDQ(inputTriggerRatio,gentype);
72+
auto seed = (gRandom->TRandom::GetSeed() % 900000000);
73+
myGen->readString("Random:setSeed on");
74+
myGen->readString("Random:seed " + std::to_string(seed));
75+
return myGen;
76+
}

‎MC/config/PWGDQ/ini/Generator_InjectedBottomoniaFwdy_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_withInjectedBottomoniaSignals_gaptriggered_dq.C
4+
funcName=GeneratorPythia8InjectedBottomoniaGapTriggered(5,0)
5+
6+
[GeneratorPythia8]
7+
config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_triggerGap.cfg

0 commit comments

Comments
 (0)