|
| 1 | +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. |
| 2 | +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. |
| 3 | +// All rights not expressly granted are reserved. |
| 4 | +// |
| 5 | +// This software is distributed under the terms of the GNU General Public |
| 6 | +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". |
| 7 | +// |
| 8 | +// In applying this license CERN does not waive the privileges and immunities |
| 9 | +// granted to it by virtue of its status as an Intergovernmental Organization |
| 10 | +// or submit itself to any jurisdiction. |
| 11 | +// |
| 12 | +/// \brief A filter task for non prompt deuterons |
| 13 | +/// \author Marta Razza marta.razza@cern.ch |
| 14 | +/// \author Francesca Ercolessi francesca.ercolessi@cern.ch |
| 15 | +/// \since Dec 17, 2025 |
| 16 | + |
| 17 | +#include "../filterTables.h" |
| 18 | + |
| 19 | +#include "PWGLF/DataModel/LFParticleIdentification.h" |
| 20 | + |
| 21 | +#include "Common/Core/RecoDecay.h" |
| 22 | +#include "Common/Core/TrackSelection.h" |
| 23 | +#include "Common/Core/trackUtilities.h" |
| 24 | +#include "Common/DataModel/Centrality.h" |
| 25 | +#include "Common/DataModel/CollisionAssociationTables.h" |
| 26 | +#include "Common/DataModel/EventSelection.h" |
| 27 | +#include "Common/DataModel/Multiplicity.h" |
| 28 | +#include "Common/DataModel/PIDResponseTOF.h" |
| 29 | +#include "Common/DataModel/TrackSelectionTables.h" |
| 30 | + |
| 31 | +#include "CCDB/BasicCCDBManager.h" |
| 32 | +#include "CommonConstants/PhysicsConstants.h" |
| 33 | +#include "DCAFitter/DCAFitterN.h" |
| 34 | +#include "DataFormatsParameters/GRPMagField.h" |
| 35 | +#include "DataFormatsParameters/GRPObject.h" |
| 36 | +#include "DetectorsBase/Propagator.h" |
| 37 | +#include "Framework/ASoAHelpers.h" |
| 38 | +#include "Framework/AnalysisDataModel.h" |
| 39 | +#include "Framework/AnalysisTask.h" |
| 40 | +#include "Framework/runDataProcessing.h" |
| 41 | +#include "ReconstructionDataFormats/Track.h" |
| 42 | +#include "ReconstructionDataFormats/TrackParametrization.h" |
| 43 | + |
| 44 | +#include "TVector3.h" |
| 45 | + |
| 46 | +#include <cmath> |
| 47 | + |
| 48 | +struct H2fromLbFilter { |
| 49 | + |
| 50 | + // Recall the output table |
| 51 | + o2::framework::Produces<o2::aod::H2fromLbFilters> table; |
| 52 | + |
| 53 | + o2::framework::Configurable<int> nBinsPt{"nBinsPt", 100, "N bins in pT histo"}; |
| 54 | + |
| 55 | + o2::base::Propagator::MatCorrType noMatCorr = o2::base::Propagator::MatCorrType::USEMatCorrNONE; |
| 56 | + |
| 57 | + // Define a histograms and registries |
| 58 | + o2::framework::HistogramRegistry QAHistos{"QAHistos", {}, o2::framework::OutputObjHandlingPolicy::AnalysisObject, false, true}; |
| 59 | + o2::framework::OutputObj<TH1D> hProcessedEvents{TH1D("hProcessedEvents", "Event filtered;; Number of events", 6, 0., 6.)}; |
| 60 | + |
| 61 | + o2::framework::ConfigurableAxis pTAxis{"pTAxis", {200, -10.f, 10.f}, "p_{T} GeV/c"}; |
| 62 | + o2::framework::ConfigurableAxis nSigmaAxis{"nSigmaAxis", {200, -10.f, 10.f}, "p_{T} GeV/c"}; |
| 63 | + o2::framework::ConfigurableAxis DCAxyAxis{"DCAxyAxis", {1000, -0.2f, 0.2f}, "DCA xy (cm)"}; |
| 64 | + o2::framework::ConfigurableAxis DCAzAxis{"DCAzAxis", {1000, -0.2f, 0.2f}, "DCA z (cm)"}; |
| 65 | + |
| 66 | + o2::framework::Configurable<float> cutzvertex{"cutzvertex", 20.0f, "Accepted z-vertex range"}; // 20 cm |
| 67 | + o2::framework::Configurable<bool> isTVX{"isTVX", true, "isTVX event selection"}; |
| 68 | + o2::framework::Configurable<bool> isNoTimeFrameBorder{"isNoTimeFrameBorder", true, "isNoTimeFrameBorder event selection"}; |
| 69 | + o2::framework::Configurable<bool> isNoITSROFrameBorder{"isNoITSROFrameBorder", true, "isNoITSROFrameBorder event selection"}; |
| 70 | + o2::framework::Configurable<float> cfgTPCNsigma{"cfgTPCNsigma", 4.0f, "TPC n sigma for deuteron PID"}; |
| 71 | + o2::framework::Configurable<float> cfgTOFNsigma{"cfgTOFNsigma", 4.0f, "TOF n sigma for deuteron PID"}; |
| 72 | + o2::framework::Configurable<float> cfgEta{"cfgEta", 0.8f, "Track eta selection"}; |
| 73 | + o2::framework::Configurable<int> cfgTPCNclsFound{"cfgTPCNclsFound", 100, "Minimum TPC clusters found"}; |
| 74 | + o2::framework::Configurable<float> cfgTPCChi2Ncl{"cfgTPCChi2Ncl", 4.0f, "Maximum TPC chi2 per N clusters"}; |
| 75 | + o2::framework::Configurable<float> cfgITSChi2Ncl{"cfgITSChi2Ncl", 36.0f, "Maximum ITS chi2 per N clusters"}; |
| 76 | + o2::framework::Configurable<int> cfgITScls{"cfgITScls", 2, "Minimum ITS clusters"}; |
| 77 | + o2::framework::Configurable<float> cfgMaxPt{"cfgMaxPt", 5.0f, "Maximum pT cut"}; |
| 78 | + o2::framework::Configurable<float> cfgMinPt{"cfgMinPt", 1.0f, "Maximum pT cut"}; |
| 79 | + o2::framework::Configurable<float> cfgDCAcut{"cfgDCAcut", 0.003f, "DCA cut for non prompt deuteron"}; |
| 80 | + |
| 81 | + o2::framework::Service<o2::ccdb::BasicCCDBManager> ccdb; |
| 82 | + |
| 83 | + void init(o2::framework::InitContext&) |
| 84 | + { |
| 85 | + ccdb->setURL("http://alice-ccdb.cern.ch"); // Set CCDB URL to get magnetic field |
| 86 | + ccdb->setCaching(true); |
| 87 | + ccdb->setLocalObjectValidityChecking(); |
| 88 | + ccdb->setCreatedNotAfter(std::chrono::duration_cast<std::chrono::milliseconds>(std::chrono::system_clock::now().time_since_epoch()).count()); |
| 89 | + |
| 90 | + o2::framework::AxisSpec ptAxis = {100, 0.0f, 10.0f, "#it{p}_{T} (GeV/#it{c})"}; |
| 91 | + |
| 92 | + // general QA histograms |
| 93 | + QAHistos.add("hVtxZ", "Z-Vertex distribution after selection;Z (cm)", o2::framework::HistType::kTH1F, {{100, -50, 50}}); |
| 94 | + QAHistos.add("hDCAxyVsPt", "DCAxy #bar{d} vs p_{T}", {o2::framework::HistType::kTH2D, {pTAxis, DCAxyAxis}}); |
| 95 | + QAHistos.add("hDCAzVsPt", "DCAz #bar{d} vs p_{T}", {o2::framework::HistType::kTH2D, {pTAxis, DCAzAxis}}); |
| 96 | + QAHistos.add("hnSigmaTPCVsPt", "n#sigma TPC vs p_{T} for #bar{d} hypothesis; p_{T} (GeV/c); n#sigma TPC", {o2::framework::HistType::kTH2D, {pTAxis, nSigmaAxis}}); |
| 97 | + QAHistos.add("hnSigmaTOFVsPt", "n#sigma TOF vs p_{T} for #bar{d} hypothesis; p_{T} (GeV/c); n#sigma TOF", {o2::framework::HistType::kTH2D, {pTAxis, nSigmaAxis}}); |
| 98 | + QAHistos.add("ptAntiDeuteron", "ptAntiDeuteron", {o2::framework::HistType::kTH1F, {ptAxis}}); |
| 99 | + QAHistos.add("etaAntideuteron", "etaAntideuteron", {o2::framework::HistType::kTH1F, {{100, -1.0f, 1.0f, "eta #bar{d}"}}}); |
| 100 | + QAHistos.add("hDCAxyVsPt-pre_selection", "DCAxy #bar{d} vs p_{T}", {o2::framework::HistType::kTH2D, {pTAxis, DCAxyAxis}}); |
| 101 | + QAHistos.add("hDCAzVsPt-pre-selection", "DCAz #bar{d} vs p_{T}", {o2::framework::HistType::kTH2D, {pTAxis, DCAzAxis}}); |
| 102 | + |
| 103 | + // processed events |
| 104 | + hProcessedEvents->GetXaxis()->SetBinLabel(1, "Events processed"); |
| 105 | + hProcessedEvents->GetXaxis()->SetBinLabel(2, "TVX"); |
| 106 | + hProcessedEvents->GetXaxis()->SetBinLabel(3, "TF border"); |
| 107 | + hProcessedEvents->GetXaxis()->SetBinLabel(4, "ITS-RO border"); |
| 108 | + hProcessedEvents->GetXaxis()->SetBinLabel(5, "z-vertex"); |
| 109 | + hProcessedEvents->GetXaxis()->SetBinLabel(6, o2::aod::filtering::H2fromLb::columnLabel()); |
| 110 | + } |
| 111 | + |
| 112 | + // Tables |
| 113 | + using CollisionCandidates = o2::soa::Join<o2::aod::Collisions, o2::aod::EvSels>; |
| 114 | + using TrackCandidates = o2::soa::Join<o2::aod::Tracks, o2::aod::TracksCov, o2::aod::TracksExtra, o2::aod::TracksDCA, o2::aod::TrackSelection, o2::aod::pidTPCFullDe, o2::aod::pidTOFFullDe>; |
| 115 | + |
| 116 | + // Single-Track Selection |
| 117 | + template <typename T1> |
| 118 | + bool passedSingleTrackSelection(const T1& track) |
| 119 | + { |
| 120 | + // Single-Track Selections |
| 121 | + if (std::abs(track.eta()) > cfgEta) |
| 122 | + return false; |
| 123 | + if (!track.hasITS()) |
| 124 | + return false; |
| 125 | + if (!track.hasTPC()) |
| 126 | + return false; |
| 127 | + if (!track.hasTOF()) |
| 128 | + return false; |
| 129 | + if (track.tpcNClsFound() < cfgTPCNclsFound) |
| 130 | + return false; |
| 131 | + if (track.tpcChi2NCl() > cfgTPCChi2Ncl) |
| 132 | + return false; |
| 133 | + if (track.itsChi2NCl() > cfgITSChi2Ncl) |
| 134 | + return false; |
| 135 | + if (track.itsNCls() < cfgITScls) |
| 136 | + return false; |
| 137 | + if (track.pt() > cfgMaxPt) |
| 138 | + return false; |
| 139 | + if (track.pt() < cfgMinPt) |
| 140 | + return false; |
| 141 | + if (track.sign() > 0) |
| 142 | + return false; |
| 143 | + |
| 144 | + return true; |
| 145 | + } |
| 146 | + |
| 147 | + void fillTriggerTable(bool keepEvent[]) |
| 148 | + { |
| 149 | + table(keepEvent[0]); |
| 150 | + } |
| 151 | + |
| 152 | + o2::framework::Preslice<o2::aod::TrackAssoc> trackIndicesPerCollision = o2::aod::track_association::collisionId; |
| 153 | + int mCurrentRun = -1; |
| 154 | + |
| 155 | + void process(CollisionCandidates const& collisions, |
| 156 | + o2::aod::TrackAssoc const& trackIndices, |
| 157 | + TrackCandidates const& tracks, |
| 158 | + o2::aod::BCsWithTimestamps const&) |
| 159 | + { |
| 160 | + for (const auto& collision : collisions) // start loop over collisions |
| 161 | + { |
| 162 | + if (mCurrentRun != collision.bc_as<o2::aod::BCsWithTimestamps>().runNumber()) { // If the run is new then we need to initialize the propagator field |
| 163 | + o2::parameters::GRPMagField* grpo = ccdb->getForTimeStamp<o2::parameters::GRPMagField>("GLO/Config/GRPMagField", collision.bc_as<o2::aod::BCsWithTimestamps>().timestamp()); |
| 164 | + o2::base::Propagator::initFieldFromGRP(grpo); |
| 165 | + mCurrentRun = collision.bc_as<o2::aod::BCsWithTimestamps>().runNumber(); |
| 166 | + } |
| 167 | + |
| 168 | + // Is event good? keepEvent[0] = non promp deuteron |
| 169 | + bool keepEvent[1]{}; // explicitly zero-initialised |
| 170 | + hProcessedEvents->Fill(0.5); |
| 171 | + // Event selection cuts |
| 172 | + if (isTVX && !collision.selection_bit(o2::aod::evsel::kIsTriggerTVX)) { |
| 173 | + fillTriggerTable(keepEvent); |
| 174 | + continue; |
| 175 | + } |
| 176 | + hProcessedEvents->Fill(1.5); |
| 177 | + if (isNoTimeFrameBorder && !collision.selection_bit(o2::aod::evsel::kNoTimeFrameBorder)) { |
| 178 | + fillTriggerTable(keepEvent); |
| 179 | + continue; |
| 180 | + } |
| 181 | + hProcessedEvents->Fill(2.5); |
| 182 | + if (isNoITSROFrameBorder && !collision.selection_bit(o2::aod::evsel::kNoITSROFrameBorder)) { |
| 183 | + fillTriggerTable(keepEvent); |
| 184 | + continue; |
| 185 | + } |
| 186 | + hProcessedEvents->Fill(3.5); |
| 187 | + if (std::abs(collision.posZ()) > cutzvertex) { |
| 188 | + fillTriggerTable(keepEvent); |
| 189 | + continue; |
| 190 | + } |
| 191 | + hProcessedEvents->Fill(4.5); |
| 192 | + QAHistos.fill(HIST("hVtxZ"), collision.posZ()); |
| 193 | + |
| 194 | + // Loop over tracks |
| 195 | + const auto& trackIdsThisCollision = trackIndices.sliceBy(trackIndicesPerCollision, collision.globalIndex()); |
| 196 | + |
| 197 | + for (const auto& trackId : trackIdsThisCollision) { // start loop over tracks |
| 198 | + |
| 199 | + const auto& track = tracks.rawIteratorAt(trackId.trackId()); |
| 200 | + |
| 201 | + std::array<float, 2> dca{track.dcaXY(), track.dcaZ()}; |
| 202 | + std::array<float, 3> pVec = track.pVector(); |
| 203 | + |
| 204 | + if (track.collisionId() != collision.globalIndex()) { |
| 205 | + auto trackPar = getTrackParCov(track); |
| 206 | + o2::base::Propagator::Instance()->propagateToDCABxByBz({collision.posX(), collision.posY(), collision.posZ()}, trackPar, 2.f, noMatCorr, &dca); |
| 207 | + getPxPyPz(trackPar, pVec); |
| 208 | + } |
| 209 | + |
| 210 | + if (!passedSingleTrackSelection(track)) { |
| 211 | + continue; |
| 212 | + } |
| 213 | + |
| 214 | + const bool isTOFDe = std::abs(track.tofNSigmaDe()) < cfgTOFNsigma; |
| 215 | + const bool isTPCDe = std::abs(track.tpcNSigmaDe()) < cfgTPCNsigma; |
| 216 | + |
| 217 | + if (isTPCDe && isTOFDe) { |
| 218 | + |
| 219 | + QAHistos.fill(HIST("hDCAxyVsPt-pre_selection"), track.pt(), dca[0]); |
| 220 | + QAHistos.fill(HIST("hDCAzVsPt-pre-selection"), track.pt(), dca[1]); |
| 221 | + if (std::abs(dca[0]) < cfgDCAcut) { |
| 222 | + continue; |
| 223 | + } |
| 224 | + keepEvent[0] = true; |
| 225 | + QAHistos.fill(HIST("ptAntiDeuteron"), track.pt()); |
| 226 | + QAHistos.fill(HIST("etaAntideuteron"), track.eta()); |
| 227 | + QAHistos.fill(HIST("hDCAxyVsPt"), track.pt(), dca[0]); |
| 228 | + QAHistos.fill(HIST("hDCAzVsPt"), track.pt(), dca[1]); |
| 229 | + QAHistos.fill(HIST("hnSigmaTPCVsPt"), track.pt(), track.tpcNSigmaDe()); |
| 230 | + QAHistos.fill(HIST("hnSigmaTOFVsPt"), track.pt(), track.tofNSigmaDe()); |
| 231 | + } |
| 232 | + |
| 233 | + } // end track loop |
| 234 | + if (keepEvent[0]) { |
| 235 | + hProcessedEvents->Fill(5.5); |
| 236 | + } |
| 237 | + |
| 238 | + fillTriggerTable(keepEvent); |
| 239 | + } |
| 240 | + } |
| 241 | +}; |
| 242 | + |
| 243 | +o2::framework::WorkflowSpec defineDataProcessing(o2::framework::ConfigContext const& cfgc) |
| 244 | +{ |
| 245 | + return o2::framework::WorkflowSpec{o2::framework::adaptAnalysisTask<H2fromLbFilter>(cfgc)}; |
| 246 | +} |
0 commit comments