Skip to content

Commit 58146a1

Browse files
lucamicheletti93lucamicheletti
lucamicheletti93
and
lucamicheletti
authored
Fix light ion generator (#2118)
Co-authored-by: Lucamicheletti93 <luca.mike93@gmail.com>
1 parent 0b2c114 commit 58146a1

6 files changed

+187
-6
lines changed

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

Lines changed: 8 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -133,11 +133,17 @@ class O2_GeneratorParamPsi : public GeneratorTGenerator
133133
static Double_t YPsipp5TeV(const Double_t* py, const Double_t* /*dummy*/)
134134
{
135135
// psi2s y in pp at 5.02 TeV, tuned on https://www.hepdata.net/record/ins1935680
136+
// WARNING! The shape extracted from data provide wired rapidity shape (low stat.), the J/psi one is used
136137
Double_t y = *py;
137138
Float_t p0, p1, p2;
139+
// Extracted from Psi(2S) Run 2 data
140+
//p0 = 1;
141+
//p1 = -17.4857;
142+
//p2 = 2.98887;
143+
// Same parametrization as J/psi
138144
p0 = 1;
139-
p1 = -17.4857;
140-
p2 = 2.98887;
145+
p1 = 0.0338222;
146+
p2 = 2.96748;
141147
return p0 * TMath::Exp(-(1. / 2.) * TMath::Power(((y - p1) / p2), 2));
142148
}
143149

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

Lines changed: 8 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -756,11 +756,17 @@ class O2_GeneratorParamPsipp5TeV : public GeneratorTGenerator
756756
static Double_t YPsipp5TeV(const Double_t* py, const Double_t* /*dummy*/)
757757
{
758758
// psi2s y in pp at 5.02 TeV, tuned on https://www.hepdata.net/record/ins1935680
759+
// WARNING! The shape extracted from data provide wired rapidity shape (low stat.), the J/psi one is used
759760
Double_t y = *py;
760761
Float_t p0, p1, p2;
762+
// Extracted from Psi(2S) Run 2 data
763+
//p0 = 1;
764+
//p1 = -17.4857;
765+
//p2 = 2.98887;
766+
// Same parametrization as J/psi
761767
p0 = 1;
762-
p1 = -17.4857;
763-
p2 = 2.98887;
768+
p1 = 0.0338222;
769+
p2 = 2.96748;
764770
return p0 * TMath::Exp(-(1. / 2.) * TMath::Power(((y - p1) / p2), 2));
765771
}
766772

‎MC/config/PWGDQ/ini/Generator_InjectedPromptCharmoniaFwdy_TriggerGap_OO5TeV.ini‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -4,4 +4,4 @@ fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/external/generator/generator_py
44
funcName=GeneratorPythia8InjectedPromptCharmoniaGapTriggered(5,11)
55

66
[GeneratorPythia8]
7-
config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_triggerGap_OO5TeV.cfg
7+
config=${O2DPG_MC_CONFIG_ROOT}/MC/config/common/pythia8/generator/pythia8_OO_536.cfg

‎MC/config/PWGDQ/ini/Generator_InjectedPromptCharmoniaMidy_TriggerGap_OO5TeV.ini‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -4,4 +4,4 @@ fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/external/generator/generator_py
44
funcName=GeneratorPythia8InjectedPromptCharmoniaGapTriggered(5,12)
55

66
[GeneratorPythia8]
7-
config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_triggerGap_OO5TeV.cfg
7+
config=${O2DPG_MC_CONFIG_ROOT}/MC/config/common/pythia8/generator/pythia8_OO_536.cfg

‎MC/config/PWGDQ/ini/tests/Generator_InjectedPromptCharmoniaFwdy_TriggerGap_pp5TeV.C‎

Lines changed: 85 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,85 @@
1+
int External()
2+
{
3+
int checkPdgSignal[] = {443,100443};
4+
int checkPdgDecay = 13;
5+
double rapiditymin = -4.3; double rapiditymax = -2.3;
6+
std::string path{"o2sim_Kine.root"};
7+
std::cout << "Check for\nsignal PDG " << checkPdgSignal << "\ndecay PDG " << checkPdgDecay << "\n";
8+
TFile file(path.c_str(), "READ");
9+
if (file.IsZombie()) {
10+
std::cerr << "Cannot open ROOT file " << path << "\n";
11+
return 1;
12+
}
13+
14+
auto tree = (TTree*)file.Get("o2sim");
15+
std::vector<o2::MCTrack>* tracks{};
16+
tree->SetBranchAddress("MCTrack", &tracks);
17+
18+
int nLeptons{};
19+
int nAntileptons{};
20+
int nLeptonPairs{};
21+
int nLeptonPairsToBeDone{};
22+
int nSignalJpsi{};
23+
int nSignalPsi2S{};
24+
int nSignalJpsiWithinAcc{};
25+
int nSignalPsi2SWithinAcc{};
26+
auto nEvents = tree->GetEntries();
27+
o2::steer::MCKinematicsReader mcreader("o2sim", o2::steer::MCKinematicsReader::Mode::kMCKine);
28+
Bool_t isInjected = kFALSE;
29+
30+
for (int i = 0; i < nEvents; i++) {
31+
tree->GetEntry(i);
32+
for (auto& track : *tracks) {
33+
auto pdg = track.GetPdgCode();
34+
auto rapidity = track.GetRapidity();
35+
auto idMoth = track.getMotherTrackId();
36+
if (pdg == checkPdgDecay) {
37+
// count leptons
38+
nLeptons++;
39+
} else if(pdg == -checkPdgDecay) {
40+
// count anti-leptons
41+
nAntileptons++;
42+
} else if (pdg == checkPdgSignal[0] || pdg == checkPdgSignal[1]) {
43+
if(idMoth < 0) {
44+
// count signal PDG
45+
pdg == checkPdgSignal[0] ? nSignalJpsi++ : nSignalPsi2S++;
46+
// count signal PDG within acceptance
47+
if(rapidity > rapiditymin && rapidity < rapiditymax) { pdg == checkPdgSignal[0] ? nSignalJpsiWithinAcc++ : nSignalPsi2SWithinAcc++;}
48+
}
49+
auto child0 = o2::mcutils::MCTrackNavigator::getDaughter0(track, *tracks);
50+
auto child1 = o2::mcutils::MCTrackNavigator::getDaughter1(track, *tracks);
51+
if (child0 != nullptr && child1 != nullptr) {
52+
// check for parent-child relations
53+
auto pdg0 = child0->GetPdgCode();
54+
auto pdg1 = child1->GetPdgCode();
55+
std::cout << "First and last children of parent " << checkPdgSignal << " are PDG0: " << pdg0 << " PDG1: " << pdg1 << "\n";
56+
if (std::abs(pdg0) == checkPdgDecay && std::abs(pdg1) == checkPdgDecay && pdg0 == -pdg1) {
57+
nLeptonPairs++;
58+
if (child0->getToBeDone() && child1->getToBeDone()) {
59+
nLeptonPairsToBeDone++;
60+
}
61+
}
62+
}
63+
}
64+
}
65+
}
66+
std::cout << "#events: " << nEvents << "\n"
67+
<< "#leptons: " << nLeptons << "\n"
68+
<< "#antileptons: " << nAntileptons << "\n"
69+
<< "#signal (prompt Jpsi): " << nSignalJpsi << "; within acceptance " << rapiditymin << " < y < " << rapiditymax << " : " << nSignalJpsiWithinAcc << "\n"
70+
<< "#signal (prompt Psi(2S)): " << nSignalPsi2S << "; within acceptance " << rapiditymin << " < y < " << rapiditymax << " : " << nSignalPsi2SWithinAcc << "\n"
71+
<< "#lepton pairs: " << nLeptonPairs << "\n"
72+
<< "#lepton pairs to be done: " << nLeptonPairs << "\n";
73+
74+
75+
if (nLeptonPairs == 0 || nLeptons == 0 || nAntileptons == 0) {
76+
std::cerr << "Number of leptons, number of anti-leptons as well as number of lepton pairs should all be greater than 1.\n";
77+
return 1;
78+
}
79+
if (nLeptonPairs != nLeptonPairsToBeDone) {
80+
std::cerr << "The number of lepton pairs should be the same as the number of lepton pairs which should be transported.\n";
81+
return 1;
82+
}
83+
84+
return 0;
85+
}

‎MC/config/PWGDQ/ini/tests/Generator_InjectedPromptCharmoniaMidy_TriggerGap_OO5TeV.C‎

Lines changed: 84 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,84 @@
1+
int External()
2+
{
3+
int checkPdgSignal[] = {443,100443};
4+
int checkPdgDecay = 11;
5+
std::string path{"o2sim_Kine.root"};
6+
std::cout << "Check for\nsignal PDG " << checkPdgSignal << "\ndecay PDG " << checkPdgDecay << "\n";
7+
TFile file(path.c_str(), "READ");
8+
if (file.IsZombie()) {
9+
std::cerr << "Cannot open ROOT file " << path << "\n";
10+
return 1;
11+
}
12+
13+
auto tree = (TTree*)file.Get("o2sim");
14+
std::vector<o2::MCTrack>* tracks{};
15+
tree->SetBranchAddress("MCTrack", &tracks);
16+
17+
int nLeptons{};
18+
int nAntileptons{};
19+
int nLeptonPairs{};
20+
int nLeptonPairsToBeDone{};
21+
int nSignalJpsi{};
22+
int nSignalPsi2S{};
23+
int nSignalJpsiWithinAcc{};
24+
int nSignalPsi2SWithinAcc{};
25+
auto nEvents = tree->GetEntries();
26+
o2::steer::MCKinematicsReader mcreader("o2sim", o2::steer::MCKinematicsReader::Mode::kMCKine);
27+
Bool_t isInjected = kFALSE;
28+
29+
for (int i = 0; i < nEvents; i++) {
30+
tree->GetEntry(i);
31+
for (auto& track : *tracks) {
32+
auto pdg = track.GetPdgCode();
33+
auto rapidity = track.GetRapidity();
34+
auto idMoth = track.getMotherTrackId();
35+
if (pdg == checkPdgDecay) {
36+
// count leptons
37+
nLeptons++;
38+
} else if(pdg == -checkPdgDecay) {
39+
// count anti-leptons
40+
nAntileptons++;
41+
} else if (pdg == checkPdgSignal[0] || pdg == checkPdgSignal[1]) {
42+
if(idMoth < 0){
43+
// count signal PDG
44+
pdg == checkPdgSignal[0] ? nSignalJpsi++ : nSignalPsi2S++;
45+
// count signal PDG within acceptance
46+
if(std::abs(rapidity) < 1.0) { pdg == checkPdgSignal[0] ? nSignalJpsiWithinAcc++ : nSignalPsi2SWithinAcc++;}
47+
}
48+
auto child0 = o2::mcutils::MCTrackNavigator::getDaughter0(track, *tracks);
49+
auto child1 = o2::mcutils::MCTrackNavigator::getDaughter1(track, *tracks);
50+
if (child0 != nullptr && child1 != nullptr) {
51+
// check for parent-child relations
52+
auto pdg0 = child0->GetPdgCode();
53+
auto pdg1 = child1->GetPdgCode();
54+
std::cout << "First and last children of parent " << checkPdgSignal << " are PDG0: " << pdg0 << " PDG1: " << pdg1 << "\n";
55+
if (std::abs(pdg0) == checkPdgDecay && std::abs(pdg1) == checkPdgDecay && pdg0 == -pdg1) {
56+
nLeptonPairs++;
57+
if (child0->getToBeDone() && child1->getToBeDone()) {
58+
nLeptonPairsToBeDone++;
59+
}
60+
}
61+
}
62+
}
63+
}
64+
}
65+
std::cout << "#events: " << nEvents << "\n"
66+
<< "#leptons: " << nLeptons << "\n"
67+
<< "#antileptons: " << nAntileptons << "\n"
68+
<< "#signal (prompt Jpsi): " << nSignalJpsi << "; within acceptance (|y| < 1): " << nSignalJpsiWithinAcc << "\n"
69+
<< "#signal (prompt Psi(2S)): " << nSignalPsi2S << "; within acceptance (|y| < 1): " << nSignalPsi2SWithinAcc << "\n"
70+
<< "#lepton pairs: " << nLeptonPairs << "\n"
71+
<< "#lepton pairs to be done: " << nLeptonPairs << "\n";
72+
73+
74+
if (nLeptonPairs == 0 || nLeptons == 0 || nAntileptons == 0) {
75+
std::cerr << "Number of leptons, number of anti-leptons as well as number of lepton pairs should all be greater than 1.\n";
76+
return 1;
77+
}
78+
if (nLeptonPairs != nLeptonPairsToBeDone) {
79+
std::cerr << "The number of lepton pairs should be the same as the number of lepton pairs which should be transported.\n";
80+
return 1;
81+
}
82+
83+
return 0;
84+
}

0 commit comments

Comments
 (0)