diff --git a/PWGJE/Tasks/jetTaggerHFQA.cxx b/PWGJE/Tasks/jetTaggerHFQA.cxx index b2d1261e335..2f51f5e787f 100644 --- a/PWGJE/Tasks/jetTaggerHFQA.cxx +++ b/PWGJE/Tasks/jetTaggerHFQA.cxx @@ -59,6 +59,7 @@ struct JetTaggerHFQA { Configurable trackPtMax{"trackPtMax", 100.0, "maximum pT acceptance for tracks"}; Configurable trackDcaXYMax{"trackDcaXYMax", 1, "minimum DCA xy acceptance for tracks [cm]"}; Configurable trackDcaZMax{"trackDcaZMax", 2, "minimum DCA z acceptance for tracks [cm]"}; + Configurable maxDeltaR{"maxDeltaR", 0.25, "maximum distance of jet axis from flavour initiating parton"}; Configurable jetEtaMin{"jetEtaMin", -99.0, "minimum jet pseudorapidity"}; Configurable jetEtaMax{"jetEtaMax", 99.0, "maximum jet pseudorapidity"}; Configurable prongChi2PCAMin{"prongChi2PCAMin", 1, "minimum Chi2 PCA of decay length of prongs"}; @@ -155,6 +156,18 @@ struct JetTaggerHFQA { registry.add("h_impact_parameter_xyz", "", {HistType::kTH1F, {{axisImpactParameterXYZ}}}); registry.add("h_impact_parameter_xyz_significance", "", {HistType::kTH1F, {{axisImpactParameterXYZSignificance}}}); } + if (doprocessValFlavourDefMCD) { + registry.add("h2_flavour_dist_quark_flavour_dist_hadron", "", {HistType::kTH2F, {{axisJetFlavour}, {axisJetFlavour}}}); + registry.add("h2_flavour_const_quark_flavour_const_hadron", "", {HistType::kTH2F, {{axisJetFlavour}, {axisJetFlavour}}}); + registry.add("h2_flavour_const_hadron_flavour_dist_hadron", "", {HistType::kTH2F, {{axisJetFlavour}, {axisJetFlavour}}}); + registry.add("h2_flavour_const_quark_flavour_dist_quark", "", {HistType::kTH2F, {{axisJetFlavour}, {axisJetFlavour}}}); + } + if (doprocessValFlavourDefMCP) { + registry.add("h2_part_flavour_dist_quark_part_flavour_dist_hadron", "", {HistType::kTH2F, {{axisJetFlavour}, {axisJetFlavour}}}); + registry.add("h2_part_flavour_const_quark_part_flavour_const_hadron", "", {HistType::kTH2F, {{axisJetFlavour}, {axisJetFlavour}}}); + registry.add("h2_part_flavour_const_hadron_part_flavour_dist_hadron", "", {HistType::kTH2F, {{axisJetFlavour}, {axisJetFlavour}}}); + registry.add("h2_part_flavour_const_quark_part_flavour_dist_quark", "", {HistType::kTH2F, {{axisJetFlavour}, {axisJetFlavour}}}); + } if (doprocessIPsData) { registry.add("h_jet_pt", "", {HistType::kTH1F, {{axisJetPt}}}); registry.add("h_jet_eta", "", {HistType::kTH1F, {{axisEta}}}); @@ -443,6 +456,47 @@ struct JetTaggerHFQA { return true; } + template + void fillValidationFlavourDefMCD(T const& mcdjet, V const& tracks, W const& particles, X const& particlesPerColl, float eventWeight = 1.0) + { + float pTHat = 10. / (std::pow(eventWeight, 1.0 / pTHatExponent)); + if (mcdjet.pt() > pTHatMaxMCD * pTHat) { + return; + } + typename V::iterator hftrack; + int jetflavourConstQuark = jettaggingutilities::mcdJetFromHFShower(mcdjet, tracks, particles, maxDeltaR, true); + int jetflavourConstHadron = jettaggingutilities::mcdJetFromHFShower(mcdjet, tracks, particles, maxDeltaR, false); + int jetflavourDistQuark = -1; + int jetflavourDistHadron = -1; + for (auto const& mcpjet : mcdjet.template matchedJetGeo_as()) { + jetflavourDistQuark = jettaggingutilities::getJetFlavor(mcpjet, particlesPerColl); + jetflavourDistHadron = jettaggingutilities::getJetFlavorHadron(mcpjet, particlesPerColl); + } + if (jetflavourDistQuark < 0 || jetflavourDistHadron < 0) + return; + registry.fill(HIST("h2_flavour_dist_quark_flavour_dist_hadron"), jetflavourDistQuark, jetflavourDistHadron, eventWeight); + registry.fill(HIST("h2_flavour_const_quark_flavour_const_hadron"), jetflavourConstQuark, jetflavourConstHadron, eventWeight); + registry.fill(HIST("h2_flavour_const_hadron_flavour_dist_hadron"), jetflavourConstHadron, jetflavourDistHadron, eventWeight); + registry.fill(HIST("h2_flavour_const_quark_flavour_dist_quark"), jetflavourConstQuark, jetflavourDistQuark, eventWeight); + } + + template + void fillValidationFlavourDefMCP(T const& mcpjet, U const& particles, V const& particlesPerColl, float eventWeight = 1.0) + { + float pTHat = 10. / (std::pow(eventWeight, 1.0 / pTHatExponent)); + if (mcpjet.pt() > pTHatMaxMCP * pTHat) { + return; + } + int jetflavourConstQuark = jettaggingutilities::mcpJetFromHFShower(mcpjet, particles, maxDeltaR, true); + int jetflavourConstHadron = jettaggingutilities::mcpJetFromHFShower(mcpjet, particles, maxDeltaR, false); + int jetflavourDistQuark = jettaggingutilities::getJetFlavor(mcpjet, particlesPerColl); + int jetflavourDistHadron = jettaggingutilities::getJetFlavorHadron(mcpjet, particlesPerColl); + registry.fill(HIST("h2_part_flavour_dist_quark_part_flavour_dist_hadron"), jetflavourDistQuark, jetflavourDistHadron, eventWeight); + registry.fill(HIST("h2_part_flavour_const_quark_part_flavour_const_hadron"), jetflavourConstQuark, jetflavourConstHadron, eventWeight); + registry.fill(HIST("h2_part_flavour_const_hadron_part_flavour_dist_hadron"), jetflavourConstHadron, jetflavourDistHadron, eventWeight); + registry.fill(HIST("h2_part_flavour_const_quark_part_flavour_dist_quark"), jetflavourConstQuark, jetflavourDistQuark, eventWeight); + } + template void fillHistogramIPsData(T const& jet, U const& /*tracks*/) { @@ -1013,6 +1067,44 @@ struct JetTaggerHFQA { } PROCESS_SWITCH(JetTaggerHFQA, processTracksDca, "Fill inclusive tracks' imformation for data", false); + void processValFlavourDefMCD(soa::Filtered>::iterator const& collision, soa::Join const& mcdjets, soa::Join const& /*mcpjets*/, JetTagTracksMCD const& tracks, aod::JetParticles const& particles) + { + if (collision.trackOccupancyInTimeRange() < trackOccupancyInTimeRangeMin || trackOccupancyInTimeRangeMax < collision.trackOccupancyInTimeRange()) { + return; + } + for (auto const& mcdjet : mcdjets) { + auto const particlesPerColl = particles.sliceBy(particlesPerCollision, collision.mcCollisionId()); + if (!jetfindingutilities::isInEtaAcceptance(mcdjet, jetEtaMin, jetEtaMax, trackEtaMin, trackEtaMax)) { + continue; + } + if (!isAcceptedJet(mcdjet)) { + continue; + } + fillValidationFlavourDefMCD>(mcdjet, tracks, particles, particlesPerColl, mcdjet.eventWeight()); + } + } + PROCESS_SWITCH(JetTaggerHFQA, processValFlavourDefMCD, "to check the validation of jet-flavour definition when compared to distance for mcd jets", false); + + void processValFlavourDefMCP(soa::Join const& mcpjets, aod::JetParticles const& particles, aod::JetMcCollisions const&) + { + for (auto const& mcpjet : mcpjets) { + auto const particlesPerColl = particles.sliceBy(particlesPerCollision, mcpjet.globalIndex()); + if (!jetfindingutilities::isInEtaAcceptance(mcpjet, jetEtaMin, jetEtaMax, trackEtaMin, trackEtaMax)) { + continue; + } + if (!isAcceptedJet(mcpjet)) { + continue; + } + int eventWeight = mcpjet.eventWeight(); + float pTHat = 10. / (std::pow(eventWeight, 1.0 / pTHatExponent)); + if (mcpjet.pt() > pTHatMaxMCD * pTHat) { + return; + } + fillValidationFlavourDefMCP(mcpjet, particles, particlesPerColl); + } + } + PROCESS_SWITCH(JetTaggerHFQA, processValFlavourDefMCP, "to check the validation of jet-flavour definition when compared to distance for mcp jets", false); + void processIPsData(soa::Filtered::iterator const& collision, soa::Join const& jets, JetTagTracksData const& tracks) { if (collision.trackOccupancyInTimeRange() < trackOccupancyInTimeRangeMin || trackOccupancyInTimeRangeMax < collision.trackOccupancyInTimeRange()) { @@ -1119,7 +1211,7 @@ struct JetTaggerHFQA { if (!mcdjet.has_matchedJetGeo()) continue; for (auto const& mcpjet : mcdjet.template matchedJetGeo_as>()) { - registry.fill(HIST("h3_jet_pt_jet_pt_part_matchedgeo_flavour"), mcpjet.pt(), mcdjet.pt(), mcdjet.origin()); + registry.fill(HIST("h3_jet_pt_jet_pt_part_matchedgeo_flavour"), mcdjet.pt(), mcpjet.pt(), mcdjet.origin()); } if (!doprocessIPsMCD) fillHistogramIPsMCD(mcdjet, tracks); @@ -1150,7 +1242,7 @@ struct JetTaggerHFQA { if (mcpjet.pt() > pTHatMaxMCP * pTHat) { continue; } - registry.fill(HIST("h3_jet_pt_jet_pt_part_matchedgeo_flavour"), mcpjet.pt(), mcdjet.pt(), mcdjet.origin(), mcdjet.eventWeight()); + registry.fill(HIST("h3_jet_pt_jet_pt_part_matchedgeo_flavour"), mcdjet.pt(), mcpjet.pt(), mcdjet.origin(), mcdjet.eventWeight()); } if (!doprocessIPsMCDWeighted) fillHistogramIPsMCD(mcdjet, tracks, mcdjet.eventWeight());