Skip to content

Commit e9d6c6c

Browse files
fcatalan92fcatalan92
committed
[PWGHF] Add configurations for D2H XiC0, XiC+ and OmegaC MCs (gen only, no gap trigger)
1 parent 73425a3 commit e9d6c6c

8 files changed

+454
-0
lines changed

‎MC/config/PWGHF/ini/GeneratorHF_D2H_bbbar_Mode2_OmegaC.ini‎

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,8 @@
1+
### The external generator derives from GeneratorPythia8.
2+
[GeneratorExternal]
3+
fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/external/generator/generator_pythia8_gaptriggered_hf.C
4+
funcName=GeneratorPythia8GapTriggeredBeauty(1, -1.5, 1.5, -1.5, 1.5, {4332})
5+
6+
[GeneratorPythia8]
7+
config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/generator/pythia8_charmhadronic_with_decays_Mode2.cfg
8+
includePartonEvent=true

‎MC/config/PWGHF/ini/GeneratorHF_D2H_bbbar_Mode2_XiC0_XiCplus.ini‎

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,8 @@
1+
### The external generator derives from GeneratorPythia8.
2+
[GeneratorExternal]
3+
fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/external/generator/generator_pythia8_gaptriggered_hf.C
4+
funcName=GeneratorPythia8GapTriggeredBeauty(1, -1.5, 1.5, -1.5, 1.5, {4132, 4232})
5+
6+
[GeneratorPythia8]
7+
config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/generator/pythia8_charmhadronic_with_decays_Mode2.cfg
8+
includePartonEvent=true

‎MC/config/PWGHF/ini/GeneratorHF_D2H_ccbar_Mode2_OmegaC.ini‎

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,8 @@
1+
### The external generator derives from GeneratorPythia8.
2+
[GeneratorExternal]
3+
fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/external/generator/generator_pythia8_gaptriggered_hf.C
4+
funcName=GeneratorPythia8GapTriggeredCharm(1, -1.5, 1.5, -1.5, 1.5, {4332})
5+
6+
[GeneratorPythia8]
7+
config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/generator/pythia8_charmhadronic_with_decays_Mode2.cfg
8+
includePartonEvent=true

‎MC/config/PWGHF/ini/GeneratorHF_D2H_ccbar_Mode2_XiC0_XiCplus.ini‎

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,8 @@
1+
### The external generator derives from GeneratorPythia8.
2+
[GeneratorExternal]
3+
fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/external/generator/generator_pythia8_gaptriggered_hf.C
4+
funcName=GeneratorPythia8GapTriggeredCharm(1, -1.5, 1.5, -1.5, 1.5, {4132, 4232})
5+
6+
[GeneratorPythia8]
7+
config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/generator/pythia8_charmhadronic_with_decays_Mode2.cfg
8+
includePartonEvent=true

‎MC/config/PWGHF/ini/tests/GeneratorHF_D2H_bbbar_Mode2_OmegaC.C‎

Lines changed: 105 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,105 @@
1+
int External() {
2+
std::string path{"o2sim_Kine.root"};
3+
4+
int checkPdgQuark{5};
5+
float ratioTrigger = 1.; // only enriched events
6+
7+
std::vector<int> checkPdgHadron{4332};
8+
std::map<int, std::vector<std::vector<int>>> checkHadronDecays{ // sorted pdg of daughters
9+
{4332, {{211, 3334}, {211, 3312}}}, // Omegac0
10+
};
11+
12+
TFile file(path.c_str(), "READ");
13+
if (file.IsZombie()) {
14+
std::cerr << "Cannot open ROOT file " << path << "\n";
15+
return 1;
16+
}
17+
18+
auto tree = (TTree *)file.Get("o2sim");
19+
std::vector<o2::MCTrack> *tracks{};
20+
tree->SetBranchAddress("MCTrack", &tracks);
21+
o2::dataformats::MCEventHeader *eventHeader = nullptr;
22+
tree->SetBranchAddress("MCEventHeader.", &eventHeader);
23+
24+
int nEventsMB{}, nEventsInj{};
25+
int nQuarks{}, nSignals{}, nSignalGoodDecay{};
26+
auto nEvents = tree->GetEntries();
27+
28+
for (int i = 0; i < nEvents; i++) {
29+
tree->GetEntry(i);
30+
31+
// check subgenerator information
32+
if (eventHeader->hasInfo(o2::mcgenid::GeneratorProperty::SUBGENERATORID)) {
33+
bool isValid = false;
34+
int subGeneratorId = eventHeader->getInfo<int>(o2::mcgenid::GeneratorProperty::SUBGENERATORID, isValid);
35+
if (subGeneratorId == 0) {
36+
nEventsMB++;
37+
} else if (subGeneratorId == checkPdgQuark) {
38+
nEventsInj++;
39+
}
40+
}
41+
42+
for (auto &track : *tracks) {
43+
auto pdg = track.GetPdgCode();
44+
if (std::abs(pdg) == checkPdgQuark) {
45+
nQuarks++;
46+
continue;
47+
}
48+
if (std::find(checkPdgHadron.begin(), checkPdgHadron.end(), std::abs(pdg)) != checkPdgHadron.end()) { // found signal
49+
nSignals++; // count signal PDG
50+
51+
std::vector<int> pdgsDecay{};
52+
std::vector<int> pdgsDecayAntiPart{};
53+
for (int j{track.getFirstDaughterTrackId()}; j <= track.getLastDaughterTrackId(); ++j) {
54+
auto pdgDau = tracks->at(j).GetPdgCode();
55+
pdgsDecay.push_back(pdgDau);
56+
if (pdgDau != 333) { // phi is antiparticle of itself
57+
pdgsDecayAntiPart.push_back(-pdgDau);
58+
} else {
59+
pdgsDecayAntiPart.push_back(pdgDau);
60+
}
61+
}
62+
63+
std::sort(pdgsDecay.begin(), pdgsDecay.end());
64+
std::sort(pdgsDecayAntiPart.begin(), pdgsDecayAntiPart.end());
65+
66+
for (auto &decay : checkHadronDecays[std::abs(pdg)]) {
67+
if (pdgsDecay == decay || pdgsDecayAntiPart == decay) {
68+
nSignalGoodDecay++;
69+
break;
70+
}
71+
}
72+
}
73+
}
74+
}
75+
76+
std::cout << "--------------------------------\n";
77+
std::cout << "# Events: " << nEvents << "\n";
78+
std::cout << "# MB events: " << nEventsMB << "\n";
79+
std::cout << Form("# events injected with %d quark pair: ", checkPdgQuark) << nEventsInj << "\n";
80+
std::cout << Form("# %d (anti)quarks: ", checkPdgQuark) << nQuarks << "\n";
81+
std::cout <<"# signal hadrons: " << nSignals << "\n";
82+
std::cout <<"# signal hadrons decaying in the correct channel: " << nSignalGoodDecay << "\n";
83+
84+
if (nEventsMB < nEvents * (1 - ratioTrigger) * 0.95 || nEventsMB > nEvents * (1 - ratioTrigger) * 1.05) { // we put some tolerance since the number of generated events is small
85+
std::cerr << "Number of generated MB events different than expected\n";
86+
return 1;
87+
}
88+
if (nEventsInj < nEvents * ratioTrigger * 0.95 || nEventsInj > nEvents * ratioTrigger * 1.05) {
89+
std::cerr << "Number of generated events injected with " << checkPdgQuark << " different than expected\n";
90+
return 1;
91+
}
92+
93+
if (nQuarks < 2 * nEvents * ratioTrigger) { // we expect anyway more because the same quark is repeated several time, after each gluon radiation
94+
std::cerr << "Number of generated (anti)quarks " << checkPdgQuark << " lower than expected\n";
95+
return 1;
96+
}
97+
98+
float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
99+
if (fracForcedDecays < 0.9) { // we put some tolerance (e.g. due to oscillations which might change the final state)
100+
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
101+
return 1;
102+
}
103+
104+
return 0;
105+
}

‎MC/config/PWGHF/ini/tests/GeneratorHF_D2H_bbbar_Mode2_XiC0_XiCplus.C‎

Lines changed: 106 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,106 @@
1+
int External() {
2+
std::string path{"o2sim_Kine.root"};
3+
4+
int checkPdgQuark{5};
5+
float ratioTrigger = 1.; // only enriched events
6+
7+
std::vector<int> checkPdgHadron{4132, 4232};
8+
std::map<int, std::vector<std::vector<int>>> checkHadronDecays{ // sorted pdg of daughters
9+
{4132, {{211, 3312}}}, // Xic0
10+
{4232, {{-313, 2212}, {-321, 3324}, {211, 211, 3312}, {-321, 211, 2212}}} // Xic+
11+
};
12+
13+
TFile file(path.c_str(), "READ");
14+
if (file.IsZombie()) {
15+
std::cerr << "Cannot open ROOT file " << path << "\n";
16+
return 1;
17+
}
18+
19+
auto tree = (TTree *)file.Get("o2sim");
20+
std::vector<o2::MCTrack> *tracks{};
21+
tree->SetBranchAddress("MCTrack", &tracks);
22+
o2::dataformats::MCEventHeader *eventHeader = nullptr;
23+
tree->SetBranchAddress("MCEventHeader.", &eventHeader);
24+
25+
int nEventsMB{}, nEventsInj{};
26+
int nQuarks{}, nSignals{}, nSignalGoodDecay{};
27+
auto nEvents = tree->GetEntries();
28+
29+
for (int i = 0; i < nEvents; i++) {
30+
tree->GetEntry(i);
31+
32+
// check subgenerator information
33+
if (eventHeader->hasInfo(o2::mcgenid::GeneratorProperty::SUBGENERATORID)) {
34+
bool isValid = false;
35+
int subGeneratorId = eventHeader->getInfo<int>(o2::mcgenid::GeneratorProperty::SUBGENERATORID, isValid);
36+
if (subGeneratorId == 0) {
37+
nEventsMB++;
38+
} else if (subGeneratorId == checkPdgQuark) {
39+
nEventsInj++;
40+
}
41+
}
42+
43+
for (auto &track : *tracks) {
44+
auto pdg = track.GetPdgCode();
45+
if (std::abs(pdg) == checkPdgQuark) {
46+
nQuarks++;
47+
continue;
48+
}
49+
if (std::find(checkPdgHadron.begin(), checkPdgHadron.end(), std::abs(pdg)) != checkPdgHadron.end()) { // found signal
50+
nSignals++; // count signal PDG
51+
52+
std::vector<int> pdgsDecay{};
53+
std::vector<int> pdgsDecayAntiPart{};
54+
for (int j{track.getFirstDaughterTrackId()}; j <= track.getLastDaughterTrackId(); ++j) {
55+
auto pdgDau = tracks->at(j).GetPdgCode();
56+
pdgsDecay.push_back(pdgDau);
57+
if (pdgDau != 333) { // phi is antiparticle of itself
58+
pdgsDecayAntiPart.push_back(-pdgDau);
59+
} else {
60+
pdgsDecayAntiPart.push_back(pdgDau);
61+
}
62+
}
63+
64+
std::sort(pdgsDecay.begin(), pdgsDecay.end());
65+
std::sort(pdgsDecayAntiPart.begin(), pdgsDecayAntiPart.end());
66+
67+
for (auto &decay : checkHadronDecays[std::abs(pdg)]) {
68+
if (pdgsDecay == decay || pdgsDecayAntiPart == decay) {
69+
nSignalGoodDecay++;
70+
break;
71+
}
72+
}
73+
}
74+
}
75+
}
76+
77+
std::cout << "--------------------------------\n";
78+
std::cout << "# Events: " << nEvents << "\n";
79+
std::cout << "# MB events: " << nEventsMB << "\n";
80+
std::cout << Form("# events injected with %d quark pair: ", checkPdgQuark) << nEventsInj << "\n";
81+
std::cout << Form("# %d (anti)quarks: ", checkPdgQuark) << nQuarks << "\n";
82+
std::cout <<"# signal hadrons: " << nSignals << "\n";
83+
std::cout <<"# signal hadrons decaying in the correct channel: " << nSignalGoodDecay << "\n";
84+
85+
if (nEventsMB < nEvents * (1 - ratioTrigger) * 0.95 || nEventsMB > nEvents * (1 - ratioTrigger) * 1.05) { // we put some tolerance since the number of generated events is small
86+
std::cerr << "Number of generated MB events different than expected\n";
87+
return 1;
88+
}
89+
if (nEventsInj < nEvents * ratioTrigger * 0.95 || nEventsInj > nEvents * ratioTrigger * 1.05) {
90+
std::cerr << "Number of generated events injected with " << checkPdgQuark << " different than expected\n";
91+
return 1;
92+
}
93+
94+
if (nQuarks < 2 * nEvents * ratioTrigger) { // we expect anyway more because the same quark is repeated several time, after each gluon radiation
95+
std::cerr << "Number of generated (anti)quarks " << checkPdgQuark << " lower than expected\n";
96+
return 1;
97+
}
98+
99+
float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
100+
if (fracForcedDecays < 0.9) { // we put some tolerance (e.g. due to oscillations which might change the final state)
101+
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
102+
return 1;
103+
}
104+
105+
return 0;
106+
}

‎MC/config/PWGHF/ini/tests/GeneratorHF_D2H_ccbar_Mode2_OmegaC.C‎

Lines changed: 105 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,105 @@
1+
int External() {
2+
std::string path{"o2sim_Kine.root"};
3+
4+
int checkPdgQuark{4};
5+
float ratioTrigger = 1.; // only enriched events
6+
7+
std::vector<int> checkPdgHadron{4332};
8+
std::map<int, std::vector<std::vector<int>>> checkHadronDecays{ // sorted pdg of daughters
9+
{4332, {{211, 3334}, {211, 3312}}}, // Omegac0
10+
};
11+
12+
TFile file(path.c_str(), "READ");
13+
if (file.IsZombie()) {
14+
std::cerr << "Cannot open ROOT file " << path << "\n";
15+
return 1;
16+
}
17+
18+
auto tree = (TTree *)file.Get("o2sim");
19+
std::vector<o2::MCTrack> *tracks{};
20+
tree->SetBranchAddress("MCTrack", &tracks);
21+
o2::dataformats::MCEventHeader *eventHeader = nullptr;
22+
tree->SetBranchAddress("MCEventHeader.", &eventHeader);
23+
24+
int nEventsMB{}, nEventsInj{};
25+
int nQuarks{}, nSignals{}, nSignalGoodDecay{};
26+
auto nEvents = tree->GetEntries();
27+
28+
for (int i = 0; i < nEvents; i++) {
29+
tree->GetEntry(i);
30+
31+
// check subgenerator information
32+
if (eventHeader->hasInfo(o2::mcgenid::GeneratorProperty::SUBGENERATORID)) {
33+
bool isValid = false;
34+
int subGeneratorId = eventHeader->getInfo<int>(o2::mcgenid::GeneratorProperty::SUBGENERATORID, isValid);
35+
if (subGeneratorId == 0) {
36+
nEventsMB++;
37+
} else if (subGeneratorId == checkPdgQuark) {
38+
nEventsInj++;
39+
}
40+
}
41+
42+
for (auto &track : *tracks) {
43+
auto pdg = track.GetPdgCode();
44+
if (std::abs(pdg) == checkPdgQuark) {
45+
nQuarks++;
46+
continue;
47+
}
48+
if (std::find(checkPdgHadron.begin(), checkPdgHadron.end(), std::abs(pdg)) != checkPdgHadron.end()) { // found signal
49+
nSignals++; // count signal PDG
50+
51+
std::vector<int> pdgsDecay{};
52+
std::vector<int> pdgsDecayAntiPart{};
53+
for (int j{track.getFirstDaughterTrackId()}; j <= track.getLastDaughterTrackId(); ++j) {
54+
auto pdgDau = tracks->at(j).GetPdgCode();
55+
pdgsDecay.push_back(pdgDau);
56+
if (pdgDau != 333) { // phi is antiparticle of itself
57+
pdgsDecayAntiPart.push_back(-pdgDau);
58+
} else {
59+
pdgsDecayAntiPart.push_back(pdgDau);
60+
}
61+
}
62+
63+
std::sort(pdgsDecay.begin(), pdgsDecay.end());
64+
std::sort(pdgsDecayAntiPart.begin(), pdgsDecayAntiPart.end());
65+
66+
for (auto &decay : checkHadronDecays[std::abs(pdg)]) {
67+
if (pdgsDecay == decay || pdgsDecayAntiPart == decay) {
68+
nSignalGoodDecay++;
69+
break;
70+
}
71+
}
72+
}
73+
}
74+
}
75+
76+
std::cout << "--------------------------------\n";
77+
std::cout << "# Events: " << nEvents << "\n";
78+
std::cout << "# MB events: " << nEventsMB << "\n";
79+
std::cout << Form("# events injected with %d quark pair: ", checkPdgQuark) << nEventsInj << "\n";
80+
std::cout << Form("# %d (anti)quarks: ", checkPdgQuark) << nQuarks << "\n";
81+
std::cout <<"# signal hadrons: " << nSignals << "\n";
82+
std::cout <<"# signal hadrons decaying in the correct channel: " << nSignalGoodDecay << "\n";
83+
84+
if (nEventsMB < nEvents * (1 - ratioTrigger) * 0.95 || nEventsMB > nEvents * (1 - ratioTrigger) * 1.05) { // we put some tolerance since the number of generated events is small
85+
std::cerr << "Number of generated MB events different than expected\n";
86+
return 1;
87+
}
88+
if (nEventsInj < nEvents * ratioTrigger * 0.95 || nEventsInj > nEvents * ratioTrigger * 1.05) {
89+
std::cerr << "Number of generated events injected with " << checkPdgQuark << " different than expected\n";
90+
return 1;
91+
}
92+
93+
if (nQuarks < 2 * nEvents * ratioTrigger) { // we expect anyway more because the same quark is repeated several time, after each gluon radiation
94+
std::cerr << "Number of generated (anti)quarks " << checkPdgQuark << " lower than expected\n";
95+
return 1;
96+
}
97+
98+
float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
99+
if (fracForcedDecays < 0.9) { // we put some tolerance (e.g. due to oscillations which might change the final state)
100+
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
101+
return 1;
102+
}
103+
104+
return 0;
105+
}

0 commit comments

Comments
 (0)