PWGHF: adding particle origin MC information to LcToK0sP tree creator (#2657)

Author: DanielSamitz

PWGHF/DataModel/CandidateReconstructionTables.h

@@ -535,6 +535,8 @@ DECLARE_SOA_DYNAMIC_COLUMN(PtV0Neg, ptV0Neg, //!
                            [](float px, float py) { return RecoDecay::pt(px, py); });
 DECLARE_SOA_COLUMN(FlagMcMatchRec, flagMcMatchRec, int8_t); //! reconstruction level
 DECLARE_SOA_COLUMN(FlagMcMatchGen, flagMcMatchGen, int8_t); //! generator level
+DECLARE_SOA_COLUMN(OriginMcRec, originMcRec, int8_t);       //! particle origin, reconstruction level
+DECLARE_SOA_COLUMN(OriginMcGen, originMcGen, int8_t);       //! particle origin, generator level
 DECLARE_SOA_COLUMN(V0X, v0x, float);
 DECLARE_SOA_COLUMN(V0Y, v0y, float);
 DECLARE_SOA_COLUMN(V0Z, v0z, float);
@@ -614,11 +616,13 @@ using HfCandCascade = HfCandCascExt;
 
 // table with results of reconstruction level MC matching for Cascade
 DECLARE_SOA_TABLE(HfCandCascadeMcRec, "AOD", "HFCANDCASCMCREC", //!
-                  hf_cand_casc::FlagMcMatchRec);
+                  hf_cand_casc::FlagMcMatchRec,
+                  hf_cand_casc::OriginMcRec);
 
 // table with results of generator level MC matching
 DECLARE_SOA_TABLE(HfCandCascadeMcGen, "AOD", "HFCANDCASCMCGEN", //!
-                  hf_cand_casc::FlagMcMatchGen);
+                  hf_cand_casc::FlagMcMatchGen,
+                  hf_cand_casc::OriginMcGen);
 
 // specific BPlus candidate properties
 namespace hf_cand_bplus

PWGHF/TableProducer/candidateCreatorCascade.cxx

@@ -264,13 +264,17 @@ struct HfCandidateCreatorCascadeMc {
                  aod::McParticles const& particlesMC)
   {
     int8_t sign = 0;
+    int8_t origin = 0;
+    int indexRec = -1;
     std::vector<int> arrDaughLcIndex;
     std::array<int, 3> arrDaughLcPDG;
     std::array<int, 3> arrDaughLcPDGRef = {2212, 211, -211};
 
     // Match reconstructed candidates.
     rowCandidateCasc->bindExternalIndices(&tracks);
-    for (auto& candidate : *rowCandidateCasc) {
+    for (const auto& candidate : *rowCandidateCasc) {
+
+      origin = 0;
 
       const auto& bach = candidate.prong0_as<aod::BigTracksMC>();
       const auto& trackV0DaughPos = candidate.posTrack_as<aod::BigTracksMC>();
@@ -299,16 +303,23 @@ struct HfCandidateCreatorCascadeMc {
         // then we check the Lc
         MY_DEBUG_MSG(sign, LOG(info) << "K0S was correct! now we check the Lc");
         MY_DEBUG_MSG(sign, LOG(info) << "index proton = " << indexBach);
-        RecoDecay::getMatchedMCRec(particlesMC, arrayDaughtersLc, pdg::Code::kLambdaCPlus, array{+kProton, +kPiPlus, -kPiPlus}, true, &sign, 3); // 3-levels Lc --> p + K0 --> p + K0s --> p + pi+ pi-
+        indexRec = RecoDecay::getMatchedMCRec(particlesMC, arrayDaughtersLc, pdg::Code::kLambdaCPlus, array{+kProton, +kPiPlus, -kPiPlus}, true, &sign, 3); // 3-levels Lc --> p + K0 --> p + K0s --> p + pi+ pi-
         MY_DEBUG_MSG(sign, LOG(info) << "Lc found with sign " << sign; printf("\n"));
       }
 
-      rowMcMatchRec(sign);
+      // Check whether the particle is non-prompt (from a b quark).
+      if (sign != 0) {
+        auto particle = particlesMC.rawIteratorAt(indexRec);
+        origin = RecoDecay::getCharmHadronOrigin(particlesMC, particle);
+      }
+
+      rowMcMatchRec(sign, origin);
     }
     //}
 
     // Match generated particles.
-    for (auto& particle : particlesMC) {
+    for (const auto& particle : particlesMC) {
+      origin = 0;
       // checking if I have a Lc --> K0S + p
       RecoDecay::isMatchedMCGen(particlesMC, particle, pdg::Code::kLambdaCPlus, array{+kProton, +kK0Short}, true, &sign, 2);
       if (sign != 0) {
@@ -329,7 +340,11 @@ struct HfCandidateCreatorCascadeMc {
           MY_DEBUG_MSG(sign == 0, LOG(info) << "Pity, the three final daughters are not p, pi+, pi-, but " << arrDaughLcPDG[0] << ", " << arrDaughLcPDG[1] << ", " << arrDaughLcPDG[2]);
         }
       }
-      rowMcMatchGen(sign);
+      // Check whether the particle is non-prompt (from a b quark).
+      if (sign != 0) {
+        origin = RecoDecay::getCharmHadronOrigin(particlesMC, particle);
+      }
+      rowMcMatchGen(sign, origin);
     }
   }
 

PWGHF/TableProducer/treeCreatorLcToK0sP.cxx

@@ -66,6 +66,8 @@ DECLARE_SOA_COLUMN(V0MAntiLambda, v0MAntiLambda, float);
 DECLARE_SOA_COLUMN(V0MK0Short, v0MK0Short, float);
 DECLARE_SOA_COLUMN(V0MGamma, v0MGamma, float);
 DECLARE_SOA_COLUMN(FlagMc, flagMc, int8_t);
+DECLARE_SOA_COLUMN(OriginMcRec, originMcRec, int8_t);
+DECLARE_SOA_COLUMN(OriginMcGen, originMcGen, int8_t);
 // Events
 DECLARE_SOA_COLUMN(IsEventReject, isEventReject, int);
 DECLARE_SOA_COLUMN(RunNumber, runNumber, int);
@@ -136,7 +138,8 @@ DECLARE_SOA_TABLE(HfCandCascFull, "AOD", "HFCANDCASCFull",
                   full::Phi,
                   full::Y,
                   full::E,
-                  full::FlagMc);
+                  full::FlagMc,
+                  full::OriginMcRec);
 
 DECLARE_SOA_TABLE(HfCandCascFullEvents, "AOD", "HFCANDCASCFullE",
                   collision::BCId,
@@ -151,7 +154,8 @@ DECLARE_SOA_TABLE(HfCandCascFullParticles, "AOD", "HFCANDCASCFullP",
                   full::Eta,
                   full::Phi,
                   full::Y,
-                  full::FlagMc);
+                  full::FlagMc,
+                  full::OriginMcGen);
 
 } // namespace o2::aod
 
@@ -168,7 +172,7 @@ struct HfTreeCreatorLcToK0sP {
   }
 
   template <typename T, typename U>
-  void fillCandidate(const T& candidate, const U& bach, int8_t flagMc)
+  void fillCandidate(const T& candidate, const U& bach, int8_t flagMc, int8_t originMcRec)
   {
     rowCandidateFull(
       bach.collision().bcId(),
@@ -235,7 +239,8 @@ struct HfTreeCreatorLcToK0sP {
       candidate.phi(),
       o2::aod::hf_cand_3prong::yLc(candidate),
       o2::aod::hf_cand_3prong::eLc(candidate),
-      flagMc);
+      flagMc,
+      originMcRec);
   }
 
   template <typename T>
@@ -268,7 +273,7 @@ struct HfTreeCreatorLcToK0sP {
       auto bach = candidate.prong0_as<aod::BigTracksPID>(); // bachelor
       double pseudoRndm = bach.pt() * 1000. - (int16_t)(bach.pt() * 1000);
       if (candidate.isSelLcToK0sP() >= 1 && pseudoRndm < downSampleBkgFactor) {
-        fillCandidate(candidate, bach, candidate.flagMcMatchRec());
+        fillCandidate(candidate, bach, candidate.flagMcMatchRec(), candidate.originMcRec());
       }
     }
 
@@ -281,8 +286,10 @@ struct HfTreeCreatorLcToK0sP {
           particle.pt(),
           particle.eta(),
           particle.phi(),
-          RecoDecay::y(array{particle.px(), particle.py(), particle.pz()}, RecoDecay::getMassPDG(particle.pdgCode())),
-          particle.flagMcMatchGen());
+          RecoDecay::y(array{particle.px(), particle.py(), particle.pz()},
+                       RecoDecay::getMassPDG(particle.pdgCode())),
+          particle.flagMcMatchGen(),
+          particle.originMcGen());
       }
     }
   }
@@ -305,7 +312,7 @@ struct HfTreeCreatorLcToK0sP {
       auto bach = candidate.prong0_as<aod::BigTracksPID>(); // bachelor
       double pseudoRndm = bach.pt() * 1000. - (int16_t)(bach.pt() * 1000);
       if (candidate.isSelLcToK0sP() >= 1 && pseudoRndm < downSampleBkgFactor) {
-        fillCandidate(candidate, bach, 0);
+        fillCandidate(candidate, bach, 0, 0);
       }
     }
   }