[PWGDQ] dqEfficiency/dilepton-track: Configurable to include only correctly associated MC

PWGDQ/Tasks/dqEfficiency_withAssoc.cxx

@@ -47,7 +47,7 @@
 #include <TString.h>
 
 #include <algorithm>
 #include <iostream>
 #include <map>
 #include <memory>
 #include <string>
@@ -253,7 +253,7 @@
 void PrintBitMap(TMap map, int nbits)
 {
   for (int i = 0; i < nbits; i++) {
     cout << ((map & (TMap(1) << i)) > 0 ? "1" : "0");
   }
 }
 
@@ -3592,6 +3592,7 @@
   Configurable<float> fConfigMCGenDileptonLegEtaAbs{"cfgMCGenDileptonLegEtaAbs", 0.9f, "eta abs range for the dilepton leg"};
   Configurable<float> fConfigMCGenHadronEtaAbs{"cfgMCGenHadronEtaAbs", 0.9f, "eta abs range for the hadron"};
   Configurable<bool> fConfigUseMCRapcut{"cfgUseMCRapcut", false, "Use Rap cut for dileptons used in the triplet vertexing(reconstructed)"};
+  Configurable<bool> fConfigUseCorrectlyAssociatedMC{"cfgUseCorrectlyAssociatedMC", true, "Use only correctly associated track to construct triplet"};
   Configurable<int> fConfigMixingDepth{"cfgMixingDepth", 5, "Event mixing pool depth"};
 
   int fCurrentRun; // needed to detect if the run changed and trigger update of calibrations etc.
@@ -3982,6 +3983,9 @@
 
     uint32_t mcDecision = static_cast<uint32_t>(0);
     size_t isig = 0;
+    bool isCorrectAssoc_lepton1 = false;
+    bool isCorrectAssoc_lepton2 = false;
+    bool isCorrectAssoc_assoc = false;
 
     for (auto dilepton : dileptons) {
       // get full track info of tracks based on the index
@@ -3998,9 +4002,15 @@
       }
       // dilepton rap cut
       float rap = dilepton.rap();
       if (fConfigUseMCRapcut && abs(rap) > fConfigDileptonRapCutAbs)
         continue;
 
+      // Correct association
+      if (fConfigUseCorrectlyAssociatedMC) {
+        isCorrectAssoc_lepton1 = (lepton1MC.reducedMCevent() == event.reducedMCevent());
+        isCorrectAssoc_lepton2 = (lepton2MC.reducedMCevent() == event.reducedMCevent());
+      }
+
       VarManager::FillTrack<fgDileptonFillMap>(dilepton, fValuesDilepton);
 
       // fill selected dilepton histograms for each specified selection
@@ -4051,11 +4061,22 @@
           if (track.globalIndex() == dilepton.index0Id() || track.globalIndex() == dilepton.index1Id()) {
             continue;
           }
+
+          auto trackMC = track.reducedMCTrack();
+          if (fConfigUseCorrectlyAssociatedMC) { // check if correctly associated
+            isCorrectAssoc_assoc = (trackMC.reducedMCevent() == event.reducedMCevent());
+          }
+          bool isCorrect_triplet = isCorrectAssoc_lepton1 && isCorrectAssoc_lepton2 && isCorrectAssoc_assoc;
+          if (fConfigUseCorrectlyAssociatedMC) {
+            if (!isCorrect_triplet) {
+              continue;
+            }
+          }
+
           // compute needed quantities
           VarManager::FillDileptonHadron(dilepton, track, fValuesHadron);
           VarManager::FillDileptonTrackVertexing<TCandidateType, TEventFillMap, TTrackFillMap>(event, lepton1, lepton2, track, fValuesHadron);
 
-          auto trackMC = track.reducedMCTrack();
           mcDecision = 0;
           isig = 0;
           for (auto sig = fRecMCSignals.begin(); sig != fRecMCSignals.end(); sig++, isig++) {
@@ -4306,12 +4327,12 @@
           mctrack.mcReducedFlags() > 0) {
         /*cout << ">>>>>>>>>>>>>>>>>>>>>>> track idx / pdg / selections: " << mctrack.globalIndex() << " / " << mctrack.pdgCode() << " / ";
         PrintBitMap(mctrack.mcReducedFlags(), 16);
         cout << endl;
         if (mctrack.has_mothers()) {
           for (auto& m : mctrack.mothersIds()) {
             if (m < mcTracks.size()) { // protect against bad mother indices
               auto aMother = mcTracks.rawIteratorAt(m);
               cout << "<<<<<< mother idx / pdg: " << m << " / " << aMother.pdgCode() << endl;
             }
           }
         }
@@ -4320,7 +4341,7 @@
           for (int d = mctrack.daughtersIds()[0]; d <= mctrack.daughtersIds()[1]; ++d) {
             if (d < mcTracks.size()) { // protect against bad daughter indices
               auto aDaughter = mcTracks.rawIteratorAt(d);
               cout << "<<<<<< daughter idx / pdg: " << d << " / " << aDaughter.pdgCode() << endl;
             }
           }
         }*/
@@ -4422,12 +4443,12 @@
       // apply kinematic cuts for signal
       if ((t1_raw.pt() < fConfigDileptonLowpTCut || t1_raw.pt() > fConfigDileptonHighpTCut))
         continue;
       if (abs(t1_raw.y()) > fConfigDileptonRapCutAbs)
         continue;
       // for the energy correlators
       for (auto& t2 : groupedMCTracks) {
         auto t2_raw = groupedMCTracks.rawIteratorAt(t2.globalIndex());
         if (TMath::Abs(t2_raw.pdgCode()) == 443 || TMath::Abs(t2_raw.pdgCode()) == 11 || TMath::Abs(t2_raw.pdgCode()) == 22)
           continue;
         if (t2_raw.pt() < fConfigMCGenHadronPtMin.value || std::abs(t2_raw.eta()) > fConfigMCGenHadronEtaAbs.value)
           continue;
@@ -4490,12 +4511,12 @@
       // apply kinematic cuts for signal
       if ((t1_raw.pt() < fConfigDileptonLowpTCut || t1_raw.pt() > fConfigDileptonHighpTCut))
         continue;
       if (abs(t1_raw.y()) > fConfigDileptonRapCutAbs)
         continue;
       // for the energy correlators
       for (auto& t2 : groupedMCTracks2) {
         auto t2_raw = groupedMCTracks2.rawIteratorAt(t2.globalIndex());
         if (TMath::Abs(t2_raw.pdgCode()) == 443 || TMath::Abs(t2_raw.pdgCode()) == 11 || TMath::Abs(t2_raw.pdgCode()) == 22)
           continue;
         if (t2_raw.pt() < fConfigMCGenHadronPtMin.value || std::abs(t2_raw.eta()) > fConfigMCGenHadronEtaAbs.value)
           continue;