diff --git a/Analysis/ALICE3/CMakeLists.txt b/Analysis/ALICE3/CMakeLists.txt
index 41a841e9e7a80..ea55031b224c9 100644
--- a/Analysis/ALICE3/CMakeLists.txt
+++ b/Analysis/ALICE3/CMakeLists.txt
@@ -26,6 +26,16 @@ o2_add_dpl_workflow(alice3-qa-multiplicity
                     PUBLIC_LINK_LIBRARIES O2::Framework O2::AnalysisDataModel O2::AnalysisCore
                     COMPONENT_NAME Analysis)
 
+o2_add_dpl_workflow(alice3-qa-singleparticle
+                    SOURCES src/alice3-qa-singleparticle.cxx
+                    PUBLIC_LINK_LIBRARIES O2::Framework O2::AnalysisDataModel O2::AnalysisCore
+                    COMPONENT_NAME Analysis)
+
+o2_add_dpl_workflow(alice3-lutmaker
+                    SOURCES src/alice3-lutmaker.cxx
+                    PUBLIC_LINK_LIBRARIES O2::Framework O2::AnalysisDataModel O2::AnalysisCore
+                    COMPONENT_NAME Analysis)
+
 o2_add_dpl_workflow(alice3-pid-rich-qa
                     SOURCES src/pidRICHqa.cxx
                     PUBLIC_LINK_LIBRARIES O2::AnalysisDataModel O2::AnalysisCore O2::ALICE3Analysis
diff --git a/Analysis/ALICE3/src/alice3-lutmaker.cxx b/Analysis/ALICE3/src/alice3-lutmaker.cxx
new file mode 100644
index 0000000000000..b17c731e915e3
--- /dev/null
+++ b/Analysis/ALICE3/src/alice3-lutmaker.cxx
@@ -0,0 +1,176 @@
+// 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 Nicolo' Jacazio <nicolo.jacazio@cern.ch>, CERN
+/// \brief Task to extract LUTs for the fast simulation from full simulation
+/// \since 27/04/2021
+
+// O2 includes
+#include "Framework/AnalysisTask.h"
+#include "AnalysisCore/MC.h"
+#include "ReconstructionDataFormats/Track.h"
+
+using namespace o2;
+using namespace framework;
+using namespace framework::expressions;
+
+void customize(std::vector<o2::framework::ConfigParamSpec>& workflowOptions)
+{
+  std::vector<ConfigParamSpec> options{
+    {"lut-el", VariantType::Int, 1, {"LUT input for the Electron PDG code"}},
+    {"lut-mu", VariantType::Int, 1, {"LUT input for the Muon PDG code"}},
+    {"lut-pi", VariantType::Int, 1, {"LUT input for the Pion PDG code"}},
+    {"lut-ka", VariantType::Int, 1, {"LUT input for the Kaon PDG code"}},
+    {"lut-pr", VariantType::Int, 1, {"LUT input for the Proton PDG code"}},
+    {"lut-tr", VariantType::Int, 0, {"LUT input for the Triton PDG code"}},
+    {"lut-de", VariantType::Int, 0, {"LUT input for the Deuteron PDG code"}},
+    {"lut-he", VariantType::Int, 0, {"LUT input for the Helium3 PDG code"}}};
+  std::swap(workflowOptions, options);
+}
+
+#include "Framework/runDataProcessing.h"
+
+template <o2::track::pid_constants::ID particle>
+struct Alice3LutMaker {
+  static constexpr int nSpecies = 8;
+  static constexpr int PDGs[nSpecies] = {kElectron, kMuonMinus, kPiPlus, kKPlus, kProton, 1000010020, 1000010030, 1000020030};
+  static_assert(particle < nSpecies && "Maximum of particles reached");
+  static constexpr int pdg = PDGs[particle];
+  Configurable<bool> selPrim{"sel-prim", false, "If true selects primaries, if not select all particles"};
+  Configurable<int> etaBins{"eta-bins", 80, "Number of eta bins"};
+  Configurable<float> etaMin{"eta-min", -4.f, "Lower limit in eta"};
+  Configurable<float> etaMax{"eta-max", 4.f, "Upper limit in eta"};
+  Configurable<int> ptBins{"pt-bins", 200, "Number of pT bins"};
+  Configurable<float> ptMin{"pt-min", -2.f, "Lower limit in pT"};
+  Configurable<float> ptMax{"pt-max", 2.f, "Upper limit in pT"};
+  Configurable<int> logPt{"log-pt", 0, "Flag to use a logarithmic pT axis, in this case the pT limits are the expontents"};
+  HistogramRegistry histos{"Histos", {}, OutputObjHandlingPolicy::AnalysisObject};
+
+  void init(InitContext&)
+  {
+    const TString commonTitle = Form(" PDG %i", pdg);
+    AxisSpec axisPt{ptBins, ptMin, ptMax, "#it{p}_{T} GeV/#it{c}"};
+    if (logPt) {
+      const double min = axisPt.binEdges[0];
+      const double width = (axisPt.binEdges[1] - axisPt.binEdges[0]) / axisPt.nBins.value();
+      axisPt.binEdges.clear();
+      axisPt.binEdges.resize(0);
+      for (int i = 0; i < axisPt.nBins.value() + 1; i++) {
+        axisPt.binEdges.push_back(std::pow(10., min + i * width));
+      }
+      axisPt.nBins = std::nullopt;
+    }
+    const AxisSpec axisEta{etaBins, etaMin, etaMax, "#it{#eta}"};
+
+    histos.add("multiplicity", "Track multiplicity;Tracks per event;Events", kTH1F, {{100, 0, 2000}});
+    histos.add("pt", "pt" + commonTitle, kTH1F, {axisPt});
+    histos.add("eta", "eta" + commonTitle, kTH1F, {axisEta});
+    histos.add("CovMat_cYY", "cYY" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_cZY", "cZY" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_cZZ", "cZZ" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_cSnpY", "cSnpY" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_cSnpZ", "cSnpZ" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_cSnpSnp", "cSnpSnp" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_cTglY", "cTglY" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_cTglZ", "cTglZ" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_cTglSnp", "cTglSnp" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_cTglTgl", "cTglTgl" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_c1PtY", "c1PtY" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_c1PtZ", "c1PtZ" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_c1PtSnp", "c1PtSnp" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_c1PtTgl", "c1PtTgl" + commonTitle, kTProfile2D, {axisPt, axisEta});
+    histos.add("CovMat_c1Pt21Pt2", "c1Pt21Pt2" + commonTitle, kTProfile2D, {axisPt, axisEta});
+
+    histos.add("Efficiency", "Efficiency" + commonTitle, kTProfile2D, {axisPt, axisEta});
+  }
+
+  void process(const soa::Join<aod::Tracks, aod::TracksCov, aod::McTrackLabels>& tracks,
+               const aod::McParticles& mcParticles)
+  {
+    std::vector<int64_t> recoTracks(tracks.size());
+    int ntrks = 0;
+
+    for (const auto& track : tracks) {
+      const auto mcParticle = track.mcParticle();
+      if (mcParticle.pdgCode() != pdg) {
+        continue;
+      }
+      if (selPrim.value && !MC::isPhysicalPrimary(mcParticles, mcParticle)) { // Requiring is physical primary
+        continue;
+      }
+
+      recoTracks[ntrks++] = mcParticle.globalIndex();
+
+      histos.fill(HIST("pt"), mcParticle.pt());
+      histos.fill(HIST("eta"), mcParticle.eta());
+      histos.fill(HIST("CovMat_cYY"), mcParticle.pt(), mcParticle.eta(), track.cYY());
+      histos.fill(HIST("CovMat_cZY"), mcParticle.pt(), mcParticle.eta(), track.cZY());
+      histos.fill(HIST("CovMat_cZZ"), mcParticle.pt(), mcParticle.eta(), track.cZZ());
+      histos.fill(HIST("CovMat_cSnpY"), mcParticle.pt(), mcParticle.eta(), track.cSnpY());
+      histos.fill(HIST("CovMat_cSnpZ"), mcParticle.pt(), mcParticle.eta(), track.cSnpZ());
+      histos.fill(HIST("CovMat_cSnpSnp"), mcParticle.pt(), mcParticle.eta(), track.cSnpSnp());
+      histos.fill(HIST("CovMat_cTglY"), mcParticle.pt(), mcParticle.eta(), track.cTglY());
+      histos.fill(HIST("CovMat_cTglZ"), mcParticle.pt(), mcParticle.eta(), track.cTglZ());
+      histos.fill(HIST("CovMat_cTglSnp"), mcParticle.pt(), mcParticle.eta(), track.cTglSnp());
+      histos.fill(HIST("CovMat_cTglTgl"), mcParticle.pt(), mcParticle.eta(), track.cTglTgl());
+      histos.fill(HIST("CovMat_c1PtY"), mcParticle.pt(), mcParticle.eta(), track.c1PtY());
+      histos.fill(HIST("CovMat_c1PtZ"), mcParticle.pt(), mcParticle.eta(), track.c1PtZ());
+      histos.fill(HIST("CovMat_c1PtSnp"), mcParticle.pt(), mcParticle.eta(), track.c1PtSnp());
+      histos.fill(HIST("CovMat_c1PtTgl"), mcParticle.pt(), mcParticle.eta(), track.c1PtTgl());
+      histos.fill(HIST("CovMat_c1Pt21Pt2"), mcParticle.pt(), mcParticle.eta(), track.c1Pt21Pt2());
+    }
+    histos.fill(HIST("multiplicity"), ntrks);
+
+    for (const auto& mcParticle : mcParticles) {
+      if (mcParticle.pdgCode() != pdg) {
+        continue;
+      }
+      if (!MC::isPhysicalPrimary(mcParticles, mcParticle)) { // Requiring is physical primary
+        continue;
+      }
+
+      if (std::find(recoTracks.begin(), recoTracks.end(), mcParticle.globalIndex()) != recoTracks.end()) {
+        histos.fill(HIST("Efficiency"), mcParticle.pt(), mcParticle.eta(), 1.);
+      } else {
+        histos.fill(HIST("Efficiency"), mcParticle.pt(), mcParticle.eta(), 0.);
+      }
+    }
+  }
+};
+
+WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
+{
+  WorkflowSpec w;
+  if (cfgc.options().get<int>("lut-el")) {
+    w.push_back(adaptAnalysisTask<Alice3LutMaker<o2::track::PID::Electron>>(cfgc, TaskName{"alice3-lutmaker-electron"}));
+  }
+  if (cfgc.options().get<int>("lut-mu")) {
+    w.push_back(adaptAnalysisTask<Alice3LutMaker<o2::track::PID::Muon>>(cfgc, TaskName{"alice3-lutmaker-muon"}));
+  }
+  if (cfgc.options().get<int>("lut-pi")) {
+    w.push_back(adaptAnalysisTask<Alice3LutMaker<o2::track::PID::Pion>>(cfgc, TaskName{"alice3-lutmaker-pion"}));
+  }
+  if (cfgc.options().get<int>("lut-ka")) {
+    w.push_back(adaptAnalysisTask<Alice3LutMaker<o2::track::PID::Kaon>>(cfgc, TaskName{"alice3-lutmaker-kaon"}));
+  }
+  if (cfgc.options().get<int>("lut-pr")) {
+    w.push_back(adaptAnalysisTask<Alice3LutMaker<o2::track::PID::Proton>>(cfgc, TaskName{"alice3-lutmaker-proton"}));
+  }
+  if (cfgc.options().get<int>("lut-de")) {
+    w.push_back(adaptAnalysisTask<Alice3LutMaker<o2::track::PID::Deuteron>>(cfgc, TaskName{"alice3-lutmaker-deuteron"}));
+  }
+  if (cfgc.options().get<int>("lut-tr")) {
+    w.push_back(adaptAnalysisTask<Alice3LutMaker<o2::track::PID::Triton>>(cfgc, TaskName{"alice3-lutmaker-triton"}));
+  }
+  if (cfgc.options().get<int>("lut-he")) {
+    w.push_back(adaptAnalysisTask<Alice3LutMaker<o2::track::PID::Helium3>>(cfgc, TaskName{"alice3-lutmaker-helium3"}));
+  }
+  return w;
+}
diff --git a/Analysis/ALICE3/src/alice3-qa-singleparticle.cxx b/Analysis/ALICE3/src/alice3-qa-singleparticle.cxx
new file mode 100644
index 0000000000000..81a085d1347a6
--- /dev/null
+++ b/Analysis/ALICE3/src/alice3-qa-singleparticle.cxx
@@ -0,0 +1,135 @@
+// 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 Nicolo' Jacazio <nicolo.jacazio@cern.ch>, CERN
+
+// O2 includes
+#include "Framework/AnalysisTask.h"
+#include "Framework/runDataProcessing.h"
+#include "AnalysisCore/MC.h"
+#include "Framework/HistogramRegistry.h"
+
+using namespace o2;
+using namespace o2::framework;
+using namespace o2::framework::expressions;
+
+struct Alice3SingleParticle {
+  Configurable<int> PDG{"PDG", 2212, "PDG code of the particle of interest"};
+  Configurable<int> IsStable{"IsStable", 0, "Flag to check stable particles"};
+  HistogramRegistry histos{"Histos", {}, OutputObjHandlingPolicy::AnalysisObject};
+  Configurable<int> ptBins{"pt-bins", 500, "Number of pT bins"};
+  Configurable<float> ptMin{"pt-min", 0.f, "Lower limit in pT"};
+  Configurable<float> ptMax{"pt-max", 5.f, "Upper limit in pT"};
+  Configurable<int> etaBins{"eta-bins", 500, "Number of eta bins"};
+  Configurable<float> etaMin{"eta-min", -3.f, "Lower limit in eta"};
+  Configurable<float> etaMax{"eta-max", 3.f, "Upper limit in eta"};
+  Configurable<float> yMin{"y-min", -3.f, "Lower limit in y"};
+  Configurable<float> yMax{"y-max", 3.f, "Upper limit in y"};
+  Configurable<int> prodBins{"prod-bins", 100, "Number of production vertex bins"};
+  Configurable<float> prodMin{"prod-min", -1.f, "Lower limit in production vertex"};
+  Configurable<float> prodMax{"prod-max", 1.f, "Upper limit in production vertex"};
+  Configurable<float> charge{"charge", 1.f, "Particle charge to scale the reconstructed momentum"};
+  Configurable<bool> doPrint{"doPrint", false, "Flag to print debug messages"};
+
+  void init(InitContext&)
+  {
+    const TString tit = Form("%i", PDG.value);
+    AxisSpec axisPt{ptBins, ptMin, ptMax};
+    AxisSpec axisEta{etaBins, etaMin, etaMax};
+    AxisSpec axisProd{prodBins, prodMin, prodMax};
+
+    histos.add("particlePt", "Particle Pt " + tit + ";#it{p}_{T} (GeV/#it{c})", kTH1D, {axisPt});
+    histos.add("prodVx", "Particle Prod. Vertex X " + tit + ";Prod. Vertex X (cm)", kTH1D, {axisProd});
+    histos.add("prodVy", "Particle Prod. Vertex Y " + tit + ";Prod. Vertex Y (cm)", kTH1D, {axisProd});
+    histos.add("prodVz", "Particle Prod. Vertex Z " + tit + ";Prod. Vertex Z (cm)", kTH1D, {axisProd});
+    histos.add("prodRadius", "Particle Prod. Vertex Radius " + tit + ";Prod. Vertex Radius (cm)", kTH1D, {axisProd});
+    histos.add("prodVxVsPt", "Particle Prod. Vertex X " + tit + ";#it{p}_{T} (GeV/#it{c});Prod. Vertex X (cm)", kTH2D, {axisPt, axisProd});
+    histos.add("prodVyVsPt", "Particle Prod. Vertex Y " + tit + ";#it{p}_{T} (GeV/#it{c});Prod. Vertex Y (cm)", kTH2D, {axisPt, axisProd});
+    histos.add("prodVzVsPt", "Particle Prod. Vertex Z " + tit + ";#it{p}_{T} (GeV/#it{c});Prod. Vertex Z (cm)", kTH2D, {axisPt, axisProd});
+    histos.add("prodRadiusVsPt", "Particle Prod. Vertex Radius " + tit + ";#it{p}_{T} (GeV/#it{c});Prod. Vertex Radius (cm)", kTH2D, {axisPt, axisProd});
+    histos.add("prodRadius3DVsPt", "Particle Prod. Vertex Radius XYZ " + tit + ";#it{p}_{T} (GeV/#it{c});Prod. Vertex Radius XYZ (cm)", kTH2D, {axisPt, axisProd});
+    histos.add("trackPt", "Track Pt " + tit + ";#it{p}_{T} (GeV/#it{c})", kTH1D, {axisPt});
+    histos.add("particleEta", "Particle Eta " + tit + ";#it{#eta}", kTH1D, {axisEta});
+    histos.add("trackEta", "Track Eta " + tit + ";#it{#eta}", kTH1D, {axisEta});
+    histos.add("particleY", "Particle Y " + tit + ";#it{y}", kTH1D, {axisEta});
+    histos.add("primaries", "Source for primaries " + tit + ";PDG Code", kTH1D, {{100, 0.f, 100.f}});
+    histos.add("secondaries", "Source for secondaries " + tit + ";PDG Code", kTH1D, {{100, 0.f, 100.f}});
+  }
+
+  void process(const soa::Join<o2::aod::Tracks, o2::aod::McTrackLabels>& tracks,
+               const aod::McParticles& mcParticles)
+  {
+
+    std::vector<int64_t> ParticlesOfInterest;
+    for (const auto& mcParticle : mcParticles) {
+      if (mcParticle.pdgCode() != PDG) {
+        continue;
+      }
+      if (mcParticle.y() < yMin || mcParticle.y() > yMax) {
+        continue;
+      }
+      histos.fill(HIST("particlePt"), mcParticle.pt());
+      histos.fill(HIST("particleEta"), mcParticle.eta());
+      histos.fill(HIST("particleY"), mcParticle.y());
+      histos.fill(HIST("prodVx"), mcParticle.vx());
+      histos.fill(HIST("prodVy"), mcParticle.vy());
+      histos.fill(HIST("prodVz"), mcParticle.vz());
+      histos.fill(HIST("prodRadius"), std::sqrt(mcParticle.vx() * mcParticle.vx() + mcParticle.vy() * mcParticle.vy()));
+      histos.fill(HIST("prodVxVsPt"), mcParticle.pt(), mcParticle.vx());
+      histos.fill(HIST("prodVyVsPt"), mcParticle.pt(), mcParticle.vy());
+      histos.fill(HIST("prodVzVsPt"), mcParticle.pt(), mcParticle.vz());
+      histos.fill(HIST("prodRadiusVsPt"), mcParticle.pt(), std::sqrt(mcParticle.vx() * mcParticle.vx() + mcParticle.vy() * mcParticle.vy()));
+      histos.fill(HIST("prodRadius3DVsPt"), mcParticle.pt(), std::sqrt(mcParticle.vx() * mcParticle.vx() + mcParticle.vy() * mcParticle.vy() + mcParticle.vz() * mcParticle.vz()));
+      ParticlesOfInterest.push_back(mcParticle.globalIndex());
+    }
+
+    for (const auto& track : tracks) {
+      const auto mcParticle = track.mcParticle();
+      if (!IsStable) {
+        if (mcParticle.mother0() < 0) {
+          continue;
+        }
+        auto mother = mcParticles.iteratorAt(mcParticle.mother0());
+        const auto ParticleIsInteresting = std::find(ParticlesOfInterest.begin(), ParticlesOfInterest.end(), mother.globalIndex()) != ParticlesOfInterest.end();
+        if (!ParticleIsInteresting) {
+          continue;
+        }
+        if (doPrint) {
+          LOG(INFO) << "Track " << track.globalIndex() << " comes from a " << mother.pdgCode() << " and is a " << mcParticle.pdgCode();
+        }
+      } else {
+        if (mcParticle.pdgCode() != PDG) {
+          continue;
+        }
+        histos.fill(HIST("trackPt"), track.pt() * charge);
+        histos.fill(HIST("trackEta"), track.eta());
+        if (mcParticle.mother0() < 0) {
+          if (doPrint) {
+            LOG(INFO) << "Track " << track.globalIndex() << " is a " << mcParticle.pdgCode();
+          }
+          continue;
+        }
+        auto mother = mcParticles.iteratorAt(mcParticle.mother0());
+        if (MC::isPhysicalPrimary(mcParticles, mcParticle)) {
+          histos.get<TH1>(HIST("primaries"))->Fill(Form("%i", mother.pdgCode()), 1.f);
+        } else {
+          histos.get<TH1>(HIST("secondaries"))->Fill(Form("%i", mother.pdgCode()), 1.f);
+        }
+        if (doPrint) {
+          LOG(INFO) << "Track " << track.globalIndex() << " is a " << mcParticle.pdgCode() << " and comes from a " << mother.pdgCode() << " and is " << (MC::isPhysicalPrimary(mcParticles, mcParticle) ? "" : "not") << " a primary";
+        }
+      }
+    }
+  }
+};
+
+WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
+{
+  return WorkflowSpec{adaptAnalysisTask<Alice3SingleParticle>(cfgc)};
+}
diff --git a/Analysis/Tasks/PWGPP/qaEfficiency.cxx b/Analysis/Tasks/PWGPP/qaEfficiency.cxx
index 4d55dd0ca67fa..dedad9adc4545 100644
--- a/Analysis/Tasks/PWGPP/qaEfficiency.cxx
+++ b/Analysis/Tasks/PWGPP/qaEfficiency.cxx
@@ -17,7 +17,7 @@
 #include "Framework/AnalysisTask.h"
 #include "Framework/HistogramRegistry.h"
 #include "ReconstructionDataFormats/DCA.h"
-#include "AnalysisCore/trackUtilities.h"
+#include "ReconstructionDataFormats/Track.h"
 #include "AnalysisCore/MC.h"
 #include "AnalysisDataModel/TrackSelectionTables.h"
 
@@ -73,8 +73,11 @@ struct QaTrackingEfficiency {
   Configurable<float> phiMax{"phi-max", 6.284f, "Upper limit in phi"};
   Configurable<float> ptMin{"pt-min", 0.f, "Lower limit in pT"};
   Configurable<float> ptMax{"pt-max", 5.f, "Upper limit in pT"};
+  Configurable<int> selPrim{"sel-prim", 1, "1 select primaries, 0 select all particles"};
+  Configurable<int> pdgSign{"pdgSign", 0, "Sign to give to the PDG. If 0 both signs are accepted."};
+  Configurable<bool> noFakes{"noFakes", false, "Flag to reject tracks that have fake hits"};
   // Event selection
-  Configurable<int> nMinNumberOfContributors{"nMinNumberOfContributors", 2, "Minimum required number of contributors to the vertex"};
+  Configurable<int> nMinNumberOfContributors{"nMinNumberOfContributors", 2, "Minimum required number of contributors to the primary vertex"};
   Configurable<float> vertexZMin{"vertex-z-min", -10.f, "Minimum position of the generated vertez in Z (cm)"};
   Configurable<float> vertexZMax{"vertex-z-max", 10.f, "Maximum position of the generated vertez in Z (cm)"};
   // Histogram configuration
@@ -82,7 +85,6 @@ struct QaTrackingEfficiency {
   Configurable<int> logPt{"log-pt", 0, "Flag to use a logarithmic pT axis"};
   Configurable<int> etaBins{"eta-bins", 500, "Number of eta bins"};
   Configurable<int> phiBins{"phi-bins", 500, "Number of phi bins"};
-  Configurable<int> selPrim{"sel-prim", 1, "1 select primaries, 0 select all particles"};
   // Task configuration
   Configurable<int> makeEff{"make-eff", 0, "Flag to produce the efficiency with TEfficiency"};
 
@@ -91,6 +93,10 @@ struct QaTrackingEfficiency {
 
   void init(InitContext&)
   {
+    if (pdgSign != 0 && pdgSign != 1 && pdgSign != -1) {
+      LOG(FATAL) << "Provide pdgSign as 0, 1, -1. Provided: " << pdgSign.value;
+      return;
+    }
     const TString tagPt = Form("%s #it{#eta} [%.2f,%.2f] #it{#varphi} [%.2f,%.2f] Prim %i",
                                o2::track::pid_constants::sNames[particle],
                                etaMin.value, etaMax.value,
@@ -138,9 +144,21 @@ struct QaTrackingEfficiency {
           list->Add(new TEfficiency(effname, efftitle, axis->GetNbins(), axis->GetXmin(), axis->GetXmax()));
         }
       };
+      auto makeEfficiency2D = [&](TString effname, TString efftitle, auto templateHistoX, auto templateHistoY) {
+        TAxis* axisX = histos.get<TH1>(templateHistoX)->GetXaxis();
+        TAxis* axisY = histos.get<TH1>(templateHistoY)->GetXaxis();
+        if (axisX->IsVariableBinSize() || axisY->IsVariableBinSize()) {
+          list->Add(new TEfficiency(effname, efftitle, axisX->GetNbins(), axisX->GetXbins()->GetArray(), axisY->GetNbins(), axisY->GetXbins()->GetArray()));
+        } else {
+          list->Add(new TEfficiency(effname, efftitle, axisX->GetNbins(), axisX->GetXmin(), axisX->GetXmax(), axisY->GetNbins(), axisY->GetXmin(), axisY->GetXmax()));
+        }
+      };
       makeEfficiency("efficiencyVsPt", "Efficiency " + tagPt + ";" + xPt + ";Efficiency", HIST("pt/num"));
+      makeEfficiency("efficiencyVsP", "Efficiency " + tagPt + ";#it{p} (GeV/#it{c});Efficiency", HIST("pt/num"));
       makeEfficiency("efficiencyVsEta", "Efficiency " + tagEta + ";" + xEta + ";Efficiency", HIST("eta/num"));
       makeEfficiency("efficiencyVsPhi", "Efficiency " + tagPhi + ";" + xPhi + ";Efficiency", HIST("phi/num"));
+
+      makeEfficiency2D("efficiencyVsPtVsEta", Form("Efficiency %s #it{#varphi} [%.2f,%.2f] Prim %i;%s;%s;Efficiency", o2::track::pid_constants::sNames[particle], phiMin.value, phiMax.value, selPrim.value, xPt.Data(), xEta.Data()), HIST("pt/num"), HIST("eta/num"));
     }
   }
 
@@ -179,12 +197,33 @@ struct QaTrackingEfficiency {
       if ((selPrim == 1) && (!MC::isPhysicalPrimary(mcParticles, mcParticle))) { // Requiring is physical primary
         continue;
       }
-      if (abs(mcParticle.pdgCode()) == particlePDG) { // Checking PDG code
-        histos.fill(HIST("pt/num"), mcParticle.pt());
-        histos.fill(HIST("eta/num"), mcParticle.eta());
-        histos.fill(HIST("phi/num"), mcParticle.phi());
-        recoTracks[ntrks++] = mcParticle.globalIndex();
+
+      if (noFakes) { // Selecting tracks with no fake hits
+        bool hasFake = false;
+        for (int i = 0; i < 10; i++) { // From ITS to TPC
+          if (track.mcMask() & 1 << i) {
+            hasFake = true;
+            break;
+          }
+        }
+        if (hasFake) {
+          continue;
+        }
+      }
+      // Selecting PDG code
+      if (pdgSign == 0 && abs(mcParticle.pdgCode()) != particlePDG) {
+        continue;
+      } else if (pdgSign == 1 && mcParticle.pdgCode() != particlePDG) {
+        continue;
+      } else if (pdgSign == -1 && mcParticle.pdgCode() != -particlePDG) {
+        continue;
+      } else {
+        continue;
       }
+      histos.fill(HIST("pt/num"), mcParticle.pt());
+      histos.fill(HIST("eta/num"), mcParticle.eta());
+      histos.fill(HIST("phi/num"), mcParticle.phi());
+      recoTracks[ntrks++] = mcParticle.globalIndex();
     }
 
     for (const auto& mcParticle : mcParticles) {
@@ -205,8 +244,10 @@ struct QaTrackingEfficiency {
         if (makeEff) {
           const auto particleReconstructed = std::find(recoTracks.begin(), recoTracks.end(), mcParticle.globalIndex()) != recoTracks.end();
           static_cast<TEfficiency*>(list->At(0))->Fill(particleReconstructed, mcParticle.pt());
-          static_cast<TEfficiency*>(list->At(1))->Fill(particleReconstructed, mcParticle.eta());
-          static_cast<TEfficiency*>(list->At(2))->Fill(particleReconstructed, mcParticle.phi());
+          static_cast<TEfficiency*>(list->At(1))->Fill(particleReconstructed, mcParticle.p());
+          static_cast<TEfficiency*>(list->At(2))->Fill(particleReconstructed, mcParticle.eta());
+          static_cast<TEfficiency*>(list->At(3))->Fill(particleReconstructed, mcParticle.phi());
+          static_cast<TEfficiency*>(list->At(4))->Fill(particleReconstructed, mcParticle.pt(), mcParticle.eta());
         }
         histos.fill(HIST("pt/den"), mcParticle.pt());
         histos.fill(HIST("eta/den"), mcParticle.eta());
diff --git a/Framework/Core/include/Framework/AnalysisDataModel.h b/Framework/Core/include/Framework/AnalysisDataModel.h
index fb9f2c7f4c7d9..8f77174c14fad 100644
--- a/Framework/Core/include/Framework/AnalysisDataModel.h
+++ b/Framework/Core/include/Framework/AnalysisDataModel.h
@@ -791,6 +791,10 @@ DECLARE_SOA_DYNAMIC_COLUMN(Eta, eta, //! Pseudorapidity
                            [](float px, float py, float pz) -> float { return 0.5f * std::log((std::sqrt(px * px + py * py + pz * pz) + pz) / (std::sqrt(px * px + py * py + pz * pz) - pz)); });
 DECLARE_SOA_DYNAMIC_COLUMN(Pt, pt, //! Transverse momentum in GeV/c
                            [](float px, float py) -> float { return std::sqrt(px * px + py * py); });
+DECLARE_SOA_DYNAMIC_COLUMN(P, p, //! Total momentum in GeV/c
+                           [](float px, float py, float pz) -> float { return std::sqrt(px * px + py * py + pz * pz); });
+DECLARE_SOA_DYNAMIC_COLUMN(Y, y, //! Particle rapidity
+                           [](float pz, float energy) -> float { return 0.5f * std::log((energy + pz) / (energy - pz)); });
 DECLARE_SOA_DYNAMIC_COLUMN(ProducedByGenerator, producedByGenerator, //! Particle produced by the generator or by the transport code
                            [](uint8_t flags) -> bool { return (flags & 0x1) == 0x0; });
 } // namespace mcparticle
@@ -805,6 +809,8 @@ DECLARE_SOA_TABLE(McParticles, "AOD", "MCPARTICLE", //! MC particle table
                   mcparticle::Phi<mcparticle::Px, mcparticle::Py>,
                   mcparticle::Eta<mcparticle::Px, mcparticle::Py, mcparticle::Pz>,
                   mcparticle::Pt<mcparticle::Px, mcparticle::Py>,
+                  mcparticle::P<mcparticle::Px, mcparticle::Py, mcparticle::Pz>,
+                  mcparticle::Y<mcparticle::Pz, mcparticle::E>,
                   mcparticle::ProducedByGenerator<mcparticle::Flags>);
 using McParticle = McParticles::iterator;