Updated derived-to-derived producer to include the bitmask check

Author: Laura Serksnyte

PWGCF/FemtoDream/TableProducer/femtoDreamProducerTaskForSpecificAnalysis.cxx

@@ -13,24 +13,24 @@
 /// \brief Tasks that reads the track tables and creates track triplets; only three identical particles can be used
 /// \author Laura Serksnyte, TU München, laura.serksnyte@tum.de
 
-#include <vector>
-#include <string>
+#include "PWGCF/DataModel/FemtoDerived.h"
+#include "PWGCF/FemtoDream/Core/femtoDreamContainerThreeBody.h"
+#include "PWGCF/FemtoDream/Core/femtoDreamDetaDphiStar.h"
+#include "PWGCF/FemtoDream/Core/femtoDreamEventHisto.h"
+#include "PWGCF/FemtoDream/Core/femtoDreamPairCleaner.h"
+#include "PWGCF/FemtoDream/Core/femtoDreamParticleHisto.h"
+#include "PWGCF/FemtoDream/Core/femtoDreamUtils.h"
+
+#include "Framework/ASoAHelpers.h"
 #include "Framework/AnalysisTask.h"
-#include "Framework/runDataProcessing.h"
 #include "Framework/HistogramRegistry.h"
-#include "Framework/ASoAHelpers.h"
+#include "Framework/O2DatabasePDGPlugin.h"
 #include "Framework/RunningWorkflowInfo.h"
 #include "Framework/StepTHn.h"
-#include "Framework/O2DatabasePDGPlugin.h"
-#include "TDatabasePDG.h"
+#include "Framework/runDataProcessing.h"
 
-#include "PWGCF/DataModel/FemtoDerived.h"
-#include "PWGCF/FemtoDream/Core/femtoDreamParticleHisto.h"
-#include "PWGCF/FemtoDream/Core/femtoDreamEventHisto.h"
-#include "PWGCF/FemtoDream/Core/femtoDreamPairCleaner.h"
-#include "PWGCF/FemtoDream/Core/femtoDreamContainerThreeBody.h"
-#include "PWGCF/FemtoDream/Core/femtoDreamDetaDphiStar.h"
-#include "PWGCF/FemtoDream/Core/femtoDreamUtils.h"
+#include <string>
+#include <vector>
 
 using namespace o2;
 using namespace o2::analysis::femtoDream;
@@ -49,6 +49,10 @@ struct femtoDreamProducerTaskForSpecificAnalysis {
   float mMassOne = -999, mMassTwo = -999, mMassThree = -999;
   int collisions = 0;
 
+  // Require bitmask selection for candidates
+  Configurable<bool> ConfRequireBitmask{"ConfRequireBitmask", false, "Require bitmask selection for candidates"};
+
+  // Number of candidates required
   Configurable<int> ConfNumberOfTracks{"ConfNumberOfTracks", 3, "Number of tracks"};
   Configurable<int> ConfNumberOfV0{"ConfNumberOfV0", 0, "Number of V0"};
   Configurable<int> ConfNumberOfCascades{"ConfNumberOfCascades", 0, "Number of Cascades"};
@@ -57,6 +61,8 @@ struct femtoDreamProducerTaskForSpecificAnalysis {
   Configurable<float> ConfPIDthrMom{"ConfPIDthrMom", 1.f, "Momentum threshold from which TPC and TOF are required for PID"};
   Configurable<o2::aod::femtodreamparticle::cutContainerType> ConfTPCPIDBit{"ConfTPCPIDBit", 16, "PID TPC bit from cutCulator "};
   Configurable<o2::aod::femtodreamparticle::cutContainerType> ConfTPCTOFPIDBit{"ConfTPCTOFPIDBit", 8, "PID TPCTOF bit from cutCulator"};
+  Configurable<o2::aod::femtodreamparticle::cutContainerType> ConfCutPart{"ConfCutPart", 0, "Track - Selection bit from cutCulator for part"};
+  Configurable<o2::aod::femtodreamparticle::cutContainerType> ConfCutPartAntiPart{"ConfCutPartAntiPart", 0, "Track - Selection bit from cutCulator for antipart"};
 
   /// Partition for selected particles
   Partition<aod::FDParticles> SelectedParts = (aod::femtodreamparticle::partType == uint8_t(aod::femtodreamparticle::ParticleType::kTrack)) &&
@@ -67,6 +73,12 @@ struct femtoDreamProducerTaskForSpecificAnalysis {
   Configurable<float> Conf_maxInvMass_V0{"Conf_maxInvMass_V0", 1.15, "Maximum invariant mass of V0 (particle)"};
   Configurable<float> Conf_minInvMassAnti_V0{"Conf_minInvMassAnti_V0", 1.08, "Minimum invariant mass of V0 (antiparticle)"};
   Configurable<float> Conf_maxInvMassAnti_V0{"Conf_maxInvMassAnti_V0", 1.15, "Maximum invariant mass of V0 (antiparticle)"};
+  Configurable<o2::aod::femtodreamparticle::cutContainerType> ConfCutV0_SameForAntipart{"ConfCutV0_SameForAntipart", 0, "V0 - Selection bit from cutCulator for part/antipart"};
+  Configurable<o2::aod::femtodreamparticle::cutContainerType> Conf_ChildPos_CutV0{"Conf_ChildPos_CutV0", 149, "Selection bit for positive child of V0"};
+  Configurable<o2::aod::femtodreamparticle::cutContainerType> Conf_ChildPos_TPCBitV0{"Conf_ChildPos_TPCBitV0", 2, "PID TPC bit for positive child of V0"};
+  Configurable<o2::aod::femtodreamparticle::cutContainerType> Conf_ChildNeg_CutV0{"Conf_ChildNeg_CutV0", 149, "Selection bit for negative child of V0"};
+  Configurable<o2::aod::femtodreamparticle::cutContainerType> Conf_ChildNeg_TPCBitV0{"Conf_ChildNeg_TPCBitV0", 2, "PID TPC bit for negative child of V0"};
+
   /// Cascade selection
   Configurable<float> Conf_minInvMass_Cascade{"Conf_minInvMass_Cascade", 1.2, "Minimum invariant mass of Cascade (particle)"};
   Configurable<float> Conf_maxInvMass_Cascade{"Conf_maxInvMass_Cascade", 1.5, "Maximum invariant mass of Cascade (particle)"};
@@ -95,8 +107,12 @@ struct femtoDreamProducerTaskForSpecificAnalysis {
   void init(InitContext&)
   {
     EventRegistry.add("hStatistiscs", ";bin;Entries", kTH1F, {{3, 0, 3}});
-    // get bit for the collision mask
+    // Never run V0s and Cascades together as this will DOUBLE the track number and induce self correlations
+    if ((doprocessCollisionsWithNTracksAndNCascades && doprocessCollisionsWithNTracksAndNV0)) {
+      LOG(fatal) << "Never run V0s and Cascades together as this will DOUBLE the track number and induce self correlations!";
+    }
   }
+
   /// This function stores accepted collisions in derived data
   /// @tparam PartitionType
   /// @tparam PartType
@@ -105,27 +121,57 @@ struct femtoDreamProducerTaskForSpecificAnalysis {
   /// @param groupSelectedV0s partition for the second particle passed by the process function
   /// @param parts femtoDreamParticles table
   template <bool isMC, typename PartitionType, typename PartType>
-  void createSpecifiedDerivedData(o2::aod::FDCollision& col, PartitionType groupSelectedTracks, PartitionType groupSelectedV0s, PartType parts)
+  void createSpecifiedDerivedData(const o2::aod::FDCollision& col, PartitionType groupSelectedTracks, PartitionType groupSelectedV0s, PartType parts)
   {
     /// check tracks
     int tracksCount = 0;
     int antitracksCount = 0;
-    for (auto& part : groupSelectedTracks) {
+    for (const auto& part : groupSelectedTracks) {
       if (part.cut() & 1) {
-        antitracksCount++;
+        if (!ConfRequireBitmask || ncheckbit(part.cut(), ConfCutPartAntiPart)) {
+          antitracksCount++;
+        }
       } else {
-        tracksCount++;
+        if (!ConfRequireBitmask || ncheckbit(part.cut(), ConfCutPart)) {
+          tracksCount++;
+        }
       }
     }
 
     /// check V0s
     int V0Count = 0;
     int antiV0Count = 0;
-    for (auto& V0 : groupSelectedV0s) {
+    for (const auto& V0 : groupSelectedV0s) {
       if ((V0.mLambda() > Conf_minInvMass_V0) && (V0.mLambda() < Conf_maxInvMass_V0)) {
-        V0Count++;
+        if (ConfRequireBitmask) {
+          if (ncheckbit(V0.cut(), ConfCutV0_SameForAntipart)) {
+            const auto& posChild = parts.iteratorAt(V0.index() - 2);
+            const auto& negChild = parts.iteratorAt(V0.index() - 1);
+            if (((posChild.cut() & Conf_ChildPos_CutV0) == Conf_ChildPos_CutV0 &&
+                 (posChild.pidcut() & Conf_ChildPos_TPCBitV0) == Conf_ChildPos_TPCBitV0 &&
+                 (negChild.cut() & Conf_ChildNeg_CutV0) == Conf_ChildNeg_CutV0 &&
+                 (negChild.pidcut() & Conf_ChildNeg_TPCBitV0) == Conf_ChildNeg_TPCBitV0)) {
+              V0Count++;
+            }
+          }
+        } else {
+          V0Count++;
+        }
       } else if ((V0.mAntiLambda() > Conf_minInvMassAnti_V0) && (V0.mAntiLambda() < Conf_maxInvMassAnti_V0)) {
-        antiV0Count++;
+        if (ConfRequireBitmask) {
+          if (ncheckbit(V0.cut(), ConfCutV0_SameForAntipart)) {
+            const auto& posChild = parts.iteratorAt(V0.index() - 2);
+            const auto& negChild = parts.iteratorAt(V0.index() - 1);
+            if (((posChild.cut() & Conf_ChildPos_CutV0) == Conf_ChildPos_CutV0 &&
+                 (posChild.pidcut() & Conf_ChildNeg_TPCBitV0) == Conf_ChildNeg_TPCBitV0 && // exchanged values because checking antiparticle daughters and pid of particles exchange
+                 (negChild.cut() & Conf_ChildNeg_CutV0) == Conf_ChildNeg_CutV0 &&
+                 (negChild.pidcut() & Conf_ChildPos_TPCBitV0) == Conf_ChildPos_TPCBitV0)) { // exchanged values because checking antiparticle daughters and pid of particles exchange
+              antiV0Count++;
+            }
+          }
+        } else {
+          antiV0Count++;
+        }
       }
     }
 
@@ -134,7 +180,7 @@ struct femtoDreamProducerTaskForSpecificAnalysis {
     if ((V0Count >= ConfNumberOfV0 && tracksCount >= ConfNumberOfTracks) || (antiV0Count >= ConfNumberOfV0 && antitracksCount >= ConfNumberOfTracks)) {
       EventRegistry.fill(HIST("hStatistiscs"), 1);
       outputCollision(col.posZ(), col.multV0M(), col.multNtr(), col.sphericity(), col.magField());
-      for (auto& femtoParticle : parts) {
+      for (const auto& femtoParticle : parts) {
         if (aod::femtodreamparticle::ParticleType::kTrack == femtoParticle.partType()) {
           std::vector<int> childIDs = {0, 0};
           outputParts(outputCollision.lastIndex(),
@@ -197,8 +243,8 @@ struct femtoDreamProducerTaskForSpecificAnalysis {
   /// process function to create derived data with only collisions containing n tracks
   /// \param col subscribe to the collision table (Data)
   /// \param parts subscribe to the femtoDreamParticleTable
-  void processCollisionsWithNTracksAndNV0(o2::aod::FDCollision& col,
-                                          o2::aod::FDParticles& parts)
+  void processCollisionsWithNTracksAndNV0(const o2::aod::FDCollision& col,
+                                          const o2::aod::FDParticles& parts)
   {
     EventRegistry.fill(HIST("hStatistiscs"), 0);
     auto thegroupSelectedParts = SelectedParts->sliceByCached(aod::femtodreamparticle::fdCollisionId, col.globalIndex(), cache);
@@ -216,13 +262,13 @@ struct femtoDreamProducerTaskForSpecificAnalysis {
   /// @param groupSelectedV0s partition for the second particle passed by the process function
   /// @param parts femtoDreamParticles table
   template <bool isMC, typename PartitionType, typename PartType>
-  void createSpecifiedDerivedData_TrkCascade(o2::aod::FDCollision& col, PartitionType groupSelectedTracks, PartitionType groupSelectedCascades, PartType parts)
+  void createSpecifiedDerivedData_TrkCascade(const o2::aod::FDCollision& col, PartitionType groupSelectedTracks, PartitionType groupSelectedCascades, PartType parts)
   {
 
     /// check tracks
     int tracksCount = 0;
     int antitracksCount = 0;
-    for (auto& part : groupSelectedTracks) {
+    for (const auto& part : groupSelectedTracks) {
       if (part.cut() & 1) {
         antitracksCount++;
       } else {
@@ -233,7 +279,7 @@ struct femtoDreamProducerTaskForSpecificAnalysis {
     /// check Cascades
     int CascadeCount = 0;
     int antiCascadeCount = 0;
-    for (auto& casc : groupSelectedCascades) {
+    for (const auto& casc : groupSelectedCascades) {
       if ((casc.cut() & kSignPlusMask) == kSignPlusMask) {
         CascadeCount++;
       } else {
@@ -247,7 +293,7 @@ struct femtoDreamProducerTaskForSpecificAnalysis {
       EventRegistry.fill(HIST("hStatistiscs"), 1);
       outputCollision(col.posZ(), col.multV0M(), col.multNtr(), col.sphericity(), col.magField());
 
-      for (auto& femtoParticle : parts) {
+      for (const auto& femtoParticle : parts) {
         if (aod::femtodreamparticle::ParticleType::kTrack == femtoParticle.partType()) {
           std::vector<int> childIDs = {0, 0};
           outputParts(outputCollision.lastIndex(),
@@ -330,8 +376,8 @@ struct femtoDreamProducerTaskForSpecificAnalysis {
   /// process function to create derived data with only collisions containing n tracks
   /// \param col subscribe to the collision table (Data)
   /// \param parts subscribe to the femtoDreamParticleTable
-  void processCollisionsWithNTracksAndNCascades(o2::aod::FDCollision& col,
-                                                o2::aod::FDParticles& parts)
+  void processCollisionsWithNTracksAndNCascades(const o2::aod::FDCollision& col,
+                                                const o2::aod::FDParticles& parts)
   {
     EventRegistry.fill(HIST("hStatistiscs"), 0);
     auto thegroupSelectedParts = SelectedParts->sliceByCached(aod::femtodreamparticle::fdCollisionId, col.globalIndex(), cache);