Skip to content

Commit 21f45c8

Browse files
nasirmehdimaliknasirmehdimalik
authored
Add files via upload
1 parent df7b3e4 commit 21f45c8

File tree

4 files changed

+568
-0
lines changed

4 files changed

+568
-0
lines changed

‎MC/config/PWGLF/ini/tests/GeneratorLF_ResonancesBaryonic_NeNe_injection.C‎

Lines changed: 142 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,142 @@
1+
int External()
2+
{
3+
std::string path{"o2sim_Kine.root"};
4+
int numberOfInjectedSignalsPerEvent{1};
5+
int numberOfGapEvents{4};
6+
int numberOfEventsProcessed{0};
7+
int numberOfEventsProcessedWithoutInjection{0};
8+
std::vector<int> injectedPDGs = {
9+
102134, // Lambda(1520)0
10+
-102134, // Lambda(1520)0bar
11+
3324, // Xi(1530)0
12+
-3324 // Xi(1530)0bar
13+
};
14+
std::vector<std::vector<int>> decayDaughters = {
15+
{2212, -321}, // Lambda(1520)0
16+
{-2212, 321}, // Lambda(1520)0bar
17+
{3312, 211}, // Xi(1530)0
18+
{-3312, -211} // Xi(1530)0bar
19+
};
20+
21+
auto nInjection = injectedPDGs.size();
22+
23+
TFile file(path.c_str(), "READ");
24+
if (file.IsZombie())
25+
{
26+
std::cerr << "Cannot open ROOT file " << path << "\n";
27+
return 1;
28+
}
29+
30+
auto tree = (TTree *)file.Get("o2sim");
31+
if (!tree)
32+
{
33+
std::cerr << "Cannot find tree o2sim in file " << path << "\n";
34+
return 1;
35+
}
36+
std::vector<o2::MCTrack> *tracks{};
37+
tree->SetBranchAddress("MCTrack", &tracks);
38+
39+
std::vector<int> nSignal;
40+
for (int i = 0; i < nInjection; i++)
41+
{
42+
nSignal.push_back(0);
43+
}
44+
std::vector<std::vector<int>> nDecays;
45+
std::vector<int> nNotDecayed;
46+
for (int i = 0; i < nInjection; i++)
47+
{
48+
std::vector<int> nDecay;
49+
for (int j = 0; j < decayDaughters[i].size(); j++)
50+
{
51+
nDecay.push_back(0);
52+
}
53+
nDecays.push_back(nDecay);
54+
nNotDecayed.push_back(0);
55+
}
56+
auto nEvents = tree->GetEntries();
57+
bool hasInjection = false;
58+
for (int i = 0; i < nEvents; i++)
59+
{
60+
hasInjection = false;
61+
numberOfEventsProcessed++;
62+
auto check = tree->GetEntry(i);
63+
for (int idxMCTrack = 0; idxMCTrack < tracks->size(); ++idxMCTrack)
64+
{
65+
auto track = tracks->at(idxMCTrack);
66+
auto pdg = track.GetPdgCode();
67+
auto it = std::find(injectedPDGs.begin(), injectedPDGs.end(), pdg);
68+
int index = std::distance(injectedPDGs.begin(), it); // index of injected PDG
69+
if (it != injectedPDGs.end()) // found
70+
{
71+
// count signal PDG
72+
nSignal[index]++;
73+
if (track.getFirstDaughterTrackId() < 0)
74+
{
75+
nNotDecayed[index]++;
76+
continue;
77+
}
78+
for (int j{track.getFirstDaughterTrackId()}; j <= track.getLastDaughterTrackId(); ++j)
79+
{
80+
auto pdgDau = tracks->at(j).GetPdgCode();
81+
bool foundDau = false;
82+
// count decay PDGs
83+
for (int idxDaughter = 0; idxDaughter < decayDaughters[index].size(); ++idxDaughter)
84+
{
85+
if (pdgDau == decayDaughters[index][idxDaughter])
86+
{
87+
nDecays[index][idxDaughter]++;
88+
foundDau = true;
89+
hasInjection = true;
90+
break;
91+
}
92+
}
93+
if (!foundDau)
94+
{
95+
std::cerr << "Decay daughter not found: " << pdg << " -> " << pdgDau << "\n";
96+
}
97+
}
98+
}
99+
}
100+
if (!hasInjection)
101+
{
102+
numberOfEventsProcessedWithoutInjection++;
103+
}
104+
}
105+
std::cout << "--------------------------------\n";
106+
std::cout << "# Events: " << nEvents << "\n";
107+
for (int i = 0; i < nInjection; i++)
108+
{
109+
std::cout << "# Mother \n";
110+
std::cout << injectedPDGs[i] << " generated: " << nSignal[i] << ", " << nNotDecayed[i] << " did not decay\n";
111+
if (nSignal[i] == 0)
112+
{
113+
std::cerr << "No generated: " << injectedPDGs[i] << "\n";
114+
// return 1; // At least one of the injected particles should be generated
115+
}
116+
for (int j = 0; j < decayDaughters[i].size(); j++)
117+
{
118+
std::cout << "# Daughter " << decayDaughters[i][j] << ": " << nDecays[i][j] << "\n";
119+
}
120+
// if (nSignal[i] != nEvents * numberOfInjectedSignalsPerEvent)
121+
// {
122+
// std::cerr << "Number of generated: " << injectedPDGs[i] << ", lower than expected\n";
123+
// // return 1; // Don't need to return 1, since the number of generated particles is not the same for each event
124+
// }
125+
}
126+
std::cout << "--------------------------------\n";
127+
std::cout << "Number of events processed: " << numberOfEventsProcessed << "\n";
128+
std::cout << "Number of input for the gap events: " << numberOfGapEvents << "\n";
129+
std::cout << "Number of events processed without injection: " << numberOfEventsProcessedWithoutInjection << "\n";
130+
// injected event + numberOfGapEvents*gap events + injected event + numberOfGapEvents*gap events + ...
131+
// total fraction of the gap event: numberOfEventsProcessedWithoutInjection/numberOfEventsProcessed
132+
float ratioOfNormalEvents = numberOfEventsProcessedWithoutInjection / numberOfEventsProcessed;
133+
if (ratioOfNormalEvents > 0.75)
134+
{
135+
std::cout << "The number of injected event is loo low!!" << std::endl;
136+
return 1;
137+
}
138+
139+
return 0;
140+
}
141+
142+
void GeneratorLF_ResonancesBaryonic_NeNe_injection() { External(); }

‎MC/config/PWGLF/ini/tests/GeneratorLF_ResonancesBaryonic_OO_injection.C‎

Lines changed: 142 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,142 @@
1+
int External()
2+
{
3+
std::string path{"o2sim_Kine.root"};
4+
int numberOfInjectedSignalsPerEvent{1};
5+
int numberOfGapEvents{4};
6+
int numberOfEventsProcessed{0};
7+
int numberOfEventsProcessedWithoutInjection{0};
8+
std::vector<int> injectedPDGs = {
9+
102134, // Lambda(1520)0
10+
-102134, // Lambda(1520)0bar
11+
3324, // Xi(1530)0
12+
-3324 // Xi(1530)0bar
13+
};
14+
std::vector<std::vector<int>> decayDaughters = {
15+
{2212, -321}, // Lambda(1520)0
16+
{-2212, 321}, // Lambda(1520)0bar
17+
{3312, 211}, // Xi(1530)0
18+
{-3312, -211} // Xi(1530)0bar
19+
};
20+
21+
auto nInjection = injectedPDGs.size();
22+
23+
TFile file(path.c_str(), "READ");
24+
if (file.IsZombie())
25+
{
26+
std::cerr << "Cannot open ROOT file " << path << "\n";
27+
return 1;
28+
}
29+
30+
auto tree = (TTree *)file.Get("o2sim");
31+
if (!tree)
32+
{
33+
std::cerr << "Cannot find tree o2sim in file " << path << "\n";
34+
return 1;
35+
}
36+
std::vector<o2::MCTrack> *tracks{};
37+
tree->SetBranchAddress("MCTrack", &tracks);
38+
39+
std::vector<int> nSignal;
40+
for (int i = 0; i < nInjection; i++)
41+
{
42+
nSignal.push_back(0);
43+
}
44+
std::vector<std::vector<int>> nDecays;
45+
std::vector<int> nNotDecayed;
46+
for (int i = 0; i < nInjection; i++)
47+
{
48+
std::vector<int> nDecay;
49+
for (int j = 0; j < decayDaughters[i].size(); j++)
50+
{
51+
nDecay.push_back(0);
52+
}
53+
nDecays.push_back(nDecay);
54+
nNotDecayed.push_back(0);
55+
}
56+
auto nEvents = tree->GetEntries();
57+
bool hasInjection = false;
58+
for (int i = 0; i < nEvents; i++)
59+
{
60+
hasInjection = false;
61+
numberOfEventsProcessed++;
62+
auto check = tree->GetEntry(i);
63+
for (int idxMCTrack = 0; idxMCTrack < tracks->size(); ++idxMCTrack)
64+
{
65+
auto track = tracks->at(idxMCTrack);
66+
auto pdg = track.GetPdgCode();
67+
auto it = std::find(injectedPDGs.begin(), injectedPDGs.end(), pdg);
68+
int index = std::distance(injectedPDGs.begin(), it); // index of injected PDG
69+
if (it != injectedPDGs.end()) // found
70+
{
71+
// count signal PDG
72+
nSignal[index]++;
73+
if (track.getFirstDaughterTrackId() < 0)
74+
{
75+
nNotDecayed[index]++;
76+
continue;
77+
}
78+
for (int j{track.getFirstDaughterTrackId()}; j <= track.getLastDaughterTrackId(); ++j)
79+
{
80+
auto pdgDau = tracks->at(j).GetPdgCode();
81+
bool foundDau = false;
82+
// count decay PDGs
83+
for (int idxDaughter = 0; idxDaughter < decayDaughters[index].size(); ++idxDaughter)
84+
{
85+
if (pdgDau == decayDaughters[index][idxDaughter])
86+
{
87+
nDecays[index][idxDaughter]++;
88+
foundDau = true;
89+
hasInjection = true;
90+
break;
91+
}
92+
}
93+
if (!foundDau)
94+
{
95+
std::cerr << "Decay daughter not found: " << pdg << " -> " << pdgDau << "\n";
96+
}
97+
}
98+
}
99+
}
100+
if (!hasInjection)
101+
{
102+
numberOfEventsProcessedWithoutInjection++;
103+
}
104+
}
105+
std::cout << "--------------------------------\n";
106+
std::cout << "# Events: " << nEvents << "\n";
107+
for (int i = 0; i < nInjection; i++)
108+
{
109+
std::cout << "# Mother \n";
110+
std::cout << injectedPDGs[i] << " generated: " << nSignal[i] << ", " << nNotDecayed[i] << " did not decay\n";
111+
if (nSignal[i] == 0)
112+
{
113+
std::cerr << "No generated: " << injectedPDGs[i] << "\n";
114+
// return 1; // At least one of the injected particles should be generated
115+
}
116+
for (int j = 0; j < decayDaughters[i].size(); j++)
117+
{
118+
std::cout << "# Daughter " << decayDaughters[i][j] << ": " << nDecays[i][j] << "\n";
119+
}
120+
// if (nSignal[i] != nEvents * numberOfInjectedSignalsPerEvent)
121+
// {
122+
// std::cerr << "Number of generated: " << injectedPDGs[i] << ", lower than expected\n";
123+
// // return 1; // Don't need to return 1, since the number of generated particles is not the same for each event
124+
// }
125+
}
126+
std::cout << "--------------------------------\n";
127+
std::cout << "Number of events processed: " << numberOfEventsProcessed << "\n";
128+
std::cout << "Number of input for the gap events: " << numberOfGapEvents << "\n";
129+
std::cout << "Number of events processed without injection: " << numberOfEventsProcessedWithoutInjection << "\n";
130+
// injected event + numberOfGapEvents*gap events + injected event + numberOfGapEvents*gap events + ...
131+
// total fraction of the gap event: numberOfEventsProcessedWithoutInjection/numberOfEventsProcessed
132+
float ratioOfNormalEvents = numberOfEventsProcessedWithoutInjection / numberOfEventsProcessed;
133+
if (ratioOfNormalEvents > 0.75)
134+
{
135+
std::cout << "The number of injected event is loo low!!" << std::endl;
136+
return 1;
137+
}
138+
139+
return 0;
140+
}
141+
142+
void GeneratorLF_ResonancesBaryonic_OO_injection() { External(); }

0 commit comments

Comments
 (0)