From 9bfc21070cb52e9476906fa2a3968dce3d269712 Mon Sep 17 00:00:00 2001 From: nzardosh Date: Wed, 5 Nov 2025 19:18:03 +0100 Subject: [PATCH 01/13] adding rct analyses to jet framework --- EventFiltering/PWGJE/jetFilter.cxx | 7 ++++--- PWGJE/Core/JetDerivedDataUtilities.h | 8 +++++++- PWGJE/DataModel/JetReducedData.h | 8 ++++++-- PWGJE/TableProducer/derivedDataProducer.cxx | 10 +++++----- PWGJE/TableProducer/derivedDataWriter.cxx | 4 ++-- PWGJE/TableProducer/luminosityProducer.cxx | 2 +- PWGJE/Tasks/jetTriggerChargedQa.cxx | 4 ++-- 7 files changed, 27 insertions(+), 16 deletions(-) diff --git a/EventFiltering/PWGJE/jetFilter.cxx b/EventFiltering/PWGJE/jetFilter.cxx index 277f749d22b..dbc971f6181 100644 --- a/EventFiltering/PWGJE/jetFilter.cxx +++ b/EventFiltering/PWGJE/jetFilter.cxx @@ -33,6 +33,7 @@ #include "PWGJE/Core/JetBkgSubUtils.h" #include "PWGJE/DataModel/EMCALClusters.h" #include "PWGJE/DataModel/Jet.h" +#include "PWGJE/Core/JetDerivedDataUtilities.h" #include "../filterTables.h" @@ -189,7 +190,7 @@ struct jetFilter { // collision process loop bool keepEvent[kTriggerObjects]{false, false, false, false}; - if (!collision.selection_bit(aod::evsel::kNoTimeFrameBorder)) { + if (!jetderiveddatautilities::selectCollision(collision, jetderiveddatautilities::initialiseEventSelectionBits(static_cast("NoTimeFrameBorder")))) { tags(keepEvent[kJetChLowPt], keepEvent[kJetChHighPt], keepEvent[kTrackLowPt], keepEvent[kTrackHighPt]); return; } @@ -287,13 +288,13 @@ struct jetFilter { tags(keepEvent[kJetChLowPt], keepEvent[kJetChHighPt], keepEvent[kTrackLowPt], keepEvent[kTrackHighPt]); } - void processWithoutRho(soa::Join::iterator const& collision, o2::aod::ChargedJets const& jets, soa::Filtered const& tracks) + void processWithoutRho(aod::JetCollision const& collision, o2::aod::ChargedJets const& jets, soa::Filtered const& tracks) { doTriggering(collision, jets, tracks); } PROCESS_SWITCH(jetFilter, processWithoutRho, "Do charged jet triggering without background estimation for filling histograms", true); - void processWithRho(soa::Join::iterator const& collision, o2::aod::ChargedJets const& jets, soa::Filtered const& tracks) + void processWithRho(soa::Join::iterator const& collision, o2::aod::ChargedJets const& jets, soa::Filtered const& tracks) { doTriggering(collision, jets, tracks); } diff --git a/PWGJE/Core/JetDerivedDataUtilities.h b/PWGJE/Core/JetDerivedDataUtilities.h index c3dc36ab18f..ef10d7bd333 100644 --- a/PWGJE/Core/JetDerivedDataUtilities.h +++ b/PWGJE/Core/JetDerivedDataUtilities.h @@ -19,6 +19,7 @@ #include "Common/CCDB/EventSelectionParams.h" #include "Common/CCDB/TriggerAliases.h" +#include "Common/CCDB/RCTSelectionFlags.h" #include @@ -53,11 +54,16 @@ enum JCollisionSubGeneratorId { }; template -bool selectCollision(T const& collision, std::vector eventSelectionMaskBits, bool skipMBGapEvents = true) +bool selectCollision(T const& collision, std::vector eventSelectionMaskBits, bool skipMBGapEvents = true, bool rctSelection = true, std::string rctLabel = "CBT_hadronPID", bool rejectLimitedAcceptanceRct = false, bool requireZDCRct = false) { if (skipMBGapEvents && collision.subGeneratorId() == JCollisionSubGeneratorId::mbGap) { return false; } + o2::aod::rctsel::RCTFlagsChecker rctChecker; + rctChecker.init(rctLabel, requireZDCRct, rejectLimitedAcceptanceRct); + if (rctSelection && !rctChecker.checkTable(collision)) { //CBT_hadronPID given as default so that TOF is included in RCT selection to benefit from better timing for tracks. Impact of this for inclusive jets should be studied + return false; + } if (eventSelectionMaskBits.size() == 0) { return true; } diff --git a/PWGJE/DataModel/JetReducedData.h b/PWGJE/DataModel/JetReducedData.h index f119d843c49..da5d616c9d1 100644 --- a/PWGJE/DataModel/JetReducedData.h +++ b/PWGJE/DataModel/JetReducedData.h @@ -38,6 +38,7 @@ DECLARE_SOA_COLUMN(GlobalBC, globalBC, uint64_t); DECLARE_SOA_COLUMN(Timestamp, timestamp, uint64_t); DECLARE_SOA_BITMAP_COLUMN(Alias, alias, 32); DECLARE_SOA_BITMAP_COLUMN(Selection, selection, 64); +DECLARE_SOA_BITMAP_COLUMN(Rct, rct, 32); DECLARE_SOA_COLUMN(ReadCounts, readCounts, std::vector); DECLARE_SOA_COLUMN(ReadCountsWithTVX, readCountsWithTVX, std::vector); DECLARE_SOA_COLUMN(ReadCountsWithTVXAndNoTFB, readCountsWithTVXAndNoTFB, std::vector); @@ -50,7 +51,8 @@ DECLARE_SOA_TABLE_STAGED(JBCs, "JBC", jbc::GlobalBC, jbc::Timestamp, jbc::Alias, - jbc::Selection); + jbc::Selection, + jbc::Rct); using JBC = JBCs::iterator; using StoredJBC = StoredJBCs::iterator; @@ -91,6 +93,7 @@ DECLARE_SOA_COLUMN(Weight, weight, float); DECLARE_SOA_COLUMN(SubGeneratorId, subGeneratorId, int); DECLARE_SOA_COLUMN(EventSel, eventSel, uint16_t); DECLARE_SOA_BITMAP_COLUMN(Alias, alias, 32); +DECLARE_SOA_BITMAP_COLUMN(Rct, rct, 32); DECLARE_SOA_COLUMN(TrackOccupancyInTimeRange, trackOccupancyInTimeRange, int); DECLARE_SOA_COLUMN(TriggerSel, triggerSel, uint64_t); DECLARE_SOA_COLUMN(ChargedTriggerSel, chargedTriggerSel, uint8_t); @@ -133,8 +136,9 @@ DECLARE_SOA_TABLE_STAGED(JCollisions, "JCOLLISION", jcollision::CentFT0CVariant1, jcollision::HadronicRate, jcollision::TrackOccupancyInTimeRange, - jcollision::EventSel, jcollision::Alias, + jcollision::EventSel, + jcollision::Rct, jcollision::TriggerSel); using JCollision = JCollisions::iterator; diff --git a/PWGJE/TableProducer/derivedDataProducer.cxx b/PWGJE/TableProducer/derivedDataProducer.cxx index 38d0bfd8acd..847fe6619e4 100644 --- a/PWGJE/TableProducer/derivedDataProducer.cxx +++ b/PWGJE/TableProducer/derivedDataProducer.cxx @@ -194,7 +194,7 @@ struct JetDerivedDataProducerTask { void processBunchCrossings(soa::Join::iterator const& bc) { - products.jBCsTable(bc.runNumber(), bc.globalBC(), bc.timestamp(), bc.alias_raw(), bc.selection_raw()); + products.jBCsTable(bc.runNumber(), bc.globalBC(), bc.timestamp(), bc.alias_raw(), bc.selection_raw(), bc.rct_raw()); products.jBCParentIndexTable(bc.globalIndex()); } PROCESS_SWITCH(JetDerivedDataProducerTask, processBunchCrossings, "produces derived bunch crossing table", false); @@ -213,7 +213,7 @@ struct JetDerivedDataProducerTask { triggerDecider.initCCDB(ccdb.service, bc.runNumber(), bc.timestamp(), jetderiveddatautilities::JTriggerMasks); triggerBit = jetderiveddatautilities::setTriggerSelectionBit(triggerDecider.getTriggerOfInterestResults(bc.globalBC())); } - products.jCollisionsTable(collision.posX(), collision.posY(), collision.posZ(), collision.multFV0A(), collision.multFV0C(), collision.multFT0A(), collision.multFT0C(), collision.centFV0A(), -1.0, collision.centFT0A(), collision.centFT0C(), collision.centFT0M(), collision.centFT0CVariant1(), hadronicRate, collision.trackOccupancyInTimeRange(), jetderiveddatautilities::setEventSelectionBit(collision), collision.alias_raw(), triggerBit); // note change multFT0C to multFT0M when problems with multFT0A are fixed + products.jCollisionsTable(collision.posX(), collision.posY(), collision.posZ(), collision.multFV0A(), collision.multFV0C(), collision.multFT0A(), collision.multFT0C(), collision.centFV0A(), -1.0, collision.centFT0A(), collision.centFT0C(), collision.centFT0M(), collision.centFT0CVariant1(), hadronicRate, collision.trackOccupancyInTimeRange(), collision.alias_raw(), jetderiveddatautilities::setEventSelectionBit(collision), collision.rct_raw(), triggerBit); // note change multFT0C to multFT0M when problems with multFT0A are fixed products.jCollisionsParentIndexTable(collision.globalIndex()); products.jCollisionsBunchCrossingIndexTable(collision.bcId()); } @@ -227,7 +227,7 @@ struct JetDerivedDataProducerTask { triggerDecider.initCCDB(ccdb.service, bc.runNumber(), bc.timestamp(), jetderiveddatautilities::JTriggerMasks); triggerBit = jetderiveddatautilities::setTriggerSelectionBit(triggerDecider.getTriggerOfInterestResults(bc.globalBC())); } - products.jCollisionsTable(collision.posX(), collision.posY(), collision.posZ(), -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1, jetderiveddatautilities::setEventSelectionBit(collision), collision.alias_raw(), triggerBit); + products.jCollisionsTable(collision.posX(), collision.posY(), collision.posZ(), -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1, collision.alias_raw(), jetderiveddatautilities::setEventSelectionBit(collision), collision.rct_raw(), triggerBit); products.jCollisionsParentIndexTable(collision.globalIndex()); products.jCollisionsBunchCrossingIndexTable(collision.bcId()); } @@ -235,7 +235,7 @@ struct JetDerivedDataProducerTask { void processCollisionsRun2(soa::Join::iterator const& collision) { - products.jCollisionsTable(collision.posX(), collision.posY(), collision.posZ(), -1.0, -1.0, -1.0, -1.0, collision.centRun2V0A(), collision.centRun2V0M(), -1.0, -1.0, -1.0, -1.0, 1.0, -1, jetderiveddatautilities::setEventSelectionBit(collision), collision.alias_raw(), 0); // note change multFT0C to multFT0M when problems with multFT0A are fixed + products.jCollisionsTable(collision.posX(), collision.posY(), collision.posZ(), -1.0, -1.0, -1.0, -1.0, collision.centRun2V0A(), collision.centRun2V0M(), -1.0, -1.0, -1.0, -1.0, 1.0, -1, collision.alias_raw(), jetderiveddatautilities::setEventSelectionBit(collision), collision.rct_raw(), 0); // note change multFT0C to multFT0M when problems with multFT0A are fixed products.jCollisionsParentIndexTable(collision.globalIndex()); products.jCollisionsBunchCrossingIndexTable(collision.bcId()); } @@ -243,7 +243,7 @@ struct JetDerivedDataProducerTask { void processCollisionsALICE3(aod::Collision const& collision) { - products.jCollisionsTable(collision.posX(), collision.posY(), collision.posZ(), -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1, -1.0, 0, 0); + products.jCollisionsTable(collision.posX(), collision.posY(), collision.posZ(), -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1, -1.0, 0, 0, 0); products.jCollisionsParentIndexTable(collision.globalIndex()); products.jCollisionsBunchCrossingIndexTable(-1); } diff --git a/PWGJE/TableProducer/derivedDataWriter.cxx b/PWGJE/TableProducer/derivedDataWriter.cxx index 3453874c25f..5e876ac03ca 100644 --- a/PWGJE/TableProducer/derivedDataWriter.cxx +++ b/PWGJE/TableProducer/derivedDataWriter.cxx @@ -445,7 +445,7 @@ struct JetDerivedDataWriter { if (collision.isCollisionSelected()) { auto bc = collision.bc_as>(); if (std::find(bcIndicies.begin(), bcIndicies.end(), bc.globalIndex()) == bcIndicies.end()) { - products.storedJBCsTable(bc.runNumber(), bc.globalBC(), bc.timestamp(), bc.alias_raw(), bc.selection_raw()); + products.storedJBCsTable(bc.runNumber(), bc.globalBC(), bc.timestamp(), bc.alias_raw(), bc.selection_raw(), bc.rct_raw()); products.storedJBCParentIndexTable(bc.bcId()); bcIndicies.push_back(bc.globalIndex()); bcMapping[bc.globalIndex()] = products.storedJBCsTable.lastIndex(); @@ -463,7 +463,7 @@ struct JetDerivedDataWriter { for (auto const& collision : collisions) { if (collision.isCollisionSelected()) { - products.storedJCollisionsTable(collision.posX(), collision.posY(), collision.posZ(), collision.multFV0A(), collision.multFV0C(), collision.multFT0A(), collision.multFT0C(), collision.centFV0A(), collision.centFV0M(), collision.centFT0A(), collision.centFT0C(), collision.centFT0M(), collision.centFT0CVariant1(), collision.hadronicRate(), collision.trackOccupancyInTimeRange(), collision.eventSel(), collision.alias_raw(), collision.triggerSel()); + products.storedJCollisionsTable(collision.posX(), collision.posY(), collision.posZ(), collision.multFV0A(), collision.multFV0C(), collision.multFT0A(), collision.multFT0C(), collision.centFV0A(), collision.centFV0M(), collision.centFT0A(), collision.centFT0C(), collision.centFT0M(), collision.centFT0CVariant1(), collision.hadronicRate(), collision.trackOccupancyInTimeRange(), collision.alias_raw(), collision.eventSel(), collision.rct_raw(), collision.triggerSel()); collisionMapping[collision.globalIndex()] = products.storedJCollisionsTable.lastIndex(); products.storedJCollisionMcInfosTable(collision.weight(), collision.subGeneratorId()); products.storedJCollisionsParentIndexTable(collision.collisionId()); diff --git a/PWGJE/TableProducer/luminosityProducer.cxx b/PWGJE/TableProducer/luminosityProducer.cxx index 455ab4e5910..ea81af2f381 100644 --- a/PWGJE/TableProducer/luminosityProducer.cxx +++ b/PWGJE/TableProducer/luminosityProducer.cxx @@ -119,7 +119,7 @@ struct LuminosityProducer { int readCollisionWithCustomCounter = 0; for (const auto& collision : collisions) { readCollisionCounter++; - if (jetderiveddatautilities::selectCollision(collision, jetderiveddatautilities::initialiseEventSelectionBits("TVX"))) { // asuumes all selections include the TVX trigger + if (jetderiveddatautilities::selectCollision(collision, jetderiveddatautilities::initialiseEventSelectionBits("TVX"))) { // asuumes all selections include the TVX trigger and also assumes the default RCT configuration is always used readCollisionWithTVXCounter++; if (std::abs(collision.posZ()) > vertexZCutForCounting) { continue; diff --git a/PWGJE/Tasks/jetTriggerChargedQa.cxx b/PWGJE/Tasks/jetTriggerChargedQa.cxx index f1cb9309cd1..d8b795aa1d6 100644 --- a/PWGJE/Tasks/jetTriggerChargedQa.cxx +++ b/PWGJE/Tasks/jetTriggerChargedQa.cxx @@ -39,7 +39,7 @@ using namespace o2; using namespace o2::framework; using namespace o2::framework::expressions; -using FilteredColl = soa::Filtered>::iterator; +using FilteredColl = soa::Filtered>::iterator; using FilteredJTracks = soa::Filtered>; using FilteredJets = soa::Filtered>; using JoinedTracks = soa::Join; @@ -162,7 +162,7 @@ struct JetTriggerChargedQa { void process(FilteredColl const& collision, FilteredJTracks const& tracks, FilteredJets const& jets, JoinedTracks const&) { - if (!collision.selection_bit(aod::evsel::kNoTimeFrameBorder)) { + if (!jetderiveddatautilities::selectCollision(collision, jetderiveddatautilities::initialiseEventSelectionBits(static_cast("NoTimeFrameBorder")))) { return; } From b9d2e1a0c45c891a8e5f2cd49468562f897f8580 Mon Sep 17 00:00:00 2001 From: nzardosh Date: Fri, 7 Nov 2025 13:07:01 +0100 Subject: [PATCH 02/13] adding process info to jet mc particle tables --- EventFiltering/PWGJE/jetFilter.cxx | 36 ++--- PWGJE/Core/JetDerivedDataUtilities.h | 4 +- PWGJE/DataModel/JetReducedData.h | 36 ++++- PWGJE/DataModel/JetSubtraction.h | 130 ++++++++++++------ PWGJE/TableProducer/derivedDataProducer.cxx | 2 +- PWGJE/TableProducer/derivedDataWriter.cxx | 2 +- .../eventwiseConstituentSubtractor.cxx | 5 +- 7 files changed, 144 insertions(+), 71 deletions(-) diff --git a/EventFiltering/PWGJE/jetFilter.cxx b/EventFiltering/PWGJE/jetFilter.cxx index dbc971f6181..47761d62c1e 100644 --- a/EventFiltering/PWGJE/jetFilter.cxx +++ b/EventFiltering/PWGJE/jetFilter.cxx @@ -11,33 +11,33 @@ // Author: Filip Krizek -#include -#include -#include -#include +#include "../filterTables.h" -#include "Framework/ASoA.h" -#include "Framework/ASoAHelpers.h" -#include "Framework/AnalysisDataModel.h" -#include "Framework/AnalysisTask.h" -#include "Framework/runDataProcessing.h" -#include "ReconstructionDataFormats/Track.h" +#include "PWGJE/Core/FastJetUtilities.h" +#include "PWGJE/Core/JetBkgSubUtils.h" +#include "PWGJE/Core/JetDerivedDataUtilities.h" +#include "PWGJE/Core/JetFinder.h" +#include "PWGJE/DataModel/EMCALClusters.h" +#include "PWGJE/DataModel/Jet.h" #include "Common/Core/TrackSelection.h" #include "Common/Core/TrackSelectionDefaults.h" #include "Common/DataModel/EventSelection.h" #include "Common/DataModel/TrackSelectionTables.h" -#include "PWGJE/Core/JetFinder.h" -#include "PWGJE/Core/FastJetUtilities.h" -#include "PWGJE/Core/JetBkgSubUtils.h" -#include "PWGJE/DataModel/EMCALClusters.h" -#include "PWGJE/DataModel/Jet.h" -#include "PWGJE/Core/JetDerivedDataUtilities.h" +#include "Framework/ASoA.h" +#include "Framework/ASoAHelpers.h" +#include "Framework/AnalysisDataModel.h" +#include "Framework/AnalysisTask.h" +#include "Framework/HistogramRegistry.h" +#include "Framework/runDataProcessing.h" +#include "ReconstructionDataFormats/Track.h" -#include "../filterTables.h" +#include -#include "Framework/HistogramRegistry.h" +#include +#include +#include using namespace o2; using namespace o2::framework; diff --git a/PWGJE/Core/JetDerivedDataUtilities.h b/PWGJE/Core/JetDerivedDataUtilities.h index ef10d7bd333..921959742c4 100644 --- a/PWGJE/Core/JetDerivedDataUtilities.h +++ b/PWGJE/Core/JetDerivedDataUtilities.h @@ -18,8 +18,8 @@ #define PWGJE_CORE_JETDERIVEDDATAUTILITIES_H_ #include "Common/CCDB/EventSelectionParams.h" -#include "Common/CCDB/TriggerAliases.h" #include "Common/CCDB/RCTSelectionFlags.h" +#include "Common/CCDB/TriggerAliases.h" #include @@ -61,7 +61,7 @@ bool selectCollision(T const& collision, std::vector eventSelectionMaskBits } o2::aod::rctsel::RCTFlagsChecker rctChecker; rctChecker.init(rctLabel, requireZDCRct, rejectLimitedAcceptanceRct); - if (rctSelection && !rctChecker.checkTable(collision)) { //CBT_hadronPID given as default so that TOF is included in RCT selection to benefit from better timing for tracks. Impact of this for inclusive jets should be studied + if (rctSelection && !rctChecker.checkTable(collision)) { // CBT_hadronPID given as default so that TOF is included in RCT selection to benefit from better timing for tracks. Impact of this for inclusive jets should be studied return false; } if (eventSelectionMaskBits.size() == 0) { diff --git a/PWGJE/DataModel/JetReducedData.h b/PWGJE/DataModel/JetReducedData.h index da5d616c9d1..846a9388900 100644 --- a/PWGJE/DataModel/JetReducedData.h +++ b/PWGJE/DataModel/JetReducedData.h @@ -340,11 +340,14 @@ DECLARE_SOA_COLUMN(Phi, phi, float); DECLARE_SOA_COLUMN(Y, y, float); DECLARE_SOA_COLUMN(E, e, float); DECLARE_SOA_COLUMN(PdgCode, pdgCode, int); -DECLARE_SOA_COLUMN(GenStatusCode, getGenStatusCode, int); // TODO : We can look at combining this with the two below -DECLARE_SOA_COLUMN(HepMCStatusCode, getHepMCStatusCode, int); -DECLARE_SOA_COLUMN(IsPhysicalPrimary, isPhysicalPrimary, bool); +DECLARE_SOA_COLUMN(StatusCode, statusCode, int); +DECLARE_SOA_COLUMN(Flags, flags, uint8_t); + DECLARE_SOA_SELF_ARRAY_INDEX_COLUMN(Mothers, mothers); DECLARE_SOA_SELF_SLICE_INDEX_COLUMN(Daughters, daughters); + +DECLARE_SOA_DYNAMIC_COLUMN(PVector, pVector, //! Momentum vector in x,y,z-directions in GeV/c + [](float px, float py, float pz) -> std::array { return std::array{px, py, pz}; }); DECLARE_SOA_DYNAMIC_COLUMN(Px, px, [](float pt, float phi) -> float { return pt * std::cos(phi); }); DECLARE_SOA_DYNAMIC_COLUMN(Py, py, @@ -355,6 +358,20 @@ DECLARE_SOA_DYNAMIC_COLUMN(P, p, [](float pt, float eta) -> float { return pt * std::cosh(eta); }); DECLARE_SOA_DYNAMIC_COLUMN(Energy, energy, [](float e) -> float { return e; }); + +DECLARE_SOA_DYNAMIC_COLUMN(ProducedByGenerator, producedByGenerator, //! True if particle produced by the generator (==TMCProcess::kPrimary); False if by the transport code + [](uint8_t flags) -> bool { return (flags & o2::aod::mcparticle::enums::ProducedByTransport) == 0x0; }); +DECLARE_SOA_DYNAMIC_COLUMN(FromBackgroundEvent, fromBackgroundEvent, //! Particle from background event + [](uint8_t flags) -> bool { return (flags & o2::aod::mcparticle::enums::FromBackgroundEvent) == o2::aod::mcparticle::enums::FromBackgroundEvent; }); +DECLARE_SOA_DYNAMIC_COLUMN(GetProcess, getProcess, //! The VMC physics code (as int) that generated this particle (see header TMCProcess.h in ROOT) + [](uint8_t flags, int statusCode) -> int { if ((flags & o2::aod::mcparticle::enums::ProducedByTransport) == 0x0) { return 0 /*TMCProcess::kPrimary*/; } else { return statusCode; } }); +DECLARE_SOA_DYNAMIC_COLUMN(GetGenStatusCode, getGenStatusCode, //! The native status code put by the generator, or -1 if a particle produced during transport + [](uint8_t flags, int statusCode) -> int { if ((flags & o2::aod::mcparticle::enums::ProducedByTransport) == 0x0) { return o2::mcgenstatus::getGenStatusCode(statusCode); } else { return -1; } }); +DECLARE_SOA_DYNAMIC_COLUMN(GetHepMCStatusCode, getHepMCStatusCode, //! The HepMC status code put by the generator, or -1 if a particle produced during transport + [](uint8_t flags, int statusCode) -> int { if ((flags & o2::aod::mcparticle::enums::ProducedByTransport) == 0x0) { return o2::mcgenstatus::getHepMCStatusCode(statusCode); } else { return -1; } }); +DECLARE_SOA_DYNAMIC_COLUMN(IsPhysicalPrimary, isPhysicalPrimary, //! True if particle is considered a physical primary according to the ALICE definition + [](uint8_t flags) -> bool { return (flags & o2::aod::mcparticle::enums::PhysicalPrimary) == o2::aod::mcparticle::enums::PhysicalPrimary; }); + } // namespace jmcparticle DECLARE_SOA_TABLE_STAGED(JMcParticles, "JMCPARTICLE", @@ -366,16 +383,21 @@ DECLARE_SOA_TABLE_STAGED(JMcParticles, "JMCPARTICLE", jmcparticle::Y, jmcparticle::E, jmcparticle::PdgCode, - jmcparticle::GenStatusCode, - jmcparticle::HepMCStatusCode, - jmcparticle::IsPhysicalPrimary, + jmcparticle::StatusCode, + jmcparticle::Flags, jmcparticle::MothersIds, jmcparticle::DaughtersIdSlice, jmcparticle::Px, jmcparticle::Py, jmcparticle::Pz, jmcparticle::P, - jmcparticle::Energy); + jmcparticle::Energy, + jmcparticle::ProducedByGenerator, + jmcparticle::FromBackgroundEvent, + jmcparticle::GetProcess, + jmcparticle::GetGenStatusCode, + jmcparticle::GetHepMCStatusCode, + jmcparticle::IsPhysicalPrimary); using JMcParticle = JMcParticles::iterator; using StoredJMcParticle = StoredJMcParticles::iterator; diff --git a/PWGJE/DataModel/JetSubtraction.h b/PWGJE/DataModel/JetSubtraction.h index 0635ff8cb77..0ca981c6e04 100644 --- a/PWGJE/DataModel/JetSubtraction.h +++ b/PWGJE/DataModel/JetSubtraction.h @@ -275,14 +275,19 @@ DECLARE_SOA_TABLE(JMcParticleSubs, "AOD", "JMcPartSubs", jmcparticle::Y, jmcparticle::E, jmcparticle::PdgCode, - jmcparticle::GenStatusCode, - jmcparticle::HepMCStatusCode, - jmcparticle::IsPhysicalPrimary, + jmcparticle::StatusCode, + jmcparticle::Flags, jmcparticle::Px, jmcparticle::Py, jmcparticle::Pz, jmcparticle::P, - jmcparticle::Energy); + jmcparticle::Energy, + jmcparticle::ProducedByGenerator, // this will give a nonsensical value and should not be used + jmcparticle::FromBackgroundEvent, // this will give a nonsensical value and should not be used + jmcparticle::GetProcess, // this will give a nonsensical value and should not be used + jmcparticle::GetGenStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::GetHepMCStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::IsPhysicalPrimary); using JMcParticleSub = JMcParticleSubs::iterator; @@ -310,14 +315,19 @@ DECLARE_SOA_TABLE(JMcParticleD0Subs, "AOD", "JMcPartD0Subs", jmcparticle::Y, jmcparticle::E, jmcparticle::PdgCode, - jmcparticle::GenStatusCode, - jmcparticle::HepMCStatusCode, - jmcparticle::IsPhysicalPrimary, + jmcparticle::StatusCode, + jmcparticle::Flags, jmcparticle::Px, jmcparticle::Py, jmcparticle::Pz, jmcparticle::P, - jmcparticle::Energy); + jmcparticle::Energy, + jmcparticle::ProducedByGenerator, // this will give a nonsensical value and should not be used + jmcparticle::FromBackgroundEvent, // this will give a nonsensical value and should not be used + jmcparticle::GetProcess, // this will give a nonsensical value and should not be used + jmcparticle::GetGenStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::GetHepMCStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::IsPhysicalPrimary); using JMcParticleD0Sub = JMcParticleD0Subs::iterator; @@ -345,14 +355,19 @@ DECLARE_SOA_TABLE(JMcParticleDplusSubs, "AOD", "JMcPartDPSubs", jmcparticle::Y, jmcparticle::E, jmcparticle::PdgCode, - jmcparticle::GenStatusCode, - jmcparticle::HepMCStatusCode, - jmcparticle::IsPhysicalPrimary, + jmcparticle::StatusCode, + jmcparticle::Flags, jmcparticle::Px, jmcparticle::Py, jmcparticle::Pz, jmcparticle::P, - jmcparticle::Energy); + jmcparticle::Energy, + jmcparticle::ProducedByGenerator, // this will give a nonsensical value and should not be used + jmcparticle::FromBackgroundEvent, // this will give a nonsensical value and should not be used + jmcparticle::GetProcess, // this will give a nonsensical value and should not be used + jmcparticle::GetGenStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::GetHepMCStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::IsPhysicalPrimary); using JMcParticleDplusSub = JMcParticleDplusSubs::iterator; @@ -380,14 +395,19 @@ DECLARE_SOA_TABLE(JMcParticleDsSubs, "AOD", "JMcPartDSSubs", jmcparticle::Y, jmcparticle::E, jmcparticle::PdgCode, - jmcparticle::GenStatusCode, - jmcparticle::HepMCStatusCode, - jmcparticle::IsPhysicalPrimary, + jmcparticle::StatusCode, + jmcparticle::Flags, jmcparticle::Px, jmcparticle::Py, jmcparticle::Pz, jmcparticle::P, - jmcparticle::Energy); + jmcparticle::Energy, + jmcparticle::ProducedByGenerator, // this will give a nonsensical value and should not be used + jmcparticle::FromBackgroundEvent, // this will give a nonsensical value and should not be used + jmcparticle::GetProcess, // this will give a nonsensical value and should not be used + jmcparticle::GetGenStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::GetHepMCStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::IsPhysicalPrimary); using JMcParticleDsSub = JMcParticleDsSubs::iterator; @@ -415,14 +435,19 @@ DECLARE_SOA_TABLE(JMcParticleDstarSubs, "AOD", "JMcPartDSTSubs", jmcparticle::Y, jmcparticle::E, jmcparticle::PdgCode, - jmcparticle::GenStatusCode, - jmcparticle::HepMCStatusCode, - jmcparticle::IsPhysicalPrimary, + jmcparticle::StatusCode, + jmcparticle::Flags, jmcparticle::Px, jmcparticle::Py, jmcparticle::Pz, jmcparticle::P, - jmcparticle::Energy); + jmcparticle::Energy, + jmcparticle::ProducedByGenerator, // this will give a nonsensical value and should not be used + jmcparticle::FromBackgroundEvent, // this will give a nonsensical value and should not be used + jmcparticle::GetProcess, // this will give a nonsensical value and should not be used + jmcparticle::GetGenStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::GetHepMCStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::IsPhysicalPrimary); using JMcParticleDstarSub = JMcParticleDstarSubs::iterator; @@ -450,14 +475,19 @@ DECLARE_SOA_TABLE(JMcParticleLcSubs, "AOD", "JMcPartLCSubs", jmcparticle::Y, jmcparticle::E, jmcparticle::PdgCode, - jmcparticle::GenStatusCode, - jmcparticle::HepMCStatusCode, - jmcparticle::IsPhysicalPrimary, + jmcparticle::StatusCode, + jmcparticle::Flags, jmcparticle::Px, jmcparticle::Py, jmcparticle::Pz, jmcparticle::P, - jmcparticle::Energy); + jmcparticle::Energy, + jmcparticle::ProducedByGenerator, // this will give a nonsensical value and should not be used + jmcparticle::FromBackgroundEvent, // this will give a nonsensical value and should not be used + jmcparticle::GetProcess, // this will give a nonsensical value and should not be used + jmcparticle::GetGenStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::GetHepMCStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::IsPhysicalPrimary); using JMcParticleLcSub = JMcParticleLcSubs::iterator; @@ -485,14 +515,19 @@ DECLARE_SOA_TABLE(JMcParticleB0Subs, "AOD", "JMcPartB0Subs", jmcparticle::Y, jmcparticle::E, jmcparticle::PdgCode, - jmcparticle::GenStatusCode, - jmcparticle::HepMCStatusCode, - jmcparticle::IsPhysicalPrimary, + jmcparticle::StatusCode, + jmcparticle::Flags, jmcparticle::Px, jmcparticle::Py, jmcparticle::Pz, jmcparticle::P, - jmcparticle::Energy); + jmcparticle::Energy, + jmcparticle::ProducedByGenerator, // this will give a nonsensical value and should not be used + jmcparticle::FromBackgroundEvent, // this will give a nonsensical value and should not be used + jmcparticle::GetProcess, // this will give a nonsensical value and should not be used + jmcparticle::GetGenStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::GetHepMCStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::IsPhysicalPrimary); using JMcParticleB0Sub = JMcParticleB0Subs::iterator; @@ -520,14 +555,19 @@ DECLARE_SOA_TABLE(JMcParticleBplusSubs, "AOD", "JMcPartBPSubs", jmcparticle::Y, jmcparticle::E, jmcparticle::PdgCode, - jmcparticle::GenStatusCode, - jmcparticle::HepMCStatusCode, - jmcparticle::IsPhysicalPrimary, + jmcparticle::StatusCode, + jmcparticle::Flags, jmcparticle::Px, jmcparticle::Py, jmcparticle::Pz, jmcparticle::P, - jmcparticle::Energy); + jmcparticle::Energy, + jmcparticle::ProducedByGenerator, // this will give a nonsensical value and should not be used + jmcparticle::FromBackgroundEvent, // this will give a nonsensical value and should not be used + jmcparticle::GetProcess, // this will give a nonsensical value and should not be used + jmcparticle::GetGenStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::GetHepMCStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::IsPhysicalPrimary); using JMcParticleBplusSub = JMcParticleBplusSubs::iterator; @@ -555,14 +595,19 @@ DECLARE_SOA_TABLE(JMcParticleXicToXiPiPiSubs, "AOD", "JMcPartXICXPPCSubs", jmcparticle::Y, jmcparticle::E, jmcparticle::PdgCode, - jmcparticle::GenStatusCode, - jmcparticle::HepMCStatusCode, - jmcparticle::IsPhysicalPrimary, + jmcparticle::StatusCode, + jmcparticle::Flags, jmcparticle::Px, jmcparticle::Py, jmcparticle::Pz, jmcparticle::P, - jmcparticle::Energy); + jmcparticle::Energy, + jmcparticle::ProducedByGenerator, // this will give a nonsensical value and should not be used + jmcparticle::FromBackgroundEvent, // this will give a nonsensical value and should not be used + jmcparticle::GetProcess, // this will give a nonsensical value and should not be used + jmcparticle::GetGenStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::GetHepMCStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::IsPhysicalPrimary); using JMcParticleXicToXiPiPiSub = JMcParticleXicToXiPiPiSubs::iterator; @@ -590,14 +635,19 @@ DECLARE_SOA_TABLE(JMcParticleDielectronSubs, "AOD", "JMcPartDIELSubs", jmcparticle::Y, jmcparticle::E, jmcparticle::PdgCode, - jmcparticle::GenStatusCode, - jmcparticle::HepMCStatusCode, - jmcparticle::IsPhysicalPrimary, + jmcparticle::StatusCode, + jmcparticle::Flags, jmcparticle::Px, jmcparticle::Py, jmcparticle::Pz, jmcparticle::P, - jmcparticle::Energy); + jmcparticle::Energy, + jmcparticle::ProducedByGenerator, // this will give a nonsensical value and should not be used + jmcparticle::FromBackgroundEvent, // this will give a nonsensical value and should not be used + jmcparticle::GetProcess, // this will give a nonsensical value and should not be used + jmcparticle::GetGenStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::GetHepMCStatusCode, // this will give a nonsensical value and should not be used + jmcparticle::IsPhysicalPrimary); using JMcParticleDielectronSub = JMcParticleDielectronSubs::iterator; diff --git a/PWGJE/TableProducer/derivedDataProducer.cxx b/PWGJE/TableProducer/derivedDataProducer.cxx index 847fe6619e4..2bf830c1a12 100644 --- a/PWGJE/TableProducer/derivedDataProducer.cxx +++ b/PWGJE/TableProducer/derivedDataProducer.cxx @@ -429,7 +429,7 @@ struct JetDerivedDataProducerTask { i++; } } - products.jMcParticlesTable(particle.mcCollisionId(), particle.pt(), particle.eta(), particle.phi(), particle.y(), particle.e(), particle.pdgCode(), particle.getGenStatusCode(), particle.getHepMCStatusCode(), particle.isPhysicalPrimary(), mothersId, daughtersId); + products.jMcParticlesTable(particle.mcCollisionId(), particle.pt(), particle.eta(), particle.phi(), particle.y(), particle.e(), particle.pdgCode(), particle.flags(), particle.statusCode(), mothersId, daughtersId); products.jParticlesParentIndexTable(particle.globalIndex()); } PROCESS_SWITCH(JetDerivedDataProducerTask, processParticles, "produces derived parrticle table", false); diff --git a/PWGJE/TableProducer/derivedDataWriter.cxx b/PWGJE/TableProducer/derivedDataWriter.cxx index 5e876ac03ca..5b98f81fd68 100644 --- a/PWGJE/TableProducer/derivedDataWriter.cxx +++ b/PWGJE/TableProducer/derivedDataWriter.cxx @@ -682,7 +682,7 @@ struct JetDerivedDataWriter { i++; } } - products.storedJMcParticlesTable(mcCollisionMapping[mcCollision.globalIndex()], o2::math_utils::detail::truncateFloatFraction(particle.pt(), precisionMomentumMask), o2::math_utils::detail::truncateFloatFraction(particle.eta(), precisionPositionMask), o2::math_utils::detail::truncateFloatFraction(particle.phi(), precisionPositionMask), o2::math_utils::detail::truncateFloatFraction(particle.y(), precisionPositionMask), o2::math_utils::detail::truncateFloatFraction(particle.e(), precisionMomentumMask), particle.pdgCode(), particle.getGenStatusCode(), particle.getHepMCStatusCode(), particle.isPhysicalPrimary(), mothersIds, daughtersIds); + products.storedJMcParticlesTable(mcCollisionMapping[mcCollision.globalIndex()], o2::math_utils::detail::truncateFloatFraction(particle.pt(), precisionMomentumMask), o2::math_utils::detail::truncateFloatFraction(particle.eta(), precisionPositionMask), o2::math_utils::detail::truncateFloatFraction(particle.phi(), precisionPositionMask), o2::math_utils::detail::truncateFloatFraction(particle.y(), precisionPositionMask), o2::math_utils::detail::truncateFloatFraction(particle.e(), precisionMomentumMask), particle.pdgCode(), particle.statusCode(), particle.flags(), mothersIds, daughtersIds); products.storedJParticlesParentIndexTable(particle.mcParticleId()); } } diff --git a/PWGJE/TableProducer/eventwiseConstituentSubtractor.cxx b/PWGJE/TableProducer/eventwiseConstituentSubtractor.cxx index 254421574ab..720d4081926 100644 --- a/PWGJE/TableProducer/eventwiseConstituentSubtractor.cxx +++ b/PWGJE/TableProducer/eventwiseConstituentSubtractor.cxx @@ -25,6 +25,7 @@ #include "Framework/O2DatabasePDGPlugin.h" #include #include +#include #include #include #include @@ -137,7 +138,7 @@ struct eventWiseConstituentSubtractorTask { tracksSubtracted = eventWiseConstituentSubtractor.JetBkgSubUtils::doEventConstSub(inputParticles, candidate.rho(), candidate.rhoM()); for (auto const& trackSubtracted : tracksSubtracted) { - particleSubTable(candidate.globalIndex(), trackSubtracted.pt(), trackSubtracted.eta(), trackSubtracted.phi(), trackSubtracted.rap(), trackSubtracted.e(), 211, 1, 1, 1); // everything after phi is artificial and should not be used for analyses + particleSubTable(candidate.globalIndex(), trackSubtracted.pt(), trackSubtracted.eta(), trackSubtracted.phi(), trackSubtracted.rap(), trackSubtracted.e(), 211, 0, static_cast(o2::aod::mcparticle::enums::PhysicalPrimary)); // everything after phi is artificial and should not be used for analyses } } } @@ -171,7 +172,7 @@ struct eventWiseConstituentSubtractorTask { tracksSubtracted = eventWiseConstituentSubtractor.JetBkgSubUtils::doEventConstSub(inputParticles, mcCollision.rho(), mcCollision.rhoM()); for (auto const& trackSubtracted : tracksSubtracted) { - particleSubtractedTable(mcCollision.globalIndex(), trackSubtracted.pt(), trackSubtracted.eta(), trackSubtracted.phi(), trackSubtracted.rap(), trackSubtracted.e(), 211, 1, 1, 1); // everything after phi is artificial and should not be used for analyses + particleSubtractedTable(mcCollision.globalIndex(), trackSubtracted.pt(), trackSubtracted.eta(), trackSubtracted.phi(), trackSubtracted.rap(), trackSubtracted.e(), 211, 0, static_cast(o2::aod::mcparticle::enums::PhysicalPrimary)); // everything after phi is artificial and should not be used for analyses } } PROCESS_SWITCH(eventWiseConstituentSubtractorTask, processMcCollisions, "Fill table of subtracted tracks for Mc collisions", false); From eeb1d7766876eca52a667839059c4ccf78fa7172 Mon Sep 17 00:00:00 2001 From: nzardosh Date: Fri, 7 Nov 2025 13:10:16 +0100 Subject: [PATCH 03/13] remove area calculation for substructure reclustering --- PWGJE/Tasks/jetSubstructure.cxx | 1 + PWGJE/Tasks/jetSubstructureHF.cxx | 1 + 2 files changed, 2 insertions(+) diff --git a/PWGJE/Tasks/jetSubstructure.cxx b/PWGJE/Tasks/jetSubstructure.cxx index 9980893e4e5..92ce238721b 100644 --- a/PWGJE/Tasks/jetSubstructure.cxx +++ b/PWGJE/Tasks/jetSubstructure.cxx @@ -119,6 +119,7 @@ struct JetSubstructureTask { jetReclusterer.isReclustering = true; jetReclusterer.algorithm = fastjet::JetAlgorithm::cambridge_algorithm; + jetReclusterer.ghostRepeatN = 0; } Preslice TracksPerCollision = aod::jtrack::collisionId; diff --git a/PWGJE/Tasks/jetSubstructureHF.cxx b/PWGJE/Tasks/jetSubstructureHF.cxx index 1117a79bd37..3693d0908ff 100644 --- a/PWGJE/Tasks/jetSubstructureHF.cxx +++ b/PWGJE/Tasks/jetSubstructureHF.cxx @@ -123,6 +123,7 @@ struct JetSubstructureHFTask { jetReclusterer.isReclustering = true; jetReclusterer.algorithm = fastjet::JetAlgorithm::cambridge_algorithm; + jetReclusterer.ghostRepeatN = 0; candMass = jetcandidateutilities::getTablePDGMass(); } From 236c5833fdcb7bba22adc111c773ee5b3799dced Mon Sep 17 00:00:00 2001 From: nzardosh Date: Fri, 7 Nov 2025 13:38:28 +0100 Subject: [PATCH 04/13] removing unfilled mc centrality tables --- PWGJE/DataModel/JetReducedData.h | 6 ------ PWGJE/TableProducer/derivedDataProducer.cxx | 12 ++++++------ PWGJE/TableProducer/derivedDataWriter.cxx | 2 +- PWGJE/Tasks/jetSpectraCharged.cxx | 3 ++- 4 files changed, 9 insertions(+), 14 deletions(-) diff --git a/PWGJE/DataModel/JetReducedData.h b/PWGJE/DataModel/JetReducedData.h index 846a9388900..900abf07396 100644 --- a/PWGJE/DataModel/JetReducedData.h +++ b/PWGJE/DataModel/JetReducedData.h @@ -196,9 +196,6 @@ DECLARE_SOA_COLUMN(PosZ, posZ, float); DECLARE_SOA_COLUMN(MultFV0A, multFV0A, float); DECLARE_SOA_COLUMN(MultFT0A, multFT0A, float); DECLARE_SOA_COLUMN(MultFT0C, multFT0C, float); -DECLARE_SOA_COLUMN(CentFV0A, centFV0A, float); -DECLARE_SOA_COLUMN(CentFT0A, centFT0A, float); -DECLARE_SOA_COLUMN(CentFT0C, centFT0C, float); DECLARE_SOA_COLUMN(CentFT0M, centFT0M, float); DECLARE_SOA_COLUMN(Weight, weight, float); DECLARE_SOA_COLUMN(SubGeneratorId, subGeneratorId, int); @@ -217,9 +214,6 @@ DECLARE_SOA_TABLE_STAGED(JMcCollisions, "JMCCOLLISION", jmccollision::MultFV0A, jmccollision::MultFT0A, jmccollision::MultFT0C, - jmccollision::CentFV0A, - jmccollision::CentFT0A, - jmccollision::CentFT0C, jmccollision::CentFT0M, jmccollision::Weight, jmccollision::SubGeneratorId, diff --git a/PWGJE/TableProducer/derivedDataProducer.cxx b/PWGJE/TableProducer/derivedDataProducer.cxx index 2bf830c1a12..4091715b441 100644 --- a/PWGJE/TableProducer/derivedDataProducer.cxx +++ b/PWGJE/TableProducer/derivedDataProducer.cxx @@ -272,30 +272,30 @@ struct JetDerivedDataProducerTask { } PROCESS_SWITCH(JetDerivedDataProducerTask, processMcCollisionLabels, "produces derived MC collision labels table", false); - void processMcCollisions(soa::Join::iterator const& mcCollision) + void processMcCollisions(soa::Join::iterator const& mcCollision) { - products.jMcCollisionsTable(mcCollision.posX(), mcCollision.posY(), mcCollision.posZ(), mcCollision.multMCFV0A(), mcCollision.multMCFT0A(), mcCollision.multMCFT0C(), mcCollision.centFV0A(), mcCollision.centFT0A(), mcCollision.centFT0C(), mcCollision.centFT0M(), mcCollision.weight(), mcCollision.getSubGeneratorId(), mcCollision.accepted(), mcCollision.attempted(), mcCollision.xsectGen(), mcCollision.xsectErr(), mcCollision.ptHard()); + products.jMcCollisionsTable(mcCollision.posX(), mcCollision.posY(), mcCollision.posZ(), mcCollision.multMCFV0A(), mcCollision.multMCFT0A(), mcCollision.multMCFT0C(), mcCollision.centFT0M(), mcCollision.weight(), mcCollision.getSubGeneratorId(), mcCollision.accepted(), mcCollision.attempted(), mcCollision.xsectGen(), mcCollision.xsectErr(), mcCollision.ptHard()); products.jMcCollisionsParentIndexTable(mcCollision.globalIndex()); } PROCESS_SWITCH(JetDerivedDataProducerTask, processMcCollisions, "produces derived MC collision table", false); void processMcCollisionsWithoutCentralityAndMultiplicity(soa::Join::iterator const& mcCollision) { - products.jMcCollisionsTable(mcCollision.posX(), mcCollision.posY(), mcCollision.posZ(), -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, mcCollision.weight(), mcCollision.getSubGeneratorId(), mcCollision.accepted(), mcCollision.attempted(), mcCollision.xsectGen(), mcCollision.xsectErr(), mcCollision.ptHard()); + products.jMcCollisionsTable(mcCollision.posX(), mcCollision.posY(), mcCollision.posZ(), -1.0, -1.0, -1.0, -1.0, mcCollision.weight(), mcCollision.getSubGeneratorId(), mcCollision.accepted(), mcCollision.attempted(), mcCollision.xsectGen(), mcCollision.xsectErr(), mcCollision.ptHard()); products.jMcCollisionsParentIndexTable(mcCollision.globalIndex()); } PROCESS_SWITCH(JetDerivedDataProducerTask, processMcCollisionsWithoutCentralityAndMultiplicity, "produces derived MC collision table without centraility and multiplicity", false); - void processMcCollisionsWithoutXsection(soa::Join::iterator const& mcCollision) + void processMcCollisionsWithoutXsection(soa::Join::iterator const& mcCollision) { - products.jMcCollisionsTable(mcCollision.posX(), mcCollision.posY(), mcCollision.posZ(), mcCollision.multMCFV0A(), mcCollision.multMCFT0A(), mcCollision.multMCFT0C(), mcCollision.centFV0A(), mcCollision.centFT0A(), mcCollision.centFT0C(), mcCollision.centFT0M(), mcCollision.weight(), mcCollision.getSubGeneratorId(), 1, 1, 1.0, 1.0, 999.0); + products.jMcCollisionsTable(mcCollision.posX(), mcCollision.posY(), mcCollision.posZ(), mcCollision.multMCFV0A(), mcCollision.multMCFT0A(), mcCollision.multMCFT0C(), mcCollision.centFT0M(), mcCollision.weight(), mcCollision.getSubGeneratorId(), 1, 1, 1.0, 1.0, 999.0); products.jMcCollisionsParentIndexTable(mcCollision.globalIndex()); } PROCESS_SWITCH(JetDerivedDataProducerTask, processMcCollisionsWithoutXsection, "produces derived MC collision table without cross section information", false); void processMcCollisionsWithoutCentralityAndMultiplicityAndXsection(aod::McCollision const& mcCollision) { - products.jMcCollisionsTable(mcCollision.posX(), mcCollision.posY(), mcCollision.posZ(), -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, -1.0, mcCollision.weight(), mcCollision.getSubGeneratorId(), 1, 1, 1.0, 1.0, 999.0); + products.jMcCollisionsTable(mcCollision.posX(), mcCollision.posY(), mcCollision.posZ(), -1.0, -1.0, -1.0, -1.0, mcCollision.weight(), mcCollision.getSubGeneratorId(), 1, 1, 1.0, 1.0, 999.0); products.jMcCollisionsParentIndexTable(mcCollision.globalIndex()); } PROCESS_SWITCH(JetDerivedDataProducerTask, processMcCollisionsWithoutCentralityAndMultiplicityAndXsection, "produces derived MC collision table without centrality, multiplicity and cross section information", false); diff --git a/PWGJE/TableProducer/derivedDataWriter.cxx b/PWGJE/TableProducer/derivedDataWriter.cxx index 5b98f81fd68..72864b14e58 100644 --- a/PWGJE/TableProducer/derivedDataWriter.cxx +++ b/PWGJE/TableProducer/derivedDataWriter.cxx @@ -640,7 +640,7 @@ struct JetDerivedDataWriter { mcCollisionMapping.resize(mcCollisions.size(), -1); for (auto const& mcCollision : mcCollisions) { if (mcCollision.isMcCollisionSelected()) { - products.storedJMcCollisionsTable(mcCollision.posX(), mcCollision.posY(), mcCollision.posZ(), mcCollision.multFV0A(), mcCollision.multFT0A(), mcCollision.multFT0C(), mcCollision.centFV0A(), mcCollision.centFT0A(), mcCollision.centFT0C(), mcCollision.centFT0M(), mcCollision.weight(), mcCollision.subGeneratorId(), mcCollision.accepted(), mcCollision.attempted(), mcCollision.xsectGen(), mcCollision.xsectErr(), mcCollision.ptHard()); + products.storedJMcCollisionsTable(mcCollision.posX(), mcCollision.posY(), mcCollision.posZ(), mcCollision.multFV0A(), mcCollision.multFT0A(), mcCollision.multFT0C(), mcCollision.centFT0M(), mcCollision.weight(), mcCollision.subGeneratorId(), mcCollision.accepted(), mcCollision.attempted(), mcCollision.xsectGen(), mcCollision.xsectErr(), mcCollision.ptHard()); products.storedJMcCollisionsParentIndexTable(mcCollision.mcCollisionId()); mcCollisionMapping[mcCollision.globalIndex()] = products.storedJMcCollisionsTable.lastIndex(); } diff --git a/PWGJE/Tasks/jetSpectraCharged.cxx b/PWGJE/Tasks/jetSpectraCharged.cxx index 45400c66ab8..90c75369ff7 100644 --- a/PWGJE/Tasks/jetSpectraCharged.cxx +++ b/PWGJE/Tasks/jetSpectraCharged.cxx @@ -303,7 +303,8 @@ struct JetSpectraCharged { bool applyMCCollisionCuts(TMCColl const& mccollision, TCollisions const& collisions, bool fillHistograms = false, bool isWeighted = false, float eventWeight = 1.0) { float centrality = -1.0; - checkCentFT0M ? centrality = mccollision.centFT0M() : centrality = mccollision.centFT0C(); + // checkCentFT0M ? centrality = mccollision.centFT0M() : centrality = mccollision.centFT0C(); + centrality = mccollision.centFT0M(); if (fillHistograms) { registry.fill(HIST("h_mccollisions"), 0.5); From 7bab3a6fff5ba83b6a14ae4073de1c14af80d3ac Mon Sep 17 00:00:00 2001 From: nzardosh Date: Sun, 9 Nov 2025 22:16:52 +0100 Subject: [PATCH 05/13] adding online mode to trigger correlation task --- PWGJE/Tasks/triggerCorrelations.cxx | 40 ++++++++++++++++++++++++----- 1 file changed, 33 insertions(+), 7 deletions(-) diff --git a/PWGJE/Tasks/triggerCorrelations.cxx b/PWGJE/Tasks/triggerCorrelations.cxx index 2d7b54707bb..73acb90145e 100644 --- a/PWGJE/Tasks/triggerCorrelations.cxx +++ b/PWGJE/Tasks/triggerCorrelations.cxx @@ -38,6 +38,9 @@ struct TriggerCorrelationsTask { HistogramRegistry registry; std::vector triggerMaskBits; + long unsigned int nChargedTriggers = 4; + long unsigned int nChargedHFTriggers = 4; + long unsigned int nFullTriggers = 13; void init(o2::framework::InitContext&) { triggerMaskBits = jetderiveddatautilities::initialiseTriggerMaskBits(jetderiveddatautilities::JTriggerMasks); @@ -52,7 +55,30 @@ struct TriggerCorrelationsTask { } template - void fillCorrelationsHistogram(T const& collision, bool fill = false, int iCurrentTrig = -1) + void fillOnlineCorrelationsHistogram(T const& collision, bool fill = false, int iCurrentTrig = -1) + { + for (std::vector::size_type iTrig = 0; iTrig < triggerMaskBits.size(); iTrig++) { + if (fill) { + if (iTrig >= 0 && iTrig < nChargedTriggers && jetderiveddatautilities::selectChargedTrigger(collision, iTrig + 1)) { + registry.fill(HIST("triggerCorrelations"), iCurrentTrig, iTrig); + } + if (iTrig >= nChargedTriggers && iTrig < (nChargedTriggers + nChargedHFTriggers) && jetderiveddatautilities::selectChargedHFTrigger(collision, iTrig - nChargedTriggers + 1)) { + registry.fill(HIST("triggerCorrelations"), iCurrentTrig, iTrig); + } + if (iTrig >= (nChargedTriggers + nChargedHFTriggers) && iTrig < (nChargedTriggers + nChargedHFTriggers + nFullTriggers) && jetderiveddatautilities::selectFullTrigger(collision, iTrig - (nChargedTriggers + nChargedHFTriggers) + 1)) { + registry.fill(HIST("triggerCorrelations"), iCurrentTrig, iTrig); + } + + } else { + if (jetderiveddatautilities::selectTrigger(collision, triggerMaskBits[iTrig])) { + fillOnlineCorrelationsHistogram(collision, true, iTrig); + } + } + } + } + + template + void fillOfflineCorrelationsHistogram(T const& collision, bool fill = false, int iCurrentTrig = -1) { for (std::vector::size_type iTrig = 0; iTrig < triggerMaskBits.size(); iTrig++) { if (fill) { @@ -61,23 +87,23 @@ struct TriggerCorrelationsTask { } } else { if (jetderiveddatautilities::selectTrigger(collision, triggerMaskBits[iTrig])) { - fillCorrelationsHistogram(collision, true, iTrig); + fillOfflineCorrelationsHistogram(collision, true, iTrig); } } } } - void processTriggeredCorrelations(soa::Join::iterator const& collision) + void processTriggeredCorrelationsOnline(soa::Join::iterator const& collision) { - fillCorrelationsHistogram(collision); + fillOnlineCorrelationsHistogram(collision); } - PROCESS_SWITCH(TriggerCorrelationsTask, processTriggeredCorrelations, "QA for trigger correlations", true); + PROCESS_SWITCH(TriggerCorrelationsTask, processTriggeredCorrelationsOnline, "QA for online trigger correlations", true); void processTriggeredCorrelationsOffline(aod::JCollision const& collision) { - fillCorrelationsHistogram(collision); + fillOfflineCorrelationsHistogram(collision); } - PROCESS_SWITCH(TriggerCorrelationsTask, processTriggeredCorrelationsOffline, "QA for trigger correlations in offline analysis", false); + PROCESS_SWITCH(TriggerCorrelationsTask, processTriggeredCorrelationsOffline, "QA for offline trigger correlations", false); }; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) From 03c4e3536e2a2f0bf33f012ab671be3dba3bcc99 Mon Sep 17 00:00:00 2001 From: nzardosh Date: Mon, 10 Nov 2025 11:53:57 +0100 Subject: [PATCH 06/13] moving perpendicular cone subtraction --- PWGJE/Core/JetBkgSubUtils.cxx | 58 +------------------ PWGJE/Core/JetBkgSubUtils.h | 10 +--- PWGJE/Core/JetCandidateUtilities.h | 11 ++-- PWGJE/Core/JetDQUtilities.h | 5 +- PWGJE/Core/JetFindingUtilities.h | 43 +++++++++----- PWGJE/Core/JetHFUtilities.h | 5 +- PWGJE/Core/JetUtilities.h | 53 +++++++++++++++++ PWGJE/Core/JetV0Utilities.h | 5 +- PWGJE/Tasks/jetSubstructure.cxx | 29 ++++++++++ PWGJE/Tasks/jetSubstructureHF.cxx | 38 ++++++++++++ .../Nuspex/angularCorrelationsInJets.cxx | 5 +- PWGLF/Tasks/Nuspex/antinucleiInJets.cxx | 13 +++-- .../Strangeness/lambdaJetpolarization.cxx | 3 +- PWGLF/Tasks/Strangeness/strangenessInJets.cxx | 7 ++- 14 files changed, 180 insertions(+), 105 deletions(-) diff --git a/PWGJE/Core/JetBkgSubUtils.cxx b/PWGJE/Core/JetBkgSubUtils.cxx index 4cc6d202130..baa8e65d727 100644 --- a/PWGJE/Core/JetBkgSubUtils.cxx +++ b/PWGJE/Core/JetBkgSubUtils.cxx @@ -48,7 +48,7 @@ JetBkgSubUtils::JetBkgSubUtils(float jetBkgR_out, float bkgEtaMin_out, float bkg void JetBkgSubUtils::initialise() { - // Note: if you are using the PerpCone method you should jetBkgR to be the same as the anit_kt jets R, otherwise use R=0.2 + // Note: recommended to use R=0.2 jetDefBkg = fastjet::JetDefinition(algorithmBkg, jetBkgR, recombSchemeBkg, fastjet::Best); areaDefBkg = fastjet::AreaDefinition(fastjet::active_area_explicit_ghosts, ghostAreaSpec); selRho = fastjet::SelectorRapRange(bkgEtaMin, bkgEtaMax) && fastjet::SelectorPhiRange(bkgPhiMin, bkgPhiMax) && !fastjet::SelectorNHardest(nHardReject); // here we have to put rap range, to be checked! @@ -104,62 +104,6 @@ std::tuple JetBkgSubUtils::estimateRhoAreaMedian(const std::vect return std::make_tuple(rho, rhoM); } -std::tuple JetBkgSubUtils::estimateRhoPerpCone(const std::vector& inputParticles, const std::vector& jets) -{ - - JetBkgSubUtils::initialise(); - if (inputParticles.size() == 0 || jets.size() == 0) { - return std::make_tuple(0.0, 0.0); - } - - double perpPtDensity1 = 0; - double perpPtDensity2 = 0; - double perpMdDensity1 = 0; - double perpMdDensity2 = 0; - - fastjet::Selector selectJet = fastjet::SelectorEtaRange(bkgEtaMin, bkgEtaMax) && fastjet::SelectorPhiRange(bkgPhiMin, bkgPhiMax); - - std::vector selectedJets = fastjet::sorted_by_pt(selectJet(jets)); - - if (selectedJets.size() == 0) { - return std::make_tuple(0.0, 0.0); - } - - fastjet::PseudoJet leadingJet = selectedJets[0]; - - double dPhi1 = 999.; - double dPhi2 = 999.; - double dEta = 999.; - double PerpendicularConeAxisPhi1 = 999., PerpendicularConeAxisPhi2 = 999.; - // build 2 perp cones in phi around the leading jet (right and left of the jet) - PerpendicularConeAxisPhi1 = RecoDecay::constrainAngle(leadingJet.phi() + (M_PI / 2.)); // This will contrain the angel between 0-2Pi - PerpendicularConeAxisPhi2 = RecoDecay::constrainAngle(leadingJet.phi() - (M_PI / 2.)); // This will contrain the angel between 0-2Pi - - for (auto& particle : inputParticles) { - // sum the momentum of all paricles that fill the two cones - dPhi1 = particle.phi() - PerpendicularConeAxisPhi1; - dPhi1 = RecoDecay::constrainAngle(dPhi1, -M_PI); // This will contrain the angel between -pi & Pi - dPhi2 = particle.phi() - PerpendicularConeAxisPhi2; - dPhi2 = RecoDecay::constrainAngle(dPhi2, -M_PI); // This will contrain the angel between -pi & Pi - dEta = leadingJet.eta() - particle.eta(); // The perp cone eta is the same as the leading jet since the cones are perpendicular only in phi - if (TMath::Sqrt(dPhi1 * dPhi1 + dEta * dEta) <= jetBkgR) { - perpPtDensity1 += particle.perp(); - perpMdDensity1 += TMath::Sqrt(particle.m() * particle.m() + particle.pt() * particle.pt()) - particle.pt(); - } - - if (TMath::Sqrt(dPhi2 * dPhi2 + dEta * dEta) <= jetBkgR) { - perpPtDensity2 += particle.perp(); - perpMdDensity2 += TMath::Sqrt(particle.m() * particle.m() + particle.pt() * particle.pt()) - particle.pt(); - } - } - - // Caculate rho as the ratio of average pT of the two cones / the cone area - double perpPtDensity = (perpPtDensity1 + perpPtDensity2) / (2 * M_PI * jetBkgR * jetBkgR); - double perpMdDensity = (perpMdDensity1 + perpMdDensity2) / (2 * M_PI * jetBkgR * jetBkgR); - - return std::make_tuple(perpPtDensity, perpMdDensity); -} - fastjet::PseudoJet JetBkgSubUtils::doRhoAreaSub(const fastjet::PseudoJet& jet, double rhoParam, double rhoMParam) { diff --git a/PWGJE/Core/JetBkgSubUtils.h b/PWGJE/Core/JetBkgSubUtils.h index 81a132c9af9..03b22433460 100644 --- a/PWGJE/Core/JetBkgSubUtils.h +++ b/PWGJE/Core/JetBkgSubUtils.h @@ -30,8 +30,8 @@ enum class BkgSubEstimator { none = 0, medianRho = 1, - medianRhoSparse = 2, - perpCone = 3 + medianRhoSparse = 2 + // perpendicular cone method is in JetUtilities }; enum class BkgSubMode { none = 0, @@ -68,12 +68,6 @@ class JetBkgSubUtils /// @return Rho, RhoM the underlying event density std::tuple estimateRhoAreaMedian(const std::vector& inputParticles, bool doSparseSub); - /// @brief Background estimator using the perpendicular cone method - /// @param inputParticles - /// @param jets (all jets in the event) - /// @return Rho, RhoM the underlying event density - std::tuple estimateRhoPerpCone(const std::vector& inputParticles, const std::vector& jets); - /// @brief method that subtracts the background from jets using the area method /// @param jet input jet to be background subtracted /// @param rhoParam the underlying evvent density vs pT (to be set) diff --git a/PWGJE/Core/JetCandidateUtilities.h b/PWGJE/Core/JetCandidateUtilities.h index fd637ca169a..de6b7faf8f9 100644 --- a/PWGJE/Core/JetCandidateUtilities.h +++ b/PWGJE/Core/JetCandidateUtilities.h @@ -118,17 +118,16 @@ constexpr bool isMatchedCandidate(T const& candidate) * * @param track track that is being checked * @param candidate candidate that is being checked - * @param tracks the track table */ -template -bool isDaughterTrack(T& track, U& candidate, V const& tracks) +template +bool isDaughterTrack(T& track, U& candidate) { if constexpr (jethfutilities::isHFCandidate()) { - return jethfutilities::isHFDaughterTrack(track, candidate, tracks); + return jethfutilities::isHFDaughterTrack(track, candidate); } else if constexpr (jetv0utilities::isV0Candidate()) { - return jetv0utilities::isV0DaughterTrack(track, candidate, tracks); + return jetv0utilities::isV0DaughterTrack(track, candidate); } else if constexpr (jetdqutilities::isDielectronCandidate()) { - return jetdqutilities::isDielectronDaughterTrack(track, candidate, tracks); + return jetdqutilities::isDielectronDaughterTrack(track, candidate); } else { return false; } diff --git a/PWGJE/Core/JetDQUtilities.h b/PWGJE/Core/JetDQUtilities.h index fbca70b178d..97749dae152 100644 --- a/PWGJE/Core/JetDQUtilities.h +++ b/PWGJE/Core/JetDQUtilities.h @@ -94,10 +94,9 @@ constexpr bool isMatchedDielectronCandidate(T const& /*candidate*/) * * @param track track that is being checked * @param candidate Dielectron candidate that is being checked - * @param tracks the track table */ -template -bool isDielectronDaughterTrack(T& track, U& candidate, V const& /*tracks*/) +template +bool isDielectronDaughterTrack(T& track, U& candidate) { if constexpr (isDielectronCandidate()) { if (candidate.prong0Id() == track.globalIndex() || candidate.prong1Id() == track.globalIndex()) { diff --git a/PWGJE/Core/JetFindingUtilities.h b/PWGJE/Core/JetFindingUtilities.h index 8c11ae2bbbb..cf435ca86e0 100644 --- a/PWGJE/Core/JetFindingUtilities.h +++ b/PWGJE/Core/JetFindingUtilities.h @@ -83,24 +83,23 @@ constexpr bool isEMCALClusterTable() } /** - * Adds tracks to a fastjet inputParticles list + * performs all track selections * - * @param inputParticles fastjet container - * @param tracks track table to be added + * @param track track to be checked * @param trackSelection track selection to be applied to tracks * @param candidate optional HF candidiate */ template -void analyseTracks(std::vector& inputParticles, T const& tracks, int trackSelection, bool applyTrackingEfficiency, std::vector trackingEfficiency, std::vector trackingEfficiencyPtBinning, const U* candidate = nullptr) +bool isTrackSelected(T const& track, int trackSelection, bool applyTrackingEfficiency, const std::vector& trackingEfficiency, const std::vector& trackingEfficiencyPtBinning, const U* candidate = nullptr) { - for (auto& track : tracks) { - if (!jetderiveddatautilities::selectTrack(track, trackSelection)) { - continue; - } + + if (!jetderiveddatautilities::selectTrack(track, trackSelection)) { + return false; + } if (candidate != nullptr) { - if (jetcandidateutilities::isDaughterTrack(track, *candidate, tracks)) { - continue; + if (jetcandidateutilities::isDaughterTrack(track, *candidate)) { + return false; } } if (applyTrackingEfficiency) { @@ -109,11 +108,29 @@ void analyseTracks(std::vector& inputParticles, T const& tra std::size_t index = std::distance(trackingEfficiencyPtBinning.begin(), iter) - 1; TRandom3 randomNumber(0); if (randomNumber.Rndm() > trackingEfficiency[index]) { - continue; + return false; } } } - fastjetutilities::fillTracks(track, inputParticles, track.globalIndex()); + return true; +} + +/** + * Adds tracks to a fastjet inputParticles list + * + * @param inputParticles fastjet container + * @param tracks track table to be added + * @param trackSelection track selection to be applied to tracks + * @param candidate optional HF candidiate + */ + +template +void analyseTracks(std::vector& inputParticles, T const& tracks, int trackSelection, bool applyTrackingEfficiency, const std::vector& trackingEfficiency, const std::vector& trackingEfficiencyPtBinning, const U* candidate = nullptr) +{ + for (auto& track : tracks) { + if (isTrackSelected(track, trackSelection, applyTrackingEfficiency, trackingEfficiency, trackingEfficiencyPtBinning, candidate)) { + fastjetutilities::fillTracks(track, inputParticles, track.globalIndex()); + } } } @@ -134,7 +151,7 @@ void analyseTracksMultipleCandidates(std::vector& inputParti continue; } for (auto& candidate : candidates) { - if (jetcandidateutilities::isDaughterTrack(track, candidate, tracks)) { + if (jetcandidateutilities::isDaughterTrack(track, candidate)) { continue; } } diff --git a/PWGJE/Core/JetHFUtilities.h b/PWGJE/Core/JetHFUtilities.h index 796285757bf..dcdf022a849 100644 --- a/PWGJE/Core/JetHFUtilities.h +++ b/PWGJE/Core/JetHFUtilities.h @@ -544,10 +544,9 @@ constexpr bool isMatchedHFCandidate(T const& candidate) * * @param track track that is being checked * @param candidate HF candidate that is being checked - * @param tracks the track table */ -template -bool isHFDaughterTrack(T& track, U& candidate, V const& /*tracks*/) +template +bool isHFDaughterTrack(T& track, U& candidate) { if constexpr (isD0Candidate()) { if (candidate.prong0Id() == track.globalIndex() || candidate.prong1Id() == track.globalIndex()) { diff --git a/PWGJE/Core/JetUtilities.h b/PWGJE/Core/JetUtilities.h index 0efaf1fd8f9..1f7ee8f9791 100644 --- a/PWGJE/Core/JetUtilities.h +++ b/PWGJE/Core/JetUtilities.h @@ -153,6 +153,59 @@ float deltaR(T const& eta1, U const& phi1, V const& eta2, W const& phi2) return std::sqrt(dEta * dEta + dPhi * dPhi); } + +/// @brief Background estimator using the perpendicular cone method +/// @param inputParticles +/// @param jet +/// @return Rho, RhoM the underlying event density + +template +std::tuple estimateRhoPerpCone(const T& inputParticles, const U& jet, V perpConeR) +{ + + if (inputParticles.size() == 0) { + return std::make_tuple(0.0, 0.0); + } + + double perpPtDensity1 = 0; + double perpPtDensity2 = 0; + double perpMdDensity1 = 0; + double perpMdDensity2 = 0; + + double dPhi1 = 999.; + double dPhi2 = 999.; + double dEta = 999.; + double PerpendicularConeAxisPhi1 = 999., PerpendicularConeAxisPhi2 = 999.; + + // build 2 perp cones in phi around the leading jet (right and left of the jet) + PerpendicularConeAxisPhi1 = RecoDecay::constrainAngle(jet.phi() + (M_PI / 2.)); // This will contrain the angel between 0-2Pi + PerpendicularConeAxisPhi2 = RecoDecay::constrainAngle(jet.phi() - (M_PI / 2.)); // This will contrain the angel between 0-2Pi + + for (auto& particle : inputParticles) { + // sum the momentum of all paricles that fill the two cones + dPhi1 = particle.phi() - PerpendicularConeAxisPhi1; + dPhi1 = RecoDecay::constrainAngle(dPhi1, -M_PI); // This will contrain the angel between -pi & Pi + dPhi2 = particle.phi() - PerpendicularConeAxisPhi2; + dPhi2 = RecoDecay::constrainAngle(dPhi2, -M_PI); // This will contrain the angel between -pi & Pi + dEta = jet.eta() - particle.eta(); // The perp cone eta is the same as the leading jet since the cones are perpendicular only in phi + if (TMath::Sqrt(dPhi1 * dPhi1 + dEta * dEta) <= static_cast(perpConeR)) { + perpPtDensity1 += particle.perp(); + perpMdDensity1 += TMath::Sqrt(particle.m() * particle.m() + particle.pt() * particle.pt()) - particle.pt(); + } + + if (TMath::Sqrt(dPhi2 * dPhi2 + dEta * dEta) <= perpConeR) { + perpPtDensity2 += particle.perp(); + perpMdDensity2 += TMath::Sqrt(particle.m() * particle.m() + particle.pt() * particle.pt()) - particle.pt(); + } + } + + // Caculate rho as the ratio of average pT of the two cones / the cone area + double perpPtDensity = (perpPtDensity1 + perpPtDensity2) / (2 * M_PI * static_cast(perpConeR) * static_cast(perpConeR)); + double perpMdDensity = (perpMdDensity1 + perpMdDensity2) / (2 * M_PI * static_cast(perpConeR) * static_cast(perpConeR)); + + return std::make_tuple(perpPtDensity, perpMdDensity); +} + }; // namespace jetutilities #endif // PWGJE_CORE_JETUTILITIES_H_ diff --git a/PWGJE/Core/JetV0Utilities.h b/PWGJE/Core/JetV0Utilities.h index 107b67c0ee8..0a8e6454294 100644 --- a/PWGJE/Core/JetV0Utilities.h +++ b/PWGJE/Core/JetV0Utilities.h @@ -77,10 +77,9 @@ constexpr bool isV0McTable() * * @param track track that is being checked * @param candidate V0 candidate that is being checked - * @param tracks the track table */ -template -bool isV0DaughterTrack(T& track, U& candidate, V const& /*tracks*/) +template +bool isV0DaughterTrack(T& track, U& candidate) { if constexpr (isV0Candidate()) { if (candidate.posTrackId() == track.globalIndex() || candidate.negTrackId() == track.globalIndex()) { diff --git a/PWGJE/Tasks/jetSubstructure.cxx b/PWGJE/Tasks/jetSubstructure.cxx index 92ce238721b..e82a923f38f 100644 --- a/PWGJE/Tasks/jetSubstructure.cxx +++ b/PWGJE/Tasks/jetSubstructure.cxx @@ -19,7 +19,9 @@ #include "RecoDecay.h" #include "PWGJE/Core/FastJetUtilities.h" +#include "PWGJE/Core/JetDerivedDataUtilities.h" #include "PWGJE/Core/JetFinder.h" +#include "PWGJE/Core/JetFindingUtilities.h" #include "PWGJE/Core/JetSubstructureUtilities.h" #include "PWGJE/Core/JetUtilities.h" #include "PWGJE/DataModel/Jet.h" @@ -74,6 +76,10 @@ struct JetSubstructureTask { Configurable alpha{"alpha", 1.0, "angularity alpha"}; Configurable doPairBkg{"doPairBkg", true, "save bkg pairs"}; Configurable pairConstituentPtMin{"pairConstituentPtMin", 1.0, "pt cut off for constituents going into pairs"}; + Configurable trackSelections{"trackSelections", "globalTracks", "set track selections"}; + Configurable applyTrackingEfficiency{"applyTrackingEfficiency", {false}, "configurable to decide whether to apply artificial tracking efficiency (discarding tracks) in jet finding"}; + Configurable> trackingEfficiencyPtBinning{"trackingEfficiencyPtBinning", {0., 10, 999.}, "pt binning of tracking efficiency array if applyTrackingEfficiency is true"}; + Configurable> trackingEfficiency{"trackingEfficiency", {1.0, 1.0}, "tracking efficiency array applied to jet finding if applyTrackingEfficiency is true"}; Service pdg; std::vector jetConstituents; @@ -103,6 +109,8 @@ struct JetSubstructureTask { HistogramRegistry registry; + int trackSelection = -1; + void init(InitContext const&) { registry.add("h2_jet_pt_jet_zg", ";#it{p}_{T,jet} (GeV/#it{c});#it{z}_{g}", {HistType::kTH2F, {{200, 0., 200.}, {22, 0.0, 1.1}}}); @@ -120,6 +128,17 @@ struct JetSubstructureTask { jetReclusterer.isReclustering = true; jetReclusterer.algorithm = fastjet::JetAlgorithm::cambridge_algorithm; jetReclusterer.ghostRepeatN = 0; + + trackSelection = jetderiveddatautilities::initialiseTrackSelection(static_cast(trackSelections)); + + if (applyTrackingEfficiency) { + if (trackingEfficiencyPtBinning->size() < 2) { + LOGP(fatal, "jetFinder workflow: trackingEfficiencyPtBinning configurable should have at least two bin edges"); + } + if (trackingEfficiency->size() + 1 != trackingEfficiencyPtBinning->size()) { + LOGP(fatal, "jetFinder workflow: trackingEfficiency configurable should have exactly one less entry than the number of bin edges set in trackingEfficiencyPtBinning configurable"); + } + } } Preslice TracksPerCollision = aod::jtrack::collisionId; @@ -247,6 +266,16 @@ struct JetSubstructureTask { std::vector tracksPerpCone1Vec; std::vector tracksPerpCone2Vec; for (auto const& track : tracksPerCollision) { + if (!jetderiveddatautilities::applyTrackKinematics(track)) { + continue; + } + + if constexpr (!std::is_same_v, aod::JetParticles>) { + if (!jetfindingutilities::isTrackSelected(track, trackSelection, applyTrackingEfficiency, trackingEfficiency, trackingEfficiencyPtBinning)) { + continue; + } + } + float deltaPhi1 = track.phi() - perpCone1Phi; deltaPhi1 = RecoDecay::constrainAngle(deltaPhi1, -M_PI); float deltaPhi2 = track.phi() - perpCone2Phi; diff --git a/PWGJE/Tasks/jetSubstructureHF.cxx b/PWGJE/Tasks/jetSubstructureHF.cxx index 3693d0908ff..7d0917e7b51 100644 --- a/PWGJE/Tasks/jetSubstructureHF.cxx +++ b/PWGJE/Tasks/jetSubstructureHF.cxx @@ -19,6 +19,7 @@ #include "PWGJE/Core/FastJetUtilities.h" #include "PWGJE/Core/JetDQUtilities.h" #include "PWGJE/Core/JetFinder.h" +#include "PWGJE/Core/JetFindingUtilities.h" #include "PWGJE/Core/JetHFUtilities.h" #include "PWGJE/Core/JetSubstructureUtilities.h" #include "PWGJE/Core/JetUtilities.h" @@ -77,6 +78,10 @@ struct JetSubstructureHFTask { Configurable alpha{"alpha", 1.0, "angularity alpha"}; Configurable doPairBkg{"doPairBkg", true, "save bkg pairs"}; Configurable pairConstituentPtMin{"pairConstituentPtMin", 1.0, "pt cut off for constituents going into pairs"}; + Configurable trackSelections{"trackSelections", "globalTracks", "set track selections"}; + Configurable applyTrackingEfficiency{"applyTrackingEfficiency", {false}, "configurable to decide whether to apply artificial tracking efficiency (discarding tracks) in jet finding"}; + Configurable> trackingEfficiencyPtBinning{"trackingEfficiencyPtBinning", {0., 10, 999.}, "pt binning of tracking efficiency array if applyTrackingEfficiency is true"}; + Configurable> trackingEfficiency{"trackingEfficiency", {1.0, 1.0}, "tracking efficiency array applied to jet finding if applyTrackingEfficiency is true"}; Service pdg; float candMass; @@ -107,6 +112,9 @@ struct JetSubstructureHFTask { float perpConeRho; HistogramRegistry registry; + + int trackSelection = -1; + void init(InitContext const&) { registry.add("h2_jet_pt_jet_zg", ";#it{p}_{T,jet} (GeV/#it{c});#it{z}_{g}", {HistType::kTH2F, {{200, 0., 200.}, {22, 0.0, 1.1}}}); @@ -126,6 +134,17 @@ struct JetSubstructureHFTask { jetReclusterer.ghostRepeatN = 0; candMass = jetcandidateutilities::getTablePDGMass(); + + trackSelection = jetderiveddatautilities::initialiseTrackSelection(static_cast(trackSelections)); + + if (applyTrackingEfficiency) { + if (trackingEfficiencyPtBinning->size() < 2) { + LOGP(fatal, "jetFinder workflow: trackingEfficiencyPtBinning configurable should have at least two bin edges"); + } + if (trackingEfficiency->size() + 1 != trackingEfficiencyPtBinning->size()) { + LOGP(fatal, "jetFinder workflow: trackingEfficiency configurable should have exactly one less entry than the number of bin edges set in trackingEfficiencyPtBinning configurable"); + } + } } Preslice TracksPerCollision = aod::jtrack::collisionId; @@ -329,6 +348,25 @@ struct JetSubstructureHFTask { std::vector tracksPerpCone1Vec; std::vector tracksPerpCone2Vec; for (auto const& track : tracksPerCollision) { + if (!jetderiveddatautilities::applyTrackKinematics(track)) { + continue; + } + + if constexpr (!std::is_same_v, aod::JetParticles>) { + if (!jetfindingutilities::isTrackSelected(track, trackSelection, applyTrackingEfficiency, trackingEfficiency, trackingEfficiencyPtBinning)) { + continue; + } + } + bool isdaughterTrack = false; + for (auto& candidate : jet.template candidates_as()) { + if (jetcandidateutilities::isDaughterTrack(track, candidate)) { + isdaughterTrack = true; + break; + } + } + if (isdaughterTrack) { + continue; + } float deltaPhi1 = track.phi() - perpCone1Phi; deltaPhi1 = RecoDecay::constrainAngle(deltaPhi1, -M_PI); float deltaPhi2 = track.phi() - perpCone2Phi; diff --git a/PWGLF/Tasks/Nuspex/angularCorrelationsInJets.cxx b/PWGLF/Tasks/Nuspex/angularCorrelationsInJets.cxx index 8119a38dd94..30d9b2125d2 100644 --- a/PWGLF/Tasks/Nuspex/angularCorrelationsInJets.cxx +++ b/PWGLF/Tasks/Nuspex/angularCorrelationsInJets.cxx @@ -15,6 +15,7 @@ /// \brief task for analysis of angular correlations in jets using Fastjet #include "PWGJE/Core/JetBkgSubUtils.h" +#include "PWGJE/Core/JetUtilities.h" #include "Common/Core/PID/PIDTOF.h" #include "Common/Core/RecoDecay.h" @@ -1110,7 +1111,7 @@ struct AngularCorrelationsInJets { registryData.fill(HIST("eventProtocol"), 3); - auto [rhoPerp, rhoMPerp] = bkgSub.estimateRhoPerpCone(jetInput, jets); + auto [rhoPerp, rhoMPerp] = jetutilities::estimateRhoPerpCone(jetInput, jets[0], jetR); for (const auto& jet : jets) { // this is where the magic happens @@ -1291,7 +1292,7 @@ struct AngularCorrelationsInJets { registryData.fill(HIST("eventProtocol"), 3); - auto [rhoPerp, rhoMPerp] = bkgSub.estimateRhoPerpCone(jetInput, jets); + auto [rhoPerp, rhoMPerp] = jetutilities::estimateRhoPerpCone(jetInput, jets[0], jetR); for (auto& jet : jets) { // o2-linter: disable=const-ref-in-for-loop (jets are modified) if (!jet.has_constituents()) diff --git a/PWGLF/Tasks/Nuspex/antinucleiInJets.cxx b/PWGLF/Tasks/Nuspex/antinucleiInJets.cxx index ab183909dfa..0a2d1bf51c1 100644 --- a/PWGLF/Tasks/Nuspex/antinucleiInJets.cxx +++ b/PWGLF/Tasks/Nuspex/antinucleiInJets.cxx @@ -17,6 +17,7 @@ #include "PWGJE/Core/JetBkgSubUtils.h" #include "PWGJE/Core/JetDerivedDataUtilities.h" +#include "PWGJE/Core/JetUtilities.h" #include "PWGJE/DataModel/Jet.h" #include "PWGJE/DataModel/JetReducedData.h" @@ -1000,7 +1001,7 @@ struct AntinucleiInJets { fastjet::AreaDefinition areaDef(fastjet::active_area, fastjet::GhostedAreaSpec(1.0)); fastjet::ClusterSequenceArea cs(fjParticles, jetDef, areaDef); std::vector jets = fastjet::sorted_by_pt(cs.inclusive_jets()); - auto [rhoPerp, rhoMPerp] = backgroundSub.estimateRhoPerpCone(fjParticles, jets); + auto [rhoPerp, rhoMPerp] = jetutilities::estimateRhoPerpCone(fjParticles, jets[0], rJet); // Loop over reconstructed jets bool isAtLeastOneJetSelected = false; @@ -1283,7 +1284,7 @@ struct AntinucleiInJets { fastjet::AreaDefinition areaDef(fastjet::active_area, fastjet::GhostedAreaSpec(1.0)); fastjet::ClusterSequenceArea cs(fjParticles, jetDef, areaDef); std::vector jets = fastjet::sorted_by_pt(cs.inclusive_jets()); - auto [rhoPerp, rhoMPerp] = backgroundSub.estimateRhoPerpCone(fjParticles, jets); + auto [rhoPerp, rhoMPerp] = jetutilities::estimateRhoPerpCone(fjParticles, jets[0], rJet); // Loop over reconstructed jets bool isAtLeastOneJetSelected = false; @@ -1352,7 +1353,7 @@ struct AntinucleiInJets { fastjet::AreaDefinition areaDef(fastjet::active_area, fastjet::GhostedAreaSpec(1.0)); fastjet::ClusterSequenceArea cs(fjParticles, jetDef, areaDef); std::vector jets = fastjet::sorted_by_pt(cs.inclusive_jets()); - auto [rhoPerp, rhoMPerp] = backgroundSub.estimateRhoPerpCone(fjParticles, jets); + auto [rhoPerp, rhoMPerp] = jetutilities::estimateRhoPerpCone(fjParticles, jets[0], rJet); // Loop over reconstructed jets int njetsInAcc(0); @@ -1820,7 +1821,7 @@ struct AntinucleiInJets { // Cluster MC particles into jets using anti-kt algorithm fastjet::ClusterSequenceArea cs(fjParticles, jetDef, areaDef); std::vector jets = fastjet::sorted_by_pt(cs.inclusive_jets()); - auto [rhoPerp, rhoMPerp] = backgroundSub.estimateRhoPerpCone(fjParticles, jets); + auto [rhoPerp, rhoMPerp] = jetutilities::estimateRhoPerpCone(fjParticles, jets[0], rJet); // Loop over clustered jets bool isAtLeastOneJetSelected = false; @@ -2032,7 +2033,7 @@ struct AntinucleiInJets { // Cluster particles using the anti-kt algorithm fastjet::ClusterSequenceArea cs(fjParticles, jetDef, areaDef); std::vector jets = fastjet::sorted_by_pt(cs.inclusive_jets()); - auto [rhoPerp, rhoMPerp] = backgroundSub.estimateRhoPerpCone(fjParticles, jets); + auto [rhoPerp, rhoMPerp] = jetutilities::estimateRhoPerpCone(fjParticles, jets[0], rJet); // Loop over reconstructed jets bool isAtLeastOneJetSelected = false; @@ -2772,7 +2773,7 @@ struct AntinucleiInJets { fastjet::AreaDefinition areaDef(fastjet::active_area, fastjet::GhostedAreaSpec(1.0)); fastjet::ClusterSequenceArea cs(fjParticles, jetDef, areaDef); std::vector jets = fastjet::sorted_by_pt(cs.inclusive_jets()); - auto [rhoPerp, rhoMPerp] = backgroundSub.estimateRhoPerpCone(fjParticles, jets); + auto [rhoPerp, rhoMPerp] = jetutilities::estimateRhoPerpCone(fjParticles, jets[0], rJet); // Loop over reconstructed jets bool isAtLeastOneJetSelected = false; diff --git a/PWGLF/Tasks/Strangeness/lambdaJetpolarization.cxx b/PWGLF/Tasks/Strangeness/lambdaJetpolarization.cxx index 710ba9361a9..7c8ca28fb05 100644 --- a/PWGLF/Tasks/Strangeness/lambdaJetpolarization.cxx +++ b/PWGLF/Tasks/Strangeness/lambdaJetpolarization.cxx @@ -15,6 +15,7 @@ #include "PWGJE/Core/JetBkgSubUtils.h" #include "PWGJE/Core/JetDerivedDataUtilities.h" +#include "PWGJE/Core/JetUtilities.h" #include "PWGJE/DataModel/Jet.h" #include "PWGJE/DataModel/JetReducedData.h" #include "PWGLF/DataModel/LFStrangenessTables.h" @@ -1158,7 +1159,7 @@ struct LfMyV0s { fastjet::AreaDefinition areaDef(fastjet::active_area, fastjet::GhostedAreaSpec(1.0)); fastjet::ClusterSequenceArea cs(fjParticles, jetDef, areaDef); std::vector jets = fastjet::sorted_by_pt(cs.inclusive_jets()); - auto [rhoPerp, rhoMPerp] = backgroundSub.estimateRhoPerpCone(fjParticles, jets); + auto [rhoPerp, rhoMPerp] = jetutilities::estimateRhoPerpCone(fjParticles, jets[0], rJet); // jet selection bool isAtLeastOneJetSelected = false; diff --git a/PWGLF/Tasks/Strangeness/strangenessInJets.cxx b/PWGLF/Tasks/Strangeness/strangenessInJets.cxx index e9b0a4a55df..404036430b5 100644 --- a/PWGLF/Tasks/Strangeness/strangenessInJets.cxx +++ b/PWGLF/Tasks/Strangeness/strangenessInJets.cxx @@ -21,6 +21,7 @@ #include "PWGJE/Core/JetBkgSubUtils.h" #include "PWGJE/Core/JetDerivedDataUtilities.h" +#include "PWGJE/Core/JetUtilities.h" #include "PWGJE/DataModel/Jet.h" #include "PWGJE/DataModel/JetReducedData.h" #include "PWGLF/DataModel/LFStrangenessTables.h" @@ -973,7 +974,7 @@ struct StrangenessInJets { fastjet::AreaDefinition areaDef(fastjet::active_area, fastjet::GhostedAreaSpec(1.0)); fastjet::ClusterSequenceArea cs(fjParticles, jetDef, areaDef); std::vector jets = fastjet::sorted_by_pt(cs.inclusive_jets()); - auto [rhoPerp, rhoMPerp] = backgroundSub.estimateRhoPerpCone(fjParticles, jets); + auto [rhoPerp, rhoMPerp] = jetutilities::estimateRhoPerpCone(fjParticles, jets[0], rJet); // Jet selection bool isAtLeastOneJetSelected = false; @@ -1253,7 +1254,7 @@ struct StrangenessInJets { std::vector jets = fastjet::sorted_by_pt(cs.inclusive_jets()); // Estimate background energy density (rho) in perpendicular cone - auto [rhoPerp, rhoMPerp] = backgroundSub.estimateRhoPerpCone(fjParticles, jets); + auto [rhoPerp, rhoMPerp] = jetutilities::estimateRhoPerpCone(fjParticles, jets[0], rJet); // Loop over clustered jets for (const auto& jet : jets) { @@ -1519,7 +1520,7 @@ struct StrangenessInJets { // Cluster particles using the anti-kt algorithm fastjet::ClusterSequenceArea cs(fjParticles, jetDef, areaDef); std::vector jets = fastjet::sorted_by_pt(cs.inclusive_jets()); - auto [rhoPerp, rhoMPerp] = backgroundSub.estimateRhoPerpCone(fjParticles, jets); + auto [rhoPerp, rhoMPerp] = jetutilities::estimateRhoPerpCone(fjParticles, jets[0], rJet); // Jet selection bool isAtLeastOneJetSelected = false; From cb08b0031222ebee8610be08773767e907be6433 Mon Sep 17 00:00:00 2001 From: nzardosh Date: Mon, 10 Nov 2025 15:49:22 +0100 Subject: [PATCH 07/13] adding the ability to select on particle level trigger tracks for derived data --- PWGJE/TableProducer/derivedDataSelector.cxx | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/PWGJE/TableProducer/derivedDataSelector.cxx b/PWGJE/TableProducer/derivedDataSelector.cxx index 215f56b00eb..64f185d6c08 100644 --- a/PWGJE/TableProducer/derivedDataSelector.cxx +++ b/PWGJE/TableProducer/derivedDataSelector.cxx @@ -82,6 +82,7 @@ struct JetDerivedDataSelector { Configurable thresholdChargedEventWiseSubtractedDielectronJetPtMin{"thresholdChargedEventWiseSubtractedDielectronJetPtMin", 0.0, "Minimum charged event-wise subtracted Dielectron jet pt to accept event"}; Configurable thresholdChargedDielectronMCPJetPtMin{"thresholdChargedDielectronMCPJetPtMin", 0.0, "Minimum charged Dielectron mcp jet pt to accept event"}; Configurable thresholdTriggerTrackPtMin{"thresholdTriggerTrackPtMin", 0.0, "Minimum trigger track pt to accept event"}; + Configurable thresholdTriggerParticlePtMin{"thresholdTriggerParticlePtMin", 0.0, "Minimum trigger particle pt to accept event"}; Configurable thresholdClusterEnergyMin{"thresholdClusterEnergyMin", 0.0, "Minimum cluster energy to accept event"}; Configurable downscaleFactor{"downscaleFactor", 1, "random downscale of selected events"}; @@ -268,6 +269,8 @@ struct JetDerivedDataSelector { selectionObjectPtMin = config.thresholdChargedDielectronMCPJetPtMin; } else if constexpr (std::is_same_v, aod::JTracks>) { selectionObjectPtMin = config.thresholdTriggerTrackPtMin; + } else if constexpr (std::is_same_v, aod::JMcParticles>) { + selectionObjectPtMin = config.thresholdTriggerParticlePtMin; } else if constexpr (std::is_same_v, aod::JClusters>) { selectionObjectPtMin = config.thresholdClusterEnergyMin; } else { @@ -293,7 +296,7 @@ struct JetDerivedDataSelector { } } if (isTriggerObject) { - if constexpr (std::is_same_v, aod::ChargedMCParticleLevelJets> || std::is_same_v, aod::NeutralMCParticleLevelJets> || std::is_same_v, aod::FullMCParticleLevelJets> || std::is_same_v, aod::D0ChargedMCParticleLevelJets> || std::is_same_v, aod::DplusChargedMCParticleLevelJets> || std::is_same_v, aod::DsChargedMCParticleLevelJets> || std::is_same_v, aod::DstarChargedMCParticleLevelJets> || std::is_same_v, aod::LcChargedMCParticleLevelJets> || std::is_same_v, aod::B0ChargedMCParticleLevelJets> || std::is_same_v, aod::BplusChargedMCParticleLevelJets> || std::is_same_v, aod::XicToXiPiPiChargedMCParticleLevelJets> || std::is_same_v, aod::DielectronChargedMCParticleLevelJets>) { + if constexpr (std::is_same_v, aod::ChargedMCParticleLevelJets> || std::is_same_v, aod::NeutralMCParticleLevelJets> || std::is_same_v, aod::FullMCParticleLevelJets> || std::is_same_v, aod::D0ChargedMCParticleLevelJets> || std::is_same_v, aod::DplusChargedMCParticleLevelJets> || std::is_same_v, aod::DsChargedMCParticleLevelJets> || std::is_same_v, aod::DstarChargedMCParticleLevelJets> || std::is_same_v, aod::LcChargedMCParticleLevelJets> || std::is_same_v, aod::B0ChargedMCParticleLevelJets> || std::is_same_v, aod::BplusChargedMCParticleLevelJets> || std::is_same_v, aod::XicToXiPiPiChargedMCParticleLevelJets> || std::is_same_v, aod::DielectronChargedMCParticleLevelJets> || std::is_same_v, aod::JMcParticles>) { if (selectionObject.mcCollisionId() >= 0) { McCollisionFlag[selectionObject.mcCollisionId()] = true; } @@ -369,6 +372,7 @@ struct JetDerivedDataSelector { PROCESS_SWITCH_FULL(JetDerivedDataSelector, processSelectionObjects, processSelectingDielectronChargedMCPJets, "process Dielectron charged mcp jets", false); PROCESS_SWITCH_FULL(JetDerivedDataSelector, processSelectionObjects, processSelectingClusters, "process EMCal clusters", false); PROCESS_SWITCH_FULL(JetDerivedDataSelector, processSelectionObjects, processSelectingTracks, "process high pt tracks", false); + PROCESS_SWITCH_FULL(JetDerivedDataSelector, processSelectionObjects, processSelectingParticles, "process high pt particles", false); PROCESS_SWITCH_FULL(JetDerivedDataSelector, processDoDownscaling, processCollisionDownscaling, "process downsaling of triggered collisions", false); PROCESS_SWITCH_FULL(JetDerivedDataSelector, processDoDownscaling, processMcCollisionDownscaling, "process downsaling of triggered mccollisions", false); PROCESS_SWITCH(JetDerivedDataSelector, processDoCollisionSelections, "process event selections for saved events", false); From 5f88f5c5bde852d7a488ec45fa59ac2b0f99de8c Mon Sep 17 00:00:00 2001 From: nzardosh Date: Mon, 10 Nov 2025 20:54:46 +0100 Subject: [PATCH 08/13] fixing cpp checks --- PWGJE/Core/JetDQUtilities.h | 2 +- PWGJE/Core/JetDerivedDataUtilities.h | 22 ++++----- PWGJE/Core/JetFindingUtilities.h | 6 +-- PWGJE/Core/JetMatchingUtilities.h | 12 ++--- PWGJE/Core/JetTaggingUtilities.h | 7 ++- PWGJE/Core/JetUtilities.h | 9 ++-- PWGJE/Core/JetV0Utilities.h | 2 +- PWGJE/JetFinders/jetFinder.cxx | 1 - PWGJE/JetFinders/jetFinderHF.cxx | 1 - PWGJE/JetFinders/jetFinderV0.cxx | 1 - PWGJE/TableProducer/derivedDataProducer.cxx | 12 ++--- PWGJE/TableProducer/derivedDataWriter.cxx | 6 +-- .../eventwiseConstituentSubtractor.cxx | 1 - PWGJE/TableProducer/jetTaggerHF.cxx | 4 +- PWGJE/Tasks/bjetTaggingML.cxx | 2 +- PWGJE/Tasks/bjetTreeCreator.cxx | 2 +- PWGJE/Tasks/dijetFinderQA.cxx | 16 +++---- PWGJE/Tasks/fullJetSpectra.cxx | 41 +++++----------- PWGJE/Tasks/fullJetTriggerQATask.cxx | 12 ++--- PWGJE/Tasks/gammaJetTreeProducer.cxx | 2 - PWGJE/Tasks/jetChCorr.cxx | 8 ++-- PWGJE/Tasks/jetChargedV2.cxx | 48 +++++++------------ PWGJE/Tasks/jetFormationTimeReclustering.cxx | 9 ++-- PWGJE/Tasks/jetFragmentation.cxx | 10 ++-- PWGJE/Tasks/jetHadronRecoil.cxx | 24 ++++------ PWGJE/Tasks/jetPlanarFlow.cxx | 2 +- PWGJE/Tasks/jetSpectraEseTask.cxx | 2 +- PWGJE/Tasks/jetSubstructure.cxx | 6 +-- PWGJE/Tasks/jetSubstructureHF.cxx | 6 +-- PWGJE/Tasks/jetValidationQA.cxx | 8 ++-- PWGJE/Tasks/mcGeneratorStudies.cxx | 2 +- PWGJE/Tasks/nucleiInJets.cxx | 17 +++---- PWGJE/Tasks/phiInJets.cxx | 10 ++-- 33 files changed, 122 insertions(+), 191 deletions(-) diff --git a/PWGJE/Core/JetDQUtilities.h b/PWGJE/Core/JetDQUtilities.h index 97749dae152..75ad0a399b9 100644 --- a/PWGJE/Core/JetDQUtilities.h +++ b/PWGJE/Core/JetDQUtilities.h @@ -282,7 +282,7 @@ bool selectDielectronParticleDecay(T const& dielectronParticle, int dielectronPa return (dielectronParticle.decayFlag() & (1 << dielectronParticleDecaySelection)); } -int initialiseDielectronParticleDecaySelection(std::string dielectronParticleDecaySelection) +int initialiseDielectronParticleDecaySelection(const std::string& dielectronParticleDecaySelection) { if (dielectronParticleDecaySelection == "JPsiToEE") { return JDielectronParticleDecays::JPsiToEE; diff --git a/PWGJE/Core/JetDerivedDataUtilities.h b/PWGJE/Core/JetDerivedDataUtilities.h index 921959742c4..1b02a380312 100644 --- a/PWGJE/Core/JetDerivedDataUtilities.h +++ b/PWGJE/Core/JetDerivedDataUtilities.h @@ -54,7 +54,7 @@ enum JCollisionSubGeneratorId { }; template -bool selectCollision(T const& collision, std::vector eventSelectionMaskBits, bool skipMBGapEvents = true, bool rctSelection = true, std::string rctLabel = "CBT_hadronPID", bool rejectLimitedAcceptanceRct = false, bool requireZDCRct = false) +bool selectCollision(T const& collision, const std::vector& eventSelectionMaskBits, bool skipMBGapEvents = true, bool rctSelection = true, std::string rctLabel = "CBT_hadronPID", bool rejectLimitedAcceptanceRct = false, bool requireZDCRct = false) { if (skipMBGapEvents && collision.subGeneratorId() == JCollisionSubGeneratorId::mbGap) { return false; @@ -75,7 +75,7 @@ bool selectCollision(T const& collision, std::vector eventSelectionMaskBits return true; } -bool eventSelectionMasksContainSelection(std::string eventSelectionMasks, std::string selection) +bool eventSelectionMasksContainSelection(const std::string& eventSelectionMasks, std::string selection) { size_t position = 0; while ((position = eventSelectionMasks.find(selection, position)) != std::string::npos) { @@ -89,7 +89,7 @@ bool eventSelectionMasksContainSelection(std::string eventSelectionMasks, std::s return false; } -std::vector initialiseEventSelectionBits(std::string eventSelectionMasks) +std::vector initialiseEventSelectionBits(const std::string& eventSelectionMasks) { std::vector eventSelectionMaskBits; if (eventSelectionMasksContainSelection(eventSelectionMasks, "sel8")) { @@ -236,7 +236,7 @@ enum JTrigSel { }; template -bool selectTrigger(T const& collision, std::vector triggerMaskBits) +bool selectTrigger(T const& collision, const std::vector& triggerMaskBits) { if (triggerMaskBits.size() == 0) { return true; @@ -258,7 +258,7 @@ bool selectTrigger(T const& collision, int triggerMaskBit) return collision.triggerSel() & (1 << triggerMaskBit); } -bool triggerMasksContainTrigger(std::string triggerMasks, std::string trigger) +bool triggerMasksContainTrigger(const std::string& triggerMasks, std::string trigger) { size_t position = 0; while ((position = triggerMasks.find(trigger, position)) != std::string::npos) { @@ -272,7 +272,7 @@ bool triggerMasksContainTrigger(std::string triggerMasks, std::string trigger) return false; } -std::vector initialiseTriggerMaskBits(std::string triggerMasks) +std::vector initialiseTriggerMaskBits(const std::string& triggerMasks) { std::vector triggerMaskBits; if (triggerMasksContainTrigger(triggerMasks, "fJetChLowPt")) { @@ -341,7 +341,7 @@ std::vector initialiseTriggerMaskBits(std::string triggerMasks) return triggerMaskBits; } -uint64_t setTriggerSelectionBit(std::vector triggerDecisions) +uint64_t setTriggerSelectionBit(const std::vector& triggerDecisions) { uint64_t bit = 0; for (std::vector::size_type i = 0; i < triggerDecisions.size(); i++) { @@ -369,7 +369,7 @@ bool selectChargedTrigger(T const& collision, int triggerSelection) return (collision.chargedTriggerSel() & (1 << triggerSelection)); } -int initialiseChargedTriggerSelection(std::string triggerSelection) +int initialiseChargedTriggerSelection(const std::string& triggerSelection) { if (triggerSelection == "jetChLowPt") { return JTrigSelCh::jetChLowPt; @@ -434,7 +434,7 @@ bool selectFullTrigger(T const& collision, int triggerSelection) return (collision.fullTriggerSel() & (1 << triggerSelection)); } -int initialiseFullTriggerSelection(std::string triggerSelection) +int initialiseFullTriggerSelection(const std::string& triggerSelection) { if (triggerSelection == "emcalReadout") { return JTrigSelFull::emcalReadout; @@ -529,7 +529,7 @@ bool selectChargedHFTrigger(T const& collision, int triggerSelection) return (collision.chargedHFTriggerSel() & (1 << triggerSelection)); } -int initialiseChargedHFTriggerSelection(std::string triggerSelection) +int initialiseChargedHFTriggerSelection(const std::string& triggerSelection) { if (triggerSelection == "jetD0ChLowPt") { return JTrigSelChHF::jetD0ChLowPt; @@ -592,7 +592,7 @@ bool selectTrack(T const& track, int trackSelection) return (track.trackSel() & (1 << trackSelection)); } -int initialiseTrackSelection(std::string trackSelection) +int initialiseTrackSelection(const std::string& trackSelection) { if (trackSelection == "globalTracks") { return JTrackSel::globalTrack; diff --git a/PWGJE/Core/JetFindingUtilities.h b/PWGJE/Core/JetFindingUtilities.h index cf435ca86e0..b082b5df4e4 100644 --- a/PWGJE/Core/JetFindingUtilities.h +++ b/PWGJE/Core/JetFindingUtilities.h @@ -144,7 +144,7 @@ void analyseTracks(std::vector& inputParticles, T const& tra */ template -void analyseTracksMultipleCandidates(std::vector& inputParticles, T const& tracks, int trackSelection, bool applyTrackingEfficiency, std::vector trackingEfficiency, std::vector trackingEfficiencyPtBinning, U const& candidates) +void analyseTracksMultipleCandidates(std::vector& inputParticles, T const& tracks, int trackSelection, bool applyTrackingEfficiency, const std::vector& trackingEfficiency, const std::vector& trackingEfficiencyPtBinning, U const& candidates) { for (auto& track : tracks) { if (!jetderiveddatautilities::selectTrack(track, trackSelection)) { @@ -312,8 +312,8 @@ void findJets(JetFinder& jetFinder, std::vector& inputPartic if (fillThnSparse) { thnSparseJet->Fill(R, jet.pt(), jet.eta(), jet.phi()); // important for normalisation in V0Jet analyses to store all jets, including those that aren't V0s } - bool isCandidateJet = false; if (doCandidateJetFinding) { + bool isCandidateJet = false; for (const auto& constituent : jet.constituents()) { auto constituentStatus = constituent.template user_info().getStatus(); if (constituentStatus == static_cast(JetConstituentStatus::candidate)) { // note currently we cannot run V0 and HF in the same jet. If we ever need to we can seperate the loops @@ -357,7 +357,7 @@ void findJets(JetFinder& jetFinder, std::vector& inputPartic * @param candidate optional hf candidiate */ template -void analyseParticles(std::vector& inputParticles, std::string particleSelection, int jetTypeParticleLevel, T const& particles, o2::framework::Service pdgDatabase, const U* candidate = nullptr) +void analyseParticles(std::vector& inputParticles, const std::string& particleSelection, int jetTypeParticleLevel, T const& particles, o2::framework::Service pdgDatabase, const U* candidate = nullptr) { for (auto& particle : particles) { if (particleSelection == "PhysicalPrimary" && !particle.isPhysicalPrimary()) { // CHECK : Does this exclude the HF hadron? diff --git a/PWGJE/Core/JetMatchingUtilities.h b/PWGJE/Core/JetMatchingUtilities.h index 16f27fa4f0f..cc972fc0a9c 100644 --- a/PWGJE/Core/JetMatchingUtilities.h +++ b/PWGJE/Core/JetMatchingUtilities.h @@ -330,27 +330,25 @@ void MatchGeo(T const& jetsBasePerCollision, U const& jetsTagPerCollision, std:: } std::tie(baseToTagMatchingGeoIndex, tagToBaseMatchingGeoIndex) = MatchJetsGeometrically(jetsBasePhi, jetsBaseEta, jetsTagPhi, jetsTagEta, maxMatchingDistance); // change max distnace to a function call int jetBaseIndex = 0; + int jetTagIndex = 0; for (const auto& jetBase : jetsBasePerCollision) { if (std::round(jetBase.r()) != std::round(jetR)) { continue; } - int jetTagIndex = baseToTagMatchingGeoIndex[jetBaseIndex]; - int jetTagGlobalIndex; + jetTagIndex = baseToTagMatchingGeoIndex[jetBaseIndex]; if (jetTagIndex > -1 && jetTagIndex < jetsTagPerCollision.size()) { - jetTagGlobalIndex = jetsTagPerCollision.iteratorAt(jetTagIndex).globalIndex(); + int jetTagGlobalIndex = jetsTagPerCollision.iteratorAt(jetTagIndex).globalIndex(); baseToTagMatchingGeo[jetBase.globalIndex()].push_back(jetTagGlobalIndex); } jetBaseIndex++; } - int jetTagIndex = 0; for (const auto& jetTag : jetsTagPerCollision) { if (std::round(jetTag.r()) != std::round(jetR)) { continue; } - int jetBaseIndex = tagToBaseMatchingGeoIndex[jetTagIndex]; - int jetBaseGlobalIndex; + jetBaseIndex = tagToBaseMatchingGeoIndex[jetTagIndex]; if (jetBaseIndex > -1 && jetBaseIndex < jetsBasePerCollision.size()) { - jetBaseGlobalIndex = jetsBasePerCollision.iteratorAt(jetBaseIndex).globalIndex(); + int jetBaseGlobalIndex = jetsBasePerCollision.iteratorAt(jetBaseIndex).globalIndex(); tagToBaseMatchingGeo[jetTag.globalIndex()].push_back(jetBaseGlobalIndex); } jetTagIndex++; diff --git a/PWGJE/Core/JetTaggingUtilities.h b/PWGJE/Core/JetTaggingUtilities.h index 81141b036f4..56c4d26b989 100644 --- a/PWGJE/Core/JetTaggingUtilities.h +++ b/PWGJE/Core/JetTaggingUtilities.h @@ -248,10 +248,9 @@ template int jetParticleFromHFShower(T const& jet, U const& particles, typename U::iterator& hfparticle, bool searchUpToQuark) { - int origin = -1; for (const auto& particle : jet.template tracks_as()) { hfparticle = particle; // for init if origin is 1 or 2, the particle is not hfparticle - origin = RecoDecay::getParticleOrigin(particles, particle, searchUpToQuark); + int origin = RecoDecay::getParticleOrigin(particles, particle, searchUpToQuark); if (origin == RecoDecay::OriginType::Prompt || origin == RecoDecay::OriginType::NonPrompt) { // 1=charm , 2=beauty hfparticle = particle; if (origin == RecoDecay::OriginType::Prompt) { @@ -1063,7 +1062,7 @@ void analyzeJetTrackInfo4ML(AnalysisJet const& analysisJet, AnyTracks const& /*a tracksParams.emplace_back(BJetTrackParams{constituent.pt(), constituent.eta(), dotProduct, dotProduct / analysisJet.p(), deltaRJetTrack, std::abs(constituent.dcaXY()) * sign, constituent.sigmadcaXY(), std::abs(constituent.dcaZ()) * sign, constituent.sigmadcaZ(), std::abs(constituent.dcaXYZ()) * sign, constituent.sigmadcaXYZ(), constituent.p() / analysisJet.p(), rClosestSV}); } - auto compare = [](BJetTrackParams& tr1, BJetTrackParams& tr2) { + auto compare = [](const BJetTrackParams& tr1, const BJetTrackParams& tr2) { return (tr1.signedIP2D / tr1.signedIP2DSign) > (tr2.signedIP2D / tr2.signedIP2DSign); }; @@ -1088,7 +1087,7 @@ void analyzeJetTrackInfo4MLnoSV(AnalysisJet const& analysisJet, AnyTracks const& tracksParams.emplace_back(BJetTrackParams{constituent.pt(), constituent.eta(), dotProduct, dotProduct / analysisJet.p(), deltaRJetTrack, std::abs(constituent.dcaXY()) * sign, constituent.sigmadcaXY(), std::abs(constituent.dcaZ()) * sign, constituent.sigmadcaZ(), std::abs(constituent.dcaXYZ()) * sign, constituent.sigmadcaXYZ(), constituent.p() / analysisJet.p(), 0.0}); } - auto compare = [](BJetTrackParams& tr1, BJetTrackParams& tr2) { + auto compare = [](const BJetTrackParams& tr1, const BJetTrackParams& tr2) { return (tr1.signedIP2D / tr1.signedIP2DSign) > (tr2.signedIP2D / tr2.signedIP2DSign); }; diff --git a/PWGJE/Core/JetUtilities.h b/PWGJE/Core/JetUtilities.h index 1f7ee8f9791..e0e7090cc9d 100644 --- a/PWGJE/Core/JetUtilities.h +++ b/PWGJE/Core/JetUtilities.h @@ -172,9 +172,6 @@ std::tuple estimateRhoPerpCone(const T& inputParticles, const U& double perpMdDensity1 = 0; double perpMdDensity2 = 0; - double dPhi1 = 999.; - double dPhi2 = 999.; - double dEta = 999.; double PerpendicularConeAxisPhi1 = 999., PerpendicularConeAxisPhi2 = 999.; // build 2 perp cones in phi around the leading jet (right and left of the jet) @@ -183,11 +180,11 @@ std::tuple estimateRhoPerpCone(const T& inputParticles, const U& for (auto& particle : inputParticles) { // sum the momentum of all paricles that fill the two cones - dPhi1 = particle.phi() - PerpendicularConeAxisPhi1; + double dPhi1 = particle.phi() - PerpendicularConeAxisPhi1; dPhi1 = RecoDecay::constrainAngle(dPhi1, -M_PI); // This will contrain the angel between -pi & Pi - dPhi2 = particle.phi() - PerpendicularConeAxisPhi2; + double dPhi2 = particle.phi() - PerpendicularConeAxisPhi2; dPhi2 = RecoDecay::constrainAngle(dPhi2, -M_PI); // This will contrain the angel between -pi & Pi - dEta = jet.eta() - particle.eta(); // The perp cone eta is the same as the leading jet since the cones are perpendicular only in phi + double dEta = jet.eta() - particle.eta(); // The perp cone eta is the same as the leading jet since the cones are perpendicular only in phi if (TMath::Sqrt(dPhi1 * dPhi1 + dEta * dEta) <= static_cast(perpConeR)) { perpPtDensity1 += particle.perp(); perpMdDensity1 += TMath::Sqrt(particle.m() * particle.m() + particle.pt() * particle.pt()) - particle.pt(); diff --git a/PWGJE/Core/JetV0Utilities.h b/PWGJE/Core/JetV0Utilities.h index 0a8e6454294..61a8dd9ec86 100644 --- a/PWGJE/Core/JetV0Utilities.h +++ b/PWGJE/Core/JetV0Utilities.h @@ -160,7 +160,7 @@ bool selectV0ParticleDecay(T const& v0Particle, int v0ParticleDecaySelection = - return (v0Particle.decayFlag() & (1 << v0ParticleDecaySelection)); } -int initialiseV0ParticleDecaySelection(std::string v0ParticleDecaySelection) +int initialiseV0ParticleDecaySelection(const std::string& v0ParticleDecaySelection) { if (v0ParticleDecaySelection == "K0sToPiPi") { return JV0ParticleDecays::K0sToPiPi; diff --git a/PWGJE/JetFinders/jetFinder.cxx b/PWGJE/JetFinders/jetFinder.cxx index 43670879a30..ff7e67966a5 100644 --- a/PWGJE/JetFinders/jetFinder.cxx +++ b/PWGJE/JetFinders/jetFinder.cxx @@ -150,7 +150,6 @@ struct JetFinderTask { } else { jetRadiiBins.push_back(jetRadiiBins[jetRadiiBins.size() - 1] + 0.1); } - std::vector jetPtBins; int jetPtMaxInt = static_cast(jetPtMax); int jetPtMinInt = static_cast(jetPtMin); jetPtMinInt = (jetPtMinInt / jetPtBinWidth) * jetPtBinWidth; diff --git a/PWGJE/JetFinders/jetFinderHF.cxx b/PWGJE/JetFinders/jetFinderHF.cxx index 2cc269d33cc..c80688aff45 100644 --- a/PWGJE/JetFinders/jetFinderHF.cxx +++ b/PWGJE/JetFinders/jetFinderHF.cxx @@ -157,7 +157,6 @@ struct JetFinderHFTask { } else { jetRadiiBins.push_back(jetRadiiBins[jetRadiiBins.size() - 1] + 0.1); } - std::vector jetPtBins; int jetPtMaxInt = static_cast(jetPtMax); int jetPtMinInt = static_cast(jetPtMin); jetPtMinInt = (jetPtMinInt / jetPtBinWidth) * jetPtBinWidth; diff --git a/PWGJE/JetFinders/jetFinderV0.cxx b/PWGJE/JetFinders/jetFinderV0.cxx index 651f05b2061..75bb1071688 100644 --- a/PWGJE/JetFinders/jetFinderV0.cxx +++ b/PWGJE/JetFinders/jetFinderV0.cxx @@ -148,7 +148,6 @@ struct JetFinderV0Task { } else { jetRadiiBins.push_back(jetRadiiBins[jetRadiiBins.size() - 1] + 0.1); } - std::vector jetPtBins; int jetPtMaxInt = static_cast(jetPtMax); int jetPtMinInt = static_cast(jetPtMin); jetPtMinInt = (jetPtMinInt / jetPtBinWidth) * jetPtBinWidth; diff --git a/PWGJE/TableProducer/derivedDataProducer.cxx b/PWGJE/TableProducer/derivedDataProducer.cxx index 4091715b441..61c9240075d 100644 --- a/PWGJE/TableProducer/derivedDataProducer.cxx +++ b/PWGJE/TableProducer/derivedDataProducer.cxx @@ -419,8 +419,8 @@ struct JetDerivedDataProducerTask { } } int daughtersId[2] = {-1, -1}; - auto i = 0; if (particle.has_daughters()) { + auto i = 0; for (auto daughterId : particle.daughtersIds()) { if (i > 1) { break; @@ -443,13 +443,11 @@ struct JetDerivedDataProducerTask { float leadingCellEnergy = -1.0; float subleadingCellEnergy = -1.0; - float cellAmplitude = -1.0; int leadingCellNumber = -1; int subleadingCellNumber = -1; - int cellNumber = -1; for (auto const& clutserCell : clusterCells) { - cellAmplitude = clutserCell.calo().amplitude(); - cellNumber = clutserCell.calo().cellNumber(); + float cellAmplitude = clutserCell.calo().amplitude(); + int cellNumber = clutserCell.calo().cellNumber(); if (cellAmplitude > subleadingCellEnergy) { subleadingCellEnergy = cellAmplitude; subleadingCellNumber = cellNumber; @@ -779,8 +777,8 @@ struct JetDerivedDataProducerTask { } } int daughtersId[2] = {-1, -1}; - auto i = 0; if (particle.has_daughters()) { + auto i = 0; for (auto daughterId : particle.daughtersIds()) { if (i > 1) { break; @@ -833,8 +831,8 @@ struct JetDerivedDataProducerTask { } } int daughtersId[2] = {-1, -1}; - auto i = 0; if (particle.has_daughters()) { + auto i = 0; for (auto daughterId : particle.daughtersIds()) { if (i > 1) { break; diff --git a/PWGJE/TableProducer/derivedDataWriter.cxx b/PWGJE/TableProducer/derivedDataWriter.cxx index 72864b14e58..1c3125a68ff 100644 --- a/PWGJE/TableProducer/derivedDataWriter.cxx +++ b/PWGJE/TableProducer/derivedDataWriter.cxx @@ -433,7 +433,7 @@ struct JetDerivedDataWriter { std::vector lcMcCollisionMapping; std::vector b0McCollisionMapping; std::vector bplusMcCollisionMapping; - std::vector dielectronMcCollisionMapping; + // std::vector dielectronMcCollisionMapping; void processBCs(soa::Join const& collisions, soa::Join const& bcs) { @@ -672,8 +672,8 @@ struct JetDerivedDataWriter { } } int daughtersIds[2] = {-1, -1}; - auto i = 0; if (particle.has_daughters()) { + auto i = 0; for (auto daughterId : particle.daughtersIds()) { if (i > 1) { break; @@ -885,8 +885,8 @@ struct JetDerivedDataWriter { DielectronMothersIds.push_back(particleMapping[DielectronMother.globalIndex()]); } } - auto i = 0; if (DielectronParticle.has_daughters()) { + auto i = 0; for (auto const& DielectronDaughter : DielectronParticle.template daughters_as()) { if (i > 1) { break; diff --git a/PWGJE/TableProducer/eventwiseConstituentSubtractor.cxx b/PWGJE/TableProducer/eventwiseConstituentSubtractor.cxx index 720d4081926..7dbfa127ee5 100644 --- a/PWGJE/TableProducer/eventwiseConstituentSubtractor.cxx +++ b/PWGJE/TableProducer/eventwiseConstituentSubtractor.cxx @@ -82,7 +82,6 @@ struct eventWiseConstituentSubtractorTask { Configurable doRhoMassSub{"doRhoMassSub", true, "perfom mass subtraction as well"}; JetBkgSubUtils eventWiseConstituentSubtractor; - float bkgPhiMax_; std::vector inputParticles; std::vector tracksSubtracted; int trackSelection = -1; diff --git a/PWGJE/TableProducer/jetTaggerHF.cxx b/PWGJE/TableProducer/jetTaggerHF.cxx index 44899e4c702..b648929de89 100644 --- a/PWGJE/TableProducer/jetTaggerHF.cxx +++ b/PWGJE/TableProducer/jetTaggerHF.cxx @@ -276,8 +276,8 @@ struct JetTaggerHFTask { } } if (doprocessAlgorithmGNN) { - float dbRange; if (jet.pt() >= jetpTMin) { + float dbRange; if (scoreML[jet.globalIndex()] < dbMin) { dbRange = 0.5; // underflow } else if (scoreML[jet.globalIndex()] < dbMax) { @@ -571,7 +571,7 @@ struct JetTaggerHFTask { } if (bMlResponse.getOutputNodes() > 1) { - auto mDb = [](std::vector scores, float fC) { + auto mDb = [](const std::vector& scores, float fC) { return std::log(scores[2] / (fC * scores[1] + (1 - fC) * scores[0])); }; diff --git a/PWGJE/Tasks/bjetTaggingML.cxx b/PWGJE/Tasks/bjetTaggingML.cxx index 64536cb0c15..6e1b135b745 100644 --- a/PWGJE/Tasks/bjetTaggingML.cxx +++ b/PWGJE/Tasks/bjetTaggingML.cxx @@ -269,7 +269,7 @@ struct BJetTaggingML { tracksParams.emplace_back(jettaggingutilities::BJetTrackParams{constituent.pt(), constituent.eta(), dotProduct, dotProduct / analysisJet.p(), deltaRJetTrack, std::abs(constituent.dcaXY()) * sign, constituent.sigmadcaXY(), std::abs(constituent.dcaXYZ()) * sign, constituent.sigmadcaXYZ(), constituent.p() / analysisJet.p(), rClosestSV}); } - auto compare = [](jettaggingutilities::BJetTrackParams& tr1, jettaggingutilities::BJetTrackParams& tr2) { + auto compare = [](const jettaggingutilities::BJetTrackParams& tr1, const jettaggingutilities::BJetTrackParams& tr2) { return (tr1.signedIP2D / tr1.signedIP2DSign) > (tr2.signedIP2D / tr2.signedIP2DSign); }; diff --git a/PWGJE/Tasks/bjetTreeCreator.cxx b/PWGJE/Tasks/bjetTreeCreator.cxx index c42453f6817..ee059aa9626 100644 --- a/PWGJE/Tasks/bjetTreeCreator.cxx +++ b/PWGJE/Tasks/bjetTreeCreator.cxx @@ -755,7 +755,6 @@ struct BJetTreeCreator { } std::vector indicesTracks; - std::vector indicesSVs; int16_t jetFlavor = analysisJet.origin(); @@ -792,6 +791,7 @@ struct BJetTreeCreator { } if (produceTree) { + std::vector indicesSVs; bjetConstituentsTable(bjetParamsTable.lastIndex() + 1, indicesTracks, indicesSVs); bjetParamsTable(analysisJet.pt(), analysisJet.eta(), analysisJet.phi(), indicesTracks.size(), nVertices, analysisJet.mass(), jetFlavor, analysisJet.r()); } diff --git a/PWGJE/Tasks/dijetFinderQA.cxx b/PWGJE/Tasks/dijetFinderQA.cxx index b5410727b7c..5c9f39f1a84 100644 --- a/PWGJE/Tasks/dijetFinderQA.cxx +++ b/PWGJE/Tasks/dijetFinderQA.cxx @@ -238,12 +238,12 @@ struct DijetFinderQATask { } if (jetPtcuts.size() >= 2) { - auto& leading_jet = jetPtcuts[0]; + const auto& leading_jet = jetPtcuts[0]; bool found_pair = false; for (size_t i = 1; i < jetPtcuts.size() && !found_pair; i++) { - auto& candidate_jet = jetPtcuts[i]; + const auto& candidate_jet = jetPtcuts[i]; Double_t dphi = fabs(leading_jet[2] - candidate_jet[2]); Double_t deta = fabs(leading_jet[1] - candidate_jet[1]); Double_t condition = fabs(dphi - M_PI); @@ -288,12 +288,12 @@ struct DijetFinderQATask { } if (jetPtcuts.size() >= 2) { - auto& leading_jet = jetPtcuts[0]; + const auto& leading_jet = jetPtcuts[0]; bool found_pair = false; for (size_t i = 1; i < jetPtcuts.size() && !found_pair; i++) { - auto& candidate_jet = jetPtcuts[i]; + const auto& candidate_jet = jetPtcuts[i]; Double_t dphi = fabs(leading_jet[2] - candidate_jet[2]); Double_t deta = fabs(leading_jet[1] - candidate_jet[1]); Double_t condition = fabs(dphi - M_PI); @@ -338,12 +338,12 @@ struct DijetFinderQATask { } if (jetPtcuts.size() >= 2) { - auto& leading_jet = jetPtcuts[0]; + const auto& leading_jet = jetPtcuts[0]; bool found_pair = false; for (size_t i = 1; i < jetPtcuts.size() && !found_pair; i++) { - auto& candidate_jet = jetPtcuts[i]; + const auto& candidate_jet = jetPtcuts[i]; Double_t dphi = fabs(leading_jet[2] - candidate_jet[2]); Double_t deta = fabs(leading_jet[1] - candidate_jet[1]); Double_t condition = fabs(dphi - M_PI); @@ -411,8 +411,8 @@ struct DijetFinderQATask { } if (jetPtcuts_D.size() >= 2 && jetPtcuts_P.size() >= 2) { - auto& leading_jet_D = jetPtcuts_D[0]; - auto& leading_jet_P = jetPtcuts_P[0]; + const auto& leading_jet_D = jetPtcuts_D[0]; + const auto& leading_jet_P = jetPtcuts_P[0]; std::array candidate_jet_D{}; std::array candidate_jet_P{}; diff --git a/PWGJE/Tasks/fullJetSpectra.cxx b/PWGJE/Tasks/fullJetSpectra.cxx index dc831567a0e..cde5bb312de 100644 --- a/PWGJE/Tasks/fullJetSpectra.cxx +++ b/PWGJE/Tasks/fullJetSpectra.cxx @@ -121,7 +121,7 @@ struct FullJetSpectra { Configurable pTHatAbsoluteMin{"pTHatAbsoluteMin", -99.0, "minimum value of pTHat"}; int trackSelection = -1; - const float kJetAreaFractionMinThreshold = -98.0f; + // const float kJetAreaFractionMinThreshold = -98.0f; const float kLeadingTrackPtMinThreshold = -98.0f; const float kLeadingTrackPtMaxThreshold = 9998.0f; const float kLeadingClusterPtMinThreshold = -98.0f; @@ -775,7 +775,6 @@ struct FullJetSpectra { void fillJetHistograms(T const& jet, float weight = 1.0) { // float neutralEnergy = 0.0; - double sumtrackE = 0.0; if (jet.r() == round(selectedJetsRadius * 100.0f)) { registry.fill(HIST("h_full_jet_pt"), jet.pt(), weight); registry.fill(HIST("h_full_jet_eta"), jet.eta(), weight); @@ -847,7 +846,7 @@ struct FullJetSpectra { registry.fill(HIST("h2_full_jet_nef_corr_oneTrack70"), jet.pt(), nef_corr_oneTrack70, weight); registry.fill(HIST("h2_full_jet_nef_corr_allTracks100"), jet.pt(), nef_corr_allTracks100, weight); registry.fill(HIST("h2_full_jet_nef_corr_allTracks70"), jet.pt(), nef_corr_allTracks70, weight); - + double sumtrackE = 0; for (const auto& jettrack : jet.template tracks_as()) { sumtrackE += jettrack.energy(); @@ -871,8 +870,8 @@ struct FullJetSpectra { template void fillRejectedJetHistograms(T const& jet, float weight = 1.0) { - float neutralEnergy = 0.0; if (jet.r() == round(selectedJetsRadius * 100.0f)) { + float neutralEnergy = 0.0; for (const auto& cluster : jet.template clusters_as()) { neutralEnergy += cluster.energy(); } @@ -884,12 +883,6 @@ struct FullJetSpectra { template void fillMCPHistograms(T const& jet, float weight = 1.0) { - float neutralEnergy = 0.0; - int neutralconsts = 0; - int chargedconsts = 0; - int mcpjetOutsideFid = 0; - int mcpjetInsideFid = 0; - auto isInFiducial = [&](auto const& jet) { return jet.eta() >= jetEtaMin && jet.eta() <= jetEtaMax && jet.phi() >= jetPhiMin && jet.phi() <= jetPhiMax; @@ -905,18 +898,19 @@ struct FullJetSpectra { if (!isInFiducial(jet)) { // jet is outside - mcpjetOutsideFid++; - registry.fill(HIST("h2_full_mcpjetOutsideFiducial_pt"), jet.pt(), mcpjetOutsideFid, weight); + registry.fill(HIST("h2_full_mcpjetOutsideFiducial_pt"), jet.pt(), 1, weight); registry.fill(HIST("h_full_mcpjetOutside_eta_part"), jet.eta(), weight); registry.fill(HIST("h_full_mcpjetOutside_phi_part"), jet.phi(), weight); } else { // jet is inside - mcpjetInsideFid++; - registry.fill(HIST("h2_full_mcpjetInsideFiducial_pt"), jet.pt(), mcpjetInsideFid, weight); + registry.fill(HIST("h2_full_mcpjetInsideFiducial_pt"), jet.pt(), 1, weight); registry.fill(HIST("h_full_mcpjetInside_eta_part"), jet.eta(), weight); registry.fill(HIST("h_full_mcpjetInside_phi_part"), jet.phi(), weight); } + float neutralEnergy = 0.0; + int neutralconsts = 0; + int chargedconsts = 0; for (const auto& constituent : jet.template tracks_as()) { auto pdgParticle = pdgDatabase->GetParticle(constituent.pdgCode()); if (pdgParticle->Charge() == 0) { @@ -2177,9 +2171,7 @@ struct FullJetSpectra { // std::cout << "******************************************** " << std::endl; if (selectedJets.size() == 0) { // no jets = no leading jet return; - } - - if (selectedJets.size() > 0) { + } else { // Jet multiplicity per event registry.fill(HIST("h_all_fulljet_Njets"), selectedJets.size(), 1.0); @@ -2346,9 +2338,7 @@ struct FullJetSpectra { if (selectedJets.size() == 0) { // no jets = no leading jet return; - } - - if (selectedJets.size() > 0) { + } else { // Jet multiplicity per event registry.fill(HIST("h_all_fulljet_Njets"), selectedJets.size(), 1.0); @@ -2506,8 +2496,7 @@ struct FullJetSpectra { if (selectedJets.size() == 0) { // no jets = no leading jet return; - } - if (selectedJets.size() > 0) { + } else { // Jet multiplicity per event registry.fill(HIST("h_all_fulljet_Njets"), selectedJets.size(), eventWeight); @@ -2669,9 +2658,7 @@ struct FullJetSpectra { if (selectedJets.size() == 0) { // no jets = no leading jet return; - } - - if (selectedJets.size() > 0) { + } else { // Jet multiplicity per event registry.fill(HIST("h_all_fulljet_Njets_part"), selectedJets.size(), 1.0); @@ -2833,9 +2820,7 @@ struct FullJetSpectra { if (selectedJets.size() == 0) { // no jets = no leading jet return; - } - - if (selectedJets.size() > 0) { + } else { // Jet multiplicity per event registry.fill(HIST("h_all_fulljet_Njets_part"), selectedJets.size(), mccollision.weight()); diff --git a/PWGJE/Tasks/fullJetTriggerQATask.cxx b/PWGJE/Tasks/fullJetTriggerQATask.cxx index 1831a3c3372..88255a2dda7 100644 --- a/PWGJE/Tasks/fullJetTriggerQATask.cxx +++ b/PWGJE/Tasks/fullJetTriggerQATask.cxx @@ -109,9 +109,6 @@ struct JetTriggerQA { Configurable cfgVertexCut{"cfgVertexCut", 10.0f, "Accepted z-vertex range"}; Configurable mClusterDefinition{"clusterDefinition", "kV3Default", "cluster definition to be selected, e.g. V3Default"}; - std::vector jetConstituents; - std::vector jetClusterConstituents; - std::vector jetReclustered; JetFinder jetReclusterer; void init(InitContext&) @@ -450,7 +447,7 @@ struct JetTriggerQA { registry.fill(HIST("hTriggerCorrelation"), mainTrigger, assocTrigger); } - std::pair fillGammaQA(const selectedClusters& clusters, std::bitset hwtrg, const std::bitset& triggerstatus) + std::pair fillGammaQA(const selectedClusters& clusters, const std::bitset& hwtrg, const std::bitset& triggerstatus) { auto isTrigger = [&triggerstatus](TriggerType_t triggertype) -> bool { return triggerstatus.test(triggertype); @@ -461,8 +458,6 @@ struct JetTriggerQA { struct ClusterData { float mTriggerObservable; - float mEta; - float mPhi; bool mEMCALcluster; }; std::vector analysedClusters; @@ -476,7 +471,7 @@ struct JetTriggerQA { } bool emcalCluster = isClusterInEmcal(cluster); double clusterObservable = (f_GammaObservable == 0) ? cluster.energy() : cluster.energy() / std::cosh(cluster.eta()); - analysedClusters.push_back({static_cast(clusterObservable), cluster.eta(), cluster.phi(), emcalCluster}); + analysedClusters.push_back({static_cast(clusterObservable), emcalCluster}); if (emcalCluster && (clusterObservable > maxClusterObservableEMCAL)) { maxClusterObservableEMCAL = clusterObservable; @@ -889,8 +884,7 @@ struct JetTriggerQA { registry.get(HIST("jetRMaxPtEtaPhiNoFiducial"))->Fill(jetR, jetPt, jetEta, jetPhi); if (maxJet.r() == std::round(f_jetR * 100)) { for (const auto& jet : jets) { - if (!b_doLightOutput) - registry.fill(HIST("hJetRMaxPtJetPtNoFiducial"), jet.r() * 1e-2, jetPt, jet.pt()); + registry.fill(HIST("hJetRMaxPtJetPtNoFiducial"), jet.r() * 1e-2, jetPt, jet.pt()); } // for jets } // if maxJet.r() == std::round(f_jetR * 100) } // for maxjet no fiducial diff --git a/PWGJE/Tasks/gammaJetTreeProducer.cxx b/PWGJE/Tasks/gammaJetTreeProducer.cxx index a054a8e97b8..3f31fbb244e 100644 --- a/PWGJE/Tasks/gammaJetTreeProducer.cxx +++ b/PWGJE/Tasks/gammaJetTreeProducer.cxx @@ -914,8 +914,6 @@ struct GammaJetTreeProducer { } if (nRecCollisions > 1) { mHistograms.fill(HIST("mcCollisionsWithRecCollisions"), 2); - } - if (nRecCollisions > 1) { mcCollisionsMultiRecCollisions.push_back(mcCollision.globalIndex()); } diff --git a/PWGJE/Tasks/jetChCorr.cxx b/PWGJE/Tasks/jetChCorr.cxx index 4311fc5c7d8..e9a4a937334 100644 --- a/PWGJE/Tasks/jetChCorr.cxx +++ b/PWGJE/Tasks/jetChCorr.cxx @@ -174,8 +174,6 @@ struct JetChCorr { fastjet::PseudoJet parentSubJet2; bool softDropped = false; auto nsd = 0.0; - auto zg = -1.0; - auto rg = -1.0; // std::vector energyMotherVec; // std::vector ptLeadingVec; @@ -219,8 +217,8 @@ struct JetChCorr { if (z >= zCut * TMath::Power(theta / (jet.r() / 100.f), beta)) { if (found1 == true && found2 == true) { // found leading and next-to-leading in seperate prongs if (!softDropped) { - zg = z; - rg = theta; + auto zg = z; + auto rg = theta; if constexpr (!isSubtracted && !isMCP) { registry.fill(HIST("h2_jet_pt_jet_zg"), jet.pt(), zg); registry.fill(HIST("h2_jet_pt_jet_rg"), jet.pt(), rg); @@ -239,7 +237,7 @@ struct JetChCorr { v2.SetXYZ(parentSubJet2.px(), parentSubJet2.py(), parentSubJet2.pz()); vR = v1 + v2; - float z = v2.Perp(vR.Orthogonal()) / (v1.Perp(vR.Orthogonal()) + v2.Perp(vR.Orthogonal())); + float z = v2.Perp(vR.Orthogonal()) / (v1.Perp(vR.Orthogonal()) + v2.Perp(vR.Orthogonal())); // you redefine z here, please check! float fT = ((2. * z * (1 - z) * vR.Mag()) / v1.Perp2(vR)) / 6.; float kt_p = v1.Perp(vR); diff --git a/PWGJE/Tasks/jetChargedV2.cxx b/PWGJE/Tasks/jetChargedV2.cxx index 938275ca5ce..bd5fbe25107 100644 --- a/PWGJE/Tasks/jetChargedV2.cxx +++ b/PWGJE/Tasks/jetChargedV2.cxx @@ -584,17 +584,14 @@ struct JetChargedV2 { fFitModulationV2v3P->SetParameter(1, 0.01); fFitModulationV2v3P->SetParameter(3, 0.01); - double ep2fix = 0.; - double ep3fix = 0.; - if (ep2 < 0) { - ep2fix = RecoDecay::constrainAngle(ep2); + double ep2fix = RecoDecay::constrainAngle(ep2); fFitModulationV2v3P->FixParameter(2, ep2fix); } else { fFitModulationV2v3P->FixParameter(2, ep2); } if (ep3 < 0) { - ep3fix = RecoDecay::constrainAngle(ep3); + double ep3fix = RecoDecay::constrainAngle(ep3); fFitModulationV2v3P->FixParameter(4, ep3fix); } else { fFitModulationV2v3P->FixParameter(4, ep3); @@ -1187,17 +1184,14 @@ struct JetChargedV2 { fFitModulationV2v3->SetParameter(1, 0.01); fFitModulationV2v3->SetParameter(3, 0.01); - double ep2fix = 0.; - double ep3fix = 0.; - if (ep2 < 0) { - ep2fix = RecoDecay::constrainAngle(ep2); + double ep2fix = RecoDecay::constrainAngle(ep2); fFitModulationV2v3->FixParameter(2, ep2fix); } else { fFitModulationV2v3->FixParameter(2, ep2); } if (ep3 < 0) { - ep3fix = RecoDecay::constrainAngle(ep3); + double ep3fix = RecoDecay::constrainAngle(ep3); fFitModulationV2v3->FixParameter(4, ep3fix); } else { fFitModulationV2v3->FixParameter(4, ep3); @@ -1241,12 +1235,12 @@ struct JetChargedV2 { chiSqr = chi2; cDF = 1. - chiSquareCDF(nDF, chiSqr); - int evtCentAreaMin = 0; - int evtCentAreaMax = 5; - int evtMidAreaMin = 30; - int evtMidAreaMax = 50; double evtcent = collision.centFT0M(); if (cfgChkFitQuality) { + int evtCentAreaMin = 0; + int evtCentAreaMax = 5; + int evtMidAreaMin = 30; + int evtMidAreaMax = 50; registry.fill(HIST("h_PvalueCDF_CombinFit"), cDF); registry.fill(HIST("h2_PvalueCDFCent_CombinFit"), collision.centFT0M(), cDF); registry.fill(HIST("h2_Chi2Cent_CombinFit"), collision.centFT0M(), chiSqr / (static_cast(nDF))); @@ -1322,7 +1316,6 @@ struct JetChargedV2 { TRandom3 randomNumber(0); float randomConeEta = randomNumber.Uniform(trackEtaMin + randomConeR, trackEtaMax - randomConeR); float randomConePhi = randomNumber.Uniform(0.0, o2::constants::math::TwoPI); - float randomConePt = 0; double integralValueRC = fFitModulationV2v3->Integral(randomConePhi - randomConeR, randomConePhi + randomConeR); double rholocalRC = collision.rho() / (2 * randomConeR * temppara[0]) * integralValueRC; @@ -1330,6 +1323,7 @@ struct JetChargedV2 { if (nmode == cfgNmodA) { double rcPhiPsi2; rcPhiPsi2 = randomConePhi - ep2; + float randomConePt = 0; for (auto const& track : tracks) { if (jetderiveddatautilities::selectTrack(track, trackSelection)) { @@ -1444,17 +1438,14 @@ struct JetChargedV2 { fFitModulationV2v3->SetParameter(1, 0.01); fFitModulationV2v3->SetParameter(3, 0.01); - double ep2fix = 0.; - double ep3fix = 0.; - if (ep2 < 0) { - ep2fix = RecoDecay::constrainAngle(ep2); + double ep2fix = RecoDecay::constrainAngle(ep2); fFitModulationV2v3->FixParameter(2, ep2fix); } else { fFitModulationV2v3->FixParameter(2, ep2); } if (ep3 < 0) { - ep3fix = RecoDecay::constrainAngle(ep3); + double ep3fix = RecoDecay::constrainAngle(ep3); fFitModulationV2v3->FixParameter(4, ep3fix); } else { fFitModulationV2v3->FixParameter(4, ep3); @@ -1498,12 +1489,12 @@ struct JetChargedV2 { chiSqr = chi2; cDF = 1. - chiSquareCDF(nDF, chiSqr); - int evtCentAreaMin = 0; - int evtCentAreaMax = 5; - int evtMidAreaMin = 30; - int evtMidAreaMax = 50; double evtcent = collision.centFT0M(); if (cfgChkFitQuality) { + int evtCentAreaMin = 0; + int evtCentAreaMax = 5; + int evtMidAreaMin = 30; + int evtMidAreaMax = 50; registry.fill(HIST("h_PvalueCDF_CombinFit"), cDF); registry.fill(HIST("h2_PvalueCDFCent_CombinFit"), collision.centFT0M(), cDF); registry.fill(HIST("h2_Chi2Cent_CombinFit"), collision.centFT0M(), chiSqr / (static_cast(nDF))); @@ -1573,7 +1564,6 @@ struct JetChargedV2 { TRandom3 randomNumber(0); float randomConeEta = randomNumber.Uniform(trackEtaMin + randomConeR, trackEtaMax - randomConeR); float randomConePhi = randomNumber.Uniform(0.0, o2::constants::math::TwoPI); - float randomConePt = 0; double integralValueRC = fFitModulationV2v3->Integral(randomConePhi - randomConeR, randomConePhi + randomConeR); double rholocalRC = collision.rho() / (2 * randomConeR * temppara[0]) * integralValueRC; @@ -1581,6 +1571,7 @@ struct JetChargedV2 { if (nmode == cfgNmodA) { double rcPhiPsi2; rcPhiPsi2 = randomConePhi - ep2; + float randomConePt = 0; for (auto const& track : tracks) { if (jetderiveddatautilities::selectTrack(track, trackSelection)) { @@ -1912,17 +1903,14 @@ struct JetChargedV2 { fFitModulationRM->SetParameter(1, 0.01); fFitModulationRM->SetParameter(3, 0.01); - double ep2fix = 0.; - double ep3fix = 0.; - if (ep2 < 0) { - ep2fix = RecoDecay::constrainAngle(ep2); + double ep2fix = RecoDecay::constrainAngle(ep2); fFitModulationRM->FixParameter(2, ep2fix); } else { fFitModulationRM->FixParameter(2, ep2); } if (ep3 < 0) { - ep3fix = RecoDecay::constrainAngle(ep3); + double ep3fix = RecoDecay::constrainAngle(ep3); fFitModulationRM->FixParameter(4, ep3fix); } else { fFitModulationRM->FixParameter(4, ep3); diff --git a/PWGJE/Tasks/jetFormationTimeReclustering.cxx b/PWGJE/Tasks/jetFormationTimeReclustering.cxx index 6e1b4fcc859..f3f2348f7d5 100644 --- a/PWGJE/Tasks/jetFormationTimeReclustering.cxx +++ b/PWGJE/Tasks/jetFormationTimeReclustering.cxx @@ -305,9 +305,6 @@ struct FormationTimeReclustering { fastjet::PseudoJet parentSubJet2; bool softDropped = false; auto nsd = 0.0; - auto zg = -1.0; - auto rg = -1.0; - auto tg = -1.0; while (daughterSubJet.has_parents(parentSubJet1, parentSubJet2)) { if (parentSubJet1.perp() < parentSubJet2.perp()) { @@ -334,9 +331,9 @@ struct FormationTimeReclustering { if (z >= zCut * TMath::Power(theta / (jet.r() / 100.f), beta)) { if (!softDropped) { - zg = z; - rg = theta; - tg = tau; + auto zg = z; + auto rg = theta; + auto tg = tau; ptgVec.push_back(jet.pt()); thetagVec.push_back(rg); taugVec.push_back(tg); diff --git a/PWGJE/Tasks/jetFragmentation.cxx b/PWGJE/Tasks/jetFragmentation.cxx index 962f962fc6a..9abaccdfaa2 100644 --- a/PWGJE/Tasks/jetFragmentation.cxx +++ b/PWGJE/Tasks/jetFragmentation.cxx @@ -1101,7 +1101,7 @@ struct JetFragmentation { // Implementation of background subtraction at runtime // --------------------------------------------------- - int getPtBin(float pt, std::vector ptBins) + int getPtBin(float pt, const std::vector& ptBins) { if (pt < ptBins.at(0)) return -1; @@ -1216,7 +1216,7 @@ struct JetFragmentation { return v; } // Return the probability associated with a given outcome - double stateWeight(std::vector state, std::vector> weights) + double stateWeight(const std::vector& state, const std::vector>& weights) { double w = 1.; for (uint32_t ip = 0; ip < state.size(); ip++) { @@ -1226,7 +1226,7 @@ struct JetFragmentation { } // Return the corrected values for z and ptjet for a given state // Scale values by the fraction of jet momentum carried by removed V0s - std::vector correctedValues(std::vector state, std::vector values) + std::vector correctedValues(const std::vector& state, const std::vector& values) { // Assumes values = (z1, z2, ..., zn, ptjet) std::vector v(values); @@ -1253,7 +1253,7 @@ struct JetFragmentation { // Return the corrected values for z and ptjet for a given state // Take into account tracks that would have been included in the jet regardless of the V0s template - std::vector correctedValuesPlusTracks(std::vector state, V const& jet) + std::vector correctedValuesPlusTracks(const std::vector& state, V const& jet) { int ip = 0; double ptToSubtract = 0., ptToAdd = 0.; @@ -1710,7 +1710,7 @@ struct JetFragmentation { registry.fill(HIST("data/V0/nV0sEventAccWeighted"), nV0s); } template - void fillDataV0sInJetWeighted(C const& coll, J const& jet, std::vector state, std::vector values, double weight) + void fillDataV0sInJetWeighted(C const& coll, J const& jet, const std::vector& state, const std::vector& values, double weight) { double jetpt = values[values.size() - 1]; int ip = 0; diff --git a/PWGJE/Tasks/jetHadronRecoil.cxx b/PWGJE/Tasks/jetHadronRecoil.cxx index c2c086cbdce..6ed91cdbbda 100644 --- a/PWGJE/Tasks/jetHadronRecoil.cxx +++ b/PWGJE/Tasks/jetHadronRecoil.cxx @@ -223,7 +223,6 @@ struct JetHadronRecoil { std::vector ptTTAr; double phiTT = 0; double ptTT = 0; - int trigNumber = 0; int nTT = 0; double leadingPT = 0; double leadingTrackPt = 0; @@ -269,7 +268,7 @@ struct JetHadronRecoil { registry.fill(HIST("hPtTrackPtHard"), track.pt() / pTHat, track.pt(), weight); } if (nTT > 0) { - trigNumber = rand->Integer(nTT); + int trigNumber = rand->Integer(nTT); phiTT = phiTTAr[trigNumber]; ptTT = ptTTAr[trigNumber]; if (isSigCol) { @@ -359,7 +358,6 @@ struct JetHadronRecoil { std::vector ptTTAr; double phiTT = 0; double ptTT = 0; - int trigNumber = 0; int nTT = 0; double leadingPT = 0; double leadingTrackPt = 0; @@ -412,7 +410,7 @@ struct JetHadronRecoil { } } if (nTT > 0) { - trigNumber = rand->Integer(nTT); + int trigNumber = rand->Integer(nTT); phiTT = phiTTAr[trigNumber]; ptTT = ptTTAr[trigNumber]; if (isSigCol) { @@ -502,7 +500,6 @@ struct JetHadronRecoil { std::vector ptTTAr; double phiTT = 0; double ptTT = 0; - int trigNumber = 0; int nTT = 0; double leadingPartPt = 0; double leadingJetPt = 0; @@ -550,7 +547,7 @@ struct JetHadronRecoil { } if (nTT > 0) { - trigNumber = rand->Integer(nTT); + int trigNumber = rand->Integer(nTT); phiTT = phiTTAr[trigNumber]; ptTT = ptTTAr[trigNumber]; if (isSigCol) { @@ -620,8 +617,6 @@ struct JetHadronRecoil { void fillMatchedHistograms(T const& jetsBase, U const&, X const& tracks, Y const& particles, float weight = 1.0, float rho = 0.0, float pTHat = 999.0) { for (const auto& jetBase : jetsBase) { - double dR = 0; - double dRp = 0; if (jetBase.pt() > pTHatMaxMCD * pTHat) { if (outlierRejectEvent) { @@ -631,7 +626,7 @@ struct JetHadronRecoil { } } - dR = getWTAaxisDifference(jetBase, tracks); + double dR = getWTAaxisDifference(jetBase, tracks); if (jetBase.has_matchedJetGeo()) { for (const auto& jetTag : jetBase.template matchedJetGeo_as>()) { @@ -643,7 +638,7 @@ struct JetHadronRecoil { } } - dRp = getWTAaxisDifference(jetTag, particles); + double dRp = getWTAaxisDifference(jetTag, particles); registry.fill(HIST("hPtMatched"), jetBase.pt() - (rho * jetBase.area()), jetTag.pt(), weight); registry.fill(HIST("hPhiMatched"), jetBase.phi(), jetTag.phi(), weight); @@ -666,7 +661,6 @@ struct JetHadronRecoil { std::vector phiTTArPart; double phiTT = 0; double phiTTPart = 0; - int trigNumber = 0; int nTT = 0; for (const auto& track : tracks) { @@ -699,7 +693,7 @@ struct JetHadronRecoil { } if (nTT > 0) { - trigNumber = rand->Integer(nTT); + int trigNumber = rand->Integer(nTT); phiTT = phiTTAr[trigNumber]; phiTTPart = phiTTArPart[trigNumber]; } else { @@ -707,8 +701,6 @@ struct JetHadronRecoil { } for (const auto& jetBase : jetsBase) { - double dR = 0; - double dRp = 0; if (jetBase.pt() > pTHatMaxMCD * pTHat) { if (outlierRejectEvent) { @@ -719,7 +711,7 @@ struct JetHadronRecoil { } float dphi = RecoDecay::constrainAngle(jetBase.phi() - phiTT); - dR = getWTAaxisDifference(jetBase, tracks); + double dR = getWTAaxisDifference(jetBase, tracks); if (jetBase.has_matchedJetGeo()) { for (const auto& jetTag : jetBase.template matchedJetGeo_as>()) { @@ -732,7 +724,7 @@ struct JetHadronRecoil { } float dphip = RecoDecay::constrainAngle(jetTag.phi() - phiTTPart); - dRp = getWTAaxisDifference(jetTag, particles); + double dRp = getWTAaxisDifference(jetTag, particles); registry.fill(HIST("hPhiMatched"), dphi, dphip, weight); registry.fill(HIST("hPhiMatched2d"), jetTag.phi(), jetTag.pt(), weight); registry.fill(HIST("hPhiResolution"), jetTag.pt(), dphip - dphi, weight); diff --git a/PWGJE/Tasks/jetPlanarFlow.cxx b/PWGJE/Tasks/jetPlanarFlow.cxx index 673d4eebffe..1edf34d8396 100644 --- a/PWGJE/Tasks/jetPlanarFlow.cxx +++ b/PWGJE/Tasks/jetPlanarFlow.cxx @@ -243,7 +243,7 @@ struct JetPlanarFlowTask { } uint8_t isInJet = 0; - for (auto& jetConstituentId : jet.tracksIds()) { + for (const auto& jetConstituentId : jet.tracksIds()) { if (track.globalIndex() == jetConstituentId) { isInJet = 1; break; diff --git a/PWGJE/Tasks/jetSpectraEseTask.cxx b/PWGJE/Tasks/jetSpectraEseTask.cxx index 007b2ea511d..7c64a6d5f44 100644 --- a/PWGJE/Tasks/jetSpectraEseTask.cxx +++ b/PWGJE/Tasks/jetSpectraEseTask.cxx @@ -538,7 +538,7 @@ struct JetSpectraEseTask { EventPlane procEP(EPCol const& vec) { constexpr std::array AmpCut{1e-8, 0.0}; - auto computeEP = [&AmpCut](std::vector vec, auto det, float n) { return vec[2] > AmpCut[det] ? (1.0 / n) * std::atan2(vec[1], vec[0]) : 999.; }; + auto computeEP = [&AmpCut](const std::vector& vec, auto det, float n) { return vec[2] > AmpCut[det] ? (1.0 / n) * std::atan2(vec[1], vec[0]) : 999.; }; std::map epMap; std::map ep3Map; auto vec1{qVecNoESE(vec)}; diff --git a/PWGJE/Tasks/jetSubstructure.cxx b/PWGJE/Tasks/jetSubstructure.cxx index e82a923f38f..4f22831fe49 100644 --- a/PWGJE/Tasks/jetSubstructure.cxx +++ b/PWGJE/Tasks/jetSubstructure.cxx @@ -160,8 +160,6 @@ struct JetSubstructureTask { fastjet::PseudoJet parentSubJet2; bool softDropped = false; auto nsd = 0.0; - auto zg = -1.0; - auto rg = -1.0; while (daughterSubJet.has_parents(parentSubJet1, parentSubJet2)) { if (parentSubJet1.perp() < parentSubJet2.perp()) { @@ -185,8 +183,8 @@ struct JetSubstructureTask { if (z >= zCut * TMath::Power(theta / (jet.r() / 100.f), beta)) { if (!softDropped) { - zg = z; - rg = theta; + auto zg = z; + auto rg = theta; if constexpr (!isSubtracted && !isMCP) { registry.fill(HIST("h2_jet_pt_jet_zg"), jet.pt(), zg); registry.fill(HIST("h2_jet_pt_jet_rg"), jet.pt(), rg); diff --git a/PWGJE/Tasks/jetSubstructureHF.cxx b/PWGJE/Tasks/jetSubstructureHF.cxx index 7d0917e7b51..ac7be7a4be3 100644 --- a/PWGJE/Tasks/jetSubstructureHF.cxx +++ b/PWGJE/Tasks/jetSubstructureHF.cxx @@ -200,8 +200,6 @@ struct JetSubstructureHFTask { fastjet::PseudoJet parentSubJet2; bool softDropped = false; auto nsd = 0.0; - auto zg = -1.0; - auto rg = -1.0; while (daughterSubJet.has_parents(parentSubJet1, parentSubJet2)) { bool isHFInSubjet1 = false; @@ -234,8 +232,8 @@ struct JetSubstructureHFTask { thetaVec.push_back(theta); if (z >= zCut * TMath::Power(theta / (jet.r() / 100.f), beta)) { if (!softDropped) { - zg = z; - rg = theta; + auto zg = z; + auto rg = theta; if constexpr (!isSubtracted && !isMCP) { registry.fill(HIST("h2_jet_pt_jet_zg"), jet.pt(), zg); registry.fill(HIST("h2_jet_pt_jet_rg"), jet.pt(), rg); diff --git a/PWGJE/Tasks/jetValidationQA.cxx b/PWGJE/Tasks/jetValidationQA.cxx index bf4c8348c10..4ccd260f47e 100644 --- a/PWGJE/Tasks/jetValidationQA.cxx +++ b/PWGJE/Tasks/jetValidationQA.cxx @@ -496,7 +496,7 @@ struct mcJetTrackCollisionQa { for (const auto& genJet : mcPartJets) { if (genJet.mcCollisionId() == collision.globalIndex()) { fillMcPartJets(genJet); - for (auto& mcParticle : genJet.tracks_as()) { + for (const auto& mcParticle : genJet.tracks_as()) { fillMcPartJetConstituents(mcParticle); } } @@ -506,7 +506,7 @@ struct mcJetTrackCollisionQa { for (const auto& detJet : mcDetJets) { if (detJet.collisionId() == collision.globalIndex()) { fillMcDetJets(detJet); - for (auto& detConst : detJet.tracks_as()) { + for (const auto& detConst : detJet.tracks_as()) { fillMcDetJetConstituents(detConst); } } @@ -530,7 +530,7 @@ struct mcJetTrackCollisionQa { for (const auto& genJet : mcPartJets) { if (genJet.mcCollisionId() == collision.globalIndex()) { fillMcPartJets(genJet); - for (auto& mcParticle : genJet.tracks_as()) { + for (const auto& mcParticle : genJet.tracks_as()) { fillMcPartJetConstituents(mcParticle); } } @@ -540,7 +540,7 @@ struct mcJetTrackCollisionQa { for (const auto& detJet : mcDetJets) { if (detJet.collisionId() == collision.globalIndex()) { fillMcDetJets(detJet); - for (auto& detConst : detJet.tracks_as()) { + for (const auto& detConst : detJet.tracks_as()) { fillMcDetJetConstituents(detConst); } } diff --git a/PWGJE/Tasks/mcGeneratorStudies.cxx b/PWGJE/Tasks/mcGeneratorStudies.cxx index a1bfe3a5500..8310b03c25a 100644 --- a/PWGJE/Tasks/mcGeneratorStudies.cxx +++ b/PWGJE/Tasks/mcGeneratorStudies.cxx @@ -248,7 +248,7 @@ struct MCGeneratorStudies { int iCellID = -1; try { iCellID = emcal::Geometry::GetInstance()->GetAbsCellIdFromEtaPhi(mcParticles.iteratorAt(daughterId).eta(), mcParticles.iteratorAt(daughterId).phi()); - } catch (emcal::InvalidPositionException& e) { + } catch (const emcal::InvalidPositionException& e) { iCellID = -1; } if (iCellID == -1) diff --git a/PWGJE/Tasks/nucleiInJets.cxx b/PWGJE/Tasks/nucleiInJets.cxx index 4afacbe3919..e496b59ec0b 100644 --- a/PWGJE/Tasks/nucleiInJets.cxx +++ b/PWGJE/Tasks/nucleiInJets.cxx @@ -1449,9 +1449,8 @@ struct nucleiInJets { jetHist.fill(HIST("tracksInc/deuteron/h3PtVsDeuteronNSigmaTPCVsPt"), trk.pt(), trk.tpcNSigmaDe(), centrality); } } - float massTOF = -999; if (addTOFplots && trk.hasTOF()) { - massTOF = trk.p() * std::sqrt(1.f / (trk.beta() * trk.beta()) - 1.f); + float massTOF = trk.p() * std::sqrt(1.f / (trk.beta() * trk.beta()) - 1.f); if (!useTPCpreSel) { jetHist.fill(HIST("tracksInc/proton/h2TOFmassProtonVsPt"), massTOF, trk.pt(), centrality); jetHist.fill(HIST("tracksInc/proton/h2TOFmass2ProtonVsPt"), massTOF * massTOF - MassProton * MassProton, trk.pt(), centrality); @@ -1508,9 +1507,8 @@ struct nucleiInJets { jetHist.fill(HIST("tracksInc/antiDeuteron/h3PtVsantiDeuteronNSigmaTPCVsPt"), trk.pt(), trk.tpcNSigmaDe(), centrality); } } - float massTOF = -999; if (addTOFplots && trk.hasTOF()) { - massTOF = trk.p() * std::sqrt(1.f / (trk.beta() * trk.beta()) - 1.f); + float massTOF = trk.p() * std::sqrt(1.f / (trk.beta() * trk.beta()) - 1.f); if (!useTPCpreSel) { jetHist.fill(HIST("tracksInc/antiProton/h2TOFmassantiProtonVsPt"), massTOF, trk.pt(), centrality); jetHist.fill(HIST("tracksInc/antiProton/h2TOFmass2antiProtonVsPt"), massTOF * massTOF - MassProton * MassProton, trk.pt(), centrality); @@ -1716,8 +1714,6 @@ struct nucleiInJets { return; // LOG(info) <<" size(mcd) "< leadingJetWithPtEtaPhi(3); for (const auto& mcdjet : mcdjets) { if (!mcdjet.has_matchedJetGeo()) @@ -1785,7 +1781,6 @@ struct nucleiInJets { if (std::fabs(mcTrack.y()) > cfgtrkMaxRap) continue; - bool isTpcPassed(true); bool isTof(completeTrack.hasTOF()); bool isTOFAndTPCPreSel(completeTrack.hasTOF() && (std::abs(completeTrack.tpcNSigmaPr()) < cfgnTPCPIDPrTOF || std::abs(completeTrack.tpcNSigmaDe()) < cfgnTPCPIDDeTOF || @@ -1820,6 +1815,7 @@ struct nucleiInJets { } // jet if (mapPDGToValue(mcTrack.pdgCode()) != 0) { + bool isTpcPassed(true); // why is this always true? jetHist.fill(HIST("eff/recmatched/pt/PtParticleType"), mcTrack.pt(), jetFlag, mapPDGToValue(mcTrack.pdgCode())); if (useMcC) { if (randUniform.Uniform(0, 1) < 0.5) @@ -1925,9 +1921,7 @@ struct nucleiInJets { if (std::abs(collision.posZ()) > 10) return; jetHist.fill(HIST("genmatched/vertexZ"), collision.posZ()); - std::vector mcdJetPt{}; - std::vector mcdJetPhi{}; - std::vector mcdJetEta{}; + std::vector mcpJetPt{}; std::vector mcpJetPhi{}; std::vector mcpJetEta{}; @@ -1957,6 +1951,9 @@ struct nucleiInJets { jetHist.fill(HIST("genmatched/leadingJet/hGenJetPt"), leadingMCPJet.pt()); if (leadingMCPJet.has_matchedJetGeo()) { jetHist.fill(HIST("genmatched/leadingJet/hGenJetPtMatched"), leadingMCPJet.pt()); + std::vector mcdJetPt{}; + std::vector mcdJetPhi{}; + std::vector mcdJetEta{}; for (const auto& mcdjet : leadingMCPJet.template matchedJetGeo_as()) { // Assuming matchedJetGeo_as returns valid MCD jets; no redundant has check needed // Store jet properties diff --git a/PWGJE/Tasks/phiInJets.cxx b/PWGJE/Tasks/phiInJets.cxx index 821169310b0..8571b76b880 100644 --- a/PWGJE/Tasks/phiInJets.cxx +++ b/PWGJE/Tasks/phiInJets.cxx @@ -768,10 +768,9 @@ struct phiInJets { JEhistos.fill(HIST("hMCRec_dEta_qa_rot_distribution"), dEta_rot_qa); lResonance = lDecayDaughter1 + lDecayDaughter2; - if (cfgIsKstar) - lRotationalResonance = lDecayDaughter1 + lRotationalTrack; if (cfgIsKstar) { + lRotationalResonance = lDecayDaughter1 + lRotationalTrack; JEhistos.fill(HIST("hMCRec_R_distribution"), dR); JEhistos.fill(HIST("hMCRec_hUSS_Rotational"), 1.0, lRotationalResonance.Pt(), lResonance.M()); JEhistos.fill(HIST("hMCRec_R_Rotation_distribution"), dR_rot); @@ -940,10 +939,9 @@ struct phiInJets { if (jetFlag) { if (cfgDaughterQAHists) { - if (lResonance.M() > 1.005 && lResonance.M() < 1.035) + if (lResonance.M() > 1.005 && lResonance.M() < 1.035) { RealPhiCandInJet++; - } - if (cfgDaughterQAHists) { + } JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_INSIDE_pt_v_eta"), lResonance.Pt(), lResonance.Eta()); // if (lResonance.Pt() > 2.0 && lResonance.Pt() < 3) // JEhistos.fill(HIST("hMCRec_nonmatch_hUSS_INSIDE_1D_2_3"), lResonance.M()); @@ -1235,11 +1233,11 @@ struct phiInJets { double TEMP_phi_dgth_pz[2] = {0}; bool good_daughter[2] = {false}; - int dgth_index = 0; // First we check for Forced BR // if we check for Phi if (!cfgIsKstar) { if (mcParticle.has_daughters()) { + int dgth_index = 0; for (auto& dgth : mcParticle.daughters_as()) { if (std::fabs(dgth.pdgCode()) != 321) { skip = true; From 187172876fe87db40c2e92aff42be2e3d3db0fc2 Mon Sep 17 00:00:00 2001 From: nzardosh Date: Tue, 11 Nov 2025 12:25:48 +0100 Subject: [PATCH 09/13] fixing angles in perpcone --- PWGJE/Core/JetUtilities.h | 31 ++++++++++++++++--------------- 1 file changed, 16 insertions(+), 15 deletions(-) diff --git a/PWGJE/Core/JetUtilities.h b/PWGJE/Core/JetUtilities.h index e0e7090cc9d..9c0df210987 100644 --- a/PWGJE/Core/JetUtilities.h +++ b/PWGJE/Core/JetUtilities.h @@ -172,33 +172,34 @@ std::tuple estimateRhoPerpCone(const T& inputParticles, const U& double perpMdDensity1 = 0; double perpMdDensity2 = 0; - double PerpendicularConeAxisPhi1 = 999., PerpendicularConeAxisPhi2 = 999.; + const double jetPhi = RecoDecay::constrainAngle(jet.phi(), -M_PI); + const double jetEta = jet.eta(); + const double radius = static_cast(perpConeR); // build 2 perp cones in phi around the leading jet (right and left of the jet) - PerpendicularConeAxisPhi1 = RecoDecay::constrainAngle(jet.phi() + (M_PI / 2.)); // This will contrain the angel between 0-2Pi - PerpendicularConeAxisPhi2 = RecoDecay::constrainAngle(jet.phi() - (M_PI / 2.)); // This will contrain the angel between 0-2Pi + double PerpendicularConeAxisPhi1 = RecoDecay::constrainAngle(jetPhi + (M_PI / 2.), -M_PI); // This will contrain the angel between -pi & Pi + double PerpendicularConeAxisPhi2 = RecoDecay::constrainAngle(jetPhi - (M_PI / 2.), -M_PI); // This will contrain the angel between -pi & Pi - for (auto& particle : inputParticles) { + for (const auto& particle : inputParticles) { // sum the momentum of all paricles that fill the two cones - double dPhi1 = particle.phi() - PerpendicularConeAxisPhi1; - dPhi1 = RecoDecay::constrainAngle(dPhi1, -M_PI); // This will contrain the angel between -pi & Pi - double dPhi2 = particle.phi() - PerpendicularConeAxisPhi2; - dPhi2 = RecoDecay::constrainAngle(dPhi2, -M_PI); // This will contrain the angel between -pi & Pi - double dEta = jet.eta() - particle.eta(); // The perp cone eta is the same as the leading jet since the cones are perpendicular only in phi - if (TMath::Sqrt(dPhi1 * dPhi1 + dEta * dEta) <= static_cast(perpConeR)) { - perpPtDensity1 += particle.perp(); + const double phi = RecoDecay::constrainAngle(particle.phi(), -M_PI); + double dPhi1 = RecoDecay::constrainAngle(phi - PerpendicularConeAxisPhi1, -M_PI); // This will contrain the angel between -pi & Pi + double dPhi2 = RecoDecay::constrainAngle(phi - PerpendicularConeAxisPhi2, -M_PI); // This will contrain the angel between -pi & Pi + double dEta = jetEta - particle.eta(); // The perp cone eta is the same as the leading jet since the cones are perpendicular only in phi + if (TMath::Sqrt(dPhi1 * dPhi1 + dEta * dEta) <= static_cast(radius)) { + perpPtDensity1 += particle.pt(); perpMdDensity1 += TMath::Sqrt(particle.m() * particle.m() + particle.pt() * particle.pt()) - particle.pt(); } - if (TMath::Sqrt(dPhi2 * dPhi2 + dEta * dEta) <= perpConeR) { - perpPtDensity2 += particle.perp(); + if (TMath::Sqrt(dPhi2 * dPhi2 + dEta * dEta) <= static_cast(radius)) { + perpPtDensity2 += particle.pt(); perpMdDensity2 += TMath::Sqrt(particle.m() * particle.m() + particle.pt() * particle.pt()) - particle.pt(); } } // Caculate rho as the ratio of average pT of the two cones / the cone area - double perpPtDensity = (perpPtDensity1 + perpPtDensity2) / (2 * M_PI * static_cast(perpConeR) * static_cast(perpConeR)); - double perpMdDensity = (perpMdDensity1 + perpMdDensity2) / (2 * M_PI * static_cast(perpConeR) * static_cast(perpConeR)); + double perpPtDensity = (perpPtDensity1 + perpPtDensity2) / (2 * M_PI * static_cast(radius) * static_cast(radius)); + double perpMdDensity = (perpMdDensity1 + perpMdDensity2) / (2 * M_PI * static_cast(radius) * static_cast(radius)); return std::make_tuple(perpPtDensity, perpMdDensity); } From 95dc3e9379058e75286b2dc32f405308c28b2ff6 Mon Sep 17 00:00:00 2001 From: nzardosh Date: Tue, 11 Nov 2025 12:30:20 +0100 Subject: [PATCH 10/13] removing pvector definition from derived data model --- PWGJE/DataModel/JetReducedData.h | 2 -- 1 file changed, 2 deletions(-) diff --git a/PWGJE/DataModel/JetReducedData.h b/PWGJE/DataModel/JetReducedData.h index 900abf07396..b18565c4254 100644 --- a/PWGJE/DataModel/JetReducedData.h +++ b/PWGJE/DataModel/JetReducedData.h @@ -340,8 +340,6 @@ DECLARE_SOA_COLUMN(Flags, flags, uint8_t); DECLARE_SOA_SELF_ARRAY_INDEX_COLUMN(Mothers, mothers); DECLARE_SOA_SELF_SLICE_INDEX_COLUMN(Daughters, daughters); -DECLARE_SOA_DYNAMIC_COLUMN(PVector, pVector, //! Momentum vector in x,y,z-directions in GeV/c - [](float px, float py, float pz) -> std::array { return std::array{px, py, pz}; }); DECLARE_SOA_DYNAMIC_COLUMN(Px, px, [](float pt, float phi) -> float { return pt * std::cos(phi); }); DECLARE_SOA_DYNAMIC_COLUMN(Py, py, From f819f96f6c38c180a98e46412583507057d728ed Mon Sep 17 00:00:00 2001 From: nzardosh Date: Tue, 11 Nov 2025 12:45:41 +0100 Subject: [PATCH 11/13] fixing status anf flags in mcparticles --- PWGJE/Core/JetDQUtilities.h | 2 +- PWGJE/DataModel/JetReducedDataDQ.h | 34 ++++++++++++++++----- PWGJE/DataModel/JetReducedDataV0.h | 34 ++++++++++++++++----- PWGJE/TableProducer/derivedDataProducer.cxx | 6 ++-- 4 files changed, 58 insertions(+), 18 deletions(-) diff --git a/PWGJE/Core/JetDQUtilities.h b/PWGJE/Core/JetDQUtilities.h index 75ad0a399b9..e51e8232c26 100644 --- a/PWGJE/Core/JetDQUtilities.h +++ b/PWGJE/Core/JetDQUtilities.h @@ -352,7 +352,7 @@ void fillDielectronCandidateTable(T const& candidate, int32_t collisionIndex, U& template void fillDielectronCandidateMcTable(T const& candidate, int32_t mcCollisionIndex, U& DielectronMcTable) { - DielectronMcTable(mcCollisionIndex, candidate.pt(), candidate.eta(), candidate.phi(), candidate.y(), candidate.e(), candidate.m(), candidate.pdgCode(), candidate.getGenStatusCode(), candidate.getHepMCStatusCode(), candidate.isPhysicalPrimary(), candidate.decayFlag(), candidate.origin()); + DielectronMcTable(mcCollisionIndex, candidate.pt(), candidate.eta(), candidate.phi(), candidate.y(), candidate.e(), candidate.m(), candidate.pdgCode(), candidate.statusCode(), candidate.flags(), candidate.decayFlag(), candidate.origin()); } }; // namespace jetdqutilities diff --git a/PWGJE/DataModel/JetReducedDataDQ.h b/PWGJE/DataModel/JetReducedDataDQ.h index 5bc14c333f3..f9395b76fba 100644 --- a/PWGJE/DataModel/JetReducedDataDQ.h +++ b/PWGJE/DataModel/JetReducedDataDQ.h @@ -72,9 +72,8 @@ DECLARE_SOA_COLUMN(Y, y, float); DECLARE_SOA_COLUMN(E, e, float); DECLARE_SOA_COLUMN(M, m, float); DECLARE_SOA_COLUMN(PdgCode, pdgCode, int); -DECLARE_SOA_COLUMN(GenStatusCode, getGenStatusCode, int); // TODO : We can look at combining this with the two below -DECLARE_SOA_COLUMN(HepMCStatusCode, getHepMCStatusCode, int); -DECLARE_SOA_COLUMN(IsPhysicalPrimary, isPhysicalPrimary, bool); +DECLARE_SOA_COLUMN(StatusCode, statusCode, int); +DECLARE_SOA_COLUMN(Flags, flags, uint8_t); DECLARE_SOA_SELF_ARRAY_INDEX_COLUMN(Mothers, mothers); DECLARE_SOA_SELF_SLICE_INDEX_COLUMN(Daughters, daughters); DECLARE_SOA_COLUMN(DecayFlag, decayFlag, int8_t); @@ -87,6 +86,21 @@ DECLARE_SOA_DYNAMIC_COLUMN(Pz, pz, [](float pt, float eta) -> float { return pt * std::sinh(eta); }); DECLARE_SOA_DYNAMIC_COLUMN(P, p, [](float pt, float eta) -> float { return pt * std::cosh(eta); }); +DECLARE_SOA_DYNAMIC_COLUMN(Energy, energy, + [](float e) -> float { return e; }); + +DECLARE_SOA_DYNAMIC_COLUMN(ProducedByGenerator, producedByGenerator, //! True if particle produced by the generator (==TMCProcess::kPrimary); False if by the transport code + [](uint8_t flags) -> bool { return (flags & o2::aod::mcparticle::enums::ProducedByTransport) == 0x0; }); +DECLARE_SOA_DYNAMIC_COLUMN(FromBackgroundEvent, fromBackgroundEvent, //! Particle from background event + [](uint8_t flags) -> bool { return (flags & o2::aod::mcparticle::enums::FromBackgroundEvent) == o2::aod::mcparticle::enums::FromBackgroundEvent; }); +DECLARE_SOA_DYNAMIC_COLUMN(GetProcess, getProcess, //! The VMC physics code (as int) that generated this particle (see header TMCProcess.h in ROOT) + [](uint8_t flags, int statusCode) -> int { if ((flags & o2::aod::mcparticle::enums::ProducedByTransport) == 0x0) { return 0 /*TMCProcess::kPrimary*/; } else { return statusCode; } }); +DECLARE_SOA_DYNAMIC_COLUMN(GetGenStatusCode, getGenStatusCode, //! The native status code put by the generator, or -1 if a particle produced during transport + [](uint8_t flags, int statusCode) -> int { if ((flags & o2::aod::mcparticle::enums::ProducedByTransport) == 0x0) { return o2::mcgenstatus::getGenStatusCode(statusCode); } else { return -1; } }); +DECLARE_SOA_DYNAMIC_COLUMN(GetHepMCStatusCode, getHepMCStatusCode, //! The HepMC status code put by the generator, or -1 if a particle produced during transport + [](uint8_t flags, int statusCode) -> int { if ((flags & o2::aod::mcparticle::enums::ProducedByTransport) == 0x0) { return o2::mcgenstatus::getHepMCStatusCode(statusCode); } else { return -1; } }); +DECLARE_SOA_DYNAMIC_COLUMN(IsPhysicalPrimary, isPhysicalPrimary, //! True if particle is considered a physical primary according to the ALICE definition + [](uint8_t flags) -> bool { return (flags & o2::aod::mcparticle::enums::PhysicalPrimary) == o2::aod::mcparticle::enums::PhysicalPrimary; }); } // namespace jdielectronmc DECLARE_SOA_TABLE_STAGED(JDielectronMcs, "JDIELMC", @@ -99,15 +113,21 @@ DECLARE_SOA_TABLE_STAGED(JDielectronMcs, "JDIELMC", jdielectronmc::E, jdielectronmc::M, jdielectronmc::PdgCode, - jdielectronmc::GenStatusCode, - jdielectronmc::HepMCStatusCode, - jdielectronmc::IsPhysicalPrimary, + jdielectronmc::StatusCode, + jdielectronmc::Flags, jdielectronmc::DecayFlag, jdielectronmc::Origin, jdielectronmc::Px, jdielectronmc::Py, jdielectronmc::Pz, - jdielectronmc::P); + jdielectronmc::P, + jdielectronmc::Energy, + jdielectronmc::ProducedByGenerator, + jdielectronmc::FromBackgroundEvent, + jdielectronmc::GetProcess, + jdielectronmc::GetGenStatusCode, + jdielectronmc::GetHepMCStatusCode, + jdielectronmc::IsPhysicalPrimary); using JDielectronMc = JDielectronMcs::iterator; using StoredJDielectronMc = StoredJDielectronMcs::iterator; diff --git a/PWGJE/DataModel/JetReducedDataV0.h b/PWGJE/DataModel/JetReducedDataV0.h index 433cd885f65..ec460d5b186 100644 --- a/PWGJE/DataModel/JetReducedDataV0.h +++ b/PWGJE/DataModel/JetReducedDataV0.h @@ -65,9 +65,8 @@ DECLARE_SOA_COLUMN(Y, y, float); DECLARE_SOA_COLUMN(E, e, float); DECLARE_SOA_COLUMN(M, m, float); DECLARE_SOA_COLUMN(PdgCode, pdgCode, int); -DECLARE_SOA_COLUMN(GenStatusCode, getGenStatusCode, int); // TODO : We can look at combining this with the two below -DECLARE_SOA_COLUMN(HepMCStatusCode, getHepMCStatusCode, int); -DECLARE_SOA_COLUMN(IsPhysicalPrimary, isPhysicalPrimary, bool); +DECLARE_SOA_COLUMN(StatusCode, statusCode, int); +DECLARE_SOA_COLUMN(Flags, flags, uint8_t); DECLARE_SOA_SELF_ARRAY_INDEX_COLUMN(Mothers, mothers); DECLARE_SOA_SELF_SLICE_INDEX_COLUMN(Daughters, daughters); DECLARE_SOA_COLUMN(DecayFlag, decayFlag, int8_t); @@ -79,6 +78,21 @@ DECLARE_SOA_DYNAMIC_COLUMN(Pz, pz, [](float pt, float eta) -> float { return pt * std::sinh(eta); }); DECLARE_SOA_DYNAMIC_COLUMN(P, p, [](float pt, float eta) -> float { return pt * std::cosh(eta); }); +DECLARE_SOA_DYNAMIC_COLUMN(Energy, energy, + [](float e) -> float { return e; }); + +DECLARE_SOA_DYNAMIC_COLUMN(ProducedByGenerator, producedByGenerator, //! True if particle produced by the generator (==TMCProcess::kPrimary); False if by the transport code + [](uint8_t flags) -> bool { return (flags & o2::aod::mcparticle::enums::ProducedByTransport) == 0x0; }); +DECLARE_SOA_DYNAMIC_COLUMN(FromBackgroundEvent, fromBackgroundEvent, //! Particle from background event + [](uint8_t flags) -> bool { return (flags & o2::aod::mcparticle::enums::FromBackgroundEvent) == o2::aod::mcparticle::enums::FromBackgroundEvent; }); +DECLARE_SOA_DYNAMIC_COLUMN(GetProcess, getProcess, //! The VMC physics code (as int) that generated this particle (see header TMCProcess.h in ROOT) + [](uint8_t flags, int statusCode) -> int { if ((flags & o2::aod::mcparticle::enums::ProducedByTransport) == 0x0) { return 0 /*TMCProcess::kPrimary*/; } else { return statusCode; } }); +DECLARE_SOA_DYNAMIC_COLUMN(GetGenStatusCode, getGenStatusCode, //! The native status code put by the generator, or -1 if a particle produced during transport + [](uint8_t flags, int statusCode) -> int { if ((flags & o2::aod::mcparticle::enums::ProducedByTransport) == 0x0) { return o2::mcgenstatus::getGenStatusCode(statusCode); } else { return -1; } }); +DECLARE_SOA_DYNAMIC_COLUMN(GetHepMCStatusCode, getHepMCStatusCode, //! The HepMC status code put by the generator, or -1 if a particle produced during transport + [](uint8_t flags, int statusCode) -> int { if ((flags & o2::aod::mcparticle::enums::ProducedByTransport) == 0x0) { return o2::mcgenstatus::getHepMCStatusCode(statusCode); } else { return -1; } }); +DECLARE_SOA_DYNAMIC_COLUMN(IsPhysicalPrimary, isPhysicalPrimary, //! True if particle is considered a physical primary according to the ALICE definition + [](uint8_t flags) -> bool { return (flags & o2::aod::mcparticle::enums::PhysicalPrimary) == o2::aod::mcparticle::enums::PhysicalPrimary; }); } // namespace jv0mc DECLARE_SOA_TABLE_STAGED(JV0Mcs, "JV0MC", @@ -91,14 +105,20 @@ DECLARE_SOA_TABLE_STAGED(JV0Mcs, "JV0MC", jv0mc::E, jv0mc::M, jv0mc::PdgCode, - jv0mc::GenStatusCode, - jv0mc::HepMCStatusCode, - jv0mc::IsPhysicalPrimary, + jv0mc::StatusCode, + jv0mc::Flags, jv0mc::DecayFlag, jv0mc::Px, jv0mc::Py, jv0mc::Pz, - jv0mc::P); + jv0mc::P, + jv0mc::Energy, + jv0mc::ProducedByGenerator, + jv0mc::FromBackgroundEvent, + jv0mc::GetProcess, + jv0mc::GetGenStatusCode, + jv0mc::GetHepMCStatusCode, + jv0mc::IsPhysicalPrimary); using JV0Mc = JV0Mcs::iterator; using StoredJV0Mc = StoredJV0Mcs::iterator; diff --git a/PWGJE/TableProducer/derivedDataProducer.cxx b/PWGJE/TableProducer/derivedDataProducer.cxx index 61c9240075d..64e0569f8f4 100644 --- a/PWGJE/TableProducer/derivedDataProducer.cxx +++ b/PWGJE/TableProducer/derivedDataProducer.cxx @@ -429,7 +429,7 @@ struct JetDerivedDataProducerTask { i++; } } - products.jMcParticlesTable(particle.mcCollisionId(), particle.pt(), particle.eta(), particle.phi(), particle.y(), particle.e(), particle.pdgCode(), particle.flags(), particle.statusCode(), mothersId, daughtersId); + products.jMcParticlesTable(particle.mcCollisionId(), particle.pt(), particle.eta(), particle.phi(), particle.y(), particle.e(), particle.pdgCode(), particle.statusCode(), particle.flags(), mothersId, daughtersId); products.jParticlesParentIndexTable(particle.globalIndex()); } PROCESS_SWITCH(JetDerivedDataProducerTask, processParticles, "produces derived parrticle table", false); @@ -788,7 +788,7 @@ struct JetDerivedDataProducerTask { } } auto pdgParticle = pdgDatabase->GetParticle(particle.pdgCode()); - products.jV0McsTable(products.jV0McCollisionsTable.lastIndex(), particle.pt(), particle.eta(), particle.phi(), particle.y(), particle.e(), pdgParticle->Mass(), particle.pdgCode(), particle.getGenStatusCode(), particle.getHepMCStatusCode(), particle.isPhysicalPrimary(), jetv0utilities::setV0ParticleDecayBit(particles, particle)); + products.jV0McsTable(products.jV0McCollisionsTable.lastIndex(), particle.pt(), particle.eta(), particle.phi(), particle.y(), particle.e(), pdgParticle->Mass(), particle.pdgCode(), particle.statusCode(), particle.flags(), jetv0utilities::setV0ParticleDecayBit(particles, particle)); products.jV0McIdsTable(mcCollision.globalIndex(), particle.globalIndex(), mothersId, daughtersId); } } @@ -842,7 +842,7 @@ struct JetDerivedDataProducerTask { } } auto pdgParticle = pdgDatabase->GetParticle(particle.pdgCode()); - products.jDielectronMcsTable(products.jDielectronMcCollisionsTable.lastIndex(), particle.pt(), particle.eta(), particle.phi(), particle.y(), particle.e(), pdgParticle->Mass(), particle.pdgCode(), particle.getGenStatusCode(), particle.getHepMCStatusCode(), particle.isPhysicalPrimary(), jetdqutilities::setDielectronParticleDecayBit(particles, particle), RecoDecay::getCharmHadronOrigin(particles, particle, false)); // Todo: should the last thing be false? + products.jDielectronMcsTable(products.jDielectronMcCollisionsTable.lastIndex(), particle.pt(), particle.eta(), particle.phi(), particle.y(), particle.e(), pdgParticle->Mass(), particle.pdgCode(), particle.statusCode(), particle.flags(), jetdqutilities::setDielectronParticleDecayBit(particles, particle), RecoDecay::getCharmHadronOrigin(particles, particle, false)); // Todo: should the last thing be false? products.jDielectronMcIdsTable(mcCollision.globalIndex(), particle.globalIndex(), mothersId, daughtersId); products.JDielectronMcRCollDummysTable(false); } From 017d284cfe640393e0a105cf3d00b175bb4254da Mon Sep 17 00:00:00 2001 From: nzardosh Date: Tue, 11 Nov 2025 14:32:06 +0100 Subject: [PATCH 12/13] fixing error which should be a warning in trigger correlations --- PWGJE/Tasks/triggerCorrelations.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PWGJE/Tasks/triggerCorrelations.cxx b/PWGJE/Tasks/triggerCorrelations.cxx index 73acb90145e..c0c272d39c7 100644 --- a/PWGJE/Tasks/triggerCorrelations.cxx +++ b/PWGJE/Tasks/triggerCorrelations.cxx @@ -59,7 +59,7 @@ struct TriggerCorrelationsTask { { for (std::vector::size_type iTrig = 0; iTrig < triggerMaskBits.size(); iTrig++) { if (fill) { - if (iTrig >= 0 && iTrig < nChargedTriggers && jetderiveddatautilities::selectChargedTrigger(collision, iTrig + 1)) { + if (iTrig < nChargedTriggers && jetderiveddatautilities::selectChargedTrigger(collision, iTrig + 1)) { registry.fill(HIST("triggerCorrelations"), iCurrentTrig, iTrig); } if (iTrig >= nChargedTriggers && iTrig < (nChargedTriggers + nChargedHFTriggers) && jetderiveddatautilities::selectChargedHFTrigger(collision, iTrig - nChargedTriggers + 1)) { From 79d298f983e22cadfeef92bf00dba216350abd6c Mon Sep 17 00:00:00 2001 From: nzardosh Date: Tue, 11 Nov 2025 14:32:31 +0100 Subject: [PATCH 13/13] adding the ability to skip parton level info saving --- PWGJE/JetFinders/jetFinder.cxx | 8 ++-- PWGJE/JetFinders/jetFinderHF.cxx | 6 +-- PWGJE/TableProducer/derivedDataWriter.cxx | 45 ++++++++++++------- .../eventwiseConstituentSubtractor.cxx | 2 +- PWGJE/TableProducer/rhoEstimator.cxx | 2 +- 5 files changed, 38 insertions(+), 25 deletions(-) diff --git a/PWGJE/JetFinders/jetFinder.cxx b/PWGJE/JetFinders/jetFinder.cxx index ff7e67966a5..5255976feb3 100644 --- a/PWGJE/JetFinders/jetFinder.cxx +++ b/PWGJE/JetFinders/jetFinder.cxx @@ -238,7 +238,7 @@ struct JetFinderTask { { // TODO: MC event selection? inputParticles.clear(); - jetfindingutilities::analyseParticles, soa::Filtered::iterator>(inputParticles, particleSelection, 1, particles, pdgDatabase); + jetfindingutilities::analyseParticles, soa::Filtered::iterator>(inputParticles, particleSelection, 1, particles, pdgDatabase); jetfindingutilities::findJets(jetFinder, inputParticles, jetPtMin, jetPtMax, jetRadius, jetAreaFractionMin, collision, jetsTable, constituentsTable, fillTHnSparse ? registry.get(HIST("hJetMCP")) : std::shared_ptr(nullptr), fillTHnSparse); } PROCESS_SWITCH(JetFinderTask, processParticleLevelChargedJets, "Particle level charged jet finding", false); @@ -247,7 +247,7 @@ struct JetFinderTask { { // TODO: MC event selection? inputParticles.clear(); - jetfindingutilities::analyseParticles, soa::Filtered::iterator>(inputParticles, particleSelection, 1, particles, pdgDatabase); + jetfindingutilities::analyseParticles, soa::Filtered::iterator>(inputParticles, particleSelection, 1, particles, pdgDatabase); jetfindingutilities::findJets(jetFinder, inputParticles, jetEWSPtMin, jetEWSPtMax, jetRadius, jetAreaFractionMin, collision, jetsEvtWiseSubTable, constituentsEvtWiseSubTable, fillTHnSparse ? registry.get(HIST("hJetEWSMCP")) : std::shared_ptr(nullptr), fillTHnSparse); } PROCESS_SWITCH(JetFinderTask, processParticleLevelChargedEvtWiseSubJets, "Particle level charged with event-wise constituent subtraction jet finding", false); @@ -256,7 +256,7 @@ struct JetFinderTask { { // TODO: MC event selection? inputParticles.clear(); - jetfindingutilities::analyseParticles, soa::Filtered::iterator>(inputParticles, particleSelection, 2, particles, pdgDatabase); + jetfindingutilities::analyseParticles, soa::Filtered::iterator>(inputParticles, particleSelection, 2, particles, pdgDatabase); jetfindingutilities::findJets(jetFinder, inputParticles, jetPtMin, jetPtMax, jetRadius, jetAreaFractionMin, collision, jetsTable, constituentsTable, fillTHnSparse ? registry.get(HIST("hJetMCP")) : std::shared_ptr(nullptr), fillTHnSparse); } PROCESS_SWITCH(JetFinderTask, processParticleLevelNeutralJets, "Particle level neutral jet finding", false); @@ -265,7 +265,7 @@ struct JetFinderTask { { // TODO: MC event selection? inputParticles.clear(); - jetfindingutilities::analyseParticles, soa::Filtered::iterator>(inputParticles, particleSelection, 0, particles, pdgDatabase); + jetfindingutilities::analyseParticles, soa::Filtered::iterator>(inputParticles, particleSelection, 0, particles, pdgDatabase); jetfindingutilities::findJets(jetFinder, inputParticles, jetPtMin, jetPtMax, jetRadius, jetAreaFractionMin, collision, jetsTable, constituentsTable, fillTHnSparse ? registry.get(HIST("hJetMCP")) : std::shared_ptr(nullptr), fillTHnSparse); } diff --git a/PWGJE/JetFinders/jetFinderHF.cxx b/PWGJE/JetFinders/jetFinderHF.cxx index c80688aff45..9e8ad7eb1e9 100644 --- a/PWGJE/JetFinders/jetFinderHF.cxx +++ b/PWGJE/JetFinders/jetFinderHF.cxx @@ -245,10 +245,10 @@ struct JetFinderHFTask { if (!jetfindingutilities::analyseCandidate(inputParticles, candidate, candPtMin, candPtMax, candYMin, candYMax)) { return; } - if constexpr (!isEvtWiseSub) { - jetfindingutilities::analyseParticles(inputParticles, particleSelection, jetTypeParticleLevel, particles, pdgDatabase, &candidate); - } else { + if constexpr (isEvtWiseSub) { jetfindingutilities::analyseParticles(inputParticles, particleSelection, jetTypeParticleLevel, particles, pdgDatabase, &candidate); + } else { + jetfindingutilities::analyseParticles(inputParticles, particleSelection, jetTypeParticleLevel, particles, pdgDatabase, &candidate); } jetfindingutilities::findJets(jetFinder, inputParticles, minJetPt, maxJetPt, jetRadius, jetAreaFractionMin, collision, jetsTableInput, constituentsTableInput, registry.get(HIST("hJetMCP")), fillTHnSparse, true); } diff --git a/PWGJE/TableProducer/derivedDataWriter.cxx b/PWGJE/TableProducer/derivedDataWriter.cxx index 1c3125a68ff..5f6ad60bf59 100644 --- a/PWGJE/TableProducer/derivedDataWriter.cxx +++ b/PWGJE/TableProducer/derivedDataWriter.cxx @@ -53,6 +53,7 @@ struct JetDerivedDataWriter { Configurable performTrackSelection{"performTrackSelection", true, "only save tracks that pass one of the track selections"}; Configurable trackPtSelectionMin{"trackPtSelectionMin", 0.15, "only save tracks that have a pT larger than this pT"}; Configurable trackEtaSelectionMax{"trackEtaSelectionMax", 0.9, "only save tracks that have an eta smaller than this eta"}; + Configurable savePartonLevelInfo{"savePartonLevelInfo", true, "save parton level info at MCP level"}; 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"}; } config; @@ -665,21 +666,27 @@ struct JetDerivedDataWriter { for (auto particle : particlesPerMcCollision) { std::vector mothersIds; - if (particle.has_mothers()) { - auto mothersIdTemps = particle.mothersIds(); - for (auto mothersIdTemp : mothersIdTemps) { - mothersIds.push_back(particleMapping[mothersIdTemp]); - } - } int daughtersIds[2] = {-1, -1}; - if (particle.has_daughters()) { - auto i = 0; - for (auto daughterId : particle.daughtersIds()) { - if (i > 1) { - break; + if (config.savePartonLevelInfo) { + if (particle.has_mothers()) { + auto mothersIdTemps = particle.mothersIds(); + for (auto mothersIdTemp : mothersIdTemps) { + mothersIds.push_back(particleMapping[mothersIdTemp]); + } + } + if (particle.has_daughters()) { + auto i = 0; + for (auto daughterId : particle.daughtersIds()) { + if (i > 1) { + break; + } + daughtersIds[i] = particleMapping[daughterId]; + i++; } - daughtersIds[i] = particleMapping[daughterId]; - i++; + } + } else { + if (!particle.isPhysicalPrimary()) { // add outgoing partons exclusion here later + continue; } } products.storedJMcParticlesTable(mcCollisionMapping[mcCollision.globalIndex()], o2::math_utils::detail::truncateFloatFraction(particle.pt(), precisionMomentumMask), o2::math_utils::detail::truncateFloatFraction(particle.eta(), precisionPositionMask), o2::math_utils::detail::truncateFloatFraction(particle.phi(), precisionPositionMask), o2::math_utils::detail::truncateFloatFraction(particle.y(), precisionPositionMask), o2::math_utils::detail::truncateFloatFraction(particle.e(), precisionMomentumMask), particle.pdgCode(), particle.statusCode(), particle.flags(), mothersIds, daughtersIds); @@ -909,14 +916,14 @@ struct JetDerivedDataWriter { } PROCESS_SWITCH(JetDerivedDataWriter, processColllisonsMcCollisionLabel, "write out collision mcCollision label output tables", false); - void processTracksMcParticleLabel(soa::Join::iterator const& collision, soa::Join const& tracks) + void processTracksMcParticleLabel(soa::Join::iterator const& collision, soa::Join const& tracks, aod::JMcParticles const&) { if (collision.isCollisionSelected()) { for (const auto& track : tracks) { if (!trackSelection(track)) { continue; } - if (track.has_mcParticle()) { + if (track.has_mcParticle() && (config.savePartonLevelInfo || track.mcParticle().isPhysicalPrimary())) { products.storedJMcTracksLabelTable(particleMapping[track.mcParticleId()]); } else { products.storedJMcTracksLabelTable(-1); @@ -926,12 +933,18 @@ struct JetDerivedDataWriter { } PROCESS_SWITCH(JetDerivedDataWriter, processTracksMcParticleLabel, "write out track mcParticle label output tables", false); - void processClusterMcLabel(soa::Join::iterator const& collision, soa::Join const& clusters) + void processClusterMcLabel(soa::Join::iterator const& collision, soa::Join const& clusters, aod::JMcParticles const& particles) { if (collision.isCollisionSelected()) { for (const auto& cluster : clusters) { std::vector clusterStoredJParticleIDs; for (const auto& clusterParticleId : cluster.mcParticlesIds()) { + if (!config.savePartonLevelInfo) { + const auto& particle = particles.iteratorAt(clusterParticleId); + if (!particle.isPhysicalPrimary()) { + continue; + } + } clusterStoredJParticleIDs.push_back(particleMapping[clusterParticleId]); } std::vector amplitudeA; diff --git a/PWGJE/TableProducer/eventwiseConstituentSubtractor.cxx b/PWGJE/TableProducer/eventwiseConstituentSubtractor.cxx index 7dbfa127ee5..53b6818cfcc 100644 --- a/PWGJE/TableProducer/eventwiseConstituentSubtractor.cxx +++ b/PWGJE/TableProducer/eventwiseConstituentSubtractor.cxx @@ -166,7 +166,7 @@ struct eventWiseConstituentSubtractorTask { } inputParticles.clear(); tracksSubtracted.clear(); - jetfindingutilities::analyseParticles, soa::Filtered::iterator>(inputParticles, particleSelection, 1, particles, pdgDatabase); + jetfindingutilities::analyseParticles, soa::Filtered::iterator>(inputParticles, particleSelection, 1, particles, pdgDatabase); tracksSubtracted = eventWiseConstituentSubtractor.JetBkgSubUtils::doEventConstSub(inputParticles, mcCollision.rho(), mcCollision.rhoM()); diff --git a/PWGJE/TableProducer/rhoEstimator.cxx b/PWGJE/TableProducer/rhoEstimator.cxx index 14ec846b74f..2ef37a290ed 100644 --- a/PWGJE/TableProducer/rhoEstimator.cxx +++ b/PWGJE/TableProducer/rhoEstimator.cxx @@ -233,7 +233,7 @@ struct RhoEstimatorTask { return; } inputParticles.clear(); - jetfindingutilities::analyseParticles, soa::Filtered::iterator>(inputParticles, particleSelection, 1, particles, pdgDatabase); + jetfindingutilities::analyseParticles, soa::Filtered::iterator>(inputParticles, particleSelection, 1, particles, pdgDatabase); auto [rho, rhoM] = bkgSub.estimateRhoAreaMedian(inputParticles, config.doSparse); rhoChargedMcTable(rho, rhoM); }