Skip to content

Commit d830cdf

Browse files
ffiondagconesab
ffionda
authored and
gconesab
committed
add configuration files for fwd-y productions + add other fixes
1 parent 9172589 commit d830cdf

16 files changed

Lines changed: 595 additions & 7 deletions

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

Lines changed: 49 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,49 @@
1+
##################
2+
## Resonant decays not included
3+
Decay B-
4+
#
5+
# B -> cc= s sum = 1.92%
6+
#
7+
#MODIFIED in order to have 50% of B+ decaying in two bodies
8+
0.001914000 J/psi K- SVS; #[Reconstructed PDG2011]
9+
#0.001430000 J/psi K*- SVV_HELAMP PKHminus PKphHminus PKHzero PKphHzero PKHplus PKphHplus; #[Reconstructed PDG2011]
10+
0.0001 J/psi K- pi0 PHSP;
11+
0.000052000 J/psi phi K- PHSP; #[Reconstructed PDG2011]
12+
#
13+
0.001070000 J/psi K- pi- pi+ PHSP; #[New mode added] #[Reconstructed PDG2011]
14+
0.000108000 J/psi eta K- PHSP; #[New mode added] #[Reconstructed PDG2011]
15+
0.000350000 J/psi omega K- PHSP; #[New mode added] #[Reconstructed PDG2011]
16+
#
17+
#0.0001 J/psi K'_1- SVV_HELAMP 0.5 0.0 1.0 0.0 0.5 0.0;
18+
#0.0005 J/psi K_2*- PHSP;
19+
#0.001800000 J/psi K_1- SVV_HELAMP 0.5 0.0 1.0 0.0 0.5 0.0; #[Reconstructed PDG2011]
20+
21+
Enddecay
22+
23+
Decay B+
24+
#
25+
# B -> cc= s sum = 1.92%
26+
#
27+
#MODIFIED in order to have 50% of B+ decaying in two bodies
28+
0.001914000 J/psi K+ SVS; #[Reconstructed PDG2011]
29+
#0.001430000 J/psi K*+ SVV_HELAMP PKHplus PKphHplus PKHzero PKphHzero PKHminus PKphHminus; #[Reconstructed PDG2011]
30+
0.0001 J/psi K+ pi0 PHSP;
31+
0.000052000 J/psi phi K+ PHSP; #[Reconstructed PDG2011]
32+
#
33+
0.001070000 J/psi K+ pi+ pi- PHSP; #[New mode added] #[Reconstructed PDG2011]
34+
0.000108000 J/psi eta K+ PHSP; #[New mode added] #[Reconstructed PDG2011]
35+
0.000350000 J/psi omega K+ PHSP; #[New mode added] #[Reconstructed PDG2011]
36+
#
37+
#0.0001 J/psi K'_1+ SVV_HELAMP 0.5 0.0 1.0 0.0 0.5 0.0;
38+
#0.0005 J/psi K_2*+ PHSP;
39+
#0.001800000 J/psi K_1+ SVV_HELAMP 0.5 0.0 1.0 0.0 0.5 0.0; #[Reconstructed PDG2011]
40+
41+
Enddecay
42+
43+
Decay J/psi
44+
1.000 e+ e- PHOTOS VLL;
45+
Enddecay
46+
47+
48+
End
49+

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

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -394,11 +394,11 @@ FairGenerator* GeneratorCocktailPromptCharmoniaToMuonEvtGen_pp13TeV()
394394
auto genCocktailEvtGen = new o2::eventgen::GeneratorEvtGen<GeneratorCocktail>();
395395

396396
auto genJpsi = new o2::eventgen::O2_GeneratorParamJpsiFwdY;
397-
genJpsi->SetNSignalPerEvent(4); // 4 J/psi generated per event by GeneratorParam
397+
genJpsi->SetNSignalPerEvent(1); // 1 J/psi generated per event by GeneratorParam
398398
auto genPsi = new o2::eventgen::O2_GeneratorParamPsiFwdY;
399-
genPsi->SetNSignalPerEvent(2); // 2 Psi(2S) generated per event by GeneratorParam
400-
genCocktailEvtGen->AddGenerator(genJpsi, 1); // 2/3 J/psi
401-
genCocktailEvtGen->AddGenerator(genPsi, 1); // 1/3 Psi(2S)
399+
genPsi->SetNSignalPerEvent(1); // 1 Psi(2S) generated per event by GeneratorParam
400+
genCocktailEvtGen->AddGenerator(genJpsi, 1); // add J/psi generator
401+
genCocktailEvtGen->AddGenerator(genPsi, 1); // add Psi(2S) generator
402402

403403
TString pdgs = "443;100443";
404404
std::string spdg;

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

Lines changed: 79 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -124,7 +124,7 @@ Bool_t Init() override
124124
/// check at hadron level if needed
125125
int ipdg=0;
126126
for (auto& pdgVal : mHadronsPDGs){
127-
if ( (TMath::Abs(event[ipa].id()) == pdgVal) && (event[ipa].y() > mHadronRapidityMin) && (event[ipa].y() < mHadronRapidityMax) ) { countH[ipdg]++; printf("ndaughters %d \n",event[ipa].daughterList().size()); }
127+
if ( (TMath::Abs(event[ipa].id()) == pdgVal) && (event[ipa].y() > mHadronRapidityMin) && (event[ipa].y() < mHadronRapidityMax) ) countH[ipdg]++;
128128
if(isOkAtPartonicLevel && countH[ipdg] >= mHadronMultiplicity) return true;
129129
ipdg++;
130130
}
@@ -258,7 +258,8 @@ FairGenerator*
258258
}
259259

260260
TString pathO2 = gSystem->ExpandPathName("$O2DPG_ROOT/MC/config/PWGDQ/EvtGen/DecayTablesEvtgen");
261-
gen->SetDecayTable(Form("%s/BPLUSTOKAONJPSITOELE.DEC", pathO2.Data()));
261+
//gen->SetDecayTable(Form("%s/BPLUSTOKAONJPSITOELE.DEC", pathO2.Data()));
262+
gen->SetDecayTable(Form("%s/BPLUSTOKAONJPSITOELEALLMODES.DEC", pathO2.Data())); // decay table including decay modes for correlated background
262263
// print debug
263264
// gen->PrintDebug();
264265
// set random seed
@@ -272,5 +273,81 @@ FairGenerator*
272273
return gen;
273274
}
274275

276+
// Predefined generators:
277+
FairGenerator*
278+
GeneratorBeautyToJpsi_EvtGenFwdY(double rapidityMin = -4.3, double rapidityMax = -2.3, bool verbose = false, TString pdgs = "511;521;531;541;5112;5122;5232;5132;5332")
279+
{
280+
auto gen = new o2::eventgen::GeneratorEvtGen<o2::eventgen::GeneratorPythia8NonPromptInjectedGapTriggeredDQ>();
281+
gen->setRapidity(rapidityMin, rapidityMax);
282+
gen->setPDG(5);
283+
gen->setRapidityHadron(rapidityMin,rapidityMax);
284+
gen->setHadronMultiplicity(1);
285+
TString pathO2table = gSystem->ExpandPathName("$O2DPG_ROOT/MC/config/PWGDQ/pythia8/decayer/switchOffBhadrons.cfg");
286+
gen->readFile(pathO2table.Data());
287+
gen->setConfigMBdecays(pathO2table);
288+
gen->setVerbose(verbose);
289+
290+
std::string spdg;
291+
TObjArray* obj = pdgs.Tokenize(";");
292+
gen->SetSizePdg(obj->GetEntriesFast());
293+
for (int i = 0; i < obj->GetEntriesFast(); i++) {
294+
spdg = obj->At(i)->GetName();
295+
gen->AddPdg(std::stoi(spdg), i);
296+
gen->addHadronPDGs(std::stoi(spdg));
297+
printf("PDG %d \n", std::stoi(spdg));
298+
}
299+
gen->SetForceDecay(kEvtBJpsiDiMuon);
300+
301+
// set random seed
302+
gen->readString("Random:setSeed on");
303+
uint random_seed;
304+
unsigned long long int random_value = 0;
305+
ifstream urandom("/dev/urandom", ios::in|ios::binary);
306+
urandom.read(reinterpret_cast<char*>(&random_value), sizeof(random_seed));
307+
gen->readString(Form("Random:seed = %d", random_value % 900000001));
308+
309+
// print debug
310+
// gen->PrintDebug();
311+
312+
return gen;
313+
}
275314

276315

316+
FairGenerator*
317+
GeneratorBeautyToPsiAndJpsi_EvtGenFwdY(double rapidityMin = -4.3, double rapidityMax = -2.3, bool verbose = false, TString pdgs = "511;521;531;541;5112;5122;5232;5132;5332")
318+
{
319+
auto gen = new o2::eventgen::GeneratorEvtGen<o2::eventgen::GeneratorPythia8NonPromptInjectedGapTriggeredDQ>();
320+
gen->setRapidity(rapidityMin, rapidityMax);
321+
gen->setPDG(5);
322+
gen->setRapidityHadron(rapidityMin,rapidityMax);
323+
gen->setHadronMultiplicity(1);
324+
TString pathO2table = gSystem->ExpandPathName("$O2DPG_ROOT/MC/config/PWGDQ/pythia8/decayer/switchOffBhadrons.cfg");
325+
gen->readFile(pathO2table.Data());
326+
gen->setConfigMBdecays(pathO2table);
327+
gen->setVerbose(verbose);
328+
329+
std::string spdg;
330+
TObjArray* obj = pdgs.Tokenize(";");
331+
gen->SetSizePdg(obj->GetEntriesFast());
332+
for (int i = 0; i < obj->GetEntriesFast(); i++) {
333+
spdg = obj->At(i)->GetName();
334+
gen->AddPdg(std::stoi(spdg), i);
335+
printf("PDG %d \n", std::stoi(spdg));
336+
gen->addHadronPDGs(std::stoi(spdg));
337+
}
338+
gen->SetForceDecay(kEvtBPsiAndJpsiDiMuon);
339+
340+
// set random seed
341+
gen->readString("Random:setSeed on");
342+
uint random_seed;
343+
unsigned long long int random_value = 0;
344+
ifstream urandom("/dev/urandom", ios::in|ios::binary);
345+
urandom.read(reinterpret_cast<char*>(&random_value), sizeof(random_seed));
346+
gen->readString(Form("Random:seed = %d", random_value % 900000001));
347+
348+
// print debug
349+
// gen->PrintDebug();
350+
351+
return gen;
352+
}
353+

‎MC/config/PWGDQ/ini/GeneratorHF_bbbar_PsiAndJpsi_fwdy_triggerGap.ini‎

Lines changed: 15 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,15 @@
1+
### The setup uses an external event generator
2+
### This part sets the path of the file and the function call to retrieve it
3+
4+
[GeneratorExternal]
5+
fileName = ${O2DPG_ROOT}/MC/config/PWGDQ/external/generator/generator_pythia8_NonPromptSignals_gaptriggered_dq.C
6+
funcName = GeneratorBeautyToPsiAndJpsi_EvtGenFwdY()
7+
8+
### The external generator derives from GeneratorPythia8.
9+
### This part configures the bits of the interface: configuration and user hooks
10+
11+
[GeneratorPythia8]
12+
config = ${O2DPG_ROOT}/MC/config/common/pythia8/generator/pythia8_hf.cfg
13+
hooksFileName = ${O2DPG_ROOT}/MC/config/PWGHF/pythia8/hooks/pythia8_userhooks_qqbar.C
14+
hooksFuncName = pythia8_userhooks_bbbar(-4.3,-2.3)
15+

‎MC/config/PWGDQ/ini/Generator_InjectedPromptCharmoniaFwdy_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_ROOT}/MC/config/PWGDQ/external/generator/generator_pythia8_withInjectedPromptSignals_gaptriggered_dq.C
4+
funcName=GeneratorPythia8InjectedPromptCharmoniaGapTriggered(5,3)
5+
6+
[GeneratorPythia8]
7+
config=${O2DPG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_triggerGap.cfg

‎MC/config/PWGDQ/ini/Generator_InjectedPromptJpsiFwdy_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_ROOT}/MC/config/PWGDQ/external/generator/generator_pythia8_withInjectedPromptSignals_gaptriggered_dq.C
4+
funcName=GeneratorPythia8InjectedPromptCharmoniaGapTriggered(5,4)
5+
6+
[GeneratorPythia8]
7+
config=${O2DPG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_triggerGap.cfg

‎MC/config/PWGDQ/ini/Generator_InjectedPromptPsi2SFwdy_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_ROOT}/MC/config/PWGDQ/external/generator/generator_pythia8_withInjectedPromptSignals_gaptriggered_dq.C
4+
funcName=GeneratorPythia8InjectedPromptCharmoniaGapTriggered(5,5)
5+
6+
[GeneratorPythia8]
7+
config=${O2DPG_ROOT}/MC/config/PWGDQ/pythia8/generator/pythia8_inel_triggerGap.cfg

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

Lines changed: 93 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,93 @@
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+
Int_t bpdgs[] = {511, 521, 531, 5112, 5122, 5232, 5132};
29+
Int_t sizePdg = sizeof(bpdgs)/sizeof(Int_t);
30+
Bool_t hasBeautyMoth = kFALSE;
31+
32+
for (int i = 0; i < nEvents; i++) {
33+
tree->GetEntry(i);
34+
for (auto& track : *tracks) {
35+
auto pdg = track.GetPdgCode();
36+
auto rapidity = track.GetRapidity();
37+
auto idMoth = track.getMotherTrackId();
38+
if (pdg == checkPdgDecay) {
39+
// count leptons
40+
nLeptons++;
41+
} else if(pdg == -checkPdgDecay) {
42+
// count anti-leptons
43+
nAntileptons++;
44+
} else if (pdg == checkPdgSignal[0] || pdg == checkPdgSignal[1]) {
45+
// check if mothers are beauty hadrons
46+
hasBeautyMoth = kFALSE;
47+
if(idMoth){ // check beauty mother
48+
auto tdM = mcreader.getTrack(i, idMoth);
49+
for(int i=0; i<sizePdg; i++){ if (TMath::Abs(tdM->GetPdgCode()) == bpdgs[i] ) hasBeautyMoth = kTRUE; }
50+
}
51+
if(hasBeautyMoth){
52+
// count signal PDG
53+
pdg == checkPdgSignal[0] ? nSignalJpsi++ : nSignalPsi2S++;
54+
// count signal PDG within acceptance
55+
if(rapidity > rapiditymin && rapidity < rapiditymax) { pdg == checkPdgSignal[0] ? nSignalJpsiWithinAcc++ : nSignalPsi2SWithinAcc++;}
56+
}
57+
auto child0 = o2::mcutils::MCTrackNavigator::getDaughter0(track, *tracks);
58+
auto child1 = o2::mcutils::MCTrackNavigator::getDaughter1(track, *tracks);
59+
if (child0 != nullptr && child1 != nullptr) {
60+
// check for parent-child relations
61+
auto pdg0 = child0->GetPdgCode();
62+
auto pdg1 = child1->GetPdgCode();
63+
std::cout << "First and last children of parent " << checkPdgSignal << " are PDG0: " << pdg0 << " PDG1: " << pdg1 << "\n";
64+
if (std::abs(pdg0) == checkPdgDecay && std::abs(pdg1) == checkPdgDecay && pdg0 == -pdg1) {
65+
nLeptonPairs++;
66+
if (child0->getToBeDone() && child1->getToBeDone()) {
67+
nLeptonPairsToBeDone++;
68+
}
69+
}
70+
}
71+
}
72+
}
73+
}
74+
std::cout << "#events: " << nEvents << "\n"
75+
<< "#leptons: " << nLeptons << "\n"
76+
<< "#antileptons: " << nAntileptons << "\n"
77+
<< "#signal (jpsi <- b): " << nSignalJpsi << "; within acceptance " << rapiditymin << " < y < " << rapiditymax << " : " << nSignalJpsiWithinAcc << "\n"
78+
<< "#signal (psi2S <- b): " << nSignalPsi2S << "; within acceptance " << rapiditymin << " < y < " << rapiditymax << " : " << nSignalPsi2SWithinAcc << "\n"
79+
<< "#lepton pairs: " << nLeptonPairs << "\n"
80+
<< "#lepton pairs to be done: " << nLeptonPairs << "\n";
81+
82+
83+
if (nLeptonPairs == 0 || nLeptons == 0 || nAntileptons == 0) {
84+
std::cerr << "Number of leptons, number of anti-leptons as well as number of lepton pairs should all be greater than 1.\n";
85+
return 1;
86+
}
87+
if (nLeptonPairs != nLeptonPairsToBeDone) {
88+
std::cerr << "The number of lepton pairs should be the same as the number of lepton pairs which should be transported.\n";
89+
return 1;
90+
}
91+
92+
return 0;
93+
}

‎MC/config/PWGDQ/ini/tests/Generator_InjectedPromptCharmoniaFwdy_TriggerGap.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+
}

0 commit comments

Comments
 (0)