Updated ppL trigger, added kaon rejection

EventFiltering/PWGCF/CFFilter.cxx

@@ -168,6 +168,8 @@ struct CFFilter {
 
       registry.add("fPtPPL", "Transverse momentum of all processed tracks", HistType::kTH1F, {{1000, 0, 10}});
       registry.add("fPtAntiPPL", "Transverse momentum of all processed antitracks", HistType::kTH1F, {{1000, 0, 10}});
+      registry.add("fMinvLambda", "Invariant mass of lambdas ", HistType::kTH1F, {{1000, 0.7, 1.5}});
+      registry.add("fMinvAntiLambda", "Invariant mass of antilambdas ", HistType::kTH1F, {{1000, 0.7, 1.5}});
     }
 
     /// Initializing CCDB
@@ -205,6 +207,7 @@ struct CFFilter {
     auto partsLambda0 = partsLambda0Part->sliceByCached(aod::femtodreamparticle::femtoDreamCollisionId, col.globalIndex());
     auto partsProton1 = partsProton1Part->sliceByCached(aod::femtodreamparticle::femtoDreamCollisionId, col.globalIndex());
     auto partsLambda1 = partsLambda1Part->sliceByCached(aod::femtodreamparticle::femtoDreamCollisionId, col.globalIndex());
+
     auto tmstamp = col.timestamp();
     auto mafneticField = getMagneticFieldTesla(tmstamp);
     registry.get<TH1>(HIST("fProcessedEvents"))->Fill(0);
@@ -302,6 +305,7 @@ struct CFFilter {
         if (partsLambda0.size() >= 1 && partsProton0.size() >= 2) {
           for (auto& partLambda : partsLambda0) {
             registry.get<TH1>(HIST("fPtPPL"))->Fill(partLambda.pt());
+            registry.get<TH1>(HIST("fMinvLambda"))->Fill(partLambda.mLambda());
             if (!pairCleanerTV.isCleanPair(partLambda, partLambda, partsFemto)) {
               continue;
             }
@@ -331,6 +335,7 @@ struct CFFilter {
           if (partsLambda1.size() >= 1 && partsProton1.size() >= 2) {
             for (auto& partLambda : partsLambda1) {
               registry.get<TH1>(HIST("fPtAntiPPL"))->Fill(partLambda.pt());
+              registry.get<TH1>(HIST("fMinvAntiLambda"))->Fill(partLambda.mAntiLambda());
               if (!pairCleanerTV.isCleanPair(partLambda, partLambda, partsFemto)) {
                 continue;
               }

PWGCF/DataModel/FemtoDerived.h

@@ -29,6 +29,7 @@ namespace femtodreamcollision
 {
 DECLARE_SOA_COLUMN(MultV0M, multV0M, float);       //! V0M multiplicity
 DECLARE_SOA_COLUMN(Sphericity, sphericity, float); //! Sphericity of the event
+
 } // namespace femtodreamcollision
 
 DECLARE_SOA_TABLE(FemtoDreamCollisions, "AOD", "FEMTODREAMCOLS",
@@ -112,6 +113,7 @@ DECLARE_SOA_COLUMN(TransRadius, transRadius, float); //! Transverse radius of th
 DECLARE_SOA_COLUMN(DecayVtxX, decayVtxX, float);     //! X position of the decay vertex
 DECLARE_SOA_COLUMN(DecayVtxY, decayVtxY, float);     //! Y position of the decay vertex
 DECLARE_SOA_COLUMN(DecayVtxZ, decayVtxZ, float);     //! Z position of the decay vertex
+DECLARE_SOA_COLUMN(MKaon, mKaon, float);             //! The invariant mass of V0 candidate, assuming kaon
 
 } // namespace femtodreamparticle
 DECLARE_SOA_TABLE(FemtoDreamParticles, "AOD", "FEMTODREAMPARTS",
@@ -161,6 +163,7 @@ DECLARE_SOA_TABLE(FemtoDreamDebugParticles, "AOD", "FEMTODEBUGPARTS",
                   femtodreamparticle::DecayVtxX,
                   femtodreamparticle::DecayVtxY,
                   femtodreamparticle::DecayVtxZ,
+                  femtodreamparticle::MKaon,
                   femtodreamparticle::TPCCrossedRowsOverFindableCls<track::TPCNClsFindable, femtodreamparticle::TPCNClsCrossedRows>,
                   pidtpc_tiny::TPCNSigmaEl<pidtpc_tiny::TPCNSigmaStoreEl>,
                   pidtpc_tiny::TPCNSigmaPi<pidtpc_tiny::TPCNSigmaStorePi>,

PWGCF/FemtoDream/FemtoDreamTrackSelection.h

@@ -18,11 +18,15 @@
 #define ANALYSIS_TASKS_PWGCF_FEMTODREAM_FEMTODREAMTRACKSELECTION_H_
 
 #include "PWGCF/DataModel/FemtoDerived.h"
+#include "Common/DataModel/TrackSelectionTables.h"
+#include "Common/Core/TrackSelection.h"
+#include "Common/Core/TrackSelectionDefaults.h"
 #include "FemtoDreamObjectSelection.h"
 
 #include "ReconstructionDataFormats/PID.h"
 #include "Framework/HistogramRegistry.h"
 #include <cmath>
+#include <iostream>
 
 using namespace o2::framework;
 
@@ -379,7 +383,8 @@ bool FemtoDreamTrackSelection::isSelectedMinimal(T const& track)
   const auto itsNClsIB = track.itsNClsInnerBarrel();
   const auto dcaXY = track.dcaXY();
   const auto dcaZ = track.dcaZ();
-  const auto dca = std::sqrt(pow(dcaXY, 2.) + pow(dcaZ, 2.));
+  const auto dca = track.dcaXY(); // Accordingly to FemtoDream in AliPhysics  as well as LF analysis,
+                                  // only dcaXY should be checked; NOT std::sqrt(pow(dcaXY, 2.) + pow(dcaZ, 2.))
   std::vector<float> pidTPC, pidTOF;
   for (auto it : mPIDspecies) {
     pidTPC.push_back(getNsigmaTPC(track, it));

PWGCF/FemtoDream/FemtoDreamV0Selection.h

@@ -26,6 +26,8 @@
 #include "Common/Core/RecoDecay.h"
 #include "Framework/HistogramRegistry.h"
 
+#include <iostream>
+
 using namespace o2::framework;
 
 namespace o2::analysis::femtoDream
@@ -75,7 +77,10 @@ class FemtoDreamV0Selection : public FemtoDreamObjectSelection<float, femtoDream
                             TranRadV0Max(-9999999.),
                             DecVtxMax(-9999999.),
                             fInvMassLowLimit(1.05),
-                            fInvMassUpLimit(1.3){};
+                            fInvMassUpLimit(1.3),
+                            fRejectKaon(false),
+                            fInvMassKaonLowLimit(0.48),
+                            fInvMassKaonUpLimit(0.515){};
   /// Initializes histograms for the task
   template <o2::aod::femtodreamparticle::ParticleType part, o2::aod::femtodreamparticle::ParticleType daugh, typename cutContainerType>
   void init(HistogramRegistry* registry);
@@ -165,6 +170,16 @@ class FemtoDreamV0Selection : public FemtoDreamObjectSelection<float, femtoDream
     fInvMassUpLimit = upLimit;
   }
 
+  /// Set limit for the kaon rejection on the invariant mass
+  /// \param lowLimit Lower limit for the invariant mass distribution
+  /// \param upLimit Upper limit for the invariant mass distribution
+  void setKaonInvMassLimits(float lowLimit, float upLimit)
+  {
+    fRejectKaon = true;
+    fInvMassKaonLowLimit = lowLimit;
+    fInvMassKaonUpLimit = upLimit;
+  }
+
   void setChildRejectNotPropagatedTracks(femtoDreamV0Selection::ChildTrackType child, bool reject)
   {
     if (child == femtoDreamV0Selection::kPosTrack) {
@@ -193,6 +208,10 @@ class FemtoDreamV0Selection : public FemtoDreamObjectSelection<float, femtoDream
   float fInvMassLowLimit;
   float fInvMassUpLimit;
 
+  bool fRejectKaon;
+  float fInvMassKaonLowLimit;
+  float fInvMassKaonUpLimit;
+
   FemtoDreamTrackSelection PosDaughTrack;
   FemtoDreamTrackSelection NegDaughTrack;
 
@@ -314,6 +333,12 @@ bool FemtoDreamV0Selection::isSelectedMinimal(C const& col, V const& v0, T const
   if ((invMassLambda < fInvMassLowLimit or invMassLambda > fInvMassUpLimit) and (invMassAntiLambda < fInvMassLowLimit or invMassAntiLambda > fInvMassUpLimit)) {
     return false;
   }
+  if (fRejectKaon) {
+    const float invMassKaon = v0.mK0Short();
+    if (invMassKaon > fInvMassKaonLowLimit && invMassKaon < fInvMassKaonUpLimit) {
+      return false;
+    }
+  }
   if (nPtV0MinSel > 0 && pT < pTV0Min) {
     return false;
   }
@@ -418,7 +443,33 @@ std::array<cutContainerType, 5> FemtoDreamV0Selection::getCutContainer(C const&
   cutContainerType output = 0;
   size_t counter = 0;
 
+  auto lambdaMassNominal = TDatabasePDG::Instance()->GetParticle(3122)->Mass();
+  auto lambdaMassHypothesis = v0.mLambda();
+  auto antiLambdaMassHypothesis = v0.mAntiLambda();
+  auto diffLambda = abs(lambdaMassNominal - lambdaMassHypothesis);
+  auto diffAntiLambda = abs(antiLambdaMassHypothesis - lambdaMassHypothesis);
+
   float sign = 0.;
+  int nSigmaPIDMax = PosDaughTrack.getSigmaPIDMax();
+  auto nSigmaPrNeg = negTrack.tpcNSigmaPr();
+  auto nSigmaPiPos = posTrack.tpcNSigmaPi();
+  auto nSigmaPiNeg = negTrack.tpcNSigmaPi();
+  auto nSigmaPrPos = posTrack.tpcNSigmaPr();
+  // check the mass and the PID of daughters
+  if (abs(nSigmaPrNeg) < nSigmaPIDMax && abs(nSigmaPiPos) < nSigmaPIDMax && diffAntiLambda > diffLambda) {
+    sign = -1.;
+  } else if (abs(nSigmaPrPos) < nSigmaPIDMax && abs(nSigmaPiNeg) < nSigmaPIDMax && diffAntiLambda < diffLambda) {
+    sign = 1.;
+  }
+  // if it happens that none of these are true, ignore the invariant mass
+  else {
+    if (abs(nSigmaPrNeg) < nSigmaPIDMax && abs(nSigmaPiPos) < nSigmaPIDMax) {
+      sign = -1.;
+    } else if (abs(nSigmaPrPos) < nSigmaPIDMax && abs(nSigmaPiNeg) < nSigmaPIDMax) {
+      sign = 1.;
+    }
+  }
+
   const auto pT = v0.pt();
   const auto tranRad = v0.v0radius();
   const auto dcaDaughv0 = v0.dcaV0daughters();
@@ -428,23 +479,8 @@ std::array<cutContainerType, 5> FemtoDreamV0Selection::getCutContainer(C const&
   float observable = 0.;
   for (auto& sel : mSelections) {
     const auto selVariable = sel.getSelectionVariable();
-    if (selVariable == femtoDreamV0Selection::kV0Sign) {
-      int nSigmaPIDMax = PosDaughTrack.getSigmaPIDMax();
-      if (sel.getSelectionValue() < 0) {
-        auto nSigmaPr = negTrack.tpcNSigmaPr();
-        auto nSigmaPi = posTrack.tpcNSigmaPi();
-        if (abs(nSigmaPr) < nSigmaPIDMax && abs(nSigmaPi) < nSigmaPIDMax) {
-          sign = -1.;
-        }
-      } else {
-        auto nSigmaPi = negTrack.tpcNSigmaPi();
-        auto nSigmaPr = posTrack.tpcNSigmaPr();
-        if (abs(nSigmaPr) < nSigmaPIDMax && abs(nSigmaPi) < nSigmaPIDMax) {
-          sign = 1.;
-        }
-      }
-      sel.checkSelectionSetBit(sign, output, counter);
-    } else if (selVariable == femtoDreamV0Selection::kDecVtxMax) {
+
+    if (selVariable == femtoDreamV0Selection::kDecVtxMax) {
       for (size_t i = 0; i < decVtx.size(); ++i) {
         auto decVtxValue = decVtx.at(i);
         sel.checkSelectionSetBit(decVtxValue, output, counter);

PWGCF/FemtoDream/femtoDreamProducerReducedTask.cxx

@@ -186,6 +186,7 @@ struct femtoDreamProducerReducedTask {
                           -999.,
                           -999.,
                           -999.,
+                          -999.,
                           -999.);
       }
     }

PWGCF/FemtoDream/femtoDreamProducerTask.cxx

@@ -93,6 +93,9 @@ struct femtoDreamProducerTask {
   Configurable<bool> ConfEvtOfflineCheck{"ConfEvtOfflineCheck", false, "Evt sel: check for offline selection"};
 
   Configurable<bool> ConfStoreV0{"ConfStoreV0", true, "True: store V0 table"};
+  // just sanity check to make sure in case there are problems in convertion or MC production it does not affect results
+  Configurable<bool> ConfRejectNotPropagatedTracks{"ConfRejectNotPropagatedTracks", false, "True: reject not propagated tracks"};
+  Configurable<bool> ConfRejectITSHitandTOFMissing{"ConfRejectITSHitandTOFMissing", false, "True: reject if neither ITS hit nor TOF timing satisfied"};
 
   FemtoDreamTrackSelection trackCuts;
   Configurable<std::vector<float>> ConfTrkCharge{FemtoDreamTrackSelection::getSelectionName(femtoDreamTrackSelection::kSign, "ConfTrk"), std::vector<float>{-1, 1}, FemtoDreamTrackSelection::getSelectionHelper(femtoDreamTrackSelection::kSign, "Track selection: ")};
@@ -111,6 +114,7 @@ struct femtoDreamProducerTask {
   Configurable<std::vector<int>> ConfTrkTPIDspecies{"ConfTrkTPIDspecies", std::vector<int>{o2::track::PID::Pion, o2::track::PID::Kaon, o2::track::PID::Proton, o2::track::PID::Deuteron}, "Trk sel: Particles species for PID"};
 
   FemtoDreamV0Selection v0Cuts;
+  TrackSelection* o2PhysicsTrackSelection;
   /// \todo Labeled array (see Track-Track task)
 
   Configurable<std::vector<float>> ConfV0Sign{FemtoDreamV0Selection::getSelectionName(femtoDreamV0Selection::kV0Sign, "ConfV0"), std::vector<float>{-1, 1}, FemtoDreamV0Selection::getSelectionHelper(femtoDreamV0Selection::kV0Sign, "V0 selection: ")};
@@ -122,6 +126,8 @@ struct femtoDreamProducerTask {
   Configurable<std::vector<float>> V0TranRadV0Max{FemtoDreamV0Selection::getSelectionName(femtoDreamV0Selection::kTranRadV0Max, "ConfV0"), std::vector<float>{100.f}, FemtoDreamV0Selection::getSelectionHelper(femtoDreamV0Selection::kTranRadV0Max, "V0 selection: ")};
   Configurable<std::vector<float>> V0DecVtxMax{FemtoDreamV0Selection::getSelectionName(femtoDreamV0Selection::kDecVtxMax, "ConfV0"), std::vector<float>{100.f}, FemtoDreamV0Selection::getSelectionHelper(femtoDreamV0Selection::kDecVtxMax, "V0 selection: ")};
 
+  Configurable<std::vector<float>> ConfV0DaughCharge{"ConfV0DaughCharge", std::vector<float>{-1, 1}, "V0 Daugh sel: Charge"};
+  Configurable<std::vector<float>> ConfDaughEta{"ConfDaughEta", std::vector<float>{0.8f}, "V0 Daugh sel: max eta"};
   Configurable<std::vector<float>> ConfV0DaughTPCnclsMin{"ConfV0DaughTPCnclsMin", std::vector<float>{80.f, 70.f, 60.f}, "V0 Daugh sel: Min. nCls TPC"};
   Configurable<std::vector<float>> ConfV0DaughDCAMin{"ConfV0DaughDCAMin", std::vector<float>{0.05f, 0.06f}, "V0 Daugh sel:  Max. DCA Daugh to PV (cm)"};
   Configurable<std::vector<float>> ConfV0DaughPIDnSigmaMax{"ConfV0DaughPIDnSigmaMax", std::vector<float>{5.f, 4.f}, "V0 Daugh sel: Max. PID nSigma TPC"};
@@ -131,6 +137,10 @@ struct femtoDreamProducerTask {
   Configurable<float> ConfInvMassLowLimit{"ConfInvMassLowLimit", 1.05, "Lower limit of the V0 invariant mass"};
   Configurable<float> ConfInvMassUpLimit{"ConfInvMassUpLimit", 1.30, "Upper limit of the V0 invariant mass"};
 
+  Configurable<bool> ConfRejectKaons{"ConfRejectKaons", false, "Switch to reject kaons"};
+  Configurable<float> ConfInvKaonMassLowLimit{"ConfInvKaonMassLowLimit", 0.48, "Lower limit of the V0 invariant mass for Kaon rejection"};
+  Configurable<float> ConfInvKaonMassUpLimit{"ConfInvKaonMassUpLimit", 0.515, "Upper limit of the V0 invariant mass for Kaon rejection"};
+
   /// \todo should we add filter on min value pT/eta of V0 and daughters?
   /*Filter v0Filter = (nabs(aod::v0data::x) < V0DecVtxMax.value) &&
                     (nabs(aod::v0data::y) < V0DecVtxMax.value) &&
@@ -167,20 +177,30 @@ struct femtoDreamProducerTask {
       v0Cuts.setSelection(ConfDCAV0DaughMax, femtoDreamV0Selection::kDCAV0DaughMax, femtoDreamSelection::kUpperLimit);
       v0Cuts.setSelection(ConfCPAV0Min, femtoDreamV0Selection::kCPAV0Min, femtoDreamSelection::kLowerLimit);
 
-      v0Cuts.setChildCuts(femtoDreamV0Selection::kPosTrack, ConfTrkCharge, femtoDreamTrackSelection::kSign, femtoDreamSelection::kEqual);
-      v0Cuts.setChildCuts(femtoDreamV0Selection::kPosTrack, ConfTrkEta, femtoDreamTrackSelection::kEtaMax, femtoDreamSelection::kAbsUpperLimit);
+      v0Cuts.setChildCuts(femtoDreamV0Selection::kPosTrack, ConfV0DaughCharge, femtoDreamTrackSelection::kSign, femtoDreamSelection::kEqual);
+      v0Cuts.setChildCuts(femtoDreamV0Selection::kPosTrack, ConfDaughEta, femtoDreamTrackSelection::kEtaMax, femtoDreamSelection::kAbsUpperLimit);
       v0Cuts.setChildCuts(femtoDreamV0Selection::kPosTrack, ConfV0DaughTPCnclsMin, femtoDreamTrackSelection::kTPCnClsMin, femtoDreamSelection::kLowerLimit);
       v0Cuts.setChildCuts(femtoDreamV0Selection::kPosTrack, ConfV0DaughDCAMin, femtoDreamTrackSelection::kDCAMin, femtoDreamSelection::kAbsLowerLimit);
       v0Cuts.setChildCuts(femtoDreamV0Selection::kPosTrack, ConfV0DaughPIDnSigmaMax, femtoDreamTrackSelection::kPIDnSigmaMax, femtoDreamSelection::kAbsUpperLimit);
-      v0Cuts.setChildCuts(femtoDreamV0Selection::kNegTrack, ConfTrkCharge, femtoDreamTrackSelection::kSign, femtoDreamSelection::kEqual);
-      v0Cuts.setChildCuts(femtoDreamV0Selection::kNegTrack, ConfTrkEta, femtoDreamTrackSelection::kEtaMax, femtoDreamSelection::kAbsUpperLimit);
+      v0Cuts.setChildCuts(femtoDreamV0Selection::kNegTrack, ConfV0DaughCharge, femtoDreamTrackSelection::kSign, femtoDreamSelection::kEqual);
+      v0Cuts.setChildCuts(femtoDreamV0Selection::kNegTrack, ConfDaughEta, femtoDreamTrackSelection::kEtaMax, femtoDreamSelection::kAbsUpperLimit);
       v0Cuts.setChildCuts(femtoDreamV0Selection::kNegTrack, ConfV0DaughTPCnclsMin, femtoDreamTrackSelection::kTPCnClsMin, femtoDreamSelection::kLowerLimit);
       v0Cuts.setChildCuts(femtoDreamV0Selection::kNegTrack, ConfV0DaughDCAMin, femtoDreamTrackSelection::kDCAMin, femtoDreamSelection::kAbsLowerLimit);
       v0Cuts.setChildCuts(femtoDreamV0Selection::kNegTrack, ConfV0DaughPIDnSigmaMax, femtoDreamTrackSelection::kPIDnSigmaMax, femtoDreamSelection::kAbsUpperLimit);
       v0Cuts.setChildPIDSpecies(femtoDreamV0Selection::kPosTrack, ConfV0DaughTPIDspecies);
       v0Cuts.setChildPIDSpecies(femtoDreamV0Selection::kNegTrack, ConfV0DaughTPIDspecies);
       v0Cuts.init<aod::femtodreamparticle::ParticleType::kV0, aod::femtodreamparticle::ParticleType::kV0Child, aod::femtodreamparticle::cutContainerType>(&qaRegistry);
       v0Cuts.setInvMassLimits(ConfInvMassLowLimit, ConfInvMassUpLimit);
+      v0Cuts.setChildRejectNotPropagatedTracks(femtoDreamV0Selection::kPosTrack, ConfRejectNotPropagatedTracks);
+      v0Cuts.setChildRejectNotPropagatedTracks(femtoDreamV0Selection::kNegTrack, ConfRejectNotPropagatedTracks);
+
+      if (ConfRejectKaons) {
+        v0Cuts.setKaonInvMassLimits(ConfInvKaonMassLowLimit, ConfInvKaonMassUpLimit);
+      }
+      if (ConfRejectITSHitandTOFMissing) {
+        o2PhysicsTrackSelection = new TrackSelection(getGlobalTrackSelection());
+        o2PhysicsTrackSelection->SetRequireHitsInITSLayers(1, {0, 1, 2, 3});
+      }
     }
   }
 
@@ -260,6 +280,7 @@ struct femtoDreamProducerTask {
                          -999.,
                          -999.,
                          -999.,
+                         -999.,
                          -999.);
       }
     }
@@ -276,6 +297,13 @@ struct femtoDreamProducerTask {
         if (!v0Cuts.isSelectedMinimal(col, v0, postrack, negtrack)) {
           continue;
         }
+
+        if (ConfRejectITSHitandTOFMissing) {
+          // Uncomment only when TOF timing is solved
+          // bool itsHit = o2PhysicsTrackSelection->IsSelected(postrack, TrackSelection::TrackCuts::kITSHits);
+          // bool itsHit = o2PhysicsTrackSelection->IsSelected(negtrack, TrackSelection::TrackCuts::kITSHits);
+        }
+
         v0Cuts.fillQA<aod::femtodreamparticle::ParticleType::kV0, aod::femtodreamparticle::ParticleType::kV0Child>(col, v0, postrack, negtrack); ///\todo fill QA also for daughters
         auto cutContainerV0 = v0Cuts.getCutContainer<aod::femtodreamparticle::cutContainerType>(col, v0, postrack, negtrack);
 
@@ -322,6 +350,7 @@ struct femtoDreamProducerTask {
                              -999.,
                              -999.,
                              -999.,
+                             -999.,
                              -999.); // QA for positive daughter
             outputDebugParts(negtrack.sign(),
                              (uint8_t)negtrack.tpcNClsFound(),
@@ -348,6 +377,7 @@ struct femtoDreamProducerTask {
                              -999.,
                              -999.,
                              -999.,
+                             -999.,
                              -999.); // QA for negative daughter
             outputDebugParts(-999.,
                              -999.,
@@ -374,7 +404,8 @@ struct femtoDreamProducerTask {
                              v0.v0radius(),
                              v0.x(),
                              v0.y(),
-                             v0.z()); // QA for V0
+                             v0.z(),
+                             v0.mK0Short()); // QA for V0
           }
         }
       }