From 7d0e582d171c9dadc3c1fa0151a6ac3622c211af Mon Sep 17 00:00:00 2001 From: Yash Patley <52608802+yashpatley@users.noreply.github.com> Date: Sat, 17 Jan 2026 22:05:15 +0530 Subject: [PATCH 1/4] [PWGCF]: Update flowEventPlane.cxx 1. Added track table to fix previous error of table size mismatch 2. Added v0 flow --- PWGCF/Flow/Tasks/flowEventPlane.cxx | 341 +++++++++++++++++++++++++--- 1 file changed, 306 insertions(+), 35 deletions(-) diff --git a/PWGCF/Flow/Tasks/flowEventPlane.cxx b/PWGCF/Flow/Tasks/flowEventPlane.cxx index 8287f3df49b..131b754dd18 100644 --- a/PWGCF/Flow/Tasks/flowEventPlane.cxx +++ b/PWGCF/Flow/Tasks/flowEventPlane.cxx @@ -13,6 +13,8 @@ /// \brief Flow calculation using event plane. /// \author Yash Patley +#include "PWGLF/DataModel/LFStrangenessTables.h" + #include "Common/Core/RecoDecay.h" #include "Common/DataModel/Centrality.h" #include "Common/DataModel/CollisionAssociationTables.h" @@ -56,16 +58,32 @@ DECLARE_SOA_TABLE(ColSPExt, "AOD", "COLSPEXT", o2::soa::Index<>, colspext::Xc, colspext::Yc); -namespace trackidext +namespace tracksid { +DECLARE_SOA_INDEX_COLUMN(Collision, collision); +DECLARE_SOA_COLUMN(Sign, sign, int8_t); +DECLARE_SOA_COLUMN(Px, px, float); +DECLARE_SOA_COLUMN(Py, py, float); +DECLARE_SOA_COLUMN(Pz, pz, float); +DECLARE_SOA_COLUMN(Pt, pt, float); +DECLARE_SOA_COLUMN(Eta, eta, float); +DECLARE_SOA_COLUMN(Phi, phi, float); DECLARE_SOA_COLUMN(IsPion, isPion, bool); DECLARE_SOA_COLUMN(IsKaon, isKaon, bool); DECLARE_SOA_COLUMN(IsProton, isProton, bool); -} // namespace trackidext -DECLARE_SOA_TABLE(TrackIdExt, "AOD", "TRACKIDEXT", o2::soa::Index<>, - trackidext::IsPion, - trackidext::IsKaon, - trackidext::IsProton); +} // namespace tracksid +DECLARE_SOA_TABLE(TracksId, "AOD", "TRACKSID", o2::soa::Index<>, + aod::tracksid::CollisionId, + tracksid::Sign, + tracksid::Px, + tracksid::Py, + tracksid::Pz, + tracksid::Pt, + tracksid::Eta, + tracksid::Phi, + tracksid::IsPion, + tracksid::IsKaon, + tracksid::IsProton); } // namespace o2::aod enum GainClibCorr { @@ -102,6 +120,11 @@ enum ParticleType { kNPart }; +enum V0Type { + kLambda = 0, + kAntiLambda +}; + struct SpectatorPlaneTableProducer { // Table producer Produces colSPExtTable; @@ -602,7 +625,7 @@ struct SpectatorPlaneTableProducer { struct IdHadronFlow { // Table producer - Produces trackIdExtTable; + Produces tracksIdTable; // Tracks Configurable cTrackMinPt{"cTrackMinPt", 0.1, "p_{T} minimum"}; @@ -614,6 +637,8 @@ struct IdHadronFlow { Configurable cTrackDcaZCut{"cTrackDcaZCut", 1., "DcaXY Cut"}; // Track PID + Configurable cTpcElRejCutMin{"cTpcElRejCutMin", -3., "Electron Rejection Cut Minimum"}; + Configurable cTpcElRejCutMax{"cTpcElRejCutMax", 5., "Electron Rejection Cut Maximum"}; Configurable cTpcNSigmaCut{"cTpcNSigmaCut", 2, "TPC NSigma Cut"}; Configurable cTpcRejCut{"cTpcRejCut", 3, "TPC Rej Cut"}; Configurable cTofNSigmaCut{"cTofNSigmaCut", 2, "TOF NSigma Cut"}; @@ -622,6 +647,9 @@ struct IdHadronFlow { Configurable cKaonPtCut{"cKaonPtCut", 0.6, "Kaon TPC pT cutoff"}; Configurable cProtonPtCut{"cProtonPtCut", 1.1, "Proton TPC pT cutoff"}; + // Flag to fill histograms + Configurable cFillIdHist{"cFillIdHist", false, "Flag to fill histograms"}; + // Histogram registry: an object to hold your histograms HistogramRegistry histos{"histos", {}, OutputObjHandlingPolicy::AnalysisObject}; @@ -659,15 +687,17 @@ struct IdHadronFlow { histos.add("DF/hCQuNeg", "u_{x}X^{C}_{1} + u_{y}Y^{C}_{1}", kTProfile2D, {axisCent, axisTrackEta}); // Identified particle - histos.add("PartId/Pion/hdEdX", "PartId/Pion/hdEdX", kTH2F, {axisTrackPt, axisTrackdEdx}); - histos.add("PartId/Pion/hTPCNSigma", "PartId/Pion/hTPCNSigma", kTH2F, {axisTrackPt, axisTrackNSigma}); - histos.add("PartId/Pion/hTOFNSigma", "PartId/Pion/hTOFNSigma", kTH2F, {axisTrackPt, axisTrackNSigma}); - histos.add("PartId/Pion/hAQuPos", "PartId/Pion/hAQuPos", kTProfile2D, {axisCent, axisTrackEta}); - histos.add("PartId/Pion/hAQuNeg", "PartId/Pion/hAQuNeg", kTProfile2D, {axisCent, axisTrackEta}); - histos.add("PartId/Pion/hCQuPos", "PartId/Pion/hCQuPos", kTProfile2D, {axisCent, axisTrackEta}); - histos.add("PartId/Pion/hCQuNeg", "PartId/Pion/hCQuNeg", kTProfile2D, {axisCent, axisTrackEta}); - histos.addClone("PartId/Pion/", "PartId/Kaon/"); - histos.addClone("PartId/Pion/", "PartId/Proton/"); + if (cFillIdHist) { + histos.add("PartId/Pion/hdEdX", "PartId/Pion/hdEdX", kTH2F, {axisTrackPt, axisTrackdEdx}); + histos.add("PartId/Pion/hTPCNSigma", "PartId/Pion/hTPCNSigma", kTH2F, {axisTrackPt, axisTrackNSigma}); + histos.add("PartId/Pion/hTOFNSigma", "PartId/Pion/hTOFNSigma", kTH2F, {axisTrackPt, axisTrackNSigma}); + histos.add("PartId/Pion/hAQuPos", "PartId/Pion/hAQuPos", kTProfile2D, {axisCent, axisTrackEta}); + histos.add("PartId/Pion/hAQuNeg", "PartId/Pion/hAQuNeg", kTProfile2D, {axisCent, axisTrackEta}); + histos.add("PartId/Pion/hCQuPos", "PartId/Pion/hCQuPos", kTProfile2D, {axisCent, axisTrackEta}); + histos.add("PartId/Pion/hCQuNeg", "PartId/Pion/hCQuNeg", kTProfile2D, {axisCent, axisTrackEta}); + histos.addClone("PartId/Pion/", "PartId/Kaon/"); + histos.addClone("PartId/Pion/", "PartId/Proton/"); + } } // Track Selection @@ -708,6 +738,12 @@ struct IdHadronFlow { template bool identifyTrack(T const& track) { + // Electron rejection + if (track.tpcNSigmaEl() > cTpcElRejCutMin && track.tpcNSigmaEl() < cTpcElRejCutMax) { + return false; + } + + // Pion, Kaon, Proton Identification std::vector vPtCut = {cPionPtCut, cKaonPtCut, cProtonPtCut}; std::vector vTpcNsig = {std::abs(track.tpcNSigmaPi()), std::abs(track.tpcNSigmaKa()), std::abs(track.tpcNSigmaPr())}; std::vector vTofNsig = {std::abs(track.tofNSigmaPi()), std::abs(track.tofNSigmaKa()), std::abs(track.tofNSigmaPr())}; @@ -765,7 +801,7 @@ struct IdHadronFlow { } using CollisionsRun3 = soa::Join; - using Tracks = soa::Join; + using Tracks = soa::Join; void processDummy(CollisionsRun3::iterator const&) {} @@ -797,7 +833,7 @@ struct IdHadronFlow { for (auto const& track : tracks) { // Select track if (!selectTrack(track)) { - trackIdExtTable(false, false, false); + tracksIdTable(collision.index(), track.sign(), track.px(), track.py(), track.pz(), track.pt(), track.eta(), track.phi(), false, false, false); continue; } @@ -821,17 +857,28 @@ struct IdHadronFlow { histos.fill(HIST("DF/hCQuNeg"), cent, track.eta(), v1c); } - // Identified directed flow - if (identifyTrack(track)) { - trackIdExtTable(true, false, false); - getIdHadronFlow(cent, track, v1a, v1c); - } else if (identifyTrack(track)) { - trackIdExtTable(false, true, false); - getIdHadronFlow(cent, track, v1a, v1c); - } else if (identifyTrack(track)) { - trackIdExtTable(false, false, true); - getIdHadronFlow(cent, track, v1a, v1c); + // Identify track + bool pionFlag = identifyTrack(track); + bool kaonFlag = identifyTrack(track); + bool protonFlag = identifyTrack(track); + + if (cFillIdHist) { + // Pion + if (pionFlag) { + getIdHadronFlow(cent, track, v1a, v1c); + } + // Kaon + if (kaonFlag) { + getIdHadronFlow(cent, track, v1a, v1c); + } + // Proton + if (protonFlag) { + getIdHadronFlow(cent, track, v1a, v1c); + } } + + // Fill track table + tracksIdTable(collision.index(), track.sign(), track.px(), track.py(), track.pz(), track.pt(), track.eta(), track.phi(), pionFlag, kaonFlag, protonFlag); } } PROCESS_SWITCH(IdHadronFlow, processIdHadronFlow, "Identified hadron flow process", false); @@ -843,6 +890,28 @@ struct FlowEventPlane { Configurable cNInvMassBins{"cNInvMassBins", 500, "# of m bins"}; Configurable cResRapCut{"cResRapCut", 0.5, "Resonance rapidity cut"}; + // V0 + // Tracks + Configurable cTrackPtCut{"cTrackPtCut", 0.1, "p_{T} minimum"}; + Configurable cTrackEtaCut{"cTrackEtaCut", 0.8, "Pseudorapidity cut"}; + Configurable cTpcNsigmaCut{"cTpcNsigmaCut", 3.0, "TPC NSigma Selection Cut"}; + + // V0s + Configurable cV0TypeSelection{"cV0TypeSelection", 1, "V0 Type Selection"}; + Configurable cMinDcaProtonToPV{"cMinDcaProtonToPV", 0.02, "Minimum Proton DCAr to PV"}; + Configurable cMinDcaPionToPV{"cMinDcaPionToPV", 0.06, "Minimum Pion DCAr to PV"}; + Configurable cDcaV0Dau{"cDcaV0Dau", 1., "DCA between V0 daughters"}; + Configurable cDcaV0ToPv{"cDcaV0ToPv", 0.1, "DCA V0 to PV"}; + Configurable cMinV0Radius{"cMinV0Radius", 0.5, "Minimum V0 radius from PV"}; + Configurable cMaxV0Radius{"cMaxV0Radius", 999.0, "Maximum V0 radius from PV"}; + Configurable cV0CTau{"cV0CTau", 30.0, "Decay length cut"}; + Configurable cV0CosPA{"cV0CosPA", 0.995, "V0 CosPA to PV"}; + Configurable cK0SMassRej{"cK0SMassRej", 0.01, "Reject K0Short Candidates"}; + Configurable cLambdaMassWin{"cLambdaMassWin", 0.007, "Lambda Mass Window"}; + Configurable cMinV0Pt{"cMinV0Pt", 0.5, "Min v0 pT"}; + Configurable cMaxV0Pt{"cMaxV0Pt", 4.5, "Max v0 pT"}; + Configurable cV0RapCut{"cV0RapCut", 0.5, "V0 rap cut"}; + // Histogram registry: an object to hold your histograms HistogramRegistry histos{"histos", {}, OutputObjHandlingPolicy::AnalysisObject}; @@ -853,13 +922,24 @@ struct FlowEventPlane { { // Define axes const AxisSpec axisCent{100, 0., 100, "FT0C%"}; - - const AxisSpec axisXYac{600, -6, 6, "Q^{t}Q^{p}"}; - const AxisSpec axisV1{400, -4, 4, "v_{1}"}; - const AxisSpec axisTrackPt{100, 0., 10., "p_{T} (GeV/#it{c})"}; const AxisSpec axisTrackRap{cNRapBins, -0.5, 0.5, "y"}; const AxisSpec axisInvMass{cNInvMassBins, 0.87, 1.12, "M_{KK} (GeV/#it{c}^{2}"}; + const AxisSpec axisMomPID(80, 0, 4, "p (GeV/#it{c})"); + const AxisSpec axisNsigma(401, -10.025, 10.025, {"n#sigma"}); + const AxisSpec axisdEdx(360, 20, 200, "#frac{dE}{dx}"); + const AxisSpec axisRadius(2000, 0, 200, "r(cm)"); + const AxisSpec axisCosPA(300, 0.97, 1.0, "cos(#theta_{PA})"); + const AxisSpec axisDcaV0PV(1000, 0., 10., "dca (cm)"); + const AxisSpec axisDcaProngPV(5000, -50., 50., "dca (cm)"); + const AxisSpec axisDcaDau(75, 0., 1.5, "Daug DCA (#sigma)"); + const AxisSpec axisCTau(2000, 0, 200, "c#tau (cm)"); + const AxisSpec axisAlpha(40, -1, 1, "#alpha"); + const AxisSpec axisQtarm(40, 0, 0.4, "q_{T}"); + const AxisSpec axisLambdaPt(50, 0, 10, "p_{T} (GeV/#it{c})"); + const AxisSpec axisLambdaInvMass(140, 1.08, 1.15, "M_{p#pi} (GeV/#it{c}^{2})"); + const AxisSpec axisXYac{600, -6, 6, "Q^{t}Q^{p}"}; + const AxisSpec axisV1{400, -4, 4, "v_{1}"}; // Create histograms // Resonance @@ -869,6 +949,112 @@ struct FlowEventPlane { histos.add("Reso/Phi/Sig/hPhiQuC", "hPhiQuC", kTProfile3D, {axisCent, axisTrackRap, axisInvMass}); histos.add("Reso/Phi/Bkg/hPhiQuA", "hPhiQuA", kTProfile3D, {axisCent, axisTrackRap, axisInvMass}); histos.add("Reso/Phi/Bkg/hPhiQuC", "hPhiQuC", kTProfile3D, {axisCent, axisTrackRap, axisInvMass}); + + // Lambda histograms + // QA Lambda + histos.add("Lambda/QA/hQtVsAlpha", "Armentros-Podolanski Plot", kTH2F, {axisAlpha, axisQtarm}); + histos.add("Lambda/QA/hDcaV0Dau", "DCA between V0 daughters", kTH1F, {axisDcaDau}); + histos.add("Lambda/QA/hDcaPosToPv", "DCA positive prong to PV", kTH1F, {axisDcaProngPV}); + histos.add("Lambda/QA/hDcaNegToPv", "DCA negative prong to PV", kTH1F, {axisDcaProngPV}); + histos.add("Lambda/QA/hDcaV0ToPv", "DCA V0 to PV", kTH1F, {axisDcaV0PV}); + histos.add("Lambda/QA/hCosPa", "cos(#theta_{PA})", kTH1F, {axisCosPA}); + histos.add("Lambda/QA/hRxy", "V_{0} Decay Radius in XY plane", kTH1F, {axisRadius}); + histos.add("Lambda/QA/hCTau", "V_{0} c#tau", kTH1F, {axisCTau}); + histos.add("Lambda/QA/hPosdEdXVsP", "TPC Signal Pos-Prong", kTH2F, {axisMomPID, axisdEdx}); + histos.add("Lambda/QA/hNegdEdXVsP", "TPC Signal Neg-Prong", kTH2F, {axisMomPID, axisdEdx}); + histos.add("Lambda/QA/hPosNsigPrVsP", "TPC n#sigma Pos Prong", kTH2F, {axisMomPID, axisNsigma}); + histos.add("Lambda/QA/hNegNsigPrVsP", "TPC n#sigma Neg Prong", kTH2F, {axisMomPID, axisNsigma}); + histos.add("Lambda/QA/hPosNsigPiVsP", "TPC n#sigma Pos Prong", kTH2F, {axisMomPID, axisNsigma}); + histos.add("Lambda/QA/hNegNsigPiVsP", "TPC n#sigma Neg Prong", kTH2F, {axisMomPID, axisNsigma}); + + histos.add("Lambda/hInvMassVsPt", "hInvMassVsPt", kTH3F, {axisCent, axisLambdaInvMass, axisLambdaPt}); + histos.add("Lambda/Flow/hQuA", "hPhiQuA", kTProfile3D, {axisCent, axisTrackRap, axisLambdaInvMass}); + histos.add("Lambda/Flow/hQuC", "hPhiQuC", kTProfile3D, {axisCent, axisTrackRap, axisLambdaInvMass}); + + histos.addClone("Lambda/", "AntiLambda/"); + } + + template + bool selV0DauTracks(V const& v0, T const& postrack, T const& negtrack) + { + // Kinematic selection + if (postrack.pt() <= cTrackPtCut || negtrack.pt() <= cTrackPtCut || std::abs(postrack.eta()) >= cTrackEtaCut || std::abs(negtrack.eta()) >= cTrackEtaCut) { + return false; + } + + // Apply DCA Selection on Daughter Tracks Based on Lambda/AntiLambda daughters + float dcaProton = 0., dcaPion = 0.; + if (part == kLambda) { + dcaProton = std::abs(v0.dcapostopv()); + dcaPion = std::abs(v0.dcanegtopv()); + } else if (part == kAntiLambda) { + dcaPion = std::abs(v0.dcapostopv()); + dcaProton = std::abs(v0.dcanegtopv()); + } else { + return false; + } + + if (dcaProton <= cMinDcaProtonToPV || dcaPion <= cMinDcaPionToPV) { + return false; + } + + // Daughter track PID + float tpcNSigmaPr = 0., tpcNSigmaPi = 0.; + + switch (part) { + // postrack = Proton, negtrack = Pion + case kLambda: + tpcNSigmaPr = postrack.tpcNSigmaPr(); + tpcNSigmaPi = negtrack.tpcNSigmaPi(); + break; + + // negtrack = Proton, postrack = Pion + case kAntiLambda: + tpcNSigmaPr = negtrack.tpcNSigmaPr(); + tpcNSigmaPi = postrack.tpcNSigmaPi(); + break; + } + + if (std::abs(tpcNSigmaPr) >= cTpcNsigmaCut || std::abs(tpcNSigmaPi) >= cTpcNsigmaCut) { + return false; + } + + return true; + } + + template + void fillV0QAHist(C const& col, V const& v0, T const&) + { + static constexpr std::string_view SubDir[] = {"Lambda/QA/", "AntiLambda/QA/"}; + + // daugthers + auto postrack = v0.template posTrack_as(); + auto negtrack = v0.template negTrack_as(); + float mass = 0.; + + if constexpr (part == kLambda) { + mass = v0.mLambda(); + } else { + mass = v0.mAntiLambda(); + } + + // ctau + float ctau = v0.distovertotmom(col.posX(), col.posY(), col.posZ()) * MassLambda0; + + histos.fill(HIST(SubDir[part]) + HIST("hQtVsAlpha"), v0.alpha(), v0.qtarm()); + histos.fill(HIST(SubDir[part]) + HIST("hDcaV0Dau"), v0.dcaV0daughters()); + histos.fill(HIST(SubDir[part]) + HIST("hDcaPosToPv"), v0.dcapostopv()); + histos.fill(HIST(SubDir[part]) + HIST("hDcaNegToPv"), v0.dcanegtopv()); + histos.fill(HIST(SubDir[part]) + HIST("hDcaV0ToPv"), v0.dcav0topv()); + histos.fill(HIST(SubDir[part]) + HIST("hCosPa"), v0.v0cosPA()); + histos.fill(HIST(SubDir[part]) + HIST("hRxy"), v0.v0radius()); + histos.fill(HIST(SubDir[part]) + HIST("hCTau"), ctau); + histos.fill(HIST(SubDir[part]) + HIST("hPosdEdXVsP"), postrack.tpcInnerParam(), postrack.tpcSignal()); + histos.fill(HIST(SubDir[part]) + HIST("hNegdEdXVsP"), negtrack.tpcInnerParam(), negtrack.tpcSignal()); + histos.fill(HIST(SubDir[part]) + HIST("hPosNsigPrVsP"), postrack.tpcInnerParam(), postrack.tpcNSigmaPr()); + histos.fill(HIST(SubDir[part]) + HIST("hNegNsigPrVsP"), negtrack.tpcInnerParam(), negtrack.tpcNSigmaPr()); + histos.fill(HIST(SubDir[part]) + HIST("hPosNsigPiVsP"), postrack.tpcInnerParam(), postrack.tpcNSigmaPi()); + histos.fill(HIST(SubDir[part]) + HIST("hNegNsigPiVsP"), negtrack.tpcInnerParam(), negtrack.tpcNSigmaPi()); } template @@ -927,10 +1113,12 @@ struct FlowEventPlane { } using CollisionsRun3 = soa::Join; - using Tracks = soa::Join; + using Tracks = aod::TracksId; + using TracksV0s = soa::Join; + // Partitions SliceCache cache; - Partition kaonTrackPartition = (aod::trackidext::isKaon == true); + Partition kaonTrackPartition = (aod::tracksid::isKaon == true); void processDummy(CollisionsRun3::iterator const&) {} @@ -954,12 +1142,95 @@ struct FlowEventPlane { vSP[kYc] = collision.yc(); // Track partitions - auto kaonTracks = kaonTrackPartition->sliceByCached(aod::track::collisionId, collision.globalIndex(), cache); + auto kaonTracks = kaonTrackPartition->sliceByCached(aod::tracksid::collisionId, collision.globalIndex(), cache); // Resonance flow getResoFlow(kaonTracks, kaonTracks, vSP); } PROCESS_SWITCH(FlowEventPlane, processResoFlow, "Resonance flow process", false); + + void processLambdaFlow(CollisionsRun3::iterator const& collision, aod::V0Datas const& V0s, TracksV0s const& tracks) + { + // Check collision + if (!collision.selColFlag()) { + return; + } + + // Set centrality + cent = collision.centFT0C(); + + // Flow vectors + std::array vSP = {0., 0., 0., 0.}; + vSP[kXa] = collision.xa(); + vSP[kYa] = collision.ya(); + vSP[kXc] = collision.xc(); + vSP[kYc] = collision.yc(); + + // Loop over v0s + for (auto const& v0 : V0s) { + // V0 kinematic selection + if (v0.pt() <= cMinV0Pt || v0.pt() >= cMaxV0Pt || std::abs(v0.yLambda()) >= cV0RapCut) { + continue; + } + + // Topological selections + float ctau = v0.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * MassLambda0; + if (v0.dcaV0daughters() >= cDcaV0Dau || v0.dcav0topv() >= cDcaV0ToPv || + v0.v0radius() <= cMinV0Radius || v0.v0radius() >= cMaxV0Radius || + v0.v0cosPA() <= cV0CosPA || ctau >= cV0CTau || v0.v0Type() != cV0TypeSelection) { + continue; + } + + // Ks Mass Rejection + if (std::abs(v0.mK0Short() - MassK0Short) <= cK0SMassRej) { + continue; + } + + // Initialize daughter tracks + auto postrack = v0.template posTrack_as(); + auto negtrack = v0.template negTrack_as(); + + // Initialize selection flags + bool lambdaFlag = false, antiLambdaFlag = false; + + // Get v0 track as lambda + if ((std::abs(v0.mLambda() - MassLambda0) < cLambdaMassWin) && (selV0DauTracks(v0, postrack, negtrack))) { + lambdaFlag = true; + } + + // Get v0 track as anti-lambda + if ((std::abs(v0.mAntiLambda() - MassLambda0) < cLambdaMassWin) && (selV0DauTracks(v0, postrack, negtrack))) { + antiLambdaFlag = true; + } + + // Lambda/Anti-Lambda selection checks + if (!lambdaFlag && !antiLambdaFlag) { // neither Lambda nor Anti-Lambda + continue; + } else if (lambdaFlag && antiLambdaFlag) { // check if the track is identified as lambda and anti-lambda both (DISCARD THIS TRACK) + continue; + } + + // We have a Lambda/Anti-Lambda + // Directed flow + float ux = std::cos(v0.phi()); + float uy = std::sin(v0.phi()); + float v1a = ux * vSP[kXa] + uy * vSP[kYa]; + float v1c = ux * vSP[kXc] + uy * vSP[kYc]; + + if (lambdaFlag) { + fillV0QAHist(collision, v0, tracks); + histos.fill(HIST("Lambda/hInvMassVsPt"), cent, v0.mLambda(), v0.pt()); + histos.fill(HIST("Lambda/Flow/hQuA"), cent, v0.yLambda(), v0.mLambda(), v1a); + histos.fill(HIST("Lambda/Flow/hQuC"), cent, v0.yLambda(), v0.mLambda(), v1c); + } else if (antiLambdaFlag) { + fillV0QAHist(collision, v0, tracks); + histos.fill(HIST("AntiLambda/hInvMassVsPt"), cent, v0.mAntiLambda(), v0.pt()); + histos.fill(HIST("AntiLambda/Flow/hQuA"), cent, v0.yLambda(), v0.mAntiLambda(), v1a); + histos.fill(HIST("AntiLambda/Flow/hQuC"), cent, v0.yLambda(), v0.mAntiLambda(), v1c); + } + } + } + PROCESS_SWITCH(FlowEventPlane, processLambdaFlow, "Lambda flow process", false); }; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) From 1e7622dfeb9fe6cd3c4774f81b36264b20699b14 Mon Sep 17 00:00:00 2001 From: Yash Patley <52608802+yashpatley@users.noreply.github.com> Date: Sun, 18 Jan 2026 00:54:34 +0530 Subject: [PATCH 2/4] Update flowEventPlane.cxx --- PWGCF/Flow/Tasks/flowEventPlane.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PWGCF/Flow/Tasks/flowEventPlane.cxx b/PWGCF/Flow/Tasks/flowEventPlane.cxx index 131b754dd18..4c5f0dfdd2b 100644 --- a/PWGCF/Flow/Tasks/flowEventPlane.cxx +++ b/PWGCF/Flow/Tasks/flowEventPlane.cxx @@ -739,7 +739,7 @@ struct IdHadronFlow { bool identifyTrack(T const& track) { // Electron rejection - if (track.tpcNSigmaEl() > cTpcElRejCutMin && track.tpcNSigmaEl() < cTpcElRejCutMax) { + if (std::abs(track.tpcNSigmaPi()) > cTpcRejCut && std::abs(track.tpcNSigmaKa()) > cTpcRejCut && std::abs(track.tpcNSigmaPr()) > cTpcRejCut && track.tpcNSigmaEl() > cTpcElRejCutMin && track.tpcNSigmaEl() < cTpcElRejCutMax) { return false; } From 82ab2eb8e196376a4e6e98d2c323cf7ae1317230 Mon Sep 17 00:00:00 2001 From: Yash Patley <52608802+yashpatley@users.noreply.github.com> Date: Sun, 18 Jan 2026 14:32:57 +0530 Subject: [PATCH 3/4] Update flowEventPlane.cxx Removed mass calculation for Lambda and AntiLambda. --- PWGCF/Flow/Tasks/flowEventPlane.cxx | 7 ------- 1 file changed, 7 deletions(-) diff --git a/PWGCF/Flow/Tasks/flowEventPlane.cxx b/PWGCF/Flow/Tasks/flowEventPlane.cxx index 4c5f0dfdd2b..512c774e531 100644 --- a/PWGCF/Flow/Tasks/flowEventPlane.cxx +++ b/PWGCF/Flow/Tasks/flowEventPlane.cxx @@ -1030,13 +1030,6 @@ struct FlowEventPlane { // daugthers auto postrack = v0.template posTrack_as(); auto negtrack = v0.template negTrack_as(); - float mass = 0.; - - if constexpr (part == kLambda) { - mass = v0.mLambda(); - } else { - mass = v0.mAntiLambda(); - } // ctau float ctau = v0.distovertotmom(col.posX(), col.posY(), col.posZ()) * MassLambda0; From 1496a8c465c228ac15da4de4275ab2b8ddc7e81a Mon Sep 17 00:00:00 2001 From: Yash Patley <52608802+yashpatley@users.noreply.github.com> Date: Thu, 22 Jan 2026 01:33:52 +0530 Subject: [PATCH 4/4] Update flowEventPlane.cxx --- PWGCF/Flow/Tasks/flowEventPlane.cxx | 544 +++++++++++++--------------- 1 file changed, 255 insertions(+), 289 deletions(-) diff --git a/PWGCF/Flow/Tasks/flowEventPlane.cxx b/PWGCF/Flow/Tasks/flowEventPlane.cxx index 512c774e531..ea895bc8232 100644 --- a/PWGCF/Flow/Tasks/flowEventPlane.cxx +++ b/PWGCF/Flow/Tasks/flowEventPlane.cxx @@ -60,27 +60,13 @@ DECLARE_SOA_TABLE(ColSPExt, "AOD", "COLSPEXT", o2::soa::Index<>, namespace tracksid { -DECLARE_SOA_INDEX_COLUMN(Collision, collision); -DECLARE_SOA_COLUMN(Sign, sign, int8_t); -DECLARE_SOA_COLUMN(Px, px, float); -DECLARE_SOA_COLUMN(Py, py, float); -DECLARE_SOA_COLUMN(Pz, pz, float); -DECLARE_SOA_COLUMN(Pt, pt, float); -DECLARE_SOA_COLUMN(Eta, eta, float); -DECLARE_SOA_COLUMN(Phi, phi, float); +DECLARE_SOA_COLUMN(IsCharged, isCharged, bool); DECLARE_SOA_COLUMN(IsPion, isPion, bool); DECLARE_SOA_COLUMN(IsKaon, isKaon, bool); DECLARE_SOA_COLUMN(IsProton, isProton, bool); } // namespace tracksid DECLARE_SOA_TABLE(TracksId, "AOD", "TRACKSID", o2::soa::Index<>, - aod::tracksid::CollisionId, - tracksid::Sign, - tracksid::Px, - tracksid::Py, - tracksid::Pz, - tracksid::Pt, - tracksid::Eta, - tracksid::Phi, + tracksid::IsCharged, tracksid::IsPion, tracksid::IsKaon, tracksid::IsProton); @@ -128,6 +114,7 @@ enum V0Type { struct SpectatorPlaneTableProducer { // Table producer Produces colSPExtTable; + Produces tracksIdTable; // Configurables // Collisions @@ -170,6 +157,26 @@ struct SpectatorPlaneTableProducer { Configurable cCcdbUrl{"cCcdbUrl", "http://ccdb-test.cern.ch:8080", "url of ccdb"}; Configurable cCcdbPath{"cCcdbPath", "Users/y/ypatley/DFOO", "Path for ccdb-object"}; + // Tracks + Configurable cTrackMinPt{"cTrackMinPt", 0.1, "p_{T} minimum"}; + Configurable cTrackMaxPt{"cTrackMaxPt", 10.0, "p_{T} maximum"}; + Configurable cNEtaBins{"cNEtaBins", 7, "# of eta bins"}; + Configurable cTrackEtaCut{"cTrackEtaCut", 0.8, "Pseudorapidity cut"}; + Configurable cTrackGlobal{"cTrackGlobal", true, "Global Track"}; + Configurable cTrackDcaXYCut{"cTrackDcaXYCut", 0.1, "DcaXY Cut"}; + Configurable cTrackDcaZCut{"cTrackDcaZCut", 1., "DcaXY Cut"}; + + // Track PID + Configurable cTpcElRejCutMin{"cTpcElRejCutMin", -3., "Electron Rejection Cut Minimum"}; + Configurable cTpcElRejCutMax{"cTpcElRejCutMax", 5., "Electron Rejection Cut Maximum"}; + Configurable cTpcNSigmaCut{"cTpcNSigmaCut", 2, "TPC NSigma Cut"}; + Configurable cTpcRejCut{"cTpcRejCut", 3, "TPC Rej Cut"}; + Configurable cTofNSigmaCut{"cTofNSigmaCut", 2, "TOF NSigma Cut"}; + Configurable cTofRejCut{"cTofRejCut", 3, "TOF Rej Cut"}; + Configurable cPionPtCut{"cPionPtCut", 0.6, "Pion TPC pT cutoff"}; + Configurable cKaonPtCut{"cKaonPtCut", 0.6, "Kaon TPC pT cutoff"}; + Configurable cProtonPtCut{"cProtonPtCut", 1.1, "Proton TPC pT cutoff"}; + // Initialize CCDB Service Service ccdbService; @@ -283,6 +290,9 @@ struct SpectatorPlaneTableProducer { histos.add("Checks/hYaYc", "Y^{A}_{1}Y^{C}_{1}", kTProfile, {axisCent}); histos.add("Checks/hXaYc", "X^{A}_{1}Y^{C}_{1}", kTProfile, {axisCent}); histos.add("Checks/hYaXc", "Y^{A}_{1}X^{C}_{1}", kTProfile, {axisCent}); + + // Directed flow QXY vector + histos.add("DF/hQaQc", "X^{A}_{1}X^{C}_{1} + Y^{A}_{1}Y^{C}_{1}", kTProfile, {axisCent}); } template @@ -591,115 +601,6 @@ struct SpectatorPlaneTableProducer { return true; } - using BCsRun3 = soa::Join; - using CollisionsRun3 = soa::Join; - - void process(CollisionsRun3::iterator const& collision, BCsRun3 const&, aod::Zdcs const&) - { - // Analyze collision and get Spectator Plane Vector - std::array vSP = {0., 0., 0., 0.}; - bool colSPExtFlag = analyzeCollision(collision, vSP); - - // Update run number - lRunNum = cRunNum; - - // Fill histograms if SP flag is true - if (colSPExtFlag) { - // Evaluate spectator plane angle and [X,Y] correlations - float psiA = std::atan2(vSP[kYa], vSP[kXa]); - float psiC = std::atan2(vSP[kYc], vSP[kXc]); - histos.fill(HIST("Checks/hPsiSPA"), cent, psiA); - histos.fill(HIST("Checks/hPsiSPC"), cent, psiC); - histos.fill(HIST("Checks/hCosPsiSPAC"), cent, std::cos(psiA - psiC)); - histos.fill(HIST("Checks/hSinPsiSPAC"), cent, std::sin(psiA - psiC)); - histos.fill(HIST("Checks/hXaXc"), cent, (vSP[kXa] * vSP[kXc])); - histos.fill(HIST("Checks/hYaYc"), cent, (vSP[kYa] * vSP[kYc])); - histos.fill(HIST("Checks/hXaYc"), cent, (vSP[kXa] * vSP[kYc])); - histos.fill(HIST("Checks/hYaXc"), cent, (vSP[kYa] * vSP[kXc])); - } - - // Fill table - colSPExtTable(colSPExtFlag, vSP[kXa], vSP[kYa], vSP[kXc], vSP[kYc]); - } -}; - -struct IdHadronFlow { - // Table producer - Produces tracksIdTable; - - // Tracks - Configurable cTrackMinPt{"cTrackMinPt", 0.1, "p_{T} minimum"}; - Configurable cTrackMaxPt{"cTrackMaxPt", 10.0, "p_{T} maximum"}; - Configurable cNEtaBins{"cNEtaBins", 7, "# of eta bins"}; - Configurable cTrackEtaCut{"cTrackEtaCut", 0.8, "Pseudorapidity cut"}; - Configurable cTrackGlobal{"cTrackGlobal", true, "Global Track"}; - Configurable cTrackDcaXYCut{"cTrackDcaXYCut", 0.1, "DcaXY Cut"}; - Configurable cTrackDcaZCut{"cTrackDcaZCut", 1., "DcaXY Cut"}; - - // Track PID - Configurable cTpcElRejCutMin{"cTpcElRejCutMin", -3., "Electron Rejection Cut Minimum"}; - Configurable cTpcElRejCutMax{"cTpcElRejCutMax", 5., "Electron Rejection Cut Maximum"}; - Configurable cTpcNSigmaCut{"cTpcNSigmaCut", 2, "TPC NSigma Cut"}; - Configurable cTpcRejCut{"cTpcRejCut", 3, "TPC Rej Cut"}; - Configurable cTofNSigmaCut{"cTofNSigmaCut", 2, "TOF NSigma Cut"}; - Configurable cTofRejCut{"cTofRejCut", 3, "TOF Rej Cut"}; - Configurable cPionPtCut{"cPionPtCut", 0.6, "Pion TPC pT cutoff"}; - Configurable cKaonPtCut{"cKaonPtCut", 0.6, "Kaon TPC pT cutoff"}; - Configurable cProtonPtCut{"cProtonPtCut", 1.1, "Proton TPC pT cutoff"}; - - // Flag to fill histograms - Configurable cFillIdHist{"cFillIdHist", false, "Flag to fill histograms"}; - - // Histogram registry: an object to hold your histograms - HistogramRegistry histos{"histos", {}, OutputObjHandlingPolicy::AnalysisObject}; - - // Global objects - float cent = 0.; - - void init(InitContext const&) - { - // Define axes - const AxisSpec axisCent{100, 0., 100, "FT0C%"}; - - const AxisSpec axisXYac{600, -6, 6, "Q^{t}Q^{p}"}; - const AxisSpec axisV1{400, -4, 4, "v_{1}"}; - - const AxisSpec axisTrackPt{100, 0., 10., "p_{T} (GeV/#it{c})"}; - const AxisSpec axisTrackEta{cNEtaBins, -0.8, 0.8, "#eta"}; - const AxisSpec axisTrackDcaXY{60, -0.15, 0.15, "DCA_{XY}"}; - const AxisSpec axisTrackDcaZ{230, -1.15, 1.15, "DCA_{XY}"}; - const AxisSpec axisTrackdEdx{360, 20, 200, "#frac{dE}{dx}"}; - const AxisSpec axisTrackNSigma{161, -4.025, 4.025, {"n#sigma"}}; - - // Create histograms - // Track QA - histos.add("TrackQA/hPtDcaXY", "DCA_{XY} vs p_{T}", kTH2F, {axisTrackPt, axisTrackDcaXY}); - histos.add("TrackQA/hPtDcaZ", "DCA_{Z} vs p_{T}", kTH2F, {axisTrackPt, axisTrackDcaZ}); - histos.add("TrackQA/hTrackTPCdEdX", "hTrackTPCdEdX", kTH2F, {axisTrackPt, axisTrackdEdx}); - - // Charged particle directed flow - histos.add("DF/hQaQc", "X^{A}_{1}X^{C}_{1} + Y^{A}_{1}Y^{C}_{1}", kTProfile, {axisCent}); - histos.add("DF/hAQu", "u_{x}X^{A}_{1} + u_{y}Y^{A}_{1}", kTProfile2D, {axisCent, axisTrackEta}); - histos.add("DF/hCQu", "u_{x}X^{C}_{1} + u_{y}Y^{C}_{1}", kTProfile2D, {axisCent, axisTrackEta}); - histos.add("DF/hAQuPos", "u_{x}X^{A}_{1} + u_{y}Y^{A}_{1}", kTProfile2D, {axisCent, axisTrackEta}); - histos.add("DF/hCQuPos", "u_{x}X^{C}_{1} + u_{y}Y^{C}_{1}", kTProfile2D, {axisCent, axisTrackEta}); - histos.add("DF/hAQuNeg", "u_{x}X^{A}_{1} + u_{y}Y^{A}_{1}", kTProfile2D, {axisCent, axisTrackEta}); - histos.add("DF/hCQuNeg", "u_{x}X^{C}_{1} + u_{y}Y^{C}_{1}", kTProfile2D, {axisCent, axisTrackEta}); - - // Identified particle - if (cFillIdHist) { - histos.add("PartId/Pion/hdEdX", "PartId/Pion/hdEdX", kTH2F, {axisTrackPt, axisTrackdEdx}); - histos.add("PartId/Pion/hTPCNSigma", "PartId/Pion/hTPCNSigma", kTH2F, {axisTrackPt, axisTrackNSigma}); - histos.add("PartId/Pion/hTOFNSigma", "PartId/Pion/hTOFNSigma", kTH2F, {axisTrackPt, axisTrackNSigma}); - histos.add("PartId/Pion/hAQuPos", "PartId/Pion/hAQuPos", kTProfile2D, {axisCent, axisTrackEta}); - histos.add("PartId/Pion/hAQuNeg", "PartId/Pion/hAQuNeg", kTProfile2D, {axisCent, axisTrackEta}); - histos.add("PartId/Pion/hCQuPos", "PartId/Pion/hCQuPos", kTProfile2D, {axisCent, axisTrackEta}); - histos.add("PartId/Pion/hCQuNeg", "PartId/Pion/hCQuNeg", kTProfile2D, {axisCent, axisTrackEta}); - histos.addClone("PartId/Pion/", "PartId/Kaon/"); - histos.addClone("PartId/Pion/", "PartId/Proton/"); - } - } - // Track Selection template bool selectTrack(T const& track) @@ -712,7 +613,7 @@ struct IdHadronFlow { return false; } - if (std::abs(track.dcaXY()) > cTrackDcaXYCut || std::abs(track.dcaZ()) > cTrackDcaZCut) { + if (std::abs(track.dcaXY()) >= cTrackDcaXYCut || std::abs(track.dcaZ()) >= cTrackDcaZCut) { return false; } @@ -762,129 +663,62 @@ struct IdHadronFlow { return retFlag; } - template - void getIdHadronFlow(float const& cent, T const& track, float const& v1a, float const& v1c) - { - static constexpr std::string_view SubDir[] = {"Pion/", "Kaon/", "Proton/"}; - float tpcNsigma = 0., tofNsigma = 0.; - if (part == kPi) { - tpcNsigma = track.tpcNSigmaPi(); - tofNsigma = track.tofNSigmaPi(); - } else if (part == kKa) { - tpcNsigma = track.tpcNSigmaKa(); - tofNsigma = track.tofNSigmaKa(); - } else if (part == kPr) { - tpcNsigma = track.tpcNSigmaPr(); - tofNsigma = track.tofNSigmaPr(); - } else { - return; - } - histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hdEdX"), track.pt(), track.tpcSignal()); - histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hTPCNSigma"), track.pt(), tpcNsigma); - if (track.hasTOF()) { - histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hTOFNSigma"), track.pt(), tofNsigma); - } - if (track.sign() > 0) { - histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hAQuPos"), cent, track.eta(), v1a); - histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hCQuPos"), cent, track.eta(), v1c); - } else { - histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hAQuNeg"), cent, track.eta(), v1a); - histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hCQuNeg"), cent, track.eta(), v1c); - } - } - - template - void fillTrackHist(T const& track) - { - histos.fill(HIST("TrackQA/hPtDcaZ"), track.pt(), track.dcaZ()); - histos.fill(HIST("TrackQA/hPtDcaXY"), track.pt(), track.dcaXY()); - } - - using CollisionsRun3 = soa::Join; + using BCsRun3 = soa::Join; + using CollisionsRun3 = soa::Join; using Tracks = soa::Join; - void processDummy(CollisionsRun3::iterator const&) {} - - PROCESS_SWITCH(IdHadronFlow, processDummy, "Dummy process", true); - - void processIdHadronFlow(CollisionsRun3::iterator const& collision, Tracks const& tracks) + void processSpectatorPlane(CollisionsRun3::iterator const& collision, BCsRun3 const&, aod::Zdcs const&) { - // Check collision - if (!collision.selColFlag()) { - return; - } - - // Set centrality - cent = collision.centFT0C(); - - // Flow vectors + // Analyze collision and get Spectator Plane Vector std::array vSP = {0., 0., 0., 0.}; - vSP[kXa] = collision.xa(); - vSP[kYa] = collision.ya(); - vSP[kXc] = collision.xc(); - vSP[kYc] = collision.yc(); + bool colSPExtFlag = analyzeCollision(collision, vSP); - // Directed flow QXY vector - float qac = (vSP[kXa] * vSP[kXc]) + (vSP[kYa] * vSP[kYc]); - histos.fill(HIST("DF/hQaQc"), cent, qac); + // Update run number + lRunNum = cRunNum; - // Loop over tracks - float ux = 0., uy = 0., v1a = 0., v1c = 0.; - for (auto const& track : tracks) { - // Select track - if (!selectTrack(track)) { - tracksIdTable(collision.index(), track.sign(), track.px(), track.py(), track.pz(), track.pt(), track.eta(), track.phi(), false, false, false); - continue; - } + // Fill histograms if SP flag is true + if (colSPExtFlag) { + // Evaluate spectator plane angle and [X,Y] correlations + float psiA = std::atan2(vSP[kYa], vSP[kXa]); + float psiC = std::atan2(vSP[kYc], vSP[kXc]); + histos.fill(HIST("Checks/hPsiSPA"), cent, psiA); + histos.fill(HIST("Checks/hPsiSPC"), cent, psiC); + histos.fill(HIST("Checks/hCosPsiSPAC"), cent, std::cos(psiA - psiC)); + histos.fill(HIST("Checks/hSinPsiSPAC"), cent, std::sin(psiA - psiC)); + histos.fill(HIST("Checks/hXaXc"), cent, (vSP[kXa] * vSP[kXc])); + histos.fill(HIST("Checks/hYaYc"), cent, (vSP[kYa] * vSP[kYc])); + histos.fill(HIST("Checks/hXaYc"), cent, (vSP[kXa] * vSP[kYc])); + histos.fill(HIST("Checks/hYaXc"), cent, (vSP[kYa] * vSP[kXc])); - // Fill track QA - fillTrackHist(track); + // Directed flow QXY vector + float qac = (vSP[kXa] * vSP[kXc]) + (vSP[kYa] * vSP[kYc]); + histos.fill(HIST("DF/hQaQc"), cent, qac); + } - // Get directed flow - ux = std::cos(track.phi()); - uy = std::sin(track.phi()); - v1a = ux * vSP[kXa] + uy * vSP[kYa]; - v1c = ux * vSP[kXc] + uy * vSP[kYc]; + // Fill table + colSPExtTable(colSPExtFlag, vSP[kXa], vSP[kYa], vSP[kXc], vSP[kYc]); + } - // Charged particle directed flow - histos.fill(HIST("DF/hAQu"), cent, track.eta(), v1a); - histos.fill(HIST("DF/hCQu"), cent, track.eta(), v1c); - if (track.sign() > 0) { - histos.fill(HIST("DF/hAQuPos"), cent, track.eta(), v1a); - histos.fill(HIST("DF/hCQuPos"), cent, track.eta(), v1c); - } else { - histos.fill(HIST("DF/hAQuNeg"), cent, track.eta(), v1a); - histos.fill(HIST("DF/hCQuNeg"), cent, track.eta(), v1c); - } + PROCESS_SWITCH(SpectatorPlaneTableProducer, processSpectatorPlane, "Spectator Plane Process", true); - // Identify track + void processIdHadrons(Tracks const& tracks) + { + for (auto const& track : tracks) { + bool chargedFlag = selectTrack(track); bool pionFlag = identifyTrack(track); bool kaonFlag = identifyTrack(track); bool protonFlag = identifyTrack(track); - - if (cFillIdHist) { - // Pion - if (pionFlag) { - getIdHadronFlow(cent, track, v1a, v1c); - } - // Kaon - if (kaonFlag) { - getIdHadronFlow(cent, track, v1a, v1c); - } - // Proton - if (protonFlag) { - getIdHadronFlow(cent, track, v1a, v1c); - } - } - - // Fill track table - tracksIdTable(collision.index(), track.sign(), track.px(), track.py(), track.pz(), track.pt(), track.eta(), track.phi(), pionFlag, kaonFlag, protonFlag); + tracksIdTable(chargedFlag, pionFlag, kaonFlag, protonFlag); } } - PROCESS_SWITCH(IdHadronFlow, processIdHadronFlow, "Identified hadron flow process", false); + + PROCESS_SWITCH(SpectatorPlaneTableProducer, processIdHadrons, "Hadron Id Process", false); }; struct FlowEventPlane { + // Tracks + Configurable cNEtaBins{"cNEtaBins", 7, "# of eta bins"}; + // Resonance Configurable cNRapBins{"cNRapBins", 5, "# of y bins"}; Configurable cNInvMassBins{"cNInvMassBins", 500, "# of m bins"}; @@ -917,17 +751,29 @@ struct FlowEventPlane { // Global objects float cent = 0.; + std::array vSP = {0., 0., 0., 0.}; void init(InitContext const&) { // Define axes const AxisSpec axisCent{100, 0., 100, "FT0C%"}; + + const AxisSpec axisXYac{600, -6, 6, "Q^{t}Q^{p}"}; + const AxisSpec axisV1{400, -4, 4, "v_{1}"}; + const AxisSpec axisTrackPt{100, 0., 10., "p_{T} (GeV/#it{c})"}; + const AxisSpec axisTrackEta{cNEtaBins, -0.8, 0.8, "#eta"}; + const AxisSpec axisTrackDcaXY{60, -0.15, 0.15, "DCA_{XY}"}; + const AxisSpec axisTrackDcaZ{230, -1.15, 1.15, "DCA_{XY}"}; + const AxisSpec axisTrackdEdx{360, 20, 200, "#frac{dE}{dx}"}; + const AxisSpec axisTrackNSigma{161, -4.025, 4.025, {"n#sigma"}}; + const AxisSpec axisTrackRap{cNRapBins, -0.5, 0.5, "y"}; const AxisSpec axisInvMass{cNInvMassBins, 0.87, 1.12, "M_{KK} (GeV/#it{c}^{2}"}; - const AxisSpec axisMomPID(80, 0, 4, "p (GeV/#it{c})"); + const AxisSpec axisMomPID(80, 0, 4, "p_{T} (GeV/#it{c})"); const AxisSpec axisNsigma(401, -10.025, 10.025, {"n#sigma"}); const AxisSpec axisdEdx(360, 20, 200, "#frac{dE}{dx}"); + const AxisSpec axisTrackTofSignal(240, 0, 1.2, "#beta"); const AxisSpec axisRadius(2000, 0, 200, "r(cm)"); const AxisSpec axisCosPA(300, 0.97, 1.0, "cos(#theta_{PA})"); const AxisSpec axisDcaV0PV(1000, 0., 10., "dca (cm)"); @@ -938,40 +784,128 @@ struct FlowEventPlane { const AxisSpec axisQtarm(40, 0, 0.4, "q_{T}"); const AxisSpec axisLambdaPt(50, 0, 10, "p_{T} (GeV/#it{c})"); const AxisSpec axisLambdaInvMass(140, 1.08, 1.15, "M_{p#pi} (GeV/#it{c}^{2})"); - const AxisSpec axisXYac{600, -6, 6, "Q^{t}Q^{p}"}; - const AxisSpec axisV1{400, -4, 4, "v_{1}"}; // Create histograms + // Charged particles + if (doprocessChargedFlow) { + histos.add("TrackQA/hPtDcaXY", "DCA_{XY} vs p_{T}", kTH2F, {axisTrackPt, axisTrackDcaXY}); + histos.add("TrackQA/hPtDcaZ", "DCA_{Z} vs p_{T}", kTH2F, {axisTrackPt, axisTrackDcaZ}); + histos.add("TrackQA/hTrackTPCdEdX", "hTrackTPCdEdX", kTH2F, {axisMomPID, axisdEdx}); + histos.add("DF/hQaQc", "X^{A}_{1}X^{C}_{1} + Y^{A}_{1}Y^{C}_{1}", kTProfile, {axisCent}); + histos.add("DF/hAQu", "u_{x}X^{A}_{1} + u_{y}Y^{A}_{1}", kTProfile2D, {axisCent, axisTrackEta}); + histos.add("DF/hCQu", "u_{x}X^{C}_{1} + u_{y}Y^{C}_{1}", kTProfile2D, {axisCent, axisTrackEta}); + histos.add("DF/hAQuPos", "u_{x}X^{A}_{1} + u_{y}Y^{A}_{1}", kTProfile2D, {axisCent, axisTrackEta}); + histos.add("DF/hCQuPos", "u_{x}X^{C}_{1} + u_{y}Y^{C}_{1}", kTProfile2D, {axisCent, axisTrackEta}); + histos.add("DF/hAQuNeg", "u_{x}X^{A}_{1} + u_{y}Y^{A}_{1}", kTProfile2D, {axisCent, axisTrackEta}); + histos.add("DF/hCQuNeg", "u_{x}X^{C}_{1} + u_{y}Y^{C}_{1}", kTProfile2D, {axisCent, axisTrackEta}); + } + + // Identified hadrons + if (doprocessIdFlow) { + histos.add("PartId/Pion/hdEdX", "dE/dx vs pT", kTH2F, {axisMomPID, axisTrackdEdx}); + histos.add("PartId/Pion/hTOFSignal", "#beta_{TOF} vs p_{T}", kTH2F, {axisMomPID, axisTrackTofSignal}); + histos.add("PartId/Pion/hTPCNSigma", "n#sigma_{TPC} vs p_{T}", kTH2F, {axisMomPID, axisTrackNSigma}); + histos.add("PartId/Pion/hTOFNSigma", "n#sigma_{TOF} vs p_{T}", kTH2F, {axisMomPID, axisTrackNSigma}); + histos.add("PartId/Pion/hAQuPos", "PartId/Pion/hAQuPos", kTProfile2D, {axisCent, axisTrackEta}); + histos.add("PartId/Pion/hAQuNeg", "PartId/Pion/hAQuNeg", kTProfile2D, {axisCent, axisTrackEta}); + histos.add("PartId/Pion/hCQuPos", "PartId/Pion/hCQuPos", kTProfile2D, {axisCent, axisTrackEta}); + histos.add("PartId/Pion/hCQuNeg", "PartId/Pion/hCQuNeg", kTProfile2D, {axisCent, axisTrackEta}); + histos.addClone("PartId/Pion/", "PartId/Kaon/"); + histos.addClone("PartId/Pion/", "PartId/Proton/"); + } + // Resonance - histos.add("Reso/Phi/hSigCentPtInvMass", "hUSCentPtInvMass", kTH3F, {axisCent, axisTrackPt, axisInvMass}); - histos.add("Reso/Phi/hBkgCentPtInvMass", "hLSCentPtInvMass", kTH3F, {axisCent, axisTrackPt, axisInvMass}); - histos.add("Reso/Phi/Sig/hPhiQuA", "hPhiQuA", kTProfile3D, {axisCent, axisTrackRap, axisInvMass}); - histos.add("Reso/Phi/Sig/hPhiQuC", "hPhiQuC", kTProfile3D, {axisCent, axisTrackRap, axisInvMass}); - histos.add("Reso/Phi/Bkg/hPhiQuA", "hPhiQuA", kTProfile3D, {axisCent, axisTrackRap, axisInvMass}); - histos.add("Reso/Phi/Bkg/hPhiQuC", "hPhiQuC", kTProfile3D, {axisCent, axisTrackRap, axisInvMass}); - - // Lambda histograms - // QA Lambda - histos.add("Lambda/QA/hQtVsAlpha", "Armentros-Podolanski Plot", kTH2F, {axisAlpha, axisQtarm}); - histos.add("Lambda/QA/hDcaV0Dau", "DCA between V0 daughters", kTH1F, {axisDcaDau}); - histos.add("Lambda/QA/hDcaPosToPv", "DCA positive prong to PV", kTH1F, {axisDcaProngPV}); - histos.add("Lambda/QA/hDcaNegToPv", "DCA negative prong to PV", kTH1F, {axisDcaProngPV}); - histos.add("Lambda/QA/hDcaV0ToPv", "DCA V0 to PV", kTH1F, {axisDcaV0PV}); - histos.add("Lambda/QA/hCosPa", "cos(#theta_{PA})", kTH1F, {axisCosPA}); - histos.add("Lambda/QA/hRxy", "V_{0} Decay Radius in XY plane", kTH1F, {axisRadius}); - histos.add("Lambda/QA/hCTau", "V_{0} c#tau", kTH1F, {axisCTau}); - histos.add("Lambda/QA/hPosdEdXVsP", "TPC Signal Pos-Prong", kTH2F, {axisMomPID, axisdEdx}); - histos.add("Lambda/QA/hNegdEdXVsP", "TPC Signal Neg-Prong", kTH2F, {axisMomPID, axisdEdx}); - histos.add("Lambda/QA/hPosNsigPrVsP", "TPC n#sigma Pos Prong", kTH2F, {axisMomPID, axisNsigma}); - histos.add("Lambda/QA/hNegNsigPrVsP", "TPC n#sigma Neg Prong", kTH2F, {axisMomPID, axisNsigma}); - histos.add("Lambda/QA/hPosNsigPiVsP", "TPC n#sigma Pos Prong", kTH2F, {axisMomPID, axisNsigma}); - histos.add("Lambda/QA/hNegNsigPiVsP", "TPC n#sigma Neg Prong", kTH2F, {axisMomPID, axisNsigma}); - - histos.add("Lambda/hInvMassVsPt", "hInvMassVsPt", kTH3F, {axisCent, axisLambdaInvMass, axisLambdaPt}); - histos.add("Lambda/Flow/hQuA", "hPhiQuA", kTProfile3D, {axisCent, axisTrackRap, axisLambdaInvMass}); - histos.add("Lambda/Flow/hQuC", "hPhiQuC", kTProfile3D, {axisCent, axisTrackRap, axisLambdaInvMass}); - - histos.addClone("Lambda/", "AntiLambda/"); + if (doprocessResoFlow) { + histos.add("Reso/Phi/hSigCentPtInvMass", "hUSCentPtInvMass", kTH3F, {axisCent, axisTrackPt, axisInvMass}); + histos.add("Reso/Phi/hBkgCentPtInvMass", "hLSCentPtInvMass", kTH3F, {axisCent, axisTrackPt, axisInvMass}); + histos.add("Reso/Phi/Sig/hPhiQuA", "hPhiQuA", kTProfile3D, {axisCent, axisTrackRap, axisInvMass}); + histos.add("Reso/Phi/Sig/hPhiQuC", "hPhiQuC", kTProfile3D, {axisCent, axisTrackRap, axisInvMass}); + histos.add("Reso/Phi/Bkg/hPhiQuA", "hPhiQuA", kTProfile3D, {axisCent, axisTrackRap, axisInvMass}); + histos.add("Reso/Phi/Bkg/hPhiQuC", "hPhiQuC", kTProfile3D, {axisCent, axisTrackRap, axisInvMass}); + } + + // Lambda + if (doprocessLambdaFlow) { + histos.add("Lambda/QA/hQtVsAlpha", "Armentros-Podolanski Plot", kTH2F, {axisAlpha, axisQtarm}); + histos.add("Lambda/QA/hDcaV0Dau", "DCA between V0 daughters", kTH1F, {axisDcaDau}); + histos.add("Lambda/QA/hDcaPosToPv", "DCA positive prong to PV", kTH1F, {axisDcaProngPV}); + histos.add("Lambda/QA/hDcaNegToPv", "DCA negative prong to PV", kTH1F, {axisDcaProngPV}); + histos.add("Lambda/QA/hDcaV0ToPv", "DCA V0 to PV", kTH1F, {axisDcaV0PV}); + histos.add("Lambda/QA/hCosPa", "cos(#theta_{PA})", kTH1F, {axisCosPA}); + histos.add("Lambda/QA/hRxy", "V_{0} Decay Radius in XY plane", kTH1F, {axisRadius}); + histos.add("Lambda/QA/hCTau", "V_{0} c#tau", kTH1F, {axisCTau}); + histos.add("Lambda/QA/hPosdEdXVsP", "TPC Signal Pos-Prong", kTH2F, {axisMomPID, axisdEdx}); + histos.add("Lambda/QA/hNegdEdXVsP", "TPC Signal Neg-Prong", kTH2F, {axisMomPID, axisdEdx}); + histos.add("Lambda/QA/hPosNsigPrVsP", "TPC n#sigma Pos Prong", kTH2F, {axisMomPID, axisNsigma}); + histos.add("Lambda/QA/hNegNsigPrVsP", "TPC n#sigma Neg Prong", kTH2F, {axisMomPID, axisNsigma}); + histos.add("Lambda/QA/hPosNsigPiVsP", "TPC n#sigma Pos Prong", kTH2F, {axisMomPID, axisNsigma}); + histos.add("Lambda/QA/hNegNsigPiVsP", "TPC n#sigma Neg Prong", kTH2F, {axisMomPID, axisNsigma}); + histos.add("Lambda/hInvMassVsPt", "hInvMassVsPt", kTH3F, {axisCent, axisLambdaInvMass, axisLambdaPt}); + histos.add("Lambda/Flow/hQuA", "hPhiQuA", kTProfile3D, {axisCent, axisTrackRap, axisLambdaInvMass}); + histos.add("Lambda/Flow/hQuC", "hPhiQuC", kTProfile3D, {axisCent, axisTrackRap, axisLambdaInvMass}); + histos.addClone("Lambda/", "AntiLambda/"); + } + } + + template + bool selCollision(C const& collision, std::array& v) + { + // Check collision + if (!collision.selColFlag()) { + return false; + } + + // Set centrality + cent = collision.centFT0C(); + + // Flow vectors + v[kXa] = collision.xa(); + v[kYa] = collision.ya(); + v[kXc] = collision.xc(); + v[kYc] = collision.yc(); + + return true; + } + + template + void getIdFlow(T const& tracks) + { + float ux = 0., uy = 0., v1a = 0., v1c = 0.; + float tpcNsigma = 0., tofNsigma = 0.; + for (auto const& track : tracks) { + static constexpr std::string_view SubDir[] = {"Pion/", "Kaon/", "Proton/"}; + if (part == kPi) { + tpcNsigma = track.tpcNSigmaPi(); + tofNsigma = track.tofNSigmaPi(); + } else if (part == kKa) { + tpcNsigma = track.tpcNSigmaKa(); + tofNsigma = track.tofNSigmaKa(); + } else if (part == kPr) { + tpcNsigma = track.tpcNSigmaPr(); + tofNsigma = track.tofNSigmaPr(); + } else { + return; + } + histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hdEdX"), track.pt(), track.tpcSignal()); + histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hTPCNSigma"), track.pt(), tpcNsigma); + if (track.hasTOF()) { + histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hTOFSignal"), track.pt(), track.beta()); + histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hTOFNSigma"), track.pt(), tofNsigma); + } + + ux = std::cos(track.phi()); + uy = std::sin(track.phi()); + v1a = ux * vSP[kXa] + uy * vSP[kYa]; + v1c = ux * vSP[kXc] + uy * vSP[kYc]; + + if (track.sign() > 0) { + histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hAQuPos"), cent, track.eta(), v1a); + histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hCQuPos"), cent, track.eta(), v1c); + } else { + histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hAQuNeg"), cent, track.eta(), v1a); + histos.fill(HIST("PartId/") + HIST(SubDir[part]) + HIST("hCQuNeg"), cent, track.eta(), v1c); + } + } } template @@ -1106,36 +1040,80 @@ struct FlowEventPlane { } using CollisionsRun3 = soa::Join; - using Tracks = aod::TracksId; + using Tracks = soa::Join; using TracksV0s = soa::Join; // Partitions SliceCache cache; - Partition kaonTrackPartition = (aod::tracksid::isKaon == true); + Partition chargedTrackPartition = (aod::tracksid::isCharged == true); + Partition pionTrackPartition = (aod::tracksid::isCharged == true) && (aod::tracksid::isPion == true); + Partition kaonTrackPartition = (aod::tracksid::isCharged == true) && (aod::tracksid::isKaon == true); + Partition protonTrackPartition = (aod::tracksid::isCharged == true) && (aod::tracksid::isProton == true); void processDummy(CollisionsRun3::iterator const&) {} PROCESS_SWITCH(FlowEventPlane, processDummy, "Dummy process", true); - void processResoFlow(CollisionsRun3::iterator const& collision, Tracks const&) + void processChargedFlow(CollisionsRun3::iterator const& collision, Tracks const&) { // Check collision - if (!collision.selColFlag()) { + if (!selCollision(collision, vSP)) { return; } + // Loop over tracks + auto chargedTracks = chargedTrackPartition->sliceByCached(aod::track::collisionId, collision.globalIndex(), cache); + float ux = 0., uy = 0., v1a = 0., v1c = 0.; + for (auto const& track : chargedTracks) { + // Fill track QA + histos.fill(HIST("TrackQA/hPtDcaZ"), track.pt(), track.dcaZ()); + histos.fill(HIST("TrackQA/hPtDcaXY"), track.pt(), track.dcaXY()); + histos.fill(HIST("TrackQA/hTrackTPCdEdX"), track.pt(), track.tpcSignal()); - // Set centrality - cent = collision.centFT0C(); + // Get directed flow + ux = std::cos(track.phi()); + uy = std::sin(track.phi()); + v1a = ux * vSP[kXa] + uy * vSP[kYa]; + v1c = ux * vSP[kXc] + uy * vSP[kYc]; - // Flow vectors - std::array vSP = {0., 0., 0., 0.}; - vSP[kXa] = collision.xa(); - vSP[kYa] = collision.ya(); - vSP[kXc] = collision.xc(); - vSP[kYc] = collision.yc(); + // Charged particle directed flow + histos.fill(HIST("DF/hAQu"), cent, track.eta(), v1a); + histos.fill(HIST("DF/hCQu"), cent, track.eta(), v1c); + if (track.sign() > 0) { + histos.fill(HIST("DF/hAQuPos"), cent, track.eta(), v1a); + histos.fill(HIST("DF/hCQuPos"), cent, track.eta(), v1c); + } else { + histos.fill(HIST("DF/hAQuNeg"), cent, track.eta(), v1a); + histos.fill(HIST("DF/hCQuNeg"), cent, track.eta(), v1c); + } + } + } + + PROCESS_SWITCH(FlowEventPlane, processChargedFlow, "Charged particle flow process", false); + + void processIdFlow(CollisionsRun3::iterator const& collision, Tracks const&) + { + if (!selCollision(collision, vSP)) { + return; + } + // Loop over tracks + auto pionTracks = pionTrackPartition->sliceByCached(aod::track::collisionId, collision.globalIndex(), cache); + auto kaonTracks = kaonTrackPartition->sliceByCached(aod::track::collisionId, collision.globalIndex(), cache); + auto protonTracks = protonTrackPartition->sliceByCached(aod::track::collisionId, collision.globalIndex(), cache); + getIdFlow(pionTracks); + getIdFlow(kaonTracks); + getIdFlow(protonTracks); + } + + PROCESS_SWITCH(FlowEventPlane, processIdFlow, "Identified particle flow process", false); + + void processResoFlow(CollisionsRun3::iterator const& collision, Tracks const&) + { + if (!selCollision(collision, vSP)) { + return; + } // Track partitions - auto kaonTracks = kaonTrackPartition->sliceByCached(aod::tracksid::collisionId, collision.globalIndex(), cache); + auto kaonTracks = kaonTrackPartition->sliceByCached(aod::track::collisionId, collision.globalIndex(), cache); // Resonance flow getResoFlow(kaonTracks, kaonTracks, vSP); @@ -1144,21 +1122,10 @@ struct FlowEventPlane { void processLambdaFlow(CollisionsRun3::iterator const& collision, aod::V0Datas const& V0s, TracksV0s const& tracks) { - // Check collision - if (!collision.selColFlag()) { + if (!selCollision(collision, vSP)) { return; } - // Set centrality - cent = collision.centFT0C(); - - // Flow vectors - std::array vSP = {0., 0., 0., 0.}; - vSP[kXa] = collision.xa(); - vSP[kYa] = collision.ya(); - vSP[kXc] = collision.xc(); - vSP[kYc] = collision.yc(); - // Loop over v0s for (auto const& v0 : V0s) { // V0 kinematic selection @@ -1230,6 +1197,5 @@ WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { return WorkflowSpec{ adaptAnalysisTask(cfgc), - adaptAnalysisTask(cfgc), adaptAnalysisTask(cfgc)}; }