From 7a7f6a932916b1787b3ea6a5b8e46adf32e88120 Mon Sep 17 00:00:00 2001 From: mhartung71 <50153519+mhartung71@users.noreply.github.com> Date: Thu, 29 May 2025 12:26:31 +0200 Subject: [PATCH 1/5] run number added --- PWGLF/DataModel/LFHypernucleiKfTables.h | 14 ++++++++------ 1 file changed, 8 insertions(+), 6 deletions(-) diff --git a/PWGLF/DataModel/LFHypernucleiKfTables.h b/PWGLF/DataModel/LFHypernucleiKfTables.h index 98f8076b2f6..0e2424f3bb5 100644 --- a/PWGLF/DataModel/LFHypernucleiKfTables.h +++ b/PWGLF/DataModel/LFHypernucleiKfTables.h @@ -45,6 +45,7 @@ DECLARE_SOA_COLUMN(Svx, svx, float); //! DECLARE_SOA_COLUMN(Svy, svy, float); //! DECLARE_SOA_COLUMN(Svz, svz, float); //! DECLARE_SOA_COLUMN(Occupancy, occupancy, int); //! +DECLARE_SOA_COLUMN(RunNumber, runNumber, int); //! DECLARE_SOA_DYNAMIC_COLUMN(Pt, pt, [](float px, float py) { return RecoDecay::pt(std::array{px, py}); }); DECLARE_SOA_DYNAMIC_COLUMN(Y, y, [](float E, float pz) { return 0.5 * std::log((E + pz) / (E - pz)); }); DECLARE_SOA_DYNAMIC_COLUMN(Mass, mass, [](float E, float px, float py, float pz) { return std::sqrt(E * E - px * px - py * py - pz * pz); }); @@ -80,7 +81,8 @@ DECLARE_SOA_TABLE(HypKfColls, "AOD", "HYPKFCOLL", cent::CentFT0A, cent::CentFT0C, cent::CentFT0M, - hykfmc::Occupancy); + hykfmc::Occupancy, + hykfmc::RunNumber); using HypKfColl = HypKfColls::iterator; namespace hykftrk @@ -102,14 +104,14 @@ DECLARE_SOA_DYNAMIC_COLUMN(Y, y, [](float pt, float eta, float mass) { return st DECLARE_SOA_DYNAMIC_COLUMN(Lambda, lambda, [](float eta) { return 1. / std::cosh(eta); }); DECLARE_SOA_DYNAMIC_COLUMN(ItsNcluster, itsNcluster, [](uint32_t itsClusterSizes) { uint8_t n = 0; - for (uint8_t i = 0; i < 7; i++) { + for (uint8_t i = 0; i < 0x08; i++) { if (itsClusterSizes >> (4 * i) & 15) n++; } return n; }); DECLARE_SOA_DYNAMIC_COLUMN(ItsFirstLayer, itsFirstLayer, [](uint32_t itsClusterSizes) { - for (int i = 0; i < 8; i++) { + for (int i = 0; i < 0x08; i++) { if (itsClusterSizes >> (4 * i) & 15) return i; } @@ -117,7 +119,7 @@ DECLARE_SOA_DYNAMIC_COLUMN(ItsFirstLayer, itsFirstLayer, [](uint32_t itsClusterS }); DECLARE_SOA_DYNAMIC_COLUMN(ItsMeanClsSize, itsMeanClsSize, [](uint32_t itsClusterSizes) { int sum = 0, n = 0; - for (int i = 0; i < 8; i++) { + for (int i = 0; i < 0x08; i++) { sum += (itsClusterSizes >> (4 * i) & 15); if (itsClusterSizes >> (4 * i) & 15) n++; @@ -194,9 +196,9 @@ DECLARE_SOA_DYNAMIC_COLUMN(Eta, eta, [](float px, float py, float pz) { return R DECLARE_SOA_DYNAMIC_COLUMN(Phi, phi, [](float px, float py) { return RecoDecay::phi(std::array{px, py}); }); DECLARE_SOA_DYNAMIC_COLUMN(P, p, [](float px, float py, float pz) { return RecoDecay::p(px, py, pz); }); // DECLARE_SOA_DYNAMIC_COLUMN(Y, y, [](float px, float py, float pz, float mass) { return RecoDecay::y(std::array{px, py, pz}, mass); }); -DECLARE_SOA_DYNAMIC_COLUMN(McTrue, mcTrue, [](int hypKfMcPartId) { return hypKfMcPartId > 0; }); +DECLARE_SOA_DYNAMIC_COLUMN(McTrue, mcTrue, [](int hypKfMcPartId) { return hypKfMcPartId >= 0; }); DECLARE_SOA_DYNAMIC_COLUMN(IsMatter, isMatter, [](int8_t species) { return species > 0; }); -DECLARE_SOA_DYNAMIC_COLUMN(Cascade, cascade, [](int hypDaughter) { return hypDaughter > 0; }); +DECLARE_SOA_DYNAMIC_COLUMN(Cascade, cascade, [](int hypDaughter) { return hypDaughter >= 0; }); } // namespace hykfhyp DECLARE_SOA_TABLE(HypKfHypNucs, "AOD", "HYPKFHYPNUC", From 3da9965b1d335b436e5cafbb828b7570fc972cf5 Mon Sep 17 00:00:00 2001 From: mhartung71 <50153519+mhartung71@users.noreply.github.com> Date: Thu, 29 May 2025 12:27:36 +0200 Subject: [PATCH 2/5] PID updated --- PWGLF/TableProducer/Nuspex/hypKfRecoTask.cxx | 413 ++++++++++++------- 1 file changed, 259 insertions(+), 154 deletions(-) diff --git a/PWGLF/TableProducer/Nuspex/hypKfRecoTask.cxx b/PWGLF/TableProducer/Nuspex/hypKfRecoTask.cxx index 68e334f5766..9bef4359979 100644 --- a/PWGLF/TableProducer/Nuspex/hypKfRecoTask.cxx +++ b/PWGLF/TableProducer/Nuspex/hypKfRecoTask.cxx @@ -39,10 +39,14 @@ #include "PWGLF/DataModel/LFHypernucleiKfTables.h" #include "TRandom3.h" #include "Common/DataModel/CollisionAssociationTables.h" +#include "Common/TableProducer/PID/pidTPCBase.h" +#include "Common/DataModel/PIDResponseTPC.h" +#include "ReconstructionDataFormats/PID.h" +#include "MetadataHelper.h" // KFParticle #ifndef HomogeneousField -#define HomogeneousField // o2-linter: disable=[name/macro] +#define HomogeneousField // o2-linter: disable=name/macro (Name is defined in KFParticle package) #endif #include "KFParticle.h" #include "KFPTrack.h" @@ -54,10 +58,11 @@ using namespace o2; using namespace o2::framework; using namespace o2::framework::expressions; -using CollisionsFull = soa::Join; +using CollisionsFull = soa::Join; using CollisionsFullMC = soa::Join; -using TracksFull = soa::Join; +using TracksFull = soa::Join; +MetadataHelper metadataInfo; // Metadata helper //---------------------------------------------------------------------------------------------------------------- namespace @@ -148,7 +153,7 @@ static const std::vector cascadeNames{"4LLH->4LHe+pi", "4XHe->4LHe+ constexpr int cascadeEnabled[nCascades][1]{{0}, {0}, {0}, {0}, {0}, {0}}; constexpr int cascadePdgCodes[nCascades][1]{ {1020010040}, - {1120010040}, + {1120020040}, {0}, {0}, {0}, @@ -167,33 +172,32 @@ constexpr double preSelectionsCascades[nCascades][nSelCas]{ {0.00, 9.90, 0, 100, -1., 100., 10., 10.}, {0.00, 9.90, 0, 100, -1., 100., 10., 10.}, {0.00, 9.90, 0, 100, -1., 100., 10., 10.}}; - //---------------------------------------------------------------------------------------------------------------- struct DaughterParticle { TString name; int pdgCode, charge; double mass, resolution; - std::vector betheParams; + std::array betheParams; bool active; DaughterParticle(std::string name_, int pdgCode_, double mass_, int charge_, LabeledArray bethe) : name(name_), pdgCode(pdgCode_), charge(charge_), mass(mass_), active(false) { resolution = bethe.get(name, "resolution"); - betheParams.clear(); - for (unsigned int i = 0; i < 5; i++) - betheParams.push_back(bethe.get(name, i)); + for (unsigned int i = 0; i < betheParams.size(); i++) + betheParams[i] = bethe.get(name, i); } + int getCentralPIDIndex() { return getPIDIndex(pdgCode); } }; // struct DaughterParticle struct HyperNucleus { TString name; int pdgCode; - bool active; + bool active, savePrimary; std::vector daughters, daughterTrackSigns; - HyperNucleus(std::string name_, int pdgCode_, bool active_, std::vector daughters_, std::vector daughterTrackSigns_) : pdgCode(pdgCode_), active(active_) + HyperNucleus(std::string name_, int pdgCode_, bool active_, std::vector daughters_, std::vector daughterTrackSigns_) : pdgCode(pdgCode_), active(active_), savePrimary(active_) { init(name_, daughters_, daughterTrackSigns_); } - HyperNucleus(std::string name_, int pdgCode_, bool active_, int hypDaughter, std::vector daughters_, std::vector daughterTrackSigns_) : pdgCode(pdgCode_), active(active_) + HyperNucleus(std::string name_, int pdgCode_, bool active_, int hypDaughter, std::vector daughters_, std::vector daughterTrackSigns_) : pdgCode(pdgCode_), active(active_), savePrimary(active_) { daughters.push_back(hypDaughter); init(name_, daughters_, daughterTrackSigns_); @@ -213,16 +217,22 @@ struct HyperNucleus { struct DaughterKf { int64_t daughterTrackId; - int hypNucId; + int species, hypNucId; KFParticle daughterKfp; - float dcaToPv, dcaToPvXY, dcaToPvZ; - DaughterKf(int64_t daughterTrackId_, KFParticle daughterKfp_, std::vector vtx) : daughterTrackId(daughterTrackId_), hypNucId(-1), daughterKfp(daughterKfp_) + float dcaToPv, dcaToPvXY, dcaToPvZ, tpcNsigma, tpcNsigmaNLP, tpcNsigmaNHP; + bool active; + std::vector vtx; + DaughterKf(int species_, int64_t daughterTrackId_, std::vector vtx_, float tpcNsigma_, float tpcNsigmaNLP_, float tpcNsigmaNHP_) : daughterTrackId(daughterTrackId_), species(species_), hypNucId(-1), tpcNsigma(tpcNsigma_), tpcNsigmaNLP(tpcNsigmaNLP_), tpcNsigmaNHP(tpcNsigmaNHP_), vtx(vtx_) {} + DaughterKf(int species_, KFParticle daughterKfp_, int hypNucId_) : daughterTrackId(-999), species(species_), hypNucId(hypNucId_), daughterKfp(daughterKfp_), dcaToPv(-999), dcaToPvXY(-999), dcaToPvZ(-999) {} + void addKfp(KFParticle daughterKfp_) { + daughterKfp = daughterKfp_; dcaToPvXY = daughterKfp.GetDistanceFromVertexXY(&vtx[0]); dcaToPv = daughterKfp.GetDistanceFromVertex(&vtx[0]); dcaToPvZ = std::sqrt(dcaToPv * dcaToPv - dcaToPvXY * dcaToPvXY); } - DaughterKf(KFParticle daughterKfp_, int hypNucId_) : daughterTrackId(-999), hypNucId(hypNucId_), daughterKfp(daughterKfp_), dcaToPv(-999), dcaToPvXY(-999), dcaToPvZ(-999) {} + + bool isTrack() { return daughterTrackId >= 0; } }; // struct DaughterKf struct HyperNucCandidate { @@ -325,7 +335,7 @@ struct HyperNucCandidate { } void calcDcaToVtx(KFPVertex& vtx) { - if (devToPvXY != 999) + if (devToPvXY != 999) // o2-linter: disable=magic-number (To be checked) return; devToPvXY = kfp.GetDeviationFromVertexXY(vtx); dcaToPvXY = kfp.GetDistanceFromVertexXY(vtx); @@ -338,11 +348,23 @@ struct HyperNucCandidate { dcaToVtxZ = getDcaMotherToVtxZ(cand.recoSV); } float getSubDaughterMass(int d1, int d2) + { + return calcSubDaughterMass(daughters.at(d1)->daughterKfp, daughters.at(d2)->daughterKfp); + } + float getSubDaughterMassCascade(int d1, int d2) + { + if (!isCascade()) { + LOGF(warning, "Primary hypernucleus has no hypernucleus daughter!"); + return -999; + } + return calcSubDaughterMass(daughters.at(d1)->daughterKfp, hypNucDaughter->daughters.at(d2)->daughterKfp); + } + float calcSubDaughterMass(KFParticle d1, KFParticle d2) { KFParticle subDaughter; subDaughter.SetConstructMethod(2); - subDaughter.AddDaughter(daughters.at(d1)->daughterKfp); - subDaughter.AddDaughter(daughters.at(d2)->daughterKfp); + subDaughter.AddDaughter(d1); + subDaughter.AddDaughter(d2); subDaughter.TransportToDecayVertex(); return subDaughter.GetMass(); } @@ -363,6 +385,14 @@ struct IndexPairs { } return false; } + bool hasIndex(int64_t a) + { + for (const auto& pair : pairs) { + if (pair.first == a) + return true; + } + return false; + } }; // struct IndexPairs struct McCollInfo { @@ -419,6 +449,12 @@ struct HypKfRecoTask { Configurable> cfgPreSelectionsSecondaries{"cfgPreSelectionsSecondaries", {preSelectionsSecondaries[0], nHyperNuclei, nSelSec, hyperNucNames, preSelectionSecNames}, "selection criteria for secondary hypernuclei"}; Configurable> cfgPreSelectionsCascades{"cfgPreSelectionsCascades", {preSelectionsCascades[0], nCascades, nSelCas, cascadeNames, preSelectionCascadeNames}, "selection criteria for cascade hypernuclei"}; + // TPC PID Response + bool usePidResponse; + o2::pid::tpc::Response* response; + std::map metadata; + std::array betheParams; + // CCDB Service ccdb; Configurable bField{"bField", -999, "bz field, -999 is automatic"}; @@ -427,10 +463,10 @@ struct HypKfRecoTask { Configurable grpmagPath{"grpmagPath", "GLO/Config/GRPMagField", "CCDB path of the GRPMagField object"}; Configurable lutPath{"lutPath", "GLO/Param/MatLUT", "Path of the Lut parametrization"}; Configurable geoPath{"geoPath", "GLO/Config/GeometryAligned", "Path of the geometry file"}; - Configurable pidPath{"pidPath", "", "Path to the PID response object"}; + Configurable pidPath{"pidPath", "Analysis/PID/TPC/Response", "Path to the PID response object"}; + std::vector activePdgs; std::vector daughterParticles; - std::vector> foundDaughters; std::vector> foundDaughterKfs, hypNucDaughterKfs; std::vector> singleHyperNucCandidates, cascadeHyperNucCandidates; std::vector singleHyperNuclei, cascadeHyperNuclei; @@ -438,7 +474,7 @@ struct HypKfRecoTask { std::vector mcCollInfos; IndexPairs trackIndices, mcPartIndices; KFPVertex kfPrimVtx; - bool collHasCandidate, collHasMcTrueCandidate, collPassedEvSel, activeCascade; + bool collHasCandidate, collHasMcTrueCandidate, collPassedEvSel, activeCascade, isMC; int64_t mcCollTableIndex; int mRunNumber, occupancy; float dBz; @@ -447,6 +483,7 @@ struct HypKfRecoTask { void init(InitContext const&) { + isMC = false; mRunNumber = 0; dBz = 0; rand.SetSeed(0); @@ -456,18 +493,37 @@ struct HypKfRecoTask { ccdb->setLocalObjectValidityChecking(); ccdb->setFatalWhenNull(false); - for (int i = 0; i < nDaughterParticles; i++) { // create daughterparticles + usePidResponse = false; + for (unsigned int i = 0; i < nDaughterParticles; i++) { // create daughterparticles daughterParticles.push_back(DaughterParticle(particleNames.at(i), particlePdgCodes.at(i), particleMasses.at(i), particleCharge.at(i), cfgBetheBlochParams)); + if (cfgTrackPIDsettings->get(i, "useBBparams") == 2 || cfgTrackPIDsettings->get(i, "useBBparams") == 0) // o2-linter: disable=magic-number (To be checked) + usePidResponse = true; } + for (unsigned int i = 0; i < nHyperNuclei; i++) { // create hypernuclei - singleHyperNuclei.push_back(HyperNucleus(hyperNucNames.at(i), cfgHyperNucPdg->get(i, 0u), cfgHyperNucsActive->get(i, 0u), getDaughterVec(i, cfgHyperNucDaughters), getDaughterSignVec(i, cfgHyperNucSigns))); + auto active = cfgHyperNucsActive->get(i, 0u); + auto pdg = cfgHyperNucPdg->get(i, 0u); + singleHyperNuclei.push_back(HyperNucleus(hyperNucNames.at(i), pdg, active, getDaughterVec(i, cfgHyperNucDaughters), getDaughterSignVec(i, cfgHyperNucSigns))); + if (active) + activePdgs.push_back(pdg); } + activeCascade = false; for (unsigned int i = 0; i < nCascades; i++) { // create cascades - cascadeHyperNuclei.push_back(HyperNucleus(cascadeNames.at(i), cfgCascadesPdg->get(i, 0u), cfgCascadesActive->get(i, 0u), getHypDaughterVec(i, cfgCascadeHypDaughter), getDaughterVec(i, cfgCascadeDaughters), getDaughterSignVec(i, cfgCascadeSigns))); - if (cfgCascadesActive->get(i, 0u)) + auto active = cfgCascadesActive->get(i, 0u); + auto pdg = cfgCascadesPdg->get(i, 0u); + auto hypDaughter = getHypDaughterVec(i, cfgCascadeHypDaughter); + cascadeHyperNuclei.push_back(HyperNucleus(cascadeNames.at(i), pdg, active, hypDaughter, getDaughterVec(i, cfgCascadeDaughters), getDaughterSignVec(i, cfgCascadeSigns))); + if (active) { + activePdgs.push_back(pdg); + if (!singleHyperNuclei.at(hypDaughter).active) { + singleHyperNuclei.at(hypDaughter).active = true; + activePdgs.push_back(singleHyperNuclei.at(hypDaughter).pdgCode); + } activeCascade = true; + } } + // define histogram axes const AxisSpec axisMagField{10, -10., 10., "magnetic field"}; const AxisSpec axisNev{3, 0., 3., "Number of events"}; @@ -506,7 +562,7 @@ struct HypKfRecoTask { } //---------------------------------------------------------------------------------------------------------------- - void findDaughterParticles(aod::TrackAssoc const& tracksByColl, TracksFull const& tracks) + void findDaughterParticles(aod::TrackAssoc const& tracksByColl, TracksFull const& tracks, CollisionsFull::iterator const& coll) { // track loop, store daughter candidates in std::vector for (const auto& trackId : tracksByColl) { @@ -525,63 +581,69 @@ struct HypKfRecoTask { continue; if (track.itsChi2NCl() > cfgTrackPIDsettings->get(i, "maxITSchi2")) continue; - if (std::abs(getTPCnSigma(track, daughterParticles.at(i))) > cfgTrackPIDsettings->get(i, "maxTPCnSigma")) + if (getRigidity(track) < cfgTrackPIDsettings->get(i, "minRigidity") || getRigidity(track) > cfgTrackPIDsettings->get(i, "maxRigidity")) + continue; + float tpcNsigma = getTPCnSigma(track, coll, daughterParticles.at(i)); + if (std::abs(tpcNsigma) > cfgTrackPIDsettings->get(i, "maxTPCnSigma")) continue; filldedx(track, i); if (getMeanItsClsSize(track) < cfgTrackPIDsettings->get(i, "minITSclsSize")) continue; if (getMeanItsClsSize(track) > cfgTrackPIDsettings->get(i, "maxITSclsSize")) continue; - if (getRigidity(track) < cfgTrackPIDsettings->get(i, "minRigidity") || getRigidity(track) > cfgTrackPIDsettings->get(i, "maxRigidity")) - continue; if (cfgTrackPIDsettings->get(i, "TOFrequiredabove") >= 0 && getRigidity(track) > cfgTrackPIDsettings->get(i, "TOFrequiredabove") && (track.mass() < cfgTrackPIDsettings->get(i, "minTOFmass") || track.mass() > cfgTrackPIDsettings->get(i, "maxTOFmass"))) continue; - foundDaughters.at(i).push_back(track.globalIndex()); + float tpcNsigmaNHP = (i == kAlpha ? -999 : getTPCnSigma(track, coll, daughterParticles.at(i + 1))); + float tpcNsigmaNLP = (i == kPion ? 999 : getTPCnSigma(track, coll, daughterParticles.at(i - 1))); + foundDaughterKfs.at(i).push_back(DaughterKf(i, track.globalIndex(), primVtx, tpcNsigma, tpcNsigmaNLP, tpcNsigmaNHP)); } } // track loop } - //---------------------------------------------------------------------------------------------------------------- - void checkMCTrueTracks(aod::McTrackLabels const& trackLabels, aod::McParticles const&) + + void checkMCTrueTracks(aod::McTrackLabels const& trackLabels, aod::McParticles const&, TracksFull const& tracks, CollisionsFull::iterator const& coll) { - std::vector activePdgs; - std::vector*> hypNucVectors = {&singleHyperNuclei, &cascadeHyperNuclei}; - for (int vec = 0; vec < 2; vec++) { - for (size_t hyperNucIter = 0; hyperNucIter < hypNucVectors.at(vec)->size(); hyperNucIter++) { - HyperNucleus* hyperNuc = &(hypNucVectors.at(vec)->at(hyperNucIter)); - if (!hyperNuc->active) - continue; - activePdgs.push_back(std::abs(hyperNuc->pdgCode)); - } - } for (int i = 0; i < nDaughterParticles; i++) { - auto& daughterVec = foundDaughters.at(i); + auto& daughterVec = foundDaughterKfs.at(i); if (!daughterVec.size()) continue; for (auto it = daughterVec.end() - 1; it >= daughterVec.begin(); it--) { - const auto& mcLab = trackLabels.rawIteratorAt(*it); + const auto& mcLab = trackLabels.rawIteratorAt(it->daughterTrackId); if (!mcLab.has_mcParticle()) { daughterVec.erase(it); continue; } const auto& mcPart = mcLab.mcParticle_as(); - if (std::abs(mcPart.pdgCode()) != daughterParticles.at(i).pdgCode) { - daughterVec.erase(it); - continue; - } - if (!mcPart.has_mothers()) { - daughterVec.erase(it); - continue; - } - bool isDaughter = false; - for (const auto& mother : mcPart.mothers_as()) { - if (std::find(activePdgs.begin(), activePdgs.end(), std::abs(mother.pdgCode())) != activePdgs.end()) { - isDaughter = true; + if (cfgSaveOnlyMcTrue) { + if (std::abs(mcPart.pdgCode()) != daughterParticles.at(i).pdgCode) { + daughterVec.erase(it); + continue; + } + if (!mcPart.has_mothers()) { + daughterVec.erase(it); + continue; + } + bool isDaughter = false; + for (const auto& mother : mcPart.mothers_as()) { + if (std::find(activePdgs.begin(), activePdgs.end(), std::abs(mother.pdgCode())) != activePdgs.end()) { + isDaughter = true; + } + } + if (!isDaughter) { + daughterVec.erase(it); + continue; } } - if (!isDaughter) { - daughterVec.erase(it); - continue; + if (cfgTrackPIDsettings->get(i, "useBBparams") == 0) { + const auto& trk = tracks.rawIteratorAt(it->daughterTrackId); + const auto tpcNsigmaMC = getTPCnSigmaMC(trk, coll, daughterParticles.at(i), daughterParticles.at(i)); + if (std::abs(tpcNsigmaMC) <= cfgTrackPIDsettings->get(i, "maxTPCnSigma")) { + it->tpcNsigma = tpcNsigmaMC; + it->tpcNsigmaNHP = (i == kAlpha ? -999 : getTPCnSigmaMC(trk, coll, daughterParticles.at(i), daughterParticles.at(i + 1))); + it->tpcNsigmaNLP = (i == kPion ? 999 : getTPCnSigmaMC(trk, coll, daughterParticles.at(i), daughterParticles.at(i - 1))); + } else { + daughterVec.erase(it); + } } } } @@ -589,33 +651,42 @@ struct HypKfRecoTask { //---------------------------------------------------------------------------------------------------------------- void createKFDaughters(TracksFull const& tracks) { + for (size_t daughterCount = 0; daughterCount < daughterParticles.size(); daughterCount++) { + daughterParticles.at(daughterCount).active = false; + } std::vector*> hypNucVectors = {&singleHyperNuclei, &cascadeHyperNuclei}; - for (size_t vec = 0; vec < 2; vec++) { + bool singleHypNucActive = false; + for (size_t vec = 0; vec < hypNucVectors.size(); vec++) { for (const auto& hyperNuc : *(hypNucVectors.at(vec))) { if (!hyperNuc.active) continue; for (size_t i = vec; i < hyperNuc.daughters.size(); i++) { - if (foundDaughters.at(hyperNuc.daughters.at(i)).size() > 0) + if (foundDaughterKfs.at(hyperNuc.daughters.at(i)).size() > 0) daughterParticles.at(hyperNuc.daughters.at(i)).active = true; else break; + if (i == hyperNuc.daughters.size() - 1) + singleHypNucActive = true; } } + if (!singleHypNucActive) + break; } + for (size_t daughterCount = 0; daughterCount < daughterParticles.size(); daughterCount++) { const auto& daughterParticle = daughterParticles.at(daughterCount); - if (!daughterParticle.active) + if (!daughterParticle.active) { + foundDaughterKfs.at(daughterCount).clear(); continue; + } const auto& daughterMass = daughterParticle.mass; const auto& daughterCharge = daughterParticle.charge; - for (const auto& daughterId : foundDaughters.at(daughterCount)) { - const auto& daughterTrack = tracks.rawIteratorAt(daughterId); - DaughterKf daughter(daughterId, createKFParticle(daughterTrack, daughterMass, daughterCharge), primVtx); - if (std::abs(daughter.dcaToPvXY) < cfgTrackPIDsettings->get(daughterCount, "minDcaToPvXY")) - continue; - if (std::abs(daughter.dcaToPvZ) < cfgTrackPIDsettings->get(daughterCount, "minDcaToPvZ")) - continue; - foundDaughterKfs.at(daughterCount).push_back(daughter); + auto& daughterVec = foundDaughterKfs.at(daughterCount); + for (auto it = daughterVec.end() - 1; it >= daughterVec.begin(); it--) { + const auto& daughterTrack = tracks.rawIteratorAt(it->daughterTrackId); + it->addKfp(createKFParticle(daughterTrack, daughterMass, daughterCharge)); + if (std::abs(it->dcaToPvXY) < cfgTrackPIDsettings->get(daughterCount, "minDcaToPvXY") || std::abs(it->dcaToPvZ) < cfgTrackPIDsettings->get(daughterCount, "minDcaToPvZ")) + daughterVec.erase(it); } } } @@ -685,7 +756,7 @@ struct HypKfRecoTask { candidate.calcDcaToVtx(kfPrimVtx); candidate.isSecondaryCandidate = true; } - if (candidate.isPrimaryCandidate || candidate.isSecondaryCandidate) + if ((candidate.isPrimaryCandidate && hyperNuc->savePrimary) || (candidate.isSecondaryCandidate && activeCascade)) singleHyperNucCandidates.at(hyperNucIter).push_back(candidate); } } @@ -711,7 +782,7 @@ struct HypKfRecoTask { for (int64_t i = 0; i < static_cast(nHypNucDaughters); i++) { if (singleHyperNucCandidates.at(hyperNuc->daughters.at(0)).at(i).isSecondaryCandidate) { auto hypNucDaughter = &(singleHyperNucCandidates.at(hyperNuc->daughters.at(0)).at(i)); - hypNucDaughterKfs.at(hyperNucIter).push_back(DaughterKf(hypNucDaughter->kfp, i)); + hypNucDaughterKfs.at(hyperNucIter).push_back(DaughterKf(hyperNuc->daughters.at(0), hypNucDaughter->kfp, i)); } } int nCombinations = hypNucDaughterKfs.at(hyperNucIter).size(); @@ -796,9 +867,8 @@ struct HypKfRecoTask { HyperNucleus* hyperNuc = &(hypNucVectors.at(vec)->at(hyperNucIter)); if (!hyperNuc->active) continue; - for (auto& hypCand : candidateVector->at(hyperNucIter)) { // o2-linter: disable=[const-ref-in-for-loop] + for (auto& hypCand : candidateVector->at(hyperNucIter)) { // o2-linter: disable=[const-ref-in-for-loop] (Object is non const and modified in loop) std::vector motherIds; - int daughterCount = 0; if (hypCand.isCascade()) { if (!hypCand.hypNucDaughter->mcTrue) continue; @@ -811,16 +881,15 @@ struct HypKfRecoTask { break; } } - daughterCount++; } for (const auto& daughter : hypCand.daughters) { - if (daughter->daughterTrackId < 0) + if (!daughter->isTrack()) continue; const auto& mcLab = trackLabels.rawIteratorAt(daughter->daughterTrackId); if (!mcLab.has_mcParticle()) continue; const auto& mcPart = mcLab.mcParticle_as(); - if (std::abs(mcPart.pdgCode()) != daughterParticles.at(hyperNuc->daughters.at(daughterCount++)).pdgCode) + if (std::abs(mcPart.pdgCode()) != daughterParticles.at(daughter->species).pdgCode) continue; if (!mcPart.has_mothers()) continue; @@ -840,6 +909,8 @@ struct HypKfRecoTask { for (auto iter = motherIds.begin(); iter != motherIds.end() - 1; iter++) if (*iter != *(iter + 1)) hypCand.mcTrue = false; + if (!mcPartIndices.hasIndex(motherIds.front())) + hypCand.mcTrue = false; if (hypCand.mcTrue) { hypCand.mcParticleId = motherIds.front(); collHasMcTrueCandidate = true; @@ -858,49 +929,39 @@ struct HypKfRecoTask { outputCollisionTable( collPassedEvSel, mcCollTableIndex, primVtx.at(0), primVtx.at(1), primVtx.at(2), - cents.at(0), cents.at(1), cents.at(2), occupancy); + cents.at(0), cents.at(1), cents.at(2), occupancy, mRunNumber); std::vector*> hypNucVectors = {&singleHyperNuclei, &cascadeHyperNuclei}; std::vector>*> candidateVectors = {&singleHyperNucCandidates, &cascadeHyperNucCandidates}; - for (int vec = 0; vec < 2; vec++) { + for (int vec = 0; vec < candidateVectors.size(); vec++) { auto candidateVector = candidateVectors.at(vec); for (size_t hyperNucIter = 0; hyperNucIter < hypNucVectors.at(vec)->size(); hyperNucIter++) { HyperNucleus* hyperNuc = &(hypNucVectors.at(vec)->at(hyperNucIter)); if (!hyperNuc->active) continue; - for (auto& hypCand : candidateVector->at(hyperNucIter)) { // o2-linter: disable=[const-ref-in-for-loop] + for (auto& hypCand : candidateVector->at(hyperNucIter)) { // o2-linter: disable=const-ref-in-for-loop (Object is non const and modified in loop) if (!hypCand.isPrimaryCandidate && !hypCand.isUsedSecondary && !hypCand.isCascade()) continue; if (saveOnlyMcTrue && !hypCand.mcTrue && !hypCand.isCascade()) continue; hInvMass[vec * nHyperNuclei + hyperNucIter]->Fill(hypCand.mass); std::vector vecDaugtherTracks, vecAddons, vecSubDaughters; - int daughterCount = 0; for (const auto& daughter : hypCand.daughters) { - const auto& daughterTrackId = daughter->daughterTrackId; - if (daughterTrackId < 0) + if (!daughter->isTrack()) continue; + const auto& daughterTrackId = daughter->daughterTrackId; int trackTableId; if (!trackIndices.getIndex(daughterTrackId, trackTableId)) { - auto daught = hyperNuc->daughters.at(daughterCount); const auto& track = tracks.rawIteratorAt(daughterTrackId); outputTrackTable( - hyperNuc->daughters.at(daughterCount) * track.sign(), - track.pt(), track.eta(), track.phi(), - daughter->dcaToPvXY, daughter->dcaToPvZ, - track.tpcNClsFound(), track.tpcChi2NCl(), - track.itsClusterSizes(), track.itsChi2NCl(), - getRigidity(track), track.tpcSignal(), getTPCnSigma(track, daughterParticles.at(daught)), - daught == kAlpha ? -999 : getTPCnSigma(track, daughterParticles.at(daught + 1)), - daught == kPion ? 999 : getTPCnSigma(track, daughterParticles.at(daught - 1)), - track.mass(), - track.isPVContributor()); + daughter->species * track.sign(), track.pt(), track.eta(), track.phi(), daughter->dcaToPvXY, daughter->dcaToPvZ, track.tpcNClsFound(), track.tpcChi2NCl(), + track.itsClusterSizes(), track.itsChi2NCl(), getRigidity(track), track.tpcSignal(), daughter->tpcNsigma, daughter->tpcNsigmaNHP, daughter->tpcNsigmaNLP, + track.mass(), track.isPVContributor()); trackTableId = outputTrackTable.lastIndex(); trackIndices.add(daughterTrackId, trackTableId); } vecDaugtherTracks.push_back(trackTableId); - daughterCount++; } for (int i = 0; i < hypCand.getNdaughters(); i++) { std::vector& posMom = hypCand.daughterPosMoms.at(i); @@ -908,7 +969,7 @@ struct HypKfRecoTask { posMom.at(0), posMom.at(1), posMom.at(2), posMom.at(3), posMom.at(4), posMom.at(5)); vecAddons.push_back(outputDaughterAddonTable.lastIndex()); } - if (hypCand.getNdaughters() > 2) { + if (hypCand.getNdaughters() > 2) { // o2-linter: disable=magic-number (To be checked) for (int i = 0; i < hypCand.getNdaughters(); i++) { for (int j = i + 1; j < hypCand.getNdaughters(); j++) { outputSubDaughterTable(hypCand.getSubDaughterMass(i, j)); @@ -916,18 +977,23 @@ struct HypKfRecoTask { } } } + if (hypCand.isCascade()) { + for (int i = 1; i < hypCand.getNdaughters(); i++) { + for (int j = 0; j < hypCand.hypNucDaughter->getNdaughters(); j++) { + outputSubDaughterTable(hypCand.getSubDaughterMassCascade(i, j)); + vecSubDaughters.push_back(outputSubDaughterTable.lastIndex()); + } + } + } hypCand.kfp.TransportToDecayVertex(); int mcPartTableId; outputHypNucTable( mcPartIndices.getIndex(hypCand.mcParticleId, mcPartTableId) ? mcPartTableId : -1, outputCollisionTable.lastIndex(), vecDaugtherTracks, vecAddons, hypCand.getDaughterTableId(), vecSubDaughters, - (vec * nHyperNuclei + hyperNucIter + 1) * hypCand.getSign(), - hypCand.isPrimaryCandidate, hypCand.mass, - hypCand.px, hypCand.py, hypCand.pz, - hypCand.dcaToPvXY, hypCand.dcaToPvZ, hypCand.devToPvXY, - hypCand.dcaToVtxXY, hypCand.dcaToVtxZ, hypCand.chi2, - hypCand.recoSV.at(0), hypCand.recoSV.at(1), hypCand.recoSV.at(2)); + (vec * nHyperNuclei + hyperNucIter + 1) * hypCand.getSign(), hypCand.isPrimaryCandidate, hypCand.mass, + hypCand.px, hypCand.py, hypCand.pz, hypCand.dcaToPvXY, hypCand.dcaToPvZ, hypCand.devToPvXY, + hypCand.dcaToVtxXY, hypCand.dcaToVtxZ, hypCand.chi2, hypCand.recoSV.at(0), hypCand.recoSV.at(1), hypCand.recoSV.at(2)); hypCand.tableId = outputHypNucTable.lastIndex(); } } @@ -935,21 +1001,23 @@ struct HypKfRecoTask { } //---------------------------------------------------------------------------------------------------------------- - void processMC(CollisionsFullMC const& collisions, aod::McCollisions const& mcColls, TracksFull const& tracks, aod::BCsWithTimestamps const&, aod::McParticles const& particlesMC, aod::McTrackLabels const& trackLabelsMC, aod::McCollisionLabels const& collLabels, aod::TrackAssoc const& tracksColl) + void processMC(CollisionsFullMC const& collisions, aod::McCollisions const& mcColls, TracksFull const& tracks, aod::BCsWithTimestamps const&, aod::McParticles const& particlesMC, aod::McTrackLabels const& trackLabelsMC, aod::McCollisionLabels const& collLabels, aod::TrackAssoc const& tracksColl, CollisionsFull const& colls) { - + isMC = true; mcCollInfos.clear(); mcCollInfos.resize(mcColls.size()); mcPartIndices.clear(); for (const auto& collision : collisions) { if (!collision.has_mcCollision()) continue; - if (collision.sel8() && std::abs(collision.posZ()) < 10) + if (collision.sel8() && std::abs(collision.posZ()) < 10) // o2-linter: disable=magic-number (To be checked) mcCollInfos.at(collision.mcCollisionId()).passedEvSel = true; } std::vector*> hypNucVectors = {&singleHyperNuclei, &cascadeHyperNuclei}; for (const auto& mcPart : particlesMC) { - for (int vec = 0; vec < 2; vec++) { + if (!mcCollInfos.at(mcPart.mcCollisionId()).passedEvSel) + continue; + for (int vec = 0; vec < hypNucVectors.size(); vec++) { for (size_t hyperNucIter = 0; hyperNucIter < hypNucVectors.at(vec)->size(); hyperNucIter++) { HyperNucleus* hyperNuc = &(hypNucVectors.at(vec)->at(hyperNucIter)); if (!hyperNuc->active) @@ -958,11 +1026,7 @@ struct HypKfRecoTask { continue; bool isDecayMode = false; float svx, svy, svz; - int daughterPdg; - if (vec == 0) - daughterPdg = daughterParticles.at(hyperNuc->daughters.at(0)).pdgCode; - else - daughterPdg = singleHyperNuclei.at(hyperNuc->daughters.at(0)).pdgCode; + int daughterPdg = vec ? singleHyperNuclei.at(hyperNuc->daughters.at(0)).pdgCode : daughterParticles.at(hyperNuc->daughters.at(0)).pdgCode; for (const auto& mcDaught : mcPart.daughters_as()) { if (std::abs(mcDaught.pdgCode()) == daughterPdg) { isDecayMode = true; @@ -999,11 +1063,14 @@ struct HypKfRecoTask { auto bc = collision.bc_as(); initCCDB(bc); initCollision(collision); + if (!collision.has_mcCollision() || !mcCollInfos.at(collision.mcCollisionId()).passedEvSel) + continue; + const uint64_t collIdx = collision.globalIndex(); auto tracksByColl = tracksColl.sliceBy(perCollision, collIdx); - findDaughterParticles(tracksByColl, tracks); - if (cfgSaveOnlyMcTrue) - checkMCTrueTracks(trackLabelsMC, particlesMC); + findDaughterParticles(tracksByColl, tracks, colls.rawIteratorAt(collision.globalIndex())); + if (cfgSaveOnlyMcTrue || usePidResponse) + checkMCTrueTracks(trackLabelsMC, particlesMC, tracks, colls.rawIteratorAt(collision.globalIndex())); createKFDaughters(tracks); createKFHypernuclei(tracks); createMCinfo(trackLabelsMC, collLabels, particlesMC, mcColls); @@ -1015,16 +1082,13 @@ struct HypKfRecoTask { if (cfgSaveOnlyMcTrue && !collHasMcTrueCandidate) continue; - mcCollTableIndex = -1; - if (collision.has_mcCollision()) { - mcCollTableIndex = mcCollInfos.at(collision.mcCollisionId()).tableIndex; - if (mcCollTableIndex < 0) { - outputMcCollisionTable( - mcCollInfos.at(collision.mcCollisionId()).passedEvSel, - collision.mcCollision().posX(), collision.mcCollision().posY(), collision.mcCollision().posZ()); - mcCollTableIndex = outputMcCollisionTable.lastIndex(); - mcCollInfos.at(collision.mcCollisionId()).tableIndex = mcCollTableIndex; - } + mcCollTableIndex = mcCollInfos.at(collision.mcCollisionId()).tableIndex; + if (mcCollTableIndex < 0) { + outputMcCollisionTable( + mcCollInfos.at(collision.mcCollisionId()).passedEvSel, + collision.mcCollision().posX(), collision.mcCollision().posY(), collision.mcCollision().posZ()); + mcCollTableIndex = outputMcCollisionTable.lastIndex(); + mcCollInfos.at(collision.mcCollisionId()).tableIndex = mcCollTableIndex; } fillTree(tracks, cfgSaveOnlyMcTrue); } @@ -1042,7 +1106,7 @@ struct HypKfRecoTask { continue; const uint64_t collIdx = collision.globalIndex(); auto tracksByColl = tracksColl.sliceBy(perCollision, collIdx); - findDaughterParticles(tracksByColl, tracks); + findDaughterParticles(tracksByColl, tracks, collision); createKFDaughters(tracks); createKFHypernuclei(tracks); createKFCascades(tracks); @@ -1065,7 +1129,7 @@ struct HypKfRecoTask { o2::parameters::GRPMagField* grpmag = 0x0; if (grpo) { o2::base::Propagator::initFieldFromGRP(grpo); - if (bField < -990) { + if (bField < -990) { // o2-linter: disable=magic-number (To be checked) // Fetch magnetic field from ccdb for current collision dBz = grpo->getNominalL3Field(); LOG(info) << "Retrieved GRP for timestamp " << run3grpTimestamp << " with magnetic field of " << dBz << " kZG"; @@ -1078,7 +1142,7 @@ struct HypKfRecoTask { LOG(fatal) << "Got nullptr from CCDB for path " << grpmagPath << " of object GRPMagField and " << grpPath << " of object GRPObject for timestamp " << run3grpTimestamp; } o2::base::Propagator::initFieldFromGRP(grpmag); - if (bField < -990) { + if (bField < -990) { // o2-linter: disable=magic-number (To be checked) // Fetch magnetic field from ccdb for current collision dBz = std::lround(5.f * grpmag->getL3Current() / 30000.f); LOG(info) << "Retrieved GRP for timestamp " << run3grpTimestamp << " with magnetic field of " << dBz << " kZG"; @@ -1088,13 +1152,35 @@ struct HypKfRecoTask { } mRunNumber = bc.runNumber(); KFParticle::SetField(dBz); + + // PID response + if (!usePidResponse) + return; + if (metadataInfo.isFullyDefined()) { + metadata["RecoPassName"] = metadataInfo.get("RecoPassName"); + LOGP(info, "Automatically setting reco pass for TPC Response to {} from AO2D", metadata["RecoPassName"]); + } else { + LOG(info) << "Setting reco pass for TPC response to default name"; + metadata["RecoPassName"] = "apass5"; + } + const std::string path = pidPath.value; + ccdb->setTimestamp(run3grpTimestamp); + response = ccdb->getSpecific(path, run3grpTimestamp, metadata); + if (!response) { + LOGF(warning, "Unable to find TPC parametrisation for specified pass name - falling back to latest object"); + response = ccdb->getForTimeStamp(path, run3grpTimestamp); + if (!response) { + LOGF(fatal, "Unable to find any TPC object corresponding to timestamp {}!", run3grpTimestamp); + } + } + LOG(info) << "Successfully retrieved TPC PID object from CCDB for timestamp " << run3grpTimestamp << ", recoPass " << metadata["RecoPassName"]; + response->PrintAll(); + betheParams = response->GetBetheBlochParams(); } //---------------------------------------------------------------------------------------------------------------- template void initCollision(const T& collision) { - foundDaughters.clear(); - foundDaughters.resize(nDaughterParticles); foundDaughterKfs.clear(); foundDaughterKfs.resize(nDaughterParticles); hypNucDaughterKfs.clear(); @@ -1108,7 +1194,7 @@ struct HypKfRecoTask { collHasMcTrueCandidate = false; histos.fill(HIST("histMagField"), dBz); histos.fill(HIST("histNev"), 0.5); - collPassedEvSel = collision.sel8() && std::abs(collision.posZ()) < 10; + collPassedEvSel = collision.sel8() && std::abs(collision.posZ()) < 10; // o2-linter: disable=magic-number (To be checked) occupancy = collision.trackOccupancyInTimeRange(); if (collPassedEvSel) { histos.fill(HIST("histNev"), 1.5); @@ -1129,45 +1215,63 @@ struct HypKfRecoTask { { const float rigidity = getRigidity(track); hDeDx[2 * species]->Fill(track.sign() * rigidity, track.tpcSignal()); - if (track.tpcNClsFound() < 100 || track.itsNCls() < 2) + if (track.tpcNClsFound() < 100 || track.itsNCls() < 2) // o2-linter: disable=magic-number (To be checked) return; hDeDx[2 * species + 1]->Fill(track.sign() * rigidity, track.tpcSignal()); } //---------------------------------------------------------------------------------------------------------------- - template - float getTPCnSigma(T const& track, DaughterParticle const& particle) + template + float getTPCnSigma(T const& track, C const& coll, DaughterParticle& particle) { const float rigidity = getRigidity(track); if (!track.hasTPC()) return -999; + float mMip = 1, chargeFactor = 1; + float* parBB; - if (particle.name == "pion" && cfgTrackPIDsettings->get("pion", "useBBparams") < 1) - return cfgTrackPIDsettings->get("pion", "useBBparams") == 0 ? track.tpcNSigmaPi() : 0; - if (particle.name == "proton" && cfgTrackPIDsettings->get("proton", "useBBparams") < 1) - return cfgTrackPIDsettings->get("proton", "useBBparams") == 0 ? track.tpcNSigmaPr() : 0; - if (particle.name == "deuteron" && cfgTrackPIDsettings->get("deuteron", "useBBparams") < 1) - return cfgTrackPIDsettings->get("deuteron", "useBBparams") == 0 ? track.tpcNSigmaDe() : 0; - if (particle.name == "triton" && cfgTrackPIDsettings->get("triton", "useBBparams") < 1) - return cfgTrackPIDsettings->get("triton", "useBBparams") == 0 ? track.tpcNSigmaTr() : 0; - if (particle.name == "helion" && cfgTrackPIDsettings->get("helion", "useBBparams") < 1) - return cfgTrackPIDsettings->get("helion", "useBBparams") == 0 ? track.tpcNSigmaHe() : 0; - if (particle.name == "alpha" && cfgTrackPIDsettings->get("alpha", "useBBparams") < 1) - return cfgTrackPIDsettings->get("alpha", "useBBparams") == 0 ? track.tpcNSigmaAl() : 0; - - double expBethe{tpc::BetheBlochAleph(static_cast(particle.charge * rigidity / particle.mass), particle.betheParams[0], particle.betheParams[1], particle.betheParams[2], particle.betheParams[3], particle.betheParams[4])}; + switch (static_cast(cfgTrackPIDsettings->get(particle.name, "useBBparams"))) { + case -1: + return 0; + case 0: + return isMC ? 0 : response->GetNumberOfSigma(coll, track, particle.getCentralPIDIndex()); + case 1: + parBB = &particle.betheParams[0]; + break; + case 2: + mMip = response->GetMIP(); + chargeFactor = std::pow(particle.charge, response->GetChargeFactor()); + parBB = &betheParams[0]; + break; + default: + return -999; + } + double expBethe{mMip * chargeFactor * o2::tpc::BetheBlochAleph(static_cast(particle.charge * rigidity / particle.mass), parBB[0], parBB[1], parBB[2], parBB[3], parBB[4])}; double expSigma{expBethe * particle.resolution}; float sigmaTPC = static_cast((track.tpcSignal() - expBethe) / expSigma); return sigmaTPC; } //---------------------------------------------------------------------------------------------------------------- + + template + float getTPCnSigmaMC(T const& trk, C const& coll, DaughterParticle& particle1, DaughterParticle& particle2) + { + const float pidval1 = particle1.getCentralPIDIndex(); + const float pidval2 = particle2.getCentralPIDIndex(); + const auto expSignal = response->GetExpectedSignal(trk, pidval2); + const auto expSigma = response->GetExpectedSigma(coll, trk, pidval2); + const auto mcTunedTPCSignal = gRandom->Gaus(expSignal, expSigma); + return response->GetNumberOfSigmaMCTuned(coll, trk, pidval1, mcTunedTPCSignal); + } + //---------------------------------------------------------------------------------------------------------------- + template float getMeanItsClsSize(T const& track) { int sum = 0, n = 0; - for (int i = 0; i < 8; i++) { - sum += (track.itsClusterSizes() >> (4 * i) & 15); - if (track.itsClusterSizes() >> (4 * i) & 15) + for (int i = 0; i < 0x08; i++) { + sum += (track.itsClusterSizes() >> (0x04 * i) & 0x0f); + if (track.itsClusterSizes() >> (0x04 * i) & 0x0f) n++; } return n > 0 ? static_cast(sum) / n : 0.f; @@ -1192,7 +1296,7 @@ struct HypKfRecoTask { trackparCov.getXYZGlo(fP); trackparCov.getPxPyPzGlo(fM); float fPM[6]; - for (int i = 0; i < 3; i++) { + for (int i = 0; i < 0x03; i++) { fPM[i] = fP[i]; fPM[i + 3] = fM[i] * std::abs(charge); } @@ -1228,7 +1332,7 @@ struct HypKfRecoTask { std::vector getDaughterVec(unsigned int hypNuc, LabeledArray cfg) { std::vector vec; - for (unsigned int i = 0; i < 4; i++) { + for (unsigned int i = 0; i < 0x04; i++) { std::string daughter = cfg.get(hypNuc, i); if (std::find(particleNames.begin(), particleNames.end(), daughter) == particleNames.end()) break; @@ -1241,7 +1345,7 @@ struct HypKfRecoTask { std::vector getDaughterSignVec(unsigned int hypNuc, LabeledArray cfg) { std::vector vec; - for (unsigned int i = 0; i < 4; i++) { + for (unsigned int i = 0; i < 0x04; i++) { std::string sign = cfg.get(hypNuc, i); if (sign != "+" && sign != "-") break; @@ -1256,6 +1360,7 @@ struct HypKfRecoTask { //---------------------------------------------------------------------------------------------------------------- WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { + metadataInfo.initMetadata(cfgc); // Parse AO2D metadata return WorkflowSpec{ adaptAnalysisTask(cfgc)}; } From c9fbc2c49e5a3122e5fd34776cc10c1f8679c53f Mon Sep 17 00:00:00 2001 From: mhartung71 <50153519+mhartung71@users.noreply.github.com> Date: Thu, 29 May 2025 12:28:07 +0200 Subject: [PATCH 3/5] Update hypKfTreeCreator.cxx --- .../TableProducer/Nuspex/hypKfTreeCreator.cxx | 100 +++++++++++------- 1 file changed, 61 insertions(+), 39 deletions(-) diff --git a/PWGLF/TableProducer/Nuspex/hypKfTreeCreator.cxx b/PWGLF/TableProducer/Nuspex/hypKfTreeCreator.cxx index da8316e3c8d..81a7569ed27 100644 --- a/PWGLF/TableProducer/Nuspex/hypKfTreeCreator.cxx +++ b/PWGLF/TableProducer/Nuspex/hypKfTreeCreator.cxx @@ -56,15 +56,15 @@ struct TrackProperties { }; struct HyperNucleus { - HyperNucleus() : pdgCode(0), isReconstructed(0), globalIndex(0), species(0), isPrimaryCandidate(0), isMatter(0), passedEvSel(0), isMatterMC(0), passedEvSelMC(0), isPhysicalPrimary(0), collisionMcTrue(0), mass(0), y(0), pt(0), ct(0), yGen(0), ptGen(0), ctGen(0), cpaPvGen(0), cpaPv(0), cpaSv(0), maxDcaTracks(0), maxDcaTracksSv(0), dcaToPvXY(0), dcaToPvZ(0), dcaToVtxXY(0), dcaToVtxZ(0), devToPvXY(0), chi2(0), pvx(0), pvy(0), pvz(0), svx(0), svy(0), svz(0), px(0), py(0), pz(0), pvxGen(0), pvyGen(0), pvzGen(0), svxGen(0), svyGen(0), svzGen(0), pxGen(0), pyGen(0), pzGen(0), nSingleDaughters(0), nCascadeDaughters(0), mcTrue(0), mcTrueVtx(0), mcPhysicalPrimary(0), hypNucDaughter(0) {} + HyperNucleus() : pdgCode(0), isReconstructed(0), globalIndex(0), species(0), speciesMC(0), isPrimaryCandidate(0), isMatter(0), isCascade(0), isCascadeMC(0), passedEvSel(0), isMatterMC(0), passedEvSelMC(0), isPhysicalPrimary(0), collisionMcTrue(0), mass(0), y(0), pt(0), ct(0), yGen(0), ptGen(0), ctGen(0), cpaPvGen(0), cpaPv(0), cpaSv(0), maxDcaTracks(0), maxDcaTracksSv(0), dcaToPvXY(0), dcaToPvZ(0), dcaToVtxXY(0), dcaToVtxZ(0), devToPvXY(0), chi2(0), pvx(0), pvy(0), pvz(0), svx(0), svy(0), svz(0), px(0), py(0), pz(0), pvxGen(0), pvyGen(0), pvzGen(0), svxGen(0), svyGen(0), svzGen(0), pxGen(0), pyGen(0), pzGen(0), nSingleDaughters(0), nCascadeDaughters(0), mcTrue(0), mcTrueVtx(0), mcPhysicalPrimary(0), hypNucDaughter(0) {} int pdgCode, isReconstructed, globalIndex; - uint8_t species; - bool isPrimaryCandidate, isMatter, passedEvSel, isMatterMC, passedEvSelMC, isPhysicalPrimary, collisionMcTrue; + uint8_t species, speciesMC; + bool isPrimaryCandidate, isMatter, isCascade, isCascadeMC, passedEvSel, isMatterMC, passedEvSelMC, isPhysicalPrimary, collisionMcTrue; float mass, y, pt, ct, yGen, ptGen, ctGen, cpaPvGen, cpaPv, cpaSv, maxDcaTracks, maxDcaTracksSv; float dcaToPvXY, dcaToPvZ, dcaToVtxXY, dcaToVtxZ, devToPvXY, chi2; float pvx, pvy, pvz, svx, svy, svz, px, py, pz; float pvxGen, pvyGen, pvzGen, svxGen, svyGen, svzGen, pxGen, pyGen, pzGen; - int nSingleDaughters, nCascadeDaughters, cent, occu; + int nSingleDaughters, nCascadeDaughters, cent, occu, runNumber; bool mcTrue, mcTrueVtx, mcPhysicalPrimary; std::vector daughterTracks; std::vector subDaughterMassVec; @@ -110,10 +110,12 @@ DECLARE_SOA_COLUMN(TvyGen, tvyGen, float); DECLARE_SOA_COLUMN(TvzGen, tvzGen, float); DECLARE_SOA_COLUMN(Centrality, centrality, int); DECLARE_SOA_COLUMN(Occupancy, occupancy, int); +DECLARE_SOA_COLUMN(RunNumber, runNumber, int); DECLARE_SOA_COLUMN(PassedEvSelMC, passedEvSelMC, bool); +DECLARE_SOA_COLUMN(SpeciesMC, speciesMC, int8_t); //! DECLARE_SOA_COLUMN(IsMatter, isMatter, bool); DECLARE_SOA_COLUMN(IsMatterGen, isMatterGen, bool); -DECLARE_SOA_COLUMN(IsReconstructed, isReconstructed, bool); +DECLARE_SOA_COLUMN(IsReconstructed, isReconstructed, int); DECLARE_SOA_COLUMN(CollMcTrue, collMcTrue, bool); DECLARE_SOA_COLUMN(D1X, d1X, float); DECLARE_SOA_COLUMN(D1Y, d1Y, float); @@ -263,15 +265,18 @@ DECLARE_SOA_COLUMN(Sd3IsPvContributor, sd3IsPvContributor, bool); DECLARE_SOA_COLUMN(Sd1sd2Mass, sd1sd2Mass, float); DECLARE_SOA_COLUMN(Sd1sd3Mass, sd1sd3Mass, float); DECLARE_SOA_COLUMN(Sd2sd3Mass, sd2sd3Mass, float); +DECLARE_SOA_COLUMN(D1sd1Mass, d1sd1Mass, float); +DECLARE_SOA_COLUMN(D1sd2Mass, d1sd2Mass, float); +DECLARE_SOA_COLUMN(D1sd3Mass, d1sd3Mass, float); } // namespace hypkftree -#define HYPKFGENBASE mcparticle::PdgCode, hypkftree::IsMatterGen, hypkftree::IsReconstructed, hykfmc::IsPhysicalPrimary, hypkftree::PassedEvSelMC, hypkftree::YGen, hypkftree::PtGen, hypkftree::CtGen +#define HYPKFGENBASE hypkftree::SpeciesMC, mcparticle::PdgCode, hypkftree::IsMatterGen, hypkftree::IsReconstructed, hykfmc::IsPhysicalPrimary, hypkftree::PassedEvSelMC, hypkftree::YGen, hypkftree::PtGen, hypkftree::CtGen #define HYPKFGENEXT hypkftree::CpaPvGen, hypkftree::PxGen, hypkftree::PyGen, hypkftree::PzGen, hypkftree::PvxGen, hypkftree::PvyGen, hypkftree::PvzGen, hypkftree::SvxGen, hypkftree::SvyGen, hypkftree::SvzGen #define HYPKFGENCAS hypkftree::TvxGen, hypkftree::TvyGen, hypkftree::TvzGen -#define HYPKFHYPNUC hykfmc::Species, hypkftree::IsMatter, hypkftree::Centrality, hypkftree::Occupancy, hykfmccoll::PassedEvSel, hykfhyp::Mass, hypkftree::Y, track::Pt, hypkftree::Ct, hypkftree::CosPa, hypkftree::DcaTracks, hypkftree::DcaTrackSv, hykfhyp::DcaToPvXY, hykfhyp::DcaToPvZ, hykfhyp::DevToPvXY, hykfhyp::Chi2, hypkftree::Pvx, hypkftree::Pvy, hypkftree::Pvz, hykfmc::Svx, hykfmc::Svy, hykfmc::Svz, hykfhyp::Px, hykfhyp::Py, hykfhyp::Pz, hypkftree::CollMcTrue +#define HYPKFHYPNUC hykfmc::Species, hypkftree::IsMatter, hypkftree::Centrality, hypkftree::Occupancy, hypkftree::RunNumber, hykfmccoll::PassedEvSel, hykfhyp::Mass, hypkftree::Y, track::Pt, hypkftree::Ct, hypkftree::CosPa, hypkftree::DcaTracks, hypkftree::DcaTrackSv, hykfhyp::DcaToPvXY, hykfhyp::DcaToPvZ, hykfhyp::DevToPvXY, hykfhyp::Chi2, hypkftree::Pvx, hypkftree::Pvy, hypkftree::Pvz, hykfmc::Svx, hykfmc::Svy, hykfmc::Svz, hykfhyp::Px, hykfhyp::Py, hykfhyp::Pz, hypkftree::CollMcTrue #define HYPKFHYPNUCMC hypkftree::McTrue, hykfmc::IsPhysicalPrimary @@ -291,6 +296,7 @@ DECLARE_SOA_COLUMN(Sd2sd3Mass, sd2sd3Mass, float); #define HYPKFSDMASS hypkftree::D1d2Mass, hypkftree::D1d3Mass, hypkftree::D2d3Mass #define HYPKFSSDMASS hypkftree::Sd1sd2Mass, hypkftree::Sd1sd3Mass, hypkftree::Sd2sd3Mass +#define HYPKFCSDMASS hypkftree::D1sd1Mass, hypkftree::D1sd2Mass, hypkftree::D1sd3Mass DECLARE_SOA_TABLE(HypKfGens, "AOD", "HYPKFGEN", HYPKFGENBASE); using HypKfGen = HypKfGens::iterator; @@ -307,10 +313,10 @@ using HypKfSingleThreeBodyCandidate = HypKfSingleThreeBodyCandidates::iterator; DECLARE_SOA_TABLE(HypKfMcSingleThreeBodyCandidates, "AOD", "HYPKFMCCAND3", HYPKFGENBASE, HYPKFGENEXT, HYPKFHYPNUC, HYPKFD1, HYPKFD2, HYPKFD3, HYPKFSDMASS); using HypKfMcSingleThreeBodyCandidate = HypKfMcSingleThreeBodyCandidates::iterator; -DECLARE_SOA_TABLE(HypKfCascadeTwoThreeCandidates, "AOD", "HYPKFCAND23", HYPKFHYPNUC, HYPKFHYPNUCMC, HYPKFD0, HYPKFD1, HYPKFSD1, HYPKFSD2, HYPKFSD3, HYPKFSSDMASS); +DECLARE_SOA_TABLE(HypKfCascadeTwoThreeCandidates, "AOD", "HYPKFCAND23", HYPKFHYPNUC, HYPKFHYPNUCMC, HYPKFD0, HYPKFD1, HYPKFSD1, HYPKFSD2, HYPKFSD3, HYPKFSSDMASS, HYPKFCSDMASS); using HypKfCascadeTwoThreeCandidate = HypKfCascadeTwoThreeCandidates::iterator; -DECLARE_SOA_TABLE(HypKfMcCascadeTwoThreeCandidates, "AOD", "HYPKFMCCAND23", HYPKFGENBASE, HYPKFGENEXT, HYPKFHYPNUC, HYPKFD0, HYPKFD1, HYPKFSD1, HYPKFSD2, HYPKFSD3, HYPKFSSDMASS); +DECLARE_SOA_TABLE(HypKfMcCascadeTwoThreeCandidates, "AOD", "HYPKFMCCAND23", HYPKFGENBASE, HYPKFGENEXT, HYPKFHYPNUC, HYPKFD0, HYPKFD1, HYPKFSD1, HYPKFSD2, HYPKFSD3, HYPKFSSDMASS, HYPKFCSDMASS); using HypKfMcCascadeTwoThreeCandidate = HypKfMcCascadeTwoThreeCandidates::iterator; DECLARE_SOA_TABLE(HypKfCascadeThreeTwoCandidates, "AOD", "HYPKFCAND32", HYPKFHYPNUC, HYPKFHYPNUCMC, HYPKFD0, HYPKFD1, HYPKFD2, HYPKFSDMASS, HYPKFSD1, HYPKFSD2); @@ -376,22 +382,22 @@ struct HypKfTreeCreator { { if (isMC && cfgMCGenerated) outputMcGenTable( - cand.pdgCode, cand.isMatterMC, cand.isReconstructed, cand.isPhysicalPrimary, cand.passedEvSelMC, cand.yGen, cand.ptGen, cand.ctGen); + cand.speciesMC, cand.pdgCode, cand.isMatterMC, cand.isReconstructed, cand.isPhysicalPrimary, cand.passedEvSelMC, cand.yGen, cand.ptGen, cand.ctGen); if (!cand.isReconstructed) { cand.daughterTracks.resize(4); cand.subDaughterMassVec.resize(4); hypDaughter.daughterTracks.resize(4); - hypDaughter.subDaughterMassVec.resize(4); + hypDaughter.subDaughterMassVec.resize(8); } - if (cfgNprimDaughters == 2 && cfgNsecDaughters == 0) { + if (cfgNprimDaughters == 2 && cfgNsecDaughters == 0) { // o2-linter: disable=magic-number (To be checked) const auto& d1 = cand.daughterTracks.at(0); const auto& d2 = cand.daughterTracks.at(1); if (!isMC || (isMC && cfgMCReconstructed && cand.isReconstructed)) outputTableTwo( - cand.species, cand.isMatter, cand.cent, cand.occu, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, cand.dcaToPvXY, - cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, cand.pz, cand.collisionMcTrue, + cand.species, cand.isMatter, cand.cent, cand.occu, cand.runNumber, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, + cand.dcaToPvXY, cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, cand.pz, cand.collisionMcTrue, cand.mcTrue, cand.mcPhysicalPrimary, d1.x, d1.y, d1.z, d1.px, d1.py, d1.pz, d1.tpcNcls, d1.tpcChi2, d1.itsNcls, d1.itsChi2, d1.itsMeanClsSizeL, d1.rigidity, d1.tpcSignal, d1.tpcNsigma, d1.tpcNsigmaNhp, d1.tpcNsigmaNlp, d1.tofMass, d1.dcaXY, d1.dcaZ, d1.isPvContributor, @@ -399,24 +405,24 @@ struct HypKfTreeCreator { d2.rigidity, d2.tpcSignal, d2.tpcNsigma, d2.tpcNsigmaNhp, d2.tpcNsigmaNlp, d2.tofMass, d2.dcaXY, d2.dcaZ, d2.isPvContributor); if (isMC && cfgMCCombined) outputTableMcTwo( - cand.pdgCode, cand.isMatterMC, cand.isReconstructed, cand.isPhysicalPrimary, cand.passedEvSelMC, cand.yGen, cand.ptGen, cand.ctGen, + cand.speciesMC, cand.pdgCode, cand.isMatterMC, cand.isReconstructed, cand.isPhysicalPrimary, cand.passedEvSelMC, cand.yGen, cand.ptGen, cand.ctGen, cand.cpaPvGen, cand.pxGen, cand.pyGen, cand.pzGen, cand.pvxGen, cand.pvyGen, cand.pvzGen, cand.svxGen, cand.svyGen, cand.svzGen, - cand.species, cand.isMatter, cand.cent, cand.occu, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, - cand.dcaToPvXY, cand.dcaToPvZ, cand.devToPvXY, + cand.species, cand.isMatter, cand.cent, cand.occu, cand.runNumber, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, + cand.maxDcaTracksSv, cand.dcaToPvXY, cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, cand.pz, cand.collisionMcTrue, d1.x, d1.y, d1.z, d1.px, d1.py, d1.pz, d1.tpcNcls, d1.tpcChi2, d1.itsNcls, d1.itsChi2, d1.itsMeanClsSizeL, d1.rigidity, d1.tpcSignal, d1.tpcNsigma, d1.tpcNsigmaNhp, d1.tpcNsigmaNlp, d1.tofMass, d1.dcaXY, d1.dcaZ, d1.isPvContributor, d2.x, d2.y, d2.z, d2.px, d2.py, d2.pz, d2.tpcNcls, d2.tpcChi2, d2.itsNcls, d2.itsChi2, d2.itsMeanClsSizeL, d2.rigidity, d2.tpcSignal, d2.tpcNsigma, d2.tpcNsigmaNhp, d2.tpcNsigmaNlp, d2.tofMass, d2.dcaXY, d2.dcaZ, d2.isPvContributor); } - if (cand.isPrimaryCandidate && ((cfgNprimDaughters == 3 && cfgNsecDaughters == 0) || (cfgNsecDaughters == 3 && cfgSpecies == 0))) { + if (((!isMC && cand.isPrimaryCandidate) || (isMC && cand.isPhysicalPrimary)) && ((cfgNprimDaughters == 3 && cfgNsecDaughters == 0) || (cfgNsecDaughters == 3 && cfgSpecies == 0))) { // o2-linter: disable=magic-number (To be checked) const auto& d1 = cand.daughterTracks.at(0); const auto& d2 = cand.daughterTracks.at(1); const auto& d3 = cand.daughterTracks.at(2); if (!isMC || (isMC && cfgMCReconstructed && cand.isReconstructed)) outputTableThree( - cand.species, cand.isMatter, cand.cent, cand.occu, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, cand.dcaToPvXY, - cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, cand.pz, cand.collisionMcTrue, + cand.species, cand.isMatter, cand.cent, cand.occu, cand.runNumber, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, + cand.dcaToPvXY, cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, cand.pz, cand.collisionMcTrue, cand.mcTrue, cand.mcPhysicalPrimary, d1.x, d1.y, d1.z, d1.px, d1.py, d1.pz, d1.tpcNcls, d1.tpcChi2, d1.itsNcls, d1.itsChi2, d1.itsMeanClsSizeL, d1.rigidity, d1.tpcSignal, d1.tpcNsigma, d1.tpcNsigmaNhp, d1.tpcNsigmaNlp, d1.tofMass, d1.dcaXY, d1.dcaZ, d1.isPvContributor, @@ -427,10 +433,10 @@ struct HypKfTreeCreator { d1.subMass, d2.subMass, d3.subMass); if (isMC && cfgMCCombined) outputTableMcThree( - cand.pdgCode, cand.isMatterMC, cand.isReconstructed, cand.isPhysicalPrimary, cand.passedEvSelMC, cand.yGen, cand.ptGen, cand.ctGen, + cand.speciesMC, cand.pdgCode, cand.isMatterMC, cand.isReconstructed, cand.isPhysicalPrimary, cand.passedEvSelMC, cand.yGen, cand.ptGen, cand.ctGen, cand.cpaPvGen, cand.pxGen, cand.pyGen, cand.pzGen, cand.pvxGen, cand.pvyGen, cand.pvzGen, cand.svxGen, cand.svyGen, cand.svzGen, - cand.species, cand.isMatter, cand.cent, cand.occu, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, - cand.dcaToPvXY, cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, + cand.species, cand.isMatter, cand.cent, cand.occu, cand.runNumber, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, + cand.maxDcaTracksSv, cand.dcaToPvXY, cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, cand.pz, cand.collisionMcTrue, d1.x, d1.y, d1.z, d1.px, d1.py, d1.pz, d1.tpcNcls, d1.tpcChi2, d1.itsNcls, d1.itsChi2, d1.itsMeanClsSizeL, d1.rigidity, d1.tpcSignal, d1.tpcNsigma, d1.tpcNsigmaNhp, d1.tpcNsigmaNlp, d1.tofMass, d1.dcaXY, d1.dcaZ, d1.isPvContributor, @@ -440,7 +446,9 @@ struct HypKfTreeCreator { d3.rigidity, d3.tpcSignal, d3.tpcNsigma, d3.tpcNsigmaNhp, d3.tpcNsigmaNlp, d3.tofMass, d3.dcaXY, d3.dcaZ, d3.isPvContributor, d1.subMass, d2.subMass, d3.subMass); } - if (cfgNprimDaughters == 2 && cfgNsecDaughters == 3) { + if ((!isMC && !cand.isCascade) || (isMC && !cand.isCascadeMC)) + return; + if (cfgNprimDaughters == 2 && cfgNsecDaughters == 3) { // o2-linter: disable=magic-number (To be checked) const auto& d0 = cand.daughterTracks.at(0); const auto& d1 = cand.daughterTracks.at(1); const auto& sd1 = hypDaughter.daughterTracks.at(0); @@ -448,8 +456,8 @@ struct HypKfTreeCreator { const auto& sd3 = hypDaughter.daughterTracks.at(2); if (!isMC || (isMC && cfgMCReconstructed && cand.isReconstructed)) outputTableTwoThree( - cand.species, cand.isMatter, cand.cent, cand.occu, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, cand.dcaToPvXY, - cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, cand.pz, cand.collisionMcTrue, + cand.species, cand.isMatter, cand.cent, cand.occu, cand.runNumber, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, + cand.dcaToPvXY, cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, cand.pz, cand.collisionMcTrue, cand.mcTrue, cand.mcPhysicalPrimary, hypDaughter.svx, hypDaughter.svy, hypDaughter.svz, d0.x, d0.y, d0.z, d0.px, d0.py, d0.pz, hypDaughter.mass, hypDaughter.ct, hypDaughter.cpaPv, hypDaughter.maxDcaTracks, hypDaughter.dcaToPvXY, hypDaughter.dcaToPvZ, hypDaughter.dcaToVtxXY, hypDaughter.dcaToVtxZ, hypDaughter.chi2, @@ -461,13 +469,13 @@ struct HypKfTreeCreator { sd2.rigidity, sd2.tpcSignal, sd2.tpcNsigma, sd2.tpcNsigmaNhp, sd2.tpcNsigmaNlp, sd2.tofMass, sd2.dcaXY, sd2.dcaZ, sd2.isPvContributor, sd3.x, sd3.y, sd3.z, sd3.px, sd3.py, sd3.pz, sd3.tpcNcls, sd3.tpcChi2, sd3.itsNcls, sd3.itsChi2, sd3.itsMeanClsSizeL, sd3.rigidity, sd3.tpcSignal, sd3.tpcNsigma, sd3.tpcNsigmaNhp, sd3.tpcNsigmaNlp, sd3.tofMass, sd3.dcaXY, sd3.dcaZ, sd3.isPvContributor, - sd1.subMass, sd2.subMass, sd3.subMass); + sd1.subMass, sd2.subMass, sd3.subMass, cand.subDaughterMassVec.at(0), cand.subDaughterMassVec.at(1), cand.subDaughterMassVec.at(2)); if (isMC && cfgMCCombined) outputTableMcTwoThree( - cand.pdgCode, cand.isMatterMC, cand.isReconstructed, cand.isPhysicalPrimary, cand.passedEvSelMC, cand.yGen, cand.ptGen, cand.ctGen, + cand.speciesMC, cand.pdgCode, cand.isMatterMC, cand.isReconstructed, cand.isPhysicalPrimary, cand.passedEvSelMC, cand.yGen, cand.ptGen, cand.ctGen, cand.cpaPvGen, cand.pxGen, cand.pyGen, cand.pzGen, cand.pvxGen, cand.pvyGen, cand.pvzGen, cand.svxGen, cand.svyGen, cand.svzGen, - cand.species, cand.isMatter, cand.cent, cand.occu, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, cand.dcaToPvXY, - cand.dcaToPvZ, cand.devToPvXY, + cand.species, cand.isMatter, cand.cent, cand.occu, cand.runNumber, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, + cand.dcaToPvXY, cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, cand.pz, cand.collisionMcTrue, hypDaughter.svx, hypDaughter.svy, hypDaughter.svz, d0.x, d0.y, d0.z, d0.px, d0.py, d0.pz, hypDaughter.mass, hypDaughter.ct, hypDaughter.cpaPv, hypDaughter.maxDcaTracks, hypDaughter.dcaToPvXY, hypDaughter.dcaToPvZ, hypDaughter.dcaToVtxXY, hypDaughter.dcaToVtxZ, hypDaughter.chi2, @@ -479,9 +487,9 @@ struct HypKfTreeCreator { sd2.rigidity, sd2.tpcSignal, sd2.tpcNsigma, sd2.tpcNsigmaNhp, sd2.tpcNsigmaNlp, sd2.tofMass, sd2.dcaXY, sd2.dcaZ, sd2.isPvContributor, sd3.x, sd3.y, sd3.z, sd3.px, sd3.py, sd3.pz, sd3.tpcNcls, sd3.tpcChi2, sd3.itsNcls, sd3.itsChi2, sd3.itsMeanClsSizeL, sd3.rigidity, sd3.tpcSignal, sd3.tpcNsigma, sd3.tpcNsigmaNhp, sd3.tpcNsigmaNlp, sd3.tofMass, sd3.dcaXY, sd3.dcaZ, sd3.isPvContributor, - sd1.subMass, sd2.subMass, sd3.subMass); + sd1.subMass, sd2.subMass, sd3.subMass, cand.subDaughterMassVec.at(0), cand.subDaughterMassVec.at(1), cand.subDaughterMassVec.at(2)); } - if (cfgNprimDaughters == 3 && cfgNsecDaughters == 1) { + if (cfgNprimDaughters == 3 && cfgNsecDaughters == 1) { // o2-linter: disable=magic-number (To be checked) const auto& d0 = cand.daughterTracks.at(0); const auto& d1 = cand.daughterTracks.at(1); const auto& d2 = cand.daughterTracks.at(2); @@ -489,8 +497,8 @@ struct HypKfTreeCreator { const auto& sd2 = hypDaughter.daughterTracks.at(1); if (!isMC || (isMC && cfgMCReconstructed && cand.isReconstructed)) outputTableThreeTwo( - cand.species, cand.isMatter, cand.cent, cand.occu, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, cand.dcaToPvXY, - cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, cand.pz, cand.collisionMcTrue, + cand.species, cand.isMatter, cand.cent, cand.occu, cand.runNumber, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, + cand.dcaToPvXY, cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, cand.pz, cand.collisionMcTrue, cand.mcTrue, cand.mcPhysicalPrimary, hypDaughter.svx, hypDaughter.svy, hypDaughter.svz, d0.x, d0.y, d0.z, d0.px, d0.py, d0.pz, hypDaughter.mass, hypDaughter.ct, hypDaughter.cpaPv, hypDaughter.maxDcaTracks, hypDaughter.dcaToPvXY, hypDaughter.dcaToPvZ, hypDaughter.dcaToVtxXY, hypDaughter.dcaToVtxZ, hypDaughter.chi2, d1.x, d1.y, d1.z, d1.px, d1.py, d1.pz, d1.tpcNcls, d1.tpcChi2, d1.itsNcls, d1.itsChi2, d1.itsMeanClsSizeL, @@ -504,10 +512,10 @@ struct HypKfTreeCreator { sd2.rigidity, sd2.tpcSignal, sd2.tpcNsigma, sd2.tpcNsigmaNhp, sd2.tpcNsigmaNlp, sd2.tofMass, sd2.dcaXY, sd2.dcaZ, sd2.isPvContributor); if (isMC && cfgMCCombined) outputTableMcThreeTwo( - cand.pdgCode, cand.isMatterMC, cand.isReconstructed, cand.isPhysicalPrimary, cand.passedEvSelMC, cand.yGen, cand.ptGen, cand.ctGen, + cand.speciesMC, cand.pdgCode, cand.isMatterMC, cand.isReconstructed, cand.isPhysicalPrimary, cand.passedEvSelMC, cand.yGen, cand.ptGen, cand.ctGen, cand.cpaPvGen, cand.pxGen, cand.pyGen, cand.pzGen, cand.pvxGen, cand.pvyGen, cand.pvzGen, cand.svxGen, cand.svyGen, cand.svzGen, - cand.species, cand.isMatter, cand.cent, cand.occu, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, cand.dcaToPvXY, - cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, cand.pz, cand.collisionMcTrue, + cand.species, cand.isMatter, cand.cent, cand.occu, cand.runNumber, cand.passedEvSel, cand.mass, cand.y, cand.pt, cand.ct, cand.cpaPv, cand.maxDcaTracks, cand.maxDcaTracksSv, + cand.dcaToPvXY, cand.dcaToPvZ, cand.devToPvXY, cand.chi2, cand.pvx, cand.pvy, cand.pvz, cand.svx, cand.svy, cand.svz, cand.px, cand.py, cand.pz, cand.collisionMcTrue, hypDaughter.svx, hypDaughter.svy, hypDaughter.svz, d0.x, d0.y, d0.z, d0.px, d0.py, d0.pz, hypDaughter.mass, hypDaughter.ct, hypDaughter.cpaPv, hypDaughter.maxDcaTracks, hypDaughter.dcaToPvXY, hypDaughter.dcaToPvZ, hypDaughter.dcaToVtxXY, hypDaughter.dcaToVtxZ, hypDaughter.chi2, d1.x, d1.y, d1.z, d1.px, d1.py, d1.pz, d1.tpcNcls, d1.tpcChi2, d1.itsNcls, d1.itsChi2, d1.itsMeanClsSizeL, @@ -547,8 +555,11 @@ struct HypKfTreeCreator { cand.species = std::abs(hypNuc.species()); cand.isPrimaryCandidate = hypNuc.primary(); cand.isMatter = hypNuc.isMatter(); + cand.mcTrue = hypNuc.mcTrue(); + cand.isCascade = std::abs(cand.species) > 10; // o2-linter: disable=magic-number (To be checked) cand.cent = coll.centFT0C(); cand.occu = coll.occupancy(); + cand.runNumber = coll.runNumber(); cand.passedEvSel = coll.passedEvSel(); cand.mass = hypNuc.mass(); cand.y = hypNuc.y(); @@ -601,11 +612,19 @@ struct HypKfTreeCreator { cand.daughterTracks.at(trackCount).z = addOn.z(); cand.daughterTracks.at(trackCount).px = addOn.px(); cand.daughterTracks.at(trackCount).py = addOn.py(); - cand.daughterTracks.at(trackCount).pz = addOn.py(); + cand.daughterTracks.at(trackCount).pz = addOn.pz(); trackCount++; } + + if (cand.isCascade) { + auto subDaughters = hypNuc.hypKfSubD_as(); + for (const auto& subDaughter : subDaughters) { + cand.subDaughterMassVec.push_back(subDaughter.subMass()); + } + } + cand.nSingleDaughters = trackCount; - if (cand.nSingleDaughters < 3) + if (cand.nSingleDaughters < 3) // o2-linter: disable=magic-number (To be checked) return; trackCount = 0; @@ -626,9 +645,12 @@ struct HypKfTreeCreator { const auto mcParticleIdx = mcHypNuc.globalIndex(); auto hypNucsByMc = hypNucs.sliceBy(perMcParticle, mcParticleIdx); HyperNucleus candidate, hypNucDaughter; + candidate.speciesMC = mcHypNuc.species(); + candidate.isCascadeMC = std::abs(candidate.speciesMC) > 10; // o2-linter: disable=magic-number (To be checked) candidate.pdgCode = mcHypNuc.pdgCode(); candidate.isMatterMC = mcHypNuc.isMatter(); candidate.isPhysicalPrimary = mcHypNuc.isPhysicalPrimary(); + candidate.mcPhysicalPrimary = mcHypNuc.isPhysicalPrimary(); candidate.passedEvSelMC = mcColl.passedEvSel(); candidate.yGen = mcHypNuc.y(); candidate.ptGen = mcHypNuc.pt(); From fecdb9ae814f12970015840075d32527a269f8e2 Mon Sep 17 00:00:00 2001 From: mhartung71 <50153519+mhartung71@users.noreply.github.com> Date: Thu, 29 May 2025 13:55:45 +0200 Subject: [PATCH 4/5] Update hypKfRecoTask.cxx --- PWGLF/TableProducer/Nuspex/hypKfRecoTask.cxx | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/PWGLF/TableProducer/Nuspex/hypKfRecoTask.cxx b/PWGLF/TableProducer/Nuspex/hypKfRecoTask.cxx index 9bef4359979..d06405f1a32 100644 --- a/PWGLF/TableProducer/Nuspex/hypKfRecoTask.cxx +++ b/PWGLF/TableProducer/Nuspex/hypKfRecoTask.cxx @@ -16,6 +16,7 @@ #include #include #include +#include #include "Framework/runDataProcessing.h" #include "Framework/AnalysisTask.h" #include "Framework/AnalysisDataModel.h" @@ -934,7 +935,7 @@ struct HypKfRecoTask { std::vector*> hypNucVectors = {&singleHyperNuclei, &cascadeHyperNuclei}; std::vector>*> candidateVectors = {&singleHyperNucCandidates, &cascadeHyperNucCandidates}; - for (int vec = 0; vec < candidateVectors.size(); vec++) { + for (unsigned int vec = 0; vec < candidateVectors.size(); vec++) { auto candidateVector = candidateVectors.at(vec); for (size_t hyperNucIter = 0; hyperNucIter < hypNucVectors.at(vec)->size(); hyperNucIter++) { HyperNucleus* hyperNuc = &(hypNucVectors.at(vec)->at(hyperNucIter)); @@ -1017,7 +1018,7 @@ struct HypKfRecoTask { for (const auto& mcPart : particlesMC) { if (!mcCollInfos.at(mcPart.mcCollisionId()).passedEvSel) continue; - for (int vec = 0; vec < hypNucVectors.size(); vec++) { + for (unsigned int vec = 0; vec < hypNucVectors.size(); vec++) { for (size_t hyperNucIter = 0; hyperNucIter < hypNucVectors.at(vec)->size(); hyperNucIter++) { HyperNucleus* hyperNuc = &(hypNucVectors.at(vec)->at(hyperNucIter)); if (!hyperNuc->active) From 18ec6fd28842a59763ff92aacafe3bdd519ec629 Mon Sep 17 00:00:00 2001 From: mhartung71 <50153519+mhartung71@users.noreply.github.com> Date: Thu, 29 May 2025 13:57:08 +0200 Subject: [PATCH 5/5] Update hypKfTreeCreator.cxx --- PWGLF/TableProducer/Nuspex/hypKfTreeCreator.cxx | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/PWGLF/TableProducer/Nuspex/hypKfTreeCreator.cxx b/PWGLF/TableProducer/Nuspex/hypKfTreeCreator.cxx index 81a7569ed27..246a67f57d6 100644 --- a/PWGLF/TableProducer/Nuspex/hypKfTreeCreator.cxx +++ b/PWGLF/TableProducer/Nuspex/hypKfTreeCreator.cxx @@ -556,7 +556,7 @@ struct HypKfTreeCreator { cand.isPrimaryCandidate = hypNuc.primary(); cand.isMatter = hypNuc.isMatter(); cand.mcTrue = hypNuc.mcTrue(); - cand.isCascade = std::abs(cand.species) > 10; // o2-linter: disable=magic-number (To be checked) + cand.isCascade = cand.species > 10; // o2-linter: disable=magic-number (To be checked) cand.cent = coll.centFT0C(); cand.occu = coll.occupancy(); cand.runNumber = coll.runNumber(); @@ -646,7 +646,7 @@ struct HypKfTreeCreator { auto hypNucsByMc = hypNucs.sliceBy(perMcParticle, mcParticleIdx); HyperNucleus candidate, hypNucDaughter; candidate.speciesMC = mcHypNuc.species(); - candidate.isCascadeMC = std::abs(candidate.speciesMC) > 10; // o2-linter: disable=magic-number (To be checked) + candidate.isCascadeMC = candidate.speciesMC > 10; // o2-linter: disable=magic-number (To be checked) candidate.pdgCode = mcHypNuc.pdgCode(); candidate.isMatterMC = mcHypNuc.isMatter(); candidate.isPhysicalPrimary = mcHypNuc.isPhysicalPrimary();