diff --git a/PWGJE/TableProducer/mcOutlierRejector.cxx b/PWGJE/TableProducer/mcOutlierRejector.cxx index 48241281565..236ca9f5a3a 100644 --- a/PWGJE/TableProducer/mcOutlierRejector.cxx +++ b/PWGJE/TableProducer/mcOutlierRejector.cxx @@ -33,6 +33,7 @@ using namespace o2::framework::expressions; struct McOutlierRejectorTask { Produces collisionOutliers; Produces mcCollisionOutliers; + Preslice> perColParticle = aod::jmccollision::mcCollisionId; Configurable checkmcCollisionForCollision{"checkmcCollisionForCollision", true, "additionally reject collision based on mcCollision"}; Configurable ptHatMax{"ptHatMax", 4.0, "maximum factor of pt hat the leading jet in the event is allowed"}; @@ -55,22 +56,36 @@ struct McOutlierRejectorTask { PROCESS_SWITCH(McOutlierRejectorTask, processSetupMcCollisionSelection, "Setup MC Collision processing", true); template - void collisionSelection(int32_t collisionIndex, T const& selectionObjects, float ptHard, std::vector& flagArray) + void collisionSelection(int32_t collisionIndex, int32_t mcCollisionId, T const& selectionObjects, float ptHard, std::vector& flagArray, std::optional>> mcCollisionsOpt = std::nullopt) { - if (selectionObjects.size() != 0) { + bool isTrueOutlier = true; float maxSelectionObjectPt = 0.0; - if constexpr (std::is_same_v, aod::JetTracks> || std::is_same_v, aod::JetParticles>) { + if constexpr (std::is_same_v, aod::JetTracksMCD> || std::is_same_v, aod::JetParticles>) { for (auto selectionObject : selectionObjects) { if (selectionObject.pt() > maxSelectionObjectPt) { maxSelectionObjectPt = selectionObject.pt(); + // may be slow - could save only MC particle then check difference only for tracks IDd as outliers? + if constexpr (std::is_same_v, aod::JetTracksMCD>) { + auto& mcCollisions = mcCollisionsOpt.value().get(); + auto mcParticle = selectionObject.template mcParticle_as>(); + auto mcCollision = mcCollisions.sliceBy(perColParticle, mcParticle.mcCollisionId()); + int subGenID = mcCollision.begin().subGeneratorId(); + int diffCollisionID = mcParticle.mcCollisionId() - mcCollisionId; + if (subGenID == jetderiveddatautilities::JCollisionSubGeneratorId::mbGap && + diffCollisionID != 0) { + isTrueOutlier = true; + } else { + isTrueOutlier = false; + } + } } } } else { maxSelectionObjectPt = selectionObjects.iteratorAt(0).pt(); } - if (maxSelectionObjectPt > ptHatMax * ptHard) { + if (maxSelectionObjectPt > ptHatMax * ptHard && isTrueOutlier) { flagArray[collisionIndex] = true; // Currently if running multiple different jet finders, then a single type of jet can veto an event for others. Decide if this is the best way } } @@ -80,13 +95,20 @@ struct McOutlierRejectorTask { void processSelectionObjects(aod::JetCollisionMCD const& collision, T const& selectionObjects, aod::JetMcCollisions const&) { auto mcCollision = collision.mcCollision_as(); - collisionSelection(collision.globalIndex(), selectionObjects, mcCollision.ptHard(), collisionFlag); + collisionSelection(collision.globalIndex(), 1., selectionObjects, mcCollision.ptHard(), collisionFlag); + } + + template + void processSelectionObjectsTracks(aod::JetCollisionMCD const& collision, T const& selectionObjects, soa::Join const& mcCollisions, soa::Join const&) + { + auto mcCollision = collision.mcCollision_as>(); + collisionSelection(collision.globalIndex(), collision.mcCollisionId(), selectionObjects, mcCollision.ptHard(), collisionFlag, mcCollisions); } template void processSelectionMcObjects(aod::JetMcCollision const& mcCollision, T const& selectionMcObjects) { - collisionSelection(mcCollision.globalIndex(), selectionMcObjects, mcCollision.ptHard(), mcCollisionFlag); + collisionSelection(mcCollision.globalIndex(), 1., selectionMcObjects, mcCollision.ptHard(), mcCollisionFlag); } PROCESS_SWITCH_FULL(McOutlierRejectorTask, processSelectionObjects, processSelectingChargedMCDetectorLevelJets, "process mc detector level charged jets", true); @@ -101,7 +123,7 @@ struct McOutlierRejectorTask { PROCESS_SWITCH_FULL(McOutlierRejectorTask, processSelectionObjects, processSelectingBplusChargedMCDetectorLevelJets, "process mc detector level Bplus charged jets", false); PROCESS_SWITCH_FULL(McOutlierRejectorTask, processSelectionObjects, processSelectingXicToXiPiPiChargedMCDetectorLevelJets, "process mc detector level XicToXiPiPi charged jets", false); PROCESS_SWITCH_FULL(McOutlierRejectorTask, processSelectionObjects, processSelectingDielectronChargedMCDetectorLevelJets, "process mc detector level Dielectron charged jets", false); - PROCESS_SWITCH_FULL(McOutlierRejectorTask, processSelectionObjects, processSelectingTracks, "process tracks", false); + PROCESS_SWITCH_FULL(McOutlierRejectorTask, processSelectionObjectsTracks, processSelectingTracks, "process tracks", false); PROCESS_SWITCH_FULL(McOutlierRejectorTask, processSelectionMcObjects, processSelectingChargedMCParticleLevelJets, "process mc particle level charged jets", true); PROCESS_SWITCH_FULL(McOutlierRejectorTask, processSelectionMcObjects, processSelectingNeutralMCParticleLevelJets, "process mc particle level neutral jets", false); PROCESS_SWITCH_FULL(McOutlierRejectorTask, processSelectionMcObjects, processSelectingFullMCParticleLevelJets, "process mc particle level full jets", false); diff --git a/PWGJE/Tasks/jetFinderQA.cxx b/PWGJE/Tasks/jetFinderQA.cxx index 4273efd027c..b0c4fad895b 100644 --- a/PWGJE/Tasks/jetFinderQA.cxx +++ b/PWGJE/Tasks/jetFinderQA.cxx @@ -329,14 +329,14 @@ struct JetFinderQATask { } if (doprocessTracks || doprocessTracksWeighted) { - registry.add("h_collisions", "event status;event status;entries", {HistType::kTH1F, {{4, 0.0, 4.0}}}); + registry.add("h_collisions", "event status;event status;entries", {HistType::kTH1F, {{5, 0.0, 5.0}}}); registry.add("h2_centrality_collisions", "centrality vs collisions; centrality; collisions", {HistType::kTH2F, {{1200, -10.0, 110.0}, {4, 0.0, 4.0}}}); registry.add("h3_centrality_track_pt_track_phi", "centrality vs track pT vs track #varphi; centrality; #it{p}_{T,track} (GeV/#it{c}); #varphi_{track}", {HistType::kTH3F, {{1200, -10.0, 110.0}, {200, 0., 200.}, {160, -1.0, 7.}}}); registry.add("h3_centrality_track_pt_track_eta", "centrality vs track pT vs track #eta; centrality; #it{p}_{T,track} (GeV/#it{c}); #eta_{track}", {HistType::kTH3F, {{1200, -10.0, 110.0}, {200, 0., 200.}, trackEtaAxis}}); registry.add("h3_centrality_track_pt_track_dcaxy", "centrality vs track pT vs track DCA_{xy}; centrality; #it{p}_{T,track} (GeV/#it{c}); track DCA_{xy}", {HistType::kTH3F, {{120, -10.0, 110.0}, {20, 0., 100.}, {200, -0.15, 0.15}}}); registry.add("h3_track_pt_track_eta_track_phi", "track pT vs track #eta vs track #varphi; #it{p}_{T,track} (GeV/#it{c}); #eta_{track}; #varphi_{track}", {HistType::kTH3F, {{200, 0., 200.}, trackEtaAxis, {160, -1.0, 7.}}}); if (doprocessTracksWeighted) { - registry.add("h_collisions_weighted", "event status;event status;entries", {HistType::kTH1F, {{4, 0.0, 4.0}}}); + registry.add("h_collisions_weighted", "event status;event status;entries", {HistType::kTH1F, {{5, 0.0, 5.0}}}); } } if (doprocessTracksSub) { @@ -904,11 +904,14 @@ struct JetFinderQATask { } PROCESS_SWITCH(JetFinderQATask, processJetsMCD, "jet finder QA mcd", false); - void processJetsMCDWeighted(soa::Filtered::iterator const& collision, soa::Join const& jets, aod::JetTracks const&) + void processJetsMCDWeighted(soa::Filtered>::iterator const& collision, soa::Join const& jets, aod::JetTracks const&) { if (collision.trackOccupancyInTimeRange() < trackOccupancyInTimeRangeMin || trackOccupancyInTimeRangeMax < collision.trackOccupancyInTimeRange()) { return; } + if (collision.isOutlier()) { + return; + } for (auto const& jet : jets) { if (!jetfindingutilities::isInEtaAcceptance(jet, jetEtaMin, jetEtaMax, trackEtaMin, trackEtaMax)) { continue; @@ -946,7 +949,7 @@ struct JetFinderQATask { } PROCESS_SWITCH(JetFinderQATask, processJetsMCP, "jet finder QA mcp", false); - void processJetsMCPWeighted(soa::Join::iterator const& jet, aod::JetParticles const&, aod::JetMcCollisions const&, soa::Filtered const& collisions) + void processJetsMCPWeighted(soa::Join::iterator const& jet, aod::JetParticles const&, aod::JetMcCollisions const&, soa::Filtered> const& collisions) { if (!jetfindingutilities::isInEtaAcceptance(jet, jetEtaMin, jetEtaMax, trackEtaMin, trackEtaMax)) { return; @@ -962,7 +965,7 @@ struct JetFinderQATask { } if (checkMcCollisionIsMatched) { auto collisionspermcpjet = collisions.sliceBy(CollisionsPerMCPCollision, jet.mcCollisionId()); - if (collisionspermcpjet.size() >= 1 && jetderiveddatautilities::selectCollision(collisionspermcpjet.begin(), eventSelectionBits)) { + if (collisionspermcpjet.size() >= 1 && jetderiveddatautilities::selectCollision(collisionspermcpjet.begin(), eventSelectionBits) && !collisionspermcpjet.begin().isOutlier()) { fillMCPHistograms(jet, jet.eventWeight()); } } else { @@ -997,7 +1000,7 @@ struct JetFinderQATask { } PROCESS_SWITCH(JetFinderQATask, processJetsMCPMCDMatched, "jet finder QA matched mcp and mcd", false); - void processJetsMCPMCDMatchedWeighted(soa::Filtered::iterator const& collision, + void processJetsMCPMCDMatchedWeighted(soa::Filtered>::iterator const& collision, soa::Join const& mcdjets, soa::Join const&, aod::JetTracks const&, aod::JetParticles const&) @@ -1005,6 +1008,9 @@ struct JetFinderQATask { if (collision.trackOccupancyInTimeRange() < trackOccupancyInTimeRangeMin || trackOccupancyInTimeRangeMax < collision.trackOccupancyInTimeRange()) { return; } + if (collision.isOutlier()) { + return; + } for (const auto& mcdjet : mcdjets) { if (!jetfindingutilities::isInEtaAcceptance(mcdjet, jetEtaMin, jetEtaMax, trackEtaMin, trackEtaMax)) { continue; @@ -1188,7 +1194,7 @@ struct JetFinderQATask { } PROCESS_SWITCH(JetFinderQATask, processTracks, "QA for charged tracks", false); - void processTracksWeighted(soa::Join::iterator const& collision, + void processTracksWeighted(soa::Join::iterator const& collision, aod::JetMcCollisions const&, soa::Filtered> const& tracks) { @@ -1208,6 +1214,11 @@ struct JetFinderQATask { } registry.fill(HIST("h_collisions"), 2.5); registry.fill(HIST("h_collisions_weighted"), 2.5, eventWeight); + if (collision.isOutlier()) { + return; + } + registry.fill(HIST("h_collisions"), 3.5); + registry.fill(HIST("h_collisions_weighted"), 3.5, eventWeight); fillTrackHistograms(collision, tracks, eventWeight); } PROCESS_SWITCH(JetFinderQATask, processTracksWeighted, "QA for charged tracks weighted", false);