Added extra QA checks

Author: gbencedi

PWGMM/Mult/Tasks/dndetaMFTPbPb.cxx

@@ -60,19 +60,10 @@ using namespace o2::aod::rctsel;
 
 auto static constexpr kMinCharge = 3.f;
 
-AxisSpec ptAxis = {1001, -0.005, 10.005};
-AxisSpec multAxis = {701, -0.5, 700.5, "N_{trk}"};
-AxisSpec zAxis = {60, -30., 30.};
-AxisSpec deltaZAxis = {61, -6.1, 6.1};
-AxisSpec dcaxyAxis = {1000, 0., 5.};
-AxisSpec dcazAxis = {1000, 0., 5.};
-AxisSpec phiAxis = {629, 0, TwoPI, "Rad", "#phi"};
-AxisSpec etaAxis = {20, -4., -2.};
-AxisSpec centAxis{100, 0, 100, "centrality"};
-AxisSpec chiSqAxis = {100, 0., 1000.};
-AxisSpec nclsAxis{10, 0.5, 10.5, "# clusters"};
-
 enum TrkSel {
+  trkSelAll,
+  trkSelWoAmbiguous,
+  trkSelNumReassoc,
   trkSelNCls,
   trkSelChi2Ncl,
   trkSelEta,
@@ -102,9 +93,30 @@ struct DndetaMFTPbPb {
   Configurable<bool> cfgDoIR{"cfgDoIR", false, "Flag to retrieve Interaction rate from CCDB"};
   Configurable<bool> cfgUseIRCut{"cfgUseIRCut", false, "Flag to cut on IR rate"};
   Configurable<bool> cfgIRCrashOnNull{"cfgIRCrashOnNull", false, "Flag to avoid CTP RateFetcher crash"};
-  Configurable<std::string> cfgIRSource{"cfgIRSource", "T0VTX", "Estimator of the interaction rate (Pb-Pb: ZNC hadronic)"};
+  Configurable<std::string> cfgIRSource{"cfgIRSource", "ZNC hadronic", "Estimator of the interaction rate (Pb-Pb: ZNC hadronic)"};
   Configurable<bool> cfgUseTrackSel{"cfgUseTrackSel", false, "Flag to apply track selection"};
 
+  struct : ConfigurableGroup {
+    ConfigurableAxis interactionRateBins{"interactionRateBins", {500, 0, 50}, "Binning for the interaction rate (kHz)"};
+    ConfigurableAxis occupancyBins{"occupancyBins", {VARIABLE_WIDTH, 0.0f, 250.0f, 500.0f, 750.0f, 1000.0f, 1500.0f, 2000.0f, 3000.0f, 4500.0f, 6000.0f, 8000.0f, 10000.0f, 50000.0f}, "Occupancy"};
+    ConfigurableAxis centralityBins{"centralityBins", {VARIABLE_WIDTH, 0, 5, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100}, "Centrality"};
+    ConfigurableAxis irBins{"irBins", {500, 0, 50}, "Interaction rate (kHz)"};
+    ConfigurableAxis pvBins{"pvBins", {501, -0.5, 500.5}, ""};
+    ConfigurableAxis fv0aMultBins{"fv0aMultBins", {501, -0.5, 500.5}, ""};
+    ConfigurableAxis ft0aMultBins{"ft0aMultBins", {501, -0.5, 500.5}, ""};
+    ConfigurableAxis ft0cMultBins{"ft0cMultBins", {501, -0.5, 500.5}, ""};
+    ConfigurableAxis PtBins{"PtBins", {101, -0.5, 10.5}, "pT binning (GeV/c)"};
+    ConfigurableAxis MultBins{"MultBins", {701, -0.5, 700.5}, "Multiplicity binning"};
+    ConfigurableAxis ZvtxBins{"ZvtxBins", {60, -30., 30.}, "Z-vtx binning (cm)"};
+    ConfigurableAxis DeltaZBins{"DeltaZBins", {120, -6., 6.}, "Delta Z-vtx binning (cm)"};
+    ConfigurableAxis DCAxyBins{"DCAxyBins", {800, -1., 1.}, "DCAxy binning (cm)"};
+    ConfigurableAxis DCAzBins{"DCAzBins", {800, -1., 1.}, "DCAz binning (cm)"};
+    ConfigurableAxis PhiBins{"PhiBins", {629, 0., TwoPI}, "#varphi binning (rad)"};
+    ConfigurableAxis EtaBins{"EtaBins", {20, -4., -2.}, "#eta binning"};
+    ConfigurableAxis ChiSqPerNclBins{"ChiSqPerNclBins", {100, 0, 100}, "#chi^{2} binning"};
+    ConfigurableAxis NclBins{"NclBins", {10, 0.5, 10.5}, "number of clusters binning"};
+  } binOpt;
+
   struct : ConfigurableGroup {
     Configurable<bool> requireRCTFlagChecker{"requireRCTFlagChecker", false, "Check event quality in run condition table"};
     Configurable<std::string> cfgEvtRCTFlagCheckerLabel{"cfgEvtRCTFlagCheckerLabel", "CBT_fw", "Evt sel: RCT flag checker label"};
@@ -123,20 +135,21 @@ struct DndetaMFTPbPb {
     Configurable<bool> useChi2Cut{"useChi2Cut", false, "use track chi2 cut"};
     Configurable<float> maxChi2NCl{"maxChi2NCl", 1000.f, "maximum chi2 per MFT clusters"};
     Configurable<bool> usePtCut{"usePtCut", false, "use track pT cut"};
-    Configurable<double> minPt{"minPt", 0., "minimum pT of the MFT tracks"};
+    Configurable<float> minPt{"minPt", 0., "minimum pT of the MFT tracks"};
     Configurable<bool> requireCA{"requireCA", false, "Use Cellular Automaton track-finding algorithm"};
+    Configurable<bool> excludeAmbiguous{"excludeAmbiguous", false, "Exclude Ambiguous tracks"};
     Configurable<float> maxDCAxy{"maxDCAxy", 0.01f, "Cut on dca XY"};
     Configurable<float> maxDCAz{"maxDCAz", 0.01f, "Cut on dca Z"};
   } trackCuts;
 
   struct : ConfigurableGroup {
-    Configurable<float> maxZvtx{"maxZvtx", 10.0f, "maximum cut on z-vtx (cm)"};
-    Configurable<float> minZvtx{"minZvtx", -10.0f, "minimum cut on z-vtx (cm)"};
+    Configurable<float> maxZvtx{"maxZvtx", 20.0f, "maximum cut on z-vtx (cm)"};
+    Configurable<float> minZvtx{"minZvtx", -20.0f, "minimum cut on z-vtx (cm)"};
     Configurable<bool> useZDiffCut{"useZDiffCut", false, "use Zvtx reco-mc diff. cut"};
     Configurable<float> maxZvtxDiff{"maxZvtxDiff", 1.0f, "max allowed Z vtx difference for reconstruced collisions (cm)"};
     Configurable<bool> requireIsGoodZvtxFT0VsPV{"requireIsGoodZvtxFT0VsPV", true, "require events with PV position along z consistent (within 1 cm) between PV reconstructed using tracks and PV using FT0 A-C time difference"};
     Configurable<bool> requireRejectSameBunchPileup{"requireRejectSameBunchPileup", true, "reject collisions in case of pileup with another collision in the same foundBC"};
-    Configurable<bool> requireNoCollInTimeRangeStrict{"requireNoCollInTimeRangeStrict", false, " requireNoCollInTimeRangeStrict"};
+    Configurable<bool> requireNoCollInTimeRangeStrict{"requireNoCollInTimeRangeStrict", true, " requireNoCollInTimeRangeStrict"};
     Configurable<bool> requireNoCollInRofStrict{"requireNoCollInRofStrict", false, "requireNoCollInRofStrict"};
     Configurable<bool> requireNoCollInRofStandard{"requireNoCollInRofStandard", false, "requireNoCollInRofStandard"};
     Configurable<bool> requireNoHighMultCollInPrevRof{"requireNoHighMultCollInPrevRof", false, "requireNoHighMultCollInPrevRof"};
@@ -149,23 +162,14 @@ struct DndetaMFTPbPb {
     Configurable<float> maxIR{"maxIR", -1, "maximum IR (kHz) collisions"};
   } eventCuts;
 
-  ConfigurableAxis occupancyBins{"occupancyBins", {VARIABLE_WIDTH, 0.0f, 250.0f, 500.0f, 750.0f, 1000.0f, 1500.0f, 2000.0f, 3000.0f, 4500.0f, 6000.0f, 8000.0f, 10000.0f, 50000.0f}, "Occupancy"};
-  ConfigurableAxis centralityBins{"centralityBins", {VARIABLE_WIDTH, 0, 5, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100}, "Centrality"};
-  ConfigurableAxis irBins{"irBins", {500, 0, 50}, "Interaction rate (kHz)"};
-  ConfigurableAxis pvBins{"pvBins", {501, -0.5, 500.5}, ""};
-  ConfigurableAxis fv0aMultBins{"fv0aMultBins", {501, -0.5, 500.5}, ""};
-  ConfigurableAxis ft0aMultBins{"ft0aMultBins", {501, -0.5, 500.5}, ""};
-  ConfigurableAxis ft0cMultBins{"ft0cMultBins", {501, -0.5, 500.5}, ""};
-
   Service<o2::framework::O2DatabasePDG> pdg;
-
   Service<ccdb::BasicCCDBManager> ccdb;
   Configurable<int64_t> ccdbNoLaterThan{"ccdbNoLaterThan", std::chrono::duration_cast<std::chrono::milliseconds>(std::chrono::system_clock::now().time_since_epoch()).count(), "latest acceptable timestamp of creation for the object"};
   Configurable<std::string> ccdbUrl{"ccdbUrl", "http://alice-ccdb.cern.ch", "url of the ccdb repository"};
 
   int mRunNumber{-1};
   uint64_t mSOR{0};
-  double mMinSeconds{-1.};
+  float mMinSeconds{-1.};
   std::unordered_map<int, TH2*> gHadronicRate;
   ctpRateFetcher rateFetcher;
   TH2* gCurrentHadronicRate;
@@ -174,12 +178,23 @@ struct DndetaMFTPbPb {
   /// @brief init function, definition of histograms
   void init(InitContext&)
   {
-    AxisSpec pvAxis = {pvBins, "PV", "PVAxis"};
-    AxisSpec multFV0aAxis = {fv0aMultBins, "fv0a", "FV0AMultAxis"};
-    AxisSpec multFT0aAxis = {ft0aMultBins, "ft0a", "FT0AMultAxis"};
-    AxisSpec multFT0cAxis = {ft0cMultBins, "ft0c", "FT0CMultAxis"};
-    AxisSpec centralityAxis = {centralityBins, "Centrality", "centralityAxis"};
-    AxisSpec occupancyAxis = {occupancyBins, "Occupancy", "occupancyAxis"};
+    const AxisSpec pvAxis = {binOpt.pvBins, "PV", "PV axis"};
+    const AxisSpec multFV0aAxis = {binOpt.fv0aMultBins, "fv0a", "FV0AMult axis"};
+    const AxisSpec multFT0aAxis = {binOpt.ft0aMultBins, "ft0a", "FT0AMult axis"};
+    const AxisSpec multFT0cAxis = {binOpt.ft0cMultBins, "ft0c", "FT0CMult axis"};
+    const AxisSpec centralityAxis = {binOpt.centralityBins, "Centrality", "centrality axis"};
+    const AxisSpec occupancyAxis = {binOpt.occupancyBins, "Occupancy", "occupancy axis"};
+    const AxisSpec irAxis = {binOpt.interactionRateBins, "Interaction Rate", "IR axis"};
+    const AxisSpec ptAxis = {binOpt.PtBins, "Pt axis (GeV/c)"};
+    const AxisSpec multAxis = {binOpt.MultBins, "N_{trk} axis"};
+    const AxisSpec zAxis = {binOpt.ZvtxBins, "Z-vtx axis"};
+    const AxisSpec deltaZAxis = {binOpt.DeltaZBins, "Delta Z-vtx axis"};
+    const AxisSpec dcaxyAxis = {binOpt.DCAxyBins, "DCA-xy axis"};
+    const AxisSpec dcazAxis = {binOpt.DCAzBins, "DCA-z axis"};
+    const AxisSpec phiAxis = {binOpt.PhiBins, "#phi axis"};
+    const AxisSpec etaAxis = {binOpt.EtaBins, "#eta axis"};
+    const AxisSpec chiSqAxis = {binOpt.ChiSqPerNclBins, "Chi2 axis"};
+    const AxisSpec nclsAxis = {binOpt.NclBins, "Number of clusters axis"};
 
     rctChecker.init(rctCuts.cfgEvtRCTFlagCheckerLabel, rctCuts.cfgEvtRCTFlagCheckerZDCCheck, rctCuts.cfgEvtRCTFlagCheckerLimitAcceptAsBad);
 
@@ -261,7 +276,10 @@ struct DndetaMFTPbPb {
     hev->GetXaxis()->SetBinLabel(13, "Above max occup.");
     hev->GetXaxis()->SetBinLabel(14, "RCT Flag Checker");
 
-    registry.add("Tracks/hTrkSel", "Number of tracks; Cut; #Tracks Passed Cut", {HistType::kTH1D, {{nTrkSel, 0, nTrkSel}}});
+    registry.add("Tracks/hTrkSel", "Number of tracks; Cut; #Tracks Passed Cut", {HistType::kTH1F, {{nTrkSel, -0.5, +nTrkSel - 0.5}}});
+    registry.get<TH1>(HIST("Tracks/hTrkSel"))->GetXaxis()->SetBinLabel(trkSelAll + 1, "All");
+    registry.get<TH1>(HIST("Tracks/hTrkSel"))->GetXaxis()->SetBinLabel(trkSelWoAmbiguous + 1, "Ambiguity cut");
+    registry.get<TH1>(HIST("Tracks/hTrkSel"))->GetXaxis()->SetBinLabel(trkSelNumReassoc + 1, "Reassociated");
     registry.get<TH1>(HIST("Tracks/hTrkSel"))->GetXaxis()->SetBinLabel(trkSelNCls + 1, "Ncl cut");
     registry.get<TH1>(HIST("Tracks/hTrkSel"))->GetXaxis()->SetBinLabel(trkSelChi2Ncl + 1, "#chi^{2}/Ncl cut");
     registry.get<TH1>(HIST("Tracks/hTrkSel"))->GetXaxis()->SetBinLabel(trkSelEta + 1, "#eta cut");
@@ -285,6 +303,8 @@ struct DndetaMFTPbPb {
       x->SetBinLabel(1, "All");
       x->SetBinLabel(2, "Selected");
 
+      registry.add("Events/hInteractionRate", "; occupancy; IR (kHz)", kTH2F, {occupancyAxis, irAxis});
+
       registry.add({"Events/NtrkZvtx",
                     "; N_{trk}; Z_{vtx} (cm); occupancy",
                     {HistType::kTHnSparseF, {multAxis, zAxis, occupancyAxis}}});
@@ -375,9 +395,11 @@ struct DndetaMFTPbPb {
       hstat->GetAxis(0)->SetBinLabel(1, "All");
       hstat->GetAxis(0)->SetBinLabel(2, "Selected");
 
+      registry.add("Events/Centrality/hInteractionRate", "; centrality; occupancy; IR (kHz)", kTHnSparseF, {centralityAxis, occupancyAxis, irAxis});
+
       qaregistry.add({"Events/Centrality/hCent",
                       "; centrality; occupancy",
-                      {HistType::kTH2F, {centAxis, occupancyAxis}},
+                      {HistType::kTH2F, {centralityAxis, occupancyAxis}},
                       true});
       qaregistry.add(
         {"Events/Centrality/hZvtxCent",
@@ -775,7 +797,6 @@ struct DndetaMFTPbPb {
   using CollGenCent = CollsGenCentFT0C::iterator;
   using MFTTracksLabeled = soa::Join<aod::MFTTracks, aod::McMFTTrackLabels>;
   using MftTracksWColls = soa::Join<aod::MFTTracks, aod::MFTTrkCompColls>;
-
   using CollsCorr = soa::Join<aod::Collisions, aod::EvSels, aod::Mults, aod::PVMults, aod::CentFT0Cs, aod::CentFV0As, aod::CentFT0CVariant1s, aod::CentFT0Ms, aod::CentNGlobals, aod::CentMFTs>;
 
   /// Filtered tables
@@ -888,7 +909,15 @@ struct DndetaMFTPbPb {
     auto nATrk = 0;
     if (besttracks.size() > 0) {
       for (auto const& atrack : besttracks) {
+        registry.fill(HIST("Tracks/hTrkSel"), trkSelAll);
+        if (trackCuts.excludeAmbiguous && atrack.ambDegree() > 1) {
+          continue;
+        }
+        registry.fill(HIST("Tracks/hTrkSel"), trkSelWoAmbiguous);
         auto itrack = atrack.template mfttrack_as<T>();
+        if (itrack.collisionId() != atrack.bestCollisionId()) {
+          registry.fill(HIST("Tracks/hTrkSel"), trkSelNumReassoc);
+        }
         if (!isTrackSelected(itrack)) {
           continue;
         }
@@ -985,8 +1014,8 @@ struct DndetaMFTPbPb {
     if (gHadronicRate.find(mRunNumber) == gHadronicRate.end()) {
       auto runDuration = ccdb->getRunDuration(mRunNumber);
       mSOR = runDuration.first;
-      mMinSeconds = std::floor(mSOR * 1.e-3);                /// round tsSOR to the highest integer lower than tsSOR
-      double maxSec = std::ceil(runDuration.second * 1.e-3); /// round tsEOR to the lowest integer higher than tsEOR
+      mMinSeconds = std::floor(mSOR * 1.e-3);               /// round tsSOR to the highest integer lower than tsSOR
+      float maxSec = std::ceil(runDuration.second * 1.e-3); /// round tsEOR to the lowest integer higher than tsEOR
       const AxisSpec axisSeconds{static_cast<int>((maxSec - mMinSeconds) / 20.f), 0, maxSec - mMinSeconds, "Seconds since SOR"};
       int hadronicRateBins = static_cast<int>(eventCuts.maxIR - eventCuts.minIR);
       gHadronicRate[mRunNumber] = registry.add<TH2>(Form("HadronicRate/%i", mRunNumber), ";Time since SOR (s);Hadronic rate (kHz)", kTH2D, {axisSeconds, {hadronicRateBins, eventCuts.minIR, eventCuts.maxIR}}).get();
@@ -1195,8 +1224,13 @@ struct DndetaMFTPbPb {
 
     if (cfgDoIR) {
       initHadronicRate(bc);
-      double ir = rateFetcher.fetch(ccdb.service, bc.timestamp(), bc.runNumber(), cfgIRSource, cfgIRCrashOnNull) * 1.e-3;
-      double seconds = bc.timestamp() * 1.e-3 - mMinSeconds;
+      float ir = !cfgIRSource.value.empty() ? rateFetcher.fetch(ccdb.service, bc.timestamp(), bc.runNumber(), cfgIRSource, cfgIRCrashOnNull) * 1.e-3 : -1;
+      if constexpr (has_reco_cent<C>) {
+        registry.fill(HIST("Events/Centrality/hInteractionRate"), c, occ, ir);
+      } else {
+        registry.fill(HIST("Events/hInteractionRate"), occ, ir);
+      }
+      float seconds = bc.timestamp() * 1.e-3 - mMinSeconds;
       if (cfgUseIRCut && (ir < eventCuts.minIR || ir > eventCuts.maxIR)) { // cut on hadronic rate
         return;
       }
@@ -1206,8 +1240,6 @@ struct DndetaMFTPbPb {
     auto z = collision.posZ();
     if constexpr (has_reco_cent<C>) {
       registry.fill(HIST("Events/Centrality/Selection"), 2., c, occ);
-      qaregistry.fill(HIST("Events/Centrality/hZvtxCent"), z, c, occ);
-      qaregistry.fill(HIST("Events/Centrality/hCent"), c, occ);
     } else {
       registry.fill(HIST("Events/Selection"), 2., occ);
     }
@@ -1237,14 +1269,19 @@ struct DndetaMFTPbPb {
       registry.fill(HIST("Events/Selection"), 1., occ);
     }
 
-    if (!isGoodEvent<false>(collision)) {
+    if (!isGoodEvent<true>(collision)) {
       return;
     }
 
     if (cfgDoIR) {
       initHadronicRate(bc);
-      double ir = rateFetcher.fetch(ccdb.service, bc.timestamp(), bc.runNumber(), cfgIRSource, cfgIRCrashOnNull) * 1.e-3;
-      double seconds = bc.timestamp() * 1.e-3 - mMinSeconds;
+      float ir = !cfgIRSource.value.empty() ? rateFetcher.fetch(ccdb.service, bc.timestamp(), bc.runNumber(), cfgIRSource, cfgIRCrashOnNull) * 1.e-3 : -1;
+      if constexpr (has_reco_cent<C>) {
+        registry.fill(HIST("Events/Centrality/hInteractionRate"), c, occ, ir);
+      } else {
+        registry.fill(HIST("Events/hInteractionRate"), occ, ir);
+      }
+      float seconds = bc.timestamp() * 1.e-3 - mMinSeconds;
       if (cfgUseIRCut && (ir < eventCuts.minIR || ir > eventCuts.maxIR)) { // cut on hadronic rate
         return;
       }
@@ -1254,6 +1291,8 @@ struct DndetaMFTPbPb {
     auto z = collision.posZ();
     if constexpr (has_reco_cent<C>) {
       registry.fill(HIST("Events/Centrality/Selection"), 2., c, occ);
+      qaregistry.fill(HIST("Events/Centrality/hZvtxCent"), z, c, occ);
+      qaregistry.fill(HIST("Events/Centrality/hCent"), c, occ);
     } else {
       registry.fill(HIST("Events/Selection"), 2., occ);
     }
@@ -1960,6 +1999,9 @@ struct DndetaMFTPbPb {
     }
 
     for (auto const& track : besttracks) {
+      if (trackCuts.excludeAmbiguous && track.ambDegree() > 1) {
+        continue;
+      }
       auto itrack = track.mfttrack_as<FiltMcMftTracks>();
       if (!isTrackSelected(itrack)) {
         continue;
@@ -2256,6 +2298,9 @@ struct DndetaMFTPbPb {
 
     auto nBestTrks = 0;
     for (auto const& atrack : besttracks) {
+      if (trackCuts.excludeAmbiguous && atrack.ambDegree() > 1) {
+        continue;
+      }
       auto itrack = atrack.template mfttrack_as<FiltMftTracks>();
       if (itrack.eta() < trackCuts.minEta || itrack.eta() > trackCuts.maxEta) {
         continue;