Skip to content

Commit 225feb3

Browse files
Yunseul-BaeYunseul-Bae
committed
[PWGLF] Add selection QA mode and update MC rec/true
- Add QA mode and histograms for track/PID/DCA selections - Add MC-true selection steps (cut-indexed GEN counting) - Add centrality table usage for MC true
1 parent 7b487e9 commit 225feb3

File tree

1 file changed

+140
-52
lines changed

1 file changed

+140
-52
lines changed

‎PWGLF/Tasks/Resonances/f0980analysis.cxx‎

Lines changed: 140 additions & 52 deletions
Original file line numberDiff line numberDiff line change
@@ -53,10 +53,9 @@ struct f0980analysis {
5353
Configurable<float> cfgMaxDCAzToPVcut{"cfgMaxDCAzToPVcut", 2.0, "Maximum longitudinal DCA"};
5454
Configurable<float> cfgMinRap{"cfgMinRap", -0.5, "Minimum rapidity for pair"};
5555
Configurable<float> cfgMaxRap{"cfgMaxRap", 0.5, "Maximum rapidity for pair"};
56-
Configurable<bool> cfgFindRT{"cfgFindRT", false, "boolean for RT analysis"};
5756

5857
// Track selections
59-
Configurable<bool> cfgPrimaryTrack{"cfgPrimaryTrack", true, "Primary track selection"}; // kGoldenChi2 | kDCAxy | kDCAz
58+
Configurable<bool> cfgPrimaryTrack{"cfgPrimaryTrack", false, "Primary track selection"}; // kGoldenChi2 | kDCAxy | kDCAz
6059
Configurable<bool> cfgGlobalTrack{"cfgGlobalTrack", false, "Global track selection"}; // kGoldenChi2 | kDCAxy | kDCAz
6160
Configurable<bool> cfgGlobalWoDCATrack{"cfgGlobalWoDCATrack", true, "Global track selection without DCA"}; // kQualityTracks (kTrackType |
6261
// kTPCNCls | kTPCCrossedRows |
@@ -84,47 +83,85 @@ struct f0980analysis {
8483
ConfigurableAxis pTAxis{"pTAxis", {VARIABLE_WIDTH, 0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.8, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0, 6.0, 7.0, 8.0, 10.0, 13.0, 20.0}, "Transverse momentum Binning"};
8584
ConfigurableAxis centAxis{"centAxis", {VARIABLE_WIDTH, 0.0, 1.0, 5.0, 10.0, 15.0, 20.0, 25.0, 30.0, 35.0, 40.0, 45.0, 50.0, 55.0, 60.0, 65.0, 70.0, 75.0, 80.0, 95.0, 100.0, 105.0, 110.0}, "Centrality Binning"};
8685

86+
Configurable<bool> cfgFindRT{"cfgFindRT", false, "boolean for RT analysis"};
87+
Configurable<bool> cfgFindEP{"cfgFindEP", false, "boolean for Event plane analysis"};
88+
Configurable<bool> cfgQAEPLT{"cfgQAEPLT", false, "Fill QA histograms for Event Plane and Leading Track"};
89+
Configurable<bool> cfgQASelection{"cfgQASelection", true, "Fill QA histograms for Selection"};
90+
Configurable<bool> cfgQAMCTrue{"cfgQAMCTrue", false, "Fill QA histograms for MC True Selection"};
91+
8792
void init(o2::framework::InitContext&)
8893
{
8994
std::vector<double> lptBinning = {0, 5.0, 13.0, 20.0, 50.0, 1000.0};
9095

9196
AxisSpec rtAxis = {3, 0, 3};
9297
AxisSpec lptAxis = {lptBinning}; // Minimum leading hadron pT selection
93-
94-
AxisSpec pidqaAxis = {60, -6, 6};
95-
AxisSpec pTqaAxis = {200, 0, 20};
96-
AxisSpec phiqaAxis = {72, 0, o2::constants::math::TwoPI}; // Azimuthal angle axis
97-
9898
AxisSpec epAxis = {10, 0, o2::constants::math::PI}; // Event Plane
9999
AxisSpec epqaAxis = {200, -o2::constants::math::PI, o2::constants::math::PI};
100100
AxisSpec epResAxis = {200, -2, 2};
101101

102+
AxisSpec pidqaAxis = {60, -6, 6, "#sigma"};
103+
AxisSpec pTqaAxis = {200, 0, 20, "#it{p}_{T} (GeV/#it{c})"};
104+
AxisSpec phiqaAxis = {72, 0, o2::constants::math::TwoPI, "#Phi"}; // Azimuthal angle axis
105+
AxisSpec etaqaAxis = {150, -2, 2, "#eta"}; // Pseudorapidity axis
106+
AxisSpec rapqaAxis = {60, -1.5, 1.5, "#it{y}"}; // Rapidity axis
107+
AxisSpec dcaxyAxis = {200, -5.0, 5.0, "DCA_{xy} (cm)"}; // DCAxy axis
108+
AxisSpec dcazAxis = {200, -5.0, 5.0, "DCA_{z} (cm)"}; // DCAz axis
109+
110+
AxisSpec collCutAxis = {4, -0.5, 3.5, "Collision cut index for MC"};
111+
102112
if (cfgFindRT) {
103-
histos.add("hInvMass_f0980_US", "unlike invariant mass", {HistType::kTHnSparseF, {massAxis, pTAxis, centAxis, rtAxis, lptAxis}});
104-
histos.add("hInvMass_f0980_LSpp", "++ invariant mass", {HistType::kTHnSparseF, {massAxis, pTAxis, centAxis, rtAxis, lptAxis}});
105-
histos.add("hInvMass_f0980_LSmm", "-- invariant mass", {HistType::kTHnSparseF, {massAxis, pTAxis, centAxis, rtAxis, lptAxis}});
106-
} else {
113+
histos.add("hInvMass_f0980_US_RT", "unlike invariant mass", {HistType::kTHnSparseF, {massAxis, pTAxis, centAxis, rtAxis, lptAxis}});
114+
histos.add("hInvMass_f0980_LSpp_RT", "++ invariant mass", {HistType::kTHnSparseF, {massAxis, pTAxis, centAxis, rtAxis, lptAxis}});
115+
histos.add("hInvMass_f0980_LSmm_RT", "-- invariant mass", {HistType::kTHnSparseF, {massAxis, pTAxis, centAxis, rtAxis, lptAxis}});
116+
} else if (cfgFindEP){
107117
histos.add("hInvMass_f0980_US_EPA", "unlike invariant mass", {HistType::kTHnSparseF, {massAxis, pTAxis, centAxis, epAxis}});
108118
histos.add("hInvMass_f0980_LSpp_EPA", "++ invariant mass", {HistType::kTHnSparseF, {massAxis, pTAxis, centAxis, epAxis}});
109119
histos.add("hInvMass_f0980_LSmm_EPA", "-- invariant mass", {HistType::kTHnSparseF, {massAxis, pTAxis, centAxis, epAxis}});
120+
} else {
121+
histos.add("hInvMass_f0980_US", "unlike invariant mass", {HistType::kTHnSparseF, {massAxis, pTAxis, centAxis}});
122+
histos.add("hInvMass_f0980_LSpp", "++ invariant mass", {HistType::kTHnSparseF, {massAxis, pTAxis, centAxis}});
123+
histos.add("hInvMass_f0980_LSmm", "-- invariant mass", {HistType::kTHnSparseF, {massAxis, pTAxis, centAxis}});
110124
}
111125

112-
histos.add("QA/EPhist", "", {HistType::kTH2F, {centAxis, epqaAxis}});
113-
histos.add("QA/hEPResAB", "", {HistType::kTH2F, {centAxis, epResAxis}});
114-
histos.add("QA/hEPResBC", "", {HistType::kTH2F, {centAxis, epResAxis}});
115-
histos.add("QA/hEPResAC", "", {HistType::kTH2F, {centAxis, epResAxis}});
116-
117-
histos.add("QA/LTpt", "", {HistType::kTH3F, {pTqaAxis, centAxis, phiqaAxis}});
118-
119-
histos.add("QA/Nsigma_TPC", "", {HistType::kTH2F, {pTqaAxis, pidqaAxis}});
120-
histos.add("QA/Nsigma_TOF", "", {HistType::kTH2F, {pTqaAxis, pidqaAxis}});
121-
histos.add("QA/Nsigma_TPC_TOF", "", {HistType::kTH2F, {pidqaAxis, pidqaAxis}});
122-
123-
if (doprocessMCLight) {
124-
histos.add("MCL/hpT_f0980_GEN", "generated f0 signals", HistType::kTH1F, {pTqaAxis});
125-
histos.add("MCL/hpT_f0980_REC", "reconstructed f0 signals", HistType::kTH3F, {massAxis, pTqaAxis, centAxis});
126+
if (cfgQAEPLT) {
127+
// Event Plane QA
128+
histos.add("QA/EPhist", "", {HistType::kTH2F, {centAxis, epqaAxis}});
129+
histos.add("QA/hEPResAB", "", {HistType::kTH2F, {centAxis, epResAxis}});
130+
histos.add("QA/hEPResBC", "", {HistType::kTH2F, {centAxis, epResAxis}});
131+
histos.add("QA/hEPResAC", "", {HistType::kTH2F, {centAxis, epResAxis}});
132+
// Leading track pT QA
133+
histos.add("QA/LTpt", "", {HistType::kTH3F, {pTqaAxis, centAxis, phiqaAxis}});
134+
}
135+
if (cfgQASelection) {
136+
// General QA
137+
histos.add("QA/TrackPt", "", {HistType::kTH1F, {pTqaAxis}});
138+
histos.add("QA/TrackEta", "", {HistType::kTH1F, {etaqaAxis}});
139+
histos.add("QA/TrackPhi", "", {HistType::kTH1F, {phiqaAxis}});
140+
// Track selection QA
141+
histos.add("QA/trkDCAxy_BC", "DCA_{xy} for pion tracks (before cuts)", HistType::kTH2F, {pTqaAxis, dcaxyAxis});
142+
histos.add("QA/trkDCAz_BC", "DCA_{z} for pion tracks (before cuts)", HistType::kTH2F, {pTqaAxis, dcazAxis});
143+
histos.add("QA/trkDCAxy", "DCA_{xy} for pion tracks (after cuts)", HistType::kTH2F, {pTqaAxis, dcaxyAxis});
144+
histos.add("QA/trkDCAz", "DCA_{z} for pion tracks (after cuts)", HistType::kTH2F, {pTqaAxis, dcazAxis});
145+
// PID QA
146+
histos.add("QA/Nsigma_TPC_BC", "TPC n#sigma^{#pi} (before PID cuts); p_{T} (GeV/c); n#sigma_{TPC}^{#pi}", {HistType::kTH2F, {pTqaAxis, pidqaAxis}});
147+
histos.add("QA/Nsigma_TOF_BC", "TOF n#sigma^{#pi} (before PID cuts); p_{T} (GeV/c); n#sigma_{TOF}^{#pi}", {HistType::kTH2F, {pTqaAxis, pidqaAxis}});
148+
histos.add("QA/Nsigma_TPC_TOF_BC", "", {HistType::kTH2F, {pidqaAxis, pidqaAxis}});
149+
histos.add("QA/Nsigma_TPC", "TPC n#sigma^{#pi} (after PID cuts); p_{T} (GeV/c); n#sigma_{TPC}^{#pi}", {HistType::kTH2F, {pTqaAxis, pidqaAxis}});
150+
histos.add("QA/Nsigma_TOF", "TOF n#sigma^{#pi} (after PID cuts); p_{T} (GeV/c); n#sigma_{TOF}^{#pi}", {HistType::kTH2F, {pTqaAxis, pidqaAxis}});
151+
histos.add("QA/Nsigma_TPC_TOF", "", {HistType::kTH2F, {pidqaAxis, pidqaAxis}});
126152
}
127153

154+
if (doprocessMCRec) {
155+
histos.add("MCL/hpT_f0980_REC", "Reconstructed f0 signals", HistType::kTH3F, {massAxis, pTqaAxis, centAxis});
156+
}
157+
if (doprocessMCTrue) {
158+
// histos.add("MCL/hpT_f0980_GEN", "Generated f0 signals", HistType::kTH2F, {pTqaAxis, centAxis});
159+
histos.add("MCL/hpT_f0980_GEN", "Generated f0 signals; selIdx; p_{T} (GeV/c); Centrality (%)", HistType::kTH3F, {collCutAxis, pTqaAxis, centAxis});
160+
if (cfgQAMCTrue) {
161+
histos.add("QAMCTrue/f0_pt_y", "Generated f0 ; #it{p}_{T} (GeV/#it{c}) ; #it{y}", HistType::kTH2F, {pTqaAxis, rapqaAxis});
162+
histos.add("QAMCTrue/f0_pt_cent", "Generated f0 ; #it{p}_{T} (GeV/#it{c}); Centrality (%)", HistType::kTH2F, {pTqaAxis, centAxis});
163+
}
164+
}
128165
histos.print();
129166
}
130167

@@ -221,41 +258,71 @@ struct f0980analysis {
221258
lhphi = trk.phi();
222259
}
223260
}
261+
if (cfgQAEPLT)
262+
histos.fill(HIST("QA/LTpt"), lhpT, collision.cent(), lhphi);
263+
} else if (cfgFindEP) {
264+
if (cfgQAEPLT) {
265+
histos.fill(HIST("QA/EPhist"), collision.cent(), collision.evtPl());
266+
histos.fill(HIST("QA/hEPResAB"), collision.cent(), collision.evtPlResAB());
267+
histos.fill(HIST("QA/hEPResBC"), collision.cent(), collision.evtPlResBC());
268+
histos.fill(HIST("QA/hEPResAC"), collision.cent(), collision.evtPlResAC());
269+
}
224270
}
225-
histos.fill(HIST("QA/EPhist"), collision.cent(), collision.evtPl());
226-
histos.fill(HIST("QA/hEPResAB"), collision.cent(), collision.evtPlResAB());
227-
histos.fill(HIST("QA/hEPResBC"), collision.cent(), collision.evtPlResBC());
228-
histos.fill(HIST("QA/hEPResAC"), collision.cent(), collision.evtPlResAC());
229-
histos.fill(HIST("QA/LTpt"), lhpT, collision.cent(), lhphi);
230271

231272
ROOT::Math::LorentzVector<ROOT::Math::PxPyPzM4D<float>> pion1, pion2, reco;
232273
for (const auto& [trk1, trk2] : combinations(CombinationsStrictlyUpperIndexPolicy(dTracks, dTracks))) {
233-
274+
if (cfgQASelection) {
275+
histos.fill(HIST("QA/trkDCAxy_BC"), trk1.pt(), trk1.dcaXY());
276+
histos.fill(HIST("QA/trkDCAz_BC"), trk1.pt(), trk1.dcaZ());
277+
}
234278
if (!selTrack(trk1) || !selTrack(trk2))
235279
continue;
236-
histos.fill(HIST("QA/Nsigma_TPC"), trk1.pt(), trk1.tpcNSigmaPi());
237-
histos.fill(HIST("QA/Nsigma_TOF"), trk1.pt(), trk1.tofNSigmaPi());
238-
histos.fill(HIST("QA/Nsigma_TPC_TOF"), trk1.tpcNSigmaPi(), trk1.tofNSigmaPi());
239-
280+
if (cfgQASelection) {
281+
histos.fill(HIST("QA/trkDCAxy"), trk1.pt(), trk1.dcaXY());
282+
histos.fill(HIST("QA/trkDCAz"), trk1.pt(), trk1.dcaZ());
283+
284+
histos.fill(HIST("QA/Nsigma_TPC_BC"), trk1.pt(), trk1.tpcNSigmaPi());
285+
if (trk1.hasTOF()) {
286+
histos.fill(HIST("QA/Nsigma_TOF_BC"), trk1.pt(), trk1.tofNSigmaPi());
287+
histos.fill(HIST("QA/Nsigma_TPC_TOF_BC"), trk1.tpcNSigmaPi(), trk1.tofNSigmaPi());
288+
}
289+
}
240290
if (!selPion(trk1) || !selPion(trk2))
241291
continue;
292+
if (cfgQASelection) {
293+
histos.fill(HIST("QA/Nsigma_TPC"), trk1.pt(), trk1.tpcNSigmaPi());
294+
if (trk1.hasTOF()) {
295+
histos.fill(HIST("QA/Nsigma_TOF"), trk1.pt(), trk1.tofNSigmaPi());
296+
histos.fill(HIST("QA/Nsigma_TPC_TOF"), trk1.tpcNSigmaPi(), trk1.tofNSigmaPi());
297+
}
298+
299+
histos.fill(HIST("QA/TrackPt"), trk1.pt());
300+
histos.fill(HIST("QA/TrackEta"), trk1.eta());
301+
histos.fill(HIST("QA/TrackPhi"), trk1.phi());
302+
}
242303

243304
pion1 = ROOT::Math::PxPyPzMVector(trk1.px(), trk1.py(), trk1.pz(), massPi);
244305
pion2 = ROOT::Math::PxPyPzMVector(trk2.px(), trk2.py(), trk2.pz(), massPi);
245306
reco = pion1 + pion2;
246307
if (reco.Rapidity() > cfgMaxRap || reco.Rapidity() < cfgMinRap)
247308
continue;
248309

249-
relphi = TVector2::Phi_0_2pi(reco.Phi() - collision.evtPl());
250-
if (relphi > o2::constants::math::PI) {
251-
relphi -= o2::constants::math::PI;
310+
if (cfgFindEP) {
311+
relphi = TVector2::Phi_0_2pi(reco.Phi() - collision.evtPl());
312+
if (relphi > o2::constants::math::PI) {
313+
relphi -= o2::constants::math::PI;
314+
}
252315
}
253316

254317
if (trk1.sign() * trk2.sign() < 0) {
255318
if (cfgFindRT) {
256-
histos.fill(HIST("hInvMass_f0980_US"), reco.M(), reco.Pt(), collision.cent(), rtIndex(reco.Phi(), lhphi), lhpT);
319+
histos.fill(HIST("hInvMass_f0980_US_RT"), reco.M(), reco.Pt(), collision.cent(), rtIndex(reco.Phi(), lhphi), lhpT);
320+
} else if (cfgFindEP) {
321+
histos.fill(HIST("hInvMass_f0980_US_EPA"), reco.M(), reco.Pt(), collision.cent(), relphi);
322+
} else {
323+
histos.fill(HIST("hInvMass_f0980_US"), reco.M(), reco.Pt(), collision.cent());
257324
}
258-
histos.fill(HIST("hInvMass_f0980_US_EPA"), reco.M(), reco.Pt(), collision.cent(), relphi);
325+
259326
if constexpr (IsMC) {
260327
if (std::abs(trk1.pdgCode()) != kPiPlus || std::abs(trk2.pdgCode()) != kPiPlus)
261328
continue;
@@ -267,14 +334,20 @@ struct f0980analysis {
267334
}
268335
} else if (trk1.sign() > 0 && trk2.sign() > 0) {
269336
if (cfgFindRT) {
270-
histos.fill(HIST("hInvMass_f0980_LSpp"), reco.M(), reco.Pt(), collision.cent(), rtIndex(reco.Phi(), lhphi), lhpT);
337+
histos.fill(HIST("hInvMass_f0980_LSpp_RT"), reco.M(), reco.Pt(), collision.cent(), rtIndex(reco.Phi(), lhphi), lhpT);
338+
} else if (cfgFindEP) {
339+
histos.fill(HIST("hInvMass_f0980_LSpp_EPA"), reco.M(), reco.Pt(), collision.cent(), relphi);
340+
} else {
341+
histos.fill(HIST("hInvMass_f0980_LSpp"), reco.M(), reco.Pt(), collision.cent());
271342
}
272-
histos.fill(HIST("hInvMass_f0980_LSpp_EPA"), reco.M(), reco.Pt(), collision.cent(), relphi);
273343
} else if (trk1.sign() < 0 && trk2.sign() < 0) {
274344
if (cfgFindRT) {
275-
histos.fill(HIST("hInvMass_f0980_LSmm"), reco.M(), reco.Pt(), collision.cent(), rtIndex(reco.Phi(), lhphi), lhpT);
345+
histos.fill(HIST("hInvMass_f0980_LSmm_RT"), reco.M(), reco.Pt(), collision.cent(), rtIndex(reco.Phi(), lhphi), lhpT);
346+
} else if (cfgFindEP) {
347+
histos.fill(HIST("hInvMass_f0980_LSmm_EPA"), reco.M(), reco.Pt(), collision.cent(), relphi);
348+
} else {
349+
histos.fill(HIST("hInvMass_f0980_LSmm"), reco.M(), reco.Pt(), collision.cent());
276350
}
277-
histos.fill(HIST("hInvMass_f0980_LSmm_EPA"), reco.M(), reco.Pt(), collision.cent(), relphi);
278351
}
279352
}
280353
}
@@ -286,15 +359,15 @@ struct f0980analysis {
286359
}
287360
PROCESS_SWITCH(f0980analysis, processData, "Process Event for data", true);
288361

289-
void processMCLight(
290-
o2::soa::Join<o2::aod::ResoCollisions, o2::aod::ResoEvtPlCollisions>::iterator const& collision,
291-
o2::soa::Join<o2::aod::ResoTracks, o2::aod::ResoMCTracks> const& resotracks)
362+
void processMCRec(o2::soa::Join<o2::aod::ResoCollisions, o2::aod::ResoEvtPlCollisions>::iterator const& collision,
363+
o2::soa::Join<o2::aod::ResoTracks, o2::aod::ResoMCTracks> const& resotracks)
292364
{
293365
fillHistograms<true>(collision, resotracks);
294366
}
295-
PROCESS_SWITCH(f0980analysis, processMCLight, "Process Event for MC", false);
367+
PROCESS_SWITCH(f0980analysis, processMCRec, "Process Event for MC", false);
296368

297-
void processMCTrue(const o2::aod::ResoMCParents& resoParents)
369+
void processMCTrue(o2::soa::Join<o2::aod::ResoCollisions, o2::aod::ResoMCCollisions>::iterator const& resoCollision,
370+
o2::aod::ResoMCParents const& resoParents)
298371
{
299372
for (const auto& part : resoParents) { // loop over all pre-filtered MC particles
300373
if (std::abs(part.pdgCode()) != 9010221)
@@ -305,14 +378,29 @@ struct f0980analysis {
305378
continue;
306379
}
307380
bool pass = false;
308-
if ((std::abs(part.daughterPDG1()) == kPiPlus &&
309-
std::abs(part.daughterPDG2()) == kPiPlus)) {
381+
if ((std::abs(part.daughterPDG1()) == kPiPlus && std::abs(part.daughterPDG2()) == kPiPlus)) {
310382
pass = true;
311383
}
312384
if (!pass) // If we have both decay products
313385
continue;
314386

315-
histos.fill(HIST("MCL/hpT_f0980_GEN"), part.pt());
387+
// no event selection
388+
histos.fill(HIST("MCL/hpT_f0980_GEN"), 0, part.pt(), resoCollision.cent());
389+
// |zvtx|<10 cm
390+
if (resoCollision.isVtxIn10()) {
391+
histos.fill(HIST("MCL/hpT_f0980_GEN"), 1, part.pt(), resoCollision.cent());
392+
}
393+
// |zvtx|<10 cm & TVX trigger
394+
if (resoCollision.isVtxIn10() && resoCollision.isTriggerTVX()) {
395+
histos.fill(HIST("MCL/hpT_f0980_GEN"), 2, part.pt(), resoCollision.cent());
396+
}
397+
if (resoCollision.isInAfterAllCuts()) {
398+
histos.fill(HIST("MCL/hpT_f0980_GEN"), 3, part.pt(), resoCollision.cent());
399+
}
400+
if (cfgQAMCTrue) {
401+
histos.fill(HIST("QAMCTrue/f0_pt_y"), part.pt(), part.y());
402+
histos.fill(HIST("QAMCTrue/f0_pt_cent"), part.pt(), resoCollision.cent());
403+
}
316404
}
317405
};
318406
PROCESS_SWITCH(f0980analysis, processMCTrue, "Process Event for MC", false);

0 commit comments

Comments
 (0)