[PWGJE] additional outlier checks on tracks (#13747)

Author: jaimenorman

PWGJE/Tasks/jetOutlierQA.cxx

@@ -84,6 +84,7 @@ struct JetOutlierQATask {
   Configurable<unsigned int> splitCollisionsDeltaBC{"splitCollisionsDeltaBC", 5, "threshold in BC to assign as split collision"};
   Configurable<int> mergeCollisionsDeltaMin{"mergeCollisionsDeltaMin", -10, "number of prior collisions to search for close Z position"};
   Configurable<int> mergeCollisionsDeltaMax{"mergeCollisionsDeltaMax", 10, "number of following collisions to search for close Z position"};
+  Configurable<int> maxNTracksJJdifferent{"maxNTracksJJdifferent", 10, "maximum number of tracks from different JJ collision to be considered for track rejection"};
 
   std::map<uint64_t, std::vector<int64_t>> fBCCollMap; // key: global BC, value: vector of reduced event global indices
 
@@ -233,6 +234,15 @@ struct JetOutlierQATask {
       registry.add("h2_outlier_collision_ID_difference", "difference in collision ID between outlier collision and analysed collision", {HistType::kTH2F, {{600, 0, 600}, {200, -100, 100}}});
       registry.add("h_DeltaZ_Outlier", "Delta Z between outlier collision and analysed collision", {HistType::kTH1F, {{1200, -30, 30}}});
       registry.add("h2_DeltaZ_Outlier_difference", "Delta Z between outlier collision and analysed collision vs difference in collision ID", {HistType::kTH2F, {{1200, -30, 30}, {200, -100, 100}}});
+
+      registry.add("h_track_pt_same_collision", "track pt from same collision or different MB collision;p_{T,track} (GeV/#it{c});entries", {HistType::kTH1F, {{200, 0, 200}}});
+      registry.add("h_track_pt_eta_same_collision", "track pt vs eta from same collision or different MB collision;p_{T,track} (GeV/#it{c});#eta_{track};entries", {HistType::kTH2F, {{200, 0, 200}, {100, -5, 5}}});
+      registry.add("h_track_pt_phi_same_collision", "track pt vs phi from same collision or different MB collision;p_{T,track} (GeV/#it{c});#varphi_{track} (rad);entries", {HistType::kTH2F, {{200, 0, 200}, {160, -1.0, 7.0}}});
+      registry.add("h2_collision_ID_difference_same_collision", "difference in collision ID between outlier collision and analysed collision", {HistType::kTH2F, {{600, 0, 600}, {200, -100, 100}}});
+      registry.add("h_track_pt_no_JJ_different", "track pt from same collision or different MB collision;p_{T,track} (GeV/#it{c});entries", {HistType::kTH1F, {{200, 0, 200}}});
+      registry.add("h_track_pt_eta_no_JJ_different", "track pt vs eta from same collision or different MB collision;p_{T,track} (GeV/#it{c});#eta_{track};entries", {HistType::kTH2F, {{200, 0, 200}, {100, -5, 5}}});
+      registry.add("h_track_pt_phi_no_JJ_different", "track pt vs phi from same collision or different MB collision;p_{T,track} (GeV/#it{c});#varphi_{track} (rad);entries", {HistType::kTH2F, {{200, 0, 200}, {160, -1.0, 7.0}}});
+      registry.add("h2_collision_ID_difference_no_JJ_different", "difference in collision ID between outlier collision and analysed collision", {HistType::kTH2F, {{600, 0, 600}, {200, -100, 100}}});
     }
   }
 
@@ -567,6 +577,25 @@ struct JetOutlierQATask {
         }
         // check outlier tracks and neighbouring collisions
         registry.fill(HIST("h2_pt_hat_track_pt"), pTHat, track.pt());
+        // get MC info about track and collision
+        auto mcParticleOutlier = track.mcParticle_as<JetParticlesWithOriginal>();
+        auto collisionMCOutlier = collisionsMC.sliceBy(perColParticle, mcParticleOutlier.mcCollisionId());
+        if (collisionMCOutlier.size() != 1) {
+          LOG(info) << "size of collision outlier not expected";
+          return;
+        }
+        int mcCollisionIDcoll = collision.mcCollisionId(); // Get the corresponding MC collision ID from the reco collision
+        int mcCollisionIDOutlier = mcParticleOutlier.mcCollisionId();
+        int subGenIDOutlier = collisionMCOutlier.begin().subGeneratorId();
+        int outlierCollisionIDDifference = mcCollisionIDOutlier - mcCollisionIDcoll;
+
+        int nMBdifferent = 0;
+        int nMBdifferentSelected = 0;
+        int nJJdifferent = 0;
+        int nJJdifferentSelected = 0;
+        int nJJsame = 0;
+        int nJJsameSelected = 0;
+        // ID outlier based on track pT relative to pTHat
         if (track.pt() > pTHatMaxMCDOutlier * pTHat) { // high weight outlier track
           registry.fill(HIST("h_track_pt_outlier"), track.pt());
           registry.fill(HIST("h2_pt_hat_track_pt_outlier"), pTHat, track.pt());
@@ -581,13 +610,6 @@ struct JetOutlierQATask {
           }
           // now match tracks to their MC particle, check the MC collision ID of this particle, and
           // check what fraction of tracks in this event are associated to this MC collision
-          int nMBdifferent = 0;
-          int nMBdifferentSelected = 0;
-          int nJJdifferent = 0;
-          int nJJdifferentSelected = 0;
-          int nJJsame = 0;
-          int nJJsameSelected = 0;
-          int mcCollisionIDcoll = collision.mcCollisionId(); // Get the corresponding MC collision ID from the reco collision
           // LOG(info) << "--- Loop over tracks in outlier event with pT/pThat = " << track.pt() / pTHat << "---";
           // LOG(info) << "N tracks in outlier event = " << tracksColl.size() << " pTHat = " << pTHat << " collisionID = " << collision.globalIndex() << " mcCollisionID = " << collision.mcCollisionId();
           for (auto const& trackOutlier : tracksColl) {
@@ -639,15 +661,6 @@ struct JetOutlierQATask {
           registry.fill(HIST("h2_outlier_event_tracks_frac_different_selected_JJ"), pTHat, float(nJJdifferentSelected) / float(nJJdifferentSelected + nJJsameSelected + nMBdifferentSelected));
           registry.fill(HIST("h2_outlier_event_tracks_frac_different_selected_MB"), pTHat, float(nMBdifferentSelected) / float(nJJdifferentSelected + nJJsameSelected + nMBdifferentSelected));
           // now check where outlier comes from
-          auto mcParticleOutlier = track.mcParticle_as<JetParticlesWithOriginal>();
-          auto collisionMCOutlier = collisionsMC.sliceBy(perColParticle, mcParticleOutlier.mcCollisionId());
-          if (collisionMCOutlier.size() != 1) {
-            LOG(info) << "size of collision outlier not expected";
-            return;
-          }
-          int mcCollisionIDOutlier = mcParticleOutlier.mcCollisionId();
-          int subGenIDOutlier = collisionMCOutlier.begin().subGeneratorId();
-          int outlierCollisionIDDifference = mcCollisionIDOutlier - mcCollisionIDcoll;
           // LOG(info) <<"outlier comes from " << (mcCollisionIDOutlier == mcCollisionIDcoll ? "same" : "different") << " event which is a " << (subGenIDOutlier == jetderiveddatautilities::JCollisionSubGeneratorId::mbGap ? " MB-gap" : " jet-jet") << " collision with mcCollisionID = " << mcCollisionIDOutlier;
           registry.fill(HIST("h2_outlier_collision_ID_difference"), pTHat, float(outlierCollisionIDDifference));
           // if outlier comes from different collision, check which type and set flags
@@ -660,7 +673,23 @@ struct JetOutlierQATask {
             isOutlierEventDifferentMBCollision = true;
           }
         }
-        // all
+        // fill for tracks from same collision or different MB collision in collisions that likely aren't fully merged
+        if (nJJdifferentSelected < maxNTracksJJdifferent &&
+            (subGenIDOutlier == jetderiveddatautilities::JCollisionSubGeneratorId::mbGap ||
+             mcCollisionIDOutlier == mcCollisionIDcoll)) {
+          registry.fill(HIST("h_track_pt_same_collision"), track.pt(), weight);
+          registry.fill(HIST("h_track_pt_eta_same_collision"), track.pt(), track.eta(), weight);
+          registry.fill(HIST("h_track_pt_phi_same_collision"), track.pt(), track.phi(), weight);
+          registry.fill(HIST("h2_collision_ID_difference_same_collision"), pTHat, float(outlierCollisionIDDifference));
+        }
+        // fill tracks for events which have no JJ outlier tracks from different events
+        if (nJJdifferentSelected == 0) {
+          registry.fill(HIST("h_track_pt_no_JJ_different"), track.pt(), weight);
+          registry.fill(HIST("h_track_pt_eta_no_JJ_different"), track.pt(), track.eta(), weight);
+          registry.fill(HIST("h_track_pt_phi_no_JJ_different"), track.pt(), track.phi(), weight);
+          registry.fill(HIST("h2_collision_ID_difference_no_JJ_different"), pTHat, float(outlierCollisionIDDifference));
+        }
+        // collision checks for all tracks
         for (auto const& collisionOutlier : collisions) { // find collisions closeby
           float eventWeightOutlier = collisionOutlier.mcCollision().weight();
           double pTHatOutlier = collisionOutlier.mcCollision().ptHard();