Skip to content

Commit dff5285

Browse files
NicoleBastidNicoleBastid
authored
Add files via upload
test MC of semimuonic decays at kinematics
1 parent db136d7 commit dff5285

File tree

2 files changed

+206
-0
lines changed

2 files changed

+206
-0
lines changed

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

Lines changed: 103 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,103 @@
1+
int External() {
2+
3+
int checkPdgDecayMuon = 13;
4+
int checkPdgQuark = 5;
5+
6+
float ratioTrigger = 1. / 5; // one event triggered out of 5
7+
8+
std::string path{"o2sim_Kine.root"};
9+
10+
TFile file(path.c_str(), "READ");
11+
if (file.IsZombie()) {
12+
std::cerr << "Cannot open ROOT file" << path << "\n";
13+
return 1;
14+
}
15+
16+
auto tree = (TTree *)file.Get("o2sim");
17+
if (!tree) {
18+
std::cerr << "Cannot find tree o2sim in file" << path << "\n";
19+
return 1;
20+
}
21+
22+
std::vector<o2::MCTrack> *tracks{};
23+
tree->SetBranchAddress("MCTrack", &tracks);
24+
25+
o2::dataformats::MCEventHeader *eventHeader = nullptr;
26+
tree->SetBranchAddress("MCEventHeader.", &eventHeader);
27+
28+
int nEventsMB{};
29+
int nEventsInj{};
30+
int nQuarks{};
31+
int nMuons{};
32+
33+
int nMuonsInAcceptance{};
34+
35+
auto nEvents = tree->GetEntries();
36+
37+
for (int i = 0; i < nEvents; i++) {
38+
tree->GetEntry(i);
39+
// check subgenerator information
40+
if (eventHeader->hasInfo(o2::mcgenid::GeneratorProperty::SUBGENERATORID)) {
41+
bool isValid = false;
42+
int subGeneratorId = eventHeader->getInfo<int>(
43+
o2::mcgenid::GeneratorProperty::SUBGENERATORID, isValid);
44+
if (subGeneratorId == 0) {
45+
nEventsMB++;
46+
} else if (subGeneratorId == checkPdgQuark) {
47+
nEventsInj++;
48+
}
49+
} // if event header
50+
51+
int nmuonsev = 0;
52+
int nmuonsevinacc = 0;
53+
54+
for (auto &track : *tracks) {
55+
auto pdg = track.GetPdgCode();
56+
if (std::abs(pdg) == checkPdgQuark) {
57+
nQuarks++;
58+
continue;
59+
} // pdgquark
60+
auto y = track.GetRapidity();
61+
if (std::abs(pdg) == checkPdgDecayMuon) {
62+
int igmother = track.getMotherTrackId();
63+
auto gmTrack = (*tracks)[igmother];
64+
int gmpdg = gmTrack.GetPdgCode();
65+
if (int(std::abs(gmpdg) / 100.) == 5 ||
66+
int(std::abs(gmpdg) / 1000.) == 5) {
67+
nMuons++;
68+
nmuonsev++;
69+
if (-4.3 < y && y < -2.2) {
70+
nMuonsInAcceptance++;
71+
nmuonsevinacc++;
72+
}
73+
} // gmpdg
74+
75+
} // pdgdecay
76+
77+
} // loop track
78+
// std::cout << "#muons per event: " << nmuonsev << "\n";
79+
// std::cout << "#muons in acceptance per event: " << nmuonsev << "\n";
80+
} // events
81+
82+
std::cout << "#events: " << nEvents << "\n";
83+
std::cout << "# MB events: " << nEventsMB << "\n";
84+
std::cout << Form("# events injected with %d quark pair: ", checkPdgQuark)
85+
<< nEventsInj << "\n";
86+
if (nEventsMB < nEvents * (1 - ratioTrigger) * 0.95 ||
87+
nEventsMB > nEvents * (1 - ratioTrigger) *
88+
1.05) { // we put some tolerance since the number of
89+
// generated events is small
90+
std::cerr << "Number of generated MB events different than expected\n";
91+
return 1;
92+
}
93+
if (nEventsInj < nEvents * ratioTrigger * 0.95 ||
94+
nEventsInj > nEvents * ratioTrigger * 1.05) {
95+
std::cerr << "Number of generated events injected with " << checkPdgQuark
96+
<< " different than expected\n";
97+
return 1;
98+
}
99+
std::cout << "#muons: " << nMuons << "\n";
100+
std::cout << "#muons in acceptance: " << nMuonsInAcceptance << "\n";
101+
102+
return 0;
103+
} // external

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

Lines changed: 103 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,103 @@
1+
int External() {
2+
3+
int checkPdgDecayMuon = 13;
4+
int checkPdgQuark = 4;
5+
6+
float ratioTrigger = 1. / 5; // one event triggered out of 5
7+
8+
std::string path{"o2sim_Kine.root"};
9+
10+
TFile file(path.c_str(), "READ");
11+
if (file.IsZombie()) {
12+
std::cerr << "Cannot open ROOT file" << path << "\n";
13+
return 1;
14+
}
15+
16+
auto tree = (TTree *)file.Get("o2sim");
17+
if (!tree) {
18+
std::cerr << "Cannot find tree o2sim in file" << path << "\n";
19+
return 1;
20+
}
21+
22+
std::vector<o2::MCTrack> *tracks{};
23+
tree->SetBranchAddress("MCTrack", &tracks);
24+
25+
o2::dataformats::MCEventHeader *eventHeader = nullptr;
26+
tree->SetBranchAddress("MCEventHeader.", &eventHeader);
27+
28+
int nEventsMB{};
29+
int nEventsInj{};
30+
int nQuarks{};
31+
int nMuons{};
32+
33+
int nMuonsInAcceptance{};
34+
35+
auto nEvents = tree->GetEntries();
36+
37+
for (int i = 0; i < nEvents; i++) {
38+
tree->GetEntry(i);
39+
// check subgenerator information
40+
if (eventHeader->hasInfo(o2::mcgenid::GeneratorProperty::SUBGENERATORID)) {
41+
bool isValid = false;
42+
int subGeneratorId = eventHeader->getInfo<int>(
43+
o2::mcgenid::GeneratorProperty::SUBGENERATORID, isValid);
44+
if (subGeneratorId == 0) {
45+
nEventsMB++;
46+
} else if (subGeneratorId == checkPdgQuark) {
47+
nEventsInj++;
48+
}
49+
} // if event header
50+
51+
int nmuonsev = 0;
52+
int nmuonsevinacc = 0;
53+
54+
for (auto &track : *tracks) {
55+
auto pdg = track.GetPdgCode();
56+
if (std::abs(pdg) == checkPdgQuark) {
57+
nQuarks++;
58+
continue;
59+
} // pdgquark
60+
auto y = track.GetRapidity();
61+
if (std::abs(pdg) == checkPdgDecayMuon) {
62+
int igmother = track.getMotherTrackId();
63+
auto gmTrack = (*tracks)[igmother];
64+
int gmpdg = gmTrack.GetPdgCode();
65+
if (int(std::abs(gmpdg) / 100.) == 5 ||
66+
int(std::abs(gmpdg) / 1000.) == 5) {
67+
nMuons++;
68+
nmuonsev++;
69+
if (-4.3 < y && y < -2.2) {
70+
nMuonsInAcceptance++;
71+
nmuonsevinacc++;
72+
}
73+
} // gmpdg
74+
75+
} // pdgdecay
76+
77+
} // loop track
78+
// std::cout << "#muons per event: " << nmuonsev << "\n";
79+
// std::cout << "#muons in acceptance per event: " << nmuonsev << "\n";
80+
} // events
81+
82+
std::cout << "#events: " << nEvents << "\n";
83+
std::cout << "# MB events: " << nEventsMB << "\n";
84+
std::cout << Form("# events injected with %d quark pair: ", checkPdgQuark)
85+
<< nEventsInj << "\n";
86+
if (nEventsMB < nEvents * (1 - ratioTrigger) * 0.95 ||
87+
nEventsMB > nEvents * (1 - ratioTrigger) *
88+
1.05) { // we put some tolerance since the number of
89+
// generated events is small
90+
std::cerr << "Number of generated MB events different than expected\n";
91+
return 1;
92+
}
93+
if (nEventsInj < nEvents * ratioTrigger * 0.95 ||
94+
nEventsInj > nEvents * ratioTrigger * 1.05) {
95+
std::cerr << "Number of generated events injected with " << checkPdgQuark
96+
<< " different than expected\n";
97+
return 1;
98+
}
99+
std::cout << "#muons: " << nMuons << "\n";
100+
std::cout << "#muons in acceptance: " << nMuonsInAcceptance << "\n";
101+
102+
return 0;
103+
} // external

0 commit comments

Comments
 (0)