Add centrality information to MC collisions

1, Added centrality for MC collisions and updated relevant event split and event loose processing functions.

Author: Tao-Fang

PWGHF/D2H/Tasks/taskXic0ToXiPi.cxx

@@ -19,6 +19,7 @@
 #include "PWGHF/Core/HfHelper.h"
 #include "PWGHF/DataModel/CandidateReconstructionTables.h"
 #include "PWGHF/DataModel/CandidateSelectionTables.h"
+#include "PWGLF/DataModel/mcCentrality.h"
 
 #include "Common/Core/RecoDecay.h"
 #include "Common/DataModel/Centrality.h"
@@ -81,6 +82,8 @@ struct HfTaskXic0ToXiPi {
   using CollisionsWithFT0M = soa::Join<aod::Collisions, aod::EvSels, aod::CentFT0Ms>;
   using CollisionsWithMcLabels = soa::Join<aod::Collisions, aod::McCollisionLabels, aod::EvSels>;
 
+  using McCollisionsCentFT0Ms = soa::Join<aod::McCollisions, aod::McCentFT0Ms>;
+
   Filter filterSelectXic0Candidates = aod::hf_sel_toxipi::resultSelections == true;
   Filter filterXicMatchedRec = nabs(aod::hf_cand_xic0_omegac0::flagMcMatchRec) == static_cast<int8_t>(BIT(aod::hf_cand_xic0_omegac0::DecayType::XiczeroToXiPi));
   Filter filterXicMatchedGen = nabs(aod::hf_cand_xic0_omegac0::flagMcMatchGen) == static_cast<int8_t>(BIT(aod::hf_cand_xic0_omegac0::DecayType::XiczeroToXiPi));
@@ -124,18 +127,25 @@ struct HfTaskXic0ToXiPi {
     const AxisSpec thnAxisGenPtD{thnConfigAxisGenPtD, "#it{p}_{T} (GeV/#it{c})"};
     const AxisSpec thnAxisGenPtB{thnConfigAxisGenPtB, "#it{p}_{T}^{B} (GeV/#it{c})"};
     const AxisSpec thnAxisNumPvContr{thnConfigAxisNumPvContr, "Number of PV contributors"};
+    const AxisSpec thnAxisCent{thnConfigAxisCent, "Centrality percentile"};
 
     if (doprocessMcWithKFParticle || doprocessMcWithKFParticleMl || doprocessMcWithDCAFitter || doprocessMcWithDCAFitterMl) {
-      std::vector<AxisSpec> const axesAcc = {thnAxisGenPtD, thnAxisGenPtB, thnAxisY, thnAxisOrigin, thnAxisNumPvContr};
+      std::vector<AxisSpec> const axesAcc = {thnAxisGenPtD, thnAxisGenPtB, thnAxisY, thnAxisOrigin, thnAxisCent, thnAxisNumPvContr};
       registry.add("hSparseAcc", "Thn for generated Xic0 from charm and beauty", HistType::kTHnSparseD, axesAcc);
       registry.get<THnSparse>(HIST("hSparseAcc"))->Sumw2();
+
+      registry.add("hSparseAccWithRecoColl", "Gen. Xic0 from charm and beauty (associated to a reco collision)", HistType::kTHnSparseD, axesAcc);
+      registry.get<THnSparse>(HIST("hSparseAccWithRecoColl"))->Sumw2();
+
+      registry.add("hNumRecoCollPerMcColl", "Number of reco collisions associated to a mc collision;Num. reco. coll. per Mc coll.;", {HistType::kTH1D, {{10, -0.5, 9.5}}});
     }
 
     std::vector<AxisSpec> axes = {thnAxisMass, thnAxisPt, thnAxisY};
-    if (doprocessMcWithKFParticle || doprocessMcWithKFParticleMl) {
+    if (doprocessMcWithKFParticle || doprocessMcWithKFParticleMl || doprocessMcWithDCAFitter || doprocessMcWithDCAFitterMl) {
       axes.push_back(thnAxisPtB);
       axes.push_back(thnAxisOrigin);
       axes.push_back(thnAxisMatchFlag);
+      axes.push_back(thnAxisCent);
       axes.push_back(thnAxisNumPvContr);
     }
     if (applyMl) {
@@ -203,12 +213,12 @@ struct HfTaskXic0ToXiPi {
     }
   }
 
-  template <bool UseKfParticle, bool ApplyMl, typename CandType, typename CollType>
+  template <bool UseKfParticle, bool ApplyMl, typename CandType, typename CollType, typename McCollisionWithCents>
   void processMc(const CandType& candidates,
                  Xic0Gen const& mcParticles,
                  TracksMc const&,
                  CollType const& collisions,
-                 aod::McCollisions const&)
+                 McCollisionWithCents const&)
   {
     // MC rec.
     for (const auto& candidate : candidates) {
@@ -225,6 +235,8 @@ struct HfTaskXic0ToXiPi {
         continue;
       }
 
+      auto collision = candidate.template collision_as<CollisionsWithMcLabels>();
+      float const mcCent = o2::hf_centrality::getCentralityColl(collision.template mcCollision_as<McCollisionWithCents>());
       auto numPvContributors = candidate.template collision_as<CollType>().numContrib();
       double const ptXic = RecoDecay::pt(candidate.pxCharmBaryon(), candidate.pyCharmBaryon());
       if constexpr (ApplyMl) {
@@ -236,6 +248,7 @@ struct HfTaskXic0ToXiPi {
                       candidate.ptBhadMotherPart(),
                       candidate.originMcRec(),
                       candidate.flagMcMatchRec(),
+                      mcCent,
                       numPvContributors);
       } else {
         registry.fill(HIST("hMassVsPtVsPtBVsYVsOriginVsXic0Type"),
@@ -245,6 +258,7 @@ struct HfTaskXic0ToXiPi {
                       candidate.ptBhadMotherPart(),
                       candidate.originMcRec(),
                       candidate.flagMcMatchRec(),
+                      mcCent,
                       numPvContributors);
       }
     }
@@ -258,6 +272,7 @@ struct HfTaskXic0ToXiPi {
       auto ptGen = particle.pt();
       auto yGen = particle.rapidityCharmBaryonGen();
 
+      float const mcCent = o2::hf_centrality::getCentralityColl(particle.template mcCollision_as<McCollisionWithCents>());
       unsigned maxNumContrib = 0;
       const auto& recoCollsPerMcColl = collisions.sliceBy(colPerMcCollision, particle.mcCollision().globalIndex());
       for (const auto& recCol : recoCollsPerMcColl) {
@@ -270,6 +285,7 @@ struct HfTaskXic0ToXiPi {
                       -1.,
                       yGen,
                       RecoDecay::OriginType::Prompt,
+                      mcCent,
                       maxNumContrib);
       } else {
         float const ptGenB = mcParticles.rawIteratorAt(particle.idxBhadMotherPart()).pt();
@@ -278,8 +294,33 @@ struct HfTaskXic0ToXiPi {
                       ptGenB,
                       yGen,
                       RecoDecay::OriginType::NonPrompt,
+                      mcCent,
                       maxNumContrib);
       }
+
+      registry.fill(HIST("hNumRecoCollPerMcColl"), recoCollsPerMcColl.size());
+
+      // fill sparse only for gen particles associated to a reconstructed collision
+      if (recoCollsPerMcColl.size() >= 1) {
+        if (particle.originMcGen() == RecoDecay::OriginType::Prompt) {
+          registry.fill(HIST("hSparseAccWithRecoColl"),
+                        ptGen,
+                        -1.,
+                        yGen,
+                        RecoDecay::OriginType::Prompt,
+                        mcCent,
+                        maxNumContrib);
+        } else {
+          float const ptGenB = mcParticles.rawIteratorAt(particle.idxBhadMotherPart()).pt();
+          registry.fill(HIST("hSparseAccWithRecoColl"),
+                        ptGen,
+                        ptGenB,
+                        yGen,
+                        RecoDecay::OriginType::NonPrompt,
+                        mcCent,
+                        maxNumContrib);
+        }
+      }
     }
   }
 
@@ -443,7 +484,7 @@ struct HfTaskXic0ToXiPi {
                               Xic0Gen const& mcParticles,
                               TracksMc const& tracks,
                               CollisionsWithMcLabels const& collisions,
-                              aod::McCollisions const& mcCollisions)
+                              McCollisionsCentFT0Ms const& mcCollisions)
   {
     processMc<false, false>(xic0CandidatesMc, mcParticles, tracks, collisions, mcCollisions);
   }
@@ -453,7 +494,7 @@ struct HfTaskXic0ToXiPi {
                                Xic0Gen const& mcParticles,
                                TracksMc const& tracks,
                                CollisionsWithMcLabels const& collisions,
-                               aod::McCollisions const& mcCollisions)
+                               McCollisionsCentFT0Ms const& mcCollisions)
   {
     processMc<true, false>(xic0CandidatesMcKf, mcParticles, tracks, collisions, mcCollisions);
   }
@@ -463,7 +504,7 @@ struct HfTaskXic0ToXiPi {
                                 Xic0Gen const& mcParticles,
                                 TracksMc const& tracks,
                                 CollisionsWithMcLabels const& collisions,
-                                aod::McCollisions const& mcCollisions)
+                                McCollisionsCentFT0Ms const& mcCollisions)
   {
     processMc<false, true>(xic0CandidatesMlMc, mcParticles, tracks, collisions, mcCollisions);
   }
@@ -473,7 +514,7 @@ struct HfTaskXic0ToXiPi {
                                  Xic0Gen const& mcParticles,
                                  TracksMc const& tracks,
                                  CollisionsWithMcLabels const& collisions,
-                                 aod::McCollisions const& mcCollisions)
+                                 McCollisionsCentFT0Ms const& mcCollisions)
   {
     processMc<true, true>(xic0CandidatesMlMcKf, mcParticles, tracks, collisions, mcCollisions);
   }