diff --git a/PWGJE/Core/JetMatchingUtilities.h b/PWGJE/Core/JetMatchingUtilities.h index cc972fc0a9c..e4fb0af8308 100644 --- a/PWGJE/Core/JetMatchingUtilities.h +++ b/PWGJE/Core/JetMatchingUtilities.h @@ -364,7 +364,6 @@ void MatchHF(T const& jetsBasePerCollision, U const& jetsTagPerCollision, std::v if (jetBase.candidatesIds().size() == 0) { continue; } - const auto candidateBase = jetBase.template candidates_first_as(); for (const auto& jetTag : jetsTagPerCollision) { if (jetTag.candidatesIds().size() == 0) { continue; @@ -372,22 +371,30 @@ void MatchHF(T const& jetsBasePerCollision, U const& jetsTagPerCollision, std::v if (std::round(jetBase.r()) != std::round(jetTag.r())) { continue; } - if constexpr (jetsBaseIsMc || jetsTagIsMc) { - if (jetcandidateutilities::isMatchedCandidate(candidateBase)) { - const auto candidateBaseMcId = jetcandidateutilities::matchedParticleId(candidateBase, tracksBase, tracksTag); - const auto candidateTag = jetTag.template candidates_first_as(); - const auto candidateTagId = candidateTag.mcParticleId(); - if (candidateBaseMcId == candidateTagId) { - baseToTagMatchingHF[jetBase.globalIndex()].push_back(jetTag.globalIndex()); - tagToBaseMatchingHF[jetTag.globalIndex()].push_back(jetBase.globalIndex()); + auto const& candidatesBase = jetBase.template candidates_as(); + std::size_t iCandidateBaseMatched = 0; + for (auto const& candidateBase : candidatesBase) { + if constexpr (jetsBaseIsMc || jetsTagIsMc) { + if (jetcandidateutilities::isMatchedCandidate(candidateBase)) { + const auto candidateBaseMcId = jetcandidateutilities::matchedParticleId(candidateBase, tracksBase, tracksTag); + for (auto const& candidateTag : jetTag.template candidates_as()) { + const auto candidateTagId = candidateTag.mcParticleId(); + if (candidateBaseMcId == candidateTagId) { + iCandidateBaseMatched++; + } + } + } + } else { + for (auto const& candidateTag : jetTag.template candidates_as()) { + if (candidateBase.globalIndex() == candidateTag.globalIndex()) { + iCandidateBaseMatched++; + } } } - } else { - const auto candidateTag = jetTag.template candidates_first_as(); - if (candidateBase.globalIndex() == candidateTag.globalIndex()) { - baseToTagMatchingHF[jetBase.globalIndex()].push_back(jetTag.globalIndex()); - tagToBaseMatchingHF[jetTag.globalIndex()].push_back(jetBase.globalIndex()); - } + } + if (iCandidateBaseMatched == candidatesBase.size()) { + baseToTagMatchingHF[jetBase.globalIndex()].push_back(jetTag.globalIndex()); + tagToBaseMatchingHF[jetTag.globalIndex()].push_back(jetBase.globalIndex()); } } } diff --git a/PWGJE/DataModel/JetSubstructure.h b/PWGJE/DataModel/JetSubstructure.h index cb9e17d2357..983a7944a80 100644 --- a/PWGJE/DataModel/JetSubstructure.h +++ b/PWGJE/DataModel/JetSubstructure.h @@ -204,9 +204,10 @@ DECLARE_SOA_DYNAMIC_COLUMN(R, r, //! { \ DECLARE_SOA_INDEX_COLUMN_CUSTOM(_jet_type_##CO, collision, _jet_description_ "COs"); \ DECLARE_SOA_INDEX_COLUMN_FULL_CUSTOM(Candidate, candidate, int, _cand_type_, _cand_description_ "S", "_0"); \ + DECLARE_SOA_ARRAY_INDEX_COLUMN_FULL_CUSTOM(Candidates, candidates, int32_t, _cand_type_, _cand_description_ "S", "_0"); \ DECLARE_SOA_INDEX_COLUMN(_jet_recoil_type_, jet); \ } \ - DECLARE_SOA_TABLE(_jet_type_##Os, "AOD", _jet_description_ "O", o2::soa::Index<>, _name_##jetoutput::_jet_type_##COId, _name_##jetoutput::CandidateId, jetoutput::JetPt, jetoutput::JetPhi, jetoutput::JetEta, jetoutput::JetY, jetoutput::JetR, jetoutput::JetArea, jetoutput::JetRho, jetoutput::JetPerpConeRho, jetoutput::JetNConstituents); \ + DECLARE_SOA_TABLE(_jet_type_##Os, "AOD", _jet_description_ "O", o2::soa::Index<>, _name_##jetoutput::_jet_type_##COId, _name_##jetoutput::CandidatesIds, jetoutput::JetPt, jetoutput::JetPhi, jetoutput::JetEta, jetoutput::JetY, jetoutput::JetR, jetoutput::JetArea, jetoutput::JetRho, jetoutput::JetPerpConeRho, jetoutput::JetNConstituents); \ DECLARE_SOA_TABLE(_jet_type_##Rs, "AOD", _jet_description_ "R", _name_##jetoutput::_jet_recoil_type_##Id, _name_##jetoutput::CandidateId, jetoutput::JetPt, jetoutput::JetPhi, jetoutput::JetEta, jetoutput::JetR, jetoutput::JetNConstituents, jetoutput::Pt, jetoutput::Phi, jetoutput::Eta, jetoutput::R); \ DECLARE_SOA_TABLE(_jet_type_##ROs, "AOD", _jet_description_ "RO", o2::soa::Index<>, _name_##jetoutput::CandidateId, jetoutput::JetPt, jetoutput::JetPhi, jetoutput::JetEta, jetoutput::JetR, jetoutput::JetNConstituents); \ using _jet_type_##O = _jet_type_##Os::iterator; \ diff --git a/PWGJE/JetFinders/CMakeLists.txt b/PWGJE/JetFinders/CMakeLists.txt index 396fe93a838..4426a43c4a4 100644 --- a/PWGJE/JetFinders/CMakeLists.txt +++ b/PWGJE/JetFinders/CMakeLists.txt @@ -209,4 +209,34 @@ o2physics_add_dpl_workflow(jet-finder-dielectron-mcp-charged PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::PWGJECore O2Physics::AnalysisCore O2::FrameworkPhysicsSupport COMPONENT_NAME Analysis) +o2physics_add_dpl_workflow(jet-finder-d0d0bar-data-charged + SOURCES jetFinderD0D0BarDataCharged.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::PWGJECore O2Physics::AnalysisCore O2::FrameworkPhysicsSupport + COMPONENT_NAME Analysis) + +o2physics_add_dpl_workflow(jet-finder-d0d0bar-mcd-charged + SOURCES jetFinderD0D0BarMCDCharged.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::PWGJECore O2Physics::AnalysisCore O2::FrameworkPhysicsSupport + COMPONENT_NAME Analysis) + +o2physics_add_dpl_workflow(jet-finder-d0d0bar-mcp-charged + SOURCES jetFinderD0D0BarMCPCharged.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::PWGJECore O2Physics::AnalysisCore O2::FrameworkPhysicsSupport + COMPONENT_NAME Analysis) + +o2physics_add_dpl_workflow(jet-finder-dplusdminus-data-charged + SOURCES jetFinderDplusDminusDataCharged.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::PWGJECore O2Physics::AnalysisCore O2::FrameworkPhysicsSupport + COMPONENT_NAME Analysis) + +o2physics_add_dpl_workflow(jet-finder-dplusdminus-mcd-charged + SOURCES jetFinderDplusDminusMCDCharged.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::PWGJECore O2Physics::AnalysisCore O2::FrameworkPhysicsSupport + COMPONENT_NAME Analysis) + +o2physics_add_dpl_workflow(jet-finder-dplusdminus-mcp-charged + SOURCES jetFinderDplusDminusMCPCharged.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::PWGJECore O2Physics::AnalysisCore O2::FrameworkPhysicsSupport + COMPONENT_NAME Analysis) + endif() diff --git a/PWGJE/JetFinders/jetFinderD0D0BarDataCharged.cxx b/PWGJE/JetFinders/jetFinderD0D0BarDataCharged.cxx new file mode 100644 index 00000000000..91505f01068 --- /dev/null +++ b/PWGJE/JetFinders/jetFinderD0D0BarDataCharged.cxx @@ -0,0 +1,41 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +// jet finder d0-d0bar data charged task +// +/// \author Nima Zardoshti + +#include "PWGJE/DataModel/Jet.h" +#include "PWGJE/JetFinders/jetFinderHFHFBar.h" + +#include +#include +#include +#include + +#include + +using namespace o2; +using namespace o2::framework; +using namespace o2::framework::expressions; + +using JetFinderD0D0BarDataCharged = JetFinderHFHFBarTask; + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + std::vector tasks; + + tasks.emplace_back(adaptAnalysisTask(cfgc, + SetDefaultProcesses{{{"processChargedJetsData", true}}}, + TaskName{"jet-finder-d0d0bar-data-charged"})); + + return WorkflowSpec{tasks}; +} diff --git a/PWGJE/JetFinders/jetFinderD0D0BarMCDCharged.cxx b/PWGJE/JetFinders/jetFinderD0D0BarMCDCharged.cxx new file mode 100644 index 00000000000..2430528207b --- /dev/null +++ b/PWGJE/JetFinders/jetFinderD0D0BarMCDCharged.cxx @@ -0,0 +1,41 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +// jet finder d0-d0bar mcd charged task +// +/// \author Nima Zardoshti + +#include "PWGJE/DataModel/Jet.h" +#include "PWGJE/JetFinders/jetFinderHFHFBar.h" + +#include +#include +#include +#include + +#include + +using namespace o2; +using namespace o2::framework; +using namespace o2::framework::expressions; + +using JetFinderD0D0BarMCDetectorLevelCharged = JetFinderHFHFBarTask; + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + std::vector tasks; + + tasks.emplace_back(adaptAnalysisTask(cfgc, + SetDefaultProcesses{{{"processChargedJetsMCD", true}}}, + TaskName{"jet-finder-d0d0bar-mcd-charged"})); + + return WorkflowSpec{tasks}; +} diff --git a/PWGJE/JetFinders/jetFinderD0D0BarMCPCharged.cxx b/PWGJE/JetFinders/jetFinderD0D0BarMCPCharged.cxx new file mode 100644 index 00000000000..0c401741136 --- /dev/null +++ b/PWGJE/JetFinders/jetFinderD0D0BarMCPCharged.cxx @@ -0,0 +1,41 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +// jet finder d0 mcp charged task +// +/// \author Nima Zardoshti + +#include "PWGJE/DataModel/Jet.h" +#include "PWGJE/JetFinders/jetFinderHFHFBar.h" + +#include +#include +#include +#include + +#include + +using namespace o2; +using namespace o2::framework; +using namespace o2::framework::expressions; + +using JetFinderD0D0BarMCParticleLevelCharged = JetFinderHFHFBarTask; + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + std::vector tasks; + + tasks.emplace_back(adaptAnalysisTask(cfgc, + SetDefaultProcesses{{{"processChargedJetsMCP", true}}}, + TaskName{"jet-finder-d0d0bar-mcp-charged"})); + + return WorkflowSpec{tasks}; +} diff --git a/PWGJE/JetFinders/jetFinderDplusDminusDataCharged.cxx b/PWGJE/JetFinders/jetFinderDplusDminusDataCharged.cxx new file mode 100644 index 00000000000..66552743530 --- /dev/null +++ b/PWGJE/JetFinders/jetFinderDplusDminusDataCharged.cxx @@ -0,0 +1,41 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +// jet finder D+ data charged task +// +/// \author Nima Zardoshti + +#include "PWGJE/DataModel/Jet.h" +#include "PWGJE/JetFinders/jetFinderHFHFBar.h" + +#include +#include +#include +#include + +#include + +using namespace o2; +using namespace o2::framework; +using namespace o2::framework::expressions; + +using JetFinderDplusDminusDataCharged = JetFinderHFHFBarTask; + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + std::vector tasks; + + tasks.emplace_back(adaptAnalysisTask(cfgc, + SetDefaultProcesses{{{"processChargedJetsData", true}}}, + TaskName{"jet-finder-dplusdminus-data-charged"})); + + return WorkflowSpec{tasks}; +} diff --git a/PWGJE/JetFinders/jetFinderDplusDminusMCDCharged.cxx b/PWGJE/JetFinders/jetFinderDplusDminusMCDCharged.cxx new file mode 100644 index 00000000000..24d47efb75a --- /dev/null +++ b/PWGJE/JetFinders/jetFinderDplusDminusMCDCharged.cxx @@ -0,0 +1,41 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +// jet finder D+ mcd charged task +// +/// \author Nima Zardoshti + +#include "PWGJE/DataModel/Jet.h" +#include "PWGJE/JetFinders/jetFinderHFHFBar.h" + +#include +#include +#include +#include + +#include + +using namespace o2; +using namespace o2::framework; +using namespace o2::framework::expressions; + +using JetFinderDplusDminusMCDetectorLevelCharged = JetFinderHFHFBarTask; + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + std::vector tasks; + + tasks.emplace_back(adaptAnalysisTask(cfgc, + SetDefaultProcesses{{{"processChargedJetsMCD", true}}}, + TaskName{"jet-finder-dplusdminus-mcd-charged"})); + + return WorkflowSpec{tasks}; +} diff --git a/PWGJE/JetFinders/jetFinderDplusDminusMCPCharged.cxx b/PWGJE/JetFinders/jetFinderDplusDminusMCPCharged.cxx new file mode 100644 index 00000000000..07a06ae4870 --- /dev/null +++ b/PWGJE/JetFinders/jetFinderDplusDminusMCPCharged.cxx @@ -0,0 +1,41 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +// jet finder D+ mcp charged task +// +/// \author Nima Zardoshti + +#include "PWGJE/DataModel/Jet.h" +#include "PWGJE/JetFinders/jetFinderHFHFBar.h" + +#include +#include +#include +#include + +#include + +using namespace o2; +using namespace o2::framework; +using namespace o2::framework::expressions; + +using JetFinderDplusDminusMCParticleLevelCharged = JetFinderHFHFBarTask; + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + std::vector tasks; + + tasks.emplace_back(adaptAnalysisTask(cfgc, + SetDefaultProcesses{{{"processChargedJetsMCP", true}}}, + TaskName{"jet-finder-dplusdminus-mcp-charged"})); + + return WorkflowSpec{tasks}; +} diff --git a/PWGJE/JetFinders/jetFinderHFHFBar.h b/PWGJE/JetFinders/jetFinderHFHFBar.h new file mode 100644 index 00000000000..9b36af83e6c --- /dev/null +++ b/PWGJE/JetFinders/jetFinderHFHFBar.h @@ -0,0 +1,320 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +// jet finder hf-hf bar task +// +/// \author Nima Zardoshti +/// \author Jochen Klein + +#ifndef PWGJE_JETFINDERS_JETFINDERHF_H_ +#define PWGJE_JETFINDERS_JETFINDERHF_H_ + +#include "PWGJE/Core/JetDerivedDataUtilities.h" +#include "PWGJE/Core/JetFinder.h" +#include "PWGJE/Core/JetFindingUtilities.h" +#include "PWGJE/DataModel/EMCALClusterDefinition.h" +#include "PWGJE/DataModel/EMCALClusters.h" +#include "PWGJE/DataModel/Jet.h" +#include "PWGJE/DataModel/JetReducedData.h" +#include "PWGJE/DataModel/JetSubtraction.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include // IWYU pragma: export + +#include +#include + +#include +#include + +#include +#include + +template +struct JetFinderHFHFBarTask { + o2::framework::Produces jetsTable; + o2::framework::Produces constituentsTable; + o2::framework::Produces jetsEvtWiseSubTable; + o2::framework::Produces constituentsEvtWiseSubTable; + + o2::framework::HistogramRegistry registry; + + // event level configurables + o2::framework::Configurable vertexZCut{"vertexZCut", 10.0f, "Accepted z-vertex range"}; + o2::framework::Configurable centralityMin{"centralityMin", -999.0, "minimum centrality"}; + o2::framework::Configurable centralityMax{"centralityMax", 999.0, "maximum centrality"}; + o2::framework::Configurable trackOccupancyInTimeRangeMax{"trackOccupancyInTimeRangeMax", 999999, "maximum occupancy of tracks in neighbouring collisions in a given time range"}; + o2::framework::Configurable eventSelections{"eventSelections", "sel8", "choose event selection"}; + o2::framework::Configurable triggerMasks{"triggerMasks", "", "possible JE Trigger masks: fJetChLowPt,fJetChHighPt,fTrackLowPt,fTrackHighPt,fJetD0ChLowPt,fJetD0ChHighPt,fJetLcChLowPt,fJetLcChHighPt,fEMCALReadout,fJetFullHighPt,fJetFullLowPt,fJetNeutralHighPt,fJetNeutralLowPt,fGammaVeryHighPtEMCAL,fGammaVeryHighPtDCAL,fGammaHighPtEMCAL,fGammaHighPtDCAL,fGammaLowPtEMCAL,fGammaLowPtDCAL,fGammaVeryLowPtEMCAL,fGammaVeryLowPtDCAL"}; + o2::framework::Configurable skipMBGapEvents{"skipMBGapEvents", true, "decide to run over MB gap events or not"}; + o2::framework::Configurable applyRCTSelections{"applyRCTSelections", true, "decide to apply RCT selections"}; + + // track level configurables + o2::framework::Configurable trackPtMin{"trackPtMin", 0.15, "minimum track pT"}; + o2::framework::Configurable trackPtMax{"trackPtMax", 1000.0, "maximum track pT"}; + o2::framework::Configurable trackEtaMin{"trackEtaMin", -0.9, "minimum track eta"}; + o2::framework::Configurable trackEtaMax{"trackEtaMax", 0.9, "maximum track eta"}; + o2::framework::Configurable trackPhiMin{"trackPhiMin", -999, "minimum track phi"}; + o2::framework::Configurable trackPhiMax{"trackPhiMax", 999, "maximum track phi"}; + o2::framework::Configurable applyTrackingEfficiency{"applyTrackingEfficiency", {false}, "configurable to decide whether to apply artificial tracking efficiency (discarding tracks) in jet finding"}; + o2::framework::Configurable> trackingEfficiencyPtBinning{"trackingEfficiencyPtBinning", {0., 10, 999.}, "pt binning of tracking efficiency array if applyTrackingEfficiency is true"}; + o2::framework::Configurable> trackingEfficiency{"trackingEfficiency", {1.0, 1.0}, "tracking efficiency array applied to jet finding if applyTrackingEfficiency is true"}; + o2::framework::Configurable trackSelections{"trackSelections", "globalTracks", "set track selections"}; + o2::framework::Configurable particleSelections{"particleSelections", "PhysicalPrimary", "set particle selections"}; + + // cluster level configurables + o2::framework::Configurable clusterDefinitionS{"clusterDefinition", "kV3Default", "cluster definition to be selected, e.g. V3Default"}; + o2::framework::Configurable clusterEtaMin{"clusterEtaMin", -0.71, "minimum cluster eta"}; // For ECMAL: |eta| < 0.7, phi = 1.40 - 3.26 + o2::framework::Configurable clusterEtaMax{"clusterEtaMax", 0.71, "maximum cluster eta"}; // For ECMAL: |eta| < 0.7, phi = 1.40 - 3.26 + o2::framework::Configurable clusterPhiMin{"clusterPhiMin", 1.39, "minimum cluster phi"}; + o2::framework::Configurable clusterPhiMax{"clusterPhiMax", 3.27, "maximum cluster phi"}; + o2::framework::Configurable clusterEnergyMin{"clusterEnergyMin", 0.5, "minimum cluster energy in EMCAL (GeV)"}; + o2::framework::Configurable clusterTimeMin{"clusterTimeMin", -25., "minimum Cluster time (ns)"}; + o2::framework::Configurable clusterTimeMax{"clusterTimeMax", 25., "maximum Cluster time (ns)"}; + o2::framework::Configurable clusterRejectExotics{"clusterRejectExotics", true, "Reject exotic clusters"}; + + // HF candidate level configurables + o2::framework::Configurable candPtMin{"candPtMin", 0.0, "minimum candidate pT"}; + o2::framework::Configurable candPtMax{"candPtMax", 100.0, "maximum candidate pT"}; + o2::framework::Configurable candYMin{"candYMin", -0.8, "minimum candidate eta"}; + o2::framework::Configurable candYMax{"candYMax", 0.8, "maximum candidate eta"}; + // HF candidiate selection configurables + o2::framework::Configurable rejectBackgroundMCDCandidates{"rejectBackgroundMCDCandidates", false, "reject background HF candidates at MC detector level"}; + o2::framework::Configurable rejectIncorrectDecaysMCP{"rejectIncorrectDecaysMCP", true, "reject HF paticles decaying to the non-analysed decay channels at MC generator level"}; + + // jet level configurables + o2::framework::Configurable> jetRadius{"jetRadius", {0.4}, "jet resolution parameters"}; + o2::framework::Configurable jetPtMin{"jetPtMin", 0.0, "minimum jet pT"}; + o2::framework::Configurable jetPtMax{"jetPtMax", 1000.0, "maximum jet pT"}; + o2::framework::Configurable jetEWSPtMin{"jetEWSPtMin", 0.0, "minimum event-wise subtracted jet pT"}; + o2::framework::Configurable jetEWSPtMax{"jetEWSPtMax", 1000.0, "maximum event-wise subtracted jet pT"}; + o2::framework::Configurable jetEtaMin{"jetEtaMin", -99.0, "minimum jet pseudorapidity"}; + o2::framework::Configurable jetEtaMax{"jetEtaMax", 99.0, "maximum jet pseudorapidity"}; + o2::framework::Configurable jetTypeParticleLevel{"jetTypeParticleLevel", 1, "Type of stored jets. 0 = full, 1 = charged, 2 = neutral"}; + o2::framework::Configurable jetAlgorithm{"jetAlgorithm", 2, "jet clustering algorithm. 0 = kT, 1 = C/A, 2 = Anti-kT"}; + o2::framework::Configurable jetRecombScheme{"jetRecombScheme", 0, "jet recombination scheme. 0 = E-scheme, 1 = pT-scheme, 2 = pT2-scheme"}; + o2::framework::Configurable jetGhostArea{"jetGhostArea", 0.005, "jet ghost area"}; + o2::framework::Configurable ghostRepeat{"ghostRepeat", 1, "set to 0 to gain speed if you dont need area calculation"}; + o2::framework::Configurable DoTriggering{"DoTriggering", false, "used for the charged jet trigger to remove the eta constraint on the jet axis"}; + o2::framework::Configurable jetAreaFractionMin{"jetAreaFractionMin", -99.0, "used to make a cut on the jet areas"}; + o2::framework::Configurable jetPtBinWidth{"jetPtBinWidth", 5, "used to define the width of the jetPt bins for the THnSparse"}; + o2::framework::Configurable fillTHnSparse{"fillTHnSparse", false, "switch to fill the THnSparse"}; + o2::framework::Configurable jetExtraParam{"jetExtraParam", -99.0, "sets the _extra_param in fastjet"}; + + o2::framework::Service pdgDatabase; + int trackSelection = -1; + std::vector eventSelectionBits; + std::string particleSelection; + + JetFinder jetFinder; + std::vector inputParticles; + + std::vector triggerMaskBits; + + void init(o2::framework::InitContext const&) + { + trackSelection = jetderiveddatautilities::initialiseTrackSelection(static_cast(trackSelections)); + triggerMaskBits = jetderiveddatautilities::initialiseTriggerMaskBits(triggerMasks); + eventSelectionBits = jetderiveddatautilities::initialiseEventSelectionBits(static_cast(eventSelections)); + particleSelection = static_cast(particleSelections); + + jetFinder.etaMin = trackEtaMin; + jetFinder.etaMax = trackEtaMax; + jetFinder.jetPtMin = jetPtMin; + jetFinder.jetPtMax = jetPtMax; + jetFinder.jetEtaMin = jetEtaMin; + jetFinder.jetEtaMax = jetEtaMax; + if (jetEtaMin < -98.0) { + jetFinder.jetEtaDefault = true; + } + jetFinder.algorithm = static_cast(static_cast(jetAlgorithm)); + jetFinder.recombScheme = static_cast(static_cast(jetRecombScheme)); + jetFinder.ghostArea = jetGhostArea; + jetFinder.ghostRepeatN = ghostRepeat; + if (DoTriggering) { + jetFinder.isTriggering = true; + } + jetFinder.fastjetExtraParam = jetExtraParam; + + auto jetRadiiBins = (std::vector)jetRadius; + if (jetRadiiBins.size() > 1) { + jetRadiiBins.push_back(jetRadiiBins[jetRadiiBins.size() - 1] + (TMath::Abs(jetRadiiBins[jetRadiiBins.size() - 1] - jetRadiiBins[jetRadiiBins.size() - 2]))); + } else { + jetRadiiBins.push_back(jetRadiiBins[jetRadiiBins.size() - 1] + 0.1); + } + int jetPtMaxInt = static_cast(jetPtMax); + int jetPtMinInt = static_cast(jetPtMin); + jetPtMinInt = (jetPtMinInt / jetPtBinWidth) * jetPtBinWidth; + jetPtMaxInt = ((jetPtMaxInt + jetPtBinWidth - 1) / jetPtBinWidth) * jetPtBinWidth; + int jetPtBinNumber = (jetPtMaxInt - jetPtMinInt) / jetPtBinWidth; + double jetPtMinDouble = static_cast(jetPtMinInt); + double jetPtMaxDouble = static_cast(jetPtMaxInt); + + registry.add("hJet", "sparse for data or mcd jets", {o2::framework::HistType::kTHnD, {{jetRadiiBins, ""}, {jetPtBinNumber, jetPtMinDouble, jetPtMaxDouble}, {40, -1.0, 1.0}, {18, 0.0, 7.0}}}); + registry.add("hJetMCP", "sparse for mcp jets", {o2::framework::HistType::kTHnD, {{jetRadiiBins, ""}, {jetPtBinNumber, jetPtMinDouble, jetPtMaxDouble}, {40, -1.0, 1.0}, {18, 0.0, 7.0}}}); + + if (applyTrackingEfficiency) { + if (trackingEfficiencyPtBinning->size() < 2) { + LOGP(fatal, "jetFinderHF workflow: trackingEfficiencyPtBinning configurable should have at least two bin edges"); + } + if (trackingEfficiency->size() + 1 != trackingEfficiencyPtBinning->size()) { + LOGP(fatal, "jetFinderHF workflow: trackingEfficiency configurable should have exactly one less entry than the number of bin edges set in trackingEfficiencyPtBinning configurable"); + } + } + } + + o2::aod::EMCALClusterDefinition clusterDefinition = o2::aod::emcalcluster::getClusterDefinitionFromString(clusterDefinitionS.value); + o2::framework::expressions::Filter collisionFilter = (nabs(o2::aod::jcollision::posZ) < vertexZCut && o2::aod::jcollision::centFT0M >= centralityMin && o2::aod::jcollision::centFT0M < centralityMax && o2::aod::jcollision::trackOccupancyInTimeRange <= trackOccupancyInTimeRangeMax); + o2::framework::expressions::Filter mcCollisionFilter = (nabs(o2::aod::jmccollision::posZ) < vertexZCut); + o2::framework::expressions::Filter trackCuts = (o2::aod::jtrack::pt >= trackPtMin && o2::aod::jtrack::pt < trackPtMax && o2::aod::jtrack::eta >= trackEtaMin && o2::aod::jtrack::eta <= trackEtaMax && o2::aod::jtrack::phi >= trackPhiMin && o2::aod::jtrack::phi <= trackPhiMax); + o2::framework::expressions::Filter partCuts = (o2::aod::jmcparticle::pt >= trackPtMin && o2::aod::jmcparticle::pt < trackPtMax && o2::aod::jmcparticle::eta >= trackEtaMin && o2::aod::jmcparticle::eta <= trackEtaMax && o2::aod::jmcparticle::phi >= trackPhiMin && o2::aod::jmcparticle::phi <= trackPhiMax); + o2::framework::expressions::Filter clusterFilter = (o2::aod::jcluster::definition == static_cast(clusterDefinition) && o2::aod::jcluster::eta >= clusterEtaMin && o2::aod::jcluster::eta <= clusterEtaMax && o2::aod::jcluster::phi >= clusterPhiMin && o2::aod::jcluster::phi <= clusterPhiMax && o2::aod::jcluster::energy >= clusterEnergyMin && o2::aod::jcluster::time > clusterTimeMin && o2::aod::jcluster::time < clusterTimeMax && (clusterRejectExotics && o2::aod::jcluster::isExotic != true)); + // o2::framework::expressions::Filter candidateCuts = (o2::aod::hfcand::pt >= candPtMin && o2::aod::hfcand::pt < candPtMax && o2::aod::hfcand::y >= candYMin && o2::aod::hfcand::y < candYMax); + + o2::framework::PresliceOptional> perD0Candidate = o2::aod::bkgd0::candidateId; + o2::framework::PresliceOptional> perD0McCandidate = o2::aod::bkgd0mc::candidateId; + o2::framework::PresliceOptional> perDplusCandidate = o2::aod::bkgdplus::candidateId; + o2::framework::PresliceOptional> perDplusMcCandidate = o2::aod::bkgdplusmc::candidateId; + o2::framework::PresliceOptional> perDsCandidate = o2::aod::bkgds::candidateId; + o2::framework::PresliceOptional> perDsMcCandidate = o2::aod::bkgdsmc::candidateId; + o2::framework::PresliceOptional> perDstarCandidate = o2::aod::bkgdstar::candidateId; + o2::framework::PresliceOptional> perDstarMcCandidate = o2::aod::bkgdstarmc::candidateId; + o2::framework::PresliceOptional> perLcCandidate = o2::aod::bkglc::candidateId; + o2::framework::PresliceOptional> perLcMcCandidate = o2::aod::bkglcmc::candidateId; + o2::framework::PresliceOptional> perB0Candidate = o2::aod::bkgb0::candidateId; + o2::framework::PresliceOptional> perB0McCandidate = o2::aod::bkgb0mc::candidateId; + o2::framework::PresliceOptional> perBplusCandidate = o2::aod::bkgbplus::candidateId; + o2::framework::PresliceOptional> perBplusMcCandidate = o2::aod::bkgbplusmc::candidateId; + o2::framework::PresliceOptional> perXicToXiPiPiCandidate = o2::aod::bkgxictoxipipi::candidateId; + o2::framework::PresliceOptional> perXicToXiPiPiMcCandidate = o2::aod::bkgxictoxipipimc::candidateId; + o2::framework::PresliceOptional> perDielectronCandidate = o2::aod::bkgdielectron::candidateId; + o2::framework::PresliceOptional> perDielectronMcCandidate = o2::aod::bkgdielectronmc::candidateId; + + // function that generalically processes Data and reco level events + template + void analyseCharged(T const& collision, U const& tracks, V const& candidate, V const& candidateBar, M& jetsTableInput, N& constituentsTableInput, O& /*originalTracks*/, float minJetPt, float maxJetPt) + { + if (candidate.globalIndex() == candidateBar.globalIndex() || candidate.candidateSelFlag() == candidateBar.candidateSelFlag()) { + return; + } + for (auto const& track : tracks) { + if (jetcandidateutilities::isDaughterTrack(track, candidate) && jetcandidateutilities::isDaughterTrack(track, candidateBar)) { + return; + } + } + if (!jetderiveddatautilities::selectCollision(collision, eventSelectionBits, skipMBGapEvents, applyRCTSelections) || !jetderiveddatautilities::selectTrigger(collision, triggerMaskBits)) { + return; + } + inputParticles.clear(); + + if constexpr (jetcandidateutilities::isCandidate()) { + if (!jetfindingutilities::analyseCandidate(inputParticles, candidate, candPtMin, candPtMax, candYMin, candYMax) || !jetfindingutilities::analyseCandidate(inputParticles, candidateBar, candPtMin, candPtMax, candYMin, candYMax)) { + return; + } + } + + if constexpr (jetcandidateutilities::isMcCandidate()) { + if (!jetfindingutilities::analyseCandidateMC(inputParticles, candidate, candPtMin, candPtMax, candYMin, candYMax, rejectBackgroundMCDCandidates) || !jetfindingutilities::analyseCandidateMC(inputParticles, candidateBar, candPtMin, candPtMax, candYMin, candYMax, rejectBackgroundMCDCandidates)) { + return; + } + } + jetfindingutilities::analyseTracks(inputParticles, tracks, trackSelection, applyTrackingEfficiency, trackingEfficiency, trackingEfficiencyPtBinning, &candidate); + + jetfindingutilities::findJets(jetFinder, inputParticles, minJetPt, maxJetPt, jetRadius, jetAreaFractionMin, collision, jetsTableInput, constituentsTableInput, registry.get(HIST("hJet")), fillTHnSparse, true); + } + + // function that generalically processes gen level events + template + void analyseMCP(T const& mcCollision, U const& particles, V const& candidate, V const& candidateBar, M& jetsTableInput, N& constituentsTableInput, int jetTypeParticleLevel, float minJetPt, float maxJetPt) + { + if (candidate.globalIndex() == candidateBar.globalIndex() || candidate.flagMcMatchGen() == candidateBar.flagMcMatchGen()) { + return; + } + for (auto const& particle : particles) { + if (jetcandidateutilities::isDaughterParticle(candidate.template mcParticle_as(), particle.globalIndex()) || jetcandidateutilities::isDaughterParticle(candidateBar.template mcParticle_as(), particle.globalIndex())) { + return; + } + } + if (!jetderiveddatautilities::selectMcCollision(mcCollision, skipMBGapEvents, applyRCTSelections)) { + return; + } + + if (rejectIncorrectDecaysMCP && (!jetcandidateutilities::isMatchedCandidate(candidate) || !jetcandidateutilities::isMatchedCandidate(candidateBar))) { // is this even needed in the new derived format? it means any simulations run have to force the decay channel + return; + } + + inputParticles.clear(); + if (!jetfindingutilities::analyseCandidate(inputParticles, candidate, candPtMin, candPtMax, candYMin, candYMax) || !jetfindingutilities::analyseCandidate(inputParticles, candidateBar, candPtMin, candPtMax, candYMin, candYMax)) { + return; + } + jetfindingutilities::analyseParticles(inputParticles, particleSelection, jetTypeParticleLevel, particles, pdgDatabase, &candidate); + + jetfindingutilities::findJets(jetFinder, inputParticles, minJetPt, maxJetPt, jetRadius, jetAreaFractionMin, mcCollision, jetsTableInput, constituentsTableInput, registry.get(HIST("hJetMCP")), fillTHnSparse, true); + } + + void processDummy(o2::aod::JetCollisions const&) + { + } + PROCESS_SWITCH(JetFinderHFHFBarTask, processDummy, "Dummy process function turned on by default", true); + + void processChargedJetsData(o2::soa::Filtered::iterator const& collision, o2::soa::Filtered const& tracks, CandidateTableData const& candidates) + { + for (auto candidateIterator = candidates.begin(); candidateIterator != candidates.end(); ++candidateIterator) { + auto candidateBarIterator = candidateIterator; + ++candidateBarIterator; + for (; candidateBarIterator != candidates.end(); ++candidateBarIterator) { + typename CandidateTableData::iterator const& candidate = *candidateIterator; + typename CandidateTableData::iterator const& candidateBar = *candidateBarIterator; + analyseCharged(collision, tracks, candidate, candidateBar, jetsTable, constituentsTable, tracks, jetPtMin, jetPtMax); + } + } + } + PROCESS_SWITCH(JetFinderHFHFBarTask, processChargedJetsData, "charged hf jet finding on data", false); + + void processChargedJetsMCD(o2::soa::Filtered::iterator const& collision, o2::soa::Filtered const& tracks, CandidateTableMCD const& candidates) + { + for (auto candidateIterator = candidates.begin(); candidateIterator != candidates.end(); ++candidateIterator) { + auto candidateBarIterator = candidateIterator; + ++candidateBarIterator; + for (; candidateBarIterator != candidates.end(); ++candidateBarIterator) { + typename CandidateTableMCD::iterator const& candidate = *candidateIterator; + typename CandidateTableMCD::iterator const& candidateBar = *candidateBarIterator; + analyseCharged(collision, tracks, candidate, candidateBar, jetsTable, constituentsTable, tracks, jetPtMin, jetPtMax); + } + } + } + PROCESS_SWITCH(JetFinderHFHFBarTask, processChargedJetsMCD, "charged hf jet finding on MC detector level", false); + + void processChargedJetsMCP(o2::soa::Filtered::iterator const& mcCollision, + o2::soa::Filtered const& particles, + CandidateTableMCP const& candidates) + { + for (auto candidateIterator = candidates.begin(); candidateIterator != candidates.end(); ++candidateIterator) { + auto candidateBarIterator = candidateIterator; + ++candidateBarIterator; + for (; candidateBarIterator != candidates.end(); ++candidateBarIterator) { + typename CandidateTableMCP::iterator const& candidate = *candidateIterator; + typename CandidateTableMCP::iterator const& candidateBar = *candidateBarIterator; + analyseMCP(mcCollision, particles, candidate, candidateBar, jetsTable, constituentsTable, 1, jetPtMin, jetPtMax); + } + } + } + PROCESS_SWITCH(JetFinderHFHFBarTask, processChargedJetsMCP, "hf jet finding on MC particle level", false); +}; + +#endif // PWGJE_JETFINDERS_JETFINDERHF_H_ diff --git a/PWGJE/Tasks/jetSubstructureHF.h b/PWGJE/Tasks/jetSubstructureHF.h index 5bb88758ea3..7aa47a466df 100644 --- a/PWGJE/Tasks/jetSubstructureHF.h +++ b/PWGJE/Tasks/jetSubstructureHF.h @@ -191,7 +191,7 @@ struct JetSubstructureHFTask { } template - void jetReclustering(T const& jet, U& splittingTable) + void jetReclustering(T const& jet, U& splittingTable, int nHFCandidates) { energyMotherVec.clear(); ptLeadingVec.clear(); @@ -207,16 +207,26 @@ struct JetSubstructureHFTask { auto nsd = 0.0; while (daughterSubJet.has_parents(parentSubJet1, parentSubJet2)) { - bool isHFInSubjet1 = false; + int nHFInSubjet1 = 0; for (auto& subjet1Constituent : parentSubJet1.constituents()) { if (subjet1Constituent.template user_info().getStatus() == static_cast(JetConstituentStatus::candidate)) { - isHFInSubjet1 = true; - break; + nHFInSubjet1++; } } - if (!isHFInSubjet1) { + if (nHFCandidates == 1 && nHFInSubjet1 == 0) { std::swap(parentSubJet1, parentSubJet2); } + if (nHFCandidates == 2) { + if (nHFInSubjet1 == 2) { + daughterSubJet = parentSubJet1; + continue; + } + if (nHFInSubjet1 == 0) { + daughterSubJet = parentSubJet2; + continue; + } + } + std::vector tracks; std::vector candidates; std::vector clusters; @@ -256,6 +266,9 @@ struct JetSubstructureHFTask { nsd++; } daughterSubJet = parentSubJet1; + if (nHFCandidates == 2 && nHFInSubjet1 == 1) { + break; + } } if constexpr (!isSubtracted && !isMCP) { registry.fill(HIST("h2_jet_pt_jet_nsd"), jet.pt(), nsd); @@ -467,11 +480,13 @@ struct JetSubstructureHFTask { for (auto& jetConstituent : jet.template tracks_as()) { fastjetutilities::fillTracks(jetConstituent, jetConstituents, jetConstituent.globalIndex()); } - for (auto& jetHFCandidate : jet.template candidates_as()) { // should only be one at the moment + int nHFCandidates = 0; + for (auto& jetHFCandidate : jet.template candidates_as()) { fastjetutilities::fillTracks(jetHFCandidate, jetConstituents, jetHFCandidate.globalIndex(), static_cast(JetConstituentStatus::candidate), candMass); + nHFCandidates++; } nSub = jetsubstructureutilities::getNSubjettiness(jet, tracks, tracks, candidates, 2, fastjet::contrib::CA_Axes(), true, zCut, beta); - jetReclustering(jet, splittingTable); + jetReclustering(jet, splittingTable, nHFCandidates); jetPairing(jet, tracks, candidates, trackSlicer, pairTable); jetSubstructureSimple(jet, tracks, candidates); outputTable(energyMotherVec, ptLeadingVec, ptSubLeadingVec, thetaVec, nSub[0], nSub[1], nSub[2], pairJetPtVec, pairJetEnergyVec, pairJetThetaVec, pairJetPerpCone1PtVec, pairJetPerpCone1EnergyVec, pairJetPerpCone1ThetaVec, pairPerpCone1PerpCone1PtVec, pairPerpCone1PerpCone1EnergyVec, pairPerpCone1PerpCone1ThetaVec, pairPerpCone1PerpCone2PtVec, pairPerpCone1PerpCone2EnergyVec, pairPerpCone1PerpCone2ThetaVec, angularity, leadingConstituentPt, perpConeRho); @@ -522,11 +537,13 @@ struct JetSubstructureHFTask { for (auto& jetConstituent : jet.template tracks_as()) { fastjetutilities::fillTracks(jetConstituent, jetConstituents, jetConstituent.globalIndex(), static_cast(JetConstituentStatus::track), pdg->Mass(jetConstituent.pdgCode())); } + int nHFCandidates = 0; for (auto& jetHFCandidate : jet.template candidates_as()) { fastjetutilities::fillTracks(jetHFCandidate, jetConstituents, jetHFCandidate.globalIndex(), static_cast(JetConstituentStatus::candidate), candMass); + nHFCandidates++; } nSub = jetsubstructureutilities::getNSubjettiness(jet, particles, particles, candidates, 2, fastjet::contrib::CA_Axes(), true, zCut, beta); - jetReclustering(jet, jetSplittingsMCPTable); + jetReclustering(jet, jetSplittingsMCPTable, nHFCandidates); jetPairing(jet, particles, candidates, ParticlesPerMcCollision, jetPairsMCPTable); jetSubstructureSimple(jet, particles, candidates); jetSubstructureMCPTable(energyMotherVec, ptLeadingVec, ptSubLeadingVec, thetaVec, nSub[0], nSub[1], nSub[2], pairJetPtVec, pairJetEnergyVec, pairJetThetaVec, pairJetPerpCone1PtVec, pairJetPerpCone1EnergyVec, pairJetPerpCone1ThetaVec, pairPerpCone1PerpCone1PtVec, pairPerpCone1PerpCone1EnergyVec, pairPerpCone1PerpCone1ThetaVec, pairPerpCone1PerpCone2PtVec, pairPerpCone1PerpCone2EnergyVec, pairPerpCone1PerpCone2ThetaVec, angularity, leadingConstituentPt, perpConeRho); diff --git a/PWGJE/Tasks/jetSubstructureHFOutput.h b/PWGJE/Tasks/jetSubstructureHFOutput.h index 377251928ad..ac35f02ab9d 100644 --- a/PWGJE/Tasks/jetSubstructureHFOutput.h +++ b/PWGJE/Tasks/jetSubstructureHFOutput.h @@ -111,6 +111,9 @@ struct JetSubstructureHFOutputTask { std::vector candidateSelectionFlagsData; std::vector candidateSelectionFlagsMCD; std::vector candidateSelectionFlagsMCP; + std::vector candidateRecoilSelectionFlagsData; + std::vector candidateRecoilSelectionFlagsMCD; + std::vector candidateRecoilSelectionFlagsMCP; std::vector> splittingMatchesGeoVecVecData; std::vector> splittingMatchesPtVecVecData; @@ -272,8 +275,8 @@ struct JetSubstructureHFOutputTask { } } - template - void fillJetTables(T const& jet, U const& /*cand*/, int32_t collisionIndex, int32_t candidateIndex, V& jetOutputTable, M& jetSubstructureOutputTable, std::vector>& splittingMatchesGeoVecVec, std::vector>& splittingMatchesPtVecVec, std::vector>& splittingMatchesHFVecVec, std::vector>& pairMatchesVecVec, float rho, std::map& jetMap) + template + void fillJetTables(T const& jet, int32_t collisionIndex, std::vector const& candidatesIndices, U& jetOutputTable, V& jetSubstructureOutputTable, std::vector>& splittingMatchesGeoVecVec, std::vector>& splittingMatchesPtVecVec, std::vector>& splittingMatchesHFVecVec, std::vector>& pairMatchesVecVec, float rho, std::map& jetMap) { std::vector energyMotherVec; std::vector ptLeadingVec; @@ -337,7 +340,7 @@ struct JetSubstructureHFOutputTask { pairMatchesVec = pairMatchesVecVec[jet.globalIndex()]; } - jetOutputTable(collisionIndex, candidateIndex, jet.pt(), jet.phi(), jet.eta(), jet.y(), jet.r(), jet.area(), rho, jet.perpConeRho(), jet.tracksIds().size() + jet.candidatesIds().size()); // here we take the decision to keep the collision index consistent with the JE framework in case it is later needed to join to other tables. The candidate Index however can be linked to the HF tables + jetOutputTable(collisionIndex, candidatesIndices, jet.pt(), jet.phi(), jet.eta(), jet.y(), jet.r(), jet.area(), rho, jet.perpConeRho(), jet.tracksIds().size() + jet.candidatesIds().size()); // here we take the decision to keep the collision index consistent with the JE framework in case it is later needed to join to other tables. The candidate Index however can be linked to the HF tables jetSubstructureOutputTable(jetOutputTable.lastIndex(), energyMotherVec, ptLeadingVec, ptSubLeadingVec, thetaVec, jet.nSub2DR(), jet.nSub1(), jet.nSub2(), pairJetPtVec, pairJetEnergyVec, pairJetThetaVec, pairJetPerpCone1PtVec, pairJetPerpCone1EnergyVec, pairJetPerpCone1ThetaVec, pairPerpCone1PerpCone1PtVec, pairPerpCone1PerpCone1EnergyVec, pairPerpCone1PerpCone1ThetaVec, pairPerpCone1PerpCone2PtVec, pairPerpCone1PerpCone2EnergyVec, pairPerpCone1PerpCone2ThetaVec, jet.angularity(), jet.ptLeadingConstituent(), splittingMatchesGeoVec, splittingMatchesPtVec, splittingMatchesHFVec, pairMatchesVec); jetMap.insert(std::make_pair(jet.globalIndex(), jetOutputTable.lastIndex())); } @@ -348,6 +351,7 @@ struct JetSubstructureHFOutputTask { int nJetInCollision = 0; int32_t collisionIndex = -1; + float rho = 0.0; for (const auto& jet : jets) { if (jet.pt() < jetPtMin) { continue; @@ -357,11 +361,13 @@ struct JetSubstructureHFOutputTask { } for (const auto& jetRadiiValue : jetRadiiValues) { if (jet.r() == round(jetRadiiValue * 100.0f)) { - auto candidate = jet.template candidates_first_as(); - int32_t candidateIndex = -1; - auto candidateTableIndex = candidateMap.find(candidate.globalIndex()); - if (candidateTableIndex != candidateMap.end()) { - candidateIndex = candidateTableIndex->second; + std::vector candidatesIndices; + for (auto const& candidate : jet.template candidates_as()) { + auto candidateTableIndex = candidateMap.find(candidate.globalIndex()); + if (candidateTableIndex != candidateMap.end()) { + candidatesIndices.push_back(candidateTableIndex->second); + } + rho = candidate.rho(); } if (nJetInCollision == 0) { float centrality = -1.0; @@ -374,7 +380,7 @@ struct JetSubstructureHFOutputTask { collisionIndex = collisionOutputTable.lastIndex(); } nJetInCollision++; - fillJetTables(jet, candidate, collisionIndex, candidateIndex, jetOutputTable, jetSubstructureOutputTable, splittingMatchesGeoVecVec, splittingMatchesPtVecVec, splittingMatchesHFVecVec, pairMatchesVecVec, candidate.rho(), jetMap); + fillJetTables(jet, collisionIndex, candidatesIndices, jetOutputTable, jetSubstructureOutputTable, splittingMatchesGeoVecVec, splittingMatchesPtVecVec, splittingMatchesHFVecVec, pairMatchesVecVec, rho, jetMap); } } } @@ -406,14 +412,17 @@ struct JetSubstructureHFOutputTask { void selectCandidates(T const& jets, U const& /*candidates*/, float jetPtMin, std::vector& candidateSelectionFlags) { for (const auto& jet : jets) { - auto candidate = jet.template candidates_first_as(); - auto candidateId = candidate.globalIndex(); + auto const& candidates = jet.template candidates_as(); if (jet.pt() < jetPtMin) { - candidateSelectionFlags[candidateId] = false; + for (const auto& candidate : candidates) { + candidateSelectionFlags[candidate.globalIndex()] = false; + } continue; } if (!jetfindingutilities::isInEtaAcceptance(jet, configs.jetEtaMin, configs.jetEtaMax, configs.trackEtaMin, configs.trackEtaMax)) { - candidateSelectionFlags[candidateId] = false; + for (const auto& candidate : candidates) { + candidateSelectionFlags[candidate.globalIndex()] = false; + } continue; } bool radiusSelected = false; @@ -423,7 +432,9 @@ struct JetSubstructureHFOutputTask { } } if (!radiusSelected) { - candidateSelectionFlags[candidateId] = false; + for (const auto& candidate : candidates) { + candidateSelectionFlags[candidate.globalIndex()] = false; + } continue; } } @@ -456,9 +467,9 @@ struct JetSubstructureHFOutputTask { } template - void analyseCandidate(T const& candidate, std::map& candidateMap, std::vector& candidateSelectionFlags) + void analyseCandidate(T const& candidate, std::map& candidateMap, std::vector const& candidateSelectionFlags, std::vector const& candidateRecoilSelectionFlags) { - if (!candidateSelectionFlags[candidate.globalIndex()]) { + if (!candidateSelectionFlags[candidate.globalIndex()] && !candidateRecoilSelectionFlags[candidate.globalIndex()]) { return; } int32_t candidateCollisionIndex = -1; @@ -603,8 +614,8 @@ struct JetSubstructureHFOutputTask { } } - template - void analyseHFCollisions(T const& collisions, U const& mcCollisions, V const& hfCollisions, M const& hfMcCollisions, N const& jets, O const& jetsMCP, P const& candidates, S const& candidatesMCP, float jetPtMin, float jetPtMinMCP = 0.0) + template + void analyseHFCollisions(T const& collisions, U const& mcCollisions, V const& hfCollisions, M const& hfMcCollisions, N const& candidates, O const& candidatesMCP, std::vector const& candidateSelectionFlags, std::vector const& candidateRecoilSelectionFlags, std::vector const& candidateMCPSelectionFlags, std::vector const& candidateMCPRecoilSelectionFlags) { collisionFlag.clear(); collisionFlag.resize(collisions.size()); @@ -615,43 +626,27 @@ struct JetSubstructureHFOutputTask { std::fill(mcCollisionFlag.begin(), mcCollisionFlag.end(), false); if constexpr (!isMCPOnly) { - for (const auto& jet : jets) { - if (jet.pt() < jetPtMin) { - continue; - } - if (!jetfindingutilities::isInEtaAcceptance(jet, configs.jetEtaMin, configs.jetEtaMax, configs.trackEtaMin, configs.trackEtaMax)) { - continue; - } - for (const auto& jetRadiiValue : jetRadiiValues) { - if (jet.r() == round(jetRadiiValue * 100.0f)) { - collisionFlag[jet.collisionId()] = true; - if constexpr (isMC) { - auto mcCollisionId = jet.template collision_as().mcCollisionId(); - if (mcCollisionId >= 0) { - mcCollisionFlag[mcCollisionId] = true; - } + for (auto const& candidate : candidates) { + + if (candidateSelectionFlags[candidate.globalIndex()] || candidateRecoilSelectionFlags[candidate.globalIndex()]) { + collisionFlag[candidate.collisionId()] = true; + if constexpr (isMC) { + auto mcCollisionId = candidate.template collision_as().mcCollisionId(); + if (mcCollisionId >= 0) { + mcCollisionFlag[mcCollisionId] = true; } } } } } if constexpr (isMC) { - for (const auto& jetMCP : jetsMCP) { - if (jetMCP.pt() < jetPtMinMCP) { - continue; - } - if (!jetfindingutilities::isInEtaAcceptance(jetMCP, configs.jetEtaMin, configs.jetEtaMax, configs.trackEtaMin, configs.trackEtaMax)) { - continue; - } - for (const auto& jetRadiiValue : jetRadiiValues) { - if (jetMCP.r() == round(jetRadiiValue * 100.0f)) { - - mcCollisionFlag[jetMCP.mcCollisionId()] = true; - if constexpr (!isMCPOnly) { - const auto collisionsPerMcCollision = collisions.sliceBy(preslices.CollisionsPerMcCollision, jetMCP.mcCollisionId()); - for (auto collision : collisionsPerMcCollision) { - collisionFlag[collision.globalIndex()] = true; - } + for (auto const& candidateMCP : candidatesMCP) { + if (candidateMCPSelectionFlags[candidateMCP.globalIndex()] || candidateMCPRecoilSelectionFlags[candidateMCP.globalIndex()]) { + mcCollisionFlag[candidateMCP.mcCollisionId()] = true; + if constexpr (!isMCPOnly) { + const auto collisionsPerMcCollision = collisions.sliceBy(preslices.CollisionsPerMcCollision, candidateMCP.mcCollisionId()); + for (auto collision : collisionsPerMcCollision) { + collisionFlag[collision.globalIndex()] = true; } } } @@ -683,7 +678,7 @@ struct JetSubstructureHFOutputTask { } jetcandidateutilities::fillCandidateMcCollisionTable(hfMcCollisionPerMcCollision, candidatesMCP, products.hfMcCollisionsTable); candidateMcCollisionMapping.insert(std::make_pair(hfMcCollisionPerMcCollision.globalIndex(), products.hfMcCollisionsTable.lastIndex())); - if constexpr (!isMCPOnly && (jethfutilities::isHFTable

() || jethfutilities::isHFMcTable())) { // the matching of mcCollision to Collision is only done for HF tables + if constexpr (!isMCPOnly && (jethfutilities::isHFTable() || jethfutilities::isHFMcTable())) { // the matching of mcCollision to Collision is only done for HF tables std::vector hfCollisionIDs; for (auto const& hfCollisionPerMcCollision : hfMcCollisionPerMcCollision.template hfCollBases_as()) { // if added for others this line needs to be templated per type auto hfCollisionIndex = candidateCollisionMapping.find(hfCollisionPerMcCollision.globalIndex()); @@ -711,10 +706,14 @@ struct JetSubstructureHFOutputTask { if (doprocessOutputCandidatesData) { candidateSelectionFlagsData.clear(); candidateSelectionFlagsData.resize(candidates.size(), true); + candidateRecoilSelectionFlagsData.clear(); + candidateRecoilSelectionFlagsData.resize(candidates.size(), true); } if (doprocessOutputCandidatesMCD) { candidateSelectionFlagsMCD.clear(); candidateSelectionFlagsMCD.resize(candidates.size(), true); + candidateRecoilSelectionFlagsMCD.clear(); + candidateRecoilSelectionFlagsMCD.resize(candidates.size(), true); } } PROCESS_SWITCH(JetSubstructureHFOutputTask, processClearMaps, "process function that clears all the non-mcp maps in each dataframe", true); @@ -730,49 +729,11 @@ struct JetSubstructureHFOutputTask { } candidateSelectionFlagsMCP.clear(); candidateSelectionFlagsMCP.resize(candidates.size(), true); + candidateRecoilSelectionFlagsMCP.clear(); + candidateRecoilSelectionFlagsMCP.resize(candidates.size(), true); } PROCESS_SWITCH(JetSubstructureHFOutputTask, processClearMapsMCP, "process function that clears all the mcp maps in each dataframe", true); - void processOutputCollisionsData(o2::aod::JetCollisions const& collisions, - JetTableData const& jets, - CandidateCollisionTable const& canidateCollisions, - CandidateTable const& candidates) - { - analyseHFCollisions(collisions, collisions, canidateCollisions, canidateCollisions, jets, jets, candidates, candidates, configs.jetPtMinData); - } - PROCESS_SWITCH(JetSubstructureHFOutputTask, processOutputCollisionsData, "hf collision output data", false); - - void processOutputCollisionsDataSub(o2::aod::JetCollisions const& collisions, - JetTableDataSub const& jets, - CandidateCollisionTable const& canidateCollisions, - CandidateTable const& candidates) - { - analyseHFCollisions(collisions, collisions, canidateCollisions, canidateCollisions, jets, jets, candidates, candidates, configs.jetPtMinDataSub); - } - PROCESS_SWITCH(JetSubstructureHFOutputTask, processOutputCollisionsDataSub, "hf collision output data eventwise constituent subtracted", false); - - void processOutputCollisionsMc(o2::soa::Join const& collisions, - o2::aod::JetMcCollisions const& mcCollisions, - JetTableMCD const& jetsMCD, - JetTableMatchedMCP const& jetsMCP, - CandidateCollisionTable const& canidateCollisions, - CandidateMcCollisionTable const& canidateMcCollisions, - CandidateTableMCD const& candidatesMCD, - CandidateTableMCP const& candidatesMCP) - { - analyseHFCollisions(collisions, mcCollisions, canidateCollisions, canidateMcCollisions, jetsMCD, jetsMCP, candidatesMCD, candidatesMCP, configs.jetPtMinMCD, configs.jetPtMinMCP); - } - PROCESS_SWITCH(JetSubstructureHFOutputTask, processOutputCollisionsMc, "hf collision output MC", false); - - void processOutputCollisionsMCPOnly(o2::aod::JetMcCollisions const& mcCollisions, - JetTableMCP const& jetsMCP, - CandidateMcOnlyCollisionTable const& canidateMcCollisions, - CandidateTableMCP const& candidatesMCP) - { - analyseHFCollisions(mcCollisions, mcCollisions, canidateMcCollisions, canidateMcCollisions, jetsMCP, jetsMCP, candidatesMCP, candidatesMCP, 0.0, configs.jetPtMinMCP); - } - PROCESS_SWITCH(JetSubstructureHFOutputTask, processOutputCollisionsMCPOnly, "hf collision output MCP only", false); - void processSelectCandidatesData(JetTableData const& jets, CandidateTable const& candidates) { selectCandidates(jets, candidates, configs.jetPtMinData, candidateSelectionFlagsData); @@ -781,7 +742,7 @@ struct JetSubstructureHFOutputTask { void processSelectRecoilCandidatesData(RecoilTableData const& jets) { - selectRecoilCandidates(jets, configs.recoilJetPtMinData, candidateSelectionFlagsData); + selectRecoilCandidates(jets, configs.recoilJetPtMinData, candidateRecoilSelectionFlagsData); } PROCESS_SWITCH(JetSubstructureHFOutputTask, processSelectRecoilCandidatesData, "select HF candidates for recoil data", false); @@ -793,7 +754,7 @@ struct JetSubstructureHFOutputTask { void processSelectRecoilCandidatesMCD(RecoilTableMCD const& jets) { - selectRecoilCandidates(jets, configs.recoilJetPtMinMCD, candidateSelectionFlagsMCD); + selectRecoilCandidates(jets, configs.recoilJetPtMinMCD, candidateRecoilSelectionFlagsMCD); } PROCESS_SWITCH(JetSubstructureHFOutputTask, processSelectRecoilCandidatesMCD, "select HF candidates for recoil mcd", false); @@ -805,26 +766,61 @@ struct JetSubstructureHFOutputTask { void processSelectRecoilCandidatesMCP(RecoilTableMCP const& jets) { - selectRecoilCandidates(jets, configs.recoilJetPtMinMCP, candidateSelectionFlagsMCP); + selectRecoilCandidates(jets, configs.recoilJetPtMinMCP, candidateRecoilSelectionFlagsMCP); } PROCESS_SWITCH(JetSubstructureHFOutputTask, processSelectRecoilCandidatesMCP, "select HF candidates for recoil MCP", false); + void processOutputCollisionsData(o2::aod::JetCollisions const& collisions, + CandidateCollisionTable const& canidateCollisions, + CandidateTable const& candidates) + { + analyseHFCollisions(collisions, collisions, canidateCollisions, canidateCollisions, candidates, candidates, candidateSelectionFlagsData, candidateRecoilSelectionFlagsData, candidateSelectionFlagsData, candidateRecoilSelectionFlagsData); + } + PROCESS_SWITCH(JetSubstructureHFOutputTask, processOutputCollisionsData, "hf collision output data", false); + + void processOutputCollisionsDataSub(o2::aod::JetCollisions const& collisions, + CandidateCollisionTable const& canidateCollisions, + CandidateTable const& candidates) + { + analyseHFCollisions(collisions, collisions, canidateCollisions, canidateCollisions, candidates, candidates, candidateSelectionFlagsData, candidateRecoilSelectionFlagsData, candidateSelectionFlagsData, candidateRecoilSelectionFlagsData); + } + PROCESS_SWITCH(JetSubstructureHFOutputTask, processOutputCollisionsDataSub, "hf collision output data eventwise constituent subtracted", false); + + void processOutputCollisionsMc(o2::soa::Join const& collisions, + o2::aod::JetMcCollisions const& mcCollisions, + CandidateCollisionTable const& canidateCollisions, + CandidateMcCollisionTable const& canidateMcCollisions, + CandidateTableMCD const& candidatesMCD, + CandidateTableMCP const& candidatesMCP) + { + analyseHFCollisions(collisions, mcCollisions, canidateCollisions, canidateMcCollisions, candidatesMCD, candidatesMCP, candidateSelectionFlagsMCD, candidateRecoilSelectionFlagsMCD, candidateSelectionFlagsMCP, candidateRecoilSelectionFlagsMCP); + } + PROCESS_SWITCH(JetSubstructureHFOutputTask, processOutputCollisionsMc, "hf collision output MC", false); + + void processOutputCollisionsMCPOnly(o2::aod::JetMcCollisions const& mcCollisions, + CandidateMcOnlyCollisionTable const& canidateMcCollisions, + CandidateTableMCP const& candidatesMCP) + { + analyseHFCollisions(mcCollisions, mcCollisions, canidateMcCollisions, canidateMcCollisions, candidatesMCP, candidatesMCP, candidateSelectionFlagsMCP, candidateRecoilSelectionFlagsMCP, candidateSelectionFlagsMCP, candidateRecoilSelectionFlagsMCP); + } + PROCESS_SWITCH(JetSubstructureHFOutputTask, processOutputCollisionsMCPOnly, "hf collision output MCP only", false); + void processOutputCandidatesData(typename CandidateTable::iterator const& candidate) { - analyseCandidate(candidate, candidateMapping, candidateSelectionFlagsData); + analyseCandidate(candidate, candidateMapping, candidateSelectionFlagsData, candidateRecoilSelectionFlagsData); } PROCESS_SWITCH(JetSubstructureHFOutputTask, processOutputCandidatesData, "hf candidate output data", false); void processOutputCandidatesMCD(typename CandidateTableMCD::iterator const& candidate) { - analyseCandidate(candidate, candidateMapping, candidateSelectionFlagsMCD); + analyseCandidate(candidate, candidateMapping, candidateSelectionFlagsMCD, candidateRecoilSelectionFlagsMCD); } PROCESS_SWITCH(JetSubstructureHFOutputTask, processOutputCandidatesMCD, "hf candidate output MCD", false); void processOutputCandidatesMCP(typename CandidateTableMCP::iterator const& candidate) { - analyseCandidate(candidate, candidateMappingMCP, candidateSelectionFlagsMCP); + analyseCandidate(candidate, candidateMappingMCP, candidateSelectionFlagsMCP, candidateRecoilSelectionFlagsMCP); } PROCESS_SWITCH(JetSubstructureHFOutputTask, processOutputCandidatesMCP, "hf candidate output MCP", false); diff --git a/PWGJE/Tasks/jetSubstructureOutput.cxx b/PWGJE/Tasks/jetSubstructureOutput.cxx index 39f1e1518cd..32a53ca6516 100644 --- a/PWGJE/Tasks/jetSubstructureOutput.cxx +++ b/PWGJE/Tasks/jetSubstructureOutput.cxx @@ -191,13 +191,14 @@ struct JetSubstructureOutputTask { std::vector splittingMatchesPtVec; std::vector splittingMatchesHFVec; std::vector pairMatchesVec; + std::vector dummyVec; if (doprocessOutputSubstructureMatchingData || doprocessOutputSubstructureMatchingMC) { splittingMatchesGeoVec = splittingMatchesGeoVecVec[jet.globalIndex()]; splittingMatchesPtVec = splittingMatchesPtVecVec[jet.globalIndex()]; splittingMatchesHFVec = splittingMatchesHFVecVec[jet.globalIndex()]; pairMatchesVec = pairMatchesVecVec[jet.globalIndex()]; } - jetOutputTable(collisionIndex, collisionIndex, jet.pt(), jet.phi(), jet.eta(), jet.y(), jet.r(), jet.area(), rho, jet.perpConeRho(), jet.tracksIds().size()); // second collision index is a dummy coloumn mirroring the hf candidate + jetOutputTable(collisionIndex, dummyVec, jet.pt(), jet.phi(), jet.eta(), jet.y(), jet.r(), jet.area(), rho, jet.perpConeRho(), jet.tracksIds().size()); // second collision index is a dummy coloumn mirroring the hf candidate jetSubstructureOutputTable(jetOutputTable.lastIndex(), energyMotherVec, ptLeadingVec, ptSubLeadingVec, thetaVec, jet.nSub2DR(), jet.nSub1(), jet.nSub2(), pairJetPtVec, pairJetEnergyVec, pairJetThetaVec, pairJetPerpCone1PtVec, pairJetPerpCone1EnergyVec, pairJetPerpCone1ThetaVec, pairPerpCone1PerpCone1PtVec, pairPerpCone1PerpCone1EnergyVec, pairPerpCone1PerpCone1ThetaVec, pairPerpCone1PerpCone2PtVec, pairPerpCone1PerpCone2EnergyVec, pairPerpCone1PerpCone2ThetaVec, jet.angularity(), jet.ptLeadingConstituent(), splittingMatchesGeoVec, splittingMatchesPtVec, splittingMatchesHFVec, pairMatchesVec); jetMapping.insert(std::make_pair(jet.globalIndex(), jetOutputTable.lastIndex())); } diff --git a/PWGJE/Tasks/substructureDebug.cxx b/PWGJE/Tasks/substructureDebug.cxx index 9a1a2fea361..c6830c3a038 100644 --- a/PWGJE/Tasks/substructureDebug.cxx +++ b/PWGJE/Tasks/substructureDebug.cxx @@ -111,8 +111,9 @@ struct SubstructureDebugTask { registry.fill(HIST("h_D0jet_eta"), jet.jetEta()); registry.fill(HIST("h_D0jet_nTracks"), jet.jetNConstituents()); registry.fill(HIST("h_D0jet_angularity"), jet.angularity()); - auto const& candidate = jet.candidate_as(); - registry.fill(HIST("h_D0injet_pt"), candidate.pt()); + for (auto const& candidate : jet.candidates_as()) { + registry.fill(HIST("h_D0injet_pt"), candidate.pt()); + } } } PROCESS_SWITCH(SubstructureDebugTask, processDataD0, "jets D0 data", true); @@ -125,8 +126,9 @@ struct SubstructureDebugTask { registry.fill(HIST("h_Lcjet_eta"), jet.jetEta()); registry.fill(HIST("h_Lcjet_nTracks"), jet.jetNConstituents()); registry.fill(HIST("h_Lcjet_angularity"), jet.angularity()); - auto const& candidate = jet.candidate_as(); - registry.fill(HIST("h_Lcinjet_pt"), candidate.pt()); + for (auto const& candidate : jet.candidates_as()) { + registry.fill(HIST("h_Lcinjet_pt"), candidate.pt()); + } } } PROCESS_SWITCH(SubstructureDebugTask, processDataLc, "jets Lc data", true); @@ -160,8 +162,9 @@ struct SubstructureDebugTask { for (auto const& mcpjet : matchedJets) { registry.fill(HIST("h_D0jet_pt_matched"), mcpjet.jetPt(), jet.jetPt()); } - auto const& candidate = jet.candidate_as(); - registry.fill(HIST("h_D0injet_pt"), candidate.pt()); + for (auto const& candidate : jet.candidates_as()) { + registry.fill(HIST("h_D0injet_pt"), candidate.pt()); + } } } PROCESS_SWITCH(SubstructureDebugTask, processMCDD0, "jets D0 mcd", false); @@ -178,8 +181,9 @@ struct SubstructureDebugTask { for (auto const& mcpjet : matchedJets) { registry.fill(HIST("h_Lcjet_pt_matched"), mcpjet.jetPt(), jet.jetPt()); } - auto const& candidate = jet.candidate_as(); - registry.fill(HIST("h_Lcinjet_pt"), candidate.pt()); + for (auto const& candidate : jet.candidates_as()) { + registry.fill(HIST("h_Lcinjet_pt"), candidate.pt()); + } } } PROCESS_SWITCH(SubstructureDebugTask, processMCDLc, "jets Lc mcd", false); @@ -205,8 +209,9 @@ struct SubstructureDebugTask { registry.fill(HIST("h_mcpD0jet_eta"), jet.jetEta()); registry.fill(HIST("h_mcpD0jet_nTracks"), jet.jetNConstituents()); registry.fill(HIST("h_mcpD0jet_angularity"), jet.angularity()); - auto const& candidate = jet.candidate_as(); - registry.fill(HIST("h_mcpD0injet_pt"), candidate.pt()); + for (auto const& candidate : jet.candidates_as()) { + registry.fill(HIST("h_mcpD0injet_pt"), candidate.pt()); + } } } PROCESS_SWITCH(SubstructureDebugTask, processMCPD0, "jets D0 mcp", false); @@ -219,8 +224,9 @@ struct SubstructureDebugTask { registry.fill(HIST("h_mcpLcjet_eta"), jet.jetEta()); registry.fill(HIST("h_mcpLcjet_nTracks"), jet.jetNConstituents()); registry.fill(HIST("h_mcpLcjet_angularity"), jet.angularity()); - auto const& candidate = jet.candidate_as(); - registry.fill(HIST("h_mcpLcinjet_pt"), candidate.pt()); + for (auto const& candidate : jet.candidates_as()) { + registry.fill(HIST("h_mcpLcinjet_pt"), candidate.pt()); + } } } PROCESS_SWITCH(SubstructureDebugTask, processMCPLc, "jets Lc mcp", false);