From 4bca7e2be5b68a983c67c8994380d9fcdf533560 Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Mon, 13 Jan 2025 21:50:17 +0100 Subject: [PATCH 01/19] Track-nucleus combinations correction --- .../femtoUniversePairTaskTrackNucleus.cxx | 417 +++++++++--------- 1 file changed, 209 insertions(+), 208 deletions(-) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx index 156d6f43292..8927326aadf 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx @@ -20,6 +20,7 @@ #include #include +#include #include "Framework/AnalysisTask.h" #include "Framework/runDataProcessing.h" #include "Framework/HistogramRegistry.h" @@ -49,14 +50,14 @@ using namespace o2::soa; namespace { -static constexpr int nPart = 2; -static constexpr int nCuts = 5; +static constexpr int Npart = 2; +static constexpr int Ncuts = 5; static const std::vector partNames{"Track", "Nucleus"}; static const std::vector cutNames{"MaxPt", "PIDthr", "nSigmaTPC", "nSigmaTPCTOF", "MaxP"}; -static const float cutsTable[nPart][nCuts]{{4.05f, 1.f, 3.f, 3.f, 100.f}, {4.05f, 1.f, 3.f, 3.f, 100.f}}; +static const float cutsTable[Npart][Ncuts]{{4.05f, 1.f, 3.f, 3.f, 100.f}, {4.05f, 1.f, 3.f, 3.f, 100.f}}; } // namespace -struct femtoUniversePairTaskTrackNucleus { +struct FemtoUniversePairTaskTrackNucleus { Service pdg; @@ -64,50 +65,50 @@ struct femtoUniversePairTaskTrackNucleus { /// Table for both objects struct : o2::framework::ConfigurableGroup { - Configurable IsKaonNsigma{"IsKaonNsigma", false, "Enable a strict cut selection for K+ and K-"}; - Configurable ConfNsigmaCombined{"ConfNsigmaCombined", 3.0f, "TPC and TOF Pion Sigma (combined) for momentum > ConfTOFpMin"}; - Configurable ConfNsigmaTPC{"ConfNsigmaTPC", 3.0f, "TPC Pion Sigma for momentum < ConfTOFpMin"}; - Configurable ConfTOFpMin{"ConfTOFpMin", 0.45f, "Min. momentum for which TOF is required for PID."}; - Configurable ConfEtaMax{"ConfEtaMax", 0.8f, "Higher limit for |Eta| (the same for both particles)"}; - - Configurable> ConfCutTable{"ConfCutTable", {cutsTable[0], nPart, nCuts, partNames, cutNames}, "Particle selections"}; - Configurable ConfNspecies{"ConfNspecies", 2, "Number of particle spieces with PID info"}; - Configurable ConfIsMC{"ConfIsMC", false, "Enable additional Histogramms in the case of a MonteCarlo Run"}; - Configurable> ConfTrkPIDnSigmaMax{"ConfTrkPIDnSigmaMax", std::vector{4.f, 3.f, 2.f}, "This configurable needs to be the same as the one used in the producer task"}; - Configurable ConfUse3D{"ConfUse3D", false, "Enable three dimensional histogramms (to be used only for analysis with high statistics): k* vs mT vs multiplicity"}; + Configurable isKaonNsigma{"isKaonNsigma", false, "Enable a strict cut selection for K+ and K-"}; + Configurable confNsigmaCombined{"confNsigmaCombined", 3.0f, "TPC and TOF Pion Sigma (combined) for momentum > confTOFpMin"}; + Configurable confNsigmaTPC{"confNsigmaTPC", 3.0f, "TPC Pion Sigma for momentum < confTOFpMin"}; + Configurable confTOFpMin{"confTOFpMin", 0.45f, "Min. momentum for which TOF is required for PID."}; + Configurable confEtaMax{"confEtaMax", 0.8f, "Higher limit for |Eta| (the same for both particles)"}; + + Configurable> confCutTable{"confCutTable", {cutsTable[0], Npart, Ncuts, partNames, cutNames}, "Particle selections"}; + Configurable confNspecies{"confNspecies", 2, "Number of particle spieces with PID info"}; + Configurable confIsMC{"confIsMC", false, "Enable additional Histogramms in the case of a MonteCarlo Run"}; + Configurable> confTrkPIDNsigmaMax{"confTrkPIDNsigmaMax", std::vector{4.f, 3.f, 2.f}, "This configurable needs to be the same as the one used in the producer task"}; + Configurable confUse3D{"confUse3D", false, "Enable three dimensional histogramms (to be used only for analysis with high statistics): k* vs mT vs multiplicity"}; } twoobjectsconfigs; /// Table for separate deuteron configurables struct : o2::framework::ConfigurableGroup { - Configurable ConfNsigmaTPCDe{"ConfNsigmaTPCDe", 2.0f, "TPC Deuteron Sigma for momentum < ConfTOFpMinDe"}; - Configurable ConfNsigmaTOFDe{"ConfNsigmaTOFDe", 2.0f, "TOF Deuteron Sigma"}; - Configurable ConfTOFpMinDe{"ConfTOFpMinDe", 1.0f, "Min. momentum for deuterons for which TOF is required for PID"}; - Configurable ConfPLowDe{"ConfPLowDe", 0.8f, "Lower limit for momentum for deuterons"}; - Configurable ConfPHighDe{"ConfPHighDe", 1.8f, "Higher limit for momentum for deuterons"}; + Configurable confNsigmaTPCDe{"confNsigmaTPCDe", 2.0f, "TPC Deuteron Sigma for momentum < confTOFpMinDe"}; + Configurable confNsigmaTOFDe{"confNsigmaTOFDe", 2.0f, "TOF Deuteron Sigma"}; + Configurable confTOFpMinDe{"confTOFpMinDe", 1.0f, "Min. momentum for deuterons for which TOF is required for PID"}; + Configurable confPlowDe{"confPlowDe", 0.8f, "Lower limit for momentum for deuterons"}; + Configurable confPhighDe{"confPhighDe", 1.8f, "Higher limit for momentum for deuterons"}; } deuteronconfigs; /// Table for linear cut for TPC Deuteron Sigma struct : o2::framework::ConfigurableGroup { - Configurable ConfIsLine{"ConfIsLine", false, "Enable a separation line for clearer TPC Deuteron Sigma"}; - Configurable LinCutpLow{"LinCutpLow", 0.0f, "Lower limit of momentum for linear cut of TPC Deuteron Sigma"}; - Configurable LinCutpHigh{"LinCutpHigh", 1.4f, "Higher limit of momentum for linear cut of TPC Deuteron Sigma"}; - Configurable LinCutParA{"LinCutParA", -167.0f, "Parameter 'A' of a linear function 'y = A * x + B'"}; - Configurable LinCutParB{"LinCutParB", 300.0f, "Parameter 'B' of a linear function 'y = A * x + B'"}; + Configurable confIsLine{"confIsLine", false, "Enable a separation line for clearer TPC Deuteron Sigma"}; + Configurable linCutPlow{"linCutPlow", 0.0f, "Lower limit of momentum for linear cut of TPC Deuteron Sigma"}; + Configurable linCutPhigh{"linCutPhigh", 1.4f, "Higher limit of momentum for linear cut of TPC Deuteron Sigma"}; + Configurable linCutParA{"linCutParA", -167.0f, "Parameter 'A' of a linear function 'y = A * x + B'"}; + Configurable linCutParB{"linCutParB", 300.0f, "Parameter 'B' of a linear function 'y = A * x + B'"}; } lincut; /// Table for polynomial 3 cut for TPC Deuteron Sigma struct : o2::framework::ConfigurableGroup { - Configurable ConfIsPol{"ConfIsPol", false, "Enable a separation polynomial 3 curve for clearer TPC Deuteron Sigma"}; - Configurable PolCutParA{"PolCutParA", -52.2f, "Parameter 'A' of a polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; - Configurable PolCutParB{"PolCutParB", 357.7f, "Parameter 'B' of a polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; - Configurable PolCutParC{"PolCutParC", -834.7f, "Parameter 'C' of a polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; - Configurable PolCutParD{"PolCutParD", 705.8f, "Parameter 'D' of a polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; + Configurable confIsPol3{"confIsPol3", false, "Enable a separation polynomial 3 curve for clearer TPC Deuteron Sigma"}; + Configurable polCutParA{"polCutParA", -52.2f, "Parameter 'A' of a polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; + Configurable polCutParB{"polCutParB", 357.7f, "Parameter 'B' of a polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; + Configurable polCutParC{"polCutParC", -834.7f, "Parameter 'C' of a polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; + Configurable polCutParD{"polCutParD", 705.8f, "Parameter 'D' of a polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; } polcut; using FemtoFullParticles = soa::Join; // Filters for selection - Filter trackAdditionalfilter = (nabs(aod::femtouniverseparticle::eta) < twoobjectsconfigs.ConfEtaMax); // example filtering on configurable + Filter trackAdditionalfilter = (nabs(aod::femtouniverseparticle::eta) < twoobjectsconfigs.confEtaMax); // example filtering on configurable using FilteredFemtoFullParticles = soa::Filtered; // using FilteredFemtoFullParticles = FemtoFullParticles; //if no filtering is applied uncomment this option @@ -116,36 +117,36 @@ struct femtoUniversePairTaskTrackNucleus { /// Track struct : o2::framework::ConfigurableGroup { - Configurable ConfPDGCodeTrack{"ConfPDGCodeTrack", 321, "Track -- PDG code"}; - // Configurable ConfCutTrack{"ConfCutTrack", 5542474, "Track -- Selection bit from cutCulator"}; - Configurable ConfPIDTrack{"ConfPIDTrack", 3, "Track -- Read from cutCulator"}; // we also need the possibility to specify whether the bit is true/false ->std>>vector>int>> - Configurable ConfPtLowTrack{"ConfPtLowTrack", 0.14, "Lower pT limit for track"}; - Configurable ConfPtHighTrack{"ConfPtHighTrack", 1.5, "Higher pT limit for track"}; - Configurable ConfChargeTrack{"ConfChargeTrack", 1, "Track sign"}; // -1 means anti-particle + Configurable confPDGCodeTrack{"confPDGCodeTrack", 321, "Track -- PDG code"}; + // Configurable confCutTrack{"confCutTrack", 5542474, "Track -- Selection bit from cutCulator"}; + Configurable confPIDTrack{"confPIDTrack", 3, "Track -- Read from cutCulator"}; // we also need the possibility to specify whether the bit is true/false ->std>>vector>int>> + Configurable confPtLowTrack{"confPtLowTrack", 0.14, "Lower pT limit for track"}; + Configurable confPtHighTrack{"confPtHighTrack", 1.5, "Higher pT limit for track"}; + Configurable confChargeTrack{"confChargeTrack", 1, "Track sign"}; // -1 means anti-particle } trackfilter; /// Partition for track - Partition partTrack = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == trackfilter.ConfChargeTrack && aod::femtouniverseparticle::pt < trackfilter.ConfPtHighTrack && aod::femtouniverseparticle::pt > trackfilter.ConfPtLowTrack; + Partition partTrack = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == trackfilter.confChargeTrack && aod::femtouniverseparticle::pt < trackfilter.confPtHighTrack && aod::femtouniverseparticle::pt > trackfilter.confPtLowTrack; - Partition> partTrackMC = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == trackfilter.ConfChargeTrack && aod::femtouniverseparticle::pt < trackfilter.ConfPtHighTrack && aod::femtouniverseparticle::pt > trackfilter.ConfPtLowTrack; + Partition> partTrackMC = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && aod::femtouniverseparticle::sign == trackfilter.confChargeTrack && aod::femtouniverseparticle::pt < trackfilter.confPtHighTrack && aod::femtouniverseparticle::pt > trackfilter.confPtLowTrack; /// Histogramming for track FemtoUniverseParticleHisto trackHistoTrack; /// Nucleus struct : o2::framework::ConfigurableGroup { - Configurable ConfPDGCodeNucleus{"ConfPDGCodeNucleus", 1000010020, "Nucleus -- PDG code"}; - // Configurable ConfCutNucleus{"ConfCutNucleus", 5542474, "Nucleus -- Selection bit"}; - Configurable ConfPIDNucleus{"ConfPIDNucleus", 5, "Nucleus -- Read from cutCulator"}; // we also need the possibility to specify whether the bit is true/false ->std>>vector> - Configurable ConfPtLowNucleus{"ConfPtLowNucleus", 0, "Lower pT limit for nucleus"}; - Configurable ConfPtHighNucleus{"ConfPtHighNucleus", 5, "Higher pT limit for nucleus"}; - Configurable ConfChargeNucleus{"ConfChargeNucleus", 1, "Nucleus sign"}; // -1 means anti-nucleus + Configurable confPDGCodeNucleus{"confPDGCodeNucleus", 1000010020, "Nucleus -- PDG code"}; + // Configurable confCutNucleus{"confCutNucleus", 5542474, "Nucleus -- Selection bit"}; + Configurable confPIDNucleus{"confPIDNucleus", 5, "Nucleus -- Read from cutCulator"}; // we also need the possibility to specify whether the bit is true/false ->std>>vector> + Configurable confPtLowNucleus{"confPtLowNucleus", 0, "Lower pT limit for nucleus"}; + Configurable confPtHighNucleus{"confPtHighNucleus", 5, "Higher pT limit for nucleus"}; + Configurable confChargeNucleus{"confChargeNucleus", 1, "Nucleus sign"}; // -1 means anti-nucleus } nucleusfilter; /// Partition for nucleus - Partition partNucleus = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && (aod::femtouniverseparticle::sign == nucleusfilter.ConfChargeNucleus) && aod::femtouniverseparticle::pt < nucleusfilter.ConfPtHighNucleus && aod::femtouniverseparticle::pt > nucleusfilter.ConfPtLowNucleus; + Partition partNucleus = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && (aod::femtouniverseparticle::sign == nucleusfilter.confChargeNucleus) && aod::femtouniverseparticle::pt < nucleusfilter.confPtHighNucleus && aod::femtouniverseparticle::pt > nucleusfilter.confPtLowNucleus; - Partition> partNucleusMC = aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack) && (aod::femtouniverseparticle::sign == nucleusfilter.ConfChargeNucleus) && aod::femtouniverseparticle::pt < nucleusfilter.ConfPtHighNucleus && aod::femtouniverseparticle::pt > nucleusfilter.ConfPtLowNucleus; + Partition> partNucleusMC = aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack) && (aod::femtouniverseparticle::sign == nucleusfilter.confChargeNucleus) && aod::femtouniverseparticle::pt < nucleusfilter.confPtHighNucleus && aod::femtouniverseparticle::pt > nucleusfilter.confPtLowNucleus; /// Histogramming for nucleus FemtoUniverseParticleHisto trackHistoNucleus; @@ -158,41 +159,41 @@ struct femtoUniversePairTaskTrackNucleus { std::vector kNsigma; /// Event part - Configurable ConfV0MLow{"ConfV0MLow", 0.0, "Lower limit for V0M multiplicity"}; - Configurable ConfV0MHigh{"ConfV0MHigh", 25000.0, "Upper limit for V0M multiplicity"}; - Configurable ConfSphericityCutMin{"ConfSphericityCutMin", 0, "Min. sphericity"}; - Configurable ConfSphericityCutMax{"ConfSphericityCutMax", 3, "Max. sphericity"}; + Configurable confV0Mlow{"confV0Mlow", 0.0, "Lower limit for V0M multiplicity"}; + Configurable confV0Mhigh{"confV0Mhigh", 25000.0, "Upper limit for V0M multiplicity"}; + Configurable confSphericityCutMin{"confSphericityCutMin", 0, "Min. sphericity"}; + Configurable confSphericityCutMax{"confSphericityCutMax", 3, "Max. sphericity"}; - Filter collV0Mfilter = ((o2::aod::femtouniversecollision::multV0M > ConfV0MLow) && (o2::aod::femtouniversecollision::multV0M < ConfV0MHigh)); - Filter colSpherfilter = ((o2::aod::femtouniversecollision::sphericity > ConfSphericityCutMin) && (o2::aod::femtouniversecollision::sphericity < ConfSphericityCutMax)); - // Filter trackAdditionalfilter = (nabs(aod::femtouniverseparticle::eta) < twoobjectsconfigs.ConfEtaMax); // example filtering on configurable + Filter collV0Mfilter = ((o2::aod::femtouniversecollision::multV0M > confV0Mlow) && (o2::aod::femtouniversecollision::multV0M < confV0Mhigh)); + Filter colSpherfilter = ((o2::aod::femtouniversecollision::sphericity > confSphericityCutMin) && (o2::aod::femtouniversecollision::sphericity < confSphericityCutMax)); + // Filter trackAdditionalfilter = (nabs(aod::femtouniverseparticle::eta) < twoobjectsconfigs.confEtaMax); // example filtering on configurable /// Particle part - ConfigurableAxis ConfTempFitVarBins{"ConfDTempFitVarBins", {300, -0.15, 0.15}, "binning of the TempFitVar in the pT vs. TempFitVar plot"}; - ConfigurableAxis ConfTempFitVarpTBins{"ConfTempFitVarpTBins", {20, 0.5, 4.05}, "pT binning of the pT vs. TempFitVar plot"}; + ConfigurableAxis confTempFitVarBins{"confTempFitVarBins", {300, -0.15, 0.15}, "Binning of the TempFitVar in the pT vs. TempFitVar plot"}; + ConfigurableAxis confTempFitVarPtBins{"confTempFitVarPtBins", {20, 0.5, 4.05}, "pT binning of the pT vs. TempFitVar plot"}; /// Correlation part - ConfigurableAxis ConfMultBins{"ConfMultBins", {VARIABLE_WIDTH, 0.0f, 10.0f, 20.0f, 30.0f, 40.0f, 50.0f, 60.0f, 70.0f, 80.0f, 90.0f, 100.0f, 99999.f}, "Mixing bins -- multiplicity or centrality"}; // \todo to be obtained from the hash task - ConfigurableAxis ConfMultKstarBins{"ConfMultKstarBins", {VARIABLE_WIDTH, 0.0f, 10.0f, 20.0f, 30.0f, 40.0f, 50.0f, 100.0f, 99999.f}, "Bins for kstar analysis in multiplicity or centrality bins (10 is maximum)"}; - ConfigurableAxis ConfKtKstarBins{"ConfKtKstarBins", {VARIABLE_WIDTH, 0.0f, 0.2f, 0.4f, 0.6f, 0.8f, 1.0f, 2.0f, 99999.f}, "Bins for kstar analysis in kT bins (10 is maximum)"}; - ConfigurableAxis ConfVtxBins{"ConfVtxBins", {VARIABLE_WIDTH, -10.0f, -5.f, 0.f, 5.f, 10.f}, "Mixing bins -- z-vertex"}; - - ConfigurableAxis ConfmTBins3D{"ConfmTBins3D", {VARIABLE_WIDTH, 1.02f, 1.14f, 1.20f, 1.26f, 1.38f, 1.56f, 1.86f, 4.50f}, "mT Binning for the 3Dimensional plot: k* vs multiplicity vs mT (set <> to true in order to use)"}; - ConfigurableAxis ConfmultBins3D{"ConfmultBins3D", {VARIABLE_WIDTH, 0.0f, 20.0f, 30.0f, 40.0f, 99999.0f}, "multiplicity Binning for the 3Dimensional plot: k* vs multiplicity vs mT (set <> to true in order to use)"}; - - ColumnBinningPolicy colBinning{{ConfVtxBins, ConfMultBins}, true}; - - ConfigurableAxis ConfkstarBins{"ConfkstarBins", {1500, 0., 6.}, "binning kstar"}; - ConfigurableAxis ConfkTBins{"ConfkTBins", {150, 0., 9.}, "binning kT"}; - ConfigurableAxis ConfmTBins{"ConfmTBins", {225, 0., 7.5}, "binning mT"}; - Configurable ConfNEventsMix{"ConfNEventsMix", 5, "Number of events for mixing"}; - Configurable ConfIsCPR{"ConfIsCPR", true, "Close Pair Rejection"}; - Configurable ConfCPRPlotPerRadii{"ConfCPRPlotPerRadii", false, "Plot CPR per radii"}; - Configurable ConfCPRdeltaPhiCutMax{"ConfCPRdeltaPhiCutMax", 0.0, "Delta Phi max cut for Close Pair Rejection"}; - Configurable ConfCPRdeltaPhiCutMin{"ConfCPRdeltaPhiCutMin", 0.0, "Delta Phi min cut for Close Pair Rejection"}; - Configurable ConfCPRdeltaEtaCutMax{"ConfCPRdeltaEtaCutMax", 0.0, "Delta Eta max cut for Close Pair Rejection"}; - Configurable ConfCPRdeltaEtaCutMin{"ConfCPRdeltaEtaCutMin", 0.0, "Delta Eta min cut for Close Pair Rejection"}; - Configurable ConfCPRChosenRadii{"ConfCPRChosenRadii", 0.80, "Delta Eta cut for Close Pair Rejection"}; + ConfigurableAxis confMultBins{"confMultBins", {VARIABLE_WIDTH, 0.0f, 10.0f, 20.0f, 30.0f, 40.0f, 50.0f, 60.0f, 70.0f, 80.0f, 90.0f, 100.0f, 99999.f}, "Mixing bins -- multiplicity or centrality"}; // \todo to be obtained from the hash task + ConfigurableAxis confMultKstarBins{"confMultKstarBins", {VARIABLE_WIDTH, 0.0f, 10.0f, 20.0f, 30.0f, 40.0f, 50.0f, 100.0f, 99999.f}, "Bins for kstar analysis in multiplicity or centrality bins (10 is maximum)"}; + ConfigurableAxis confKtKstarBins{"confKtKstarBins", {VARIABLE_WIDTH, 0.0f, 0.2f, 0.4f, 0.6f, 0.8f, 1.0f, 2.0f, 99999.f}, "Bins for kstar analysis in kT bins (10 is maximum)"}; + ConfigurableAxis confVtxBins{"confVtxBins", {VARIABLE_WIDTH, -10.0f, -5.f, 0.f, 5.f, 10.f}, "Mixing bins -- z-vertex"}; + + ConfigurableAxis confMtBins3D{"confMtBins3D", {VARIABLE_WIDTH, 1.02f, 1.14f, 1.20f, 1.26f, 1.38f, 1.56f, 1.86f, 4.50f}, "mT Binning for the 3Dimensional plot: k* vs multiplicity vs mT (set <> to true in order to use)"}; + ConfigurableAxis confMultBins3D{"confMultBins3D", {VARIABLE_WIDTH, 0.0f, 20.0f, 30.0f, 40.0f, 99999.0f}, "Multiplicity binning for the 3Dimensional plot: k* vs multiplicity vs mT (set <> to true in order to use)"}; + + ColumnBinningPolicy colBinning{{confVtxBins, confMultBins}, true}; + + ConfigurableAxis confKstarBins{"confKstarBins", {1500, 0., 6.}, "Binning kstar"}; + ConfigurableAxis confKtBins{"confKtBins", {150, 0., 9.}, "Binning kT"}; + ConfigurableAxis confMtBins{"confMtBins", {225, 0., 7.5}, "Binning mT"}; + Configurable confNEventsMix{"confNEventsMix", 5, "Number of events for mixing"}; + Configurable confIsCPR{"confIsCPR", true, "Close Pair Rejection"}; + Configurable confCPRPlotPerRadii{"confCPRPlotPerRadii", false, "Plot CPR per radii"}; + Configurable confCPRDeltaPhiCutMax{"confCPRDeltaPhiCutMax", 0.0, "Delta Phi max cut for Close Pair Rejection"}; + Configurable confCPRDeltaPhiCutMin{"confCPRDeltaPhiCutMin", 0.0, "Delta Phi min cut for Close Pair Rejection"}; + Configurable confCPRDeltaEtaCutMax{"confCPRDeltaEtaCutMax", 0.0, "Delta Eta max cut for Close Pair Rejection"}; + Configurable confCPRDeltaEtaCutMin{"confCPRDeltaEtaCutMin", 0.0, "Delta Eta min cut for Close Pair Rejection"}; + Configurable confCPRChosenRadii{"confCPRChosenRadii", 0.80, "Delta Eta cut for Close Pair Rejection"}; Configurable cfgProcessPP{"cfgProcessPP", true, "Process positively charged particles (plus-plus)"}; Configurable cfgProcessMM{"cfgProcessMM", true, "Process negatively charged particles (minus-minus)"}; @@ -240,68 +241,68 @@ struct femtoUniversePairTaskTrackNucleus { HistogramRegistry resultRegistryMM{"CorrelationsMM", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; HistogramRegistry resultRegistryPM{"CorrelationsPM", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; HistogramRegistry resultRegistryMP{"CorrelationsMP", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; - HistogramRegistry MixQaRegistry{"MixQaRegistry", {}, OutputObjHandlingPolicy::AnalysisObject}; + HistogramRegistry mixQARegistry{"mixQARegistry", {}, OutputObjHandlingPolicy::AnalysisObject}; - HistogramRegistry SameMultRegistryPP{"SameMultRegistryPP", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; - HistogramRegistry MixedMultRegistryPP{"MixedMultRegistryPP", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; + HistogramRegistry sameMultRegistryPP{"sameMultRegistryPP", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; + HistogramRegistry mixedMultRegistryPP{"mixedMultRegistryPP", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; - HistogramRegistry SameMultRegistryMM{"SameMultRegistryMM", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; - HistogramRegistry MixedMultRegistryMM{"MixedMultRegistryMM", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; + HistogramRegistry sameMultRegistryMM{"sameMultRegistryMM", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; + HistogramRegistry mixedMultRegistryMM{"mixedMultRegistryMM", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; - HistogramRegistry SameMultRegistryPM{"SameMultRegistryPM", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; - HistogramRegistry MixedMultRegistryPM{"MixedMultRegistryPM", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; + HistogramRegistry sameMultRegistryPM{"sameMultRegistryPM", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; + HistogramRegistry mixedMultRegistryPM{"mixedMultRegistryPM", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; - HistogramRegistry SameMultRegistryMP{"SameMultRegistryMP", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; - HistogramRegistry MixedMultRegistryMP{"MixedMultRegistryMP", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; + HistogramRegistry sameMultRegistryMP{"sameMultRegistryMP", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; + HistogramRegistry mixedMultRegistryMP{"mixedMultRegistryMP", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; - HistogramRegistry sphericityRegistry{"SphericityHisto", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; + HistogramRegistry sphericityRegistry{"sphericityRegistry", {}, OutputObjHandlingPolicy::AnalysisObject, true, true}; /// TPC Pion/Kaon/Proton Sigma selection (general) - bool IsNSigma(float mom, float nsigmaTPC, float nsigmaTOF) + bool isNsigma(float mom, float nsigmaTPC, float nsigmaTOF) { // |nsigma_TPC| < 3 for p < 0.5 GeV/c // |nsigma_combined| < 3 for p > 0.5 // using configurables: - // ConfTOFpMin -- momentum value when we start using TOF; set to 1000 if TOF not needed - // ConfNsigmaTPC -> TPC Sigma for momentum < ConfTOFpMin - // ConfNsigmaCombined -> TPC and TOF Sigma (combined) for momentum > ConfTOFpMin + // confTOFpMin -- momentum value when we start using TOF; set to 1000 if TOF not needed + // confNsigmaTPC -> TPC Sigma for momentum < confTOFpMin + // confNsigmaCombined -> TPC and TOF Sigma (combined) for momentum > confTOFpMin - if (mom < twoobjectsconfigs.ConfTOFpMin) { - return TMath::Abs(nsigmaTPC) < twoobjectsconfigs.ConfNsigmaTPC; + if (mom < twoobjectsconfigs.confTOFpMin) { + return abs(nsigmaTPC) < twoobjectsconfigs.confNsigmaTPC; } else { - return TMath::Hypot(nsigmaTOF, nsigmaTPC) < twoobjectsconfigs.ConfNsigmaCombined; + return hypot(nsigmaTOF, nsigmaTPC) < twoobjectsconfigs.confNsigmaCombined; } } /// TPC Kaon Sigma selection (stricter cuts for K+ and K-) -- based on Run2 results - bool IsKaonNSigma(float mom, float nsigmaTPCK, float nsigmaTOFK) + bool isKaonNsigma(float mom, float nsigmaTPCK, float nsigmaTOFK) { - if (twoobjectsconfigs.IsKaonNsigma == true) { + if (twoobjectsconfigs.isKaonNsigma == true) { if (mom < 0.4) { - return TMath::Abs(nsigmaTPCK) < 2; + return abs(nsigmaTPCK) < 2; } else if (mom > 0.4 && mom < 0.45) { - return TMath::Abs(nsigmaTPCK) < 1; + return abs(nsigmaTPCK) < 1; } else if (mom > 0.45 && mom < 0.8) { - return (TMath::Abs(nsigmaTPCK) < 3 && TMath::Abs(nsigmaTOFK) < 2); + return (abs(nsigmaTPCK) < 3 && abs(nsigmaTOFK) < 2); } else if (mom > 0.8 && mom < 1.5) { - return (TMath::Abs(nsigmaTPCK) < 3 && TMath::Abs(nsigmaTOFK) < 1.5); + return (abs(nsigmaTPCK) < 3 && abs(nsigmaTOFK) < 1.5); } else { return false; } } else { - return IsNSigma(mom, nsigmaTPCK, nsigmaTOFK); + return isNsigma(mom, nsigmaTPCK, nsigmaTOFK); } } /// TPC Deuteron Sigma selection - bool IsDeuteronNSigma(float mom, float nsigmaTPCDe, float nsigmaTOFDe) + bool isDeuteronNsigma(float mom, float nsigmaTPCDe, float nsigmaTOFDe) { - if (mom > deuteronconfigs.ConfPLowDe && mom < deuteronconfigs.ConfPHighDe) { - if (mom < deuteronconfigs.ConfTOFpMinDe) { - return (TMath::Abs(nsigmaTPCDe) < deuteronconfigs.ConfNsigmaTPCDe); + if (mom > deuteronconfigs.confPlowDe && mom < deuteronconfigs.confPhighDe) { + if (mom < deuteronconfigs.confTOFpMinDe) { + return (abs(nsigmaTPCDe) < deuteronconfigs.confNsigmaTPCDe); } else { - return (TMath::Abs(nsigmaTOFDe) < deuteronconfigs.ConfNsigmaTOFDe && (TMath::Abs(nsigmaTPCDe) < deuteronconfigs.ConfNsigmaTPCDe)); + return (abs(nsigmaTOFDe) < deuteronconfigs.confNsigmaTOFDe && (abs(nsigmaTPCDe) < deuteronconfigs.confNsigmaTPCDe)); } } else { return false; @@ -309,88 +310,88 @@ struct femtoUniversePairTaskTrackNucleus { } /// Linear cut for clearer TPC Deuteron Sigma - bool IsDeuteronNSigmaLinearCut(float mom, float nsigmaTPCDe, float nsigmaTOFDe, float tpcSignal) + bool isDeuteronNsigmaLinearCut(float mom, float nsigmaTPCDe, float nsigmaTOFDe, float tpcSignal) { - if (lincut.ConfIsLine == true) { - if (mom > lincut.LinCutpLow && mom < lincut.LinCutpHigh) { - if (tpcSignal > lincut.LinCutParA * mom + lincut.LinCutParB) { - return IsDeuteronNSigma(mom, nsigmaTPCDe, nsigmaTOFDe); + if (lincut.confIsLine == true) { + if (mom > lincut.linCutPlow && mom < lincut.linCutPhigh) { + if (tpcSignal > lincut.linCutParA * mom + lincut.linCutParB) { + return isDeuteronNsigma(mom, nsigmaTPCDe, nsigmaTOFDe); } else { return false; } } else { - return IsDeuteronNSigma(mom, nsigmaTPCDe, nsigmaTOFDe); + return isDeuteronNsigma(mom, nsigmaTPCDe, nsigmaTOFDe); } } else { - return IsDeuteronNSigma(mom, nsigmaTPCDe, nsigmaTOFDe); + return isDeuteronNsigma(mom, nsigmaTPCDe, nsigmaTOFDe); } } /// Polynomial 3 cut for clearer TPC Deuteron Sigma - bool IsDeuteronNSigmaPolCut(float mom, float nsigmaTPCDe, float nsigmaTOFDe, float tpcSignal) + bool isDeuteronNsigmaPol3Cut(float mom, float nsigmaTPCDe, float nsigmaTOFDe, float tpcSignal) { - if (polcut.ConfIsPol == true) { - if (tpcSignal > polcut.PolCutParA * TMath::Power(mom, 3) + polcut.PolCutParB * TMath::Power(mom, 2) + polcut.PolCutParC * mom + polcut.PolCutParD) { - return IsDeuteronNSigma(mom, nsigmaTPCDe, nsigmaTOFDe); + if (polcut.confIsPol3 == true) { + if (tpcSignal > polcut.polCutParA * pow(mom, 3) + polcut.polCutParB * pow(mom, 2) + polcut.polCutParC * mom + polcut.polCutParD) { + return isDeuteronNsigma(mom, nsigmaTPCDe, nsigmaTOFDe); } else { return false; } } else { - return IsDeuteronNSigma(mom, nsigmaTPCDe, nsigmaTOFDe); + return isDeuteronNsigma(mom, nsigmaTPCDe, nsigmaTOFDe); } } - bool IsParticleNSigma(int8_t object_number, float mom, float nsigmaTPCPr, float nsigmaTOFPr, float nsigmaTPCPi, float nsigmaTOFPi, float nsigmaTPCK, float nsigmaTOFK, float nsigmaTPCDe, float nsigmaTOFDe, float tpcSignal) + bool isParticleNsigma(int8_t object_number, float mom, float nsigmaTPCPr, float nsigmaTOFPr, float nsigmaTPCPi, float nsigmaTOFPi, float nsigmaTPCK, float nsigmaTOFK, float nsigmaTPCDe, float nsigmaTOFDe, float tpcSignal) { if (object_number == 1) { - switch (trackfilter.ConfPDGCodeTrack) { + switch (trackfilter.confPDGCodeTrack) { case 2212: // Proton case -2212: // Anti-proton - return IsNSigma(mom, nsigmaTPCPr, nsigmaTOFPr); + return isNsigma(mom, nsigmaTPCPr, nsigmaTOFPr); break; case 211: // Pion+ case -211: // Pion- case 111: // Pion 0 - return IsNSigma(mom, nsigmaTPCPi, nsigmaTOFPi); + return isNsigma(mom, nsigmaTPCPi, nsigmaTOFPi); break; case 321: // Kaon+ case -321: // Kaon- - return IsKaonNSigma(mom, nsigmaTPCK, nsigmaTOFK); + return isKaonNsigma(mom, nsigmaTPCK, nsigmaTOFK); break; case 130: // Kaon 0 LONG case 310: // Kaon 0 SHORT - return IsNSigma(mom, nsigmaTPCK, nsigmaTOFK); + return isNsigma(mom, nsigmaTPCK, nsigmaTOFK); break; case 1000010020: // Deuteron case -1000010020: // Anti-deuteron - return IsDeuteronNSigmaPolCut(mom, nsigmaTPCDe, nsigmaTOFDe, tpcSignal); + return isDeuteronNsigmaPol3Cut(mom, nsigmaTPCDe, nsigmaTOFDe, tpcSignal); break; default: return false; } return false; } else if (object_number == 2) { - switch (nucleusfilter.ConfPDGCodeNucleus) { + switch (nucleusfilter.confPDGCodeNucleus) { case 2212: // Proton case -2212: // Anti-proton - return IsNSigma(mom, nsigmaTPCPr, nsigmaTOFPr); + return isNsigma(mom, nsigmaTPCPr, nsigmaTOFPr); break; case 211: // Pion+ case -211: // Pion- case 111: // Pion 0 - return IsNSigma(mom, nsigmaTPCPi, nsigmaTOFPi); + return isNsigma(mom, nsigmaTPCPi, nsigmaTOFPi); break; case 321: // Kaon+ case -321: // Kaon- - return IsKaonNSigma(mom, nsigmaTPCK, nsigmaTOFK); + return isKaonNsigma(mom, nsigmaTPCK, nsigmaTOFK); break; case 130: // Kaon 0 LONG case 310: // Kaon 0 SHORT - return IsNSigma(mom, nsigmaTPCK, nsigmaTOFK); + return isNsigma(mom, nsigmaTPCK, nsigmaTOFK); break; case 1000010020: // Deuteron case -1000010020: // Anti-deuteron - return IsDeuteronNSigmaPolCut(mom, nsigmaTPCDe, nsigmaTOFDe, tpcSignal); + return isDeuteronNsigmaPol3Cut(mom, nsigmaTPCDe, nsigmaTOFDe, tpcSignal); break; default: return false; @@ -407,80 +408,80 @@ struct femtoUniversePairTaskTrackNucleus { eventHisto.init(&qaRegistry); sphericityRegistry.add("sphericity", ";Sphericity;Entries", kTH1F, {{150, 0.0, 3, "Sphericity"}}); - trackHistoTrack.init(&qaRegistry, ConfTempFitVarpTBins, ConfTempFitVarBins, twoobjectsconfigs.ConfIsMC, trackfilter.ConfPDGCodeTrack, true); + trackHistoTrack.init(&qaRegistry, confTempFitVarPtBins, confTempFitVarBins, twoobjectsconfigs.confIsMC, trackfilter.confPDGCodeTrack, true); - trackHistoNucleus.init(&qaRegistry, ConfTempFitVarpTBins, ConfTempFitVarBins, twoobjectsconfigs.ConfIsMC, nucleusfilter.ConfPDGCodeNucleus, true); + trackHistoNucleus.init(&qaRegistry, confTempFitVarPtBins, confTempFitVarBins, twoobjectsconfigs.confIsMC, nucleusfilter.confPDGCodeNucleus, true); - MixQaRegistry.add("MixingQA/hSECollisionBins", ";bin;Entries", kTH1F, {{120, -0.5, 119.5}}); - MixQaRegistry.add("MixingQA/hMECollisionBins", ";bin;Entries", kTH1F, {{120, -0.5, 119.5}}); + mixQARegistry.add("MixingQA/hSECollisionBins", ";bin;Entries", kTH1F, {{120, -0.5, 119.5}}); + mixQARegistry.add("MixingQA/hMECollisionBins", ";bin;Entries", kTH1F, {{120, -0.5, 119.5}}); - mass1 = pdg->Mass(trackfilter.ConfPDGCodeTrack); - mass2 = pdg->Mass(nucleusfilter.ConfPDGCodeNucleus); + mass1 = pdg->Mass(trackfilter.confPDGCodeTrack); + mass2 = pdg->Mass(nucleusfilter.confPDGCodeNucleus); if (cfgProcessPP) { - sameEventContPP.init(&resultRegistryPP, ConfkstarBins, ConfMultBins, ConfkTBins, ConfmTBins, ConfmultBins3D, ConfmTBins3D, twoobjectsconfigs.ConfIsMC, twoobjectsconfigs.ConfUse3D); - mixedEventContPP.init(&resultRegistryPP, ConfkstarBins, ConfMultBins, ConfkTBins, ConfmTBins, ConfmultBins3D, ConfmTBins3D, twoobjectsconfigs.ConfIsMC, twoobjectsconfigs.ConfUse3D); - sameEventContPP.setPDGCodes(trackfilter.ConfPDGCodeTrack, nucleusfilter.ConfPDGCodeNucleus); - mixedEventContPP.setPDGCodes(trackfilter.ConfPDGCodeTrack, nucleusfilter.ConfPDGCodeNucleus); + sameEventContPP.init(&resultRegistryPP, confKstarBins, confMultBins, confKtBins, confMtBins, confMultBins3D, confMtBins3D, twoobjectsconfigs.confIsMC, twoobjectsconfigs.confUse3D); + mixedEventContPP.init(&resultRegistryPP, confKstarBins, confMultBins, confKtBins, confMtBins, confMultBins3D, confMtBins3D, twoobjectsconfigs.confIsMC, twoobjectsconfigs.confUse3D); + sameEventContPP.setPDGCodes(trackfilter.confPDGCodeTrack, nucleusfilter.confPDGCodeNucleus); + mixedEventContPP.setPDGCodes(trackfilter.confPDGCodeTrack, nucleusfilter.confPDGCodeNucleus); if (cfgProcessMultBins) { - sameEventMultContPP.init(&SameMultRegistryPP, ConfkstarBins, ConfMultKstarBins, ConfKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); - mixedEventMultContPP.init(&MixedMultRegistryPP, ConfkstarBins, ConfMultKstarBins, ConfKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); + sameEventMultContPP.init(&sameMultRegistryPP, confKstarBins, confMultKstarBins, confKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); + mixedEventMultContPP.init(&mixedMultRegistryPP, confKstarBins, confMultKstarBins, confKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); } } if (cfgProcessMM) { - sameEventContMM.init(&resultRegistryMM, ConfkstarBins, ConfMultBins, ConfkTBins, ConfmTBins, ConfmultBins3D, ConfmTBins3D, twoobjectsconfigs.ConfIsMC, twoobjectsconfigs.ConfUse3D); - mixedEventContMM.init(&resultRegistryMM, ConfkstarBins, ConfMultBins, ConfkTBins, ConfmTBins, ConfmultBins3D, ConfmTBins3D, twoobjectsconfigs.ConfIsMC, twoobjectsconfigs.ConfUse3D); - sameEventContMM.setPDGCodes(trackfilter.ConfPDGCodeTrack, nucleusfilter.ConfPDGCodeNucleus); - mixedEventContMM.setPDGCodes(trackfilter.ConfPDGCodeTrack, nucleusfilter.ConfPDGCodeNucleus); + sameEventContMM.init(&resultRegistryMM, confKstarBins, confMultBins, confKtBins, confMtBins, confMultBins3D, confMtBins3D, twoobjectsconfigs.confIsMC, twoobjectsconfigs.confUse3D); + mixedEventContMM.init(&resultRegistryMM, confKstarBins, confMultBins, confKtBins, confMtBins, confMultBins3D, confMtBins3D, twoobjectsconfigs.confIsMC, twoobjectsconfigs.confUse3D); + sameEventContMM.setPDGCodes(trackfilter.confPDGCodeTrack, nucleusfilter.confPDGCodeNucleus); + mixedEventContMM.setPDGCodes(trackfilter.confPDGCodeTrack, nucleusfilter.confPDGCodeNucleus); if (cfgProcessMultBins) { - sameEventMultContMM.init(&SameMultRegistryMM, ConfkstarBins, ConfMultKstarBins, ConfKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); - mixedEventMultContMM.init(&MixedMultRegistryMM, ConfkstarBins, ConfMultKstarBins, ConfKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); + sameEventMultContMM.init(&sameMultRegistryMM, confKstarBins, confMultKstarBins, confKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); + mixedEventMultContMM.init(&mixedMultRegistryMM, confKstarBins, confMultKstarBins, confKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); } } if (cfgProcessPM) { - sameEventContPM.init(&resultRegistryPM, ConfkstarBins, ConfMultBins, ConfkTBins, ConfmTBins, ConfmultBins3D, ConfmTBins3D, twoobjectsconfigs.ConfIsMC, twoobjectsconfigs.ConfUse3D); - mixedEventContPM.init(&resultRegistryPM, ConfkstarBins, ConfMultBins, ConfkTBins, ConfmTBins, ConfmultBins3D, ConfmTBins3D, twoobjectsconfigs.ConfIsMC, twoobjectsconfigs.ConfUse3D); + sameEventContPM.init(&resultRegistryPM, confKstarBins, confMultBins, confKtBins, confMtBins, confMultBins3D, confMtBins3D, twoobjectsconfigs.confIsMC, twoobjectsconfigs.confUse3D); + mixedEventContPM.init(&resultRegistryPM, confKstarBins, confMultBins, confKtBins, confMtBins, confMultBins3D, confMtBins3D, twoobjectsconfigs.confIsMC, twoobjectsconfigs.confUse3D); - sameEventContPM.setPDGCodes(trackfilter.ConfPDGCodeTrack, nucleusfilter.ConfPDGCodeNucleus); - mixedEventContPM.setPDGCodes(trackfilter.ConfPDGCodeTrack, nucleusfilter.ConfPDGCodeNucleus); + sameEventContPM.setPDGCodes(trackfilter.confPDGCodeTrack, nucleusfilter.confPDGCodeNucleus); + mixedEventContPM.setPDGCodes(trackfilter.confPDGCodeTrack, nucleusfilter.confPDGCodeNucleus); if (cfgProcessMultBins) { - sameEventMultContPM.init(&SameMultRegistryPM, ConfkstarBins, ConfMultKstarBins, ConfKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); - mixedEventMultContPM.init(&MixedMultRegistryPM, ConfkstarBins, ConfMultKstarBins, ConfKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); + sameEventMultContPM.init(&sameMultRegistryPM, confKstarBins, confMultKstarBins, confKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); + mixedEventMultContPM.init(&mixedMultRegistryPM, confKstarBins, confMultKstarBins, confKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); } } if (cfgProcessMP) { - sameEventContMP.init(&resultRegistryMP, ConfkstarBins, ConfMultBins, ConfkTBins, ConfmTBins, ConfmultBins3D, ConfmTBins3D, twoobjectsconfigs.ConfIsMC, twoobjectsconfigs.ConfUse3D); - mixedEventContMP.init(&resultRegistryMP, ConfkstarBins, ConfMultBins, ConfkTBins, ConfmTBins, ConfmultBins3D, ConfmTBins3D, twoobjectsconfigs.ConfIsMC, twoobjectsconfigs.ConfUse3D); + sameEventContMP.init(&resultRegistryMP, confKstarBins, confMultBins, confKtBins, confMtBins, confMultBins3D, confMtBins3D, twoobjectsconfigs.confIsMC, twoobjectsconfigs.confUse3D); + mixedEventContMP.init(&resultRegistryMP, confKstarBins, confMultBins, confKtBins, confMtBins, confMultBins3D, confMtBins3D, twoobjectsconfigs.confIsMC, twoobjectsconfigs.confUse3D); - sameEventContMP.setPDGCodes(trackfilter.ConfPDGCodeTrack, nucleusfilter.ConfPDGCodeNucleus); - mixedEventContMP.setPDGCodes(trackfilter.ConfPDGCodeTrack, nucleusfilter.ConfPDGCodeNucleus); + sameEventContMP.setPDGCodes(trackfilter.confPDGCodeTrack, nucleusfilter.confPDGCodeNucleus); + mixedEventContMP.setPDGCodes(trackfilter.confPDGCodeTrack, nucleusfilter.confPDGCodeNucleus); if (cfgProcessMultBins) { - sameEventMultContMP.init(&SameMultRegistryMP, ConfkstarBins, ConfMultKstarBins, ConfKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); - mixedEventMultContMP.init(&MixedMultRegistryMP, ConfkstarBins, ConfMultKstarBins, ConfKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); + sameEventMultContMP.init(&sameMultRegistryMP, confKstarBins, confMultKstarBins, confKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); + mixedEventMultContMP.init(&mixedMultRegistryMP, confKstarBins, confMultKstarBins, confKtKstarBins, cfgProcessKtBins, cfgProcessKtMt3DCF); } } pairCleaner.init(&qaRegistry); - if (ConfIsCPR.value) { - pairCloseRejection.init(&resultRegistry, &qaRegistry, ConfCPRdeltaPhiCutMin.value, ConfCPRdeltaPhiCutMax.value, ConfCPRdeltaEtaCutMin.value, ConfCPRdeltaEtaCutMax.value, ConfCPRChosenRadii.value, ConfCPRPlotPerRadii.value); + if (confIsCPR.value) { + pairCloseRejection.init(&resultRegistry, &qaRegistry, confCPRDeltaPhiCutMin.value, confCPRDeltaPhiCutMax.value, confCPRDeltaEtaCutMin.value, confCPRDeltaEtaCutMax.value, confCPRChosenRadii.value, confCPRPlotPerRadii.value); } - vPIDTrack = trackfilter.ConfPIDTrack.value; - vPIDNucleus = nucleusfilter.ConfPIDNucleus.value; - kNsigma = twoobjectsconfigs.ConfTrkPIDnSigmaMax.value; + vPIDTrack = trackfilter.confPIDTrack.value; + vPIDNucleus = nucleusfilter.confPIDNucleus.value; + kNsigma = twoobjectsconfigs.confTrkPIDNsigmaMax.value; } template void fillCollision(CollisionType col) { - MixQaRegistry.fill(HIST("MixingQA/hSECollisionBins"), colBinning.getBin({col.posZ(), col.multV0M()})); + mixQARegistry.fill(HIST("MixingQA/hSECollisionBins"), colBinning.getBin({col.posZ(), col.multV0M()})); eventHisto.fillQA(col); } @@ -499,32 +500,32 @@ struct femtoUniversePairTaskTrackNucleus { template void doSameEvent(PartitionType groupTrack, PartitionType groupNucleus, PartType parts, float magFieldTesla, int multCol, int pairType, bool /*fillQA*/) { - for (auto& part : groupTrack) { - if (!IsParticleNSigma((int8_t)1, part.p(), trackCuts.getNsigmaTPC(part, o2::track::PID::Proton), trackCuts.getNsigmaTOF(part, o2::track::PID::Proton), trackCuts.getNsigmaTPC(part, o2::track::PID::Pion), trackCuts.getNsigmaTOF(part, o2::track::PID::Pion), trackCuts.getNsigmaTPC(part, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(part, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(part, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(part, o2::track::PID::Deuteron), part.tpcSignal())) { + for (const auto& part : groupTrack) { + if (!isParticleNsigma((int8_t)1, part.p(), trackCuts.getNsigmaTPC(part, o2::track::PID::Proton), trackCuts.getNsigmaTOF(part, o2::track::PID::Proton), trackCuts.getNsigmaTPC(part, o2::track::PID::Pion), trackCuts.getNsigmaTOF(part, o2::track::PID::Pion), trackCuts.getNsigmaTPC(part, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(part, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(part, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(part, o2::track::PID::Deuteron), part.tpcSignal())) { continue; } trackHistoTrack.fillQA(part); } - for (auto& part : groupNucleus) { - if (!IsParticleNSigma((int8_t)2, part.p(), trackCuts.getNsigmaTPC(part, o2::track::PID::Proton), trackCuts.getNsigmaTOF(part, o2::track::PID::Proton), trackCuts.getNsigmaTPC(part, o2::track::PID::Pion), trackCuts.getNsigmaTOF(part, o2::track::PID::Pion), trackCuts.getNsigmaTPC(part, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(part, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(part, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(part, o2::track::PID::Deuteron), part.tpcSignal())) { + for (const auto& part : groupNucleus) { + if (!isParticleNsigma((int8_t)2, part.p(), trackCuts.getNsigmaTPC(part, o2::track::PID::Proton), trackCuts.getNsigmaTOF(part, o2::track::PID::Proton), trackCuts.getNsigmaTPC(part, o2::track::PID::Pion), trackCuts.getNsigmaTOF(part, o2::track::PID::Pion), trackCuts.getNsigmaTPC(part, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(part, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(part, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(part, o2::track::PID::Deuteron), part.tpcSignal())) { continue; } trackHistoNucleus.fillQA(part); } - /// Combinations for identical pairs - for (auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupTrack, groupTrack))) { + /// Combinations creation + for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupTrack, groupNucleus))) { - if (!IsParticleNSigma((int8_t)2, p1.p(), trackCuts.getNsigmaTPC(p1, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p1, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p1, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p1, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(p1, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(p1, o2::track::PID::Deuteron), p1.tpcSignal())) { + if (!isParticleNsigma((int8_t)1, p1.p(), trackCuts.getNsigmaTPC(p1, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p1, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p1, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p1, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(p1, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(p1, o2::track::PID::Deuteron), p1.tpcSignal())) { continue; } - if (!IsParticleNSigma((int8_t)2, p2.p(), trackCuts.getNsigmaTPC(p2, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p2, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p2, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p2, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(p2, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(p2, o2::track::PID::Deuteron), p2.tpcSignal())) { + if (!isParticleNsigma((int8_t)2, p2.p(), trackCuts.getNsigmaTPC(p2, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p2, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p2, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p2, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(p2, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(p2, o2::track::PID::Deuteron), p2.tpcSignal())) { continue; } - if (ConfIsCPR.value) { + if (confIsCPR.value) { if (pairCloseRejection.isClosePair(p1, p2, parts, magFieldTesla, femto_universe_container::EventType::same)) { continue; } @@ -540,7 +541,7 @@ struct femtoUniversePairTaskTrackNucleus { float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2); float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2); - sameEventContPP.setPair(p1, p2, multCol, twoobjectsconfigs.ConfUse3D); + sameEventContPP.setPair(p1, p2, multCol, twoobjectsconfigs.confUse3D); if (cfgProcessMultBins) sameEventMultContPP.fill(kstar, multCol, kT); @@ -552,7 +553,7 @@ struct femtoUniversePairTaskTrackNucleus { float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2); float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2); - sameEventContMM.setPair(p1, p2, multCol, twoobjectsconfigs.ConfUse3D); + sameEventContMM.setPair(p1, p2, multCol, twoobjectsconfigs.confUse3D); if (cfgProcessMultBins) sameEventMultContMM.fill(kstar, multCol, kT); @@ -564,7 +565,7 @@ struct femtoUniversePairTaskTrackNucleus { float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2); float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2); - sameEventContPM.setPair(p1, p2, multCol, twoobjectsconfigs.ConfUse3D); + sameEventContPM.setPair(p1, p2, multCol, twoobjectsconfigs.confUse3D); if (cfgProcessMultBins) sameEventMultContPM.fill(kstar, multCol, kT); @@ -576,7 +577,7 @@ struct femtoUniversePairTaskTrackNucleus { float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2); float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2); - sameEventContMP.setPair(p1, p2, multCol, twoobjectsconfigs.ConfUse3D); + sameEventContMP.setPair(p1, p2, multCol, twoobjectsconfigs.confUse3D); if (cfgProcessMultBins) sameEventMultContMP.fill(kstar, multCol, kT); @@ -593,8 +594,8 @@ struct femtoUniversePairTaskTrackNucleus { /// Process function to call doSameEvent with Data /// \param col subscribe to the collision table (Data) /// \param parts subscribe to the femtoUniverseParticleTable - void processSameEvent(soa::Filtered::iterator& col, - FilteredFemtoFullParticles& parts) + void processSameEvent(const soa::Filtered::iterator& col, + const FilteredFemtoFullParticles& parts) { fillCollision(col); sphericityRegistry.fill(HIST("sphericity"), col.sphericity()); @@ -620,15 +621,15 @@ struct femtoUniversePairTaskTrackNucleus { fillQA = false; } } - PROCESS_SWITCH(femtoUniversePairTaskTrackNucleus, processSameEvent, "Enable processing same event", true); + PROCESS_SWITCH(FemtoUniversePairTaskTrackNucleus, processSameEvent, "Enable processing same event", true); /// Process function to call doSameEvent with Monte Carlo /// \param col subscribe to the collision table (Monte Carlo Reconstructed reconstructed) /// \param parts subscribe to joined table FemtoUniverseParticles and FemtoUniverseMCLables to access Monte Carlo truth /// \param FemtoUniverseMCParticles subscribe to the Monte Carlo Truth table - void processSameEventMC(o2::aod::FdCollision& col, - soa::Join& parts, - o2::aod::FdMCParticles&) + void processSameEventMC(const o2::aod::FdCollision& col, + const soa::Join& parts, + const o2::aod::FdMCParticles&) { fillCollision(col); @@ -653,7 +654,7 @@ struct femtoUniversePairTaskTrackNucleus { fillQA = false; } } - PROCESS_SWITCH(femtoUniversePairTaskTrackNucleus, processSameEventMC, "Enable processing same event for Monte Carlo", false); + PROCESS_SWITCH(FemtoUniversePairTaskTrackNucleus, processSameEventMC, "Enable processing same event for Monte Carlo", false); /// This function processes 'mixed event' /// \todo the trivial loops over the collisions and tracks should be factored out since they will be common to all combinations of T-T, T-V0, V0-V0, ... @@ -669,17 +670,17 @@ struct femtoUniversePairTaskTrackNucleus { template void doMixedEvent(PartitionType groupTrack, PartitionType groupNucleus, PartType parts, float magFieldTesla, int multCol, int pairType) { - for (auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupTrack, groupNucleus))) { + for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupTrack, groupNucleus))) { - if (!IsParticleNSigma((int8_t)2, p1.p(), trackCuts.getNsigmaTPC(p1, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p1, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p1, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p1, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(p1, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(p1, o2::track::PID::Deuteron), p1.tpcSignal())) { + if (!isParticleNsigma((int8_t)2, p1.p(), trackCuts.getNsigmaTPC(p1, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p1, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p1, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p1, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(p1, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(p1, o2::track::PID::Deuteron), p1.tpcSignal())) { continue; } - if (!IsParticleNSigma((int8_t)2, p2.p(), trackCuts.getNsigmaTPC(p2, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p2, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p2, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p2, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(p2, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(p2, o2::track::PID::Deuteron), p2.tpcSignal())) { + if (!isParticleNsigma((int8_t)2, p2.p(), trackCuts.getNsigmaTPC(p2, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p2, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p2, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p2, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(p2, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(p2, o2::track::PID::Deuteron), p2.tpcSignal())) { continue; } - if (ConfIsCPR.value) { + if (confIsCPR.value) { if (pairCloseRejection.isClosePair(p1, p2, parts, magFieldTesla, femto_universe_container::EventType::mixed)) { continue; } @@ -690,7 +691,7 @@ struct femtoUniversePairTaskTrackNucleus { float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2); float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2); - mixedEventContPP.setPair(p1, p2, multCol, twoobjectsconfigs.ConfUse3D); + mixedEventContPP.setPair(p1, p2, multCol, twoobjectsconfigs.confUse3D); if (cfgProcessMultBins) mixedEventMultContPP.fill(kstar, multCol, kT); @@ -702,7 +703,7 @@ struct femtoUniversePairTaskTrackNucleus { float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2); float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2); - mixedEventContMM.setPair(p1, p2, multCol, twoobjectsconfigs.ConfUse3D); + mixedEventContMM.setPair(p1, p2, multCol, twoobjectsconfigs.confUse3D); if (cfgProcessMultBins) mixedEventMultContMM.fill(kstar, multCol, kT); @@ -714,7 +715,7 @@ struct femtoUniversePairTaskTrackNucleus { float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2); float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2); - mixedEventContPM.setPair(p1, p2, multCol, twoobjectsconfigs.ConfUse3D); + mixedEventContPM.setPair(p1, p2, multCol, twoobjectsconfigs.confUse3D); if (cfgProcessMultBins) mixedEventMultContPM.fill(kstar, multCol, kT); @@ -726,7 +727,7 @@ struct femtoUniversePairTaskTrackNucleus { float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2); float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2); - mixedEventContMP.setPair(p1, p2, multCol, twoobjectsconfigs.ConfUse3D); + mixedEventContMP.setPair(p1, p2, multCol, twoobjectsconfigs.confUse3D); if (cfgProcessMultBins) mixedEventMultContMP.fill(kstar, multCol, kT); @@ -743,13 +744,13 @@ struct femtoUniversePairTaskTrackNucleus { /// Process function to call doMixedEvent with Data /// \param cols subscribe to the collisions table (Data) /// \param parts subscribe to the femtoUniverseParticleTable - void processMixedEvent(soa::Filtered& cols, - FilteredFemtoFullParticles& parts) + void processMixedEvent(const soa::Filtered& cols, + const FilteredFemtoFullParticles& parts) { - for (auto& [collision1, collision2] : soa::selfCombinations(colBinning, 5, -1, cols, cols)) { + for (const auto& [collision1, collision2] : soa::selfCombinations(colBinning, 5, -1, cols, cols)) { const int multiplicityCol = collision1.multV0M(); - MixQaRegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinning.getBin({collision1.posZ(), multiplicityCol})); + mixQARegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinning.getBin({collision1.posZ(), multiplicityCol})); const auto& magFieldTesla1 = collision1.magField(); const auto& magFieldTesla2 = collision2.magField(); @@ -783,20 +784,20 @@ struct femtoUniversePairTaskTrackNucleus { } } } - PROCESS_SWITCH(femtoUniversePairTaskTrackNucleus, processMixedEvent, "Enable processing mixed events", true); + PROCESS_SWITCH(FemtoUniversePairTaskTrackNucleus, processMixedEvent, "Enable processing mixed events", true); /// Process function to call doMixedEvent with Monte Carlo /// \param cols subscribe to the collisions table (Monte Carlo Reconstructed reconstructed) /// \param parts subscribe to joined table FemtoUniverseParticles and FemtoUniverseMCLables to access Monte Carlo truth /// \param FemtoUniverseMCParticles subscribe to the Monte Carlo truth table - void processMixedEventMC(o2::aod::FdCollisions& cols, - soa::Join& parts, - o2::aod::FdMCParticles&) + void processMixedEventMC(const o2::aod::FdCollisions& cols, + const soa::Join& parts, + const o2::aod::FdMCParticles&) { - for (auto& [collision1, collision2] : soa::selfCombinations(colBinning, 5, -1, cols, cols)) { + for (const auto& [collision1, collision2] : soa::selfCombinations(colBinning, 5, -1, cols, cols)) { const int multiplicityCol = collision1.multV0M(); - MixQaRegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinning.getBin({collision1.posZ(), multiplicityCol})); + mixQARegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinning.getBin({collision1.posZ(), multiplicityCol})); const auto& magFieldTesla1 = collision1.magField(); const auto& magFieldTesla2 = collision2.magField(); @@ -805,7 +806,7 @@ struct femtoUniversePairTaskTrackNucleus { continue; } /// \todo before mixing we should check whether both collisions contain a pair of particles! - // if partTrack.size() == 0 || nPart2Evt1 == 0 || nPart1Evt2 == 0 || partNucleus.size() == 0 ) continue; + // if partTrack.size() == 0 || Npart2Evt1 == 0 || Npart1Evt2 == 0 || partNucleus.size() == 0 ) continue; if (cfgProcessPP) { auto groupTrack = partTrack->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision1.globalIndex(), cache); @@ -832,13 +833,13 @@ struct femtoUniversePairTaskTrackNucleus { } } } - PROCESS_SWITCH(femtoUniversePairTaskTrackNucleus, processMixedEventMC, "Enable processing mixed events MC", false); + PROCESS_SWITCH(FemtoUniversePairTaskTrackNucleus, processMixedEventMC, "Enable processing mixed events MC", false); }; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { WorkflowSpec workflow{ - adaptAnalysisTask(cfgc), + adaptAnalysisTask(cfgc), }; return workflow; From 781919c704e76823b0d10a7bc8bd9379f1bebd36 Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Mon, 13 Jan 2025 22:04:12 +0100 Subject: [PATCH 02/19] Track-nucleus combinations correction --- .../femtoUniversePairTaskTrackNucleus.cxx | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx index 8927326aadf..e672fc99af3 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx @@ -269,9 +269,9 @@ struct FemtoUniversePairTaskTrackNucleus { // confNsigmaCombined -> TPC and TOF Sigma (combined) for momentum > confTOFpMin if (mom < twoobjectsconfigs.confTOFpMin) { - return abs(nsigmaTPC) < twoobjectsconfigs.confNsigmaTPC; + return std::abs(nsigmaTPC) < twoobjectsconfigs.confNsigmaTPC; } else { - return hypot(nsigmaTOF, nsigmaTPC) < twoobjectsconfigs.confNsigmaCombined; + return std::hypot(nsigmaTOF, nsigmaTPC) < twoobjectsconfigs.confNsigmaCombined; } } @@ -280,13 +280,13 @@ struct FemtoUniversePairTaskTrackNucleus { { if (twoobjectsconfigs.isKaonNsigma == true) { if (mom < 0.4) { - return abs(nsigmaTPCK) < 2; + return std::abs(nsigmaTPCK) < 2; } else if (mom > 0.4 && mom < 0.45) { - return abs(nsigmaTPCK) < 1; + return std::abs(nsigmaTPCK) < 1; } else if (mom > 0.45 && mom < 0.8) { - return (abs(nsigmaTPCK) < 3 && abs(nsigmaTOFK) < 2); + return (std::abs(nsigmaTPCK) < 3 && std::abs(nsigmaTOFK) < 2); } else if (mom > 0.8 && mom < 1.5) { - return (abs(nsigmaTPCK) < 3 && abs(nsigmaTOFK) < 1.5); + return (std::abs(nsigmaTPCK) < 3 && std::abs(nsigmaTOFK) < 1.5); } else { return false; } @@ -300,9 +300,9 @@ struct FemtoUniversePairTaskTrackNucleus { { if (mom > deuteronconfigs.confPlowDe && mom < deuteronconfigs.confPhighDe) { if (mom < deuteronconfigs.confTOFpMinDe) { - return (abs(nsigmaTPCDe) < deuteronconfigs.confNsigmaTPCDe); + return (std::abs(nsigmaTPCDe) < deuteronconfigs.confNsigmaTPCDe); } else { - return (abs(nsigmaTOFDe) < deuteronconfigs.confNsigmaTOFDe && (abs(nsigmaTPCDe) < deuteronconfigs.confNsigmaTPCDe)); + return (std::abs(nsigmaTOFDe) < deuteronconfigs.confNsigmaTOFDe && (std::abs(nsigmaTPCDe) < deuteronconfigs.confNsigmaTPCDe)); } } else { return false; @@ -331,7 +331,7 @@ struct FemtoUniversePairTaskTrackNucleus { bool isDeuteronNsigmaPol3Cut(float mom, float nsigmaTPCDe, float nsigmaTOFDe, float tpcSignal) { if (polcut.confIsPol3 == true) { - if (tpcSignal > polcut.polCutParA * pow(mom, 3) + polcut.polCutParB * pow(mom, 2) + polcut.polCutParC * mom + polcut.polCutParD) { + if (tpcSignal > polcut.polCutParA * std::pow(mom, 3) + polcut.polCutParB * std::pow(mom, 2) + polcut.polCutParC * mom + polcut.polCutParD) { return isDeuteronNsigma(mom, nsigmaTPCDe, nsigmaTOFDe); } else { return false; From fcea1f75df22641b89d9ba8275a49acb9acc2b76 Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Tue, 14 Jan 2025 14:27:49 +0100 Subject: [PATCH 03/19] Track-nucleus combinations correction --- .../Tasks/femtoUniversePairTaskTrackNucleus.cxx | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx index e672fc99af3..3c9151940d6 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx @@ -92,17 +92,17 @@ struct FemtoUniversePairTaskTrackNucleus { Configurable confIsLine{"confIsLine", false, "Enable a separation line for clearer TPC Deuteron Sigma"}; Configurable linCutPlow{"linCutPlow", 0.0f, "Lower limit of momentum for linear cut of TPC Deuteron Sigma"}; Configurable linCutPhigh{"linCutPhigh", 1.4f, "Higher limit of momentum for linear cut of TPC Deuteron Sigma"}; - Configurable linCutParA{"linCutParA", -167.0f, "Parameter 'A' of a linear function 'y = A * x + B'"}; - Configurable linCutParB{"linCutParB", 300.0f, "Parameter 'B' of a linear function 'y = A * x + B'"}; + Configurable linCutParA{"linCutParA", -167.0f, "Parameter 'A' of the linear function 'y = A * x + B'"}; + Configurable linCutParB{"linCutParB", 300.0f, "Parameter 'B' of the linear function 'y = A * x + B'"}; } lincut; /// Table for polynomial 3 cut for TPC Deuteron Sigma struct : o2::framework::ConfigurableGroup { Configurable confIsPol3{"confIsPol3", false, "Enable a separation polynomial 3 curve for clearer TPC Deuteron Sigma"}; - Configurable polCutParA{"polCutParA", -52.2f, "Parameter 'A' of a polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; - Configurable polCutParB{"polCutParB", 357.7f, "Parameter 'B' of a polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; - Configurable polCutParC{"polCutParC", -834.7f, "Parameter 'C' of a polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; - Configurable polCutParD{"polCutParD", 705.8f, "Parameter 'D' of a polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; + Configurable polCutParA{"polCutParA", -52.2f, "Parameter 'A' of the polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; + Configurable polCutParB{"polCutParB", 357.7f, "Parameter 'B' of the polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; + Configurable polCutParC{"polCutParC", -834.7f, "Parameter 'C' of the polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; + Configurable polCutParD{"polCutParD", 705.8f, "Parameter 'D' of the polynomial function 'y = A * x^3 + B * x^2 + C * x + D'"}; } polcut; using FemtoFullParticles = soa::Join; From 18521071411e61433f717e974c31e00f58f783a8 Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Thu, 16 Jan 2025 13:37:37 +0100 Subject: [PATCH 04/19] Correct combinations mixing --- .../FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx index 3c9151940d6..2fc7a7cf2fa 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackNucleus.cxx @@ -176,7 +176,7 @@ struct FemtoUniversePairTaskTrackNucleus { ConfigurableAxis confMultBins{"confMultBins", {VARIABLE_WIDTH, 0.0f, 10.0f, 20.0f, 30.0f, 40.0f, 50.0f, 60.0f, 70.0f, 80.0f, 90.0f, 100.0f, 99999.f}, "Mixing bins -- multiplicity or centrality"}; // \todo to be obtained from the hash task ConfigurableAxis confMultKstarBins{"confMultKstarBins", {VARIABLE_WIDTH, 0.0f, 10.0f, 20.0f, 30.0f, 40.0f, 50.0f, 100.0f, 99999.f}, "Bins for kstar analysis in multiplicity or centrality bins (10 is maximum)"}; ConfigurableAxis confKtKstarBins{"confKtKstarBins", {VARIABLE_WIDTH, 0.0f, 0.2f, 0.4f, 0.6f, 0.8f, 1.0f, 2.0f, 99999.f}, "Bins for kstar analysis in kT bins (10 is maximum)"}; - ConfigurableAxis confVtxBins{"confVtxBins", {VARIABLE_WIDTH, -10.0f, -5.f, 0.f, 5.f, 10.f}, "Mixing bins -- z-vertex"}; + ConfigurableAxis confVtxBins{"confVtxBins", {VARIABLE_WIDTH, -10.0f, 10.f}, "Mixing bins -- z-vertex"}; ConfigurableAxis confMtBins3D{"confMtBins3D", {VARIABLE_WIDTH, 1.02f, 1.14f, 1.20f, 1.26f, 1.38f, 1.56f, 1.86f, 4.50f}, "mT Binning for the 3Dimensional plot: k* vs multiplicity vs mT (set <> to true in order to use)"}; ConfigurableAxis confMultBins3D{"confMultBins3D", {VARIABLE_WIDTH, 0.0f, 20.0f, 30.0f, 40.0f, 99999.0f}, "Multiplicity binning for the 3Dimensional plot: k* vs multiplicity vs mT (set <> to true in order to use)"}; @@ -672,7 +672,7 @@ struct FemtoUniversePairTaskTrackNucleus { { for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupTrack, groupNucleus))) { - if (!isParticleNsigma((int8_t)2, p1.p(), trackCuts.getNsigmaTPC(p1, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p1, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p1, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p1, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(p1, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(p1, o2::track::PID::Deuteron), p1.tpcSignal())) { + if (!isParticleNsigma((int8_t)1, p1.p(), trackCuts.getNsigmaTPC(p1, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p1, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p1, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p1, o2::track::PID::Pion), trackCuts.getNsigmaTPC(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTPC(p1, o2::track::PID::Deuteron), trackCuts.getNsigmaTOF(p1, o2::track::PID::Deuteron), p1.tpcSignal())) { continue; } From 7c7177f270853b72bd05ea2a0da349227261d744 Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Mon, 19 May 2025 19:18:51 +0200 Subject: [PATCH 05/19] Add helicity angle analysis --- PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h | 8 +++--- .../Tasks/femtoUniverseDebugV0.cxx | 25 +++++++++++++++---- 2 files changed, 23 insertions(+), 10 deletions(-) diff --git a/PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h b/PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h index 556f075d819..58fa2e943d9 100644 --- a/PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h +++ b/PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h @@ -64,13 +64,14 @@ class FemtoUniverseMath return 0.5 * trackRelK.P(); } + /// Boost particles from LAB Frame to Pair Rest Frame (for lambda daughters) /// \tparam T type of tracks /// \param part1 Particle 1 /// \param mass1 Mass of particle 1 /// \param part2 Particle 2 /// \param mass2 Mass of particle 2 template - static float getthetastar(const T& part1, const float mass1, const T& part2, const float mass2) + static ROOT::Math::PxPyPzMVector boostPRF(const T& part1, const float mass1, const T& part2, const float mass2) { const ROOT::Math::PtEtaPhiMVector vecpart1(part1.pt(), part1.eta(), part1.phi(), mass1); const ROOT::Math::PtEtaPhiMVector vecpart2(part2.pt(), part2.eta(), part2.phi(), mass2); @@ -88,10 +89,7 @@ class FemtoUniverseMath partOneCMS = boostPRF(partOneCMS); partTwoCMS = boostPRF(partTwoCMS); - const ROOT::Math::PtEtaPhiMVector partOneCMSGeo(partOneCMS); - const ROOT::Math::PtEtaPhiMVector partTwoCMSGeo(partTwoCMS); - - return (partOneCMSGeo.Theta() - partTwoCMSGeo.Theta()); + return partOneCMS; } /// Compute the qij of a pair of particles diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx index 219c43472dd..4a94c7e1f3a 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx @@ -24,12 +24,14 @@ #include "Framework/ASoAHelpers.h" #include "Framework/RunningWorkflowInfo.h" #include "Framework/StepTHn.h" +#include "Framework/O2DatabasePDGPlugin.h" #include "DataFormatsParameters/GRPObject.h" #include "PWGCF/FemtoUniverse/DataModel/FemtoDerived.h" #include "PWGCF/FemtoUniverse/Core/FemtoUniverseParticleHisto.h" #include "PWGCF/FemtoUniverse/Core/FemtoUniverseEventHisto.h" #include "PWGCF/FemtoUniverse/Core/femtoUtils.h" +#include "PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h" using namespace o2; using namespace o2::analysis::femto_universe; @@ -38,6 +40,9 @@ using namespace o2::framework::expressions; using namespace o2::soa; struct femtoUniverseDebugV0 { + + Service pdg; + SliceCache cache; Configurable ConfPDGCodeV0{"ConfPDGCodePartOne", 3122, "V0 - PDG code"}; @@ -71,6 +76,7 @@ struct femtoUniverseDebugV0 { /// Histogram output HistogramRegistry EventRegistry{"Event", {}, OutputObjHandlingPolicy::AnalysisObject}; HistogramRegistry V0Registry{"FullV0QA", {}, OutputObjHandlingPolicy::AnalysisObject}; + HistogramRegistry thetaRegistry{"ThetaQA", {}, OutputObjHandlingPolicy::AnalysisObject}; void init(InitContext&) { @@ -78,14 +84,16 @@ struct femtoUniverseDebugV0 { posChildHistos.init(&V0Registry, ConfChildTempFitVarpTBins, ConfChildTempFitVarBins, false, ConfPDGCodeChildPos.value, true); negChildHistos.init(&V0Registry, ConfChildTempFitVarpTBins, ConfChildTempFitVarBins, false, ConfPDGCodeChildNeg, true); V0Histos.init(&V0Registry, ConfV0TempFitVarpTBins, ConfV0TempFitVarBins, false, ConfPDGCodeV0.value, true); + + thetaRegistry.add("Theta/hTheta", " ; p (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {50, -5, 5}}); } - /// Porduce QA plots for V0 selection in FemtoUniverse framework + /// Produce QA plots for V0 selection in FemtoUniverse framework void process(o2::aod::FdCollision const& col, FemtoFullParticles const& parts) { auto groupPartsOne = partsOne->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); eventHisto.fillQA(col); - for (auto& part : groupPartsOne) { + for (const auto& part : groupPartsOne) { if (!part.has_children()) { continue; } @@ -95,14 +103,21 @@ struct femtoUniverseDebugV0 { LOG(warn) << "Indices of V0 children do not match"; continue; } - // check cuts on V0 children - if ((posChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && (posChild.cut() & ConfCutChildPos) == ConfCutChildPos) && - (negChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && (negChild.cut() & ConfCutChildNeg) == ConfCutChildNeg) && + + // Check cuts on V0 children + if (posChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && + negChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && isFullPIDSelected(posChild.pidCut(), posChild.p(), 999.f, ConfChildPosIndex.value, ConfChildnSpecies.value, ConfChildPIDnSigmaMax.value, ConfChildPosPidnSigmaMax.value, 1.f) && isFullPIDSelected(negChild.pidCut(), negChild.p(), 999.f, ConfChildNegIndex.value, ConfChildnSpecies.value, ConfChildPIDnSigmaMax.value, ConfChildNegPidnSigmaMax.value, 1.f)) { + auto protonMass = pdg->Mass(2212); + auto pionMass = pdg->Mass(211); + auto protonBoosted = FemtoUniverseMath::boostPRF(posChild, protonMass, negChild, pionMass); + auto cosineTheta = (protonBoosted.Px() * part.px() + protonBoosted.Py() * part.py() + protonBoosted.Pz() * part.pz()) / (protonBoosted.P() * part.p()); + V0Histos.fillQA(part); posChildHistos.fillQA(posChild); negChildHistos.fillQA(negChild); + thetaRegistry.fill(HIST("Theta/hTheta"), part.p(), cosineTheta); } } } From e7fbbdc8265660136b81f062011d6e9d8dbaa8a0 Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Mon, 26 May 2025 16:14:43 +0200 Subject: [PATCH 06/19] Correction --- .../Tasks/femtoUniverseDebugV0.cxx | 57 ++++++++++--------- 1 file changed, 29 insertions(+), 28 deletions(-) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx index 4a94c7e1f3a..73fc43c00e8 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx @@ -16,8 +16,8 @@ #include #include -#include #include +#include #include "Framework/AnalysisTask.h" #include "Framework/runDataProcessing.h" #include "Framework/HistogramRegistry.h" @@ -39,32 +39,32 @@ using namespace o2::framework; using namespace o2::framework::expressions; using namespace o2::soa; -struct femtoUniverseDebugV0 { +struct FemtoUniverseDebugV0 { Service pdg; SliceCache cache; - Configurable ConfPDGCodeV0{"ConfPDGCodePartOne", 3122, "V0 - PDG code"}; - Configurable ConfPDGCodeChildPos{"ConfPDGCodeChildPos", 2212, "Positive Child - PDG code"}; - Configurable ConfPDGCodeChildNeg{"ConfPDGCodeChildNeg", 211, "Negative Child- PDG code"}; - Configurable ConfCutV0{"ConfCutV0", 338, "V0 - Selection bit from cutCulator"}; - ConfigurableAxis ConfV0TempFitVarBins{"ConfV0TempFitVarBins", {300, 0.95, 1.}, "V0: binning of the TempFitVar in the pT vs. TempFitVar plot"}; - ConfigurableAxis ConfV0TempFitVarpTBins{"ConfV0TempFitVarpTBins", {20, 0.5, 4.05}, "V0: pT binning of the pT vs. TempFitVar plot"}; - - Configurable ConfCutChildPos{"ConfCutChildPos", 150, "Positive Child of V0 - Selection bit from cutCulator"}; - Configurable ConfCutChildNeg{"ConfCutChildNeg", 149, "Negative Child of V0 - Selection bit from cutCulator"}; - Configurable ConfChildPosPidnSigmaMax{"ConfChildPosPidnSigmaMax", 3.f, "Positive Child of V0 - Selection bit from cutCulator"}; - Configurable ConfChildNegPidnSigmaMax{"ConfChildNegPidnSigmaMax", 3.f, "Negative Child of V0 - Selection bit from cutCulator"}; - Configurable ConfChildPosIndex{"ConfChildPosIndex", 1, "Positive Child of V0 - Index from cutCulator"}; - Configurable ConfChildNegIndex{"ConfChildNegIndex", 0, "Negative Child of V0 - Index from cutCulator"}; - Configurable> ConfChildPIDnSigmaMax{"ConfChildPIDnSigmaMax", std::vector{4.f, 3.f}, "V0 child sel: Max. PID nSigma TPC"}; - Configurable ConfChildnSpecies{"ConfChildnSpecies", 2, "Number of particle spieces (for V0 children) with PID info"}; - ConfigurableAxis ConfChildTempFitVarBins{"ConfChildTempFitVarBins", {300, -0.15, 0.15}, "V0 child: binning of the TempFitVar in the pT vs. TempFitVar plot"}; - ConfigurableAxis ConfChildTempFitVarpTBins{"ConfChildTempFitVarpTBins", {20, 0.5, 4.05}, "V0 child: pT binning of the pT vs. TempFitVar plot"}; + Configurable confPDGCodeV0{"confPDGCodeV0", 3122, "V0 -- PDG code"}; + Configurable confPDGCodePositiveChild{"confPDGCodePositiveChild", 2212, "Positive Child -- PDG code"}; + Configurable confPDGCodeNegativeChild{"confPDGCodeNegativeChild", 211, "Negative Child -- PDG code"}; + Configurable confCutV0{"confCutV0", 338, "V0 -- Selection bit from cutCulator"}; + ConfigurableAxis confV0TempFitVarBins{"confV0TempFitVarBins", {300, 0.95, 1.}, "V0: binning of the TempFitVar in the pT vs. TempFitVar plot"}; + ConfigurableAxis confV0TempFitVarpTBins{"confV0TempFitVarpTBins", {20, 0.5, 4.05}, "V0: pT binning of the pT vs. TempFitVar plot"}; + + Configurable confCutPositiveChild{"confCutPositiveChild", 150, "Positive Child of V0 -- Selection bit from cutCulator"}; + Configurable confCutNegativeChild{"confCutNegativeChild", 149, "Negative Child of V0 -- Selection bit from cutCulator"}; + Configurable confPositiveChildPIDnSigmaMax{"confPositiveChildPIDnSigmaMax", 3.f, "Positive Child of V0 -- Selection bit from cutCulator"}; + Configurable confNegativeChildPIDnSigmaMax{"confNegativeChildPIDnSigmaMax", 3.f, "Negative Child of V0 -- Selection bit from cutCulator"}; + Configurable confPositiveChildIndex{"confPositiveChildIndex", 1, "Positive Child of V0 -- Index from cutCulator"}; + Configurable confNegativeChildIndex{"confNegativeChildIndex", 0, "Negative Child of V0 -- Index from cutCulator"}; + Configurable> confChildPIDnSigmaMax{"confChildPIDnSigmaMax", std::vector{4.f, 3.f}, "V0 child selection: max. PID nSigma TPC"}; + Configurable confChildnSpecies{"confChildnSpecies", 2, "Number of particle spieces (for V0 children) with PID info"}; + ConfigurableAxis confChildTempFitVarBins{"confChildTempFitVarBins", {300, -0.15, 0.15}, "V0 child: binning of the TempFitVar in the pT vs. TempFitVar plot"}; + ConfigurableAxis confChildTempFitVarpTBins{"confChildTempFitVarpTBins", {20, 0.5, 4.05}, "V0 child: pT binning of the pT vs. TempFitVar plot"}; using FemtoFullParticles = soa::Join; - Partition partsOne = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kV0)) && ((aod::femtouniverseparticle::cut & ConfCutV0) == ConfCutV0); + Partition partsOne = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kV0)) && ((aod::femtouniverseparticle::cut & confCutV0) == confCutV0); Preslice perCol = aod::femtouniverseparticle::fdCollisionId; /// Histogramming @@ -81,9 +81,9 @@ struct femtoUniverseDebugV0 { void init(InitContext&) { eventHisto.init(&EventRegistry); - posChildHistos.init(&V0Registry, ConfChildTempFitVarpTBins, ConfChildTempFitVarBins, false, ConfPDGCodeChildPos.value, true); - negChildHistos.init(&V0Registry, ConfChildTempFitVarpTBins, ConfChildTempFitVarBins, false, ConfPDGCodeChildNeg, true); - V0Histos.init(&V0Registry, ConfV0TempFitVarpTBins, ConfV0TempFitVarBins, false, ConfPDGCodeV0.value, true); + posChildHistos.init(&V0Registry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, confPDGCodePositiveChild.value, true); + negChildHistos.init(&V0Registry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, confPDGCodeNegativeChild, true); + V0Histos.init(&V0Registry, confV0TempFitVarpTBins, confV0TempFitVarBins, false, confPDGCodeV0.value, true); thetaRegistry.add("Theta/hTheta", " ; p (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {50, -5, 5}}); } @@ -107,10 +107,11 @@ struct femtoUniverseDebugV0 { // Check cuts on V0 children if (posChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && negChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && - isFullPIDSelected(posChild.pidCut(), posChild.p(), 999.f, ConfChildPosIndex.value, ConfChildnSpecies.value, ConfChildPIDnSigmaMax.value, ConfChildPosPidnSigmaMax.value, 1.f) && - isFullPIDSelected(negChild.pidCut(), negChild.p(), 999.f, ConfChildNegIndex.value, ConfChildnSpecies.value, ConfChildPIDnSigmaMax.value, ConfChildNegPidnSigmaMax.value, 1.f)) { - auto protonMass = pdg->Mass(2212); - auto pionMass = pdg->Mass(211); + isFullPIDSelected(posChild.pidCut(), posChild.p(), 999.f, confPositiveChildIndex.value, confChildnSpecies.value, confChildPIDnSigmaMax.value, confPositiveChildPIDnSigmaMax.value, 1.f) && + isFullPIDSelected(negChild.pidCut(), negChild.p(), 999.f, confNegativeChildIndex.value, confChildnSpecies.value, confChildPIDnSigmaMax.value, confNegativeChildPIDnSigmaMax.value, 1.f)) { + auto pdgDB = TDatabasePDG::Instance(); + auto protonMass = pdgDB->GetParticle(confPDGCodePositiveChild)->Mass(); + auto pionMass = pdgDB->GetParticle(confPDGCodeNegativeChild)->Mass(); auto protonBoosted = FemtoUniverseMath::boostPRF(posChild, protonMass, negChild, pionMass); auto cosineTheta = (protonBoosted.Px() * part.px() + protonBoosted.Py() * part.py() + protonBoosted.Pz() * part.pz()) / (protonBoosted.P() * part.p()); @@ -127,7 +128,7 @@ WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { WorkflowSpec workflow{ - adaptAnalysisTask(cfgc), + adaptAnalysisTask(cfgc), }; return workflow; } From 117b5e7a741795bd9d953665630586960d0cc15f Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Fri, 6 Jun 2025 17:13:04 +0200 Subject: [PATCH 07/19] Fix mass getters --- .../Tasks/femtoUniverseDebugV0.cxx | 30 +++++++++++++++---- 1 file changed, 25 insertions(+), 5 deletions(-) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx index 73fc43c00e8..6c3f4520fe4 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx @@ -69,8 +69,8 @@ struct FemtoUniverseDebugV0 { /// Histogramming FemtoUniverseEventHisto eventHisto; - FemtoUniverseParticleHisto posChildHistos; - FemtoUniverseParticleHisto negChildHistos; + FemtoUniverseParticleHisto positiveChildHistos; + FemtoUniverseParticleHisto negativeChildHistos; FemtoUniverseParticleHisto V0Histos; /// Histogram output @@ -81,8 +81,13 @@ struct FemtoUniverseDebugV0 { void init(InitContext&) { eventHisto.init(&EventRegistry); +<<<<<<< Updated upstream posChildHistos.init(&V0Registry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, confPDGCodePositiveChild.value, true); negChildHistos.init(&V0Registry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, confPDGCodeNegativeChild, true); +======= + positiveChildHistos.init(&V0Registry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, confPDGCodePositiveChild.value, true); + negativeChildHistos.init(&V0Registry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, confPDGCodeNegativeChild, true); +>>>>>>> Stashed changes V0Histos.init(&V0Registry, confV0TempFitVarpTBins, confV0TempFitVarBins, false, confPDGCodeV0.value, true); thetaRegistry.add("Theta/hTheta", " ; p (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {50, -5, 5}}); @@ -97,14 +102,15 @@ struct FemtoUniverseDebugV0 { if (!part.has_children()) { continue; } - const auto& posChild = parts.iteratorAt(part.index() - 2); - const auto& negChild = parts.iteratorAt(part.index() - 1); - if (posChild.globalIndex() != part.childrenIds()[0] || negChild.globalIndex() != part.childrenIds()[1]) { + const auto& positiveChild = parts.iteratorAt(part.index() - 2); + const auto& negativeChild = parts.iteratorAt(part.index() - 1); + if (positiveChild.globalIndex() != part.childrenIds()[0] || negativeChild.globalIndex() != part.childrenIds()[1]) { LOG(warn) << "Indices of V0 children do not match"; continue; } // Check cuts on V0 children +<<<<<<< Updated upstream if (posChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && negChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && isFullPIDSelected(posChild.pidCut(), posChild.p(), 999.f, confPositiveChildIndex.value, confChildnSpecies.value, confChildPIDnSigmaMax.value, confPositiveChildPIDnSigmaMax.value, 1.f) && @@ -118,6 +124,20 @@ struct FemtoUniverseDebugV0 { V0Histos.fillQA(part); posChildHistos.fillQA(posChild); negChildHistos.fillQA(negChild); +======= + if (positiveChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && + negativeChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && + isFullPIDSelected(positiveChild.pidCut(), positiveChild.p(), 999.f, confPositiveChildIndex.value, confChildnSpecies.value, confChildPIDnSigmaMax.value, confPositiveChildPIDnSigmaMax.value, 1.f) && + isFullPIDSelected(negativeChild.pidCut(), negativeChild.p(), 999.f, confNegativeChildIndex.value, confChildnSpecies.value, confChildPIDnSigmaMax.value, confNegativeChildPIDnSigmaMax.value, 1.f)) { + auto positiveChildMass = pdg->Mass(confPDGCodePositiveChild); + auto negativeChildMass = pdg->Mass(confPDGCodeNegativeChild); + auto positiveChildBoosted = FemtoUniverseMath::boostPRF(positiveChild, positiveChildMass, negativeChild, negativeChildMass); + auto cosineTheta = (positiveChildBoosted.Px() * part.px() + positiveChildBoosted.Py() * part.py() + positiveChildBoosted.Pz() * part.pz()) / (positiveChildBoosted.P() * part.p()); + + V0Histos.fillQA(part); + positiveChildHistos.fillQA(positiveChild); + negativeChildHistos.fillQA(negativeChild); +>>>>>>> Stashed changes thetaRegistry.fill(HIST("Theta/hTheta"), part.p(), cosineTheta); } } From cbbe847aa2a09f63fc42a36b5e9fd7497be65545 Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Fri, 6 Jun 2025 17:51:00 +0200 Subject: [PATCH 08/19] Improve mass getters --- .../Tasks/femtoUniverseDebugV0.cxx | 21 ------------------- 1 file changed, 21 deletions(-) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx index 6c3f4520fe4..9b1ff267b6f 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx @@ -81,13 +81,8 @@ struct FemtoUniverseDebugV0 { void init(InitContext&) { eventHisto.init(&EventRegistry); -<<<<<<< Updated upstream - posChildHistos.init(&V0Registry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, confPDGCodePositiveChild.value, true); - negChildHistos.init(&V0Registry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, confPDGCodeNegativeChild, true); -======= positiveChildHistos.init(&V0Registry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, confPDGCodePositiveChild.value, true); negativeChildHistos.init(&V0Registry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, confPDGCodeNegativeChild, true); ->>>>>>> Stashed changes V0Histos.init(&V0Registry, confV0TempFitVarpTBins, confV0TempFitVarBins, false, confPDGCodeV0.value, true); thetaRegistry.add("Theta/hTheta", " ; p (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {50, -5, 5}}); @@ -110,21 +105,6 @@ struct FemtoUniverseDebugV0 { } // Check cuts on V0 children -<<<<<<< Updated upstream - if (posChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && - negChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && - isFullPIDSelected(posChild.pidCut(), posChild.p(), 999.f, confPositiveChildIndex.value, confChildnSpecies.value, confChildPIDnSigmaMax.value, confPositiveChildPIDnSigmaMax.value, 1.f) && - isFullPIDSelected(negChild.pidCut(), negChild.p(), 999.f, confNegativeChildIndex.value, confChildnSpecies.value, confChildPIDnSigmaMax.value, confNegativeChildPIDnSigmaMax.value, 1.f)) { - auto pdgDB = TDatabasePDG::Instance(); - auto protonMass = pdgDB->GetParticle(confPDGCodePositiveChild)->Mass(); - auto pionMass = pdgDB->GetParticle(confPDGCodeNegativeChild)->Mass(); - auto protonBoosted = FemtoUniverseMath::boostPRF(posChild, protonMass, negChild, pionMass); - auto cosineTheta = (protonBoosted.Px() * part.px() + protonBoosted.Py() * part.py() + protonBoosted.Pz() * part.pz()) / (protonBoosted.P() * part.p()); - - V0Histos.fillQA(part); - posChildHistos.fillQA(posChild); - negChildHistos.fillQA(negChild); -======= if (positiveChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && negativeChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && isFullPIDSelected(positiveChild.pidCut(), positiveChild.p(), 999.f, confPositiveChildIndex.value, confChildnSpecies.value, confChildPIDnSigmaMax.value, confPositiveChildPIDnSigmaMax.value, 1.f) && @@ -137,7 +117,6 @@ struct FemtoUniverseDebugV0 { V0Histos.fillQA(part); positiveChildHistos.fillQA(positiveChild); negativeChildHistos.fillQA(negativeChild); ->>>>>>> Stashed changes thetaRegistry.fill(HIST("Theta/hTheta"), part.p(), cosineTheta); } } From 17f5d872d839c3ddd0df3faaf2794b2ee2fd25c4 Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Tue, 29 Jul 2025 12:51:06 +0200 Subject: [PATCH 09/19] Add MC analysis to debug-V0 task --- .../Tasks/femtoUniverseDebugV0.cxx | 146 +++++++++++++----- 1 file changed, 111 insertions(+), 35 deletions(-) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx index 9b1ff267b6f..a5ebb4a0055 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx @@ -13,11 +13,11 @@ /// \brief Tasks that reads the particle tables and fills QA histograms for V0s /// \author Luca Barioglio, TU München, luca.barioglio@cern.ch /// \author Zuzanna Chochulska, WUT Warsaw & CTU Prague, zchochul@cern.ch +/// \author Anna-Mariia Andrushko, WUT Warsaw, anna-mariia.andrushko@cern.ch #include #include #include -#include #include "Framework/AnalysisTask.h" #include "Framework/runDataProcessing.h" #include "Framework/HistogramRegistry.h" @@ -26,6 +26,8 @@ #include "Framework/StepTHn.h" #include "Framework/O2DatabasePDGPlugin.h" #include "DataFormatsParameters/GRPObject.h" +#include "ReconstructionDataFormats/PID.h" +#include "Common/DataModel/PIDResponse.h" #include "PWGCF/FemtoUniverse/DataModel/FemtoDerived.h" #include "PWGCF/FemtoUniverse/Core/FemtoUniverseParticleHisto.h" @@ -45,51 +47,73 @@ struct FemtoUniverseDebugV0 { SliceCache cache; - Configurable confPDGCodeV0{"confPDGCodeV0", 3122, "V0 -- PDG code"}; - Configurable confPDGCodePositiveChild{"confPDGCodePositiveChild", 2212, "Positive Child -- PDG code"}; - Configurable confPDGCodeNegativeChild{"confPDGCodeNegativeChild", 211, "Negative Child -- PDG code"}; - Configurable confCutV0{"confCutV0", 338, "V0 -- Selection bit from cutCulator"}; - ConfigurableAxis confV0TempFitVarBins{"confV0TempFitVarBins", {300, 0.95, 1.}, "V0: binning of the TempFitVar in the pT vs. TempFitVar plot"}; - ConfigurableAxis confV0TempFitVarpTBins{"confV0TempFitVarpTBins", {20, 0.5, 4.05}, "V0: pT binning of the pT vs. TempFitVar plot"}; - - Configurable confCutPositiveChild{"confCutPositiveChild", 150, "Positive Child of V0 -- Selection bit from cutCulator"}; - Configurable confCutNegativeChild{"confCutNegativeChild", 149, "Negative Child of V0 -- Selection bit from cutCulator"}; - Configurable confPositiveChildPIDnSigmaMax{"confPositiveChildPIDnSigmaMax", 3.f, "Positive Child of V0 -- Selection bit from cutCulator"}; - Configurable confNegativeChildPIDnSigmaMax{"confNegativeChildPIDnSigmaMax", 3.f, "Negative Child of V0 -- Selection bit from cutCulator"}; - Configurable confPositiveChildIndex{"confPositiveChildIndex", 1, "Positive Child of V0 -- Index from cutCulator"}; - Configurable confNegativeChildIndex{"confNegativeChildIndex", 0, "Negative Child of V0 -- Index from cutCulator"}; - Configurable> confChildPIDnSigmaMax{"confChildPIDnSigmaMax", std::vector{4.f, 3.f}, "V0 child selection: max. PID nSigma TPC"}; - Configurable confChildnSpecies{"confChildnSpecies", 2, "Number of particle spieces (for V0 children) with PID info"}; - ConfigurableAxis confChildTempFitVarBins{"confChildTempFitVarBins", {300, -0.15, 0.15}, "V0 child: binning of the TempFitVar in the pT vs. TempFitVar plot"}; - ConfigurableAxis confChildTempFitVarpTBins{"confChildTempFitVarpTBins", {20, 0.5, 4.05}, "V0 child: pT binning of the pT vs. TempFitVar plot"}; - + /// V0 configurables + struct : o2::framework::ConfigurableGroup { + Configurable confPDGCodeV0{"confPDGCodeV0", 3122, "V0 -- PDG code"}; + Configurable confPDGCodePositiveChild{"confPDGCodePositiveChild", 2212, "Positive Child -- PDG code"}; + Configurable confPDGCodeNegativeChild{"confPDGCodeNegativeChild", 211, "Negative Child -- PDG code"}; + Configurable confCutV0{"confCutV0", 338, "V0 -- Selection bit from cutCulator"}; + ConfigurableAxis confV0TempFitVarBins{"confV0TempFitVarBins", {300, 0.95, 1.}, "V0: binning of the TempFitVar in the pT vs. TempFitVar plot"}; + ConfigurableAxis confV0TempFitVarpTBins{"confV0TempFitVarpTBins", {20, 0.5, 4.05}, "V0: pT binning of the pT vs. TempFitVar plot"}; + } V0configs; // o2-linter: disable=name/function-variable + + /// Children configurables + struct : o2::framework::ConfigurableGroup { + Configurable confCutPositiveChild{"confCutPositiveChild", 150, "Positive Child of V0 -- Selection bit from cutCulator"}; + Configurable confCutNegativeChild{"confCutNegativeChild", 149, "Negative Child of V0 -- Selection bit from cutCulator"}; + Configurable confPositiveChildPIDnSigmaMax{"confPositiveChildPIDnSigmaMax", 3.f, "Positive Child of V0 -- Selection bit from cutCulator"}; + Configurable confNegativeChildPIDnSigmaMax{"confNegativeChildPIDnSigmaMax", 3.f, "Negative Child of V0 -- Selection bit from cutCulator"}; + Configurable confPositiveChildIndex{"confPositiveChildIndex", 1, "Positive Child of V0 -- Index from cutCulator"}; + Configurable confNegativeChildIndex{"confNegativeChildIndex", 0, "Negative Child of V0 -- Index from cutCulator"}; + Configurable> confChildPIDnSigmaMax{"confChildPIDnSigmaMax", std::vector{4.f, 3.f}, "V0 child selection: max. PID nSigma TPC"}; + Configurable confChildnSpecies{"confChildnSpecies", 2, "Number of particle spieces (for V0 children) with PID info"}; + ConfigurableAxis confChildTempFitVarBins{"confChildTempFitVarBins", {300, -0.15, 0.15}, "V0 child: binning of the TempFitVar in the pT vs. TempFitVar plot"}; + ConfigurableAxis confChildTempFitVarpTBins{"confChildTempFitVarpTBins", {20, 0.5, 4.05}, "V0 child: pT binning of the pT vs. TempFitVar plot"}; + } childconfigs; + + Configurable confIsMC{"confIsMC", false, "Enable additional histograms in the case of a Monte Carlo run"}; + + + /// Partitioning using FemtoFullParticles = soa::Join; - Partition partsOne = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kV0)) && ((aod::femtouniverseparticle::cut & confCutV0) == confCutV0); + Partition partsOne = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kV0)) && ((aod::femtouniverseparticle::cut & V0configs.confCutV0) == V0configs.confCutV0); + + Partition> partsOneMC = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kV0)) && ((aod::femtouniverseparticle::cut & V0configs.confCutV0) == V0configs.confCutV0); + Preslice perCol = aod::femtouniverseparticle::fdCollisionId; + /// Histogramming FemtoUniverseEventHisto eventHisto; FemtoUniverseParticleHisto positiveChildHistos; FemtoUniverseParticleHisto negativeChildHistos; - FemtoUniverseParticleHisto V0Histos; + FemtoUniverseParticleHisto V0Histos; // o2-linter: disable=name/function-variable /// Histogram output - HistogramRegistry EventRegistry{"Event", {}, OutputObjHandlingPolicy::AnalysisObject}; - HistogramRegistry V0Registry{"FullV0QA", {}, OutputObjHandlingPolicy::AnalysisObject}; + HistogramRegistry eventRegistry{"Event", {}, OutputObjHandlingPolicy::AnalysisObject}; + HistogramRegistry V0Registry{"FullV0QA", {}, OutputObjHandlingPolicy::AnalysisObject}; // o2-linter: disable=name/function-variable HistogramRegistry thetaRegistry{"ThetaQA", {}, OutputObjHandlingPolicy::AnalysisObject}; + void init(InitContext&) { - eventHisto.init(&EventRegistry); - positiveChildHistos.init(&V0Registry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, confPDGCodePositiveChild.value, true); - negativeChildHistos.init(&V0Registry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, confPDGCodeNegativeChild, true); - V0Histos.init(&V0Registry, confV0TempFitVarpTBins, confV0TempFitVarBins, false, confPDGCodeV0.value, true); - - thetaRegistry.add("Theta/hTheta", " ; p (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {50, -5, 5}}); + eventHisto.init(&eventRegistry); + positiveChildHistos.init(&V0Registry, childconfigs.confChildTempFitVarpTBins, childconfigs.confChildTempFitVarBins, confIsMC, V0configs.confPDGCodePositiveChild.value, true); + negativeChildHistos.init(&V0Registry, childconfigs.confChildTempFitVarpTBins, childconfigs.confChildTempFitVarBins, confIsMC, V0configs.confPDGCodeNegativeChild.value, true); + V0Histos.init(&V0Registry, V0configs.confV0TempFitVarpTBins, V0configs.confV0TempFitVarBins, confIsMC, V0configs.confPDGCodeV0.value, true); + + thetaRegistry.add("Theta/hTheta", " ; p (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/PositiveChild/hThetaPt", " ; p_{T} (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/PositiveChild/hThetaEta", " ; #eta; cos(#theta)", kTH2F, {{100, -1, 1}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/PositiveChild/hThetaPhi", " ; #phi; cos(#theta)", kTH2F, {{100, -1, 7}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/NegativeChild/hThetaPt", " ; p_{T} (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/NegativeChild/hThetaEta", " ; #eta; cos(#theta)", kTH2F, {{100, -1, 1}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/NegativeChild/hThetaPhi", " ; #phi; cos(#theta)", kTH2F, {{100, -1, 7}, {110, -1.1, 1.1}}); } - /// Produce QA plots for V0 selection in FemtoUniverse framework - void process(o2::aod::FdCollision const& col, FemtoFullParticles const& parts) + + /// Produce QA plots for V0 and its children on real data + void processData(o2::aod::FdCollision const& col, FemtoFullParticles const& parts) { auto groupPartsOne = partsOne->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); eventHisto.fillQA(col); @@ -107,22 +131,74 @@ struct FemtoUniverseDebugV0 { // Check cuts on V0 children if (positiveChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && negativeChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && - isFullPIDSelected(positiveChild.pidCut(), positiveChild.p(), 999.f, confPositiveChildIndex.value, confChildnSpecies.value, confChildPIDnSigmaMax.value, confPositiveChildPIDnSigmaMax.value, 1.f) && - isFullPIDSelected(negativeChild.pidCut(), negativeChild.p(), 999.f, confNegativeChildIndex.value, confChildnSpecies.value, confChildPIDnSigmaMax.value, confNegativeChildPIDnSigmaMax.value, 1.f)) { - auto positiveChildMass = pdg->Mass(confPDGCodePositiveChild); - auto negativeChildMass = pdg->Mass(confPDGCodeNegativeChild); + isFullPIDSelected(positiveChild.pidCut(), positiveChild.p(), 999.f, childconfigs.confPositiveChildIndex.value, childconfigs.confChildnSpecies.value, childconfigs.confChildPIDnSigmaMax.value, childconfigs.confPositiveChildPIDnSigmaMax.value, 1.f) && + isFullPIDSelected(negativeChild.pidCut(), negativeChild.p(), 999.f, childconfigs.confNegativeChildIndex.value, childconfigs.confChildnSpecies.value, childconfigs.confChildPIDnSigmaMax.value, childconfigs.confNegativeChildPIDnSigmaMax.value, 1.f)) { + auto positiveChildMass = pdg->Mass(V0configs.confPDGCodePositiveChild); + auto negativeChildMass = pdg->Mass(V0configs.confPDGCodeNegativeChild); + auto positiveChildBoosted = FemtoUniverseMath::boostPRF(positiveChild, positiveChildMass, negativeChild, negativeChildMass); + auto cosineTheta = (positiveChildBoosted.Px() * part.px() + positiveChildBoosted.Py() * part.py() + positiveChildBoosted.Pz() * part.pz()) / (positiveChildBoosted.P() * part.p()); + + V0Histos.fillQA(part); + positiveChildHistos.fillQA(positiveChild); + negativeChildHistos.fillQA(negativeChild); + + thetaRegistry.fill(HIST("Theta/hTheta"), part.p(), cosineTheta); + thetaRegistry.fill(HIST("Theta/PositiveChild/hThetaPt"), positiveChild.pt(), cosineTheta); + thetaRegistry.fill(HIST("Theta/PositiveChild/hThetaEta"), positiveChild.eta(), cosineTheta); + thetaRegistry.fill(HIST("Theta/PositiveChild/hThetaPhi"), positiveChild.phi(), cosineTheta); + thetaRegistry.fill(HIST("Theta/NegativeChild/hThetaPt"), negativeChild.pt(), cosineTheta); + thetaRegistry.fill(HIST("Theta/NegativeChild/hThetaEta"), negativeChild.eta(), cosineTheta); + thetaRegistry.fill(HIST("Theta/NegativeChild/hThetaPhi"), negativeChild.phi(), cosineTheta); + } + } + } + PROCESS_SWITCH(FemtoUniverseDebugV0, processData, "Enable processing on real data", true); + + + /// Produce QA plots for V0 and its children on MonteCarlo data + void processMC(o2::aod::FdCollision const& col, soa::Join const& parts, o2::aod::FdMCParticles const&) + { + auto groupPartsOne = partsOneMC->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); + eventHisto.fillQA(col); + for (const auto& part : groupPartsOne) { + if (!part.has_children()) { + continue; + } + const auto& positiveChild = parts.iteratorAt(part.index() - 2); + const auto& negativeChild = parts.iteratorAt(part.index() - 1); + if (positiveChild.globalIndex() != part.childrenIds()[0] || negativeChild.globalIndex() != part.childrenIds()[1]) { + LOG(warn) << "Indices of V0 children do not match"; + continue; + } + + // Check cuts on V0 children + if (positiveChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && + negativeChild.partType() == uint8_t(aod::femtouniverseparticle::ParticleType::kV0Child) && + isFullPIDSelected(positiveChild.pidCut(), positiveChild.p(), 999.f, childconfigs.confPositiveChildIndex.value, childconfigs.confChildnSpecies.value, childconfigs.confChildPIDnSigmaMax.value, childconfigs.confPositiveChildPIDnSigmaMax.value, 1.f) && + isFullPIDSelected(negativeChild.pidCut(), negativeChild.p(), 999.f, childconfigs.confNegativeChildIndex.value, childconfigs.confChildnSpecies.value, childconfigs.confChildPIDnSigmaMax.value, childconfigs.confNegativeChildPIDnSigmaMax.value, 1.f)) { + auto positiveChildMass = pdg->Mass(V0configs.confPDGCodePositiveChild); + auto negativeChildMass = pdg->Mass(V0configs.confPDGCodeNegativeChild); auto positiveChildBoosted = FemtoUniverseMath::boostPRF(positiveChild, positiveChildMass, negativeChild, negativeChildMass); auto cosineTheta = (positiveChildBoosted.Px() * part.px() + positiveChildBoosted.Py() * part.py() + positiveChildBoosted.Pz() * part.pz()) / (positiveChildBoosted.P() * part.p()); V0Histos.fillQA(part); positiveChildHistos.fillQA(positiveChild); negativeChildHistos.fillQA(negativeChild); + thetaRegistry.fill(HIST("Theta/hTheta"), part.p(), cosineTheta); + thetaRegistry.fill(HIST("Theta/PositiveChild/hThetaPt"), positiveChild.pt(), cosineTheta); + thetaRegistry.fill(HIST("Theta/PositiveChild/hThetaEta"), positiveChild.eta(), cosineTheta); + thetaRegistry.fill(HIST("Theta/PositiveChild/hThetaPhi"), positiveChild.phi(), cosineTheta); + thetaRegistry.fill(HIST("Theta/NegativeChild/hThetaPt"), negativeChild.pt(), cosineTheta); + thetaRegistry.fill(HIST("Theta/NegativeChild/hThetaEta"), negativeChild.eta(), cosineTheta); + thetaRegistry.fill(HIST("Theta/NegativeChild/hThetaPhi"), negativeChild.phi(), cosineTheta); } } } + PROCESS_SWITCH(FemtoUniverseDebugV0, processMC, "Enable processing on Monte Carlo", false); }; + WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { From 44af01ea8d65a2199dc562ed3b7d7d99fcd9be0d Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Tue, 29 Jul 2025 14:09:11 +0200 Subject: [PATCH 10/19] Add MC analysis to debug-V0 task --- PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx | 6 ------ 1 file changed, 6 deletions(-) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx index a5ebb4a0055..48018d1b95f 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx @@ -73,7 +73,6 @@ struct FemtoUniverseDebugV0 { Configurable confIsMC{"confIsMC", false, "Enable additional histograms in the case of a Monte Carlo run"}; - /// Partitioning using FemtoFullParticles = soa::Join; Partition partsOne = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kV0)) && ((aod::femtouniverseparticle::cut & V0configs.confCutV0) == V0configs.confCutV0); @@ -82,7 +81,6 @@ struct FemtoUniverseDebugV0 { Preslice perCol = aod::femtouniverseparticle::fdCollisionId; - /// Histogramming FemtoUniverseEventHisto eventHisto; FemtoUniverseParticleHisto positiveChildHistos; @@ -94,7 +92,6 @@ struct FemtoUniverseDebugV0 { HistogramRegistry V0Registry{"FullV0QA", {}, OutputObjHandlingPolicy::AnalysisObject}; // o2-linter: disable=name/function-variable HistogramRegistry thetaRegistry{"ThetaQA", {}, OutputObjHandlingPolicy::AnalysisObject}; - void init(InitContext&) { eventHisto.init(&eventRegistry); @@ -111,7 +108,6 @@ struct FemtoUniverseDebugV0 { thetaRegistry.add("Theta/NegativeChild/hThetaPhi", " ; #phi; cos(#theta)", kTH2F, {{100, -1, 7}, {110, -1.1, 1.1}}); } - /// Produce QA plots for V0 and its children on real data void processData(o2::aod::FdCollision const& col, FemtoFullParticles const& parts) { @@ -154,7 +150,6 @@ struct FemtoUniverseDebugV0 { } PROCESS_SWITCH(FemtoUniverseDebugV0, processData, "Enable processing on real data", true); - /// Produce QA plots for V0 and its children on MonteCarlo data void processMC(o2::aod::FdCollision const& col, soa::Join const& parts, o2::aod::FdMCParticles const&) { @@ -198,7 +193,6 @@ struct FemtoUniverseDebugV0 { PROCESS_SWITCH(FemtoUniverseDebugV0, processMC, "Enable processing on Monte Carlo", false); }; - WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { From d8409fa6eba19731634a95c807c16efab9ff1c3c Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Tue, 29 Jul 2025 14:24:28 +0200 Subject: [PATCH 11/19] Add MC analysis to debug-V0 task --- PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx index 48018d1b95f..48f8c496eee 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx @@ -150,7 +150,7 @@ struct FemtoUniverseDebugV0 { } PROCESS_SWITCH(FemtoUniverseDebugV0, processData, "Enable processing on real data", true); - /// Produce QA plots for V0 and its children on MonteCarlo data + /// Produce QA plots for V0 and its children on Monte Carlo void processMC(o2::aod::FdCollision const& col, soa::Join const& parts, o2::aod::FdMCParticles const&) { auto groupPartsOne = partsOneMC->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); From 9544464c2c1ad91a7bb410f5896f419f70da68fa Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Tue, 29 Jul 2025 14:47:22 +0200 Subject: [PATCH 12/19] Add MC analysis to debug-V0 task --- .../Tasks/femtoUniverseDebugV0.cxx | 32 +++++++++++-------- 1 file changed, 19 insertions(+), 13 deletions(-) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx index 48f8c496eee..4812f8329c1 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx @@ -15,25 +15,30 @@ /// \author Zuzanna Chochulska, WUT Warsaw & CTU Prague, zchochul@cern.ch /// \author Anna-Mariia Andrushko, WUT Warsaw, anna-mariia.andrushko@cern.ch -#include -#include -#include +#include "PWGCF/FemtoUniverse/Core/FemtoUniverseEventHisto.h" +#include "PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h" +#include "PWGCF/FemtoUniverse/Core/FemtoUniverseParticleHisto.h" +#include "PWGCF/FemtoUniverse/Core/femtoUtils.h" +#include "PWGCF/FemtoUniverse/DataModel/FemtoDerived.h" + +#include "Common/DataModel/PIDResponse.h" + +#include "DataFormatsParameters/GRPObject.h" +#include "Framework/ASoAHelpers.h" #include "Framework/AnalysisTask.h" -#include "Framework/runDataProcessing.h" + #include "Framework/HistogramRegistry.h" -#include "Framework/ASoAHelpers.h" +#include "Framework/O2DatabasePDGPlugin.h" #include "Framework/RunningWorkflowInfo.h" #include "Framework/StepTHn.h" -#include "Framework/O2DatabasePDGPlugin.h" -#include "DataFormatsParameters/GRPObject.h" +#include "Framework/runDataProcessing.h" + #include "ReconstructionDataFormats/PID.h" -#include "Common/DataModel/PIDResponse.h" -#include "PWGCF/FemtoUniverse/DataModel/FemtoDerived.h" -#include "PWGCF/FemtoUniverse/Core/FemtoUniverseParticleHisto.h" -#include "PWGCF/FemtoUniverse/Core/FemtoUniverseEventHisto.h" -#include "PWGCF/FemtoUniverse/Core/femtoUtils.h" -#include "PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h" +#include + +#include +#include using namespace o2; using namespace o2::analysis::femto_universe; @@ -75,6 +80,7 @@ struct FemtoUniverseDebugV0 { /// Partitioning using FemtoFullParticles = soa::Join; + Partition partsOne = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kV0)) && ((aod::femtouniverseparticle::cut & V0configs.confCutV0) == V0configs.confCutV0); Partition> partsOneMC = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kV0)) && ((aod::femtouniverseparticle::cut & V0configs.confCutV0) == V0configs.confCutV0); From 0f68e9fe0a3d8469e07acd1d2d02ad7504c461bb Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Tue, 29 Jul 2025 14:55:49 +0200 Subject: [PATCH 13/19] Add MC analysis to debug-V0 task --- PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx | 2 -- 1 file changed, 2 deletions(-) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx index 4812f8329c1..b0700b5acae 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniverseDebugV0.cxx @@ -26,13 +26,11 @@ #include "DataFormatsParameters/GRPObject.h" #include "Framework/ASoAHelpers.h" #include "Framework/AnalysisTask.h" - #include "Framework/HistogramRegistry.h" #include "Framework/O2DatabasePDGPlugin.h" #include "Framework/RunningWorkflowInfo.h" #include "Framework/StepTHn.h" #include "Framework/runDataProcessing.h" - #include "ReconstructionDataFormats/PID.h" #include From e0ff0c34971598613d7b64e8260e46b89e056337 Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Sat, 2 Aug 2025 15:16:53 +0200 Subject: [PATCH 14/19] Add helicity task --- PWGCF/FemtoUniverse/Tasks/CMakeLists.txt | 7 +- .../femtoUniversePairTaskTrackV0Helicity.cxx | 948 ++++++++++++++++++ 2 files changed, 954 insertions(+), 1 deletion(-) create mode 100644 PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx diff --git a/PWGCF/FemtoUniverse/Tasks/CMakeLists.txt b/PWGCF/FemtoUniverse/Tasks/CMakeLists.txt index 0a4682a7f7a..387a2ed1958 100644 --- a/PWGCF/FemtoUniverse/Tasks/CMakeLists.txt +++ b/PWGCF/FemtoUniverse/Tasks/CMakeLists.txt @@ -33,6 +33,11 @@ o2physics_add_dpl_workflow(femtouniverse-pair-track-nucleus SOURCES femtoUniversePairTaskTrackNucleus.cxx PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore COMPONENT_NAME Analysis) + +o2physics_add_dpl_workflow(femtouniverse-pair-track-v0-helicity + SOURCES femtoUniversePairTaskTrackV0Helicity.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore + COMPONENT_NAME Analysis) o2physics_add_dpl_workflow(femtouniverse-pair-track-track-threedrelmom-mult-kt-extended SOURCES femtoUniversePairTaskTrackTrack3DMultKtExtended.cxx @@ -97,4 +102,4 @@ o2physics_add_dpl_workflow(femtouniverse-efficiency-base o2physics_add_executable(femtouniverse-cutculator SOURCES femtoUniverseCutCulator.cxx PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore - COMPONENT_NAME Analysis) \ No newline at end of file + COMPONENT_NAME Analysis) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx new file mode 100644 index 00000000000..6e7052e5f3a --- /dev/null +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx @@ -0,0 +1,948 @@ +// Copyright 2019-2022 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// \file femtoUniversePairTaskTrackV0Helicity.cxx +/// \brief Tasks that build pairs of track particles and v0s +/// \author Andi Mathis, TU München, andreas.mathis@ph.tum.de +/// \author Zuzanna Chochulska, WUT Warsaw & CTU Prague, zchochul@cern.ch +/// \author Shirajum Monira, WUT Warsaw, shirajum.monira.dokt@pw.edu.pl + +#include +#include +#include +#include "Framework/AnalysisTask.h" +#include "Framework/runDataProcessing.h" +#include "Framework/HistogramRegistry.h" +#include "Framework/ASoAHelpers.h" +#include "Framework/RunningWorkflowInfo.h" +#include "Framework/StepTHn.h" +#include "PWGCF/FemtoUniverse/DataModel/FemtoDerived.h" +#include "PWGCF/FemtoUniverse/Core/FemtoUniverseParticleHisto.h" +#include "PWGCF/FemtoUniverse/Core/FemtoUniverseEventHisto.h" +#include "PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h" +#include "PWGCF/FemtoUniverse/Core/FemtoUniversePairCleaner.h" +#include "PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h" +#include "PWGCF/FemtoUniverse/Core/FemtoUniverseDetaDphiStar.h" +#include "PWGCF/FemtoUniverse/Core/femtoUtils.h" +#include "Framework/O2DatabasePDGPlugin.h" +#include +#include + +using namespace o2; +using namespace o2::soa; +using namespace o2::framework; +using namespace o2::framework::expressions; +using namespace o2::analysis::femto_universe; +using namespace o2::aod::pidutils; +using namespace o2::track; + +struct FemtoUniversePairTaskTrackV0Helicity { + + Service pdg; + Service pdgMC; + + SliceCache cache; + using FemtoFullParticles = soa::Join; + Preslice perCol = aod::femtouniverseparticle::fdCollisionId; + + using FemtoRecoParticles = soa::Join; + Preslice perColMC = aod::femtouniverseparticle::fdCollisionId; + + /// To apply narrow cut + Configurable confZVertexCut{"confZVertexCut", 10.f, "Event sel: Maximum z-Vertex (cm)"}; + Configurable confEta{"confEta", 0.8, "Eta cut for the global track"}; + + /// Particle 1 (track) + Configurable confTrkPDGCodePartOne{"confTrkPDGCodePartOne", 211, "Particle 1 (Track) - PDG code"}; + Configurable confTrackChoicePartOne{"confTrackChoicePartOne", 1, "0:Proton, 1:Pion, 2:Kaon"}; + ConfigurableAxis confTrkTempFitVarBins{"confTrkTempFitVarBins", {300, -0.15, 0.15}, "binning of the TempFitVar in the pT vs. TempFitVar plot"}; + ConfigurableAxis confTrkTempFitVarpTBins{"confTrkTempFitVarpTBins", {20, 0.5, 4.05}, "pT binning of the pT vs. TempFitVar plot"}; + Configurable confChargePart1{"confChargePart1", 0, "sign of particle 1"}; + Configurable confHPtPart1{"confHPtPart1", 4.0f, "higher limit for pt of particle 1"}; + Configurable confLPtPart1{"confLPtPart1", 0.3f, "lower limit for pt of particle 1"}; + Configurable confmom{"confmom", 0.5, "momentum threshold for particle identification using TOF"}; + Configurable confNsigmaTPCParticle{"confNsigmaTPCParticle", 3.0, "TPC Sigma for particle momentum < confmom"}; + Configurable confNsigmaCombinedParticle{"confNsigmaCombinedParticle", 3.0, "TPC and TOF Sigma (combined) for particle momentum > confmom"}; + + Filter collisionFilter = (nabs(aod::collision::posZ) < confZVertexCut); + using FilteredFDCollisions = soa::Filtered; + using FilteredFDCollision = FilteredFDCollisions::iterator; + + /// Partition for particle 1 + Partition partsOne = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && (aod::femtouniverseparticle::sign == confChargePart1) && (nabs(aod::femtouniverseparticle::eta) < confEta) && (aod::femtouniverseparticle::pt < confHPtPart1) && (aod::femtouniverseparticle::pt > confLPtPart1); + Partition partsOneMC = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kMCTruthTrack)) && (nabs(aod::femtouniverseparticle::eta) < confEta) && (aod::femtouniverseparticle::pt < confHPtPart1) && (aod::femtouniverseparticle::pt > confLPtPart1); + Partition partsOneMCReco = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kTrack)) && (aod::femtouniverseparticle::sign == confChargePart1) && (nabs(aod::femtouniverseparticle::eta) < confEta) && (aod::femtouniverseparticle::pt < confHPtPart1) && (aod::femtouniverseparticle::pt > confLPtPart1); + + /// Histogramming for particle 1 + FemtoUniverseParticleHisto trackHistoPartOnePos; + FemtoUniverseParticleHisto trackHistoPartOneNeg; + + /// Particle 2 (V0) + Configurable confV0PDGCodePartTwo{"confV0PDGCodePartTwo", 3122, "Particle 2 (V0) - PDG code"}; + ConfigurableAxis confV0TempFitVarBins{"confV0TempFitVarBins", {300, 0.95, 1.}, "V0: binning of the TempFitVar in the pT vs. TempFitVar plot"}; + ConfigurableAxis confV0TempFitVarpTBins{"confV0TempFitVarpTBins", {20, 0.5, 4.05}, "V0: pT binning of the pT vs. TempFitVar plot"}; + Configurable confV0Type1{"confV0Type1", 0, "select one of the V0s (lambda = 0, anti-lambda = 1, k0 = 2) for v0-v0 and Track-v0 combination"}; + Configurable confV0Type2{"confV0Type2", 0, "select one of the V0s (lambda = 0, anti-lambda = 1, k0 = 2) for v0-v0 combination"}; + Configurable confV0InvMassLowLimit{"confV0InvMassLowLimit", 1.10, "Lower limit of the V0 invariant mass"}; + Configurable confV0InvMassUpLimit{"confV0InvMassUpLimit", 1.13, "Upper limit of the V0 invariant mass"}; + ConfigurableAxis confChildTempFitVarBins{"confChildTempFitVarBins", {300, -0.15, 0.15}, "V0 child: binning of the TempFitVar in the pT vs. TempFitVar plot"}; + ConfigurableAxis confChildTempFitVarpTBins{"confChildTempFitVarpTBins", {20, 0.5, 4.05}, "V0 child: pT binning of the pT vs. TempFitVar plot"}; + Configurable confHPtPart2{"confHPtPart2", 4.0f, "higher limit for pt of particle 2"}; + Configurable confLPtPart2{"confLPtPart2", 0.3f, "lower limit for pt of particle 2"}; + Configurable confPDGCodePosChild{"confPDGCodePosChild", 2212, "Positive Child -- PDG code"}; + Configurable confPDGCodeNegChild{"confPDGCodeNegChild", 211, "Negative Child -- PDG code"}; + + /// Partition for particle 2 + Partition partsTwo = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kV0)) && (aod::femtouniverseparticle::pt < confHPtPart2) && (aod::femtouniverseparticle::pt > confLPtPart2); + Partition partsTwoMC = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kMCTruthTrack)) && (aod::femtouniverseparticle::pt < confHPtPart2) && (aod::femtouniverseparticle::pt > confLPtPart2); + Partition partsTwoMCReco = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kV0)) && (aod::femtouniverseparticle::pt < confHPtPart2) && (aod::femtouniverseparticle::pt > confLPtPart2); + + /// Histogramming for particle 2 + FemtoUniverseParticleHisto trackHistoPartTwo; + FemtoUniverseParticleHisto posChildHistos; + FemtoUniverseParticleHisto negChildHistos; + + FemtoUniverseParticleHisto trackHistoV0Type1; + FemtoUniverseParticleHisto posChildV0Type1; + FemtoUniverseParticleHisto negChildV0Type1; + FemtoUniverseParticleHisto trackHistoV0Type2; + FemtoUniverseParticleHisto posChildV0Type2; + FemtoUniverseParticleHisto negChildV0Type2; + + /// Histogramming for Event + FemtoUniverseEventHisto eventHisto; + + /// Correlation part + // Configurable confTrackChoicePartTwo{"confTrackChoicePartTwo", 1, "0:Proton, 1:Pion, 2:Kaon"}; //not used + Configurable confIsMC{"confIsMC", false, "Enable additional Histograms in the case of a MonteCarlo Run"}; + Configurable confUse3D{"confUse3D", false, "Enable three dimensional histogramms (to be used only for analysis with high statistics): k* vs mT vs multiplicity"}; + Configurable confUseCent{"confUseCent", false, "Use centrality in place of multiplicity"}; + ConfigurableAxis confMultBins{"confMultBins", {VARIABLE_WIDTH, 0.0f, 20.0f, 40.0f, 60.0f, 80.0f, 100.0f, 200.0f, 99999.f}, "Mixing bins - multiplicity"}; + ConfigurableAxis confVtxBins{"confVtxBins", {VARIABLE_WIDTH, -10.0f, -8.f, -6.f, -4.f, -2.f, 0.f, 2.f, 4.f, 6.f, 8.f, 10.f}, "Mixing bins - z-vertex"}; + Configurable confNEventsMix{"confNEventsMix", 5, "Number of events for mixing"}; + ConfigurableAxis confkstarBins{"confkstarBins", {1500, 0., 6.}, "binning kstar"}; + ConfigurableAxis confkTBins{"confkTBins", {150, 0., 9.}, "binning kT"}; + ConfigurableAxis confmTBins{"confmTBins", {225, 0., 7.5}, "binning mT"}; + Configurable confIsCPR{"confIsCPR", true, "Close Pair Rejection"}; + Configurable confCPRPlotPerRadii{"confCPRPlotPerRadii", false, "Plot CPR per radii"}; + Configurable confCPRdeltaPhiCutMax{"confCPRdeltaPhiCutMax", 0.0, "Delta Phi max cut for Close Pair Rejection"}; + Configurable confCPRdeltaPhiCutMin{"confCPRdeltaPhiCutMin", 0.0, "Delta Phi min cut for Close Pair Rejection"}; + Configurable confCPRdeltaEtaCutMax{"confCPRdeltaEtaCutMax", 0.0, "Delta Eta max cut for Close Pair Rejection"}; + Configurable confCPRdeltaEtaCutMin{"confCPRdeltaEtaCutMin", 0.0, "Delta Eta min cut for Close Pair Rejection"}; + Configurable confCPRChosenRadii{"confCPRChosenRadii", 0.80, "Delta Eta cut for Close Pair Rejection"}; + Configurable confPhiBins{"confPhiBins", 29, "Number of phi bins in deta dphi"}; + Configurable confEtaBins{"confEtaBins", 29, "Number of eta bins in deta dphi"}; + ConfigurableAxis confmTBins3D{"confmTBins3D", {VARIABLE_WIDTH, 1.02f, 1.14f, 1.20f, 1.26f, 1.38f, 1.56f, 1.86f, 4.50f}, "mT Binning for the 3Dimensional plot: k* vs multiplicity vs mT (set <> to true in order to use)"}; + ConfigurableAxis confMultBins3D{"confMultBins3D", {VARIABLE_WIDTH, 0.0f, 20.0f, 30.0f, 40.0f, 99999.0f}, "multiplicity Binning for the 3Dimensional plot: k* vs multiplicity vs mT (set <> to true in order to use)"}; + + // Efficiency + Configurable confLocalEfficiency{"confLocalEfficiency", "", "Local path to efficiency .root file"}; + + static constexpr unsigned int V0ChildTable[][2] = {{0, 1}, {1, 0}, {1, 1}}; // Table to select the V0 children + + FemtoUniverseContainer sameEventCont; + FemtoUniverseContainer mixedEventCont; + FemtoUniversePairCleaner pairCleaner; + FemtoUniversePairCleaner pairCleanerV0; + FemtoUniverseDetaDphiStar pairCloseRejection; + FemtoUniverseDetaDphiStar pairCloseRejectionV0; + + /// Histogram output + HistogramRegistry qaRegistry{"TrackQA", {}, OutputObjHandlingPolicy::AnalysisObject}; + HistogramRegistry resultRegistry{"Correlations", {}, OutputObjHandlingPolicy::AnalysisObject}; + HistogramRegistry registryMCtruth{"MCtruthHistos", {}, OutputObjHandlingPolicy::AnalysisObject, false, true}; + HistogramRegistry registryMCreco{"MCrecoHistos", {}, OutputObjHandlingPolicy::AnalysisObject, false, true}; + HistogramRegistry mixQaRegistry{"mixQaRegistry", {}, OutputObjHandlingPolicy::AnalysisObject}; + + HistogramRegistry thetaRegistry{"ThetaQA", {}, OutputObjHandlingPolicy::AnalysisObject}; + + std::unique_ptr plocalEffFile; + std::unique_ptr plocalEffp1; + std::unique_ptr plocalEffp2; + + bool isNSigmaCombined(float mom, float nsigmaTPCParticle, float nsigmaTOFParticle) + { + if (mom <= confmom) { + return (std::abs(nsigmaTPCParticle) < confNsigmaTPCParticle); + } else { + return (std::hypot(nsigmaTOFParticle, nsigmaTPCParticle) < confNsigmaCombinedParticle); + } + } + + bool invMLambda(float invMassLambda, float invMassAntiLambda) + { + if ((invMassLambda < confV0InvMassLowLimit || invMassLambda > confV0InvMassUpLimit) && (invMassAntiLambda < confV0InvMassLowLimit || invMassAntiLambda > confV0InvMassUpLimit)) { + return false; + } + return true; + } + + bool isNSigmaTPC(float nsigmaTPCParticle) + { + if (std::abs(nsigmaTPCParticle) < confNsigmaTPCParticle) { + return true; + } else { + return false; + } + } + + template + bool isParticleCombined(const T& part, int id) + { + const float tpcNSigmas[3] = {unPackInTable(part.tpcNSigmaStorePr()), unPackInTable(part.tpcNSigmaStorePi()), unPackInTable(part.tpcNSigmaStoreKa())}; + // const float tofNSigmas[3] = {part.tofNSigmaPr(), part.tofNSigmaPi(), part.tofNSigmaKa()}; + const float tofNSigmas[3] = {unPackInTable(part.tofNSigmaStorePr()), unPackInTable(part.tofNSigmaStorePi()), unPackInTable(part.tofNSigmaStoreKa())}; + + return isNSigmaCombined(part.p(), tpcNSigmas[id], tofNSigmas[id]); + } + + template + bool isParticleTPC(const T& part, int id) + { + const float tpcNSigmas[3] = {unPackInTable(part.tpcNSigmaStorePr()), unPackInTable(part.tpcNSigmaStorePi()), unPackInTable(part.tpcNSigmaStoreKa())}; + + return isNSigmaTPC(tpcNSigmas[id]); + } + + void init(InitContext&) + { + eventHisto.init(&qaRegistry); + qaRegistry.add("Tracks_pos/nSigmaTPC", "; #it{p} (GeV/#it{c}); n#sigma_{TPC}", kTH2F, {{100, 0, 10}, {200, -4.975, 5.025}}); + qaRegistry.add("Tracks_pos/nSigmaTOF", "; #it{p} (GeV/#it{c}); n#sigma_{TOF}", kTH2F, {{100, 0, 10}, {200, -4.975, 5.025}}); + qaRegistry.add("Tracks_neg/nSigmaTPC", "; #it{p} (GeV/#it{c}); n#sigma_{TPC}", kTH2F, {{100, 0, 10}, {200, -4.975, 5.025}}); + qaRegistry.add("Tracks_neg/nSigmaTOF", "; #it{p} (GeV/#it{c}); n#sigma_{TOF}", kTH2F, {{100, 0, 10}, {200, -4.975, 5.025}}); + trackHistoPartOnePos.init(&qaRegistry, confTrkTempFitVarpTBins, confTrkTempFitVarBins, confIsMC, confTrkPDGCodePartOne); + trackHistoPartOneNeg.init(&qaRegistry, confTrkTempFitVarpTBins, confTrkTempFitVarBins, confIsMC, confTrkPDGCodePartOne); + trackHistoPartTwo.init(&qaRegistry, confV0TempFitVarpTBins, confV0TempFitVarBins, confIsMC, confV0PDGCodePartTwo, true); + posChildHistos.init(&qaRegistry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, 0, true); + negChildHistos.init(&qaRegistry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, 0, true); + + trackHistoV0Type1.init(&qaRegistry, confV0TempFitVarpTBins, confV0TempFitVarBins, confIsMC, confV0PDGCodePartTwo, true, "V0Type1"); + posChildV0Type1.init(&qaRegistry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, 0, true, "posChildV0Type1"); + negChildV0Type1.init(&qaRegistry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, 0, true, "negChildV0Type1"); + trackHistoV0Type2.init(&qaRegistry, confV0TempFitVarpTBins, confV0TempFitVarBins, confIsMC, confV0PDGCodePartTwo, true, "V0Type2"); + posChildV0Type2.init(&qaRegistry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, 0, true, "posChildV0Type2"); + negChildV0Type2.init(&qaRegistry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, 0, true, "negChildV0Type2"); + + mixQaRegistry.add("MixingQA/hMECollisionBins", ";bin;Entries", kTH1F, {{120, -0.5, 119.5}}); + + // MC truth + registryMCtruth.add("plus/MCtruthLambda", "MC truth Lambdas;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + registryMCtruth.add("minus/MCtruthLambda", "MC truth Lambdas;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + + registryMCtruth.add("plus/MCtruthAllPt", "MC truth all;#it{p}_{T} (GeV/c); #eta", {HistType::kTH1F, {{500, 0, 5}}}); + registryMCtruth.add("minus/MCtruthAllPt", "MC truth all;#it{p}_{T} (GeV/c); #eta", {HistType::kTH1F, {{500, 0, 5}}}); + + registryMCtruth.add("plus/MCtruthPi", "MC truth pions;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + registryMCtruth.add("plus/MCtruthPr", "MC truth protons;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + + registryMCtruth.add("minus/MCtruthPi", "MC truth pions;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + registryMCtruth.add("minus/MCtruthPr", "MC truth protons;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + + registryMCtruth.add("plus/MCtruthPiPt", "MC truth pions;#it{p}_{T} (GeV/c)", {HistType::kTH1F, {{500, 0, 5}}}); + registryMCtruth.add("plus/MCtruthPrPt", "MC truth protons;#it{p}_{T} (GeV/c)", {HistType::kTH1F, {{500, 0, 5}}}); + registryMCtruth.add("minus/MCtruthPiPt", "MC truth pions;#it{p}_{T} (GeV/c)", {HistType::kTH1F, {{500, 0, 5}}}); + registryMCtruth.add("minus/MCtruthPrPt", "MC truth protons;#it{p}_{T} (GeV/c)", {HistType::kTH1F, {{500, 0, 5}}}); + + // MC reco + registryMCreco.add("plus/MCrecoLambda", "MC reco Lambdas;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + registryMCreco.add("plus/MCrecoLambdaChildPr", "MC reco Lambdas;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + registryMCreco.add("plus/MCrecoLambdaChildPi", "MC reco Lambdas;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + registryMCreco.add("minus/MCrecoLambda", "MC reco Lambdas;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + registryMCreco.add("minus/MCrecoLambdaChildPr", "MC reco Lambdas;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + registryMCreco.add("minus/MCrecoLambdaChildPi", "MC reco Lambdas;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + + registryMCreco.add("plus/MCrecoAllPt", "MC reco all;#it{p}_{T} (GeV/c); #eta", {HistType::kTH1F, {{500, 0, 5}}}); + registryMCreco.add("minus/MCrecoAllPt", "MC reco all;#it{p}_{T} (GeV/c); #eta", {HistType::kTH1F, {{500, 0, 5}}}); + + registryMCreco.add("plus/MCrecoPi", "MC reco pions;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + registryMCreco.add("plus/MCrecoPr", "MC reco protons;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + + registryMCreco.add("minus/MCrecoPi", "MC reco pions;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + registryMCreco.add("minus/MCrecoPr", "MC reco protons;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); + + registryMCreco.add("plus/MCrecoPiPt", "MC reco pions;#it{p}_{T} (GeV/c)", {HistType::kTH1F, {{500, 0, 5}}}); + registryMCreco.add("plus/MCrecoPrPt", "MC reco protons;#it{p}_{T} (GeV/c)", {HistType::kTH1F, {{500, 0, 5}}}); + registryMCreco.add("minus/MCrecoPiPt", "MC reco pions;#it{p}_{T} (GeV/c)", {HistType::kTH1F, {{500, 0, 5}}}); + registryMCreco.add("minus/MCrecoPrPt", "MC reco protons;#it{p}_{T} (GeV/c)", {HistType::kTH1F, {{500, 0, 5}}}); + + // Helicity angle + thetaRegistry.add("Theta/hTheta", " ; p (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/PositiveChild/hThetaPt", " ; p_{T} (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/PositiveChild/hThetaEta", " ; #eta; cos(#theta)", kTH2F, {{100, -1, 1}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/PositiveChild/hThetaPhi", " ; #phi; cos(#theta)", kTH2F, {{100, -1, 7}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/NegativeChild/hThetaPt", " ; p_{T} (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/NegativeChild/hThetaEta", " ; #eta; cos(#theta)", kTH2F, {{100, -1, 1}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/NegativeChild/hThetaPhi", " ; #phi; cos(#theta)", kTH2F, {{100, -1, 7}, {110, -1.1, 1.1}}); + + sameEventCont.init(&resultRegistry, confkstarBins, confMultBins, confkTBins, confmTBins, confMultBins3D, confmTBins3D, confEtaBins, confPhiBins, confIsMC, confUse3D); + sameEventCont.setPDGCodes(confTrkPDGCodePartOne, confV0PDGCodePartTwo); + mixedEventCont.init(&resultRegistry, confkstarBins, confMultBins, confkTBins, confmTBins, confMultBins3D, confmTBins3D, confEtaBins, confPhiBins, confIsMC, confUse3D); + mixedEventCont.setPDGCodes(confTrkPDGCodePartOne, confV0PDGCodePartTwo); + + pairCleaner.init(&qaRegistry); + pairCleanerV0.init(&qaRegistry); + if (confIsCPR.value) { + pairCloseRejection.init(&resultRegistry, &qaRegistry, confCPRdeltaPhiCutMin.value, confCPRdeltaPhiCutMax.value, confCPRdeltaEtaCutMin.value, confCPRdeltaEtaCutMax.value, confCPRChosenRadii.value, confCPRPlotPerRadii.value); + pairCloseRejectionV0.init(&resultRegistry, &qaRegistry, confCPRdeltaPhiCutMin.value, confCPRdeltaPhiCutMax.value, confCPRdeltaEtaCutMin.value, confCPRdeltaEtaCutMax.value, confCPRChosenRadii.value, confCPRPlotPerRadii.value); + } + + if (!confLocalEfficiency.value.empty()) { + plocalEffFile = std::unique_ptr(TFile::Open(confLocalEfficiency.value.c_str(), "read")); + if (!plocalEffFile || plocalEffFile.get()->IsZombie()) + LOGF(fatal, "Could not load efficiency histogram from %s", confLocalEfficiency.value.c_str()); + if (doprocessSameEvent || doprocessMixedEvent) { + plocalEffp1 = (confChargePart1 > 0) ? std::unique_ptr(plocalEffFile.get()->Get("PrPlus")) : std::unique_ptr(plocalEffFile.get()->Get("PrMinus")); // note: works only for protons for now + plocalEffp2 = (confV0Type1 == 0) ? std::unique_ptr(plocalEffFile.get()->Get("Lambda")) : std::unique_ptr(plocalEffFile.get()->Get("AntiLambda")); + LOGF(info, "Loaded efficiency histograms for track-V0."); + } else if (doprocessSameEventV0 || doprocessMixedEventV0) { + plocalEffp1 = (confV0Type1 == 0) ? std::unique_ptr(plocalEffFile.get()->Get("Lambda")) : std::unique_ptr(plocalEffFile.get()->Get("AntiLambda")); + plocalEffp2 = (confV0Type2 == 0) ? std::unique_ptr(plocalEffFile.get()->Get("Lambda")) : std::unique_ptr(plocalEffFile.get()->Get("AntiLambda")); + LOGF(info, "Loaded efficiency histograms for V0-V0."); + } + } + } + /// This function processes the same event for track - V0 + template + void doSameEvent(FilteredFDCollision const& col, PartType const& parts, PartitionType& groupPartsOne, PartitionType& groupPartsTwo, [[maybe_unused]] MCParticles mcParts = nullptr) + { + const auto& magFieldTesla = col.magField(); + + const int multCol = confUseCent ? col.multV0M() : col.multNtr(); + + eventHisto.fillQA(col); + + /// Histogramming same event + for (const auto& part : groupPartsTwo) { + if (!invMLambda(part.mLambda(), part.mAntiLambda())) + continue; + const auto& posChild = parts.iteratorAt(part.index() - 2); + const auto& negChild = parts.iteratorAt(part.index() - 1); + /// Daughters that do not pass this condition are not selected + if (!isParticleTPC(posChild, V0ChildTable[confV0Type1][0]) || !isParticleTPC(negChild, V0ChildTable[confV0Type1][1])) + continue; + + auto posChildMass = pdg->Mass(confPDGCodePosChild); + auto negChildMass = pdg->Mass(confPDGCodeNegChild); + auto posChildBoosted = FemtoUniverseMath::boostPRF(posChild, posChildMass, negChild, negChildMass); + auto cosineTheta = (posChildBoosted.Px() * part.px() + posChildBoosted.Py() * part.py() + posChildBoosted.Pz() * part.pz()) / (posChildBoosted.P() * part.p()); + + trackHistoPartTwo.fillQA(part); + posChildHistos.fillQA(posChild); + negChildHistos.fillQA(negChild); + + thetaRegistry.fill(HIST("Theta/hTheta"), part.p(), cosineTheta); + thetaRegistry.fill(HIST("Theta/PositiveChild/hThetaPt"), posChild.pt(), cosineTheta); + thetaRegistry.fill(HIST("Theta/PositiveChild/hThetaEta"), posChild.eta(), cosineTheta); + thetaRegistry.fill(HIST("Theta/PositiveChild/hThetaPhi"), posChild.phi(), cosineTheta); + thetaRegistry.fill(HIST("Theta/NegativeChild/hThetaPt"), negChild.pt(), cosineTheta); + thetaRegistry.fill(HIST("Theta/NegativeChild/hThetaEta"), negChild.eta(), cosineTheta); + thetaRegistry.fill(HIST("Theta/NegativeChild/hThetaPhi"), negChild.phi(), cosineTheta); + } + + for (const auto& part : groupPartsOne) { + /// PID plot for particle 1 + const float tpcNSigmas[3] = {unPackInTable(part.tpcNSigmaStorePr()), unPackInTable(part.tpcNSigmaStorePi()), unPackInTable(part.tpcNSigmaStoreKa())}; + const float tofNSigmas[3] = {unPackInTable(part.tofNSigmaStorePr()), unPackInTable(part.tofNSigmaStorePi()), unPackInTable(part.tofNSigmaStoreKa())}; + + if (!isNSigmaCombined(part.p(), tpcNSigmas[confTrackChoicePartOne], tofNSigmas[confTrackChoicePartOne])) + continue; + if (part.sign() > 0) { + qaRegistry.fill(HIST("Tracks_pos/nSigmaTPC"), part.p(), tpcNSigmas[confTrackChoicePartOne]); + qaRegistry.fill(HIST("Tracks_pos/nSigmaTOF"), part.p(), tofNSigmas[confTrackChoicePartOne]); + trackHistoPartOnePos.fillQA(part); + } else if (part.sign() < 0) { + qaRegistry.fill(HIST("Tracks_neg/nSigmaTPC"), part.p(), tpcNSigmas[confTrackChoicePartOne]); + qaRegistry.fill(HIST("Tracks_neg/nSigmaTOF"), part.p(), tofNSigmas[confTrackChoicePartOne]); + trackHistoPartOneNeg.fillQA(part); + } + } + + /// Now build the combinations + for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsTwo))) { + // Lambda invariant mass cut + if (!invMLambda(p2.mLambda(), p2.mAntiLambda())) + continue; + /// PID using stored binned nsigma + if (!isParticleCombined(p1, confTrackChoicePartOne)) + continue; + // track cleaning + if (!pairCleaner.isCleanPair(p1, p2, parts)) { + continue; + } + if (confIsCPR.value) { + if (pairCloseRejection.isClosePair(p1, p2, parts, magFieldTesla, femto_universe_container::EventType::same)) { + continue; + } + } + const auto& posChild = parts.iteratorAt(p2.index() - 2); + const auto& negChild = parts.iteratorAt(p2.index() - 1); + + /// Daughters that do not pass this condition are not selected + if (!isParticleTPC(posChild, V0ChildTable[confV0Type1][0]) || !isParticleTPC(negChild, V0ChildTable[confV0Type1][1])) + continue; + + float weight = 1.0f; + if (plocalEffp1) + weight = plocalEffp1.get()->GetBinContent(plocalEffp1->FindBin(p1.pt(), p1.eta())) * plocalEffp2.get()->GetBinContent(plocalEffp2->FindBin(p2.pt(), p2.eta())); + if constexpr (std::is_same::value) + sameEventCont.setPair(p1, p2, multCol, confUse3D, weight); + else + sameEventCont.setPair(p1, p2, multCol, confUse3D, weight); + } + } + + void processSameEvent(FilteredFDCollision const& col, FemtoFullParticles const& parts) + { + auto groupPartsOne = partsOne->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); + auto groupPartsTwo = partsTwo->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); + doSameEvent(col, parts, groupPartsOne, groupPartsTwo); + } + PROCESS_SWITCH(FemtoUniversePairTaskTrackV0Helicity, processSameEvent, "Enable processing same event for track - V0", false); + + void processSameEventMCReco(FilteredFDCollision const& col, FemtoRecoParticles const& parts, aod::FdMCParticles const& mcparts) + { + auto groupPartsOne = partsOneMCReco->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); + auto groupPartsTwo = partsTwoMCReco->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); + doSameEvent(col, parts, groupPartsOne, groupPartsTwo, mcparts); + } + PROCESS_SWITCH(FemtoUniversePairTaskTrackV0Helicity, processSameEventMCReco, "Enable processing same event for track - V0 MC Reco", false); + + /// This function processes the same event for V0 - V0 + void processSameEventV0(FilteredFDCollision const& col, FemtoFullParticles const& parts) + { + const auto& magFieldTesla = col.magField(); + + auto groupPartsTwo = partsTwo->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); + const int multCol = confUseCent ? col.multV0M() : col.multNtr(); + + eventHisto.fillQA(col); + + /// Histogramming same event + for (const auto& part : groupPartsTwo) { + if (!invMLambda(part.mLambda(), part.mAntiLambda())) + continue; + const auto& posChild = parts.iteratorAt(part.index() - 2); + const auto& negChild = parts.iteratorAt(part.index() - 1); + + /// Check daughters of first V0 particle + if (isParticleTPC(posChild, V0ChildTable[confV0Type1][0]) && isParticleTPC(negChild, V0ChildTable[confV0Type1][1])) { + trackHistoV0Type1.fillQABase(part, HIST("V0Type1")); + posChildV0Type1.fillQABase(posChild, HIST("posChildV0Type1")); + negChildV0Type1.fillQABase(negChild, HIST("negChildV0Type1")); + } + /// Check daughters of second V0 particle + if (isParticleTPC(posChild, V0ChildTable[confV0Type2][0]) && isParticleTPC(negChild, V0ChildTable[confV0Type2][1])) { + trackHistoV0Type2.fillQABase(part, HIST("V0Type2")); + posChildV0Type2.fillQABase(posChild, HIST("posChildV0Type2")); + negChildV0Type2.fillQABase(negChild, HIST("negChildV0Type2")); + } + } + + auto pairProcessFunc = [&](auto& p1, auto& p2) -> void { + // Lambda invariant mass cut for p1 + if (!invMLambda(p1.mLambda(), p1.mAntiLambda())) + return; + // Lambda invariant mass cut for p2 + if (!invMLambda(p2.mLambda(), p2.mAntiLambda())) + return; + // track cleaning + if (!pairCleanerV0.isCleanPair(p1, p2, parts)) { + return; + } + if (confIsCPR.value) { + if (pairCloseRejectionV0.isClosePair(p1, p2, parts, magFieldTesla, femto_universe_container::EventType::same)) { + return; + } + } + const auto& posChild1 = parts.iteratorAt(p1.index() - 2); + const auto& negChild1 = parts.iteratorAt(p1.index() - 1); + /// Daughters that do not pass this condition are not selected + if (!isParticleTPC(posChild1, V0ChildTable[confV0Type1][0]) || !isParticleTPC(negChild1, V0ChildTable[confV0Type1][1])) + return; + + const auto& posChild2 = parts.iteratorAt(p2.index() - 2); + const auto& negChild2 = parts.iteratorAt(p2.index() - 1); + /// Daughters that do not pass this condition are not selected + if (!isParticleTPC(posChild2, V0ChildTable[confV0Type2][0]) || !isParticleTPC(negChild2, V0ChildTable[confV0Type2][1])) + return; + + sameEventCont.setPair(p1, p2, multCol, confUse3D); + }; + if (confV0Type1 == confV0Type2) { + /// Now build the combinations for identical V0s + for (const auto& [p1, p2] : combinations(CombinationsStrictlyUpperIndexPolicy(groupPartsTwo, groupPartsTwo))) { + pairProcessFunc(p1, p2); + } + } else { + /// Now build the combinations for not identical identical V0s + for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsTwo, groupPartsTwo))) { + pairProcessFunc(p1, p2); + } + } + } + + PROCESS_SWITCH(FemtoUniversePairTaskTrackV0Helicity, processSameEventV0, "Enable processing same event for V0 - V0", false); + + /// This function processes MC same events for Track - V0 + void processMCSameEvent(FilteredFDCollision const& col, FemtoFullParticles const& parts) + { + const auto& magFieldTesla = col.magField(); + + auto groupPartsOne = partsOneMC->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); + auto groupPartsTwo = partsTwoMC->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); + const int multCol = confUseCent ? col.multV0M() : col.multNtr(); + + eventHisto.fillQA(col); + + /// Histogramming same event + for (const auto& part : groupPartsTwo) { + int pdgCode = static_cast(part.pidCut()); + if ((confV0Type1 == 0 && pdgCode != 3122) || (confV0Type1 == 1 && pdgCode != -3122)) + continue; + trackHistoPartTwo.fillQA(part); + } + + for (const auto& part : groupPartsOne) { + int pdgCode = static_cast(part.pidCut()); + if (pdgCode != confTrkPDGCodePartOne) + continue; + const auto& pdgParticle = pdgMC->GetParticle(pdgCode); + if (!pdgParticle) { + continue; + } + /// PID plot for particle 1 + if (pdgParticle->Charge() > 0.0) { + trackHistoPartOnePos.fillQA(part); + } else if (pdgParticle->Charge() < 0.0) { + trackHistoPartOneNeg.fillQA(part); + } + } + + /// Now build the combinations + for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsTwo))) { + if (static_cast(p1.pidCut()) != confTrkPDGCodePartOne) + continue; + int pdgCode2 = static_cast(p2.pidCut()); + if ((confV0Type1 == 0 && pdgCode2 != 3122) || (confV0Type1 == 1 && pdgCode2 != -3122)) + continue; + // track cleaning + if (confIsCPR.value) { + if (pairCloseRejection.isClosePair(p1, p2, parts, magFieldTesla, femto_universe_container::EventType::same)) { + continue; + } + } + sameEventCont.setPair(p1, p2, multCol, confUse3D); + } + } + + PROCESS_SWITCH(FemtoUniversePairTaskTrackV0Helicity, processMCSameEvent, "Enable processing same event for MC truth track - V0", false); + + /// This function processes MC same events for V0 - V0 + void processMCSameEventV0(FilteredFDCollision const& col, FemtoFullParticles const& /*parts*/) + { + auto groupPartsTwo = partsTwoMC->sliceByCached(aod::femtouniverseparticle::fdCollisionId, col.globalIndex(), cache); + const int multCol = confUseCent ? col.multV0M() : col.multNtr(); + + eventHisto.fillQA(col); + + /// Histogramming same event + for (const auto& part : groupPartsTwo) { + int pdgCode = static_cast(part.pidCut()); + if ((confV0Type1 == 0 && pdgCode != 3122) || (confV0Type1 == 1 && pdgCode != -3122)) + continue; + trackHistoPartTwo.fillQA(part); + } + + auto pairProcessFunc = [&](auto& p1, auto& p2) -> void { + int pdgCode1 = static_cast(p1.pidCut()); + if ((confV0Type1 == 0 && pdgCode1 != 3122) || (confV0Type1 == 1 && pdgCode1 != -3122)) + return; + int pdgCode2 = static_cast(p2.pidCut()); + if ((confV0Type2 == 0 && pdgCode2 != 3122) || (confV0Type2 == 1 && pdgCode2 != -3122)) + return; + sameEventCont.setPair(p1, p2, multCol, confUse3D); + }; + /// Now build the combinations + if (confV0Type1 == confV0Type2) { + /// Now build the combinations for identical V0s + for (const auto& [p1, p2] : combinations(CombinationsStrictlyUpperIndexPolicy(groupPartsTwo, groupPartsTwo))) { + pairProcessFunc(p1, p2); + } + } else { + /// Now build the combinations for not identical identical V0s + for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsTwo, groupPartsTwo))) { + pairProcessFunc(p1, p2); + } + } + } + + PROCESS_SWITCH(FemtoUniversePairTaskTrackV0Helicity, processMCSameEventV0, "Enable processing same event for MC truth V0 - V0", false); + + /// This function processes the mixed event for track - V0 + template + void doMixedEvent(FilteredFDCollisions const& cols, PartType const& parts, PartitionType& partitionOne, PartitionType& partitionTwo, [[maybe_unused]] MCParticles mcParts = nullptr) + { + ColumnBinningPolicy colBinningMult{{confVtxBins, confMultBins}, true}; + ColumnBinningPolicy colBinningCent{{confVtxBins, confMultBins}, true}; + + auto mixedCollProcessFunc = [&](auto& collision1, auto& collision2) -> void { + const int multCol = confUseCent ? collision1.multV0M() : collision1.multNtr(); + + auto groupPartsOne = partitionOne->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision1.globalIndex(), cache); + auto groupPartsTwo = partitionTwo->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision2.globalIndex(), cache); + + const auto& magFieldTesla1 = collision1.magField(); + const auto& magFieldTesla2 = collision2.magField(); + + if (magFieldTesla1 != magFieldTesla2) { + return; + } + + for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsTwo))) { + // Lambda invariant mass cut + if (!invMLambda(p2.mLambda(), p2.mAntiLambda())) + continue; + /// PID using stored binned nsigma + if (!isParticleCombined(p1, confTrackChoicePartOne)) + continue; + + const auto& posChild = parts.iteratorAt(p2.globalIndex() - 2); + const auto& negChild = parts.iteratorAt(p2.globalIndex() - 1); + /// Daughters that do not pass this condition are not selected + if (!isParticleTPC(posChild, V0ChildTable[confV0Type1][0]) || !isParticleTPC(negChild, V0ChildTable[confV0Type1][1])) + continue; + + // track cleaning + if (!pairCleaner.isCleanPair(p1, p2, parts)) { + continue; + } + if (confIsCPR.value) { + if (pairCloseRejection.isClosePair(p1, p2, parts, magFieldTesla1, femto_universe_container::EventType::mixed)) { + continue; + } + } + float weight = 1.0f; + if (plocalEffp1) + weight = plocalEffp1.get()->GetBinContent(plocalEffp1->FindBin(p1.pt(), p1.eta())) * plocalEffp2.get()->GetBinContent(plocalEffp2->FindBin(p2.pt(), p2.eta())); + + if constexpr (std::is_same::value) + mixedEventCont.setPair(p1, p2, multCol, confUse3D, weight); + else + mixedEventCont.setPair(p1, p2, multCol, confUse3D, weight); + } + }; + + if (confUseCent) { + for (const auto& [collision1, collision2] : soa::selfCombinations(colBinningCent, confNEventsMix, -1, cols, cols)) { + mixedCollProcessFunc(collision1, collision2); + mixQaRegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinningCent.getBin({collision1.posZ(), collision1.multV0M()})); + } + } else { + for (const auto& [collision1, collision2] : soa::selfCombinations(colBinningMult, confNEventsMix, -1, cols, cols)) { + mixedCollProcessFunc(collision1, collision2); + mixQaRegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinningMult.getBin({collision1.posZ(), collision1.multNtr()})); + } + } + } + + void processMixedEvent(FilteredFDCollisions const& cols, FemtoFullParticles const& parts) + { + doMixedEvent(cols, parts, partsOne, partsTwo); + } + PROCESS_SWITCH(FemtoUniversePairTaskTrackV0Helicity, processMixedEvent, "Enable processing mixed event for track - V0", false); + + void processMixedEventMCReco(FilteredFDCollisions const& cols, FemtoRecoParticles const& parts, aod::FdMCParticles const& mcparts) + { + doMixedEvent(cols, parts, partsOneMCReco, partsTwoMCReco, mcparts); + } + PROCESS_SWITCH(FemtoUniversePairTaskTrackV0Helicity, processMixedEventMCReco, "Enable processing mixed event for track - V0 for MC Reco", false); + + /// This function processes the mixed event for V0 - V0 + void processMixedEventV0(FilteredFDCollisions const& cols, FemtoFullParticles const& parts) + { + ColumnBinningPolicy colBinningMult{{confVtxBins, confMultBins}, true}; + ColumnBinningPolicy colBinningCent{{confVtxBins, confMultBins}, true}; + + auto mixedCollProcessFunc = [&](auto& collision1, auto& collision2) -> void { + const int multCol = confUseCent ? collision1.multV0M() : collision1.multNtr(); + + auto groupPartsOne = partsTwo->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision1.globalIndex(), cache); + auto groupPartsTwo = partsTwo->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision2.globalIndex(), cache); + + const auto& magFieldTesla1 = collision1.magField(); + const auto& magFieldTesla2 = collision2.magField(); + + if (magFieldTesla1 != magFieldTesla2) { + return; + } + + for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsTwo))) { + // Lambda invariant mass cut for p1 + if (!invMLambda(p1.mLambda(), p1.mAntiLambda())) { + continue; + } + // Lambda invariant mass cut for p2 + if (!invMLambda(p2.mLambda(), p2.mAntiLambda())) { + continue; + } + + const auto& posChild1 = parts.iteratorAt(p1.globalIndex() - 2); + const auto& negChild1 = parts.iteratorAt(p1.globalIndex() - 1); + /// Daughters that do not pass this condition are not selected + if (!isParticleTPC(posChild1, V0ChildTable[confV0Type1][0]) || !isParticleTPC(negChild1, V0ChildTable[confV0Type1][1])) + continue; + + const auto& posChild2 = parts.iteratorAt(p2.globalIndex() - 2); + const auto& negChild2 = parts.iteratorAt(p2.globalIndex() - 1); + /// Daughters that do not pass this condition are not selected + if (!isParticleTPC(posChild2, V0ChildTable[confV0Type2][0]) || !isParticleTPC(negChild2, V0ChildTable[confV0Type2][1])) + continue; + + // track cleaning + if (!pairCleanerV0.isCleanPair(p1, p2, parts)) { + continue; + } + if (confIsCPR.value) { + if (pairCloseRejectionV0.isClosePair(p1, p2, parts, magFieldTesla1, femto_universe_container::EventType::mixed)) { + continue; + } + } + mixedEventCont.setPair(p1, p2, multCol, confUse3D); + } + }; + + if (confUseCent) { + for (const auto& [collision1, collision2] : soa::selfCombinations(colBinningCent, confNEventsMix, -1, cols, cols)) { + mixedCollProcessFunc(collision1, collision2); + mixQaRegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinningCent.getBin({collision1.posZ(), collision1.multV0M()})); + } + } else { + for (const auto& [collision1, collision2] : soa::selfCombinations(colBinningMult, confNEventsMix, -1, cols, cols)) { + mixedCollProcessFunc(collision1, collision2); + mixQaRegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinningMult.getBin({collision1.posZ(), collision1.multNtr()})); + } + } + } + PROCESS_SWITCH(FemtoUniversePairTaskTrackV0Helicity, processMixedEventV0, "Enable processing mixed events for V0 - V0", false); + + /// This function processes MC mixed events for Track - V0 + void processMCMixedEvent(FilteredFDCollisions const& cols, FemtoFullParticles const& parts) + { + ColumnBinningPolicy colBinningMult{{confVtxBins, confMultBins}, true}; + ColumnBinningPolicy colBinningCent{{confVtxBins, confMultBins}, true}; + + auto mixedCollProcessFunc = [&](auto& collision1, auto& collision2) -> void { + const int multCol = confUseCent ? collision1.multV0M() : collision1.multNtr(); + + auto groupPartsOne = partsOneMC->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision1.globalIndex(), cache); + auto groupPartsTwo = partsTwoMC->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision2.globalIndex(), cache); + + const auto& magFieldTesla1 = collision1.magField(); + const auto& magFieldTesla2 = collision2.magField(); + + if (magFieldTesla1 != magFieldTesla2) { + return; + } + for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsTwo))) { + if (static_cast(p1.pidCut()) != confTrkPDGCodePartOne) + continue; + int pdgCode2 = static_cast(p2.pidCut()); + if ((confV0Type1 == 0 && pdgCode2 != 3122) || (confV0Type1 == 1 && pdgCode2 != -3122)) + continue; + if (confIsCPR.value) { + if (pairCloseRejection.isClosePair(p1, p2, parts, magFieldTesla1, femto_universe_container::EventType::mixed)) { + continue; + } + } + mixedEventCont.setPair(p1, p2, multCol, confUse3D); + } + }; + + if (confUseCent) { + for (const auto& [collision1, collision2] : soa::selfCombinations(colBinningCent, confNEventsMix, -1, cols, cols)) { + mixedCollProcessFunc(collision1, collision2); + mixQaRegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinningCent.getBin({collision1.posZ(), collision1.multV0M()})); + } + } else { + for (const auto& [collision1, collision2] : soa::selfCombinations(colBinningMult, confNEventsMix, -1, cols, cols)) { + mixedCollProcessFunc(collision1, collision2); + mixQaRegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinningMult.getBin({collision1.posZ(), collision1.multNtr()})); + } + } + } + + PROCESS_SWITCH(FemtoUniversePairTaskTrackV0Helicity, processMCMixedEvent, "Enable processing mixed events for MC truth track - V0", false); + + /// This function processes MC mixed events for V0 - V0 + void processMCMixedEventV0(FilteredFDCollisions const& cols, FemtoFullParticles const& /*parts*/) + { + ColumnBinningPolicy colBinningMult{{confVtxBins, confMultBins}, true}; + ColumnBinningPolicy colBinningCent{{confVtxBins, confMultBins}, true}; + + auto mixedCollProcessFunc = [&](auto& collision1, auto& collision2) -> void { + const int multCol = confUseCent ? collision1.multV0M() : collision1.multNtr(); + + auto groupPartsOne = partsTwoMC->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision1.globalIndex(), cache); + auto groupPartsTwo = partsTwoMC->sliceByCached(aod::femtouniverseparticle::fdCollisionId, collision2.globalIndex(), cache); + + for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsTwo))) { + int pdgCode1 = static_cast(p1.pidCut()); + if ((confV0Type1 == 0 && pdgCode1 != 3122) || (confV0Type1 == 1 && pdgCode1 != -3122)) + continue; + int pdgCode2 = static_cast(p2.pidCut()); + if ((confV0Type2 == 0 && pdgCode2 != 3122) || (confV0Type2 == 1 && pdgCode2 != -3122)) + continue; + mixedEventCont.setPair(p1, p2, multCol, confUse3D); + } + }; + + if (confUseCent) { + for (const auto& [collision1, collision2] : soa::selfCombinations(colBinningCent, confNEventsMix, -1, cols, cols)) { + mixedCollProcessFunc(collision1, collision2); + mixQaRegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinningCent.getBin({collision1.posZ(), collision1.multV0M()})); + } + } else { + for (const auto& [collision1, collision2] : soa::selfCombinations(colBinningMult, confNEventsMix, -1, cols, cols)) { + mixedCollProcessFunc(collision1, collision2); + mixQaRegistry.fill(HIST("MixingQA/hMECollisionBins"), colBinningMult.getBin({collision1.posZ(), collision1.multNtr()})); + } + } + } + + PROCESS_SWITCH(FemtoUniversePairTaskTrackV0Helicity, processMCMixedEventV0, "Enable processing mixed events for MC truth V0 - V0", false); + ///--------------------------------------------MC-------------------------------------------------/// + + /// This function fills MC truth particles from derived MC table + void processMCTruth(aod::FDParticles const& parts) + { + for (const auto& part : parts) { + if (part.partType() != uint8_t(aod::femtouniverseparticle::ParticleType::kMCTruthTrack)) + continue; + + int pdgCode = static_cast(part.pidCut()); + const auto& pdgParticle = pdgMC->GetParticle(pdgCode); + if (!pdgParticle) { + continue; + } + + if (pdgCode == 3122) { + registryMCtruth.fill(HIST("plus/MCtruthLambda"), part.pt(), part.eta()); + continue; + } else if (pdgCode == -3122) { + registryMCtruth.fill(HIST("minus/MCtruthLambda"), part.pt(), part.eta()); + continue; + } + + if (pdgParticle->Charge() > 0.0) { + registryMCtruth.fill(HIST("plus/MCtruthAllPt"), part.pt()); + } + if (pdgCode == 211) { + registryMCtruth.fill(HIST("plus/MCtruthPi"), part.pt(), part.eta()); + registryMCtruth.fill(HIST("plus/MCtruthPiPt"), part.pt()); + } + if (pdgCode == 2212) { + registryMCtruth.fill(HIST("plus/MCtruthPr"), part.pt(), part.eta()); + registryMCtruth.fill(HIST("plus/MCtruthPrPt"), part.pt()); + } + + if (pdgParticle->Charge() < 0.0) { + registryMCtruth.fill(HIST("minus/MCtruthAllPt"), part.pt()); + } + if (pdgCode == -211) { + registryMCtruth.fill(HIST("minus/MCtruthPi"), part.pt(), part.eta()); + registryMCtruth.fill(HIST("minus/MCtruthPiPt"), part.pt()); + } + if (pdgCode == -2212) { + registryMCtruth.fill(HIST("minus/MCtruthPr"), part.pt(), part.eta()); + registryMCtruth.fill(HIST("minus/MCtruthPrPt"), part.pt()); + } + } + } + + PROCESS_SWITCH(FemtoUniversePairTaskTrackV0Helicity, processMCTruth, "Process MC truth data", false); + + void processMCReco(FemtoRecoParticles const& parts, aod::FdMCParticles const& mcparts) + { + for (const auto& part : parts) { + auto mcPartId = part.fdMCParticleId(); + if (mcPartId == -1) + continue; // no MC particle + const auto& mcpart = mcparts.iteratorAt(mcPartId); + // + if (part.partType() == aod::femtouniverseparticle::ParticleType::kV0) { + if (mcpart.pdgMCTruth() == 3122) { + const auto& posChild = parts.iteratorAt(part.globalIndex() - 2); + const auto& negChild = parts.iteratorAt(part.globalIndex() - 1); + /// Daughters that do not pass this condition are not selected + if (isParticleTPC(posChild, 0) && isParticleTPC(negChild, 1)) { + registryMCreco.fill(HIST("plus/MCrecoLambda"), mcpart.pt(), mcpart.eta()); // lambda + if (auto mcpartIdChild = posChild.fdMCParticleId(); mcpartIdChild != -1) { + const auto& mcpartChild = mcparts.iteratorAt(mcpartIdChild); + registryMCreco.fill(HIST("plus/MCrecoLambdaChildPr"), mcpartChild.pt(), mcpartChild.eta()); // lambda proton child + } + if (auto mcpartIdChild = negChild.fdMCParticleId(); mcpartIdChild != -1) { + const auto& mcpartChild = mcparts.iteratorAt(mcpartIdChild); + registryMCreco.fill(HIST("plus/MCrecoLambdaChildPi"), mcpartChild.pt(), mcpartChild.eta()); // lambda pion child + } + } + } else if (mcpart.pdgMCTruth() == -3122) { + const auto& posChild = parts.iteratorAt(part.globalIndex() - 2); + const auto& negChild = parts.iteratorAt(part.globalIndex() - 1); + /// Daughters that do not pass this condition are not selected + if (isParticleTPC(posChild, 1) && isParticleTPC(negChild, 0)) { + registryMCreco.fill(HIST("minus/MCrecoLambda"), mcpart.pt(), mcpart.eta()); // anti-lambda + if (auto mcpartIdChild = posChild.fdMCParticleId(); mcpartIdChild != -1) { + const auto& mcpartChild = mcparts.iteratorAt(mcpartIdChild); + registryMCreco.fill(HIST("minus/MCrecoLambdaChildPi"), mcpartChild.pt(), mcpartChild.eta()); // anti-lambda pion child + } + if (auto mcpartIdChild = negChild.fdMCParticleId(); mcpartIdChild != -1) { + const auto& mcpartChild = mcparts.iteratorAt(mcpartIdChild); + registryMCreco.fill(HIST("minus/MCrecoLambdaChildPr"), mcpartChild.pt(), mcpartChild.eta()); // anti-lambda proton child + } + } + } + } else if (part.partType() == aod::femtouniverseparticle::ParticleType::kTrack) { + if (part.sign() > 0) { + registryMCreco.fill(HIST("plus/MCrecoAllPt"), mcpart.pt()); + if (mcpart.pdgMCTruth() == 211 && isNSigmaCombined(part.p(), unPackInTable(part.tpcNSigmaStorePi()), unPackInTable(part.tofNSigmaStorePi()))) { + registryMCreco.fill(HIST("plus/MCrecoPi"), mcpart.pt(), mcpart.eta()); + registryMCreco.fill(HIST("plus/MCrecoPiPt"), mcpart.pt()); + } else if (mcpart.pdgMCTruth() == 2212 && isNSigmaCombined(part.p(), unPackInTable(part.tpcNSigmaStorePr()), unPackInTable(part.tofNSigmaStorePr()))) { + registryMCreco.fill(HIST("plus/MCrecoPr"), mcpart.pt(), mcpart.eta()); + registryMCreco.fill(HIST("plus/MCrecoPrPt"), mcpart.pt()); + } + } + + if (part.sign() < 0) { + registryMCreco.fill(HIST("minus/MCrecoAllPt"), mcpart.pt()); + if (mcpart.pdgMCTruth() == -211 && isNSigmaCombined(part.p(), unPackInTable(part.tpcNSigmaStorePi()), unPackInTable(part.tofNSigmaStorePi()))) { + registryMCreco.fill(HIST("minus/MCrecoPi"), mcpart.pt(), mcpart.eta()); + registryMCreco.fill(HIST("minus/MCrecoPiPt"), mcpart.pt()); + } else if (mcpart.pdgMCTruth() == -2212 && isNSigmaCombined(part.p(), unPackInTable(part.tpcNSigmaStorePr()), unPackInTable(part.tofNSigmaStorePr()))) { + registryMCreco.fill(HIST("minus/MCrecoPr"), mcpart.pt(), mcpart.eta()); + registryMCreco.fill(HIST("minus/MCrecoPrPt"), mcpart.pt()); + } + } + } // partType + } + } + + PROCESS_SWITCH(FemtoUniversePairTaskTrackV0Helicity, processMCReco, "Process MC reco data", false); +}; + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + WorkflowSpec workflow{ + adaptAnalysisTask(cfgc), + }; + return workflow; +} From 436e61b280e5c12cea2363b99c44fddb33d66b28 Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Sat, 2 Aug 2025 15:46:47 +0200 Subject: [PATCH 15/19] Add helicity task --- PWGCF/FemtoUniverse/Tasks/CMakeLists.txt | 2 +- .../femtoUniversePairTaskTrackV0Helicity.cxx | 17 +++++++++-------- 2 files changed, 10 insertions(+), 9 deletions(-) diff --git a/PWGCF/FemtoUniverse/Tasks/CMakeLists.txt b/PWGCF/FemtoUniverse/Tasks/CMakeLists.txt index 387a2ed1958..a568e8b1252 100644 --- a/PWGCF/FemtoUniverse/Tasks/CMakeLists.txt +++ b/PWGCF/FemtoUniverse/Tasks/CMakeLists.txt @@ -33,7 +33,7 @@ o2physics_add_dpl_workflow(femtouniverse-pair-track-nucleus SOURCES femtoUniversePairTaskTrackNucleus.cxx PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore COMPONENT_NAME Analysis) - + o2physics_add_dpl_workflow(femtouniverse-pair-track-v0-helicity SOURCES femtoUniversePairTaskTrackV0Helicity.cxx PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx index 6e7052e5f3a..d2119528229 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx @@ -97,6 +97,7 @@ struct FemtoUniversePairTaskTrackV0Helicity { ConfigurableAxis confChildTempFitVarpTBins{"confChildTempFitVarpTBins", {20, 0.5, 4.05}, "V0 child: pT binning of the pT vs. TempFitVar plot"}; Configurable confHPtPart2{"confHPtPart2", 4.0f, "higher limit for pt of particle 2"}; Configurable confLPtPart2{"confLPtPart2", 0.3f, "lower limit for pt of particle 2"}; + Configurable confPDGCodeV0{"confPDGCodeV0", 3122, "V0 -- PDG code"}; Configurable confPDGCodePosChild{"confPDGCodePosChild", 2212, "Positive Child -- PDG code"}; Configurable confPDGCodeNegChild{"confPDGCodeNegChild", 211, "Negative Child -- PDG code"}; @@ -506,7 +507,7 @@ struct FemtoUniversePairTaskTrackV0Helicity { /// Histogramming same event for (const auto& part : groupPartsTwo) { int pdgCode = static_cast(part.pidCut()); - if ((confV0Type1 == 0 && pdgCode != 3122) || (confV0Type1 == 1 && pdgCode != -3122)) + if ((confV0Type1 == 0 && pdgCode != confPDGCodeV0) || (confV0Type1 == 1 && pdgCode != -confPDGCodeV0)) continue; trackHistoPartTwo.fillQA(part); } @@ -532,7 +533,7 @@ struct FemtoUniversePairTaskTrackV0Helicity { if (static_cast(p1.pidCut()) != confTrkPDGCodePartOne) continue; int pdgCode2 = static_cast(p2.pidCut()); - if ((confV0Type1 == 0 && pdgCode2 != 3122) || (confV0Type1 == 1 && pdgCode2 != -3122)) + if ((confV0Type1 == 0 && pdgCode2 != confPDGCodeV0) || (confV0Type1 == 1 && pdgCode2 != -confPDGCodeV0)) continue; // track cleaning if (confIsCPR.value) { @@ -557,17 +558,17 @@ struct FemtoUniversePairTaskTrackV0Helicity { /// Histogramming same event for (const auto& part : groupPartsTwo) { int pdgCode = static_cast(part.pidCut()); - if ((confV0Type1 == 0 && pdgCode != 3122) || (confV0Type1 == 1 && pdgCode != -3122)) + if ((confV0Type1 == 0 && pdgCode != confPDGCodeV0) || (confV0Type1 == 1 && pdgCode != -confPDGCodeV0)) continue; trackHistoPartTwo.fillQA(part); } auto pairProcessFunc = [&](auto& p1, auto& p2) -> void { int pdgCode1 = static_cast(p1.pidCut()); - if ((confV0Type1 == 0 && pdgCode1 != 3122) || (confV0Type1 == 1 && pdgCode1 != -3122)) + if ((confV0Type1 == 0 && pdgCode1 != confPDGCodeV0) || (confV0Type1 == 1 && pdgCode1 != -confPDGCodeV0)) return; int pdgCode2 = static_cast(p2.pidCut()); - if ((confV0Type2 == 0 && pdgCode2 != 3122) || (confV0Type2 == 1 && pdgCode2 != -3122)) + if ((confV0Type2 == 0 && pdgCode2 != confPDGCodeV0) || (confV0Type2 == 1 && pdgCode2 != -confPDGCodeV0)) return; sameEventCont.setPair(p1, p2, multCol, confUse3D); }; @@ -756,7 +757,7 @@ struct FemtoUniversePairTaskTrackV0Helicity { if (static_cast(p1.pidCut()) != confTrkPDGCodePartOne) continue; int pdgCode2 = static_cast(p2.pidCut()); - if ((confV0Type1 == 0 && pdgCode2 != 3122) || (confV0Type1 == 1 && pdgCode2 != -3122)) + if ((confV0Type1 == 0 && pdgCode2 != confPDGCodeV0) || (confV0Type1 == 1 && pdgCode2 != -confPDGCodeV0)) continue; if (confIsCPR.value) { if (pairCloseRejection.isClosePair(p1, p2, parts, magFieldTesla1, femto_universe_container::EventType::mixed)) { @@ -796,10 +797,10 @@ struct FemtoUniversePairTaskTrackV0Helicity { for (const auto& [p1, p2] : combinations(CombinationsFullIndexPolicy(groupPartsOne, groupPartsTwo))) { int pdgCode1 = static_cast(p1.pidCut()); - if ((confV0Type1 == 0 && pdgCode1 != 3122) || (confV0Type1 == 1 && pdgCode1 != -3122)) + if ((confV0Type1 == 0 && pdgCode1 != confPDGCodeV0) || (confV0Type1 == 1 && pdgCode1 != -confPDGCodeV0)) continue; int pdgCode2 = static_cast(p2.pidCut()); - if ((confV0Type2 == 0 && pdgCode2 != 3122) || (confV0Type2 == 1 && pdgCode2 != -3122)) + if ((confV0Type2 == 0 && pdgCode2 != confPDGCodeV0) || (confV0Type2 == 1 && pdgCode2 != -confPDGCodeV0)) continue; mixedEventCont.setPair(p1, p2, multCol, confUse3D); } From 5709ca399ebc429d55d46f53e89275e90f08b99d Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Sun, 3 Aug 2025 14:43:15 +0200 Subject: [PATCH 16/19] Add a new helicity task --- .../Core/FemtoUniverseContainer.h | 1 + .../femtoUniversePairTaskTrackV0Helicity.cxx | 29 ++++++++++--------- 2 files changed, 17 insertions(+), 13 deletions(-) diff --git a/PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h b/PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h index ceacc7fbf3e..e034ec66ba7 100644 --- a/PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h +++ b/PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h @@ -26,6 +26,7 @@ #include "Framework/HistogramRegistry.h" #include "Common/Core/RecoDecay.h" +#include "PWGCF/FemtoUniverse/DataModel/FemtoDerived.h" #include "PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h" #include "Math/Vector4D.h" diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx index d2119528229..46def7bba0e 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx @@ -15,27 +15,30 @@ /// \author Zuzanna Chochulska, WUT Warsaw & CTU Prague, zchochul@cern.ch /// \author Shirajum Monira, WUT Warsaw, shirajum.monira.dokt@pw.edu.pl -#include -#include -#include -#include "Framework/AnalysisTask.h" -#include "Framework/runDataProcessing.h" -#include "Framework/HistogramRegistry.h" -#include "Framework/ASoAHelpers.h" -#include "Framework/RunningWorkflowInfo.h" -#include "Framework/StepTHn.h" -#include "PWGCF/FemtoUniverse/DataModel/FemtoDerived.h" -#include "PWGCF/FemtoUniverse/Core/FemtoUniverseParticleHisto.h" +#include "PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h" +#include "PWGCF/FemtoUniverse/Core/FemtoUniverseDetaDphiStar.h" #include "PWGCF/FemtoUniverse/Core/FemtoUniverseEventHisto.h" #include "PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h" #include "PWGCF/FemtoUniverse/Core/FemtoUniversePairCleaner.h" -#include "PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h" -#include "PWGCF/FemtoUniverse/Core/FemtoUniverseDetaDphiStar.h" +#include "PWGCF/FemtoUniverse/Core/FemtoUniverseParticleHisto.h" #include "PWGCF/FemtoUniverse/Core/femtoUtils.h" +#include "PWGCF/FemtoUniverse/DataModel/FemtoDerived.h" + +#include "Framework/ASoAHelpers.h" +#include "Framework/AnalysisTask.h" +#include "Framework/HistogramRegistry.h" #include "Framework/O2DatabasePDGPlugin.h" +#include "Framework/RunningWorkflowInfo.h" +#include "Framework/StepTHn.h" +#include "Framework/runDataProcessing.h" + #include #include +#include +#include +#include + using namespace o2; using namespace o2::soa; using namespace o2::framework; From e3d01ca27b9754eac61d32888793367c8e7fa2ea Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Sun, 3 Aug 2025 14:52:25 +0200 Subject: [PATCH 17/19] Add a new helicity task --- .../FemtoUniverse/Core/FemtoUniverseContainer.h | 17 ++++++++++------- 1 file changed, 10 insertions(+), 7 deletions(-) diff --git a/PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h b/PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h index e034ec66ba7..71cd731a4f2 100644 --- a/PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h +++ b/PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h @@ -20,18 +20,21 @@ #ifndef PWGCF_FEMTOUNIVERSE_CORE_FEMTOUNIVERSECONTAINER_H_ #define PWGCF_FEMTOUNIVERSE_CORE_FEMTOUNIVERSECONTAINER_H_ -#include -#include -#include +#include "PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h" +#include "PWGCF/FemtoUniverse/DataModel/FemtoDerived.h" -#include "Framework/HistogramRegistry.h" #include "Common/Core/RecoDecay.h" -#include "PWGCF/FemtoUniverse/DataModel/FemtoDerived.h" -#include "PWGCF/FemtoUniverse/Core/FemtoUniverseMath.h" + +#include "Framework/HistogramRegistry.h" #include "Math/Vector4D.h" -#include "TMath.h" #include "TDatabasePDG.h" +#include "TMath.h" + +#include + +#include +#include using namespace o2::framework; From ad359ba4495bc779092ad5a03ce9cddb39adefc3 Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Mon, 4 Aug 2025 16:39:04 +0200 Subject: [PATCH 18/19] Add a new helicity task --- PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h b/PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h index 71cd731a4f2..9bd60ccc34f 100644 --- a/PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h +++ b/PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h @@ -57,7 +57,7 @@ enum EventType { same, ///< Pair from same event /// \brief Container for all histogramming related to the correlation function. The two /// particles of the pair are passed here, and the correlation function and QA histograms /// are filled according to the specified observable -/// \tparam eventType Type of the event (same/mixed) +/// \tparam eventType Type of the event (same or mixed) /// \tparam obs Observable to be computed (k*/Q_inv/...) template class FemtoUniverseContainer From a117e4258d9d22d27299e3678539407f9c05f414 Mon Sep 17 00:00:00 2001 From: Anna-Mariia Andrushko <01175703@pw.edu.pl> Date: Tue, 5 Aug 2025 14:03:04 +0200 Subject: [PATCH 19/19] Add MC Truth to helicity angle analysis --- .../femtoUniversePairTaskTrackV0Helicity.cxx | 46 ++++++++++++++----- 1 file changed, 35 insertions(+), 11 deletions(-) diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx index 46def7bba0e..b5fe5fbdc5a 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackV0Helicity.cxx @@ -161,13 +161,12 @@ struct FemtoUniversePairTaskTrackV0Helicity { /// Histogram output HistogramRegistry qaRegistry{"TrackQA", {}, OutputObjHandlingPolicy::AnalysisObject}; + HistogramRegistry thetaRegistry{"ThetaQA", {}, OutputObjHandlingPolicy::AnalysisObject}; HistogramRegistry resultRegistry{"Correlations", {}, OutputObjHandlingPolicy::AnalysisObject}; HistogramRegistry registryMCtruth{"MCtruthHistos", {}, OutputObjHandlingPolicy::AnalysisObject, false, true}; HistogramRegistry registryMCreco{"MCrecoHistos", {}, OutputObjHandlingPolicy::AnalysisObject, false, true}; HistogramRegistry mixQaRegistry{"mixQaRegistry", {}, OutputObjHandlingPolicy::AnalysisObject}; - HistogramRegistry thetaRegistry{"ThetaQA", {}, OutputObjHandlingPolicy::AnalysisObject}; - std::unique_ptr plocalEffFile; std::unique_ptr plocalEffp1; std::unique_ptr plocalEffp2; @@ -236,6 +235,15 @@ struct FemtoUniversePairTaskTrackV0Helicity { posChildV0Type2.init(&qaRegistry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, 0, true, "posChildV0Type2"); negChildV0Type2.init(&qaRegistry, confChildTempFitVarpTBins, confChildTempFitVarBins, false, 0, true, "negChildV0Type2"); + // Helicity angle + thetaRegistry.add("Theta/hTheta", " ; p (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/PositiveChild/hThetaPt", " ; p_{T} (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/PositiveChild/hThetaEta", " ; #eta; cos(#theta)", kTH2F, {{100, -1, 1}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/PositiveChild/hThetaPhi", " ; #phi; cos(#theta)", kTH2F, {{100, -1, 7}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/NegativeChild/hThetaPt", " ; p_{T} (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/NegativeChild/hThetaEta", " ; #eta; cos(#theta)", kTH2F, {{100, -1, 1}, {110, -1.1, 1.1}}); + thetaRegistry.add("Theta/NegativeChild/hThetaPhi", " ; #phi; cos(#theta)", kTH2F, {{100, -1, 7}, {110, -1.1, 1.1}}); + mixQaRegistry.add("MixingQA/hMECollisionBins", ";bin;Entries", kTH1F, {{120, -0.5, 119.5}}); // MC truth @@ -256,6 +264,14 @@ struct FemtoUniversePairTaskTrackV0Helicity { registryMCtruth.add("minus/MCtruthPiPt", "MC truth pions;#it{p}_{T} (GeV/c)", {HistType::kTH1F, {{500, 0, 5}}}); registryMCtruth.add("minus/MCtruthPrPt", "MC truth protons;#it{p}_{T} (GeV/c)", {HistType::kTH1F, {{500, 0, 5}}}); + registryMCtruth.add("ThetaMCTruth/hTheta", " ; p (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); + registryMCtruth.add("ThetaMCTruth/PositiveChild/hThetaPt", " ; p_{T} (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); + registryMCtruth.add("ThetaMCTruth/PositiveChild/hThetaEta", " ; #eta; cos(#theta)", kTH2F, {{100, -1, 1}, {110, -1.1, 1.1}}); + registryMCtruth.add("ThetaMCTruth/PositiveChild/hThetaPhi", " ; #phi; cos(#theta)", kTH2F, {{100, -1, 7}, {110, -1.1, 1.1}}); + registryMCtruth.add("ThetaMCTruth/NegativeChild/hThetaPt", " ; p_{T} (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); + registryMCtruth.add("ThetaMCTruth/NegativeChild/hThetaEta", " ; #eta; cos(#theta)", kTH2F, {{100, -1, 1}, {110, -1.1, 1.1}}); + registryMCtruth.add("ThetaMCTruth/NegativeChild/hThetaPhi", " ; #phi; cos(#theta)", kTH2F, {{100, -1, 7}, {110, -1.1, 1.1}}); + // MC reco registryMCreco.add("plus/MCrecoLambda", "MC reco Lambdas;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); registryMCreco.add("plus/MCrecoLambdaChildPr", "MC reco Lambdas;#it{p}_{T} (GeV/c); #eta", {HistType::kTH2F, {{500, 0, 5}, {400, -1.0, 1.0}}}); @@ -278,15 +294,6 @@ struct FemtoUniversePairTaskTrackV0Helicity { registryMCreco.add("minus/MCrecoPiPt", "MC reco pions;#it{p}_{T} (GeV/c)", {HistType::kTH1F, {{500, 0, 5}}}); registryMCreco.add("minus/MCrecoPrPt", "MC reco protons;#it{p}_{T} (GeV/c)", {HistType::kTH1F, {{500, 0, 5}}}); - // Helicity angle - thetaRegistry.add("Theta/hTheta", " ; p (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); - thetaRegistry.add("Theta/PositiveChild/hThetaPt", " ; p_{T} (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); - thetaRegistry.add("Theta/PositiveChild/hThetaEta", " ; #eta; cos(#theta)", kTH2F, {{100, -1, 1}, {110, -1.1, 1.1}}); - thetaRegistry.add("Theta/PositiveChild/hThetaPhi", " ; #phi; cos(#theta)", kTH2F, {{100, -1, 7}, {110, -1.1, 1.1}}); - thetaRegistry.add("Theta/NegativeChild/hThetaPt", " ; p_{T} (GeV/#it{c}); cos(#theta)", kTH2F, {{100, 0, 10}, {110, -1.1, 1.1}}); - thetaRegistry.add("Theta/NegativeChild/hThetaEta", " ; #eta; cos(#theta)", kTH2F, {{100, -1, 1}, {110, -1.1, 1.1}}); - thetaRegistry.add("Theta/NegativeChild/hThetaPhi", " ; #phi; cos(#theta)", kTH2F, {{100, -1, 7}, {110, -1.1, 1.1}}); - sameEventCont.init(&resultRegistry, confkstarBins, confMultBins, confkTBins, confmTBins, confMultBins3D, confmTBins3D, confEtaBins, confPhiBins, confIsMC, confUse3D); sameEventCont.setPDGCodes(confTrkPDGCodePartOne, confV0PDGCodePartTwo); mixedEventCont.init(&resultRegistry, confkstarBins, confMultBins, confkTBins, confmTBins, confMultBins3D, confmTBins3D, confEtaBins, confPhiBins, confIsMC, confUse3D); @@ -869,6 +876,23 @@ struct FemtoUniversePairTaskTrackV0Helicity { registryMCtruth.fill(HIST("minus/MCtruthPr"), part.pt(), part.eta()); registryMCtruth.fill(HIST("minus/MCtruthPrPt"), part.pt()); } + + // Helicity angle + const auto& posChild = parts.iteratorAt(part.index() - 2); + const auto& negChild = parts.iteratorAt(part.index() - 1); + + auto posChildMass = pdg->Mass(confPDGCodePosChild); + auto negChildMass = pdg->Mass(confPDGCodeNegChild); + auto posChildBoosted = FemtoUniverseMath::boostPRF(posChild, posChildMass, negChild, negChildMass); + auto cosineTheta = (posChildBoosted.Px() * part.px() + posChildBoosted.Py() * part.py() + posChildBoosted.Pz() * part.pz()) / (posChildBoosted.P() * part.p()); + + registryMCtruth.fill(HIST("ThetaMCTruth/hTheta"), part.p(), cosineTheta); + registryMCtruth.fill(HIST("ThetaMCTruth/PositiveChild/hThetaPt"), posChild.pt(), cosineTheta); + registryMCtruth.fill(HIST("ThetaMCTruth/PositiveChild/hThetaEta"), posChild.eta(), cosineTheta); + registryMCtruth.fill(HIST("ThetaMCTruth/PositiveChild/hThetaPhi"), posChild.phi(), cosineTheta); + registryMCtruth.fill(HIST("ThetaMCTruth/NegativeChild/hThetaPt"), negChild.pt(), cosineTheta); + registryMCtruth.fill(HIST("ThetaMCTruth/NegativeChild/hThetaEta"), negChild.eta(), cosineTheta); + registryMCtruth.fill(HIST("ThetaMCTruth/NegativeChild/hThetaPhi"), negChild.phi(), cosineTheta); } }