Setting up Skimming Tasks - spectraTPC

Analysis/Tasks/CMakeLists.txt

@@ -15,6 +15,7 @@ add_subdirectory(PWGJE)
 add_subdirectory(PWGLF)
 add_subdirectory(PWGUD)
 add_subdirectory(ALICE3)
+add_subdirectory(SkimmingTutorials)
 
 
 o2_add_dpl_workflow(trackextension

Analysis/Tasks/SkimmingTutorials/CMakeLists.txt

@@ -0,0 +1,25 @@
+# Copyright CERN and copyright holders of ALICE O2. This software is distributed
+# under the terms of the GNU General Public License v3 (GPL Version 3), copied
+# verbatim in the file "COPYING".
+#
+# See http://alice-o2.web.cern.ch/license for full licensing information.
+#
+# 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.
+
+
+o2_add_dpl_workflow(tpcspectra-task-skim-reference
+                    SOURCES spectraTPCReference.cxx
+                    PUBLIC_LINK_LIBRARIES O2::Framework O2::DetectorsBase O2::AnalysisDataModel
+                    COMPONENT_NAME Analysis)
+
+o2_add_dpl_workflow(tpcspectra-task-skim-provider
+                    SOURCES spectraTPCProvider.cxx
+                    PUBLIC_LINK_LIBRARIES O2::Framework O2::DetectorsBase O2::AnalysisDataModel
+                    COMPONENT_NAME Analysis)
+
+o2_add_dpl_workflow(tpcspectra-task-skim-analyser
+                    SOURCES spectraTPCAnalyser.cxx
+                    PUBLIC_LINK_LIBRARIES O2::Framework O2::DetectorsBase O2::AnalysisDataModel
+                    COMPONENT_NAME Analysis)

Analysis/Tasks/SkimmingTutorials/DataModel/LFDerived.h

@@ -0,0 +1,37 @@
+// Copyright CERN and copyright holders of ALICE O2. This software is
+// distributed under the terms of the GNU General Public License v3 (GPL
+// Version 3), copied verbatim in the file "COPYING".
+//
+// See http://alice-o2.web.cern.ch/license for full licensing information.
+//
+// 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.
+#ifndef O2_ANALYSIS_LFDERIVED_H
+#define O2_ANALYSIS_LFDERIVED_H
+
+#include "Framework/ASoA.h"
+#include "Framework/AnalysisDataModel.h"
+
+namespace o2::aod
+{
+//DECLARE_SOA_TABLE(LFCollisions, "AOD", "LFCOLLISION", o2::soa::Index<>,
+//               o2::aod::bc::RunNumber, o2::aod::collision::PosZ);
+//using LFCollision = LFCollisions::iterator;
+
+namespace lftrack
+{
+//DECLARE_SOA_INDEX_COLUMN(LFCollision, lfCollision);
+DECLARE_SOA_COLUMN(Pt, pt, float);
+DECLARE_SOA_COLUMN(P, p, float);
+DECLARE_SOA_COLUMN(Eta, eta, float);
+DECLARE_SOA_COLUMN(TPCNSigma, tpcNSigma, float[9]);
+} // namespace lftrack
+DECLARE_SOA_TABLE(LFTracks, "AOD", "LFTRACK", o2::soa::Index<>,
+                  //lftrack::LFCollisionId,
+                  lftrack::Pt, lftrack::P, lftrack::Eta,
+                  lftrack::TPCNSigma);
+using LFTrack = LFTracks::iterator;
+} // namespace o2::aod
+
+#endif // O2_ANALYSIS_LFDERIVED_H

Analysis/Tasks/SkimmingTutorials/LFDerived.h

@@ -0,0 +1,38 @@
+// Copyright CERN and copyright holders of ALICE O2. This software is
+// distributed under the terms of the GNU General Public License v3 (GPL
+// Version 3), copied verbatim in the file "COPYING".
+//
+// See http://alice-o2.web.cern.ch/license for full licensing information.
+//
+// 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.
+/// \author Nima Zardoshti <nima.zardoshti@cern.ch>, CERN
+#ifndef O2_ANALYSIS_LFDERIVED_H
+#define O2_ANALYSIS_LFDERIVED_H
+
+#include "Framework/ASoA.h"
+#include "Framework/AnalysisDataModel.h"
+
+namespace o2::aod
+{
+//DECLARE_SOA_TABLE(LFCollisions, "AOD", "LFCOLLISION", o2::soa::Index<>,
+//               o2::aod::bc::RunNumber, o2::aod::collision::PosZ);
+//using LFCollision = LFCollisions::iterator;
+
+namespace lftrack
+{
+//DECLARE_SOA_INDEX_COLUMN(LFCollision, lfCollision);
+DECLARE_SOA_COLUMN(Pt, pt, float);
+DECLARE_SOA_COLUMN(P, p, float);
+DECLARE_SOA_COLUMN(Eta, eta, float);
+DECLARE_SOA_COLUMN(TPCNSigma, tpcNSigma, float[9]);
+} // namespace lftrack
+DECLARE_SOA_TABLE(LFTracks, "AOD", "LFTRACK", o2::soa::Index<>,
+                  //lftrack::LFCollisionId,
+                  lftrack::Pt, lftrack::P, lftrack::Eta,
+                  lftrack::TPCNSigma);
+using LFTrack = LFTracks::iterator;
+} // namespace o2::aod
+
+#endif // O2_ANALYSIS_LFDERIVED_H

Analysis/Tasks/SkimmingTutorials/spectraTPCAnalyser.cxx

@@ -0,0 +1,110 @@
+// Copyright CERN and copyright holders of ALICE O2. This software is
+// distributed under the terms of the GNU General Public License v3 (GPL
+// Version 3), copied verbatim in the file "COPYING".
+//
+// See http://alice-o2.web.cern.ch/license for full licensing information.
+//
+// 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.
+
+// O2 includes
+#include "ReconstructionDataFormats/Track.h"
+#include "Framework/AnalysisTask.h"
+#include "Framework/AnalysisDataModel.h"
+#include "Framework/ASoAHelpers.h"
+#include "AnalysisDataModel/PID/PIDResponse.h"
+#include "AnalysisDataModel/TrackSelectionTables.h"
+#include "DataModel/LFDerived.h"
+
+using namespace o2;
+using namespace o2::framework;
+using namespace o2::framework::expressions;
+
+void customize(std::vector<o2::framework::ConfigParamSpec>& workflowOptions)
+{
+  std::vector<ConfigParamSpec> options{
+    {"add-tof-histos", VariantType::Int, 0, {"Generate TPC with TOF histograms"}}};
+  std::swap(workflowOptions, options);
+}
+
+#include "Framework/runDataProcessing.h"
+
+// FIXME: we should put this function in some common header so it has to be defined only once
+template <typename T>
+void makelogaxis(T h)
+{
+  const int nbins = h->GetNbinsX();
+  double binp[nbins + 1];
+  double max = h->GetXaxis()->GetBinUpEdge(nbins);
+  double min = h->GetXaxis()->GetBinLowEdge(1);
+  if (min <= 0) {
+    min = 0.00001;
+  }
+  double lmin = TMath::Log10(min);
+  double ldelta = (TMath::Log10(max) - lmin) / ((double)nbins);
+  for (int i = 0; i < nbins; i++) {
+    binp[i] = TMath::Exp(TMath::Log(10) * (lmin + i * ldelta));
+  }
+  binp[nbins] = max + 1;
+  h->GetXaxis()->Set(nbins, binp);
+}
+
+constexpr int Np = 9;
+struct TPCSpectraAnalyserTask {
+
+  static constexpr const char* pT[Np] = {"e", "#mu", "#pi", "K", "p", "d", "t", "^{3}He", "#alpha"};
+  static constexpr std::string_view hp[Np] = {"p/El", "p/Mu", "p/Pi", "p/Ka", "p/Pr", "p/De", "p/Tr", "p/He", "p/Al"};
+  static constexpr std::string_view hpt[Np] = {"pt/El", "pt/Mu", "pt/Pi", "pt/Ka", "pt/Pr", "pt/De", "pt/Tr", "pt/He", "pt/Al"};
+  HistogramRegistry histos{"Histos", {}, OutputObjHandlingPolicy::AnalysisObject};
+
+  void init(o2::framework::InitContext&)
+  {
+    histos.add("p/Unselected", "Unselected;#it{p} (GeV/#it{c})", kTH1F, {{100, 0, 20}});
+    histos.add("pt/Unselected", "Unselected;#it{p}_{T} (GeV/#it{c})", kTH1F, {{100, 0, 20}});
+    for (int i = 0; i < Np; i++) {
+      histos.add(hp[i].data(), Form("%s;#it{p} (GeV/#it{c})", pT[i]), kTH1F, {{100, 0, 20}});
+      histos.add(hpt[i].data(), Form("%s;#it{p}_{T} (GeV/#it{c})", pT[i]), kTH1F, {{100, 0, 20}});
+    }
+  }
+
+  Configurable<float> nsigmacut{"nsigmacut", 3, "Value of the Nsigma cut"};
+  Configurable<float> cfgCutVertex{"cfgCutVertex", 10.0f, "Accepted z-vertex range"};
+  Configurable<float> cfgCutEta{"cfgCutEta", 0.8f, "Eta range for tracks"};
+
+  template <std::size_t i, typename T>
+  void fillParticleHistos(const T& track, const float nsigma[])
+  {
+    if (abs(nsigma[i]) > nsigmacut.value) {
+      return;
+    }
+    histos.fill(HIST(hp[i]), track.p());
+    histos.fill(HIST(hpt[i]), track.pt());
+  }
+
+  Filter collisionFilter = nabs(aod::collision::posZ) < cfgCutVertex; //collision filters not doing anything now?
+  Filter trackFilter = nabs(aod::lftrack::eta) < cfgCutEta;
+
+  void process(soa::Filtered<aod::LFTracks>::iterator const& track)
+  {
+    auto nsigma = track.tpcNSigma();
+    histos.fill(HIST("p/Unselected"), track.p());
+    histos.fill(HIST("pt/Unselected"), track.pt());
+
+    fillParticleHistos<0>(track, nsigma);
+    fillParticleHistos<1>(track, nsigma);
+    fillParticleHistos<2>(track, nsigma);
+    fillParticleHistos<3>(track, nsigma);
+    fillParticleHistos<4>(track, nsigma);
+    fillParticleHistos<5>(track, nsigma);
+    fillParticleHistos<6>(track, nsigma);
+    fillParticleHistos<7>(track, nsigma);
+    fillParticleHistos<8>(track, nsigma);
+  }
+};
+
+WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
+{
+  WorkflowSpec workflow{adaptAnalysisTask<TPCSpectraAnalyserTask>("tpcspectra-task-skim-analyser")};
+  return workflow;
+}

Analysis/Tasks/SkimmingTutorials/spectraTPCProvider.cxx

@@ -0,0 +1,61 @@
+// Copyright CERN and copyright holders of ALICE O2. This software is
+// distributed under the terms of the GNU General Public License v3 (GPL
+// Version 3), copied verbatim in the file "COPYING".
+//
+// See http://alice-o2.web.cern.ch/license for full licensing information.
+//
+// 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.
+/// \author Nima Zardoshti <nima.zardoshti@cern.ch>, CERN
+
+// O2 includes
+#include "ReconstructionDataFormats/Track.h"
+#include "Framework/AnalysisTask.h"
+#include "Framework/AnalysisDataModel.h"
+#include "Framework/ASoAHelpers.h"
+#include "AnalysisDataModel/PID/PIDResponse.h"
+#include "AnalysisDataModel/TrackSelectionTables.h"
+#include "MathUtils/Utils.h"
+#include "DataModel/LFDerived.h"
+
+using namespace o2;
+using namespace o2::framework;
+using namespace o2::framework::expressions;
+using namespace o2::math_utils::detail;
+
+#include "Framework/runDataProcessing.h"
+
+struct TPCSpectraProviderTask {
+
+  //Produces<aod::LFCollisions> outputCollisions; //currently it seems in the spectraTPC task no loop over the collision is made. Leave this here in case it will be added
+  Produces<aod::LFTracks> outputTracks;
+
+  //Configurable<float> cfgCutVertex{"cfgCutVertex", 10.0f, "Accepted z-vertex range"};
+  Configurable<float> trackEtaCut{"trackEtaCut", 0.9f, "Eta range for tracks"};
+
+  //Filter collisionFilter = nabs(aod::collision::posZ) < cfgCutVertex;
+  Filter trackFilter = (nabs(aod::track::eta) < trackEtaCut) && (aod::track::isGlobalTrack == (uint8_t) true);
+
+  Configurable<float> nsigmacut{"nsigmacut", 3, "Value of the Nsigma cut"}; //can we add an upper limit?
+
+  using TrackCandidates = soa::Filtered<soa::Join<aod::Tracks, aod::TracksExtra, aod::pidRespTPC, aod::TrackSelection>>;
+  void process(TrackCandidates::iterator const& track)
+  {
+    uint32_t pNsigma = 0xFFFFFF00; //15 bit precision for Nsigma
+    float nsigma[9] = {truncateFloatFraction(track.tpcNSigmaEl(), pNsigma), truncateFloatFraction(track.tpcNSigmaMu(), pNsigma),
+                       truncateFloatFraction(track.tpcNSigmaPi(), pNsigma), truncateFloatFraction(track.tpcNSigmaKa(), pNsigma),
+                       truncateFloatFraction(track.tpcNSigmaPr(), pNsigma), truncateFloatFraction(track.tpcNSigmaDe(), pNsigma),
+                       truncateFloatFraction(track.tpcNSigmaTr(), pNsigma), truncateFloatFraction(track.tpcNSigmaHe(), pNsigma),
+                       truncateFloatFraction(track.tpcNSigmaAl(), pNsigma)}; //the significance needs to be discussed
+
+    //outputTracks(outputCollisions.lastIndex(), track.pt(), track.p(), track.eta(), nsigma);
+    outputTracks(track.pt(), track.p(), track.eta(), nsigma);
+  }
+};
+
+WorkflowSpec defineDataProcessing(ConfigContext const&)
+{
+  WorkflowSpec workflow{adaptAnalysisTask<TPCSpectraProviderTask>("tpcspectra-task-skim-provider")};
+  return workflow;
+}

Analysis/Tasks/SkimmingTutorials/spectraTPCReference.cxx

@@ -0,0 +1,112 @@
+// Copyright CERN and copyright holders of ALICE O2. This software is
+// distributed under the terms of the GNU General Public License v3 (GPL
+// Version 3), copied verbatim in the file "COPYING".
+//
+// See http://alice-o2.web.cern.ch/license for full licensing information.
+//
+// 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.
+
+// O2 includes
+#include "ReconstructionDataFormats/Track.h"
+#include "Framework/AnalysisTask.h"
+#include "Framework/AnalysisDataModel.h"
+#include "Framework/ASoAHelpers.h"
+#include "AnalysisDataModel/PID/PIDResponse.h"
+#include "AnalysisDataModel/TrackSelectionTables.h"
+
+using namespace o2;
+using namespace o2::framework;
+using namespace o2::framework::expressions;
+
+void customize(std::vector<o2::framework::ConfigParamSpec>& workflowOptions)
+{
+  std::vector<ConfigParamSpec> options{
+    {"add-tof-histos", VariantType::Int, 0, {"Generate TPC with TOF histograms"}}};
+  std::swap(workflowOptions, options);
+}
+
+#include "Framework/runDataProcessing.h"
+
+// FIXME: we should put this function in some common header so it has to be defined only once
+template <typename T>
+void makelogaxis(T h)
+{
+  const int nbins = h->GetNbinsX();
+  double binp[nbins + 1];
+  double max = h->GetXaxis()->GetBinUpEdge(nbins);
+  double min = h->GetXaxis()->GetBinLowEdge(1);
+  if (min <= 0) {
+    min = 0.00001;
+  }
+  double lmin = TMath::Log10(min);
+  double ldelta = (TMath::Log10(max) - lmin) / ((double)nbins);
+  for (int i = 0; i < nbins; i++) {
+    binp[i] = TMath::Exp(TMath::Log(10) * (lmin + i * ldelta));
+  }
+  binp[nbins] = max + 1;
+  h->GetXaxis()->Set(nbins, binp);
+}
+
+constexpr int Np = 9;
+struct TPCSpectraReferenceTask {
+  static constexpr const char* pT[Np] = {"e", "#mu", "#pi", "K", "p", "d", "t", "^{3}He", "#alpha"};
+  static constexpr std::string_view hp[Np] = {"p/El", "p/Mu", "p/Pi", "p/Ka", "p/Pr", "p/De", "p/Tr", "p/He", "p/Al"};
+  static constexpr std::string_view hpt[Np] = {"pt/El", "pt/Mu", "pt/Pi", "pt/Ka", "pt/Pr", "pt/De", "pt/Tr", "pt/He", "pt/Al"};
+  HistogramRegistry histos{"Histos", {}, OutputObjHandlingPolicy::AnalysisObject};
+
+  void init(o2::framework::InitContext&)
+  {
+    histos.add("p/Unselected", "Unselected;#it{p} (GeV/#it{c})", kTH1F, {{100, 0, 20}});
+    histos.add("pt/Unselected", "Unselected;#it{p}_{T} (GeV/#it{c})", kTH1F, {{100, 0, 20}});
+    for (int i = 0; i < Np; i++) {
+      histos.add(hp[i].data(), Form("%s;#it{p} (GeV/#it{c})", pT[i]), kTH1F, {{100, 0, 20}});
+      histos.add(hpt[i].data(), Form("%s;#it{p}_{T} (GeV/#it{c})", pT[i]), kTH1F, {{100, 0, 20}});
+    }
+  }
+
+  //Defining filters and input
+  Configurable<float> cfgCutVertex{"cfgCutVertex", 10.0f, "Accepted z-vertex range"};
+  Configurable<float> cfgCutEta{"cfgCutEta", 0.8f, "Eta range for tracks"};
+  Filter collisionFilter = nabs(aod::collision::posZ) < cfgCutVertex;
+  Filter trackFilter = (nabs(aod::track::eta) < cfgCutEta) && (aod::track::isGlobalTrack == (uint8_t) true);
+
+  Configurable<float> nsigmacut{"nsigmacut", 3, "Value of the Nsigma cut"};
+
+  template <std::size_t i, typename T>
+  void fillParticleHistos(const T& track, const float nsigma[])
+  {
+    if (abs(nsigma[i]) > nsigmacut.value) {
+      return;
+    }
+    histos.fill(HIST(hp[i]), track.p());
+    histos.fill(HIST(hpt[i]), track.pt());
+  }
+
+  using TrackCandidates = soa::Filtered<soa::Join<aod::Tracks, aod::TracksExtra, aod::pidRespTPC, aod::TrackSelection>>;
+  void process(TrackCandidates::iterator const& track)
+  {
+    const float nsigma[Np] = {track.tpcNSigmaEl(), track.tpcNSigmaMu(), track.tpcNSigmaPi(),
+                              track.tpcNSigmaKa(), track.tpcNSigmaPr(), track.tpcNSigmaDe(),
+                              track.tpcNSigmaTr(), track.tpcNSigmaHe(), track.tpcNSigmaAl()};
+    histos.fill(HIST("p/Unselected"), track.p());
+    histos.fill(HIST("pt/Unselected"), track.pt());
+
+    fillParticleHistos<0>(track, nsigma);
+    fillParticleHistos<1>(track, nsigma);
+    fillParticleHistos<2>(track, nsigma);
+    fillParticleHistos<3>(track, nsigma);
+    fillParticleHistos<4>(track, nsigma);
+    fillParticleHistos<5>(track, nsigma);
+    fillParticleHistos<6>(track, nsigma);
+    fillParticleHistos<7>(track, nsigma);
+    fillParticleHistos<8>(track, nsigma);
+  }
+};
+
+WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
+{
+  WorkflowSpec workflow{adaptAnalysisTask<TPCSpectraReferenceTask>("tpcspectra-task-skim-reference")};
+  return workflow;
+}