Skip to content

Commit fb65caa

Browse files
sawankumawatsawankumawat
committed
modified the test script to print relevant info regarding injected particles and their daughters
1 parent ecbfc88 commit fb65caa

1 file changed

Lines changed: 111 additions & 99 deletions

File tree

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

Lines changed: 111 additions & 99 deletions
Original file line numberDiff line numberDiff line change
@@ -1,61 +1,6 @@
1-
#if !defined(__CLING__) || defined(__ROOTCLING__)
2-
#include "FairGenerator.h"
3-
#include "TDatabasePDG.h"
4-
#include "TFile.h"
5-
#include "TMath.h"
6-
#include "TSystem.h"
7-
#include "TTree.h"
8-
#include "SimulationDataFormat/MCTrack.h"
9-
#include <iostream>
10-
#include <vector>
11-
#include <algorithm>
12-
#endif
13-
14-
// Include the underlying rapidity generator header/macro
15-
#include "generator_pythia8_LF_rapidity.C"
16-
17-
/// Entry point to configure the particle gun generator for resonance simulation
18-
FairGenerator *generatePhiResonanceGun(int pdg = 999999, // Custom PDG or specific target PDG
19-
float ptMin = 0.0,
20-
float ptMax = 50.0,
21-
float yMin = -1.0,
22-
float yMax = 1.0,
23-
std::string pythiaCfg = "${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator/pythia8_inel_136tev.cfg",
24-
int nInject = 3)
25-
{
26-
// Configure particle parameters (PDG, count, ptMin, ptMax, yMin, yMax)
27-
GeneratorPythia8LFRapidity::ConfigContainer cfg(pdg, nInject, ptMin, ptMax, yMin, yMax);
28-
29-
std::vector<GeneratorPythia8LFRapidity::ConfigContainer> cfgVec;
30-
std::vector<GeneratorPythia8LFRapidity::ConfigContainer> cfgVecGenDecayed;
31-
32-
// Let Pythia generator handle decays internally
33-
cfgVecGenDecayed.push_back(cfg);
34-
35-
return generateLFRapidity(cfgVec, cfgVecGenDecayed,
36-
/*injectOnePDGPerEvent=*/true,
37-
/*gapBetweenInjection=*/0,
38-
/*useTrigger=*/false,
39-
/*useRapidity=*/true,
40-
/*pythiaCfgMb=*/pythiaCfg,
41-
/*pythiaCfgSignal=*/"");
42-
}
43-
44-
/// Validation function to analyze o2sim_Kine.root post-simulation
451
int External()
462
{
47-
std::string path{"o2sim_Kine.root"};
48-
int numberOfGapEvents{0};
49-
int numberOfEventsProcessed{0};
50-
int numberOfEventsProcessedWithoutInjection{0};
51-
52-
// Target PDG state and decaying daughters (e.g. Phi -> K+ K-)
53-
std::vector<int> injectedPDGs = {999999};
54-
std::vector<std::vector<int>> decayDaughters = {
55-
{333, 333} // Decaying into phi-phi (PDG 333, 333)
56-
};
57-
58-
auto nInjection = injectedPDGs.size();
3+
const std::string path{"/home/sawan/alice/practice/testMC/PhiPhi/o2sim_Kine.root"};
594

605
TFile file(path.c_str(), "READ");
616
if (file.IsZombie())
@@ -74,83 +19,150 @@ int External()
7419
std::vector<o2::MCTrack> *tracks{};
7520
tree->SetBranchAddress("MCTrack", &tracks);
7621

77-
std::vector<int> nSignal(nInjection, 0);
78-
std::vector<std::vector<int>> nDecays;
79-
std::vector<int> nNotDecayed(nInjection, 0);
22+
// Counters
23+
int nResonance999999 = 0;
24+
int nNotDecayed999999 = 0;
25+
int nPhiFromResonance = 0;
8026

81-
for (size_t i = 0; i < nInjection; i++)
82-
{
83-
nDecays.push_back(std::vector<int>(decayDaughters[i].size(), 0));
84-
}
27+
// Decay counts into K+ K- (PDG 321, -321)
28+
int nKPlusFromResonancePhi = 0;
29+
int nKMinusFromResonancePhi = 0;
30+
31+
int nDirectInjectedPhi = 0;
32+
int nMBPhi = 0;
33+
34+
int nKPlusFromDirectPhi = 0;
35+
int nKMinusFromDirectPhi = 0;
36+
int nKPlusFromMBPhi = 0;
37+
int nKMinusFromMBPhi = 0;
8538

86-
auto nEvents = tree->GetEntries();
87-
bool hasInjection = false;
39+
int numberOfEventsProcessed = 0;
40+
int numberOfEventsProcessedWithoutInjection = 0;
8841

89-
for (int i = 0; i < nEvents; i++)
42+
for (Long64_t i = 0; i < tree->GetEntries(); ++i)
9043
{
91-
hasInjection = false;
92-
numberOfEventsProcessed++;
9344
tree->GetEntry(i);
45+
++numberOfEventsProcessed;
46+
bool hasInjection = false;
9447

95-
for (size_t idxMCTrack = 0; idxMCTrack < tracks->size(); ++idxMCTrack)
48+
for (size_t idx = 0; idx < tracks->size(); ++idx)
9649
{
97-
auto track = tracks->at(idxMCTrack);
98-
auto pdg = track.GetPdgCode();
99-
auto it = std::find(injectedPDGs.begin(), injectedPDGs.end(), pdg);
50+
const auto &track = tracks->at(idx);
51+
const auto pdg = track.GetPdgCode();
10052

101-
if (it != injectedPDGs.end())
53+
// 1. Process Custom Resonance 999999
54+
if (pdg == 999999)
10255
{
103-
int index = std::distance(injectedPDGs.begin(), it);
104-
nSignal[index]++;
56+
++nResonance999999;
57+
hasInjection = true;
10558

10659
if (track.getFirstDaughterTrackId() < 0)
10760
{
108-
nNotDecayed[index]++;
61+
++nNotDecayed999999;
62+
continue;
63+
}
64+
65+
// Loop through daughters of 999999 (Phi mesons)
66+
for (int j = track.getFirstDaughterTrackId(); j <= track.getLastDaughterTrackId(); ++j)
67+
{
68+
const auto &phiTrack = tracks->at(j);
69+
if (phiTrack.GetPdgCode() == 333)
70+
{
71+
++nPhiFromResonance;
72+
73+
// Check daughters of this Phi (granddaughters of 999999)
74+
if (phiTrack.getFirstDaughterTrackId() >= 0)
75+
{
76+
for (int k = phiTrack.getFirstDaughterTrackId(); k <= phiTrack.getLastDaughterTrackId(); ++k)
77+
{
78+
auto grandDauPdg = tracks->at(k).GetPdgCode();
79+
if (grandDauPdg == 321) ++nKPlusFromResonancePhi;
80+
if (grandDauPdg == -321) ++nKMinusFromResonancePhi;
81+
}
82+
}
83+
}
84+
}
85+
}
86+
87+
// 2. Process Phi (333) Mesons
88+
else if (pdg == 333)
89+
{
90+
int motherId = track.getMotherTrackId();
91+
int motherPdg = (motherId >= 0 && motherId < (int)tracks->size()) ? tracks->at(motherId).GetPdgCode() : 0;
92+
93+
// Skip Phi from 999999 here as it was handled above
94+
if (motherPdg == 999999)
95+
{
10996
continue;
11097
}
11198

112-
for (int j{track.getFirstDaughterTrackId()}; j <= track.getLastDaughterTrackId(); ++j)
99+
bool isDirectInjected = (motherId < 0);
100+
101+
if (isDirectInjected)
113102
{
114-
auto pdgDau = tracks->at(j).GetPdgCode();
115-
bool foundDau = false;
103+
++nDirectInjectedPhi;
104+
hasInjection = true;
116105

117-
for (size_t idxDaughter = 0; idxDaughter < decayDaughters[index].size(); ++idxDaughter)
106+
if (track.getFirstDaughterTrackId() >= 0)
118107
{
119-
if (pdgDau == decayDaughters[index][idxDaughter])
108+
for (int j = track.getFirstDaughterTrackId(); j <= track.getLastDaughterTrackId(); ++j)
120109
{
121-
nDecays[index][idxDaughter]++;
122-
foundDau = true;
123-
hasInjection = true;
124-
break;
110+
auto dauPdg = tracks->at(j).GetPdgCode();
111+
if (dauPdg == 321) ++nKPlusFromDirectPhi;
112+
if (dauPdg == -321) ++nKMinusFromDirectPhi;
125113
}
126114
}
127-
if (!foundDau)
115+
}
116+
else
117+
{
118+
// Minimum Bias Phi
119+
++nMBPhi;
120+
121+
if (track.getFirstDaughterTrackId() >= 0)
128122
{
129-
std::cerr << "Decay daughter not found: " << pdg << " -> " << pdgDau << "\n";
123+
for (int j = track.getFirstDaughterTrackId(); j <= track.getLastDaughterTrackId(); ++j)
124+
{
125+
auto dauPdg = tracks->at(j).GetPdgCode();
126+
if (dauPdg == 321) ++nKPlusFromMBPhi;
127+
if (dauPdg == -321) ++nKMinusFromMBPhi;
128+
}
130129
}
131130
}
132131
}
133132
}
133+
134134
if (!hasInjection)
135135
{
136-
numberOfEventsProcessedWithoutInjection++;
136+
++numberOfEventsProcessedWithoutInjection;
137137
}
138138
}
139139

140140
std::cout << "--------------------------------\n";
141-
std::cout << "# Events: " << nEvents << "\n";
142-
for (size_t i = 0; i < nInjection; i++)
143-
{
144-
std::cout << "# Mother PDG " << injectedPDGs[i] << " generated: "
145-
<< nSignal[i] << ", " << nNotDecayed[i] << " did not decay\n";
146-
for (size_t j = 0; j < decayDaughters[i].size(); j++)
147-
{
148-
std::cout << "# Daughter PDG " << decayDaughters[i][j] << ": " << nDecays[i][j] << "\n";
149-
}
150-
}
141+
std::cout << "Total Events Processed: " << tree->GetEntries() << "\n\n";
142+
143+
std::cout << "--- 1. INJECTED RESONANCE (999999) ---\n";
144+
std::cout << "Total Injected Resonance 999999: " << nResonance999999 << "\n";
145+
std::cout << "Resonances not decayed: " << nNotDecayed999999 << "\n";
146+
std::cout << "Daughter Phi (333) produced from 999999: " << nPhiFromResonance << "\n";
147+
std::cout << " -> Decayed to K+: " << nKPlusFromResonancePhi << "\n";
148+
std::cout << " -> Decayed to K-: " << nKMinusFromResonancePhi << "\n\n";
149+
150+
std::cout << "--- 2. DIRECTLY INJECTED PHI (333) ---\n";
151+
std::cout << "Total Directly Injected Phi (333): " << nDirectInjectedPhi << "\n";
152+
std::cout << " -> Decayed to K+: " << nKPlusFromDirectPhi << "\n";
153+
std::cout << " -> Decayed to K-: " << nKMinusFromDirectPhi << "\n\n";
154+
155+
std::cout << "--- 3. MINIMUM BIAS PHI (333) ---\n";
156+
std::cout << "Total Minimum Bias Phi (333): " << nMBPhi << "\n";
157+
std::cout << " -> Decayed to K+: " << nKPlusFromMBPhi << "\n";
158+
std::cout << " -> Decayed to K-: " << nKMinusFromMBPhi << "\n";
151159
std::cout << "--------------------------------\n";
160+
std::cout << "Events processed without signal injection: " << numberOfEventsProcessedWithoutInjection << "\n";
152161

153162
return 0;
154163
}
155164

156-
void GeneratorLF_phiphi2() { External(); }
165+
void GeneratorLF_phiphi2()
166+
{
167+
External();
168+
}

0 commit comments

Comments
 (0)