From 63d8c1f19d9d5ca9f1cf8bd101518493a6106f91 Mon Sep 17 00:00:00 2001 From: Swati Saha Date: Fri, 24 Apr 2026 23:56:22 +0530 Subject: [PATCH] Update event selection and PID cuts --- .../Tasks/v0ptHadPiKaProt.cxx | 244 ++++++++++++++---- 1 file changed, 198 insertions(+), 46 deletions(-) diff --git a/PWGCF/EbyEFluctuations/Tasks/v0ptHadPiKaProt.cxx b/PWGCF/EbyEFluctuations/Tasks/v0ptHadPiKaProt.cxx index cffbe84540c..ff59948570e 100644 --- a/PWGCF/EbyEFluctuations/Tasks/v0ptHadPiKaProt.cxx +++ b/PWGCF/EbyEFluctuations/Tasks/v0ptHadPiKaProt.cxx @@ -13,7 +13,8 @@ /// \brief Task for analyzing v0(pT) of inclusive hadrons, pions, kaons, and, protons /// \author Swati Saha -#include "Common/CCDB/EventSelectionParams.h" +#include "Common/Core/TrackSelection.h" +#include "Common/Core/trackUtilities.h" #include "Common/DataModel/Centrality.h" #include "Common/DataModel/EventSelection.h" #include "Common/DataModel/Multiplicity.h" @@ -22,37 +23,43 @@ #include "Common/DataModel/PIDResponseTPC.h" #include "Common/DataModel/TrackSelectionTables.h" +#include "CommonConstants/MathConstants.h" +#include "CommonConstants/PhysicsConstants.h" +#include "DataFormatsParameters/GRPMagField.h" +#include "DataFormatsParameters/GRPObject.h" +#include "Framework/ASoAHelpers.h" +#include "Framework/AnalysisDataModel.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 "ReconstructionDataFormats/Track.h" #include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include #include +#include #include -#include +#include +#include +#include +#include +#include +#include #include #include #include +#include +#include #include #include #include -#include #include -#include #include -#include #include -#include #include #include @@ -73,6 +80,8 @@ struct V0ptHadPiKaProt { Configurable ccdbNoLaterThan{"ccdbNoLaterThan", std::chrono::duration_cast(std::chrono::system_clock::now().time_since_epoch()).count(), "latest acceptable timestamp of creation for the object"}; Configurable ccdbUrl{"ccdbUrl", "https://alice-ccdb.cern.ch", "url of the ccdb repository"}; Configurable ccdbPath{"ccdbPath", "Users/s/swati/PhiWeight", "CCDB path to ccdb object containing phi weight in a 3D histogram"}; + Configurable ccdbNoLaterThanPtEff{"ccdbNoLaterThanPtEff", std::chrono::duration_cast(std::chrono::system_clock::now().time_since_epoch()).count(), "latest acceptable timestamp of creation for the object"}; + Configurable ccdbPathPtEff{"ccdbPathPtEff", "Users/s/swati/EfficiencyWeight", "CCDB path to ccdb object containing pt-dependent efficiency"}; enum Particles { PIONS = 0, @@ -110,6 +119,10 @@ struct V0ptHadPiKaProt { Configurable cfgnSigmaOtherParticles{"cfgnSigmaOtherParticles", 3.0f, "PID nSigma cut to remove other particles (default:3)"}; Configurable cfgnSigmaCutTPC{"cfgnSigmaCutTPC", 2.0f, "PID nSigma cut for TPC"}; Configurable cfgnSigmaCutTOF{"cfgnSigmaCutTOF", 2.0f, "PID nSigma cut for TOF"}; + Configurable cfgUseNewSeperationPid{"cfgUseNewSeperationPid", true, "Use seperation based PID cuts (NEW)"}; + Configurable cfgnSigmaCutTPCHigherPt{"cfgnSigmaCutTPCHigherPt", 2.0f, "PID nSigma cut for TPC at higher pt"}; + Configurable cfgnSigmaCutTOFHigherPt{"cfgnSigmaCutTOFHigherPt", 2.0f, "PID nSigma cut for TOF at higher pt"}; + Configurable cfgnSigmaSeperationCut{"cfgnSigmaSeperationCut", 3.5f, "PID nSigma of other species must be greater than the vale"}; Configurable cfgnSigmaCutCombTPCTOF{"cfgnSigmaCutCombTPCTOF", 2.0f, "PID nSigma combined cut for TPC and TOF"}; ConfigurableAxis nchAxis{"nchAxis", {5000, 0.5, 5000.5}, ""}; ConfigurableAxis centAxis{"centAxis", {90, 0., 90.}, "Centrality/Multiplicity percentile bining"}; @@ -122,7 +135,7 @@ struct V0ptHadPiKaProt { Configurable cfgCutPtUpperPID{"cfgCutPtUpperPID", 6.0f, "Higher pT cut for identified particle analysis"}; Configurable cfgCutEta{"cfgCutEta", 0.8f, "absolute Eta cut"}; Configurable cfgCutEtaLeft{"cfgCutEtaLeft", 0.8f, "Left end of eta gap"}; - Configurable cfgCutEtaRight{"cfgCutEtaRight", 0.8f, "Right end of eta gap"}; + Configurable cfgCutEtaRight{"cfgCutEtaRight", 0.8f, "Right end of eta g//ap"}; Configurable cfgNSubsample{"cfgNSubsample", 20, "Number of subsamples"}; Configurable cfgCentralityChoice{"cfgCentralityChoice", 0, "Which centrality estimator? 0-->FT0C, 1-->FT0A, 2-->FT0M, 3-->FV0A"}; Configurable cfgEvSelkNoSameBunchPileup{"cfgEvSelkNoSameBunchPileup", true, "Pileup removal"}; @@ -130,6 +143,12 @@ struct V0ptHadPiKaProt { Configurable cfgEvSelkNoITSROFrameBorder{"cfgEvSelkNoITSROFrameBorder", true, "ITSROFrame border event selection cut"}; Configurable cfgEvSelkNoTimeFrameBorder{"cfgEvSelkNoTimeFrameBorder", true, "TimeFrame border event selection cut"}; Configurable cfgEvSelUseGoodZvtxFT0vsPV{"cfgEvSelUseGoodZvtxFT0vsPV", true, "GoodZvertex and FT0 vs PV cut"}; + + Configurable cfgEvSelUseOcuppancyTimeCut{"cfgEvSelUseOcuppancyTimeCut", true, "Occupancy Time pattern cut"}; + Configurable cfgEvSelSetOcuppancyRange{"cfgEvSelSetOcuppancyRange", true, "Use cut on occupancy range"}; + Configurable cfgMinOccupancy{"cfgMinOccupancy", 0, "min. value of occupancy"}; + Configurable cfgMaxOccupancy{"cfgMaxOccupancy", 3000, "max. value of occupancy"}; + Configurable cfgUseItsPID{"cfgUseItsPID", false, "Use ITS PID for particle identification"}; Configurable cfgPtCutTOF{"cfgPtCutTOF", 0.3f, "Minimum pt to use TOF N-sigma"}; Configurable> nSigmas{"nSigmas", {LongArrayFloat[0], 3, 6, {"TPC", "TOF", "ITS"}, {"pos_pi", "pos_ka", "pos_pr", "neg_pi", "neg_ka", "neg_pr"}}, "Labeled array for n-sigma values for TPC, TOF, ITS for pions, kaons, protons (positive and negative)"}; @@ -138,6 +157,7 @@ struct V0ptHadPiKaProt { Configurable cfgCutPtMaxForV02{"cfgCutPtMaxForV02", 3.0f, "Max. pT for v02(pT)"}; Configurable cfgCutEtaWindowB{"cfgCutEtaWindowB", 0.4f, "value of x in |eta| cfgLoadPhiWeights{"cfgLoadPhiWeights", false, "Load phi weights from CCDB to take care of non-uniform acceptance"}; + Configurable cfgLoadPtEffWeights{"cfgLoadPtEffWeights", false, "Load pt-dependent efficiency weights from CCDB to take care of detector inefficiency"}; // pT dep DCAxy and DCAz cuts Configurable cfgUsePtDepDCAxy{"cfgUsePtDepDCAxy", true, "Use pt-dependent DCAxy cut"}; @@ -193,7 +213,12 @@ struct V0ptHadPiKaProt { TRandom3* funRndm = new TRandom3(0); // Phi weight histograms initialization - TH2F* hWeightPhiFunctionVzEtaPhi = nullptr; + TH3D* hWeightPhiFunctionVzEtaPhi = nullptr; + // Efficiency of diff. particle histograms initialization + TH1D* hEffAllCharged = nullptr; + TH1D* hEffPion = nullptr; + TH1D* hEffKaon = nullptr; + TH1D* hEffProton = nullptr; // Filter command*********** Filter collisionFilter = nabs(aod::collision::posZ) < cfgCutVertex; @@ -235,7 +260,7 @@ struct V0ptHadPiKaProt { // Loading phi weight histograms from CCDB if (cfgLoadPhiWeights) { - // Accessing eff histograms + // Accessing phi weight histograms ccdb->setURL(ccdbUrl.value); // Enabling object caching, otherwise each call goes to the CCDB server ccdb->setCaching(true); @@ -245,10 +270,26 @@ struct V0ptHadPiKaProt { ccdb->setCreatedNotAfter(ccdbNoLaterThan.value); LOGF(info, "Getting object %s", ccdbPath.value.data()); TList* lst = ccdb->getForTimeStamp(ccdbPath.value, ccdbNoLaterThan.value); - hWeightPhiFunctionVzEtaPhi = reinterpret_cast(lst->FindObject("hWeightPhiFunctionVzEtaPhi")); + hWeightPhiFunctionVzEtaPhi = reinterpret_cast(lst->FindObject("hWeightPhiFunctionVzEtaPhi")); if (!hWeightPhiFunctionVzEtaPhi) LOGF(info, "FATAL!! could not get phi weights---------> check"); } + // Loading pT-dependent efficiency histograms from CCDB + if (cfgLoadPtEffWeights) { + + ccdb->setURL(ccdbUrl.value); + ccdb->setCaching(true); + ccdb->setLocalObjectValidityChecking(); + ccdb->setCreatedNotAfter(ccdbNoLaterThanPtEff.value); + LOGF(info, "Getting object %s", ccdbPathPtEff.value.data()); + TList* lst = ccdb->getForTimeStamp(ccdbPathPtEff.value, ccdbNoLaterThanPtEff.value); + hEffAllCharged = reinterpret_cast(lst->FindObject("hEffAllCharged")); + hEffPion = reinterpret_cast(lst->FindObject("hEffPion")); + hEffKaon = reinterpret_cast(lst->FindObject("hEffKaon")); + hEffProton = reinterpret_cast(lst->FindObject("hEffProton")); + if (!hEffAllCharged || !hEffPion || !hEffKaon || !hEffProton) + LOGF(info, "FATAL!! could not get efficiency files---------> !!! check !!!"); + } // Define axes std::vector ptBin = {0.2, 0.4, 0.6, 0.8, 1.0, 1.2, 1.4, 1.6, 1.8, 2.0, 2.2, 2.4, 2.6, 2.8, 3.0, 3.5, 4.0, 5.0, 6.0, 8.0, 10.0}; @@ -458,9 +499,17 @@ struct V0ptHadPiKaProt { flag2 += 1; if (combNSigmaKa < cfgnSigmaOtherParticles) flag2 += 1; - if (!(flag2 > 1) && !(combNSigmaPr > combNSigmaPi) && !(combNSigmaPr > combNSigmaKa)) { - if (combNSigmaPr < cfgnSigmaCutCombTPCTOF) { - flag = 1; + + if (cfgUseNewSeperationPid) { + if (std::abs(candidate.tpcNSigmaPr()) < cfgnSigmaCutTPCHigherPt && std::abs(candidate.tofNSigmaPr()) < cfgnSigmaCutTOFHigherPt) { + if (!(flag2 > 1) && std::abs(candidate.tpcNSigmaPi()) > cfgnSigmaSeperationCut && std::abs(candidate.tofNSigmaPi()) > cfgnSigmaSeperationCut && std::abs(candidate.tpcNSigmaKa()) > cfgnSigmaSeperationCut && std::abs(candidate.tofNSigmaKa()) > cfgnSigmaSeperationCut) + flag = 1; + } + } else { + if (!(flag2 > 1) && !(combNSigmaPr > combNSigmaPi) && !(combNSigmaPr > combNSigmaKa)) { + if (combNSigmaPr < cfgnSigmaCutCombTPCTOF) { + flag = 1; + } } } } @@ -497,9 +546,17 @@ struct V0ptHadPiKaProt { flag2 += 1; if (combNSigmaKa < cfgnSigmaOtherParticles) flag2 += 1; - if (!(flag2 > 1) && !(combNSigmaPi > combNSigmaPr) && !(combNSigmaPi > combNSigmaKa)) { - if (combNSigmaPi < cfgnSigmaCutCombTPCTOF) { - flag = 1; + + if (cfgUseNewSeperationPid) { + if (std::abs(candidate.tpcNSigmaPi()) < cfgnSigmaCutTPCHigherPt && std::abs(candidate.tofNSigmaPi()) < cfgnSigmaCutTOFHigherPt) { + if (!(flag2 > 1) && std::abs(candidate.tpcNSigmaKa()) > cfgnSigmaSeperationCut && std::abs(candidate.tofNSigmaKa()) > cfgnSigmaSeperationCut && std::abs(candidate.tpcNSigmaPr()) > cfgnSigmaSeperationCut && std::abs(candidate.tofNSigmaPr()) > cfgnSigmaSeperationCut) + flag = 1; + } + } else { + if (!(flag2 > 1) && !(combNSigmaPi > combNSigmaPr) && !(combNSigmaPi > combNSigmaKa)) { + if (combNSigmaPi < cfgnSigmaCutCombTPCTOF) { + flag = 1; + } } } } @@ -536,9 +593,17 @@ struct V0ptHadPiKaProt { flag2 += 1; if (combNSigmaKa < cfgnSigmaOtherParticles) flag2 += 1; - if (!(flag2 > 1) && !(combNSigmaKa > combNSigmaPi) && !(combNSigmaKa > combNSigmaPr)) { - if (combNSigmaKa < cfgnSigmaCutCombTPCTOF) { - flag = 1; + + if (cfgUseNewSeperationPid) { + if (std::abs(candidate.tpcNSigmaKa()) < cfgnSigmaCutTPCHigherPt && std::abs(candidate.tofNSigmaKa()) < cfgnSigmaCutTOFHigherPt) { + if (!(flag2 > 1) && std::abs(candidate.tpcNSigmaPi()) > cfgnSigmaSeperationCut && std::abs(candidate.tofNSigmaPi()) > cfgnSigmaSeperationCut && std::abs(candidate.tpcNSigmaPr()) > cfgnSigmaSeperationCut && std::abs(candidate.tofNSigmaPr()) > cfgnSigmaSeperationCut) + flag = 1; + } + } else { + if (!(flag2 > 1) && !(combNSigmaKa > combNSigmaPi) && !(combNSigmaKa > combNSigmaPr)) { + if (combNSigmaKa < cfgnSigmaCutCombTPCTOF) { + flag = 1; + } } } } @@ -675,6 +740,18 @@ struct V0ptHadPiKaProt { } histos.fill(HIST("hEventStatData"), 6.5); + // events with selection bits based on occupancy time pattern + if (cfgEvSelUseOcuppancyTimeCut && !(coll.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard))) { + return 0; + } + + histos.fill(HIST("hEventStatData"), 7.5); + int occupancy = coll.trackOccupancyInTimeRange(); + if (cfgEvSelSetOcuppancyRange && (occupancy < cfgMinOccupancy || occupancy > cfgMaxOccupancy)) { + return 0; + } + + histos.fill(HIST("hEventStatData"), 8.5); return 1; } @@ -714,6 +791,66 @@ struct V0ptHadPiKaProt { return weight; } + template + float getEffAllCharged(const T& candidate) + { + if (!cfgLoadPtEffWeights || !hEffAllCharged) { + return 1.0; + } + int bin = hEffAllCharged->FindBin(candidate.pt()); + float eff = hEffAllCharged->GetBinContent(bin); + float ptweight = 1.0 / eff; + if (!std::isfinite(ptweight) || ptweight <= 0) { + return 1.0; + } + return ptweight; + } + + template + float getEffPion(const T& candidate) + { + if (!cfgLoadPtEffWeights || !hEffPion) { + return 1.0; + } + int bin = hEffPion->FindBin(candidate.pt()); + float eff = hEffPion->GetBinContent(bin); + float ptweight = 1.0 / eff; + if (!std::isfinite(ptweight) || ptweight <= 0) { + return 1.0; + } + return ptweight; + } + + template + float getEffKaon(const T& candidate) + { + if (!cfgLoadPtEffWeights || !hEffKaon) { + return 1.0; + } + int bin = hEffKaon->FindBin(candidate.pt()); + float eff = hEffKaon->GetBinContent(bin); + float ptweight = 1.0 / eff; + if (!std::isfinite(ptweight) || ptweight <= 0) { + return 1.0; + } + return ptweight; + } + + template + float getEffProton(const T& candidate) + { + if (!cfgLoadPtEffWeights || !hEffProton) { + return 1.0; + } + int bin = hEffProton->FindBin(candidate.pt()); + float eff = hEffProton->GetBinContent(bin); + float ptweight = 1.0 / eff; + if (!std::isfinite(ptweight) || ptweight <= 0) { + return 1.0; + } + return ptweight; + } + // process Data void process(AodCollisions::iterator const& coll, aod::BCsWithTimestamps const&, AodTracks const& inputTracks) { @@ -820,30 +957,37 @@ struct V0ptHadPiKaProt { } } - // fill subevent B for f(pT) in v02(pT) - if (track.sign() != 0 && trkPt < cfgCutPtMaxForV02) { - if (std::abs(trkEta) < cfgCutEtaWindowB) { - fPtProfileHadInWinB->Fill(trkPt); - nSumInWinB += 1.0; - } - } double phiweight = 1.0; if (cfgLoadPhiWeights) { - phiweight = getPhiWeight(track, coll.posZ()); + phiweight = getPhiWeight(track, coll.posZ()); // NUA weight + } + double effweight = 1.0; + if (cfgLoadPtEffWeights) { + effweight = 1.0 / getEffAllCharged(track); // NUE weight } + double weight = phiweight * effweight; + // fill subevent C for v2^2 in v02(pT) if (track.sign() != 0 && trkPt < cfgCutPtMaxForV02) { histos.fill(HIST("h3DVtxZetaPhi"), coll.posZ(), trkEta, trkPhi); if (cfgCutEtaWindowB < trkEta && trkEta < 0.8) { - vecQInWinC += phiweight * TComplex(TMath::Cos(2. * trkPhi), TMath::Sin(2. * trkPhi)); - nSumInWinC += phiweight; + vecQInWinC += weight * TComplex(TMath::Cos(2. * trkPhi), TMath::Sin(2. * trkPhi)); + nSumInWinC += weight; } } // fill subevent A for v2^2 in v02(pT) if (track.sign() != 0 && trkPt < cfgCutPtMaxForV02) { if (-0.8 < trkEta && trkEta < -1.0 * cfgCutEtaWindowB) { - vecQInWinA += phiweight * TComplex(TMath::Cos(2. * trkPhi), TMath::Sin(2. * trkPhi)); - nSumInWinA += phiweight; + vecQInWinA += weight * TComplex(TMath::Cos(2. * trkPhi), TMath::Sin(2. * trkPhi)); + nSumInWinA += weight; + } + } + + // fill subevent B for f(pT) in v02(pT) + if (track.sign() != 0 && trkPt < cfgCutPtMaxForV02) { + if (std::abs(trkEta) < cfgCutEtaWindowB) { + fPtProfileHadInWinB->Fill(trkPt, effweight); + nSumInWinB += effweight; } } @@ -919,17 +1063,25 @@ struct V0ptHadPiKaProt { } } + double effweightPion = 1.0; + double effweightKaon = 1.0; + double effweightProton = 1.0; + if (cfgLoadPtEffWeights) { + effweightPion = 1.0 / getEffPion(track); // NUE weight for pion + effweightKaon = 1.0 / getEffKaon(track); // NUE weight for kaon + effweightProton = 1.0 / getEffProton(track); // NUE weight for proton + } // fill subevent B for ***identified particles'*** f(pT) in v02(pT) if (track.sign() != 0 && trkPt < cfgCutPtMaxForV02) { if (std::abs(trkEta) < cfgCutEtaWindowB) { if (isPion) { - fPtProfilePiInWinB->Fill(trkPt); + fPtProfilePiInWinB->Fill(trkPt, effweightPion); } if (isKaon) { - fPtProfileKaInWinB->Fill(trkPt); + fPtProfileKaInWinB->Fill(trkPt, effweightKaon); } if (isProton && trkPt > cfgCutPtLowerProt) { - fPtProfileProtInWinB->Fill(trkPt); + fPtProfileProtInWinB->Fill(trkPt, effweightProton); } } }