Skip to content

Commit 902b6f9

Browse files
[PWGLF] modify MC truth tagging part (#9465)
Co-authored-by: ALICE Action Bot <alibuild@cern.ch>
1 parent 9498b2d commit 902b6f9

File tree

1 file changed

+32
-91
lines changed

1 file changed

+32
-91
lines changed

PWGLF/Tasks/Resonances/rho770analysis.cxx

Lines changed: 32 additions & 91 deletions
Original file line numberDiff line numberDiff line change
@@ -12,7 +12,7 @@
1212
/// \file rho770analysis.cxx
1313
/// \brief rho(770)0 analysis in pp 13 & 13.6 TeV
1414
/// \author Hyunji Lim (hyunji.lim@cern.ch)
15-
/// \since 12/04/2024
15+
/// \since 01/23/2025
1616

1717
#include <Framework/Configurable.h>
1818
#include <TLorentzVector.h>
@@ -38,6 +38,8 @@ struct rho770analysis {
3838
SliceCache cache;
3939
HistogramRegistry histos{"histos", {}, OutputObjHandlingPolicy::AnalysisObject};
4040

41+
using ResoMCCols = soa::Join<aod::ResoCollisions, aod::ResoMCCollisions>;
42+
4143
Configurable<float> cfgMinPt{"cfgMinPt", 0.15, "Minimum transverse momentum for charged track"};
4244
Configurable<float> cfgMaxEta{"cfgMaxEta", 0.8, "Maximum pseudorapidiy for charged track"};
4345
Configurable<float> cfgMaxDCArToPVcut{"cfgMaxDCArToPVcut", 0.15, "Maximum transverse DCA"};
@@ -87,6 +89,7 @@ struct rho770analysis {
8789
{
8890
AxisSpec pidqaAxis = {120, -6, 6};
8991
AxisSpec pTqaAxis = {200, 0, 20};
92+
AxisSpec mcLabelAxis = {5, -0.5, 4.5, "MC Label"};
9093

9194
histos.add("hInvMass_rho770_US", "unlike invariant mass", {HistType::kTHnSparseF, {massAxis, ptAxis, centAxis}});
9295
histos.add("hInvMass_rho770_LSpp", "++ invariant mass", {HistType::kTHnSparseF, {massAxis, ptAxis, centAxis}});
@@ -105,7 +108,7 @@ struct rho770analysis {
105108
histos.add("QA/TPC_TOF", "", {HistType::kTH2F, {pidqaAxis, pidqaAxis}});
106109

107110
if (doprocessMCLight) {
108-
histos.add("MCL/hpT_rho770_GEN", "generated rho770 signals", HistType::kTH1F, {ptAxis});
111+
histos.add("MCL/hpT_rho770_GEN", "generated rho770 signals", HistType::kTHnSparseF, {mcLabelAxis, massAxis, ptAxis, centAxis});
109112
histos.add("MCL/hpT_rho770_REC", "reconstructed rho770 signals", HistType::kTHnSparseF, {massAxis, ptAxis, centAxis});
110113
histos.add("MCL/hpT_omega_REC", "reconstructed omega signals", HistType::kTHnSparseF, {massAxis, ptAxis, centAxis});
111114
histos.add("MCL/hpT_K0s_REC", "reconstructed K0s signals", HistType::kTHnSparseF, {massAxis, ptAxis, centAxis});
@@ -204,6 +207,7 @@ struct rho770analysis {
204207
{
205208
TLorentzVector part1, part2, reco;
206209
for (const auto& [trk1, trk2] : combinations(CombinationsUpperIndexPolicy(dTracks, dTracks))) {
210+
207211
if (trk1.index() == trk2.index()) {
208212
if (!selTrack(trk1))
209213
continue;
@@ -242,84 +246,8 @@ struct rho770analysis {
242246
histos.fill(HIST("MCL/hpT_K0s_pipi_REC"), reco.M(), reco.Pt(), collision.cent());
243247
}
244248
} else if ((std::abs(trk1.pdgCode()) == 211 && std::abs(trk2.pdgCode()) == 321) || (std::abs(trk1.pdgCode()) == 321 && std::abs(trk2.pdgCode()) == 211)) {
245-
if (std::abs(trk1.motherPDG()) == 313)
249+
if (std::abs(trk1.motherPDG()) == 313) {
246250
histos.fill(HIST("MCL/hpT_Kstar_REC"), reco.M(), reco.Pt(), collision.cent());
247-
histos.fill(HIST("QA/Nsigma_TPC"), trk1.pt(), trk1.tpcNSigmaPi());
248-
histos.fill(HIST("QA/Nsigma_TOF"), trk1.pt(), trk1.tofNSigmaPi());
249-
histos.fill(HIST("QA/TPC_TOF"), trk1.tpcNSigmaPi(), trk1.tofNSigmaPi());
250-
continue;
251-
}
252-
253-
if (!selTrack(trk1) || !selTrack(trk2))
254-
continue;
255-
256-
if (selPion(trk1) && selPion(trk2)) {
257-
part1.SetXYZM(trk1.px(), trk1.py(), trk1.pz(), massPi);
258-
part2.SetXYZM(trk2.px(), trk2.py(), trk2.pz(), massPi);
259-
reco = part1 + part2;
260-
261-
if (reco.Rapidity() > cfgMaxRap || reco.Rapidity() < cfgMinRap)
262-
continue;
263-
264-
if (trk1.sign() * trk2.sign() < 0) {
265-
histos.fill(HIST("hInvMass_rho770_US"), reco.M(), reco.Pt(), collision.cent());
266-
histos.fill(HIST("hInvMass_K0s_US"), reco.M(), reco.Pt(), collision.cent());
267-
268-
if constexpr (IsMC) {
269-
if (trk1.motherId() != trk2.motherId())
270-
continue;
271-
if (std::abs(trk1.pdgCode()) == 211 && std::abs(trk2.pdgCode()) == 211) {
272-
if (std::abs(trk1.motherPDG()) == 113) {
273-
histos.fill(HIST("MCL/hpT_rho770_REC"), reco.M(), reco.Pt(), collision.cent());
274-
} else if (std::abs(trk1.motherPDG()) == 223) {
275-
histos.fill(HIST("MCL/hpT_omega_REC"), reco.M(), reco.Pt(), collision.cent());
276-
} else if (std::abs(trk1.motherPDG()) == 310) {
277-
histos.fill(HIST("MCL/hpT_K0s_REC"), reco.M(), reco.Pt(), collision.cent());
278-
histos.fill(HIST("MCL/hpT_K0s_pipi_REC"), reco.M(), reco.Pt(), collision.cent());
279-
}
280-
} else if ((std::abs(trk1.pdgCode()) == 211 && std::abs(trk2.pdgCode()) == 321) || (std::abs(trk1.pdgCode()) == 321 && std::abs(trk2.pdgCode()) == 211)) {
281-
if (std::abs(trk1.motherPDG()) == 313)
282-
histos.fill(HIST("MCL/hpT_Kstar_REC"), reco.M(), reco.Pt(), collision.cent());
283-
}
284-
}
285-
} else if (trk1.sign() > 0 && trk2.sign() > 0) {
286-
histos.fill(HIST("hInvMass_rho770_LSpp"), reco.M(), reco.Pt(), collision.cent());
287-
histos.fill(HIST("hInvMass_K0s_LSpp"), reco.M(), reco.Pt(), collision.cent());
288-
} else if (trk1.sign() < 0 && trk2.sign() < 0) {
289-
histos.fill(HIST("hInvMass_rho770_LSmm"), reco.M(), reco.Pt(), collision.cent());
290-
histos.fill(HIST("hInvMass_K0s_LSmm"), reco.M(), reco.Pt(), collision.cent());
291-
}
292-
}
293-
294-
if ((selPion(trk1) && selKaon(trk2)) || (selKaon(trk1) && selPion(trk2))) {
295-
if (selPion(trk1)) {
296-
part1.SetXYZM(trk1.px(), trk1.py(), trk1.pz(), massPi);
297-
part2.SetXYZM(trk2.px(), trk2.py(), trk2.pz(), massKa);
298-
} else if (selPion(trk2)) {
299-
part1.SetXYZM(trk1.px(), trk1.py(), trk1.pz(), massKa);
300-
part2.SetXYZM(trk2.px(), trk2.py(), trk2.pz(), massPi);
301-
}
302-
reco = part1 + part2;
303-
304-
if (reco.Rapidity() > cfgMaxRap || reco.Rapidity() < cfgMinRap)
305-
continue;
306-
307-
if (trk1.sign() * trk2.sign() < 0) {
308-
histos.fill(HIST("hInvMass_Kstar_US"), reco.M(), reco.Pt(), collision.cent());
309-
310-
if constexpr (IsMC) {
311-
if (trk1.motherId() != trk2.motherId())
312-
continue;
313-
314-
if ((std::abs(trk1.pdgCode()) == 211 && std::abs(trk2.pdgCode()) == 321) || (std::abs(trk1.pdgCode()) == 321 && std::abs(trk2.pdgCode()) == 211)) {
315-
if (std::abs(trk1.motherPDG()) == 313)
316-
histos.fill(HIST("MCL/hpT_Kstar_Kpi_REC"), reco.M(), reco.Pt(), collision.cent());
317-
}
318-
}
319-
} else if (trk1.sign() > 0 && trk2.sign() > 0) {
320-
histos.fill(HIST("hInvMass_Kstar_LSpp"), reco.M(), reco.Pt(), collision.cent());
321-
} else if (trk1.sign() < 0 && trk2.sign() < 0) {
322-
histos.fill(HIST("hInvMass_Kstar_LSmm"), reco.M(), reco.Pt(), collision.cent());
323251
}
324252
}
325253
}
@@ -333,10 +261,10 @@ struct rho770analysis {
333261
}
334262

335263
if ((selPion(trk1) && selKaon(trk2)) || (selKaon(trk1) && selPion(trk2))) {
336-
if (selPion(trk1)) {
264+
if (selPion(trk1) && selKaon(trk2)) {
337265
part1.SetXYZM(trk1.px(), trk1.py(), trk1.pz(), massPi);
338266
part2.SetXYZM(trk2.px(), trk2.py(), trk2.pz(), massKa);
339-
} else if (selPion(trk2)) {
267+
} else if (selKaon(trk1) && selPion(trk2)) {
340268
part1.SetXYZM(trk1.px(), trk1.py(), trk1.pz(), massKa);
341269
part2.SetXYZM(trk2.px(), trk2.py(), trk2.pz(), massPi);
342270
}
@@ -351,10 +279,10 @@ struct rho770analysis {
351279
if constexpr (IsMC) {
352280
if (trk1.motherId() != trk2.motherId())
353281
continue;
354-
355282
if ((std::abs(trk1.pdgCode()) == 211 && std::abs(trk2.pdgCode()) == 321) || (std::abs(trk1.pdgCode()) == 321 && std::abs(trk2.pdgCode()) == 211)) {
356-
if (std::abs(trk1.motherPDG()) == 313)
283+
if (std::abs(trk1.motherPDG()) == 313) {
357284
histos.fill(HIST("MCL/hpT_Kstar_Kpi_REC"), reco.M(), reco.Pt(), collision.cent());
285+
}
358286
}
359287
}
360288
} else if (trk1.sign() > 0 && trk2.sign() > 0) {
@@ -378,29 +306,42 @@ struct rho770analysis {
378306
}
379307
PROCESS_SWITCH(rho770analysis, processMCLight, "Process Event for MC", false);
380308

381-
void processMCTrue(aod::ResoMCParents const& resoParents)
309+
void processMCTrue(ResoMCCols::iterator const& collision, aod::ResoMCParents const& resoParents)
382310
{
311+
auto multiplicity = collision.cent();
312+
383313
for (const auto& part : resoParents) { // loop over all pre-filtered MC particles
384314
if (std::abs(part.pdgCode()) != 113)
385315
continue;
386316
if (!part.producedByGenerator())
387317
continue;
388318
if (part.y() < cfgMinRap || part.y() > cfgMaxRap)
389319
continue;
390-
bool pass = false;
391-
if ((std::abs(part.daughterPDG1()) == 211 && std::abs(part.daughterPDG2()) == 211))
392-
pass = true;
393-
if (!pass) // If we have both decay products
320+
if (!(std::abs(part.daughterPDG1()) == 211 && std::abs(part.daughterPDG2()) == 211))
394321
continue;
395322

396-
histos.fill(HIST("MCL/hpT_rho770_GEN"), part.pt());
323+
TLorentzVector truthpar;
324+
truthpar.SetPxPyPzE(part.px(), part.py(), part.pz(), part.e());
325+
auto mass = truthpar.M();
326+
327+
if (collision.isVtxIn10()) {
328+
histos.fill(HIST("MCL/hpT_rho770_GEN"), 0, mass, part.pt(), multiplicity);
329+
}
330+
if (collision.isVtxIn10() && collision.isInSel8()) {
331+
histos.fill(HIST("MCL/hpT_rho770_GEN"), 1, mass, part.pt(), multiplicity);
332+
}
333+
if (collision.isVtxIn10() && collision.isTriggerTVX()) {
334+
histos.fill(HIST("MCL/hpT_rho770_GEN"), 2, mass, part.pt(), multiplicity);
335+
}
336+
if (collision.isInAfterAllCuts()) {
337+
histos.fill(HIST("MCL/hpT_rho770_GEN"), 3, mass, part.pt(), multiplicity);
338+
}
397339
}
398340
};
399341
PROCESS_SWITCH(rho770analysis, processMCTrue, "Process Event for MC", false);
400342
};
401343

402344
WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
403345
{
404-
return WorkflowSpec{
405-
adaptAnalysisTask<rho770analysis>(cfgc)};
346+
return WorkflowSpec{adaptAnalysisTask<rho770analysis>(cfgc)};
406347
}

0 commit comments

Comments
 (0)