From cc5bcbe43a953b99ffb7d832851c20bd3a204b5e Mon Sep 17 00:00:00 2001 From: sandor-lokos Date: Mon, 9 Feb 2026 13:15:21 +0100 Subject: [PATCH 01/11] Initial commit for FT0 and FV0 amplitude bits per channel propagation throught SGCandProducer.cxx --- PWGUD/Core/UDHelpers.h | 113 +++++++++++++++++++++++++ PWGUD/DataModel/UDTables.h | 39 +++++++++ PWGUD/TableProducer/SGCandProducer.cxx | 75 ++++++++++++++++ PWGUD/Tasks/CMakeLists.txt | 5 ++ PWGUD/Tasks/upcTestFITBitMapping.cxx | 108 +++++++++++++++++++++++ 5 files changed, 340 insertions(+) create mode 100644 PWGUD/Tasks/upcTestFITBitMapping.cxx diff --git a/PWGUD/Core/UDHelpers.h b/PWGUD/Core/UDHelpers.h index a4137784ac5..824857a6343 100644 --- a/PWGUD/Core/UDHelpers.h +++ b/PWGUD/Core/UDHelpers.h @@ -27,6 +27,7 @@ #include "DataFormatsFIT/Triggers.h" #include "DataFormatsFT0/Digit.h" #include "Framework/Logger.h" +#include "Common/Core/RecoDecay.h" #include "TLorentzVector.h" @@ -534,6 +535,118 @@ bool FITveto(T const& bc, DGCutparHolder const& diffCuts) return false; } +// ----------------------------------------------------------------------------- +// return eta and phi of a given FIT channel based on the bitset +// Bit layout contract: +constexpr int kFT0Bits = 208; // FT0 total channels +constexpr int kFV0Bits = 48; // FV0A channels +constexpr int kTotalBits = 256; // 4*64 + +// FT0 side split (adjust if your channelization differs in your tag) +constexpr int kFT0AChannels = 96; // FT0A channels are [0..95] +constexpr int kFT0CChannels = 112; // FT0C channels are [96..207] +static_assert(kFT0AChannels + kFT0CChannels == kFT0Bits); + +using Bits256 = std::array; + +inline Bits256 makeBits256(uint64_t w0, uint64_t w1, uint64_t w2, uint64_t w3) +{ + return {w0, w1, w2, w3}; +} + +inline bool testBit(Bits256 const& w, int bit) +{ + if (bit < 0 || bit >= kTotalBits) { + return false; + } + return (w[bit >> 6] >> (bit & 63)) & 1ULL; +} + +struct FitBitRef { + enum class Det : uint8_t { FT0, FV0, Unknown }; + Det det = Det::Unknown; + int ch = -1; // FT0: 0..207, FV0: 0..47 + bool isC = false; // only meaningful for FT0 +}; + +inline FitBitRef decodeFitBit(int bit) +{ + FitBitRef out; + if (bit >= 0 && bit < kFT0Bits) { + out.det = FitBitRef::Det::FT0; + out.ch = bit; // FT0 channel id + out.isC = (bit >= kFT0AChannels); // C side if in upper range + return out; + } + if (bit >= kFT0Bits && bit < kTotalBits) { + out.det = FitBitRef::Det::FV0; + out.ch = bit - kFT0Bits; // FV0A channel id 0..47 + return out; + } + return out; +} + +// You can pass whatever type offsetFT0 is (vector/array/span of objects with getX/Y/Z). +template +inline double getPhiFT0_fromChannel(FT0DetT& ft0Det, int ft0Ch, OffsetsT const& offsetFT0, int i) +{ + ft0Det.calculateChannelCenter(); + auto chPos = ft0Det.getChannelCenter(ft0Ch); + + const double x = chPos.X() + offsetFT0[i].getX(); + const double y = chPos.Y() + offsetFT0[i].getY(); + + return RecoDecay::phi(x, y); +} + +template +inline double getEtaFT0_fromChannel(FT0DetT& ft0Det, int ft0Ch, OffsetsT const& offsetFT0, int i) +{ + ft0Det.calculateChannelCenter(); + auto chPos = ft0Det.getChannelCenter(ft0Ch); + + double x = chPos.X() + offsetFT0[i].getX(); + double y = chPos.Y() + offsetFT0[i].getY(); + double z = chPos.Z() + offsetFT0[i].getZ(); + + // If this channel belongs to FT0C, flip z (matches your original intent) + const bool isC = (ft0Ch >= kFT0AChannels); + if (isC) { + z = -z; + } + + const double r = std::sqrt(x * x + y * y); + const double theta = std::atan2(r, z); + return -std::log(std::tan(0.5 * theta)); +} + +template +inline bool getPhiEtaFromFitBit(FT0DetT& ft0Det, + int bit, + OffsetsT const& offsetFT0, + int iRunOffset, + double& phi, + double& eta) +{ + auto ref = decodeFitBit(bit); + if (ref.det != FitBitRef::Det::FT0) { + return false; + } + + // packed mapping: bit index == ft0Ch in the detector numbering you used for getChannelCenter + // FT0A: 0..95, FT0C: 96..207 + const int ft0Ch = bit; + + // (Optional) if ft0Det expects FT0C channels numbered 0..111 separately, then map: + // const int ft0Ch = (bit < kFT0AChannels) ? bit : (bit - kFT0AChannels); + + phi = getPhiFT0_fromChannel(ft0Det, ft0Ch, offsetFT0, iRunOffset); + eta = getEtaFT0_fromChannel(ft0Det, ft0Ch, offsetFT0, iRunOffset); + return true; + +} + + // ----------------------------------------------------------------------------- template diff --git a/PWGUD/DataModel/UDTables.h b/PWGUD/DataModel/UDTables.h index 8a4251ec74d..fa03b4b7111 100644 --- a/PWGUD/DataModel/UDTables.h +++ b/PWGUD/DataModel/UDTables.h @@ -406,6 +406,45 @@ DECLARE_SOA_TABLE(UDTracksFlags, "AOD", "UDTRACKFLAG", DECLARE_SOA_TABLE(UDTracksLabels, "AOD", "UDTRACKLABEL", udtrack::TrackId); +/* +using Bits256 = std::array; + +DECLARE_SOA_COLUMN(FT0BitsThr1MIP, ft0BitsThr1MIP, Bits256); +DECLARE_SOA_COLUMN(FT0BitsThr2MIP, ft0BitsThr2MIP, Bits256); + +DECLARE_SOA_TABLE(UDCollisionFT0Bits, "AOD", "UDCOLLFT0BITS", + o2::soa::Index<>, + FT0BitsThr1MIP, + FT0BitsThr2MIP); +*/ + +namespace udcollfitbits +{ + DECLARE_SOA_COLUMN(Thr1W0, thr1W0, uint64_t); /// 1 MIP thresholds for FT0A ch 0 - ch 63 + DECLARE_SOA_COLUMN(Thr1W1, thr1W1, uint64_t); /// 1 MIP thresholds for FT0A ch 64 - ch 96 & FT0C ch 0 - ch 31 + DECLARE_SOA_COLUMN(Thr1W2, thr1W2, uint64_t); /// 1 MIP thresholds for FT0C ch 32 - ch 96 + DECLARE_SOA_COLUMN(Thr1W3, thr1W3, uint64_t); /// 1 MIP thresholds for FT0C ch 97 - 112 & FV0 0 - 47 + + DECLARE_SOA_COLUMN(Thr2W0, thr2W0, uint64_t); /// 2 MIP thresholds for FT0A ch 0 - ch 63 + DECLARE_SOA_COLUMN(Thr2W1, thr2W1, uint64_t); /// 2 MIP thresholds for FT0A ch 63 - ch 96 & FT0C ch 0 - ch 31 + DECLARE_SOA_COLUMN(Thr2W2, thr2W2, uint64_t); /// 2 MIP thresholds for FT0C ch 32 - ch 96 + DECLARE_SOA_COLUMN(Thr2W3, thr2W3, uint64_t); /// 2 MIP thresholds for FT0C ch 97 - 112 & FV0 0 - 47 +} // namespace udcollfitbits + +DECLARE_SOA_TABLE(UDCollisionFT0Bits, "AOD", "UDCOLLFITBITS", + o2::soa::Index<>, + udcollfitbits::Thr1W0, /// 1 MIP thresholds for FT0A ch 0 - ch 63 + udcollfitbits::Thr1W1, /// 1 MIP thresholds for FT0A ch 63 - ch 96 & FT0C ch 0 - ch 31 + udcollfitbits::Thr1W2, /// 1 MIP thresholds for FT0C ch 32 - ch 96 + udcollfitbits::Thr1W3, /// 1 MIP thresholds for FT0C ch 97 - 112 & FV0 0 - 47 + udcollfitbits::Thr2W0, /// 2 MIP thresholds for FT0A ch 0 - ch 63 + udcollfitbits::Thr2W1, /// 2 MIP thresholds for FT0A ch 63 - ch 96 & FT0C ch 0 - ch 31 + udcollfitbits::Thr2W2, /// 2 MIP thresholds for FT0C ch 32 - ch 96 + udcollfitbits::Thr2W3 /// 2 MIP thresholds for FT0C ch 97 - 112 & FV0 0 - 47 + ); + + + using UDTrack = UDTracks::iterator; using UDTrackCov = UDTracksCov::iterator; using UDTrackExtra = UDTracksExtra::iterator; diff --git a/PWGUD/TableProducer/SGCandProducer.cxx b/PWGUD/TableProducer/SGCandProducer.cxx index c3cb33cb43e..c2cbbba0fae 100644 --- a/PWGUD/TableProducer/SGCandProducer.cxx +++ b/PWGUD/TableProducer/SGCandProducer.cxx @@ -45,6 +45,7 @@ #include #include #include +#include using namespace o2; using namespace o2::framework; @@ -87,6 +88,14 @@ struct SGCandProducer { Configurable storeSG{"storeSG", true, "Store SG events in the output"}; Configurable storeDG{"storeDG", true, "Store DG events in the output"}; + Configurable saveFITbitsets{"saveFITbitsets", true, "Write FT0 and FV0 bitset tables to output"}; + Configurable thr1_FV0A{"thr1_FV0A", 50., "Threshold for 1 MIP in FV0A"}; + Configurable thr1_FT0A{"thr1_FT0A", 50., "Threshold for 1 MIP in FT0A"}; + Configurable thr1_FT0C{"thr1_FT0C", 50., "Threshold for 1 MIP in FT0C"}; + Configurable thr2_FV0A{"thr2_FV0A", 100., "Threshold for 2 MIP in FV0A"}; + Configurable thr2_FT0A{"thr2_FT0A", 100., "Threshold for 2 MIP in FT0A"}; + Configurable thr2_FT0C{"thr2_FT0C", 100., "Threshold for 2 MIP in FT0C"}; + Configurable isGoodRCTCollision{"isGoodRCTCollision", true, "Check RCT flags for FT0,ITS,TPC and tracking"}; Configurable isGoodRCTZdc{"isGoodRCTZdc", false, "Check RCT flags for ZDC if present in run"}; @@ -119,6 +128,7 @@ struct SGCandProducer { Produces outputFwdTracks; Produces outputFwdTracksExtra; Produces outputTracksLabel; + Produces outputFT0Bits; // initialize histogram registry HistogramRegistry registry{ @@ -274,6 +284,60 @@ struct SGCandProducer { } } } + + static inline void setBit(uint64_t w[4], int bit, bool val) + { + if (!val) { + return; + } + const int word = bit >> 6; // /64 + const int offs = bit & 63; // %64 + w[word] |= (uint64_t(1) << offs); + } + + template + static inline void buildFT0FV0Words(TFT0 const& ft0, TFV0A const& fv0a, + uint64_t thr1[4], uint64_t thr2[4], + float thr1_FT0A = 25., float thr1_FT0C = 50., float thr1_FV0A = 50., + float thr2_FT0A = 50., float thr2_FT0C = 100., float thr2_FV0A = 100. ) + { + thr1[0]=thr1[1]=thr1[2]=thr1[3]=0ull; + thr2[0]=thr2[1]=thr2[2]=thr2[3]=0ull; + + constexpr int kFT0AOffset = 0; + constexpr int kFT0COffset = 96; + constexpr int kFV0Offset = 208; + + // Setting bits for FT0A + auto ampsA = ft0.amplitudeA(); + const int nA = std::min(ampsA.size(), 96); + for (int i = 0; i < nA; ++i) { + const auto a = ampsA[i]; + setBit(thr1, kFT0AOffset + i, a >= thr1_FT0A); + setBit(thr2, kFT0AOffset + i, a >= thr2_FT0A); + } + + // Setting bits for FT0C + auto ampsC = ft0.amplitudeC(); + const int nC = std::min(ampsC.size(), 112); + for (int i = 0; i < nC; ++i) { + const auto a = ampsC[i]; + setBit(thr1, kFT0COffset + i, a >= thr1_FT0C); + setBit(thr2, kFT0COffset + i, a >= thr2_FT0C); + } + + // Setting bits for FV0A + auto ampsV = fv0a.amplitude(); + const int nV = std::min(ampsV.size(), 48); + for (int i = 0; i < nV; ++i) { + const auto a = ampsV[i]; + setBit(thr1, kFV0Offset + i, a >= thr1_FV0A); + setBit(thr2, kFV0Offset + i, a >= thr2_FV0A); + } + } + + + PROCESS_SWITCH(SGCandProducer, processCountersTrg, "Produce trigger counters and luminosity histograms", true); @@ -390,6 +454,17 @@ struct SGCandProducer { fitInfo.BBFDDApf, fitInfo.BBFDDCpf, fitInfo.BGFDDApf, fitInfo.BGFDDCpf); outputCollisionSelExtras(chFT0A, chFT0C, chFDDA, chFDDC, chFV0A, occ, ir, trs, trofs, hmpr, tfb, itsROFb, sbp, zVtxFT0vPv, vtxITSTPC, collision.rct_raw()); outputCollsLabels(collision.globalIndex()); + + uint64_t w1[4] = {0ull, 0ull, 0ull, 0ull}; + uint64_t w2[4] = {0ull, 0ull, 0ull, 0ull}; + + if ( saveFITbitsets && newbc.has_foundFT0() && newbc.has_fv0a() ) { + buildFT0FV0Words(newbc.ft0(), newbc.fv0a(), w1, w2, thr1_FT0A, thr1_FT0C, thr1_FV0A, thr2_FT0A, thr2_FT0C, thr2_FV0A); + } + + outputFT0Bits(w1[0], w1[1], w1[2], w1[3], + w2[0], w2[1], w2[2], w2[3]); + if (newbc.has_zdc()) { auto zdc = newbc.zdc(); udZdcsReduced(outputCollisions.lastIndex(), zdc.timeZNA(), zdc.timeZNC(), zdc.energyCommonZNA(), zdc.energyCommonZNC()); diff --git a/PWGUD/Tasks/CMakeLists.txt b/PWGUD/Tasks/CMakeLists.txt index 62c6abc732f..bd33c051da1 100644 --- a/PWGUD/Tasks/CMakeLists.txt +++ b/PWGUD/Tasks/CMakeLists.txt @@ -268,3 +268,8 @@ o2physics_add_dpl_workflow(sg-exclusive-jpsi-midrapidity SOURCES sgExclusiveJpsiMidrapidity.cxx PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore COMPONENT_NAME Analysis) + +o2physics_add_dpl_workflow(fitbit-mapping + SOURCES upcTestFITBitMapping.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore + COMPONENT_NAME PWGUD) \ No newline at end of file diff --git a/PWGUD/Tasks/upcTestFITBitMapping.cxx b/PWGUD/Tasks/upcTestFITBitMapping.cxx new file mode 100644 index 00000000000..98189bfd9e5 --- /dev/null +++ b/PWGUD/Tasks/upcTestFITBitMapping.cxx @@ -0,0 +1,108 @@ +// UpcTestFITBitMapping.cxx +// Minimal task to test FT0 bitset -> (eta,phi) mapping using udhelpers::getPhiEtaFromFitBit +// +// Run example (from O2Physics build dir): +// ./stage/bin/o2-pwgud-ud-fitbit-mapping --aod-file AO2D_UD.root -b +// +// Note: Add this file to PWGUD/Tasks/CMakeLists.txt. + +#include "Framework/AnalysisTask.h" +#include "Framework/runDataProcessing.h" +#include "Framework/HistogramRegistry.h" + +#include "PWGUD/DataModel/UDTables.h" // aod::UDCollisionFT0Bits +#include "PWGUD/Core/UDHelpers.h" // udhelpers::Bits256, makeBits256, testBit, getPhiEtaFromFitBit +#include "FT0Base/Geometry.h" // o2::ft0::Geometry + +#include +#include + +using namespace o2; +using namespace o2::framework; + +struct UpcTestFITBitMapping +{ + Configurable whichThr{"whichThr", 1, "Use 1=Thr1 bits or 2=Thr2 bits"}; + Configurable maxEvents{"maxEvents", -1, "Process at most this many rows (-1 = all)"}; + + // Minimal offset container compatible with UDHelpers.h expectations: getX/getY/getZ + struct OffsetXYZ { + double x{0}, y{0}, z{0}; + double getX() const { return x; } + double getY() const { return y; } + double getZ() const { return z; } + }; + + std::array offsetFT0{}; // iRunOffset = 0 for now + int iRunOffset = 0; + + o2::ft0::Geometry ft0Det{}; + + HistogramRegistry registry{ + "registry", + { + {"hPhiA", "FT0A #varphi;#varphi;counts", {HistType::kTH1F, {{18, 0.0, 2.*M_PI}}}}, + {"hEtaA", "FT0A #eta;#eta;counts", {HistType::kTH1F, {{8, 3.5, 5.0}}}}, + {"hEtaPhiA", "FT0A #eta vs #varphi;#eta;#varphi", {HistType::kTH2F, {{8, 3.5, 5.0}, {18, 0.0, 2.*M_PI}}}}, + + {"hPhiC", "FT0C #varphi;#varphi;counts", {HistType::kTH1F, {{18, 0.0, 2.*M_PI}}}}, + {"hEtaC", "FT0C #eta;#eta;counts", {HistType::kTH1F, {{8, -3.5, -2.0}}}}, + {"hEtaPhiC", "FT0C #eta vs #varphi;#eta;#varphi", {HistType::kTH2F, {{8, -3.5, -2.0}, {18, 0.0, 2.*M_PI}}}}, + }}; + + void init(InitContext&) + { + // UDHelpers calls calculateChannelCenter() inside, but doing it once here is fine. + ft0Det.calculateChannelCenter(); + } + + void process(aod::UDCollisionFT0Bits const& bitsTable) + { + int64_t nProcessed = 0; + + for (auto const& row : bitsTable) { + if (maxEvents >= 0 && nProcessed >= maxEvents) { + break; + } + ++nProcessed; + + // Use udhelpers' canonical packed type + builder + udhelpers::Bits256 w{}; + if (whichThr == 2) { + w = udhelpers::makeBits256(row.thr2W0(), row.thr2W1(), row.thr2W2(), row.thr2W3()); + } else { + w = udhelpers::makeBits256(row.thr1W0(), row.thr1W1(), row.thr1W2(), row.thr1W3()); + } + + // Loop FT0 bits only (0..207). FV0 starts at 208 but ignored here. + for (int bit = 0; bit < udhelpers::kFT0Bits; ++bit) { + if (!udhelpers::testBit(w, bit)) { + continue; + } + + double phi = 0., eta = 0.; + const bool ok = udhelpers::getPhiEtaFromFitBit(ft0Det, bit, offsetFT0, iRunOffset, phi, eta); + if (!ok) { + continue; + } + + if (bit < udhelpers::kFT0AChannels) { + registry.fill(HIST("hPhiA"), phi); + registry.fill(HIST("hEtaA"), eta); + registry.fill(HIST("hEtaPhiA"), eta, phi); + } else { + registry.fill(HIST("hPhiC"), phi); + registry.fill(HIST("hEtaC"), eta); + registry.fill(HIST("hEtaPhiC"), eta, phi); + } + } + } + } +}; + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + return WorkflowSpec{ + adaptAnalysisTask(cfgc, TaskName{"fitbit-mapping"}) + }; +} From 18a90acc521a2a0f2198a7e19c64916695250ce4 Mon Sep 17 00:00:00 2001 From: sandor-lokos Date: Mon, 23 Mar 2026 17:52:57 +0100 Subject: [PATCH 02/11] Refactored FIT bit implementation. NOTE: no functionality changed! --- PWGUD/Core/CMakeLists.txt | 8 +++ PWGUD/Core/FITCutParHolder.cxx | 40 +++++++++++++ PWGUD/Core/FITCutParHolder.h | 60 ++++++++++++++++++++ PWGUD/Core/FITCutParHolderLinkDef.h | 19 +++++++ PWGUD/Core/UDHelpers.h | 58 ++++++++++++++++--- PWGUD/DataModel/UDTables.h | 14 +---- PWGUD/TableProducer/CMakeLists.txt | 2 +- PWGUD/TableProducer/SGCandProducer.cxx | 78 +++++--------------------- PWGUD/Tasks/upcTestFITBitMapping.cxx | 4 +- 9 files changed, 196 insertions(+), 87 deletions(-) create mode 100644 PWGUD/Core/FITCutParHolder.cxx create mode 100644 PWGUD/Core/FITCutParHolder.h create mode 100644 PWGUD/Core/FITCutParHolderLinkDef.h diff --git a/PWGUD/Core/CMakeLists.txt b/PWGUD/Core/CMakeLists.txt index 7686edd8eea..4a44b937160 100644 --- a/PWGUD/Core/CMakeLists.txt +++ b/PWGUD/Core/CMakeLists.txt @@ -64,3 +64,11 @@ o2physics_add_library(decayTree o2physics_target_root_dictionary(decayTree HEADERS decayTree.h LINKDEF decayTreeLinkDef.h) + +o2physics_add_library(FITCutParHolder + SOURCES FITCutParHolder.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore) + +o2physics_target_root_dictionary(FITCutParHolder + HEADERS FITCutParHolder.h + LINKDEF FITCutParHolderLinkDef.h) \ No newline at end of file diff --git a/PWGUD/Core/FITCutParHolder.cxx b/PWGUD/Core/FITCutParHolder.cxx new file mode 100644 index 00000000000..401f560c968 --- /dev/null +++ b/PWGUD/Core/FITCutParHolder.cxx @@ -0,0 +1,40 @@ +#include "FITCutParHolder.h" + +// setter +void FITCutParHolder::SetSaveFITbitsets(bool saveFITbitsets) +{ + mSaveFITbitsets = saveFITbitsets; +} +void FITCutParHolder::SetThr1FV0A(float thr1_FV0A) +{ + mThr1FV0A = thr1_FV0A; +} +void FITCutParHolder::SetThr1FT0A(float thr1_FT0A) +{ + mThr1FT0A = thr1_FT0A; +} +void FITCutParHolder::SetThr1FT0C(float thr1_FT0C) +{ + mThr1FT0C = thr1_FT0C; +} +void FITCutParHolder::SetThr2FV0A(float thr2_FV0A) +{ + mThr2FV0A = thr2_FV0A; +} +void FITCutParHolder::SetThr2FT0A(float thr2_FT0A) +{ + mThr2FT0A = thr2_FT0A; +} +void FITCutParHolder::SetThr2FT0C(float thr2_FT0C) +{ + mThr2FT0C = thr2_FT0C; +} + +// getter +bool FITCutParHolder::saveFITbitsets() const { return mSaveFITbitsets; } +float FITCutParHolder::thr1_FV0A() const { return mThr1FV0A; } +float FITCutParHolder::thr1_FT0A() const { return mThr1FT0A; } +float FITCutParHolder::thr1_FT0C() const { return mThr1FT0C; } +float FITCutParHolder::thr2_FV0A() const { return mThr2FV0A; } +float FITCutParHolder::thr2_FT0A() const { return mThr2FT0A; } +float FITCutParHolder::thr2_FT0C() const { return mThr2FT0C; } \ No newline at end of file diff --git a/PWGUD/Core/FITCutParHolder.h b/PWGUD/Core/FITCutParHolder.h new file mode 100644 index 00000000000..d87afea1ecd --- /dev/null +++ b/PWGUD/Core/FITCutParHolder.h @@ -0,0 +1,60 @@ +#ifndef PWGUD_CORE_FITCUTPARHOLDER_H_ +#define PWGUD_CORE_FITCUTPARHOLDER_H_ + +#include + +// object to hold customizable FIT bit thresholds +class FITCutParHolder +{ + public: + // constructor + FITCutParHolder(bool saveFITbitsets = true, + float thr1_FV0A = 50., + float thr1_FT0A = 50., + float thr1_FT0C = 50., + float thr2_FV0A = 100., + float thr2_FT0A = 100., + float thr2_FT0C = 100.) + : mSaveFITbitsets{saveFITbitsets}, + mThr1FV0A{thr1_FV0A}, + mThr1FT0A{thr1_FT0A}, + mThr1FT0C{thr1_FT0C}, + mThr2FV0A{thr2_FV0A}, + mThr2FT0A{thr2_FT0A}, + mThr2FT0C{thr2_FT0C} + { + } + + // setters + void SetSaveFITbitsets(bool); + void SetThr1FV0A(float); + void SetThr1FT0A(float); + void SetThr1FT0C(float); + void SetThr2FV0A(float); + void SetThr2FT0A(float); + void SetThr2FT0C(float); + + // getters + bool saveFITbitsets() const; + float thr1_FV0A() const; + float thr1_FT0A() const; + float thr1_FT0C() const; + float thr2_FV0A() const; + float thr2_FT0A() const; + float thr2_FT0C() const; + + private: + bool mSaveFITbitsets; + + float mThr1FV0A; + float mThr1FT0A; + float mThr1FT0C; + + float mThr2FV0A; + float mThr2FT0A; + float mThr2FT0C; + + ClassDefNV(FITCutParHolder, 1); +}; + +#endif // PWGUD_CORE_FITCUTPARHOLDER_H_ \ No newline at end of file diff --git a/PWGUD/Core/FITCutParHolderLinkDef.h b/PWGUD/Core/FITCutParHolderLinkDef.h new file mode 100644 index 00000000000..70341fc2ee4 --- /dev/null +++ b/PWGUD/Core/FITCutParHolderLinkDef.h @@ -0,0 +1,19 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. +#ifndef PWGUD_CORE_FITCUTPARHOLDERLINKDEF_H_ +#define PWGUD_CORE_FITCUTPARHOLDERLINKDEF_H_ + +#pragma link off all globals; +#pragma link off all classes; +#pragma link off all functions; +#pragma link C++ class FITCutParHolder + ; + +#endif // PWGUD_CORE_FITCUTPARHOLDERLINKDEF_H_ \ No newline at end of file diff --git a/PWGUD/Core/UDHelpers.h b/PWGUD/Core/UDHelpers.h index 824857a6343..cd0d10e2626 100644 --- a/PWGUD/Core/UDHelpers.h +++ b/PWGUD/Core/UDHelpers.h @@ -535,6 +535,56 @@ bool FITveto(T const& bc, DGCutparHolder const& diffCuts) return false; } +inline void setBit(uint64_t w[4], int bit, bool val) +{ + if (!val) { + return; + } + const int word = bit >> 6; + const int offs = bit & 63; + w[word] |= (uint64_t(1) << offs); +} + +template +inline void buildFT0FV0Words(TFT0 const& ft0, TFV0A const& fv0a, + uint64_t thr1[4], uint64_t thr2[4], + float thr1_FT0A = 25., float thr1_FT0C = 50., float thr1_FV0A = 50., + float thr2_FT0A = 50., float thr2_FT0C = 100., float thr2_FV0A = 100.) +{ + thr1[0] = thr1[1] = thr1[2] = thr1[3] = 0ull; + thr2[0] = thr2[1] = thr2[2] = thr2[3] = 0ull; + + constexpr int kFT0AOffset = 0; + constexpr int kFT0COffset = 96; + constexpr int kFV0Offset = 208; + + auto ampsA = ft0.amplitudeA(); + const int nA = std::min(ampsA.size(), 96); + for (int i = 0; i < nA; ++i) { + const auto a = ampsA[i]; + setBit(thr1, kFT0AOffset + i, a >= thr1_FT0A); + setBit(thr2, kFT0AOffset + i, a >= thr2_FT0A); + } + + auto ampsC = ft0.amplitudeC(); + const int nC = std::min(ampsC.size(), 112); + for (int i = 0; i < nC; ++i) { + const auto a = ampsC[i]; + setBit(thr1, kFT0COffset + i, a >= thr1_FT0C); + setBit(thr2, kFT0COffset + i, a >= thr2_FT0C); + } + + auto ampsV = fv0a.amplitude(); + const int nV = std::min(ampsV.size(), 48); + for (int i = 0; i < nV; ++i) { + const auto a = ampsV[i]; + setBit(thr1, kFV0Offset + i, a >= thr1_FV0A); + setBit(thr2, kFV0Offset + i, a >= thr2_FV0A); + } +} + + + // ----------------------------------------------------------------------------- // return eta and phi of a given FIT channel based on the bitset // Bit layout contract: @@ -542,7 +592,7 @@ constexpr int kFT0Bits = 208; // FT0 total channels constexpr int kFV0Bits = 48; // FV0A channels constexpr int kTotalBits = 256; // 4*64 -// FT0 side split (adjust if your channelization differs in your tag) +// FT0 side split constexpr int kFT0AChannels = 96; // FT0A channels are [0..95] constexpr int kFT0CChannels = 112; // FT0C channels are [96..207] static_assert(kFT0AChannels + kFT0CChannels == kFT0Bits); @@ -586,7 +636,6 @@ inline FitBitRef decodeFitBit(int bit) return out; } -// You can pass whatever type offsetFT0 is (vector/array/span of objects with getX/Y/Z). template inline double getPhiFT0_fromChannel(FT0DetT& ft0Det, int ft0Ch, OffsetsT const& offsetFT0, int i) { @@ -633,13 +682,8 @@ inline bool getPhiEtaFromFitBit(FT0DetT& ft0Det, return false; } - // packed mapping: bit index == ft0Ch in the detector numbering you used for getChannelCenter // FT0A: 0..95, FT0C: 96..207 const int ft0Ch = bit; - - // (Optional) if ft0Det expects FT0C channels numbered 0..111 separately, then map: - // const int ft0Ch = (bit < kFT0AChannels) ? bit : (bit - kFT0AChannels); - phi = getPhiFT0_fromChannel(ft0Det, ft0Ch, offsetFT0, iRunOffset); eta = getEtaFT0_fromChannel(ft0Det, ft0Ch, offsetFT0, iRunOffset); return true; diff --git a/PWGUD/DataModel/UDTables.h b/PWGUD/DataModel/UDTables.h index fa03b4b7111..b6978027df7 100644 --- a/PWGUD/DataModel/UDTables.h +++ b/PWGUD/DataModel/UDTables.h @@ -406,18 +406,6 @@ DECLARE_SOA_TABLE(UDTracksFlags, "AOD", "UDTRACKFLAG", DECLARE_SOA_TABLE(UDTracksLabels, "AOD", "UDTRACKLABEL", udtrack::TrackId); -/* -using Bits256 = std::array; - -DECLARE_SOA_COLUMN(FT0BitsThr1MIP, ft0BitsThr1MIP, Bits256); -DECLARE_SOA_COLUMN(FT0BitsThr2MIP, ft0BitsThr2MIP, Bits256); - -DECLARE_SOA_TABLE(UDCollisionFT0Bits, "AOD", "UDCOLLFT0BITS", - o2::soa::Index<>, - FT0BitsThr1MIP, - FT0BitsThr2MIP); -*/ - namespace udcollfitbits { DECLARE_SOA_COLUMN(Thr1W0, thr1W0, uint64_t); /// 1 MIP thresholds for FT0A ch 0 - ch 63 @@ -431,7 +419,7 @@ namespace udcollfitbits DECLARE_SOA_COLUMN(Thr2W3, thr2W3, uint64_t); /// 2 MIP thresholds for FT0C ch 97 - 112 & FV0 0 - 47 } // namespace udcollfitbits -DECLARE_SOA_TABLE(UDCollisionFT0Bits, "AOD", "UDCOLLFITBITS", +DECLARE_SOA_TABLE(UDCollisionFITBits, "AOD", "UDCOLLFITBITS", o2::soa::Index<>, udcollfitbits::Thr1W0, /// 1 MIP thresholds for FT0A ch 0 - ch 63 udcollfitbits::Thr1W1, /// 1 MIP thresholds for FT0A ch 63 - ch 96 & FT0C ch 0 - ch 31 diff --git a/PWGUD/TableProducer/CMakeLists.txt b/PWGUD/TableProducer/CMakeLists.txt index c2fa718ddd0..6e2b62968db 100644 --- a/PWGUD/TableProducer/CMakeLists.txt +++ b/PWGUD/TableProducer/CMakeLists.txt @@ -18,7 +18,7 @@ o2physics_add_dpl_workflow(dgcand-producer o2physics_add_dpl_workflow(sgcand-producer SOURCES SGCandProducer.cxx - PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2Physics::SGCutParHolder O2Physics::AnalysisCCDB + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2Physics::SGCutParHolder O2Physics::AnalysisCCDB O2Physics::FITCutParHolder COMPONENT_NAME Analysis) o2physics_add_dpl_workflow(dgbccand-producer diff --git a/PWGUD/TableProducer/SGCandProducer.cxx b/PWGUD/TableProducer/SGCandProducer.cxx index c2cbbba0fae..0a94da3d15b 100644 --- a/PWGUD/TableProducer/SGCandProducer.cxx +++ b/PWGUD/TableProducer/SGCandProducer.cxx @@ -20,6 +20,7 @@ #include "PWGUD/Core/SGSelector.h" #include "PWGUD/Core/UPCHelpers.h" +#include "PWGUD/Core/FITCutParHolder.h" #include "PWGUD/DataModel/UDTables.h" #include "Common/CCDB/EventSelectionParams.h" @@ -76,6 +77,8 @@ struct SGCandProducer { // get an SGCutparHolder SGCutParHolder sameCuts = SGCutParHolder(); // SGCutparHolder Configurable SGCuts{"SGCuts", {}, "SG event cuts"}; + FITCutParHolder fitCuts = FITCutParHolder(); + Configurable FITCuts{"FITCuts", {}, "FIT bitset cuts"}; Configurable verboseInfo{"verboseInfo", false, "Print general info to terminal; default it false."}; Configurable saveAllTracks{"saveAllTracks", true, "save only PV contributors or all tracks associated to a collision"}; Configurable savenonPVCITSOnlyTracks{"savenonPVCITSOnlyTracks", false, "save non PV contributors with ITS only information"}; @@ -89,12 +92,6 @@ struct SGCandProducer { Configurable storeDG{"storeDG", true, "Store DG events in the output"}; Configurable saveFITbitsets{"saveFITbitsets", true, "Write FT0 and FV0 bitset tables to output"}; - Configurable thr1_FV0A{"thr1_FV0A", 50., "Threshold for 1 MIP in FV0A"}; - Configurable thr1_FT0A{"thr1_FT0A", 50., "Threshold for 1 MIP in FT0A"}; - Configurable thr1_FT0C{"thr1_FT0C", 50., "Threshold for 1 MIP in FT0C"}; - Configurable thr2_FV0A{"thr2_FV0A", 100., "Threshold for 2 MIP in FV0A"}; - Configurable thr2_FT0A{"thr2_FT0A", 100., "Threshold for 2 MIP in FT0A"}; - Configurable thr2_FT0C{"thr2_FT0C", 100., "Threshold for 2 MIP in FT0C"}; Configurable isGoodRCTCollision{"isGoodRCTCollision", true, "Check RCT flags for FT0,ITS,TPC and tracking"}; Configurable isGoodRCTZdc{"isGoodRCTZdc", false, "Check RCT flags for ZDC if present in run"}; @@ -128,7 +125,7 @@ struct SGCandProducer { Produces outputFwdTracks; Produces outputFwdTracksExtra; Produces outputTracksLabel; - Produces outputFT0Bits; + Produces outputFITBits; // initialize histogram registry HistogramRegistry registry{ @@ -284,60 +281,6 @@ struct SGCandProducer { } } } - - static inline void setBit(uint64_t w[4], int bit, bool val) - { - if (!val) { - return; - } - const int word = bit >> 6; // /64 - const int offs = bit & 63; // %64 - w[word] |= (uint64_t(1) << offs); - } - - template - static inline void buildFT0FV0Words(TFT0 const& ft0, TFV0A const& fv0a, - uint64_t thr1[4], uint64_t thr2[4], - float thr1_FT0A = 25., float thr1_FT0C = 50., float thr1_FV0A = 50., - float thr2_FT0A = 50., float thr2_FT0C = 100., float thr2_FV0A = 100. ) - { - thr1[0]=thr1[1]=thr1[2]=thr1[3]=0ull; - thr2[0]=thr2[1]=thr2[2]=thr2[3]=0ull; - - constexpr int kFT0AOffset = 0; - constexpr int kFT0COffset = 96; - constexpr int kFV0Offset = 208; - - // Setting bits for FT0A - auto ampsA = ft0.amplitudeA(); - const int nA = std::min(ampsA.size(), 96); - for (int i = 0; i < nA; ++i) { - const auto a = ampsA[i]; - setBit(thr1, kFT0AOffset + i, a >= thr1_FT0A); - setBit(thr2, kFT0AOffset + i, a >= thr2_FT0A); - } - - // Setting bits for FT0C - auto ampsC = ft0.amplitudeC(); - const int nC = std::min(ampsC.size(), 112); - for (int i = 0; i < nC; ++i) { - const auto a = ampsC[i]; - setBit(thr1, kFT0COffset + i, a >= thr1_FT0C); - setBit(thr2, kFT0COffset + i, a >= thr2_FT0C); - } - - // Setting bits for FV0A - auto ampsV = fv0a.amplitude(); - const int nV = std::min(ampsV.size(), 48); - for (int i = 0; i < nV; ++i) { - const auto a = ampsV[i]; - setBit(thr1, kFV0Offset + i, a >= thr1_FV0A); - setBit(thr2, kFV0Offset + i, a >= thr2_FV0A); - } - } - - - PROCESS_SWITCH(SGCandProducer, processCountersTrg, "Produce trigger counters and luminosity histograms", true); @@ -458,11 +401,17 @@ struct SGCandProducer { uint64_t w1[4] = {0ull, 0ull, 0ull, 0ull}; uint64_t w2[4] = {0ull, 0ull, 0ull, 0ull}; - if ( saveFITbitsets && newbc.has_foundFT0() && newbc.has_fv0a() ) { - buildFT0FV0Words(newbc.ft0(), newbc.fv0a(), w1, w2, thr1_FT0A, thr1_FT0C, thr1_FV0A, thr2_FT0A, thr2_FT0C, thr2_FV0A); + if (fitCuts.saveFITbitsets() && newbc.has_foundFT0() && newbc.has_fv0a()) { + udhelpers::buildFT0FV0Words(newbc.ft0(), newbc.fv0a(), w1, w2, + fitCuts.thr1_FT0A(), + fitCuts.thr1_FT0C(), + fitCuts.thr1_FV0A(), + fitCuts.thr2_FT0A(), + fitCuts.thr2_FT0C(), + fitCuts.thr2_FV0A()); } - outputFT0Bits(w1[0], w1[1], w1[2], w1[3], + outputFITBits(w1[0], w1[1], w1[2], w1[3], w2[0], w2[1], w2[2], w2[3]); if (newbc.has_zdc()) { @@ -503,6 +452,7 @@ struct SGCandProducer { ccdb->setCaching(true); ccdb->setFatalWhenNull(false); sameCuts = (SGCutParHolder)SGCuts; + fitCuts = (FITCutParHolder)FITCuts; // add histograms for the different process functions histPointers.clear(); diff --git a/PWGUD/Tasks/upcTestFITBitMapping.cxx b/PWGUD/Tasks/upcTestFITBitMapping.cxx index 98189bfd9e5..8469b9c127b 100644 --- a/PWGUD/Tasks/upcTestFITBitMapping.cxx +++ b/PWGUD/Tasks/upcTestFITBitMapping.cxx @@ -10,7 +10,7 @@ #include "Framework/runDataProcessing.h" #include "Framework/HistogramRegistry.h" -#include "PWGUD/DataModel/UDTables.h" // aod::UDCollisionFT0Bits +#include "PWGUD/DataModel/UDTables.h" // aod::UDCollisionFITBits #include "PWGUD/Core/UDHelpers.h" // udhelpers::Bits256, makeBits256, testBit, getPhiEtaFromFitBit #include "FT0Base/Geometry.h" // o2::ft0::Geometry @@ -56,7 +56,7 @@ struct UpcTestFITBitMapping ft0Det.calculateChannelCenter(); } - void process(aod::UDCollisionFT0Bits const& bitsTable) + void process(aod::UDCollisionFITBits const& bitsTable) { int64_t nProcessed = 0; From cb615e96ab845e4645e66e7b2bcba28cc5fcb156 Mon Sep 17 00:00:00 2001 From: sandor-lokos Date: Tue, 24 Mar 2026 13:59:10 +0100 Subject: [PATCH 03/11] deleted processCountersTrg function from SGCandProducer --- PWGUD/TableProducer/SGCandProducer.cxx | 72 -------------------------- 1 file changed, 72 deletions(-) diff --git a/PWGUD/TableProducer/SGCandProducer.cxx b/PWGUD/TableProducer/SGCandProducer.cxx index 0a94da3d15b..25629e475e4 100644 --- a/PWGUD/TableProducer/SGCandProducer.cxx +++ b/PWGUD/TableProducer/SGCandProducer.cxx @@ -212,78 +212,6 @@ struct SGCandProducer { outputTracksLabel(track.globalIndex()); } - // function to process trigger counters, accounting for BC selection bits - void processCountersTrg(BCs const& bcs, aod::FT0s const&, aod::Zdcs const&) - { - const auto& firstBc = bcs.iteratorAt(0); - if (runNumber != firstBc.runNumber()) - runNumber = firstBc.runNumber(); - - auto hCountersTrg = getHist(TH1, "reco/hCountersTrg"); - auto hCountersTrgBcSel = getHist(TH1, "reco/hCountersTrgBcSel"); - auto hLumi = getHist(TH1, "reco/hLumi"); - auto hLumiBcSel = getHist(TH1, "reco/hLumiBcSel"); - - // Cross sections in ub. Using dummy -1 if lumi estimator is not reliable - float csTCE = 10.36e6; - const float csZEM = 415.2e6; // see AN: https://alice-notes.web.cern.ch/node/1515 - const float csZNC = 214.5e6; // see AN: https://alice-notes.web.cern.ch/node/1515 - if (runNumber > 543437 && runNumber < 543514) { - csTCE = 8.3e6; - } - if (runNumber >= 543514) { - csTCE = 4.10e6; // see AN: https://alice-notes.web.cern.ch/node/1515 - } - - for (const auto& bc : bcs) { - bool hasFT0 = bc.has_foundFT0(); - bool hasZDC = bc.has_foundZDC(); - if (!hasFT0 && !hasZDC) - continue; - bool isSelectedBc = true; - if (rejectAtTFBoundary && !bc.selection_bit(aod::evsel::kNoTimeFrameBorder)) - isSelectedBc = false; - if (noITSROFrameBorder && !bc.selection_bit(aod::evsel::kNoITSROFrameBorder)) - isSelectedBc = false; - if (hasFT0) { - auto ft0TrgMask = bc.ft0().triggerMask(); - if (TESTBIT(ft0TrgMask, o2::fit::Triggers::bitVertex)) { - hCountersTrg->Fill("TVX", 1); - if (isSelectedBc) - hCountersTrgBcSel->Fill("TVX", 1); - } - if (TESTBIT(ft0TrgMask, o2::fit::Triggers::bitVertex) && TESTBIT(ft0TrgMask, o2::fit::Triggers::bitCen)) { - hCountersTrg->Fill("TCE", 1); - hLumi->Fill("TCE", 1. / csTCE); - if (isSelectedBc) { - hCountersTrgBcSel->Fill("TCE", 1); - hLumiBcSel->Fill("TCE", 1. / csTCE); - } - } - } - if (hasZDC) { - if (bc.selection_bit(aod::evsel::kIsBBZNA) || bc.selection_bit(aod::evsel::kIsBBZNC)) { - hCountersTrg->Fill("ZEM", 1); - hLumi->Fill("ZEM", 1. / csZEM); - if (isSelectedBc) { - hCountersTrgBcSel->Fill("ZEM", 1); - hLumiBcSel->Fill("ZEM", 1. / csZEM); - } - } - if (bc.selection_bit(aod::evsel::kIsBBZNC)) { - hCountersTrg->Fill("ZNC", 1); - hLumi->Fill("ZNC", 1. / csZNC); - if (isSelectedBc) { - hCountersTrgBcSel->Fill("ZNC", 1); - hLumiBcSel->Fill("ZNC", 1. / csZNC); - } - } - } - } - } - - PROCESS_SWITCH(SGCandProducer, processCountersTrg, "Produce trigger counters and luminosity histograms", true); - // function to process reconstructed data template void processReco(std::string histdir, TCol const& collision, BCs const& bcs, From 3cd66583c57f94dcd3b1147f5e60055302a8c22f Mon Sep 17 00:00:00 2001 From: sandor-lokos Date: Tue, 24 Mar 2026 19:11:23 +0100 Subject: [PATCH 04/11] Set default thresholds for the FIT bits --- PWGUD/Core/FITCutParHolder.h | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/PWGUD/Core/FITCutParHolder.h b/PWGUD/Core/FITCutParHolder.h index d87afea1ecd..d19ff01c97f 100644 --- a/PWGUD/Core/FITCutParHolder.h +++ b/PWGUD/Core/FITCutParHolder.h @@ -9,12 +9,12 @@ class FITCutParHolder public: // constructor FITCutParHolder(bool saveFITbitsets = true, - float thr1_FV0A = 50., - float thr1_FT0A = 50., - float thr1_FT0C = 50., - float thr2_FV0A = 100., - float thr2_FT0A = 100., - float thr2_FT0C = 100.) + float thr1_FV0A = 8., + float thr1_FT0A = 8., + float thr1_FT0C = 8., + float thr2_FV0A = 20., + float thr2_FT0A = 20., + float thr2_FT0C = 20.) : mSaveFITbitsets{saveFITbitsets}, mThr1FV0A{thr1_FV0A}, mThr1FT0A{thr1_FT0A}, From d0b324b32949d187d6de30772d7f06bf13efbe95 Mon Sep 17 00:00:00 2001 From: sandor-lokos Date: Tue, 24 Mar 2026 20:22:19 +0100 Subject: [PATCH 05/11] Formatting fixes --- PWGUD/Core/CMakeLists.txt | 2 +- PWGUD/Core/FITCutParHolder.cxx | 95 ++++++----- PWGUD/Core/FITCutParHolder.h | 135 ++++++++------- PWGUD/Core/FITCutParHolderLinkDef.h | 43 ++--- PWGUD/Core/SGCutParHolderLinkDef.h | 1 + PWGUD/Core/UDHelpers.h | 32 ++-- PWGUD/DataModel/UDTables.h | 22 ++- PWGUD/TableProducer/SGCandProducer.cxx | 30 ++-- PWGUD/Tasks/upcTestFITBitMapping.cxx | 221 +++++++++++++------------ 9 files changed, 309 insertions(+), 272 deletions(-) diff --git a/PWGUD/Core/CMakeLists.txt b/PWGUD/Core/CMakeLists.txt index 4a44b937160..3a7fe0950cf 100644 --- a/PWGUD/Core/CMakeLists.txt +++ b/PWGUD/Core/CMakeLists.txt @@ -64,7 +64,7 @@ o2physics_add_library(decayTree o2physics_target_root_dictionary(decayTree HEADERS decayTree.h LINKDEF decayTreeLinkDef.h) - + o2physics_add_library(FITCutParHolder SOURCES FITCutParHolder.cxx PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore) diff --git a/PWGUD/Core/FITCutParHolder.cxx b/PWGUD/Core/FITCutParHolder.cxx index 401f560c968..d57e2c8e12f 100644 --- a/PWGUD/Core/FITCutParHolder.cxx +++ b/PWGUD/Core/FITCutParHolder.cxx @@ -1,40 +1,55 @@ -#include "FITCutParHolder.h" - -// setter -void FITCutParHolder::SetSaveFITbitsets(bool saveFITbitsets) -{ - mSaveFITbitsets = saveFITbitsets; -} -void FITCutParHolder::SetThr1FV0A(float thr1_FV0A) -{ - mThr1FV0A = thr1_FV0A; -} -void FITCutParHolder::SetThr1FT0A(float thr1_FT0A) -{ - mThr1FT0A = thr1_FT0A; -} -void FITCutParHolder::SetThr1FT0C(float thr1_FT0C) -{ - mThr1FT0C = thr1_FT0C; -} -void FITCutParHolder::SetThr2FV0A(float thr2_FV0A) -{ - mThr2FV0A = thr2_FV0A; -} -void FITCutParHolder::SetThr2FT0A(float thr2_FT0A) -{ - mThr2FT0A = thr2_FT0A; -} -void FITCutParHolder::SetThr2FT0C(float thr2_FT0C) -{ - mThr2FT0C = thr2_FT0C; -} - -// getter -bool FITCutParHolder::saveFITbitsets() const { return mSaveFITbitsets; } -float FITCutParHolder::thr1_FV0A() const { return mThr1FV0A; } -float FITCutParHolder::thr1_FT0A() const { return mThr1FT0A; } -float FITCutParHolder::thr1_FT0C() const { return mThr1FT0C; } -float FITCutParHolder::thr2_FV0A() const { return mThr2FV0A; } -float FITCutParHolder::thr2_FT0A() const { return mThr2FT0A; } -float FITCutParHolder::thr2_FT0C() const { return mThr2FT0C; } \ No newline at end of file +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. +// +// \FIT bit thresholds +// \author Sandor Lokos, sandor.lokos@cern.ch +// \since March 2026 + +#include "FITCutParHolder.h" + +// setter +void FITCutParHolder::SetSaveFITbitsets(bool saveFITbitsets) +{ + mSaveFITbitsets = saveFITbitsets; +} +void FITCutParHolder::SetThr1FV0A(float thr1_FV0A) +{ + mThr1FV0A = thr1_FV0A; +} +void FITCutParHolder::SetThr1FT0A(float thr1_FT0A) +{ + mThr1FT0A = thr1_FT0A; +} +void FITCutParHolder::SetThr1FT0C(float thr1_FT0C) +{ + mThr1FT0C = thr1_FT0C; +} +void FITCutParHolder::SetThr2FV0A(float thr2_FV0A) +{ + mThr2FV0A = thr2_FV0A; +} +void FITCutParHolder::SetThr2FT0A(float thr2_FT0A) +{ + mThr2FT0A = thr2_FT0A; +} +void FITCutParHolder::SetThr2FT0C(float thr2_FT0C) +{ + mThr2FT0C = thr2_FT0C; +} + +// getter +bool FITCutParHolder::saveFITbitsets() const { return mSaveFITbitsets; } +float FITCutParHolder::thr1_FV0A() const { return mThr1FV0A; } +float FITCutParHolder::thr1_FT0A() const { return mThr1FT0A; } +float FITCutParHolder::thr1_FT0C() const { return mThr1FT0C; } +float FITCutParHolder::thr2_FV0A() const { return mThr2FV0A; } +float FITCutParHolder::thr2_FT0A() const { return mThr2FT0A; } +float FITCutParHolder::thr2_FT0C() const { return mThr2FT0C; } diff --git a/PWGUD/Core/FITCutParHolder.h b/PWGUD/Core/FITCutParHolder.h index d19ff01c97f..e03c4161c27 100644 --- a/PWGUD/Core/FITCutParHolder.h +++ b/PWGUD/Core/FITCutParHolder.h @@ -1,60 +1,75 @@ -#ifndef PWGUD_CORE_FITCUTPARHOLDER_H_ -#define PWGUD_CORE_FITCUTPARHOLDER_H_ - -#include - -// object to hold customizable FIT bit thresholds -class FITCutParHolder -{ - public: - // constructor - FITCutParHolder(bool saveFITbitsets = true, - float thr1_FV0A = 8., - float thr1_FT0A = 8., - float thr1_FT0C = 8., - float thr2_FV0A = 20., - float thr2_FT0A = 20., - float thr2_FT0C = 20.) - : mSaveFITbitsets{saveFITbitsets}, - mThr1FV0A{thr1_FV0A}, - mThr1FT0A{thr1_FT0A}, - mThr1FT0C{thr1_FT0C}, - mThr2FV0A{thr2_FV0A}, - mThr2FT0A{thr2_FT0A}, - mThr2FT0C{thr2_FT0C} - { - } - - // setters - void SetSaveFITbitsets(bool); - void SetThr1FV0A(float); - void SetThr1FT0A(float); - void SetThr1FT0C(float); - void SetThr2FV0A(float); - void SetThr2FT0A(float); - void SetThr2FT0C(float); - - // getters - bool saveFITbitsets() const; - float thr1_FV0A() const; - float thr1_FT0A() const; - float thr1_FT0C() const; - float thr2_FV0A() const; - float thr2_FT0A() const; - float thr2_FT0C() const; - - private: - bool mSaveFITbitsets; - - float mThr1FV0A; - float mThr1FT0A; - float mThr1FT0C; - - float mThr2FV0A; - float mThr2FT0A; - float mThr2FT0C; - - ClassDefNV(FITCutParHolder, 1); -}; - -#endif // PWGUD_CORE_FITCUTPARHOLDER_H_ \ No newline at end of file +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. +// +// \FIT bit thresholds +// \author Sandor Lokos, sandor.lokos@cern.ch +// \since March 2026 + +#ifndef PWGUD_CORE_FITCUTPARHOLDER_H_ +#define PWGUD_CORE_FITCUTPARHOLDER_H_ + +#include + +// object to hold customizable FIT bit thresholds +class FITCutParHolder +{ + public: + // constructor + FITCutParHolder(bool saveFITbitsets = true, + float thr1_FV0A = 8., + float thr1_FT0A = 8., + float thr1_FT0C = 8., + float thr2_FV0A = 20., + float thr2_FT0A = 20., + float thr2_FT0C = 20.) + : mSaveFITbitsets{saveFITbitsets}, + mThr1FV0A{thr1_FV0A}, + mThr1FT0A{thr1_FT0A}, + mThr1FT0C{thr1_FT0C}, + mThr2FV0A{thr2_FV0A}, + mThr2FT0A{thr2_FT0A}, + mThr2FT0C{thr2_FT0C} + { + } + + // setters + void SetSaveFITbitsets(bool); + void SetThr1FV0A(float); + void SetThr1FT0A(float); + void SetThr1FT0C(float); + void SetThr2FV0A(float); + void SetThr2FT0A(float); + void SetThr2FT0C(float); + + // getters + bool saveFITbitsets() const; + float thr1_FV0A() const; + float thr1_FT0A() const; + float thr1_FT0C() const; + float thr2_FV0A() const; + float thr2_FT0A() const; + float thr2_FT0C() const; + + private: + bool mSaveFITbitsets; + + float mThr1FV0A; + float mThr1FT0A; + float mThr1FT0C; + + float mThr2FV0A; + float mThr2FT0A; + float mThr2FT0C; + + ClassDefNV(FITCutParHolder, 1); +}; + +#endif // PWGUD_CORE_FITCUTPARHOLDER_H_ diff --git a/PWGUD/Core/FITCutParHolderLinkDef.h b/PWGUD/Core/FITCutParHolderLinkDef.h index 70341fc2ee4..3194a6754b2 100644 --- a/PWGUD/Core/FITCutParHolderLinkDef.h +++ b/PWGUD/Core/FITCutParHolderLinkDef.h @@ -1,19 +1,24 @@ -// Copyright 2019-2020 CERN and copyright holders of ALICE O2. -// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. -// All rights not expressly granted are reserved. -// -// This software is distributed under the terms of the GNU General Public -// License v3 (GPL Version 3), copied verbatim in the file "COPYING". -// -// In applying this license CERN does not waive the privileges and immunities -// granted to it by virtue of its status as an Intergovernmental Organization -// or submit itself to any jurisdiction. -#ifndef PWGUD_CORE_FITCUTPARHOLDERLINKDEF_H_ -#define PWGUD_CORE_FITCUTPARHOLDERLINKDEF_H_ - -#pragma link off all globals; -#pragma link off all classes; -#pragma link off all functions; -#pragma link C++ class FITCutParHolder + ; - -#endif // PWGUD_CORE_FITCUTPARHOLDERLINKDEF_H_ \ No newline at end of file +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. +// +// \FIT bit thresholds +// \author Sandor Lokos, sandor.lokos@cern.ch +// \since March 2026 + +#ifndef PWGUD_CORE_FITCUTPARHOLDERLINKDEF_H_ +#define PWGUD_CORE_FITCUTPARHOLDERLINKDEF_H_ + +#pragma link off all globals; +#pragma link off all classes; +#pragma link off all functions; +#pragma link C++ class FITCutParHolder + ; + +#endif // PWGUD_CORE_FITCUTPARHOLDERLINKDEF_H_ diff --git a/PWGUD/Core/SGCutParHolderLinkDef.h b/PWGUD/Core/SGCutParHolderLinkDef.h index f9216fa0c37..5dcadadc69e 100644 --- a/PWGUD/Core/SGCutParHolderLinkDef.h +++ b/PWGUD/Core/SGCutParHolderLinkDef.h @@ -8,6 +8,7 @@ // In applying this license CERN does not waive the privileges and immunities // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. + #ifndef PWGUD_CORE_SGCUTPARHOLDERLINKDEF_H_ #define PWGUD_CORE_SGCUTPARHOLDERLINKDEF_H_ diff --git a/PWGUD/Core/UDHelpers.h b/PWGUD/Core/UDHelpers.h index cd0d10e2626..ea137fbc17a 100644 --- a/PWGUD/Core/UDHelpers.h +++ b/PWGUD/Core/UDHelpers.h @@ -19,6 +19,7 @@ #include "PWGUD/Core/DGCutparHolder.h" #include "PWGUD/Core/UPCHelpers.h" +#include "Common/Core/RecoDecay.h" #include "Common/DataModel/EventSelection.h" #include "Common/DataModel/PIDResponse.h" #include "Common/DataModel/TrackSelectionTables.h" @@ -27,7 +28,6 @@ #include "DataFormatsFIT/Triggers.h" #include "DataFormatsFT0/Digit.h" #include "Framework/Logger.h" -#include "Common/Core/RecoDecay.h" #include "TLorentzVector.h" @@ -583,18 +583,16 @@ inline void buildFT0FV0Words(TFT0 const& ft0, TFV0A const& fv0a, } } - - // ----------------------------------------------------------------------------- // return eta and phi of a given FIT channel based on the bitset // Bit layout contract: -constexpr int kFT0Bits = 208; // FT0 total channels -constexpr int kFV0Bits = 48; // FV0A channels -constexpr int kTotalBits = 256; // 4*64 +constexpr int kFT0Bits = 208; // FT0 total channels +constexpr int kFV0Bits = 48; // FV0A channels +constexpr int kTotalBits = 256; // 4*64 // FT0 side split -constexpr int kFT0AChannels = 96; // FT0A channels are [0..95] -constexpr int kFT0CChannels = 112; // FT0C channels are [96..207] +constexpr int kFT0AChannels = 96; // FT0A channels are [0..95] +constexpr int kFT0CChannels = 112; // FT0C channels are [96..207] static_assert(kFT0AChannels + kFT0CChannels == kFT0Bits); using Bits256 = std::array; @@ -606,14 +604,16 @@ inline Bits256 makeBits256(uint64_t w0, uint64_t w1, uint64_t w2, uint64_t w3) inline bool testBit(Bits256 const& w, int bit) { - if (bit < 0 || bit >= kTotalBits) { + if (bit < 0 || bit >= kTotalBits) { return false; } return (w[bit >> 6] >> (bit & 63)) & 1ULL; } struct FitBitRef { - enum class Det : uint8_t { FT0, FV0, Unknown }; + enum class Det : uint8_t { FT0, + FV0, + Unknown }; Det det = Det::Unknown; int ch = -1; // FT0: 0..207, FV0: 0..47 bool isC = false; // only meaningful for FT0 @@ -624,13 +624,13 @@ inline FitBitRef decodeFitBit(int bit) FitBitRef out; if (bit >= 0 && bit < kFT0Bits) { out.det = FitBitRef::Det::FT0; - out.ch = bit; // FT0 channel id - out.isC = (bit >= kFT0AChannels); // C side if in upper range + out.ch = bit; // FT0 channel id + out.isC = (bit >= kFT0AChannels); // C side if in upper range return out; } if (bit >= kFT0Bits && bit < kTotalBits) { out.det = FitBitRef::Det::FV0; - out.ch = bit - kFT0Bits; // FV0A channel id 0..47 + out.ch = bit - kFT0Bits; // FV0A channel id 0..47 return out; } return out; @@ -677,20 +677,18 @@ inline bool getPhiEtaFromFitBit(FT0DetT& ft0Det, double& phi, double& eta) { - auto ref = decodeFitBit(bit); + auto ref = decodeFitBit(bit); if (ref.det != FitBitRef::Det::FT0) { return false; } - + // FT0A: 0..95, FT0C: 96..207 const int ft0Ch = bit; phi = getPhiFT0_fromChannel(ft0Det, ft0Ch, offsetFT0, iRunOffset); eta = getEtaFT0_fromChannel(ft0Det, ft0Ch, offsetFT0, iRunOffset); return true; - } - // ----------------------------------------------------------------------------- template diff --git a/PWGUD/DataModel/UDTables.h b/PWGUD/DataModel/UDTables.h index b6978027df7..6789342e81a 100644 --- a/PWGUD/DataModel/UDTables.h +++ b/PWGUD/DataModel/UDTables.h @@ -408,15 +408,15 @@ DECLARE_SOA_TABLE(UDTracksLabels, "AOD", "UDTRACKLABEL", namespace udcollfitbits { - DECLARE_SOA_COLUMN(Thr1W0, thr1W0, uint64_t); /// 1 MIP thresholds for FT0A ch 0 - ch 63 - DECLARE_SOA_COLUMN(Thr1W1, thr1W1, uint64_t); /// 1 MIP thresholds for FT0A ch 64 - ch 96 & FT0C ch 0 - ch 31 - DECLARE_SOA_COLUMN(Thr1W2, thr1W2, uint64_t); /// 1 MIP thresholds for FT0C ch 32 - ch 96 - DECLARE_SOA_COLUMN(Thr1W3, thr1W3, uint64_t); /// 1 MIP thresholds for FT0C ch 97 - 112 & FV0 0 - 47 - - DECLARE_SOA_COLUMN(Thr2W0, thr2W0, uint64_t); /// 2 MIP thresholds for FT0A ch 0 - ch 63 - DECLARE_SOA_COLUMN(Thr2W1, thr2W1, uint64_t); /// 2 MIP thresholds for FT0A ch 63 - ch 96 & FT0C ch 0 - ch 31 - DECLARE_SOA_COLUMN(Thr2W2, thr2W2, uint64_t); /// 2 MIP thresholds for FT0C ch 32 - ch 96 - DECLARE_SOA_COLUMN(Thr2W3, thr2W3, uint64_t); /// 2 MIP thresholds for FT0C ch 97 - 112 & FV0 0 - 47 +DECLARE_SOA_COLUMN(Thr1W0, thr1W0, uint64_t); /// 1 MIP thresholds for FT0A ch 0 - ch 63 +DECLARE_SOA_COLUMN(Thr1W1, thr1W1, uint64_t); /// 1 MIP thresholds for FT0A ch 64 - ch 96 & FT0C ch 0 - ch 31 +DECLARE_SOA_COLUMN(Thr1W2, thr1W2, uint64_t); /// 1 MIP thresholds for FT0C ch 32 - ch 96 +DECLARE_SOA_COLUMN(Thr1W3, thr1W3, uint64_t); /// 1 MIP thresholds for FT0C ch 97 - 112 & FV0 0 - 47 + +DECLARE_SOA_COLUMN(Thr2W0, thr2W0, uint64_t); /// 2 MIP thresholds for FT0A ch 0 - ch 63 +DECLARE_SOA_COLUMN(Thr2W1, thr2W1, uint64_t); /// 2 MIP thresholds for FT0A ch 63 - ch 96 & FT0C ch 0 - ch 31 +DECLARE_SOA_COLUMN(Thr2W2, thr2W2, uint64_t); /// 2 MIP thresholds for FT0C ch 32 - ch 96 +DECLARE_SOA_COLUMN(Thr2W3, thr2W3, uint64_t); /// 2 MIP thresholds for FT0C ch 97 - 112 & FV0 0 - 47 } // namespace udcollfitbits DECLARE_SOA_TABLE(UDCollisionFITBits, "AOD", "UDCOLLFITBITS", @@ -429,9 +429,7 @@ DECLARE_SOA_TABLE(UDCollisionFITBits, "AOD", "UDCOLLFITBITS", udcollfitbits::Thr2W1, /// 2 MIP thresholds for FT0A ch 63 - ch 96 & FT0C ch 0 - ch 31 udcollfitbits::Thr2W2, /// 2 MIP thresholds for FT0C ch 32 - ch 96 udcollfitbits::Thr2W3 /// 2 MIP thresholds for FT0C ch 97 - 112 & FV0 0 - 47 - ); - - +); using UDTrack = UDTracks::iterator; using UDTrackCov = UDTracksCov::iterator; diff --git a/PWGUD/TableProducer/SGCandProducer.cxx b/PWGUD/TableProducer/SGCandProducer.cxx index 25629e475e4..73e8531557e 100644 --- a/PWGUD/TableProducer/SGCandProducer.cxx +++ b/PWGUD/TableProducer/SGCandProducer.cxx @@ -12,15 +12,15 @@ /// \file SGCandProducer.cxx /// \brief Produces PWGUD derived table from standard tables /// -/// \author Alexander Bylinkin , Uniersity of Bergen +/// \author Alexander Bylinkin , University of Bergen /// \since 23.11.2023 /// \author Adam Matyja , INP PAN Krakow, Poland /// \since May 2025 // +#include "PWGUD/Core/FITCutParHolder.h" #include "PWGUD/Core/SGSelector.h" #include "PWGUD/Core/UPCHelpers.h" -#include "PWGUD/Core/FITCutParHolder.h" #include "PWGUD/DataModel/UDTables.h" #include "Common/CCDB/EventSelectionParams.h" @@ -42,11 +42,11 @@ #include "Framework/runDataProcessing.h" #include "ReconstructionDataFormats/Vertex.h" +#include #include #include #include #include -#include using namespace o2; using namespace o2::framework; @@ -325,23 +325,23 @@ struct SGCandProducer { fitInfo.BBFDDApf, fitInfo.BBFDDCpf, fitInfo.BGFDDApf, fitInfo.BGFDDCpf); outputCollisionSelExtras(chFT0A, chFT0C, chFDDA, chFDDC, chFV0A, occ, ir, trs, trofs, hmpr, tfb, itsROFb, sbp, zVtxFT0vPv, vtxITSTPC, collision.rct_raw()); outputCollsLabels(collision.globalIndex()); - + uint64_t w1[4] = {0ull, 0ull, 0ull, 0ull}; uint64_t w2[4] = {0ull, 0ull, 0ull, 0ull}; if (fitCuts.saveFITbitsets() && newbc.has_foundFT0() && newbc.has_fv0a()) { udhelpers::buildFT0FV0Words(newbc.ft0(), newbc.fv0a(), w1, w2, - fitCuts.thr1_FT0A(), - fitCuts.thr1_FT0C(), - fitCuts.thr1_FV0A(), - fitCuts.thr2_FT0A(), - fitCuts.thr2_FT0C(), - fitCuts.thr2_FV0A()); - } - - outputFITBits(w1[0], w1[1], w1[2], w1[3], - w2[0], w2[1], w2[2], w2[3]); - + fitCuts.thr1_FT0A(), + fitCuts.thr1_FT0C(), + fitCuts.thr1_FV0A(), + fitCuts.thr2_FT0A(), + fitCuts.thr2_FT0C(), + fitCuts.thr2_FV0A()); + } + + outputFITBits(w1[0], w1[1], w1[2], w1[3], + w2[0], w2[1], w2[2], w2[3]); + if (newbc.has_zdc()) { auto zdc = newbc.zdc(); udZdcsReduced(outputCollisions.lastIndex(), zdc.timeZNA(), zdc.timeZNC(), zdc.energyCommonZNA(), zdc.energyCommonZNC()); diff --git a/PWGUD/Tasks/upcTestFITBitMapping.cxx b/PWGUD/Tasks/upcTestFITBitMapping.cxx index 8469b9c127b..8456c46dc85 100644 --- a/PWGUD/Tasks/upcTestFITBitMapping.cxx +++ b/PWGUD/Tasks/upcTestFITBitMapping.cxx @@ -1,108 +1,113 @@ -// UpcTestFITBitMapping.cxx -// Minimal task to test FT0 bitset -> (eta,phi) mapping using udhelpers::getPhiEtaFromFitBit -// -// Run example (from O2Physics build dir): -// ./stage/bin/o2-pwgud-ud-fitbit-mapping --aod-file AO2D_UD.root -b -// -// Note: Add this file to PWGUD/Tasks/CMakeLists.txt. - -#include "Framework/AnalysisTask.h" -#include "Framework/runDataProcessing.h" -#include "Framework/HistogramRegistry.h" - -#include "PWGUD/DataModel/UDTables.h" // aod::UDCollisionFITBits -#include "PWGUD/Core/UDHelpers.h" // udhelpers::Bits256, makeBits256, testBit, getPhiEtaFromFitBit -#include "FT0Base/Geometry.h" // o2::ft0::Geometry - -#include -#include - -using namespace o2; -using namespace o2::framework; - -struct UpcTestFITBitMapping -{ - Configurable whichThr{"whichThr", 1, "Use 1=Thr1 bits or 2=Thr2 bits"}; - Configurable maxEvents{"maxEvents", -1, "Process at most this many rows (-1 = all)"}; - - // Minimal offset container compatible with UDHelpers.h expectations: getX/getY/getZ - struct OffsetXYZ { - double x{0}, y{0}, z{0}; - double getX() const { return x; } - double getY() const { return y; } - double getZ() const { return z; } - }; - - std::array offsetFT0{}; // iRunOffset = 0 for now - int iRunOffset = 0; - - o2::ft0::Geometry ft0Det{}; - - HistogramRegistry registry{ - "registry", - { - {"hPhiA", "FT0A #varphi;#varphi;counts", {HistType::kTH1F, {{18, 0.0, 2.*M_PI}}}}, - {"hEtaA", "FT0A #eta;#eta;counts", {HistType::kTH1F, {{8, 3.5, 5.0}}}}, - {"hEtaPhiA", "FT0A #eta vs #varphi;#eta;#varphi", {HistType::kTH2F, {{8, 3.5, 5.0}, {18, 0.0, 2.*M_PI}}}}, - - {"hPhiC", "FT0C #varphi;#varphi;counts", {HistType::kTH1F, {{18, 0.0, 2.*M_PI}}}}, - {"hEtaC", "FT0C #eta;#eta;counts", {HistType::kTH1F, {{8, -3.5, -2.0}}}}, - {"hEtaPhiC", "FT0C #eta vs #varphi;#eta;#varphi", {HistType::kTH2F, {{8, -3.5, -2.0}, {18, 0.0, 2.*M_PI}}}}, - }}; - - void init(InitContext&) - { - // UDHelpers calls calculateChannelCenter() inside, but doing it once here is fine. - ft0Det.calculateChannelCenter(); - } - - void process(aod::UDCollisionFITBits const& bitsTable) - { - int64_t nProcessed = 0; - - for (auto const& row : bitsTable) { - if (maxEvents >= 0 && nProcessed >= maxEvents) { - break; - } - ++nProcessed; - - // Use udhelpers' canonical packed type + builder - udhelpers::Bits256 w{}; - if (whichThr == 2) { - w = udhelpers::makeBits256(row.thr2W0(), row.thr2W1(), row.thr2W2(), row.thr2W3()); - } else { - w = udhelpers::makeBits256(row.thr1W0(), row.thr1W1(), row.thr1W2(), row.thr1W3()); - } - - // Loop FT0 bits only (0..207). FV0 starts at 208 but ignored here. - for (int bit = 0; bit < udhelpers::kFT0Bits; ++bit) { - if (!udhelpers::testBit(w, bit)) { - continue; - } - - double phi = 0., eta = 0.; - const bool ok = udhelpers::getPhiEtaFromFitBit(ft0Det, bit, offsetFT0, iRunOffset, phi, eta); - if (!ok) { - continue; - } - - if (bit < udhelpers::kFT0AChannels) { - registry.fill(HIST("hPhiA"), phi); - registry.fill(HIST("hEtaA"), eta); - registry.fill(HIST("hEtaPhiA"), eta, phi); - } else { - registry.fill(HIST("hPhiC"), phi); - registry.fill(HIST("hEtaC"), eta); - registry.fill(HIST("hEtaPhiC"), eta, phi); - } - } - } - } -}; - -WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) -{ - return WorkflowSpec{ - adaptAnalysisTask(cfgc, TaskName{"fitbit-mapping"}) - }; -} +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. +// +// \FIT bits to phi, eta mapping +// \author Sandor Lokos, sandor.lokos@cern.ch +// \since March 2026 + +#include "PWGUD/Core/UDHelpers.h" // udhelpers::Bits256, makeBits256, testBit, getPhiEtaFromFitBit +#include "PWGUD/DataModel/UDTables.h" // aod::UDCollisionFITBits + +#include "FT0Base/Geometry.h" // o2::ft0::Geometry +#include "Framework/AnalysisTask.h" +#include "Framework/HistogramRegistry.h" +#include "Framework/runDataProcessing.h" + +#include +#include + +using namespace o2; +using namespace o2::framework; + +struct UpcTestFITBitMapping { + Configurable whichThr{"whichThr", 1, "Use 1=Thr1 bits or 2=Thr2 bits"}; + Configurable maxEvents{"maxEvents", -1, "Process at most this many rows (-1 = all)"}; + + // Minimal offset container compatible with UDHelpers.h expectations: getX/getY/getZ + struct OffsetXYZ { + double x{0}, y{0}, z{0}; + double getX() const { return x; } + double getY() const { return y; } + double getZ() const { return z; } + }; + + std::array offsetFT0{}; // iRunOffset = 0 for now + int iRunOffset = 0; + + o2::ft0::Geometry ft0Det{}; + + HistogramRegistry registry{ + "registry", + { + {"hPhiA", "FT0A #varphi;#varphi;counts", {HistType::kTH1F, {{18, 0.0, 2. * M_PI}}}}, + {"hEtaA", "FT0A #eta;#eta;counts", {HistType::kTH1F, {{8, 3.5, 5.0}}}}, + {"hEtaPhiA", "FT0A #eta vs #varphi;#eta;#varphi", {HistType::kTH2F, {{8, 3.5, 5.0}, {18, 0.0, 2. * M_PI}}}}, + + {"hPhiC", "FT0C #varphi;#varphi;counts", {HistType::kTH1F, {{18, 0.0, 2. * M_PI}}}}, + {"hEtaC", "FT0C #eta;#eta;counts", {HistType::kTH1F, {{8, -3.5, -2.0}}}}, + {"hEtaPhiC", "FT0C #eta vs #varphi;#eta;#varphi", {HistType::kTH2F, {{8, -3.5, -2.0}, {18, 0.0, 2. * M_PI}}}}, + }}; + + void init(InitContext&) + { + // UDHelpers calls calculateChannelCenter() inside, but doing it once here is fine. + ft0Det.calculateChannelCenter(); + } + + void process(aod::UDCollisionFITBits const& bitsTable) + { + int64_t nProcessed = 0; + + for (auto const& row : bitsTable) { + if (maxEvents >= 0 && nProcessed >= maxEvents) { + break; + } + ++nProcessed; + + // Use udhelpers' canonical packed type + builder + udhelpers::Bits256 w{}; + if (whichThr == 2) { + w = udhelpers::makeBits256(row.thr2W0(), row.thr2W1(), row.thr2W2(), row.thr2W3()); + } else { + w = udhelpers::makeBits256(row.thr1W0(), row.thr1W1(), row.thr1W2(), row.thr1W3()); + } + + // Loop FT0 bits only (0..207). FV0 starts at 208 but ignored here. + for (int bit = 0; bit < udhelpers::kFT0Bits; ++bit) { + if (!udhelpers::testBit(w, bit)) { + continue; + } + + double phi = 0., eta = 0.; + const bool ok = udhelpers::getPhiEtaFromFitBit(ft0Det, bit, offsetFT0, iRunOffset, phi, eta); + if (!ok) { + continue; + } + + if (bit < udhelpers::kFT0AChannels) { + registry.fill(HIST("hPhiA"), phi); + registry.fill(HIST("hEtaA"), eta); + registry.fill(HIST("hEtaPhiA"), eta, phi); + } else { + registry.fill(HIST("hPhiC"), phi); + registry.fill(HIST("hEtaC"), eta); + registry.fill(HIST("hEtaPhiC"), eta, phi); + } + } + } + } +}; + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + return WorkflowSpec{ + adaptAnalysisTask(cfgc, TaskName{"fitbit-mapping"})}; +} From 30d11c328a2103bba549e551e05c1e5925bea351 Mon Sep 17 00:00:00 2001 From: sandor-lokos Date: Tue, 24 Mar 2026 20:28:06 +0100 Subject: [PATCH 06/11] Formatting fixes --- PWGUD/TableProducer/SGCandProducer.cxx | 1542 ++++++++++++------------ 1 file changed, 771 insertions(+), 771 deletions(-) diff --git a/PWGUD/TableProducer/SGCandProducer.cxx b/PWGUD/TableProducer/SGCandProducer.cxx index 73e8531557e..80ec14a0b1f 100644 --- a/PWGUD/TableProducer/SGCandProducer.cxx +++ b/PWGUD/TableProducer/SGCandProducer.cxx @@ -13,9 +13,9 @@ /// \brief Produces PWGUD derived table from standard tables /// /// \author Alexander Bylinkin , University of Bergen -/// \since 23.11.2023 +/// \since 23.11.2023 /// \author Adam Matyja , INP PAN Krakow, Poland -/// \since May 2025 +/// \since May 2025 // #include "PWGUD/Core/FITCutParHolder.h" @@ -57,274 +57,274 @@ using namespace o2::aod::rctsel; #define getHist(type, name) std::get>(histPointers[name]) struct SGCandProducer { - Service ccdb; - // data inputs - using CCs = soa::Join; - using CC = CCs::iterator; - using BCs = soa::Join; - using BC = BCs::iterator; - using TCs = soa::Join; - using FWs = aod::FwdTracks; - - using MCCCs = soa::Join; - using MCCC = MCCCs::iterator; - - // get an SGCutparHolder - SGCutParHolder sameCuts = SGCutParHolder(); // SGCutparHolder - Configurable SGCuts{"SGCuts", {}, "SG event cuts"}; + Service ccdb; + // data inputs + using CCs = soa::Join; + using CC = CCs::iterator; + using BCs = soa::Join; + using BC = BCs::iterator; + using TCs = soa::Join; + using FWs = aod::FwdTracks; + + using MCCCs = soa::Join; + using MCCC = MCCCs::iterator; + + // get an SGCutparHolder + SGCutParHolder sameCuts = SGCutParHolder(); // SGCutparHolder + Configurable SGCuts{"SGCuts", {}, "SG event cuts"}; FITCutParHolder fitCuts = FITCutParHolder(); Configurable FITCuts{"FITCuts", {}, "FIT bitset cuts"}; - Configurable verboseInfo{"verboseInfo", false, "Print general info to terminal; default it false."}; - Configurable saveAllTracks{"saveAllTracks", true, "save only PV contributors or all tracks associated to a collision"}; - Configurable savenonPVCITSOnlyTracks{"savenonPVCITSOnlyTracks", false, "save non PV contributors with ITS only information"}; - Configurable rejectAtTFBoundary{"rejectAtTFBoundary", true, "reject collisions at a TF boundary"}; - Configurable noITSROFrameBorder{"noITSROFrameBorder", true, "reject ITS RO Frame Border"}; - Configurable noSameBunchPileUp{"noSameBunchPileUp", true, "reject SameBunchPileUp"}; - Configurable IsGoodVertex{"IsGoodVertex", false, "Select FT0 PV vertex matching"}; - Configurable ITSTPCVertex{"ITSTPCVertex", true, "reject ITS-only vertex"}; // if one wants to look at Single Gap pp events - Configurable> generatorIds{"generatorIds", std::vector{-1}, "MC generatorIds to process"}; - Configurable storeSG{"storeSG", true, "Store SG events in the output"}; - Configurable storeDG{"storeDG", true, "Store DG events in the output"}; - - Configurable saveFITbitsets{"saveFITbitsets", true, "Write FT0 and FV0 bitset tables to output"}; - - Configurable isGoodRCTCollision{"isGoodRCTCollision", true, "Check RCT flags for FT0,ITS,TPC and tracking"}; - Configurable isGoodRCTZdc{"isGoodRCTZdc", false, "Check RCT flags for ZDC if present in run"}; - - // Configurables to decide which tables are filled - Configurable fillTrackTables{"fillTrackTables", true, "Fill track tables"}; - Configurable fillFwdTrackTables{"fillFwdTrackTables", true, "Fill forward track tables"}; - - // SG selector - SGSelector sgSelector; - ctpRateFetcher mRateFetcher; - - // initialize RCT flag checker - RCTFlagsChecker myRCTChecker{"CBT"}; - - // data tables - Produces outputSGCollisions; - Produces outputCollisions; - Produces outputCollisionsSels; - Produces outputCollisionSelExtras; - Produces outputCollsLabels; - Produces outputZdcs; - Produces udZdcsReduced; - Produces outputTracks; - Produces outputTracksCov; - Produces outputTracksDCA; - Produces outputTracksPID; - Produces outputTracksPIDExtra; - Produces outputTracksExtra; - Produces outputTracksFlag; - Produces outputFwdTracks; - Produces outputFwdTracksExtra; - Produces outputTracksLabel; + Configurable verboseInfo{"verboseInfo", false, "Print general info to terminal; default it false."}; + Configurable saveAllTracks{"saveAllTracks", true, "save only PV contributors or all tracks associated to a collision"}; + Configurable savenonPVCITSOnlyTracks{"savenonPVCITSOnlyTracks", false, "save non PV contributors with ITS only information"}; + Configurable rejectAtTFBoundary{"rejectAtTFBoundary", true, "reject collisions at a TF boundary"}; + Configurable noITSROFrameBorder{"noITSROFrameBorder", true, "reject ITS RO Frame Border"}; + Configurable noSameBunchPileUp{"noSameBunchPileUp", true, "reject SameBunchPileUp"}; + Configurable IsGoodVertex{"IsGoodVertex", false, "Select FT0 PV vertex matching"}; + Configurable ITSTPCVertex{"ITSTPCVertex", true, "reject ITS-only vertex"}; // if one wants to look at Single Gap pp events + Configurable> generatorIds{"generatorIds", std::vector{-1}, "MC generatorIds to process"}; + Configurable storeSG{"storeSG", true, "Store SG events in the output"}; + Configurable storeDG{"storeDG", true, "Store DG events in the output"}; + + Configurable saveFITbitsets{"saveFITbitsets", true, "Write FT0 and FV0 bitset tables to output"}; + + Configurable isGoodRCTCollision{"isGoodRCTCollision", true, "Check RCT flags for FT0,ITS,TPC and tracking"}; + Configurable isGoodRCTZdc{"isGoodRCTZdc", false, "Check RCT flags for ZDC if present in run"}; + + // Configurables to decide which tables are filled + Configurable fillTrackTables{"fillTrackTables", true, "Fill track tables"}; + Configurable fillFwdTrackTables{"fillFwdTrackTables", true, "Fill forward track tables"}; + + // SG selector + SGSelector sgSelector; + ctpRateFetcher mRateFetcher; + + // initialize RCT flag checker + RCTFlagsChecker myRCTChecker{"CBT"}; + + // data tables + Produces outputSGCollisions; + Produces outputCollisions; + Produces outputCollisionsSels; + Produces outputCollisionSelExtras; + Produces outputCollsLabels; + Produces outputZdcs; + Produces udZdcsReduced; + Produces outputTracks; + Produces outputTracksCov; + Produces outputTracksDCA; + Produces outputTracksPID; + Produces outputTracksPIDExtra; + Produces outputTracksExtra; + Produces outputTracksFlag; + Produces outputFwdTracks; + Produces outputFwdTracksExtra; + Produces outputTracksLabel; Produces outputFITBits; - // initialize histogram registry - HistogramRegistry registry{ - "registry", - {}}; - std::map histPointers; - - int runNumber = -1; - - // function to update UDFwdTracks, UDFwdTracksExtra - template - void updateUDFwdTrackTables(TFwdTrack const& fwdtrack, uint64_t const& bcnum) - { - outputFwdTracks(outputCollisions.lastIndex(), - fwdtrack.px(), fwdtrack.py(), fwdtrack.pz(), fwdtrack.sign(), - bcnum, fwdtrack.trackTime(), fwdtrack.trackTimeRes()); - outputFwdTracksExtra(fwdtrack.trackType(), - fwdtrack.nClusters(), - fwdtrack.pDca(), - fwdtrack.rAtAbsorberEnd(), - fwdtrack.chi2(), - fwdtrack.chi2MatchMCHMID(), - fwdtrack.chi2MatchMCHMFT(), - fwdtrack.mchBitMap(), - fwdtrack.midBitMap(), - fwdtrack.midBoards()); - } - - // function to update UDTracks, UDTracksCov, UDTracksDCA, UDTracksPID, UDTracksExtra, UDTracksFlag, - // and UDTrackCollisionIDs - template - void updateUDTrackTables(int64_t lastIndex, TTrack const& track, uint64_t const& bcnum) - { - outputTracks(lastIndex, - track.px(), track.py(), track.pz(), track.sign(), - bcnum, track.trackTime(), track.trackTimeRes()); - - // float sigmaY = track.sigmaY(); - // float sigmaZ = track.sigmaZ(); - float sigmaY = -1.; - float sigmaZ = -1.; - outputTracksCov(track.x(), track.y(), track.z(), sigmaY, sigmaZ); - - outputTracksDCA(track.dcaZ(), track.dcaXY()); - outputTracksPID(track.tpcNSigmaEl(), - track.tpcNSigmaMu(), - track.tpcNSigmaPi(), - track.tpcNSigmaKa(), - track.tpcNSigmaPr(), - track.beta(), - track.betaerror(), - track.tofNSigmaEl(), - track.tofNSigmaMu(), - track.tofNSigmaPi(), - track.tofNSigmaKa(), - track.tofNSigmaPr()); - outputTracksPIDExtra(track.tpcNSigmaDe(), - track.tpcNSigmaTr(), - track.tpcNSigmaHe(), - track.tpcNSigmaAl(), - track.tofNSigmaDe(), - track.tofNSigmaTr(), - track.tofNSigmaHe(), - track.tofNSigmaAl()); - outputTracksExtra(track.tpcInnerParam(), - track.itsClusterSizes(), - track.tpcNClsFindable(), - track.tpcNClsFindableMinusFound(), - track.tpcNClsFindableMinusCrossedRows(), - track.tpcNClsShared(), - track.trdPattern(), - track.itsChi2NCl(), - track.tpcChi2NCl(), - track.trdChi2(), - track.tofChi2(), - track.tpcSignal(), - track.tofSignal(), - track.trdSignal(), - track.length(), - track.tofExpMom(), - track.detectorMap()); - outputTracksFlag(track.has_collision(), - track.isPVContributor()); - outputTracksLabel(track.globalIndex()); - } - - // function to process reconstructed data - template - void processReco(std::string histdir, TCol const& collision, BCs const& bcs, - TCs const& tracks, FWs const& fwdtracks, - aod::FV0As const& fv0as, aod::FT0s const& ft0s, aod::FDDs const& fdds) - { - if (verboseInfo) - LOGF(debug, " collision %d", collision.globalIndex()); - getHist(TH1, histdir + "/Stat")->Fill(0., 1.); - // reject collisions at TF boundaries - if (rejectAtTFBoundary && !collision.selection_bit(aod::evsel::kNoTimeFrameBorder)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(1., 1.); - // reject collisions at ITS RO TF boundaries - if (noITSROFrameBorder && !collision.selection_bit(aod::evsel::kNoITSROFrameBorder)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(2., 1.); - // reject Same Bunch PileUp - if (noSameBunchPileUp && !collision.selection_bit(aod::evsel::kNoSameBunchPileup)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(3., 1.); - // check vertex matching to FT0 - if (IsGoodVertex && !collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(4., 1.); - // reject ITS Only vertices - if (ITSTPCVertex && !collision.selection_bit(aod::evsel::kIsVertexITSTPC)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(5., 1.); - // nominal BC - if (!collision.has_foundBC()) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(6., 1.); - // RCT CBT for collision check - if (isGoodRCTCollision && !myRCTChecker(collision)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(7., 1.); - // RCT CBT+ZDC for collision check - if (isGoodRCTZdc && !myRCTChecker(collision)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(8., 1.); - - // - const int trs = collision.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard) ? 1 : 0; - const int trofs = collision.selection_bit(o2::aod::evsel::kNoCollInRofStandard) ? 1 : 0; - const int hmpr = collision.selection_bit(o2::aod::evsel::kNoHighMultCollInPrevRof) ? 1 : 0; - const int tfb = collision.selection_bit(o2::aod::evsel::kNoTimeFrameBorder) ? 1 : 0; - const int itsROFb = collision.selection_bit(o2::aod::evsel::kNoITSROFrameBorder) ? 1 : 0; - const int sbp = collision.selection_bit(o2::aod::evsel::kNoSameBunchPileup) ? 1 : 0; - const int zVtxFT0vPv = collision.selection_bit(o2::aod::evsel::kIsGoodZvtxFT0vsPV) ? 1 : 0; - const int vtxITSTPC = collision.selection_bit(o2::aod::evsel::kIsVertexITSTPC) ? 1 : 0; - auto bc = collision.template foundBC_as(); - double ir = 0.; - const uint64_t ts = bc.timestamp(); - const int runnumber = bc.runNumber(); - if (bc.has_zdc()) { - ir = mRateFetcher.fetch(ccdb.service, ts, runnumber, "ZNC hadronic") * 1.e-3; - } - auto newbc = bc; - - // obtain slice of compatible BCs - auto bcRange = udhelpers::compatibleBCs(collision, sameCuts.NDtcoll(), bcs, sameCuts.minNBCs()); - auto isSGEvent = sgSelector.IsSelected(sameCuts, collision, bcRange, bc); - // auto isSGEvent = sgSelector.IsSelected(sameCuts, collision, bcRange, tracks); - int issgevent = isSGEvent.value; - if (isSGEvent.bc && issgevent < 2) { - newbc = *(isSGEvent.bc); - } else { - if (verboseInfo) - LOGF(info, "No Newbc %i", bc.globalBC()); - } - getHist(TH1, histdir + "/Stat")->Fill(issgevent + 10, 1.); - if ((storeDG && issgevent == o2::aod::sgselector::DoubleGap) || (storeSG && (issgevent == o2::aod::sgselector::SingleGapA || issgevent == o2::aod::sgselector::SingleGapC))) { - if (verboseInfo) - LOGF(info, "Current BC: %i, %i, %i", bc.globalBC(), newbc.globalBC(), issgevent); - if (sameCuts.minRgtrwTOF()) { - if (udhelpers::rPVtrwTOF(tracks, collision.numContrib()) < sameCuts.minRgtrwTOF()) - return; - } - upchelpers::FITInfo fitInfo{}; - const uint8_t chFT0A = 0; - const uint8_t chFT0C = 0; - const uint8_t chFDDA = 0; - const uint8_t chFDDC = 0; - const uint8_t chFV0A = 0; - const int occ = collision.trackOccupancyInTimeRange(); - udhelpers::getFITinfo(fitInfo, newbc, bcs, ft0s, fv0as, fdds); - const int upc_flag = (collision.flags() & dataformats::Vertex>::Flags::UPCMode) ? 1 : 0; - // update SG candidates tables - outputCollisions(bc.globalBC(), bc.runNumber(), - collision.posX(), collision.posY(), collision.posZ(), upc_flag, - collision.numContrib(), udhelpers::netCharge(tracks), - 1.); // rtrwTOF); //omit the calculation to speed up the things while skimming - - outputSGCollisions(issgevent); - outputCollisionsSels(fitInfo.ampFT0A, fitInfo.ampFT0C, fitInfo.timeFT0A, fitInfo.timeFT0C, - fitInfo.triggerMaskFT0, - fitInfo.ampFDDA, fitInfo.ampFDDC, fitInfo.timeFDDA, fitInfo.timeFDDC, - fitInfo.triggerMaskFDD, - fitInfo.ampFV0A, fitInfo.timeFV0A, fitInfo.triggerMaskFV0A, - fitInfo.BBFT0Apf, fitInfo.BBFT0Cpf, fitInfo.BGFT0Apf, fitInfo.BGFT0Cpf, - fitInfo.BBFV0Apf, fitInfo.BGFV0Apf, - fitInfo.BBFDDApf, fitInfo.BBFDDCpf, fitInfo.BGFDDApf, fitInfo.BGFDDCpf); - outputCollisionSelExtras(chFT0A, chFT0C, chFDDA, chFDDC, chFV0A, occ, ir, trs, trofs, hmpr, tfb, itsROFb, sbp, zVtxFT0vPv, vtxITSTPC, collision.rct_raw()); - outputCollsLabels(collision.globalIndex()); + // initialize histogram registry + HistogramRegistry registry{ + "registry", + {}}; + std::map histPointers; + + int runNumber = -1; + + // function to update UDFwdTracks, UDFwdTracksExtra + template + void updateUDFwdTrackTables(TFwdTrack const& fwdtrack, uint64_t const& bcnum) + { + outputFwdTracks(outputCollisions.lastIndex(), + fwdtrack.px(), fwdtrack.py(), fwdtrack.pz(), fwdtrack.sign(), + bcnum, fwdtrack.trackTime(), fwdtrack.trackTimeRes()); + outputFwdTracksExtra(fwdtrack.trackType(), + fwdtrack.nClusters(), + fwdtrack.pDca(), + fwdtrack.rAtAbsorberEnd(), + fwdtrack.chi2(), + fwdtrack.chi2MatchMCHMID(), + fwdtrack.chi2MatchMCHMFT(), + fwdtrack.mchBitMap(), + fwdtrack.midBitMap(), + fwdtrack.midBoards()); + } + + // function to update UDTracks, UDTracksCov, UDTracksDCA, UDTracksPID, UDTracksExtra, UDTracksFlag, + // and UDTrackCollisionIDs + template + void updateUDTrackTables(int64_t lastIndex, TTrack const& track, uint64_t const& bcnum) + { + outputTracks(lastIndex, + track.px(), track.py(), track.pz(), track.sign(), + bcnum, track.trackTime(), track.trackTimeRes()); + + // float sigmaY = track.sigmaY(); + // float sigmaZ = track.sigmaZ(); + float sigmaY = -1.; + float sigmaZ = -1.; + outputTracksCov(track.x(), track.y(), track.z(), sigmaY, sigmaZ); + + outputTracksDCA(track.dcaZ(), track.dcaXY()); + outputTracksPID(track.tpcNSigmaEl(), + track.tpcNSigmaMu(), + track.tpcNSigmaPi(), + track.tpcNSigmaKa(), + track.tpcNSigmaPr(), + track.beta(), + track.betaerror(), + track.tofNSigmaEl(), + track.tofNSigmaMu(), + track.tofNSigmaPi(), + track.tofNSigmaKa(), + track.tofNSigmaPr()); + outputTracksPIDExtra(track.tpcNSigmaDe(), + track.tpcNSigmaTr(), + track.tpcNSigmaHe(), + track.tpcNSigmaAl(), + track.tofNSigmaDe(), + track.tofNSigmaTr(), + track.tofNSigmaHe(), + track.tofNSigmaAl()); + outputTracksExtra(track.tpcInnerParam(), + track.itsClusterSizes(), + track.tpcNClsFindable(), + track.tpcNClsFindableMinusFound(), + track.tpcNClsFindableMinusCrossedRows(), + track.tpcNClsShared(), + track.trdPattern(), + track.itsChi2NCl(), + track.tpcChi2NCl(), + track.trdChi2(), + track.tofChi2(), + track.tpcSignal(), + track.tofSignal(), + track.trdSignal(), + track.length(), + track.tofExpMom(), + track.detectorMap()); + outputTracksFlag(track.has_collision(), + track.isPVContributor()); + outputTracksLabel(track.globalIndex()); + } + + // function to process reconstructed data + template + void processReco(std::string histdir, TCol const& collision, BCs const& bcs, + TCs const& tracks, FWs const& fwdtracks, + aod::FV0As const& fv0as, aod::FT0s const& ft0s, aod::FDDs const& fdds) + { + if (verboseInfo) + LOGF(debug, " collision %d", collision.globalIndex()); + getHist(TH1, histdir + "/Stat")->Fill(0., 1.); + // reject collisions at TF boundaries + if (rejectAtTFBoundary && !collision.selection_bit(aod::evsel::kNoTimeFrameBorder)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(1., 1.); + // reject collisions at ITS RO TF boundaries + if (noITSROFrameBorder && !collision.selection_bit(aod::evsel::kNoITSROFrameBorder)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(2., 1.); + // reject Same Bunch PileUp + if (noSameBunchPileUp && !collision.selection_bit(aod::evsel::kNoSameBunchPileup)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(3., 1.); + // check vertex matching to FT0 + if (IsGoodVertex && !collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(4., 1.); + // reject ITS Only vertices + if (ITSTPCVertex && !collision.selection_bit(aod::evsel::kIsVertexITSTPC)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(5., 1.); + // nominal BC + if (!collision.has_foundBC()) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(6., 1.); + // RCT CBT for collision check + if (isGoodRCTCollision && !myRCTChecker(collision)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(7., 1.); + // RCT CBT+ZDC for collision check + if (isGoodRCTZdc && !myRCTChecker(collision)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(8., 1.); + + // + const int trs = collision.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard) ? 1 : 0; + const int trofs = collision.selection_bit(o2::aod::evsel::kNoCollInRofStandard) ? 1 : 0; + const int hmpr = collision.selection_bit(o2::aod::evsel::kNoHighMultCollInPrevRof) ? 1 : 0; + const int tfb = collision.selection_bit(o2::aod::evsel::kNoTimeFrameBorder) ? 1 : 0; + const int itsROFb = collision.selection_bit(o2::aod::evsel::kNoITSROFrameBorder) ? 1 : 0; + const int sbp = collision.selection_bit(o2::aod::evsel::kNoSameBunchPileup) ? 1 : 0; + const int zVtxFT0vPv = collision.selection_bit(o2::aod::evsel::kIsGoodZvtxFT0vsPV) ? 1 : 0; + const int vtxITSTPC = collision.selection_bit(o2::aod::evsel::kIsVertexITSTPC) ? 1 : 0; + auto bc = collision.template foundBC_as(); + double ir = 0.; + const uint64_t ts = bc.timestamp(); + const int runnumber = bc.runNumber(); + if (bc.has_zdc()) { + ir = mRateFetcher.fetch(ccdb.service, ts, runnumber, "ZNC hadronic") * 1.e-3; + } + auto newbc = bc; + + // obtain slice of compatible BCs + auto bcRange = udhelpers::compatibleBCs(collision, sameCuts.NDtcoll(), bcs, sameCuts.minNBCs()); + auto isSGEvent = sgSelector.IsSelected(sameCuts, collision, bcRange, bc); + // auto isSGEvent = sgSelector.IsSelected(sameCuts, collision, bcRange, tracks); + int issgevent = isSGEvent.value; + if (isSGEvent.bc && issgevent < 2) { + newbc = *(isSGEvent.bc); + } else { + if (verboseInfo) + LOGF(info, "No Newbc %i", bc.globalBC()); + } + getHist(TH1, histdir + "/Stat")->Fill(issgevent + 10, 1.); + if ((storeDG && issgevent == o2::aod::sgselector::DoubleGap) || (storeSG && (issgevent == o2::aod::sgselector::SingleGapA || issgevent == o2::aod::sgselector::SingleGapC))) { + if (verboseInfo) + LOGF(info, "Current BC: %i, %i, %i", bc.globalBC(), newbc.globalBC(), issgevent); + if (sameCuts.minRgtrwTOF()) { + if (udhelpers::rPVtrwTOF(tracks, collision.numContrib()) < sameCuts.minRgtrwTOF()) + return; + } + upchelpers::FITInfo fitInfo{}; + const uint8_t chFT0A = 0; + const uint8_t chFT0C = 0; + const uint8_t chFDDA = 0; + const uint8_t chFDDC = 0; + const uint8_t chFV0A = 0; + const int occ = collision.trackOccupancyInTimeRange(); + udhelpers::getFITinfo(fitInfo, newbc, bcs, ft0s, fv0as, fdds); + const int upc_flag = (collision.flags() & dataformats::Vertex>::Flags::UPCMode) ? 1 : 0; + // update SG candidates tables + outputCollisions(bc.globalBC(), bc.runNumber(), + collision.posX(), collision.posY(), collision.posZ(), upc_flag, + collision.numContrib(), udhelpers::netCharge(tracks), + 1.); // rtrwTOF); //omit the calculation to speed up the things while skimming + + outputSGCollisions(issgevent); + outputCollisionsSels(fitInfo.ampFT0A, fitInfo.ampFT0C, fitInfo.timeFT0A, fitInfo.timeFT0C, + fitInfo.triggerMaskFT0, + fitInfo.ampFDDA, fitInfo.ampFDDC, fitInfo.timeFDDA, fitInfo.timeFDDC, + fitInfo.triggerMaskFDD, + fitInfo.ampFV0A, fitInfo.timeFV0A, fitInfo.triggerMaskFV0A, + fitInfo.BBFT0Apf, fitInfo.BBFT0Cpf, fitInfo.BGFT0Apf, fitInfo.BGFT0Cpf, + fitInfo.BBFV0Apf, fitInfo.BGFV0Apf, + fitInfo.BBFDDApf, fitInfo.BBFDDCpf, fitInfo.BGFDDApf, fitInfo.BGFDDCpf); + outputCollisionSelExtras(chFT0A, chFT0C, chFDDA, chFDDC, chFV0A, occ, ir, trs, trofs, hmpr, tfb, itsROFb, sbp, zVtxFT0vPv, vtxITSTPC, collision.rct_raw()); + outputCollsLabels(collision.globalIndex()); uint64_t w1[4] = {0ull, 0ull, 0ull, 0ull}; uint64_t w2[4] = {0ull, 0ull, 0ull, 0ull}; @@ -342,518 +342,518 @@ struct SGCandProducer { outputFITBits(w1[0], w1[1], w1[2], w1[3], w2[0], w2[1], w2[2], w2[3]); - if (newbc.has_zdc()) { - auto zdc = newbc.zdc(); - udZdcsReduced(outputCollisions.lastIndex(), zdc.timeZNA(), zdc.timeZNC(), zdc.energyCommonZNA(), zdc.energyCommonZNC()); - } else { - udZdcsReduced(outputCollisions.lastIndex(), -999, -999, -999, -999); - } - // update SGTracks tables - if (fillTrackTables) { - for (const auto& track : tracks) { - if (track.pt() > sameCuts.minPt() && track.eta() > sameCuts.minEta() && track.eta() < sameCuts.maxEta()) { - if (track.isPVContributor()) { - updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); - } else if (saveAllTracks) { - if (track.itsClusterSizes() && track.itsChi2NCl() > 0 && ((track.tpcNClsFindable() == 0 && savenonPVCITSOnlyTracks) || track.tpcNClsFindable() > 50)) - updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); - // if (track.isPVContributor()) updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); - } - } - } - } - // update SGFwdTracks tables - if (fillFwdTrackTables) { - if (sameCuts.withFwdTracks()) { - for (const auto& fwdtrack : fwdtracks) { - if (!sgSelector.FwdTrkSelector(fwdtrack)) - updateUDFwdTrackTables(fwdtrack, bc.globalBC()); - } - } - } - } - } - - void init(InitContext& context) - { - ccdb->setURL("http://alice-ccdb.cern.ch"); - ccdb->setCaching(true); - ccdb->setFatalWhenNull(false); - sameCuts = (SGCutParHolder)SGCuts; + if (newbc.has_zdc()) { + auto zdc = newbc.zdc(); + udZdcsReduced(outputCollisions.lastIndex(), zdc.timeZNA(), zdc.timeZNC(), zdc.energyCommonZNA(), zdc.energyCommonZNC()); + } else { + udZdcsReduced(outputCollisions.lastIndex(), -999, -999, -999, -999); + } + // update SGTracks tables + if (fillTrackTables) { + for (const auto& track : tracks) { + if (track.pt() > sameCuts.minPt() && track.eta() > sameCuts.minEta() && track.eta() < sameCuts.maxEta()) { + if (track.isPVContributor()) { + updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); + } else if (saveAllTracks) { + if (track.itsClusterSizes() && track.itsChi2NCl() > 0 && ((track.tpcNClsFindable() == 0 && savenonPVCITSOnlyTracks) || track.tpcNClsFindable() > 50)) + updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); + // if (track.isPVContributor()) updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); + } + } + } + } + // update SGFwdTracks tables + if (fillFwdTrackTables) { + if (sameCuts.withFwdTracks()) { + for (const auto& fwdtrack : fwdtracks) { + if (!sgSelector.FwdTrkSelector(fwdtrack)) + updateUDFwdTrackTables(fwdtrack, bc.globalBC()); + } + } + } + } + } + + void init(InitContext& context) + { + ccdb->setURL("http://alice-ccdb.cern.ch"); + ccdb->setCaching(true); + ccdb->setFatalWhenNull(false); + sameCuts = (SGCutParHolder)SGCuts; fitCuts = (FITCutParHolder)FITCuts; - // add histograms for the different process functions - histPointers.clear(); - if (context.mOptions.get("processData")) { - histPointers.insert({"reco/Stat", registry.add("reco/Stat", "Cut statistics; Selection criterion; Collisions", {HistType::kTH1F, {{14, -0.5, 13.5}}})}); - - const AxisSpec axisCountersTrg{10, 0.5, 10.5, ""}; - histPointers.insert({"reco/hCountersTrg", registry.add("reco/hCountersTrg", "Trigger counts before selections; Trigger; Counts", {HistType::kTH1F, {axisCountersTrg}})}); - histPointers.insert({"reco/hCountersTrgBcSel", registry.add("reco/hCountersTrgSel", "Trigger counts after BC selections; Trigger; Counts", {HistType::kTH1F, {axisCountersTrg}})}); - histPointers.insert({"reco/hLumi", registry.add("reco/hLumi", "Integrated luminosity before selections; Trigger; Luminosity, 1/#mub", {HistType::kTH1F, {axisCountersTrg}})}); - histPointers.insert({"reco/hLumiBcSel", registry.add("reco/hLumiBcSel", "Integrated luminosity before selections; Trigger; Luminosity, 1/#mub", {HistType::kTH1F, {axisCountersTrg}})}); - auto hCountersTrg = getHist(TH1, "reco/hCountersTrg"); - auto hCountersTrgBcSel = getHist(TH1, "reco/hCountersTrgBcSel"); - auto hLumi = getHist(TH1, "reco/hLumi"); - auto hLumiBcSel = getHist(TH1, "reco/hLumiBcSel"); - for (const auto& h : {hCountersTrg, hCountersTrgBcSel, hLumi, hLumiBcSel}) { - h->GetXaxis()->SetBinLabel(1, "TVX"); - h->GetXaxis()->SetBinLabel(2, "TCE"); - h->GetXaxis()->SetBinLabel(3, "ZEM"); - h->GetXaxis()->SetBinLabel(4, "ZNC"); - } - } - if (context.mOptions.get("processMcData")) { - histPointers.insert({"MCreco/Stat", registry.add("MCreco/Stat", "Cut statistics; Selection criterion; Collisions", {HistType::kTH1F, {{14, -0.5, 13.5}}})}); - } - - if (isGoodRCTZdc) { - myRCTChecker.init("CBT", true); - } - } - - // process function for reconstructed data - void processData(CC const& collision, BCs const& bcs, TCs const& tracks, FWs const& fwdtracks, - aod::Zdcs const& /*zdcs*/, aod::FV0As const& fv0as, aod::FT0s const& ft0s, aod::FDDs const& fdds) - { - processReco(std::string("reco"), collision, bcs, tracks, fwdtracks, fv0as, ft0s, fdds); - } - PROCESS_SWITCH(SGCandProducer, processData, "Produce UD table with data", true); - - // process function for reconstructed MC data - void processMcData(MCCC const& collision, aod::McCollisions const& /*mccollisions*/, BCs const& bcs, - TCs const& tracks, FWs const& fwdtracks, aod::Zdcs const& /*zdcs*/, aod::FV0As const& fv0as, - aod::FT0s const& ft0s, aod::FDDs const& fdds) - { - // select specific processes with the GeneratorID - if (!collision.has_mcCollision()) - return; - auto mccol = collision.mcCollision(); - if (verboseInfo) - LOGF(info, "GeneratorId %d (%d)", mccol.getGeneratorId(), generatorIds->size()); - - if (std::find(generatorIds->begin(), generatorIds->end(), mccol.getGeneratorId()) != generatorIds->end()) { - if (verboseInfo) - LOGF(info, "Event with good generatorId"); - processReco(std::string("MCreco"), collision, bcs, tracks, fwdtracks, fv0as, ft0s, fdds); - } - } - PROCESS_SWITCH(SGCandProducer, processMcData, "Produce UD tables with MC data", false); + // add histograms for the different process functions + histPointers.clear(); + if (context.mOptions.get("processData")) { + histPointers.insert({"reco/Stat", registry.add("reco/Stat", "Cut statistics; Selection criterion; Collisions", {HistType::kTH1F, {{14, -0.5, 13.5}}})}); + + const AxisSpec axisCountersTrg{10, 0.5, 10.5, ""}; + histPointers.insert({"reco/hCountersTrg", registry.add("reco/hCountersTrg", "Trigger counts before selections; Trigger; Counts", {HistType::kTH1F, {axisCountersTrg}})}); + histPointers.insert({"reco/hCountersTrgBcSel", registry.add("reco/hCountersTrgSel", "Trigger counts after BC selections; Trigger; Counts", {HistType::kTH1F, {axisCountersTrg}})}); + histPointers.insert({"reco/hLumi", registry.add("reco/hLumi", "Integrated luminosity before selections; Trigger; Luminosity, 1/#mub", {HistType::kTH1F, {axisCountersTrg}})}); + histPointers.insert({"reco/hLumiBcSel", registry.add("reco/hLumiBcSel", "Integrated luminosity before selections; Trigger; Luminosity, 1/#mub", {HistType::kTH1F, {axisCountersTrg}})}); + auto hCountersTrg = getHist(TH1, "reco/hCountersTrg"); + auto hCountersTrgBcSel = getHist(TH1, "reco/hCountersTrgBcSel"); + auto hLumi = getHist(TH1, "reco/hLumi"); + auto hLumiBcSel = getHist(TH1, "reco/hLumiBcSel"); + for (const auto& h : {hCountersTrg, hCountersTrgBcSel, hLumi, hLumiBcSel}) { + h->GetXaxis()->SetBinLabel(1, "TVX"); + h->GetXaxis()->SetBinLabel(2, "TCE"); + h->GetXaxis()->SetBinLabel(3, "ZEM"); + h->GetXaxis()->SetBinLabel(4, "ZNC"); + } + } + if (context.mOptions.get("processMcData")) { + histPointers.insert({"MCreco/Stat", registry.add("MCreco/Stat", "Cut statistics; Selection criterion; Collisions", {HistType::kTH1F, {{14, -0.5, 13.5}}})}); + } + + if (isGoodRCTZdc) { + myRCTChecker.init("CBT", true); + } + } + + // process function for reconstructed data + void processData(CC const& collision, BCs const& bcs, TCs const& tracks, FWs const& fwdtracks, + aod::Zdcs const& /*zdcs*/, aod::FV0As const& fv0as, aod::FT0s const& ft0s, aod::FDDs const& fdds) + { + processReco(std::string("reco"), collision, bcs, tracks, fwdtracks, fv0as, ft0s, fdds); + } + PROCESS_SWITCH(SGCandProducer, processData, "Produce UD table with data", true); + + // process function for reconstructed MC data + void processMcData(MCCC const& collision, aod::McCollisions const& /*mccollisions*/, BCs const& bcs, + TCs const& tracks, FWs const& fwdtracks, aod::Zdcs const& /*zdcs*/, aod::FV0As const& fv0as, + aod::FT0s const& ft0s, aod::FDDs const& fdds) + { + // select specific processes with the GeneratorID + if (!collision.has_mcCollision()) + return; + auto mccol = collision.mcCollision(); + if (verboseInfo) + LOGF(info, "GeneratorId %d (%d)", mccol.getGeneratorId(), generatorIds->size()); + + if (std::find(generatorIds->begin(), generatorIds->end(), mccol.getGeneratorId()) != generatorIds->end()) { + if (verboseInfo) + LOGF(info, "Event with good generatorId"); + processReco(std::string("MCreco"), collision, bcs, tracks, fwdtracks, fv0as, ft0s, fdds); + } + } + PROCESS_SWITCH(SGCandProducer, processMcData, "Produce UD tables with MC data", false); }; struct McSGCandProducer { - // MC tables - Configurable verboseInfoMC{"verboseInfoMC", false, "Print general info to terminal; default it false."}; - Produces outputMcCollisions; - Produces outputMcParticles; - Produces outputMcCollsLabels; - Produces outputMcTrackLabels; - - // save all McTruth, even if the collisions is not reconstructed - Configurable> generatorIds{"generatorIds", std::vector{-1}, "MC generatorIds to process"}; - Configurable saveAllMcCollisions{"saveAllMcCollisions", true, "save all McCollisions"}; - - using CCs = soa::Join; - using BCs = soa::Join; - using TCs = soa::Join; - using UDCCs = soa::Join; - using UDTCs = soa::Join; - - // prepare slices - SliceCache cache; - PresliceUnsorted mcPartsPerMcCollision = aod::mcparticle::mcCollisionId; - Preslice udtracksPerUDCollision = aod::udtrack::udCollisionId; - - // initialize histogram registry - HistogramRegistry registry{ - "registry", - {}}; - - template - void updateUDMcCollisions(TMcCollision const& mccol, uint64_t globBC) - { - // save mccol - outputMcCollisions(globBC, - mccol.generatorsID(), - mccol.posX(), - mccol.posY(), - mccol.posZ(), - mccol.t(), - mccol.weight(), - mccol.impactParameter()); - } - - template - void updateUDMcParticle(TMcParticle const& McPart, int64_t McCollisionId, std::map& mcPartIsSaved) - { - // save McPart - // mother and daughter indices are set to -1 - // ATTENTION: this can be improved to also include mother and daughter indices - std::vector newmids; - int32_t newdids[2] = {-1, -1}; - - // update UDMcParticles - if (mcPartIsSaved.find(McPart.globalIndex()) == mcPartIsSaved.end()) { - outputMcParticles(McCollisionId, - McPart.pdgCode(), - McPart.statusCode(), - McPart.flags(), - newmids, - newdids, - McPart.weight(), - McPart.px(), - McPart.py(), - McPart.pz(), - McPart.e()); - mcPartIsSaved[McPart.globalIndex()] = outputMcParticles.lastIndex(); - } - } - - template - void updateUDMcParticles(TMcParticles const& McParts, int64_t McCollisionId, std::map& mcPartIsSaved) - { - // save McParts - // new mother and daughter ids - std::vector newmids; - int32_t newdids[2] = {-1, -1}; - int64_t newval = -1; - - // Determine the particle indices within the UDMcParticles table - // before filling the table - // This is needed to be able to assign the new daughter indices - std::map oldnew; - auto lastId = outputMcParticles.lastIndex(); - for (const auto& mcpart : McParts) { - auto oldId = mcpart.globalIndex(); - if (mcPartIsSaved.find(oldId) != mcPartIsSaved.end()) { - oldnew[oldId] = mcPartIsSaved[oldId]; - } else { - lastId++; - oldnew[oldId] = lastId; - } - } - - // all particles of the McCollision are saved - for (const auto& mcpart : McParts) { - if (mcPartIsSaved.find(mcpart.globalIndex()) == mcPartIsSaved.end()) { - // mothers - newmids.clear(); - auto oldmids = mcpart.mothersIds(); - for (const auto& oldmid : oldmids) { - auto m = McParts.rawIteratorAt(oldmid); - if (verboseInfoMC) - LOGF(debug, " m %d", m.globalIndex()); - if (mcPartIsSaved.find(oldmid) != mcPartIsSaved.end()) { - newval = mcPartIsSaved[oldmid]; - } else { - newval = -1; - } - newmids.push_back(newval); - } - // daughters - auto olddids = mcpart.daughtersIds(); - for (uint ii = 0; ii < olddids.size(); ii++) { - if (oldnew.find(olddids[ii]) != oldnew.end()) { - newval = oldnew[olddids[ii]]; - } else { - newval = -1; - } - newdids[ii] = newval; - } - if (verboseInfoMC) - LOGF(debug, " ms %i ds %i", oldmids.size(), olddids.size()); - - // update UDMcParticles - outputMcParticles(McCollisionId, - mcpart.pdgCode(), - mcpart.statusCode(), - mcpart.flags(), - newmids, - newdids, - mcpart.weight(), - mcpart.px(), - mcpart.py(), - mcpart.pz(), - mcpart.e()); - mcPartIsSaved[mcpart.globalIndex()] = outputMcParticles.lastIndex(); - } - } - } - - template - void updateUDMcTrackLabel(TTrack const& udtrack, std::map& mcPartIsSaved) - { - // udtrack (UDTCs) -> track (TCs) -> mcTrack (McParticles) -> udMcTrack (UDMcParticles) - auto trackId = udtrack.trackId(); - if (trackId >= 0) { - auto track = udtrack.template track_as(); - auto mcTrackId = track.mcParticleId(); - if (mcTrackId >= 0) { - if (mcPartIsSaved.find(mcTrackId) != mcPartIsSaved.end()) { - outputMcTrackLabels(mcPartIsSaved[mcTrackId], track.mcMask()); - } else { - outputMcTrackLabels(-1, track.mcMask()); - } - } else { - outputMcTrackLabels(-1, track.mcMask()); - } - } else { - outputMcTrackLabels(-1, -1); - } - } - - template - void updateUDMcTrackLabels(TTrack const& udtracks, std::map& mcPartIsSaved) - { - // loop over all tracks - for (const auto& udtrack : udtracks) { - // udtrack (UDTCs) -> track (TCs) -> mcTrack (McParticles) -> udMcTrack (UDMcParticles) - auto trackId = udtrack.trackId(); - if (trackId >= 0) { - auto track = udtrack.template track_as(); - auto mcTrackId = track.mcParticleId(); - if (mcTrackId >= 0) { - if (mcPartIsSaved.find(mcTrackId) != mcPartIsSaved.end()) { - outputMcTrackLabels(mcPartIsSaved[mcTrackId], track.mcMask()); - } else { - outputMcTrackLabels(-1, track.mcMask()); - } - } else { - outputMcTrackLabels(-1, track.mcMask()); - } - } else { - outputMcTrackLabels(-1, -1); - } - } - } - - // updating McTruth data and links to reconstructed data - void procWithSgCand(aod::McCollisions const& mccols, aod::McParticles const& mcparts, - UDCCs const& sgcands, UDTCs const& udtracks) - { - // use a hash table to keep track of the McCollisions which have been added to the UDMcCollision table - // {McCollisionId : udMcCollisionId} - // similar for the McParticles which have been added to the UDMcParticle table - // {McParticleId : udMcParticleId} - std::map mcColIsSaved; - std::map mcPartIsSaved; - - // loop over McCollisions and UDCCs simultaneously - auto mccol = mccols.iteratorAt(0); - auto mcOfInterest = std::find(generatorIds->begin(), generatorIds->end(), mccol.getGeneratorId()) != generatorIds->end(); - auto lastmccol = mccols.iteratorAt(mccols.size() - 1); - auto mccolAtEnd = false; - - auto sgcand = sgcands.iteratorAt(0); - auto lastsgcand = sgcands.iteratorAt(sgcands.size() - 1); - auto sgcandAtEnd = false; - - // advance dgcand and mccol until both are AtEnd - int64_t mccolId = mccol.globalIndex(); - int64_t mcsgId = -1; - bool goon = true; - while (goon) { - auto globBC = mccol.bc_as().globalBC(); - // check if dgcand has an associated McCollision - if (sgcand.has_collision()) { - auto sgcandCol = sgcand.collision_as(); - // colId = sgcandCol.globalIndex(); - if (sgcandCol.has_mcCollision()) { - mcsgId = sgcandCol.mcCollision().globalIndex(); - } else { - mcsgId = -1; - } - } else { - mcsgId = -1; - } - if (verboseInfoMC) - LOGF(info, "\nStart of loop mcsgId %d mccolId %d", mcsgId, mccolId); - - // two cases to consider - // 1. mcdgId <= mccolId: the event to process is a dgcand. In this case the Mc tables as well as the McLabel tables are updated - // 2. mccolId < mcdgId: the event to process is an MC event of interest without reconstructed dgcand. In this case only the Mc tables are updated - if ((!sgcandAtEnd && !mccolAtEnd && (mcsgId <= mccolId)) || mccolAtEnd) { - // this is case 1. - if (verboseInfoMC) - LOGF(info, "Doing case 1 with mcsgId %d", mcsgId); - - // update UDMcCollisions and UDMcColsLabels (for each UDCollision -> UDMcCollisions) - // update UDMcParticles and UDMcTrackLabels (for each UDTrack -> UDMcParticles) - // get dgcand tracks - auto sgTracks = udtracks.sliceByCached(aod::udtrack::udCollisionId, sgcand.globalIndex(), cache); - - // If the sgcand has an associated McCollision then the McCollision and all associated - // McParticles are saved - // but only consider generated events of interest - if (mcsgId >= 0 && mcOfInterest) { - if (mcColIsSaved.find(mcsgId) == mcColIsSaved.end()) { - if (verboseInfoMC) - LOGF(info, " Saving McCollision %d", mcsgId); - // update UDMcCollisions - auto sgcandMcCol = sgcand.collision_as().mcCollision(); - updateUDMcCollisions(sgcandMcCol, globBC); - mcColIsSaved[mcsgId] = outputMcCollisions.lastIndex(); - } - - // update UDMcColsLabels (for each UDCollision -> UDMcCollisions) - outputMcCollsLabels(mcColIsSaved[mcsgId]); - - // update UDMcParticles - auto mcPartsSlice = mcparts.sliceBy(mcPartsPerMcCollision, mcsgId); - updateUDMcParticles(mcPartsSlice, mcColIsSaved[mcsgId], mcPartIsSaved); - - // update UDMcTrackLabels (for each UDTrack -> UDMcParticles) - updateUDMcTrackLabels(sgTracks, mcPartIsSaved); - - } else { - // If the sgcand has no associated McCollision then only the McParticles which are associated - // with the tracks of the sgcand are saved - if (verboseInfoMC) - LOGF(info, " Saving McCollision %d", -1); - - // update UDMcColsLabels (for each UDCollision -> UDMcCollisions) - outputMcCollsLabels(-1); - - // update UDMcParticles and UDMcTrackLabels (for each UDTrack -> UDMcParticles) - // loop over tracks of dgcand - for (const auto& sgtrack : sgTracks) { - if (sgtrack.has_track()) { - auto track = sgtrack.track_as(); - if (track.has_mcParticle()) { - auto mcPart = track.mcParticle(); - auto mcCol = mcPart.mcCollision(); - if (mcColIsSaved.find(mcCol.globalIndex()) == mcColIsSaved.end()) { - updateUDMcCollisions(mcCol, globBC); - mcColIsSaved[mcCol.globalIndex()] = outputMcCollisions.lastIndex(); - } - updateUDMcParticle(mcPart, mcColIsSaved[mcCol.globalIndex()], mcPartIsSaved); - updateUDMcTrackLabel(sgtrack, mcPartIsSaved); - } else { - outputMcTrackLabels(-1, track.mcMask()); - } - } else { - outputMcTrackLabels(-1, -1); - } - } - } - // advance sgcand - if (sgcand != lastsgcand) { - sgcand++; - } else { - sgcandAtEnd = true; - } - } else { - // this is case 2. - if (verboseInfoMC) - LOGF(info, "Doing case 2"); - - // update UDMcCollisions and UDMcParticles - // but only consider generated events of interest - if (mcOfInterest && mcColIsSaved.find(mccolId) == mcColIsSaved.end()) { - if (verboseInfoMC) - LOGF(info, " Saving McCollision %d", mccolId); - // update UDMcCollisions - updateUDMcCollisions(mccol, globBC); - mcColIsSaved[mccolId] = outputMcCollisions.lastIndex(); - - // update UDMcParticles - auto mcPartsSlice = mcparts.sliceBy(mcPartsPerMcCollision, mccolId); - updateUDMcParticles(mcPartsSlice, mcColIsSaved[mccolId], mcPartIsSaved); - } - - // advance mccol - if (mccol != lastmccol) { - mccol++; - mcOfInterest = std::find(generatorIds->begin(), generatorIds->end(), mccol.getGeneratorId()) != generatorIds->end(); - mccolId = mccol.globalIndex(); - } else { - mccolAtEnd = true; - } - } - - goon = !sgcandAtEnd || !mccolAtEnd; - if (verboseInfoMC) - LOGF(info, "End of loop mcsgId %d mccolId %d", mcsgId, mccolId); - } - } - - // updating McTruth data only - void procWithoutSgCand(aod::McCollisions const& mccols, aod::McParticles const& mcparts) - { - // use a hash table to keep track of the McCollisions which have been added to the UDMcCollision table - // {McCollisionId : udMcCollisionId} - // similar for the McParticles which have been added to the UDMcParticle table - // {McParticleId : udMcParticleId} - std::map mcColIsSaved; - std::map mcPartIsSaved; - - // loop over McCollisions - for (auto const& mccol : mccols) { - int64_t mccolId = mccol.globalIndex(); - uint64_t globBC = mccol.bc_as().globalBC(); - - // update UDMcCollisions and UDMcParticles - if (mcColIsSaved.find(mccolId) == mcColIsSaved.end()) { - if (verboseInfoMC) - LOGF(info, " Saving McCollision %d", mccolId); - - // update UDMcCollisions - updateUDMcCollisions(mccol, globBC); - mcColIsSaved[mccolId] = outputMcCollisions.lastIndex(); - - // update UDMcParticles - auto mcPartsSlice = mcparts.sliceBy(mcPartsPerMcCollision, mccolId); - updateUDMcParticles(mcPartsSlice, mcColIsSaved[mccolId], mcPartIsSaved); - } - } - } - - void init(InitContext& context) - { - // add histograms for the different process functions - if (context.mOptions.get("processMC")) { - registry.add("mcTruth/collisions", "Number of associated collisions", {HistType::kTH1F, {{11, -0.5, 10.5}}}); - registry.add("mcTruth/collType", "Collision type", {HistType::kTH1F, {{5, -0.5, 4.5}}}); - registry.add("mcTruth/IVMpt", "Invariant mass versus p_{T}", {HistType::kTH2F, {{150, 0.0, 3.0}, {150, 0.0, 3.0}}}); - } - } - - // process function for MC data - // save the MC truth of all events of interest and of the DG events - void processMC(aod::McCollisions const& mccols, aod::McParticles const& mcparts, - UDCCs const& sgcands, UDTCs const& udtracks, - CCs const& /*collisions*/, BCs const& /*bcs*/, TCs const& /*tracks*/) - { - if (verboseInfoMC) { - LOGF(info, "Number of McCollisions %d", mccols.size()); - LOGF(info, "Number of SG candidates %d", sgcands.size()); - LOGF(info, "Number of UD tracks %d", udtracks.size()); - } - if (mccols.size() > 0) { - if (sgcands.size() > 0) { - procWithSgCand(mccols, mcparts, sgcands, udtracks); - } else { - if (saveAllMcCollisions) { - procWithoutSgCand(mccols, mcparts); - } - } - } - } - PROCESS_SWITCH(McSGCandProducer, processMC, "Produce MC tables", false); - - void processDummy(aod::Collisions const& /*collisions*/) - { - // do nothing - if (verboseInfoMC) - LOGF(info, "Running dummy process function!"); - } - PROCESS_SWITCH(McSGCandProducer, processDummy, "Dummy function", true); + // MC tables + Configurable verboseInfoMC{"verboseInfoMC", false, "Print general info to terminal; default it false."}; + Produces outputMcCollisions; + Produces outputMcParticles; + Produces outputMcCollsLabels; + Produces outputMcTrackLabels; + + // save all McTruth, even if the collisions is not reconstructed + Configurable> generatorIds{"generatorIds", std::vector{-1}, "MC generatorIds to process"}; + Configurable saveAllMcCollisions{"saveAllMcCollisions", true, "save all McCollisions"}; + + using CCs = soa::Join; + using BCs = soa::Join; + using TCs = soa::Join; + using UDCCs = soa::Join; + using UDTCs = soa::Join; + + // prepare slices + SliceCache cache; + PresliceUnsorted mcPartsPerMcCollision = aod::mcparticle::mcCollisionId; + Preslice udtracksPerUDCollision = aod::udtrack::udCollisionId; + + // initialize histogram registry + HistogramRegistry registry{ + "registry", + {}}; + + template + void updateUDMcCollisions(TMcCollision const& mccol, uint64_t globBC) + { + // save mccol + outputMcCollisions(globBC, + mccol.generatorsID(), + mccol.posX(), + mccol.posY(), + mccol.posZ(), + mccol.t(), + mccol.weight(), + mccol.impactParameter()); + } + + template + void updateUDMcParticle(TMcParticle const& McPart, int64_t McCollisionId, std::map& mcPartIsSaved) + { + // save McPart + // mother and daughter indices are set to -1 + // ATTENTION: this can be improved to also include mother and daughter indices + std::vector newmids; + int32_t newdids[2] = {-1, -1}; + + // update UDMcParticles + if (mcPartIsSaved.find(McPart.globalIndex()) == mcPartIsSaved.end()) { + outputMcParticles(McCollisionId, + McPart.pdgCode(), + McPart.statusCode(), + McPart.flags(), + newmids, + newdids, + McPart.weight(), + McPart.px(), + McPart.py(), + McPart.pz(), + McPart.e()); + mcPartIsSaved[McPart.globalIndex()] = outputMcParticles.lastIndex(); + } + } + + template + void updateUDMcParticles(TMcParticles const& McParts, int64_t McCollisionId, std::map& mcPartIsSaved) + { + // save McParts + // new mother and daughter ids + std::vector newmids; + int32_t newdids[2] = {-1, -1}; + int64_t newval = -1; + + // Determine the particle indices within the UDMcParticles table + // before filling the table + // This is needed to be able to assign the new daughter indices + std::map oldnew; + auto lastId = outputMcParticles.lastIndex(); + for (const auto& mcpart : McParts) { + auto oldId = mcpart.globalIndex(); + if (mcPartIsSaved.find(oldId) != mcPartIsSaved.end()) { + oldnew[oldId] = mcPartIsSaved[oldId]; + } else { + lastId++; + oldnew[oldId] = lastId; + } + } + + // all particles of the McCollision are saved + for (const auto& mcpart : McParts) { + if (mcPartIsSaved.find(mcpart.globalIndex()) == mcPartIsSaved.end()) { + // mothers + newmids.clear(); + auto oldmids = mcpart.mothersIds(); + for (const auto& oldmid : oldmids) { + auto m = McParts.rawIteratorAt(oldmid); + if (verboseInfoMC) + LOGF(debug, " m %d", m.globalIndex()); + if (mcPartIsSaved.find(oldmid) != mcPartIsSaved.end()) { + newval = mcPartIsSaved[oldmid]; + } else { + newval = -1; + } + newmids.push_back(newval); + } + // daughters + auto olddids = mcpart.daughtersIds(); + for (uint ii = 0; ii < olddids.size(); ii++) { + if (oldnew.find(olddids[ii]) != oldnew.end()) { + newval = oldnew[olddids[ii]]; + } else { + newval = -1; + } + newdids[ii] = newval; + } + if (verboseInfoMC) + LOGF(debug, " ms %i ds %i", oldmids.size(), olddids.size()); + + // update UDMcParticles + outputMcParticles(McCollisionId, + mcpart.pdgCode(), + mcpart.statusCode(), + mcpart.flags(), + newmids, + newdids, + mcpart.weight(), + mcpart.px(), + mcpart.py(), + mcpart.pz(), + mcpart.e()); + mcPartIsSaved[mcpart.globalIndex()] = outputMcParticles.lastIndex(); + } + } + } + + template + void updateUDMcTrackLabel(TTrack const& udtrack, std::map& mcPartIsSaved) + { + // udtrack (UDTCs) -> track (TCs) -> mcTrack (McParticles) -> udMcTrack (UDMcParticles) + auto trackId = udtrack.trackId(); + if (trackId >= 0) { + auto track = udtrack.template track_as(); + auto mcTrackId = track.mcParticleId(); + if (mcTrackId >= 0) { + if (mcPartIsSaved.find(mcTrackId) != mcPartIsSaved.end()) { + outputMcTrackLabels(mcPartIsSaved[mcTrackId], track.mcMask()); + } else { + outputMcTrackLabels(-1, track.mcMask()); + } + } else { + outputMcTrackLabels(-1, track.mcMask()); + } + } else { + outputMcTrackLabels(-1, -1); + } + } + + template + void updateUDMcTrackLabels(TTrack const& udtracks, std::map& mcPartIsSaved) + { + // loop over all tracks + for (const auto& udtrack : udtracks) { + // udtrack (UDTCs) -> track (TCs) -> mcTrack (McParticles) -> udMcTrack (UDMcParticles) + auto trackId = udtrack.trackId(); + if (trackId >= 0) { + auto track = udtrack.template track_as(); + auto mcTrackId = track.mcParticleId(); + if (mcTrackId >= 0) { + if (mcPartIsSaved.find(mcTrackId) != mcPartIsSaved.end()) { + outputMcTrackLabels(mcPartIsSaved[mcTrackId], track.mcMask()); + } else { + outputMcTrackLabels(-1, track.mcMask()); + } + } else { + outputMcTrackLabels(-1, track.mcMask()); + } + } else { + outputMcTrackLabels(-1, -1); + } + } + } + + // updating McTruth data and links to reconstructed data + void procWithSgCand(aod::McCollisions const& mccols, aod::McParticles const& mcparts, + UDCCs const& sgcands, UDTCs const& udtracks) + { + // use a hash table to keep track of the McCollisions which have been added to the UDMcCollision table + // {McCollisionId : udMcCollisionId} + // similar for the McParticles which have been added to the UDMcParticle table + // {McParticleId : udMcParticleId} + std::map mcColIsSaved; + std::map mcPartIsSaved; + + // loop over McCollisions and UDCCs simultaneously + auto mccol = mccols.iteratorAt(0); + auto mcOfInterest = std::find(generatorIds->begin(), generatorIds->end(), mccol.getGeneratorId()) != generatorIds->end(); + auto lastmccol = mccols.iteratorAt(mccols.size() - 1); + auto mccolAtEnd = false; + + auto sgcand = sgcands.iteratorAt(0); + auto lastsgcand = sgcands.iteratorAt(sgcands.size() - 1); + auto sgcandAtEnd = false; + + // advance dgcand and mccol until both are AtEnd + int64_t mccolId = mccol.globalIndex(); + int64_t mcsgId = -1; + bool goon = true; + while (goon) { + auto globBC = mccol.bc_as().globalBC(); + // check if dgcand has an associated McCollision + if (sgcand.has_collision()) { + auto sgcandCol = sgcand.collision_as(); + // colId = sgcandCol.globalIndex(); + if (sgcandCol.has_mcCollision()) { + mcsgId = sgcandCol.mcCollision().globalIndex(); + } else { + mcsgId = -1; + } + } else { + mcsgId = -1; + } + if (verboseInfoMC) + LOGF(info, "\nStart of loop mcsgId %d mccolId %d", mcsgId, mccolId); + + // two cases to consider + // 1. mcdgId <= mccolId: the event to process is a dgcand. In this case the Mc tables as well as the McLabel tables are updated + // 2. mccolId < mcdgId: the event to process is an MC event of interest without reconstructed dgcand. In this case only the Mc tables are updated + if ((!sgcandAtEnd && !mccolAtEnd && (mcsgId <= mccolId)) || mccolAtEnd) { + // this is case 1. + if (verboseInfoMC) + LOGF(info, "Doing case 1 with mcsgId %d", mcsgId); + + // update UDMcCollisions and UDMcColsLabels (for each UDCollision -> UDMcCollisions) + // update UDMcParticles and UDMcTrackLabels (for each UDTrack -> UDMcParticles) + // get dgcand tracks + auto sgTracks = udtracks.sliceByCached(aod::udtrack::udCollisionId, sgcand.globalIndex(), cache); + + // If the sgcand has an associated McCollision then the McCollision and all associated + // McParticles are saved + // but only consider generated events of interest + if (mcsgId >= 0 && mcOfInterest) { + if (mcColIsSaved.find(mcsgId) == mcColIsSaved.end()) { + if (verboseInfoMC) + LOGF(info, " Saving McCollision %d", mcsgId); + // update UDMcCollisions + auto sgcandMcCol = sgcand.collision_as().mcCollision(); + updateUDMcCollisions(sgcandMcCol, globBC); + mcColIsSaved[mcsgId] = outputMcCollisions.lastIndex(); + } + + // update UDMcColsLabels (for each UDCollision -> UDMcCollisions) + outputMcCollsLabels(mcColIsSaved[mcsgId]); + + // update UDMcParticles + auto mcPartsSlice = mcparts.sliceBy(mcPartsPerMcCollision, mcsgId); + updateUDMcParticles(mcPartsSlice, mcColIsSaved[mcsgId], mcPartIsSaved); + + // update UDMcTrackLabels (for each UDTrack -> UDMcParticles) + updateUDMcTrackLabels(sgTracks, mcPartIsSaved); + + } else { + // If the sgcand has no associated McCollision then only the McParticles which are associated + // with the tracks of the sgcand are saved + if (verboseInfoMC) + LOGF(info, " Saving McCollision %d", -1); + + // update UDMcColsLabels (for each UDCollision -> UDMcCollisions) + outputMcCollsLabels(-1); + + // update UDMcParticles and UDMcTrackLabels (for each UDTrack -> UDMcParticles) + // loop over tracks of dgcand + for (const auto& sgtrack : sgTracks) { + if (sgtrack.has_track()) { + auto track = sgtrack.track_as(); + if (track.has_mcParticle()) { + auto mcPart = track.mcParticle(); + auto mcCol = mcPart.mcCollision(); + if (mcColIsSaved.find(mcCol.globalIndex()) == mcColIsSaved.end()) { + updateUDMcCollisions(mcCol, globBC); + mcColIsSaved[mcCol.globalIndex()] = outputMcCollisions.lastIndex(); + } + updateUDMcParticle(mcPart, mcColIsSaved[mcCol.globalIndex()], mcPartIsSaved); + updateUDMcTrackLabel(sgtrack, mcPartIsSaved); + } else { + outputMcTrackLabels(-1, track.mcMask()); + } + } else { + outputMcTrackLabels(-1, -1); + } + } + } + // advance sgcand + if (sgcand != lastsgcand) { + sgcand++; + } else { + sgcandAtEnd = true; + } + } else { + // this is case 2. + if (verboseInfoMC) + LOGF(info, "Doing case 2"); + + // update UDMcCollisions and UDMcParticles + // but only consider generated events of interest + if (mcOfInterest && mcColIsSaved.find(mccolId) == mcColIsSaved.end()) { + if (verboseInfoMC) + LOGF(info, " Saving McCollision %d", mccolId); + // update UDMcCollisions + updateUDMcCollisions(mccol, globBC); + mcColIsSaved[mccolId] = outputMcCollisions.lastIndex(); + + // update UDMcParticles + auto mcPartsSlice = mcparts.sliceBy(mcPartsPerMcCollision, mccolId); + updateUDMcParticles(mcPartsSlice, mcColIsSaved[mccolId], mcPartIsSaved); + } + + // advance mccol + if (mccol != lastmccol) { + mccol++; + mcOfInterest = std::find(generatorIds->begin(), generatorIds->end(), mccol.getGeneratorId()) != generatorIds->end(); + mccolId = mccol.globalIndex(); + } else { + mccolAtEnd = true; + } + } + + goon = !sgcandAtEnd || !mccolAtEnd; + if (verboseInfoMC) + LOGF(info, "End of loop mcsgId %d mccolId %d", mcsgId, mccolId); + } + } + + // updating McTruth data only + void procWithoutSgCand(aod::McCollisions const& mccols, aod::McParticles const& mcparts) + { + // use a hash table to keep track of the McCollisions which have been added to the UDMcCollision table + // {McCollisionId : udMcCollisionId} + // similar for the McParticles which have been added to the UDMcParticle table + // {McParticleId : udMcParticleId} + std::map mcColIsSaved; + std::map mcPartIsSaved; + + // loop over McCollisions + for (auto const& mccol : mccols) { + int64_t mccolId = mccol.globalIndex(); + uint64_t globBC = mccol.bc_as().globalBC(); + + // update UDMcCollisions and UDMcParticles + if (mcColIsSaved.find(mccolId) == mcColIsSaved.end()) { + if (verboseInfoMC) + LOGF(info, " Saving McCollision %d", mccolId); + + // update UDMcCollisions + updateUDMcCollisions(mccol, globBC); + mcColIsSaved[mccolId] = outputMcCollisions.lastIndex(); + + // update UDMcParticles + auto mcPartsSlice = mcparts.sliceBy(mcPartsPerMcCollision, mccolId); + updateUDMcParticles(mcPartsSlice, mcColIsSaved[mccolId], mcPartIsSaved); + } + } + } + + void init(InitContext& context) + { + // add histograms for the different process functions + if (context.mOptions.get("processMC")) { + registry.add("mcTruth/collisions", "Number of associated collisions", {HistType::kTH1F, {{11, -0.5, 10.5}}}); + registry.add("mcTruth/collType", "Collision type", {HistType::kTH1F, {{5, -0.5, 4.5}}}); + registry.add("mcTruth/IVMpt", "Invariant mass versus p_{T}", {HistType::kTH2F, {{150, 0.0, 3.0}, {150, 0.0, 3.0}}}); + } + } + + // process function for MC data + // save the MC truth of all events of interest and of the DG events + void processMC(aod::McCollisions const& mccols, aod::McParticles const& mcparts, + UDCCs const& sgcands, UDTCs const& udtracks, + CCs const& /*collisions*/, BCs const& /*bcs*/, TCs const& /*tracks*/) + { + if (verboseInfoMC) { + LOGF(info, "Number of McCollisions %d", mccols.size()); + LOGF(info, "Number of SG candidates %d", sgcands.size()); + LOGF(info, "Number of UD tracks %d", udtracks.size()); + } + if (mccols.size() > 0) { + if (sgcands.size() > 0) { + procWithSgCand(mccols, mcparts, sgcands, udtracks); + } else { + if (saveAllMcCollisions) { + procWithoutSgCand(mccols, mcparts); + } + } + } + } + PROCESS_SWITCH(McSGCandProducer, processMC, "Produce MC tables", false); + + void processDummy(aod::Collisions const& /*collisions*/) + { + // do nothing + if (verboseInfoMC) + LOGF(info, "Running dummy process function!"); + } + PROCESS_SWITCH(McSGCandProducer, processDummy, "Dummy function", true); }; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { - WorkflowSpec workflow{ - adaptAnalysisTask(cfgc, TaskName{"sgcandproducer"}), - adaptAnalysisTask(cfgc, TaskName{"mcsgcandproducer"})}; - return workflow; + WorkflowSpec workflow{ + adaptAnalysisTask(cfgc, TaskName{"sgcandproducer"}), + adaptAnalysisTask(cfgc, TaskName{"mcsgcandproducer"})}; + return workflow; } From 79d59ac43088e91b074439a40a285135a46731f2 Mon Sep 17 00:00:00 2001 From: sandor-lokos Date: Tue, 24 Mar 2026 20:31:33 +0100 Subject: [PATCH 07/11] Formatting fixes --- PWGUD/Core/UDHelpers.h | 2 +- PWGUD/TableProducer/SGCandProducer.cxx | 1584 ++++++++++++------------ 2 files changed, 793 insertions(+), 793 deletions(-) diff --git a/PWGUD/Core/UDHelpers.h b/PWGUD/Core/UDHelpers.h index ea137fbc17a..9dbd70c6897 100644 --- a/PWGUD/Core/UDHelpers.h +++ b/PWGUD/Core/UDHelpers.h @@ -677,7 +677,7 @@ inline bool getPhiEtaFromFitBit(FT0DetT& ft0Det, double& phi, double& eta) { - auto ref = decodeFitBit(bit); + auto ref = decodeFitBit(bit); if (ref.det != FitBitRef::Det::FT0) { return false; } diff --git a/PWGUD/TableProducer/SGCandProducer.cxx b/PWGUD/TableProducer/SGCandProducer.cxx index 80ec14a0b1f..da3c052c036 100644 --- a/PWGUD/TableProducer/SGCandProducer.cxx +++ b/PWGUD/TableProducer/SGCandProducer.cxx @@ -57,803 +57,803 @@ using namespace o2::aod::rctsel; #define getHist(type, name) std::get>(histPointers[name]) struct SGCandProducer { - Service ccdb; - // data inputs - using CCs = soa::Join; - using CC = CCs::iterator; - using BCs = soa::Join; - using BC = BCs::iterator; - using TCs = soa::Join; - using FWs = aod::FwdTracks; - - using MCCCs = soa::Join; - using MCCC = MCCCs::iterator; - - // get an SGCutparHolder - SGCutParHolder sameCuts = SGCutParHolder(); // SGCutparHolder - Configurable SGCuts{"SGCuts", {}, "SG event cuts"}; - FITCutParHolder fitCuts = FITCutParHolder(); - Configurable FITCuts{"FITCuts", {}, "FIT bitset cuts"}; - Configurable verboseInfo{"verboseInfo", false, "Print general info to terminal; default it false."}; - Configurable saveAllTracks{"saveAllTracks", true, "save only PV contributors or all tracks associated to a collision"}; - Configurable savenonPVCITSOnlyTracks{"savenonPVCITSOnlyTracks", false, "save non PV contributors with ITS only information"}; - Configurable rejectAtTFBoundary{"rejectAtTFBoundary", true, "reject collisions at a TF boundary"}; - Configurable noITSROFrameBorder{"noITSROFrameBorder", true, "reject ITS RO Frame Border"}; - Configurable noSameBunchPileUp{"noSameBunchPileUp", true, "reject SameBunchPileUp"}; - Configurable IsGoodVertex{"IsGoodVertex", false, "Select FT0 PV vertex matching"}; - Configurable ITSTPCVertex{"ITSTPCVertex", true, "reject ITS-only vertex"}; // if one wants to look at Single Gap pp events - Configurable> generatorIds{"generatorIds", std::vector{-1}, "MC generatorIds to process"}; - Configurable storeSG{"storeSG", true, "Store SG events in the output"}; - Configurable storeDG{"storeDG", true, "Store DG events in the output"}; - - Configurable saveFITbitsets{"saveFITbitsets", true, "Write FT0 and FV0 bitset tables to output"}; - - Configurable isGoodRCTCollision{"isGoodRCTCollision", true, "Check RCT flags for FT0,ITS,TPC and tracking"}; - Configurable isGoodRCTZdc{"isGoodRCTZdc", false, "Check RCT flags for ZDC if present in run"}; - - // Configurables to decide which tables are filled - Configurable fillTrackTables{"fillTrackTables", true, "Fill track tables"}; - Configurable fillFwdTrackTables{"fillFwdTrackTables", true, "Fill forward track tables"}; - - // SG selector - SGSelector sgSelector; - ctpRateFetcher mRateFetcher; - - // initialize RCT flag checker - RCTFlagsChecker myRCTChecker{"CBT"}; - - // data tables - Produces outputSGCollisions; - Produces outputCollisions; - Produces outputCollisionsSels; - Produces outputCollisionSelExtras; - Produces outputCollsLabels; - Produces outputZdcs; - Produces udZdcsReduced; - Produces outputTracks; - Produces outputTracksCov; - Produces outputTracksDCA; - Produces outputTracksPID; - Produces outputTracksPIDExtra; - Produces outputTracksExtra; - Produces outputTracksFlag; - Produces outputFwdTracks; - Produces outputFwdTracksExtra; - Produces outputTracksLabel; - Produces outputFITBits; - - // initialize histogram registry - HistogramRegistry registry{ - "registry", - {}}; - std::map histPointers; - - int runNumber = -1; - - // function to update UDFwdTracks, UDFwdTracksExtra - template - void updateUDFwdTrackTables(TFwdTrack const& fwdtrack, uint64_t const& bcnum) - { - outputFwdTracks(outputCollisions.lastIndex(), - fwdtrack.px(), fwdtrack.py(), fwdtrack.pz(), fwdtrack.sign(), - bcnum, fwdtrack.trackTime(), fwdtrack.trackTimeRes()); - outputFwdTracksExtra(fwdtrack.trackType(), - fwdtrack.nClusters(), - fwdtrack.pDca(), - fwdtrack.rAtAbsorberEnd(), - fwdtrack.chi2(), - fwdtrack.chi2MatchMCHMID(), - fwdtrack.chi2MatchMCHMFT(), - fwdtrack.mchBitMap(), - fwdtrack.midBitMap(), - fwdtrack.midBoards()); - } - - // function to update UDTracks, UDTracksCov, UDTracksDCA, UDTracksPID, UDTracksExtra, UDTracksFlag, - // and UDTrackCollisionIDs - template - void updateUDTrackTables(int64_t lastIndex, TTrack const& track, uint64_t const& bcnum) - { - outputTracks(lastIndex, - track.px(), track.py(), track.pz(), track.sign(), - bcnum, track.trackTime(), track.trackTimeRes()); - - // float sigmaY = track.sigmaY(); - // float sigmaZ = track.sigmaZ(); - float sigmaY = -1.; - float sigmaZ = -1.; - outputTracksCov(track.x(), track.y(), track.z(), sigmaY, sigmaZ); - - outputTracksDCA(track.dcaZ(), track.dcaXY()); - outputTracksPID(track.tpcNSigmaEl(), - track.tpcNSigmaMu(), - track.tpcNSigmaPi(), - track.tpcNSigmaKa(), - track.tpcNSigmaPr(), - track.beta(), - track.betaerror(), - track.tofNSigmaEl(), - track.tofNSigmaMu(), - track.tofNSigmaPi(), - track.tofNSigmaKa(), - track.tofNSigmaPr()); - outputTracksPIDExtra(track.tpcNSigmaDe(), - track.tpcNSigmaTr(), - track.tpcNSigmaHe(), - track.tpcNSigmaAl(), - track.tofNSigmaDe(), - track.tofNSigmaTr(), - track.tofNSigmaHe(), - track.tofNSigmaAl()); - outputTracksExtra(track.tpcInnerParam(), - track.itsClusterSizes(), - track.tpcNClsFindable(), - track.tpcNClsFindableMinusFound(), - track.tpcNClsFindableMinusCrossedRows(), - track.tpcNClsShared(), - track.trdPattern(), - track.itsChi2NCl(), - track.tpcChi2NCl(), - track.trdChi2(), - track.tofChi2(), - track.tpcSignal(), - track.tofSignal(), - track.trdSignal(), - track.length(), - track.tofExpMom(), - track.detectorMap()); - outputTracksFlag(track.has_collision(), - track.isPVContributor()); - outputTracksLabel(track.globalIndex()); - } - - // function to process reconstructed data - template - void processReco(std::string histdir, TCol const& collision, BCs const& bcs, - TCs const& tracks, FWs const& fwdtracks, - aod::FV0As const& fv0as, aod::FT0s const& ft0s, aod::FDDs const& fdds) - { - if (verboseInfo) - LOGF(debug, " collision %d", collision.globalIndex()); - getHist(TH1, histdir + "/Stat")->Fill(0., 1.); - // reject collisions at TF boundaries - if (rejectAtTFBoundary && !collision.selection_bit(aod::evsel::kNoTimeFrameBorder)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(1., 1.); - // reject collisions at ITS RO TF boundaries - if (noITSROFrameBorder && !collision.selection_bit(aod::evsel::kNoITSROFrameBorder)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(2., 1.); - // reject Same Bunch PileUp - if (noSameBunchPileUp && !collision.selection_bit(aod::evsel::kNoSameBunchPileup)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(3., 1.); - // check vertex matching to FT0 - if (IsGoodVertex && !collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(4., 1.); - // reject ITS Only vertices - if (ITSTPCVertex && !collision.selection_bit(aod::evsel::kIsVertexITSTPC)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(5., 1.); - // nominal BC - if (!collision.has_foundBC()) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(6., 1.); - // RCT CBT for collision check - if (isGoodRCTCollision && !myRCTChecker(collision)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(7., 1.); - // RCT CBT+ZDC for collision check - if (isGoodRCTZdc && !myRCTChecker(collision)) { - return; - } - getHist(TH1, histdir + "/Stat")->Fill(8., 1.); - - // - const int trs = collision.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard) ? 1 : 0; - const int trofs = collision.selection_bit(o2::aod::evsel::kNoCollInRofStandard) ? 1 : 0; - const int hmpr = collision.selection_bit(o2::aod::evsel::kNoHighMultCollInPrevRof) ? 1 : 0; - const int tfb = collision.selection_bit(o2::aod::evsel::kNoTimeFrameBorder) ? 1 : 0; - const int itsROFb = collision.selection_bit(o2::aod::evsel::kNoITSROFrameBorder) ? 1 : 0; - const int sbp = collision.selection_bit(o2::aod::evsel::kNoSameBunchPileup) ? 1 : 0; - const int zVtxFT0vPv = collision.selection_bit(o2::aod::evsel::kIsGoodZvtxFT0vsPV) ? 1 : 0; - const int vtxITSTPC = collision.selection_bit(o2::aod::evsel::kIsVertexITSTPC) ? 1 : 0; - auto bc = collision.template foundBC_as(); - double ir = 0.; - const uint64_t ts = bc.timestamp(); - const int runnumber = bc.runNumber(); - if (bc.has_zdc()) { - ir = mRateFetcher.fetch(ccdb.service, ts, runnumber, "ZNC hadronic") * 1.e-3; - } - auto newbc = bc; - - // obtain slice of compatible BCs - auto bcRange = udhelpers::compatibleBCs(collision, sameCuts.NDtcoll(), bcs, sameCuts.minNBCs()); - auto isSGEvent = sgSelector.IsSelected(sameCuts, collision, bcRange, bc); - // auto isSGEvent = sgSelector.IsSelected(sameCuts, collision, bcRange, tracks); - int issgevent = isSGEvent.value; - if (isSGEvent.bc && issgevent < 2) { - newbc = *(isSGEvent.bc); - } else { - if (verboseInfo) - LOGF(info, "No Newbc %i", bc.globalBC()); - } - getHist(TH1, histdir + "/Stat")->Fill(issgevent + 10, 1.); - if ((storeDG && issgevent == o2::aod::sgselector::DoubleGap) || (storeSG && (issgevent == o2::aod::sgselector::SingleGapA || issgevent == o2::aod::sgselector::SingleGapC))) { - if (verboseInfo) - LOGF(info, "Current BC: %i, %i, %i", bc.globalBC(), newbc.globalBC(), issgevent); - if (sameCuts.minRgtrwTOF()) { - if (udhelpers::rPVtrwTOF(tracks, collision.numContrib()) < sameCuts.minRgtrwTOF()) - return; - } - upchelpers::FITInfo fitInfo{}; - const uint8_t chFT0A = 0; - const uint8_t chFT0C = 0; - const uint8_t chFDDA = 0; - const uint8_t chFDDC = 0; - const uint8_t chFV0A = 0; - const int occ = collision.trackOccupancyInTimeRange(); - udhelpers::getFITinfo(fitInfo, newbc, bcs, ft0s, fv0as, fdds); - const int upc_flag = (collision.flags() & dataformats::Vertex>::Flags::UPCMode) ? 1 : 0; - // update SG candidates tables - outputCollisions(bc.globalBC(), bc.runNumber(), - collision.posX(), collision.posY(), collision.posZ(), upc_flag, - collision.numContrib(), udhelpers::netCharge(tracks), - 1.); // rtrwTOF); //omit the calculation to speed up the things while skimming - - outputSGCollisions(issgevent); - outputCollisionsSels(fitInfo.ampFT0A, fitInfo.ampFT0C, fitInfo.timeFT0A, fitInfo.timeFT0C, - fitInfo.triggerMaskFT0, - fitInfo.ampFDDA, fitInfo.ampFDDC, fitInfo.timeFDDA, fitInfo.timeFDDC, - fitInfo.triggerMaskFDD, - fitInfo.ampFV0A, fitInfo.timeFV0A, fitInfo.triggerMaskFV0A, - fitInfo.BBFT0Apf, fitInfo.BBFT0Cpf, fitInfo.BGFT0Apf, fitInfo.BGFT0Cpf, - fitInfo.BBFV0Apf, fitInfo.BGFV0Apf, - fitInfo.BBFDDApf, fitInfo.BBFDDCpf, fitInfo.BGFDDApf, fitInfo.BGFDDCpf); - outputCollisionSelExtras(chFT0A, chFT0C, chFDDA, chFDDC, chFV0A, occ, ir, trs, trofs, hmpr, tfb, itsROFb, sbp, zVtxFT0vPv, vtxITSTPC, collision.rct_raw()); - outputCollsLabels(collision.globalIndex()); - - uint64_t w1[4] = {0ull, 0ull, 0ull, 0ull}; - uint64_t w2[4] = {0ull, 0ull, 0ull, 0ull}; - - if (fitCuts.saveFITbitsets() && newbc.has_foundFT0() && newbc.has_fv0a()) { - udhelpers::buildFT0FV0Words(newbc.ft0(), newbc.fv0a(), w1, w2, - fitCuts.thr1_FT0A(), - fitCuts.thr1_FT0C(), - fitCuts.thr1_FV0A(), - fitCuts.thr2_FT0A(), - fitCuts.thr2_FT0C(), - fitCuts.thr2_FV0A()); - } - - outputFITBits(w1[0], w1[1], w1[2], w1[3], - w2[0], w2[1], w2[2], w2[3]); - - if (newbc.has_zdc()) { - auto zdc = newbc.zdc(); - udZdcsReduced(outputCollisions.lastIndex(), zdc.timeZNA(), zdc.timeZNC(), zdc.energyCommonZNA(), zdc.energyCommonZNC()); - } else { - udZdcsReduced(outputCollisions.lastIndex(), -999, -999, -999, -999); - } - // update SGTracks tables - if (fillTrackTables) { - for (const auto& track : tracks) { - if (track.pt() > sameCuts.minPt() && track.eta() > sameCuts.minEta() && track.eta() < sameCuts.maxEta()) { - if (track.isPVContributor()) { - updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); - } else if (saveAllTracks) { - if (track.itsClusterSizes() && track.itsChi2NCl() > 0 && ((track.tpcNClsFindable() == 0 && savenonPVCITSOnlyTracks) || track.tpcNClsFindable() > 50)) - updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); - // if (track.isPVContributor()) updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); - } - } - } - } - // update SGFwdTracks tables - if (fillFwdTrackTables) { - if (sameCuts.withFwdTracks()) { - for (const auto& fwdtrack : fwdtracks) { - if (!sgSelector.FwdTrkSelector(fwdtrack)) - updateUDFwdTrackTables(fwdtrack, bc.globalBC()); - } - } - } - } - } - - void init(InitContext& context) - { - ccdb->setURL("http://alice-ccdb.cern.ch"); - ccdb->setCaching(true); - ccdb->setFatalWhenNull(false); - sameCuts = (SGCutParHolder)SGCuts; - fitCuts = (FITCutParHolder)FITCuts; - - // add histograms for the different process functions - histPointers.clear(); - if (context.mOptions.get("processData")) { - histPointers.insert({"reco/Stat", registry.add("reco/Stat", "Cut statistics; Selection criterion; Collisions", {HistType::kTH1F, {{14, -0.5, 13.5}}})}); - - const AxisSpec axisCountersTrg{10, 0.5, 10.5, ""}; - histPointers.insert({"reco/hCountersTrg", registry.add("reco/hCountersTrg", "Trigger counts before selections; Trigger; Counts", {HistType::kTH1F, {axisCountersTrg}})}); - histPointers.insert({"reco/hCountersTrgBcSel", registry.add("reco/hCountersTrgSel", "Trigger counts after BC selections; Trigger; Counts", {HistType::kTH1F, {axisCountersTrg}})}); - histPointers.insert({"reco/hLumi", registry.add("reco/hLumi", "Integrated luminosity before selections; Trigger; Luminosity, 1/#mub", {HistType::kTH1F, {axisCountersTrg}})}); - histPointers.insert({"reco/hLumiBcSel", registry.add("reco/hLumiBcSel", "Integrated luminosity before selections; Trigger; Luminosity, 1/#mub", {HistType::kTH1F, {axisCountersTrg}})}); - auto hCountersTrg = getHist(TH1, "reco/hCountersTrg"); - auto hCountersTrgBcSel = getHist(TH1, "reco/hCountersTrgBcSel"); - auto hLumi = getHist(TH1, "reco/hLumi"); - auto hLumiBcSel = getHist(TH1, "reco/hLumiBcSel"); - for (const auto& h : {hCountersTrg, hCountersTrgBcSel, hLumi, hLumiBcSel}) { - h->GetXaxis()->SetBinLabel(1, "TVX"); - h->GetXaxis()->SetBinLabel(2, "TCE"); - h->GetXaxis()->SetBinLabel(3, "ZEM"); - h->GetXaxis()->SetBinLabel(4, "ZNC"); - } - } - if (context.mOptions.get("processMcData")) { - histPointers.insert({"MCreco/Stat", registry.add("MCreco/Stat", "Cut statistics; Selection criterion; Collisions", {HistType::kTH1F, {{14, -0.5, 13.5}}})}); - } - - if (isGoodRCTZdc) { - myRCTChecker.init("CBT", true); - } - } - - // process function for reconstructed data - void processData(CC const& collision, BCs const& bcs, TCs const& tracks, FWs const& fwdtracks, - aod::Zdcs const& /*zdcs*/, aod::FV0As const& fv0as, aod::FT0s const& ft0s, aod::FDDs const& fdds) - { - processReco(std::string("reco"), collision, bcs, tracks, fwdtracks, fv0as, ft0s, fdds); - } - PROCESS_SWITCH(SGCandProducer, processData, "Produce UD table with data", true); - - // process function for reconstructed MC data - void processMcData(MCCC const& collision, aod::McCollisions const& /*mccollisions*/, BCs const& bcs, - TCs const& tracks, FWs const& fwdtracks, aod::Zdcs const& /*zdcs*/, aod::FV0As const& fv0as, - aod::FT0s const& ft0s, aod::FDDs const& fdds) - { - // select specific processes with the GeneratorID - if (!collision.has_mcCollision()) - return; - auto mccol = collision.mcCollision(); - if (verboseInfo) - LOGF(info, "GeneratorId %d (%d)", mccol.getGeneratorId(), generatorIds->size()); - - if (std::find(generatorIds->begin(), generatorIds->end(), mccol.getGeneratorId()) != generatorIds->end()) { - if (verboseInfo) - LOGF(info, "Event with good generatorId"); - processReco(std::string("MCreco"), collision, bcs, tracks, fwdtracks, fv0as, ft0s, fdds); - } - } - PROCESS_SWITCH(SGCandProducer, processMcData, "Produce UD tables with MC data", false); + Service ccdb; + // data inputs + using CCs = soa::Join; + using CC = CCs::iterator; + using BCs = soa::Join; + using BC = BCs::iterator; + using TCs = soa::Join; + using FWs = aod::FwdTracks; + + using MCCCs = soa::Join; + using MCCC = MCCCs::iterator; + + // get an SGCutparHolder + SGCutParHolder sameCuts = SGCutParHolder(); // SGCutparHolder + Configurable SGCuts{"SGCuts", {}, "SG event cuts"}; + FITCutParHolder fitCuts = FITCutParHolder(); + Configurable FITCuts{"FITCuts", {}, "FIT bitset cuts"}; + Configurable verboseInfo{"verboseInfo", false, "Print general info to terminal; default it false."}; + Configurable saveAllTracks{"saveAllTracks", true, "save only PV contributors or all tracks associated to a collision"}; + Configurable savenonPVCITSOnlyTracks{"savenonPVCITSOnlyTracks", false, "save non PV contributors with ITS only information"}; + Configurable rejectAtTFBoundary{"rejectAtTFBoundary", true, "reject collisions at a TF boundary"}; + Configurable noITSROFrameBorder{"noITSROFrameBorder", true, "reject ITS RO Frame Border"}; + Configurable noSameBunchPileUp{"noSameBunchPileUp", true, "reject SameBunchPileUp"}; + Configurable IsGoodVertex{"IsGoodVertex", false, "Select FT0 PV vertex matching"}; + Configurable ITSTPCVertex{"ITSTPCVertex", true, "reject ITS-only vertex"}; // if one wants to look at Single Gap pp events + Configurable> generatorIds{"generatorIds", std::vector{-1}, "MC generatorIds to process"}; + Configurable storeSG{"storeSG", true, "Store SG events in the output"}; + Configurable storeDG{"storeDG", true, "Store DG events in the output"}; + + Configurable saveFITbitsets{"saveFITbitsets", true, "Write FT0 and FV0 bitset tables to output"}; + + Configurable isGoodRCTCollision{"isGoodRCTCollision", true, "Check RCT flags for FT0,ITS,TPC and tracking"}; + Configurable isGoodRCTZdc{"isGoodRCTZdc", false, "Check RCT flags for ZDC if present in run"}; + + // Configurables to decide which tables are filled + Configurable fillTrackTables{"fillTrackTables", true, "Fill track tables"}; + Configurable fillFwdTrackTables{"fillFwdTrackTables", true, "Fill forward track tables"}; + + // SG selector + SGSelector sgSelector; + ctpRateFetcher mRateFetcher; + + // initialize RCT flag checker + RCTFlagsChecker myRCTChecker{"CBT"}; + + // data tables + Produces outputSGCollisions; + Produces outputCollisions; + Produces outputCollisionsSels; + Produces outputCollisionSelExtras; + Produces outputCollsLabels; + Produces outputZdcs; + Produces udZdcsReduced; + Produces outputTracks; + Produces outputTracksCov; + Produces outputTracksDCA; + Produces outputTracksPID; + Produces outputTracksPIDExtra; + Produces outputTracksExtra; + Produces outputTracksFlag; + Produces outputFwdTracks; + Produces outputFwdTracksExtra; + Produces outputTracksLabel; + Produces outputFITBits; + + // initialize histogram registry + HistogramRegistry registry{ + "registry", + {}}; + std::map histPointers; + + int runNumber = -1; + + // function to update UDFwdTracks, UDFwdTracksExtra + template + void updateUDFwdTrackTables(TFwdTrack const& fwdtrack, uint64_t const& bcnum) + { + outputFwdTracks(outputCollisions.lastIndex(), + fwdtrack.px(), fwdtrack.py(), fwdtrack.pz(), fwdtrack.sign(), + bcnum, fwdtrack.trackTime(), fwdtrack.trackTimeRes()); + outputFwdTracksExtra(fwdtrack.trackType(), + fwdtrack.nClusters(), + fwdtrack.pDca(), + fwdtrack.rAtAbsorberEnd(), + fwdtrack.chi2(), + fwdtrack.chi2MatchMCHMID(), + fwdtrack.chi2MatchMCHMFT(), + fwdtrack.mchBitMap(), + fwdtrack.midBitMap(), + fwdtrack.midBoards()); + } + + // function to update UDTracks, UDTracksCov, UDTracksDCA, UDTracksPID, UDTracksExtra, UDTracksFlag, + // and UDTrackCollisionIDs + template + void updateUDTrackTables(int64_t lastIndex, TTrack const& track, uint64_t const& bcnum) + { + outputTracks(lastIndex, + track.px(), track.py(), track.pz(), track.sign(), + bcnum, track.trackTime(), track.trackTimeRes()); + + // float sigmaY = track.sigmaY(); + // float sigmaZ = track.sigmaZ(); + float sigmaY = -1.; + float sigmaZ = -1.; + outputTracksCov(track.x(), track.y(), track.z(), sigmaY, sigmaZ); + + outputTracksDCA(track.dcaZ(), track.dcaXY()); + outputTracksPID(track.tpcNSigmaEl(), + track.tpcNSigmaMu(), + track.tpcNSigmaPi(), + track.tpcNSigmaKa(), + track.tpcNSigmaPr(), + track.beta(), + track.betaerror(), + track.tofNSigmaEl(), + track.tofNSigmaMu(), + track.tofNSigmaPi(), + track.tofNSigmaKa(), + track.tofNSigmaPr()); + outputTracksPIDExtra(track.tpcNSigmaDe(), + track.tpcNSigmaTr(), + track.tpcNSigmaHe(), + track.tpcNSigmaAl(), + track.tofNSigmaDe(), + track.tofNSigmaTr(), + track.tofNSigmaHe(), + track.tofNSigmaAl()); + outputTracksExtra(track.tpcInnerParam(), + track.itsClusterSizes(), + track.tpcNClsFindable(), + track.tpcNClsFindableMinusFound(), + track.tpcNClsFindableMinusCrossedRows(), + track.tpcNClsShared(), + track.trdPattern(), + track.itsChi2NCl(), + track.tpcChi2NCl(), + track.trdChi2(), + track.tofChi2(), + track.tpcSignal(), + track.tofSignal(), + track.trdSignal(), + track.length(), + track.tofExpMom(), + track.detectorMap()); + outputTracksFlag(track.has_collision(), + track.isPVContributor()); + outputTracksLabel(track.globalIndex()); + } + + // function to process reconstructed data + template + void processReco(std::string histdir, TCol const& collision, BCs const& bcs, + TCs const& tracks, FWs const& fwdtracks, + aod::FV0As const& fv0as, aod::FT0s const& ft0s, aod::FDDs const& fdds) + { + if (verboseInfo) + LOGF(debug, " collision %d", collision.globalIndex()); + getHist(TH1, histdir + "/Stat")->Fill(0., 1.); + // reject collisions at TF boundaries + if (rejectAtTFBoundary && !collision.selection_bit(aod::evsel::kNoTimeFrameBorder)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(1., 1.); + // reject collisions at ITS RO TF boundaries + if (noITSROFrameBorder && !collision.selection_bit(aod::evsel::kNoITSROFrameBorder)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(2., 1.); + // reject Same Bunch PileUp + if (noSameBunchPileUp && !collision.selection_bit(aod::evsel::kNoSameBunchPileup)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(3., 1.); + // check vertex matching to FT0 + if (IsGoodVertex && !collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(4., 1.); + // reject ITS Only vertices + if (ITSTPCVertex && !collision.selection_bit(aod::evsel::kIsVertexITSTPC)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(5., 1.); + // nominal BC + if (!collision.has_foundBC()) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(6., 1.); + // RCT CBT for collision check + if (isGoodRCTCollision && !myRCTChecker(collision)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(7., 1.); + // RCT CBT+ZDC for collision check + if (isGoodRCTZdc && !myRCTChecker(collision)) { + return; + } + getHist(TH1, histdir + "/Stat")->Fill(8., 1.); + + // + const int trs = collision.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard) ? 1 : 0; + const int trofs = collision.selection_bit(o2::aod::evsel::kNoCollInRofStandard) ? 1 : 0; + const int hmpr = collision.selection_bit(o2::aod::evsel::kNoHighMultCollInPrevRof) ? 1 : 0; + const int tfb = collision.selection_bit(o2::aod::evsel::kNoTimeFrameBorder) ? 1 : 0; + const int itsROFb = collision.selection_bit(o2::aod::evsel::kNoITSROFrameBorder) ? 1 : 0; + const int sbp = collision.selection_bit(o2::aod::evsel::kNoSameBunchPileup) ? 1 : 0; + const int zVtxFT0vPv = collision.selection_bit(o2::aod::evsel::kIsGoodZvtxFT0vsPV) ? 1 : 0; + const int vtxITSTPC = collision.selection_bit(o2::aod::evsel::kIsVertexITSTPC) ? 1 : 0; + auto bc = collision.template foundBC_as(); + double ir = 0.; + const uint64_t ts = bc.timestamp(); + const int runnumber = bc.runNumber(); + if (bc.has_zdc()) { + ir = mRateFetcher.fetch(ccdb.service, ts, runnumber, "ZNC hadronic") * 1.e-3; + } + auto newbc = bc; + + // obtain slice of compatible BCs + auto bcRange = udhelpers::compatibleBCs(collision, sameCuts.NDtcoll(), bcs, sameCuts.minNBCs()); + auto isSGEvent = sgSelector.IsSelected(sameCuts, collision, bcRange, bc); + // auto isSGEvent = sgSelector.IsSelected(sameCuts, collision, bcRange, tracks); + int issgevent = isSGEvent.value; + if (isSGEvent.bc && issgevent < 2) { + newbc = *(isSGEvent.bc); + } else { + if (verboseInfo) + LOGF(info, "No Newbc %i", bc.globalBC()); + } + getHist(TH1, histdir + "/Stat")->Fill(issgevent + 10, 1.); + if ((storeDG && issgevent == o2::aod::sgselector::DoubleGap) || (storeSG && (issgevent == o2::aod::sgselector::SingleGapA || issgevent == o2::aod::sgselector::SingleGapC))) { + if (verboseInfo) + LOGF(info, "Current BC: %i, %i, %i", bc.globalBC(), newbc.globalBC(), issgevent); + if (sameCuts.minRgtrwTOF()) { + if (udhelpers::rPVtrwTOF(tracks, collision.numContrib()) < sameCuts.minRgtrwTOF()) + return; + } + upchelpers::FITInfo fitInfo{}; + const uint8_t chFT0A = 0; + const uint8_t chFT0C = 0; + const uint8_t chFDDA = 0; + const uint8_t chFDDC = 0; + const uint8_t chFV0A = 0; + const int occ = collision.trackOccupancyInTimeRange(); + udhelpers::getFITinfo(fitInfo, newbc, bcs, ft0s, fv0as, fdds); + const int upc_flag = (collision.flags() & dataformats::Vertex>::Flags::UPCMode) ? 1 : 0; + // update SG candidates tables + outputCollisions(bc.globalBC(), bc.runNumber(), + collision.posX(), collision.posY(), collision.posZ(), upc_flag, + collision.numContrib(), udhelpers::netCharge(tracks), + 1.); // rtrwTOF); //omit the calculation to speed up the things while skimming + + outputSGCollisions(issgevent); + outputCollisionsSels(fitInfo.ampFT0A, fitInfo.ampFT0C, fitInfo.timeFT0A, fitInfo.timeFT0C, + fitInfo.triggerMaskFT0, + fitInfo.ampFDDA, fitInfo.ampFDDC, fitInfo.timeFDDA, fitInfo.timeFDDC, + fitInfo.triggerMaskFDD, + fitInfo.ampFV0A, fitInfo.timeFV0A, fitInfo.triggerMaskFV0A, + fitInfo.BBFT0Apf, fitInfo.BBFT0Cpf, fitInfo.BGFT0Apf, fitInfo.BGFT0Cpf, + fitInfo.BBFV0Apf, fitInfo.BGFV0Apf, + fitInfo.BBFDDApf, fitInfo.BBFDDCpf, fitInfo.BGFDDApf, fitInfo.BGFDDCpf); + outputCollisionSelExtras(chFT0A, chFT0C, chFDDA, chFDDC, chFV0A, occ, ir, trs, trofs, hmpr, tfb, itsROFb, sbp, zVtxFT0vPv, vtxITSTPC, collision.rct_raw()); + outputCollsLabels(collision.globalIndex()); + + uint64_t w1[4] = {0ull, 0ull, 0ull, 0ull}; + uint64_t w2[4] = {0ull, 0ull, 0ull, 0ull}; + + if (fitCuts.saveFITbitsets() && newbc.has_foundFT0() && newbc.has_fv0a()) { + udhelpers::buildFT0FV0Words(newbc.ft0(), newbc.fv0a(), w1, w2, + fitCuts.thr1_FT0A(), + fitCuts.thr1_FT0C(), + fitCuts.thr1_FV0A(), + fitCuts.thr2_FT0A(), + fitCuts.thr2_FT0C(), + fitCuts.thr2_FV0A()); + } + + outputFITBits(w1[0], w1[1], w1[2], w1[3], + w2[0], w2[1], w2[2], w2[3]); + + if (newbc.has_zdc()) { + auto zdc = newbc.zdc(); + udZdcsReduced(outputCollisions.lastIndex(), zdc.timeZNA(), zdc.timeZNC(), zdc.energyCommonZNA(), zdc.energyCommonZNC()); + } else { + udZdcsReduced(outputCollisions.lastIndex(), -999, -999, -999, -999); + } + // update SGTracks tables + if (fillTrackTables) { + for (const auto& track : tracks) { + if (track.pt() > sameCuts.minPt() && track.eta() > sameCuts.minEta() && track.eta() < sameCuts.maxEta()) { + if (track.isPVContributor()) { + updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); + } else if (saveAllTracks) { + if (track.itsClusterSizes() && track.itsChi2NCl() > 0 && ((track.tpcNClsFindable() == 0 && savenonPVCITSOnlyTracks) || track.tpcNClsFindable() > 50)) + updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); + // if (track.isPVContributor()) updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); + } + } + } + } + // update SGFwdTracks tables + if (fillFwdTrackTables) { + if (sameCuts.withFwdTracks()) { + for (const auto& fwdtrack : fwdtracks) { + if (!sgSelector.FwdTrkSelector(fwdtrack)) + updateUDFwdTrackTables(fwdtrack, bc.globalBC()); + } + } + } + } + } + + void init(InitContext& context) + { + ccdb->setURL("http://alice-ccdb.cern.ch"); + ccdb->setCaching(true); + ccdb->setFatalWhenNull(false); + sameCuts = (SGCutParHolder)SGCuts; + fitCuts = (FITCutParHolder)FITCuts; + + // add histograms for the different process functions + histPointers.clear(); + if (context.mOptions.get("processData")) { + histPointers.insert({"reco/Stat", registry.add("reco/Stat", "Cut statistics; Selection criterion; Collisions", {HistType::kTH1F, {{14, -0.5, 13.5}}})}); + + const AxisSpec axisCountersTrg{10, 0.5, 10.5, ""}; + histPointers.insert({"reco/hCountersTrg", registry.add("reco/hCountersTrg", "Trigger counts before selections; Trigger; Counts", {HistType::kTH1F, {axisCountersTrg}})}); + histPointers.insert({"reco/hCountersTrgBcSel", registry.add("reco/hCountersTrgSel", "Trigger counts after BC selections; Trigger; Counts", {HistType::kTH1F, {axisCountersTrg}})}); + histPointers.insert({"reco/hLumi", registry.add("reco/hLumi", "Integrated luminosity before selections; Trigger; Luminosity, 1/#mub", {HistType::kTH1F, {axisCountersTrg}})}); + histPointers.insert({"reco/hLumiBcSel", registry.add("reco/hLumiBcSel", "Integrated luminosity before selections; Trigger; Luminosity, 1/#mub", {HistType::kTH1F, {axisCountersTrg}})}); + auto hCountersTrg = getHist(TH1, "reco/hCountersTrg"); + auto hCountersTrgBcSel = getHist(TH1, "reco/hCountersTrgBcSel"); + auto hLumi = getHist(TH1, "reco/hLumi"); + auto hLumiBcSel = getHist(TH1, "reco/hLumiBcSel"); + for (const auto& h : {hCountersTrg, hCountersTrgBcSel, hLumi, hLumiBcSel}) { + h->GetXaxis()->SetBinLabel(1, "TVX"); + h->GetXaxis()->SetBinLabel(2, "TCE"); + h->GetXaxis()->SetBinLabel(3, "ZEM"); + h->GetXaxis()->SetBinLabel(4, "ZNC"); + } + } + if (context.mOptions.get("processMcData")) { + histPointers.insert({"MCreco/Stat", registry.add("MCreco/Stat", "Cut statistics; Selection criterion; Collisions", {HistType::kTH1F, {{14, -0.5, 13.5}}})}); + } + + if (isGoodRCTZdc) { + myRCTChecker.init("CBT", true); + } + } + + // process function for reconstructed data + void processData(CC const& collision, BCs const& bcs, TCs const& tracks, FWs const& fwdtracks, + aod::Zdcs const& /*zdcs*/, aod::FV0As const& fv0as, aod::FT0s const& ft0s, aod::FDDs const& fdds) + { + processReco(std::string("reco"), collision, bcs, tracks, fwdtracks, fv0as, ft0s, fdds); + } + PROCESS_SWITCH(SGCandProducer, processData, "Produce UD table with data", true); + + // process function for reconstructed MC data + void processMcData(MCCC const& collision, aod::McCollisions const& /*mccollisions*/, BCs const& bcs, + TCs const& tracks, FWs const& fwdtracks, aod::Zdcs const& /*zdcs*/, aod::FV0As const& fv0as, + aod::FT0s const& ft0s, aod::FDDs const& fdds) + { + // select specific processes with the GeneratorID + if (!collision.has_mcCollision()) + return; + auto mccol = collision.mcCollision(); + if (verboseInfo) + LOGF(info, "GeneratorId %d (%d)", mccol.getGeneratorId(), generatorIds->size()); + + if (std::find(generatorIds->begin(), generatorIds->end(), mccol.getGeneratorId()) != generatorIds->end()) { + if (verboseInfo) + LOGF(info, "Event with good generatorId"); + processReco(std::string("MCreco"), collision, bcs, tracks, fwdtracks, fv0as, ft0s, fdds); + } + } + PROCESS_SWITCH(SGCandProducer, processMcData, "Produce UD tables with MC data", false); }; struct McSGCandProducer { - // MC tables - Configurable verboseInfoMC{"verboseInfoMC", false, "Print general info to terminal; default it false."}; - Produces outputMcCollisions; - Produces outputMcParticles; - Produces outputMcCollsLabels; - Produces outputMcTrackLabels; - - // save all McTruth, even if the collisions is not reconstructed - Configurable> generatorIds{"generatorIds", std::vector{-1}, "MC generatorIds to process"}; - Configurable saveAllMcCollisions{"saveAllMcCollisions", true, "save all McCollisions"}; - - using CCs = soa::Join; - using BCs = soa::Join; - using TCs = soa::Join; - using UDCCs = soa::Join; - using UDTCs = soa::Join; - - // prepare slices - SliceCache cache; - PresliceUnsorted mcPartsPerMcCollision = aod::mcparticle::mcCollisionId; - Preslice udtracksPerUDCollision = aod::udtrack::udCollisionId; - - // initialize histogram registry - HistogramRegistry registry{ - "registry", - {}}; - - template - void updateUDMcCollisions(TMcCollision const& mccol, uint64_t globBC) - { - // save mccol - outputMcCollisions(globBC, - mccol.generatorsID(), - mccol.posX(), - mccol.posY(), - mccol.posZ(), - mccol.t(), - mccol.weight(), - mccol.impactParameter()); - } - - template - void updateUDMcParticle(TMcParticle const& McPart, int64_t McCollisionId, std::map& mcPartIsSaved) - { - // save McPart - // mother and daughter indices are set to -1 - // ATTENTION: this can be improved to also include mother and daughter indices - std::vector newmids; - int32_t newdids[2] = {-1, -1}; - - // update UDMcParticles - if (mcPartIsSaved.find(McPart.globalIndex()) == mcPartIsSaved.end()) { - outputMcParticles(McCollisionId, - McPart.pdgCode(), - McPart.statusCode(), - McPart.flags(), - newmids, - newdids, - McPart.weight(), - McPart.px(), - McPart.py(), - McPart.pz(), - McPart.e()); - mcPartIsSaved[McPart.globalIndex()] = outputMcParticles.lastIndex(); - } - } - - template - void updateUDMcParticles(TMcParticles const& McParts, int64_t McCollisionId, std::map& mcPartIsSaved) - { - // save McParts - // new mother and daughter ids - std::vector newmids; - int32_t newdids[2] = {-1, -1}; - int64_t newval = -1; - - // Determine the particle indices within the UDMcParticles table - // before filling the table - // This is needed to be able to assign the new daughter indices - std::map oldnew; - auto lastId = outputMcParticles.lastIndex(); - for (const auto& mcpart : McParts) { - auto oldId = mcpart.globalIndex(); - if (mcPartIsSaved.find(oldId) != mcPartIsSaved.end()) { - oldnew[oldId] = mcPartIsSaved[oldId]; - } else { - lastId++; - oldnew[oldId] = lastId; - } - } - - // all particles of the McCollision are saved - for (const auto& mcpart : McParts) { - if (mcPartIsSaved.find(mcpart.globalIndex()) == mcPartIsSaved.end()) { - // mothers - newmids.clear(); - auto oldmids = mcpart.mothersIds(); - for (const auto& oldmid : oldmids) { - auto m = McParts.rawIteratorAt(oldmid); - if (verboseInfoMC) - LOGF(debug, " m %d", m.globalIndex()); - if (mcPartIsSaved.find(oldmid) != mcPartIsSaved.end()) { - newval = mcPartIsSaved[oldmid]; - } else { - newval = -1; - } - newmids.push_back(newval); - } - // daughters - auto olddids = mcpart.daughtersIds(); - for (uint ii = 0; ii < olddids.size(); ii++) { - if (oldnew.find(olddids[ii]) != oldnew.end()) { - newval = oldnew[olddids[ii]]; - } else { - newval = -1; - } - newdids[ii] = newval; - } - if (verboseInfoMC) - LOGF(debug, " ms %i ds %i", oldmids.size(), olddids.size()); - - // update UDMcParticles - outputMcParticles(McCollisionId, - mcpart.pdgCode(), - mcpart.statusCode(), - mcpart.flags(), - newmids, - newdids, - mcpart.weight(), - mcpart.px(), - mcpart.py(), - mcpart.pz(), - mcpart.e()); - mcPartIsSaved[mcpart.globalIndex()] = outputMcParticles.lastIndex(); - } - } - } - - template - void updateUDMcTrackLabel(TTrack const& udtrack, std::map& mcPartIsSaved) - { - // udtrack (UDTCs) -> track (TCs) -> mcTrack (McParticles) -> udMcTrack (UDMcParticles) - auto trackId = udtrack.trackId(); - if (trackId >= 0) { - auto track = udtrack.template track_as(); - auto mcTrackId = track.mcParticleId(); - if (mcTrackId >= 0) { - if (mcPartIsSaved.find(mcTrackId) != mcPartIsSaved.end()) { - outputMcTrackLabels(mcPartIsSaved[mcTrackId], track.mcMask()); - } else { - outputMcTrackLabels(-1, track.mcMask()); - } - } else { - outputMcTrackLabels(-1, track.mcMask()); - } - } else { - outputMcTrackLabels(-1, -1); - } - } - - template - void updateUDMcTrackLabels(TTrack const& udtracks, std::map& mcPartIsSaved) - { - // loop over all tracks - for (const auto& udtrack : udtracks) { - // udtrack (UDTCs) -> track (TCs) -> mcTrack (McParticles) -> udMcTrack (UDMcParticles) - auto trackId = udtrack.trackId(); - if (trackId >= 0) { - auto track = udtrack.template track_as(); - auto mcTrackId = track.mcParticleId(); - if (mcTrackId >= 0) { - if (mcPartIsSaved.find(mcTrackId) != mcPartIsSaved.end()) { - outputMcTrackLabels(mcPartIsSaved[mcTrackId], track.mcMask()); - } else { - outputMcTrackLabels(-1, track.mcMask()); - } - } else { - outputMcTrackLabels(-1, track.mcMask()); - } - } else { - outputMcTrackLabels(-1, -1); - } - } - } - - // updating McTruth data and links to reconstructed data - void procWithSgCand(aod::McCollisions const& mccols, aod::McParticles const& mcparts, - UDCCs const& sgcands, UDTCs const& udtracks) - { - // use a hash table to keep track of the McCollisions which have been added to the UDMcCollision table - // {McCollisionId : udMcCollisionId} - // similar for the McParticles which have been added to the UDMcParticle table - // {McParticleId : udMcParticleId} - std::map mcColIsSaved; - std::map mcPartIsSaved; - - // loop over McCollisions and UDCCs simultaneously - auto mccol = mccols.iteratorAt(0); - auto mcOfInterest = std::find(generatorIds->begin(), generatorIds->end(), mccol.getGeneratorId()) != generatorIds->end(); - auto lastmccol = mccols.iteratorAt(mccols.size() - 1); - auto mccolAtEnd = false; - - auto sgcand = sgcands.iteratorAt(0); - auto lastsgcand = sgcands.iteratorAt(sgcands.size() - 1); - auto sgcandAtEnd = false; - - // advance dgcand and mccol until both are AtEnd - int64_t mccolId = mccol.globalIndex(); - int64_t mcsgId = -1; - bool goon = true; - while (goon) { - auto globBC = mccol.bc_as().globalBC(); - // check if dgcand has an associated McCollision - if (sgcand.has_collision()) { - auto sgcandCol = sgcand.collision_as(); - // colId = sgcandCol.globalIndex(); - if (sgcandCol.has_mcCollision()) { - mcsgId = sgcandCol.mcCollision().globalIndex(); - } else { - mcsgId = -1; - } - } else { - mcsgId = -1; - } - if (verboseInfoMC) - LOGF(info, "\nStart of loop mcsgId %d mccolId %d", mcsgId, mccolId); - - // two cases to consider - // 1. mcdgId <= mccolId: the event to process is a dgcand. In this case the Mc tables as well as the McLabel tables are updated - // 2. mccolId < mcdgId: the event to process is an MC event of interest without reconstructed dgcand. In this case only the Mc tables are updated - if ((!sgcandAtEnd && !mccolAtEnd && (mcsgId <= mccolId)) || mccolAtEnd) { - // this is case 1. - if (verboseInfoMC) - LOGF(info, "Doing case 1 with mcsgId %d", mcsgId); - - // update UDMcCollisions and UDMcColsLabels (for each UDCollision -> UDMcCollisions) - // update UDMcParticles and UDMcTrackLabels (for each UDTrack -> UDMcParticles) - // get dgcand tracks - auto sgTracks = udtracks.sliceByCached(aod::udtrack::udCollisionId, sgcand.globalIndex(), cache); - - // If the sgcand has an associated McCollision then the McCollision and all associated - // McParticles are saved - // but only consider generated events of interest - if (mcsgId >= 0 && mcOfInterest) { - if (mcColIsSaved.find(mcsgId) == mcColIsSaved.end()) { - if (verboseInfoMC) - LOGF(info, " Saving McCollision %d", mcsgId); - // update UDMcCollisions - auto sgcandMcCol = sgcand.collision_as().mcCollision(); - updateUDMcCollisions(sgcandMcCol, globBC); - mcColIsSaved[mcsgId] = outputMcCollisions.lastIndex(); - } - - // update UDMcColsLabels (for each UDCollision -> UDMcCollisions) - outputMcCollsLabels(mcColIsSaved[mcsgId]); - - // update UDMcParticles - auto mcPartsSlice = mcparts.sliceBy(mcPartsPerMcCollision, mcsgId); - updateUDMcParticles(mcPartsSlice, mcColIsSaved[mcsgId], mcPartIsSaved); - - // update UDMcTrackLabels (for each UDTrack -> UDMcParticles) - updateUDMcTrackLabels(sgTracks, mcPartIsSaved); - - } else { - // If the sgcand has no associated McCollision then only the McParticles which are associated - // with the tracks of the sgcand are saved - if (verboseInfoMC) - LOGF(info, " Saving McCollision %d", -1); - - // update UDMcColsLabels (for each UDCollision -> UDMcCollisions) - outputMcCollsLabels(-1); - - // update UDMcParticles and UDMcTrackLabels (for each UDTrack -> UDMcParticles) - // loop over tracks of dgcand - for (const auto& sgtrack : sgTracks) { - if (sgtrack.has_track()) { - auto track = sgtrack.track_as(); - if (track.has_mcParticle()) { - auto mcPart = track.mcParticle(); - auto mcCol = mcPart.mcCollision(); - if (mcColIsSaved.find(mcCol.globalIndex()) == mcColIsSaved.end()) { - updateUDMcCollisions(mcCol, globBC); - mcColIsSaved[mcCol.globalIndex()] = outputMcCollisions.lastIndex(); - } - updateUDMcParticle(mcPart, mcColIsSaved[mcCol.globalIndex()], mcPartIsSaved); - updateUDMcTrackLabel(sgtrack, mcPartIsSaved); - } else { - outputMcTrackLabels(-1, track.mcMask()); - } - } else { - outputMcTrackLabels(-1, -1); - } - } - } - // advance sgcand - if (sgcand != lastsgcand) { - sgcand++; - } else { - sgcandAtEnd = true; - } - } else { - // this is case 2. - if (verboseInfoMC) - LOGF(info, "Doing case 2"); - - // update UDMcCollisions and UDMcParticles - // but only consider generated events of interest - if (mcOfInterest && mcColIsSaved.find(mccolId) == mcColIsSaved.end()) { - if (verboseInfoMC) - LOGF(info, " Saving McCollision %d", mccolId); - // update UDMcCollisions - updateUDMcCollisions(mccol, globBC); - mcColIsSaved[mccolId] = outputMcCollisions.lastIndex(); - - // update UDMcParticles - auto mcPartsSlice = mcparts.sliceBy(mcPartsPerMcCollision, mccolId); - updateUDMcParticles(mcPartsSlice, mcColIsSaved[mccolId], mcPartIsSaved); - } - - // advance mccol - if (mccol != lastmccol) { - mccol++; - mcOfInterest = std::find(generatorIds->begin(), generatorIds->end(), mccol.getGeneratorId()) != generatorIds->end(); - mccolId = mccol.globalIndex(); - } else { - mccolAtEnd = true; - } - } - - goon = !sgcandAtEnd || !mccolAtEnd; - if (verboseInfoMC) - LOGF(info, "End of loop mcsgId %d mccolId %d", mcsgId, mccolId); - } - } - - // updating McTruth data only - void procWithoutSgCand(aod::McCollisions const& mccols, aod::McParticles const& mcparts) - { - // use a hash table to keep track of the McCollisions which have been added to the UDMcCollision table - // {McCollisionId : udMcCollisionId} - // similar for the McParticles which have been added to the UDMcParticle table - // {McParticleId : udMcParticleId} - std::map mcColIsSaved; - std::map mcPartIsSaved; - - // loop over McCollisions - for (auto const& mccol : mccols) { - int64_t mccolId = mccol.globalIndex(); - uint64_t globBC = mccol.bc_as().globalBC(); - - // update UDMcCollisions and UDMcParticles - if (mcColIsSaved.find(mccolId) == mcColIsSaved.end()) { - if (verboseInfoMC) - LOGF(info, " Saving McCollision %d", mccolId); - - // update UDMcCollisions - updateUDMcCollisions(mccol, globBC); - mcColIsSaved[mccolId] = outputMcCollisions.lastIndex(); - - // update UDMcParticles - auto mcPartsSlice = mcparts.sliceBy(mcPartsPerMcCollision, mccolId); - updateUDMcParticles(mcPartsSlice, mcColIsSaved[mccolId], mcPartIsSaved); - } - } - } - - void init(InitContext& context) - { - // add histograms for the different process functions - if (context.mOptions.get("processMC")) { - registry.add("mcTruth/collisions", "Number of associated collisions", {HistType::kTH1F, {{11, -0.5, 10.5}}}); - registry.add("mcTruth/collType", "Collision type", {HistType::kTH1F, {{5, -0.5, 4.5}}}); - registry.add("mcTruth/IVMpt", "Invariant mass versus p_{T}", {HistType::kTH2F, {{150, 0.0, 3.0}, {150, 0.0, 3.0}}}); - } - } - - // process function for MC data - // save the MC truth of all events of interest and of the DG events - void processMC(aod::McCollisions const& mccols, aod::McParticles const& mcparts, - UDCCs const& sgcands, UDTCs const& udtracks, - CCs const& /*collisions*/, BCs const& /*bcs*/, TCs const& /*tracks*/) - { - if (verboseInfoMC) { - LOGF(info, "Number of McCollisions %d", mccols.size()); - LOGF(info, "Number of SG candidates %d", sgcands.size()); - LOGF(info, "Number of UD tracks %d", udtracks.size()); - } - if (mccols.size() > 0) { - if (sgcands.size() > 0) { - procWithSgCand(mccols, mcparts, sgcands, udtracks); - } else { - if (saveAllMcCollisions) { - procWithoutSgCand(mccols, mcparts); - } - } - } - } - PROCESS_SWITCH(McSGCandProducer, processMC, "Produce MC tables", false); - - void processDummy(aod::Collisions const& /*collisions*/) - { - // do nothing - if (verboseInfoMC) - LOGF(info, "Running dummy process function!"); - } - PROCESS_SWITCH(McSGCandProducer, processDummy, "Dummy function", true); + // MC tables + Configurable verboseInfoMC{"verboseInfoMC", false, "Print general info to terminal; default it false."}; + Produces outputMcCollisions; + Produces outputMcParticles; + Produces outputMcCollsLabels; + Produces outputMcTrackLabels; + + // save all McTruth, even if the collisions is not reconstructed + Configurable> generatorIds{"generatorIds", std::vector{-1}, "MC generatorIds to process"}; + Configurable saveAllMcCollisions{"saveAllMcCollisions", true, "save all McCollisions"}; + + using CCs = soa::Join; + using BCs = soa::Join; + using TCs = soa::Join; + using UDCCs = soa::Join; + using UDTCs = soa::Join; + + // prepare slices + SliceCache cache; + PresliceUnsorted mcPartsPerMcCollision = aod::mcparticle::mcCollisionId; + Preslice udtracksPerUDCollision = aod::udtrack::udCollisionId; + + // initialize histogram registry + HistogramRegistry registry{ + "registry", + {}}; + + template + void updateUDMcCollisions(TMcCollision const& mccol, uint64_t globBC) + { + // save mccol + outputMcCollisions(globBC, + mccol.generatorsID(), + mccol.posX(), + mccol.posY(), + mccol.posZ(), + mccol.t(), + mccol.weight(), + mccol.impactParameter()); + } + + template + void updateUDMcParticle(TMcParticle const& McPart, int64_t McCollisionId, std::map& mcPartIsSaved) + { + // save McPart + // mother and daughter indices are set to -1 + // ATTENTION: this can be improved to also include mother and daughter indices + std::vector newmids; + int32_t newdids[2] = {-1, -1}; + + // update UDMcParticles + if (mcPartIsSaved.find(McPart.globalIndex()) == mcPartIsSaved.end()) { + outputMcParticles(McCollisionId, + McPart.pdgCode(), + McPart.statusCode(), + McPart.flags(), + newmids, + newdids, + McPart.weight(), + McPart.px(), + McPart.py(), + McPart.pz(), + McPart.e()); + mcPartIsSaved[McPart.globalIndex()] = outputMcParticles.lastIndex(); + } + } + + template + void updateUDMcParticles(TMcParticles const& McParts, int64_t McCollisionId, std::map& mcPartIsSaved) + { + // save McParts + // new mother and daughter ids + std::vector newmids; + int32_t newdids[2] = {-1, -1}; + int64_t newval = -1; + + // Determine the particle indices within the UDMcParticles table + // before filling the table + // This is needed to be able to assign the new daughter indices + std::map oldnew; + auto lastId = outputMcParticles.lastIndex(); + for (const auto& mcpart : McParts) { + auto oldId = mcpart.globalIndex(); + if (mcPartIsSaved.find(oldId) != mcPartIsSaved.end()) { + oldnew[oldId] = mcPartIsSaved[oldId]; + } else { + lastId++; + oldnew[oldId] = lastId; + } + } + + // all particles of the McCollision are saved + for (const auto& mcpart : McParts) { + if (mcPartIsSaved.find(mcpart.globalIndex()) == mcPartIsSaved.end()) { + // mothers + newmids.clear(); + auto oldmids = mcpart.mothersIds(); + for (const auto& oldmid : oldmids) { + auto m = McParts.rawIteratorAt(oldmid); + if (verboseInfoMC) + LOGF(debug, " m %d", m.globalIndex()); + if (mcPartIsSaved.find(oldmid) != mcPartIsSaved.end()) { + newval = mcPartIsSaved[oldmid]; + } else { + newval = -1; + } + newmids.push_back(newval); + } + // daughters + auto olddids = mcpart.daughtersIds(); + for (uint ii = 0; ii < olddids.size(); ii++) { + if (oldnew.find(olddids[ii]) != oldnew.end()) { + newval = oldnew[olddids[ii]]; + } else { + newval = -1; + } + newdids[ii] = newval; + } + if (verboseInfoMC) + LOGF(debug, " ms %i ds %i", oldmids.size(), olddids.size()); + + // update UDMcParticles + outputMcParticles(McCollisionId, + mcpart.pdgCode(), + mcpart.statusCode(), + mcpart.flags(), + newmids, + newdids, + mcpart.weight(), + mcpart.px(), + mcpart.py(), + mcpart.pz(), + mcpart.e()); + mcPartIsSaved[mcpart.globalIndex()] = outputMcParticles.lastIndex(); + } + } + } + + template + void updateUDMcTrackLabel(TTrack const& udtrack, std::map& mcPartIsSaved) + { + // udtrack (UDTCs) -> track (TCs) -> mcTrack (McParticles) -> udMcTrack (UDMcParticles) + auto trackId = udtrack.trackId(); + if (trackId >= 0) { + auto track = udtrack.template track_as(); + auto mcTrackId = track.mcParticleId(); + if (mcTrackId >= 0) { + if (mcPartIsSaved.find(mcTrackId) != mcPartIsSaved.end()) { + outputMcTrackLabels(mcPartIsSaved[mcTrackId], track.mcMask()); + } else { + outputMcTrackLabels(-1, track.mcMask()); + } + } else { + outputMcTrackLabels(-1, track.mcMask()); + } + } else { + outputMcTrackLabels(-1, -1); + } + } + + template + void updateUDMcTrackLabels(TTrack const& udtracks, std::map& mcPartIsSaved) + { + // loop over all tracks + for (const auto& udtrack : udtracks) { + // udtrack (UDTCs) -> track (TCs) -> mcTrack (McParticles) -> udMcTrack (UDMcParticles) + auto trackId = udtrack.trackId(); + if (trackId >= 0) { + auto track = udtrack.template track_as(); + auto mcTrackId = track.mcParticleId(); + if (mcTrackId >= 0) { + if (mcPartIsSaved.find(mcTrackId) != mcPartIsSaved.end()) { + outputMcTrackLabels(mcPartIsSaved[mcTrackId], track.mcMask()); + } else { + outputMcTrackLabels(-1, track.mcMask()); + } + } else { + outputMcTrackLabels(-1, track.mcMask()); + } + } else { + outputMcTrackLabels(-1, -1); + } + } + } + + // updating McTruth data and links to reconstructed data + void procWithSgCand(aod::McCollisions const& mccols, aod::McParticles const& mcparts, + UDCCs const& sgcands, UDTCs const& udtracks) + { + // use a hash table to keep track of the McCollisions which have been added to the UDMcCollision table + // {McCollisionId : udMcCollisionId} + // similar for the McParticles which have been added to the UDMcParticle table + // {McParticleId : udMcParticleId} + std::map mcColIsSaved; + std::map mcPartIsSaved; + + // loop over McCollisions and UDCCs simultaneously + auto mccol = mccols.iteratorAt(0); + auto mcOfInterest = std::find(generatorIds->begin(), generatorIds->end(), mccol.getGeneratorId()) != generatorIds->end(); + auto lastmccol = mccols.iteratorAt(mccols.size() - 1); + auto mccolAtEnd = false; + + auto sgcand = sgcands.iteratorAt(0); + auto lastsgcand = sgcands.iteratorAt(sgcands.size() - 1); + auto sgcandAtEnd = false; + + // advance dgcand and mccol until both are AtEnd + int64_t mccolId = mccol.globalIndex(); + int64_t mcsgId = -1; + bool goon = true; + while (goon) { + auto globBC = mccol.bc_as().globalBC(); + // check if dgcand has an associated McCollision + if (sgcand.has_collision()) { + auto sgcandCol = sgcand.collision_as(); + // colId = sgcandCol.globalIndex(); + if (sgcandCol.has_mcCollision()) { + mcsgId = sgcandCol.mcCollision().globalIndex(); + } else { + mcsgId = -1; + } + } else { + mcsgId = -1; + } + if (verboseInfoMC) + LOGF(info, "\nStart of loop mcsgId %d mccolId %d", mcsgId, mccolId); + + // two cases to consider + // 1. mcdgId <= mccolId: the event to process is a dgcand. In this case the Mc tables as well as the McLabel tables are updated + // 2. mccolId < mcdgId: the event to process is an MC event of interest without reconstructed dgcand. In this case only the Mc tables are updated + if ((!sgcandAtEnd && !mccolAtEnd && (mcsgId <= mccolId)) || mccolAtEnd) { + // this is case 1. + if (verboseInfoMC) + LOGF(info, "Doing case 1 with mcsgId %d", mcsgId); + + // update UDMcCollisions and UDMcColsLabels (for each UDCollision -> UDMcCollisions) + // update UDMcParticles and UDMcTrackLabels (for each UDTrack -> UDMcParticles) + // get dgcand tracks + auto sgTracks = udtracks.sliceByCached(aod::udtrack::udCollisionId, sgcand.globalIndex(), cache); + + // If the sgcand has an associated McCollision then the McCollision and all associated + // McParticles are saved + // but only consider generated events of interest + if (mcsgId >= 0 && mcOfInterest) { + if (mcColIsSaved.find(mcsgId) == mcColIsSaved.end()) { + if (verboseInfoMC) + LOGF(info, " Saving McCollision %d", mcsgId); + // update UDMcCollisions + auto sgcandMcCol = sgcand.collision_as().mcCollision(); + updateUDMcCollisions(sgcandMcCol, globBC); + mcColIsSaved[mcsgId] = outputMcCollisions.lastIndex(); + } + + // update UDMcColsLabels (for each UDCollision -> UDMcCollisions) + outputMcCollsLabels(mcColIsSaved[mcsgId]); + + // update UDMcParticles + auto mcPartsSlice = mcparts.sliceBy(mcPartsPerMcCollision, mcsgId); + updateUDMcParticles(mcPartsSlice, mcColIsSaved[mcsgId], mcPartIsSaved); + + // update UDMcTrackLabels (for each UDTrack -> UDMcParticles) + updateUDMcTrackLabels(sgTracks, mcPartIsSaved); + + } else { + // If the sgcand has no associated McCollision then only the McParticles which are associated + // with the tracks of the sgcand are saved + if (verboseInfoMC) + LOGF(info, " Saving McCollision %d", -1); + + // update UDMcColsLabels (for each UDCollision -> UDMcCollisions) + outputMcCollsLabels(-1); + + // update UDMcParticles and UDMcTrackLabels (for each UDTrack -> UDMcParticles) + // loop over tracks of dgcand + for (const auto& sgtrack : sgTracks) { + if (sgtrack.has_track()) { + auto track = sgtrack.track_as(); + if (track.has_mcParticle()) { + auto mcPart = track.mcParticle(); + auto mcCol = mcPart.mcCollision(); + if (mcColIsSaved.find(mcCol.globalIndex()) == mcColIsSaved.end()) { + updateUDMcCollisions(mcCol, globBC); + mcColIsSaved[mcCol.globalIndex()] = outputMcCollisions.lastIndex(); + } + updateUDMcParticle(mcPart, mcColIsSaved[mcCol.globalIndex()], mcPartIsSaved); + updateUDMcTrackLabel(sgtrack, mcPartIsSaved); + } else { + outputMcTrackLabels(-1, track.mcMask()); + } + } else { + outputMcTrackLabels(-1, -1); + } + } + } + // advance sgcand + if (sgcand != lastsgcand) { + sgcand++; + } else { + sgcandAtEnd = true; + } + } else { + // this is case 2. + if (verboseInfoMC) + LOGF(info, "Doing case 2"); + + // update UDMcCollisions and UDMcParticles + // but only consider generated events of interest + if (mcOfInterest && mcColIsSaved.find(mccolId) == mcColIsSaved.end()) { + if (verboseInfoMC) + LOGF(info, " Saving McCollision %d", mccolId); + // update UDMcCollisions + updateUDMcCollisions(mccol, globBC); + mcColIsSaved[mccolId] = outputMcCollisions.lastIndex(); + + // update UDMcParticles + auto mcPartsSlice = mcparts.sliceBy(mcPartsPerMcCollision, mccolId); + updateUDMcParticles(mcPartsSlice, mcColIsSaved[mccolId], mcPartIsSaved); + } + + // advance mccol + if (mccol != lastmccol) { + mccol++; + mcOfInterest = std::find(generatorIds->begin(), generatorIds->end(), mccol.getGeneratorId()) != generatorIds->end(); + mccolId = mccol.globalIndex(); + } else { + mccolAtEnd = true; + } + } + + goon = !sgcandAtEnd || !mccolAtEnd; + if (verboseInfoMC) + LOGF(info, "End of loop mcsgId %d mccolId %d", mcsgId, mccolId); + } + } + + // updating McTruth data only + void procWithoutSgCand(aod::McCollisions const& mccols, aod::McParticles const& mcparts) + { + // use a hash table to keep track of the McCollisions which have been added to the UDMcCollision table + // {McCollisionId : udMcCollisionId} + // similar for the McParticles which have been added to the UDMcParticle table + // {McParticleId : udMcParticleId} + std::map mcColIsSaved; + std::map mcPartIsSaved; + + // loop over McCollisions + for (auto const& mccol : mccols) { + int64_t mccolId = mccol.globalIndex(); + uint64_t globBC = mccol.bc_as().globalBC(); + + // update UDMcCollisions and UDMcParticles + if (mcColIsSaved.find(mccolId) == mcColIsSaved.end()) { + if (verboseInfoMC) + LOGF(info, " Saving McCollision %d", mccolId); + + // update UDMcCollisions + updateUDMcCollisions(mccol, globBC); + mcColIsSaved[mccolId] = outputMcCollisions.lastIndex(); + + // update UDMcParticles + auto mcPartsSlice = mcparts.sliceBy(mcPartsPerMcCollision, mccolId); + updateUDMcParticles(mcPartsSlice, mcColIsSaved[mccolId], mcPartIsSaved); + } + } + } + + void init(InitContext& context) + { + // add histograms for the different process functions + if (context.mOptions.get("processMC")) { + registry.add("mcTruth/collisions", "Number of associated collisions", {HistType::kTH1F, {{11, -0.5, 10.5}}}); + registry.add("mcTruth/collType", "Collision type", {HistType::kTH1F, {{5, -0.5, 4.5}}}); + registry.add("mcTruth/IVMpt", "Invariant mass versus p_{T}", {HistType::kTH2F, {{150, 0.0, 3.0}, {150, 0.0, 3.0}}}); + } + } + + // process function for MC data + // save the MC truth of all events of interest and of the DG events + void processMC(aod::McCollisions const& mccols, aod::McParticles const& mcparts, + UDCCs const& sgcands, UDTCs const& udtracks, + CCs const& /*collisions*/, BCs const& /*bcs*/, TCs const& /*tracks*/) + { + if (verboseInfoMC) { + LOGF(info, "Number of McCollisions %d", mccols.size()); + LOGF(info, "Number of SG candidates %d", sgcands.size()); + LOGF(info, "Number of UD tracks %d", udtracks.size()); + } + if (mccols.size() > 0) { + if (sgcands.size() > 0) { + procWithSgCand(mccols, mcparts, sgcands, udtracks); + } else { + if (saveAllMcCollisions) { + procWithoutSgCand(mccols, mcparts); + } + } + } + } + PROCESS_SWITCH(McSGCandProducer, processMC, "Produce MC tables", false); + + void processDummy(aod::Collisions const& /*collisions*/) + { + // do nothing + if (verboseInfoMC) + LOGF(info, "Running dummy process function!"); + } + PROCESS_SWITCH(McSGCandProducer, processDummy, "Dummy function", true); }; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { - WorkflowSpec workflow{ - adaptAnalysisTask(cfgc, TaskName{"sgcandproducer"}), - adaptAnalysisTask(cfgc, TaskName{"mcsgcandproducer"})}; - return workflow; + WorkflowSpec workflow{ + adaptAnalysisTask(cfgc, TaskName{"sgcandproducer"}), + adaptAnalysisTask(cfgc, TaskName{"mcsgcandproducer"})}; + return workflow; } From c4d6505ee811f8006389692f60b8cf35cf6899df Mon Sep 17 00:00:00 2001 From: sandor-lokos Date: Tue, 24 Mar 2026 20:35:11 +0100 Subject: [PATCH 08/11] Formatting fixes --- PWGUD/TableProducer/SGCandProducer.cxx | 20 ++++++++++---------- 1 file changed, 10 insertions(+), 10 deletions(-) diff --git a/PWGUD/TableProducer/SGCandProducer.cxx b/PWGUD/TableProducer/SGCandProducer.cxx index da3c052c036..2182197f595 100644 --- a/PWGUD/TableProducer/SGCandProducer.cxx +++ b/PWGUD/TableProducer/SGCandProducer.cxx @@ -13,9 +13,9 @@ /// \brief Produces PWGUD derived table from standard tables /// /// \author Alexander Bylinkin , University of Bergen -/// \since 23.11.2023 +/// \since 23.11.2023 /// \author Adam Matyja , INP PAN Krakow, Poland -/// \since May 2025 +/// \since May 2025 // #include "PWGUD/Core/FITCutParHolder.h" @@ -100,7 +100,7 @@ struct SGCandProducer { Configurable fillTrackTables{"fillTrackTables", true, "Fill track tables"}; Configurable fillFwdTrackTables{"fillFwdTrackTables", true, "Fill forward track tables"}; - // SG selector + // SG selector SGSelector sgSelector; ctpRateFetcher mRateFetcher; @@ -219,7 +219,7 @@ struct SGCandProducer { aod::FV0As const& fv0as, aod::FT0s const& ft0s, aod::FDDs const& fdds) { if (verboseInfo) - LOGF(debug, " collision %d", collision.globalIndex()); + LOGF(debug, " collision %d", collision.globalIndex()); getHist(TH1, histdir + "/Stat")->Fill(0., 1.); // reject collisions at TF boundaries if (rejectAtTFBoundary && !collision.selection_bit(aod::evsel::kNoTimeFrameBorder)) { @@ -357,7 +357,7 @@ struct SGCandProducer { } else if (saveAllTracks) { if (track.itsClusterSizes() && track.itsChi2NCl() > 0 && ((track.tpcNClsFindable() == 0 && savenonPVCITSOnlyTracks) || track.tpcNClsFindable() > 50)) updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); - // if (track.isPVContributor()) updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); + // if (track.isPVContributor()) updateUDTrackTables(outputCollisions.lastIndex(), track, bc.globalBC()); } } } @@ -542,7 +542,7 @@ struct McSGCandProducer { for (const auto& oldmid : oldmids) { auto m = McParts.rawIteratorAt(oldmid); if (verboseInfoMC) - LOGF(debug, " m %d", m.globalIndex()); + LOGF(debug, "m %d", m.globalIndex()); if (mcPartIsSaved.find(oldmid) != mcPartIsSaved.end()) { newval = mcPartIsSaved[oldmid]; } else { @@ -688,7 +688,7 @@ struct McSGCandProducer { if (mcsgId >= 0 && mcOfInterest) { if (mcColIsSaved.find(mcsgId) == mcColIsSaved.end()) { if (verboseInfoMC) - LOGF(info, " Saving McCollision %d", mcsgId); + LOGF(info, "Saving McCollision %d", mcsgId); // update UDMcCollisions auto sgcandMcCol = sgcand.collision_as().mcCollision(); updateUDMcCollisions(sgcandMcCol, globBC); @@ -709,7 +709,7 @@ struct McSGCandProducer { // If the sgcand has no associated McCollision then only the McParticles which are associated // with the tracks of the sgcand are saved if (verboseInfoMC) - LOGF(info, " Saving McCollision %d", -1); + LOGF(info, "Saving McCollision %d", -1); // update UDMcColsLabels (for each UDCollision -> UDMcCollisions) outputMcCollsLabels(-1); @@ -751,7 +751,7 @@ struct McSGCandProducer { // but only consider generated events of interest if (mcOfInterest && mcColIsSaved.find(mccolId) == mcColIsSaved.end()) { if (verboseInfoMC) - LOGF(info, " Saving McCollision %d", mccolId); + LOGF(info, "Saving McCollision %d", mccolId); // update UDMcCollisions updateUDMcCollisions(mccol, globBC); mcColIsSaved[mccolId] = outputMcCollisions.lastIndex(); @@ -795,7 +795,7 @@ struct McSGCandProducer { // update UDMcCollisions and UDMcParticles if (mcColIsSaved.find(mccolId) == mcColIsSaved.end()) { if (verboseInfoMC) - LOGF(info, " Saving McCollision %d", mccolId); + LOGF(info, "Saving McCollision %d", mccolId); // update UDMcCollisions updateUDMcCollisions(mccol, globBC); From 2ac21b02489f158c17d43ec0c1af3567cfaeb2f2 Mon Sep 17 00:00:00 2001 From: sandor-lokos Date: Tue, 24 Mar 2026 20:43:26 +0100 Subject: [PATCH 09/11] Formatting fixes --- PWGUD/Core/UDHelpers.h | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/PWGUD/Core/UDHelpers.h b/PWGUD/Core/UDHelpers.h index 9dbd70c6897..aefcd7739af 100644 --- a/PWGUD/Core/UDHelpers.h +++ b/PWGUD/Core/UDHelpers.h @@ -33,6 +33,7 @@ #include #include +#include // namespace with helpers for UD framework namespace udhelpers @@ -542,7 +543,7 @@ inline void setBit(uint64_t w[4], int bit, bool val) } const int word = bit >> 6; const int offs = bit & 63; - w[word] |= (uint64_t(1) << offs); + w[word] |= (static_cast(1) << offs); } template From 068cdb3035c56740fea74e953bd8de7e0ef03cf5 Mon Sep 17 00:00:00 2001 From: sandor-lokos Date: Tue, 24 Mar 2026 20:45:21 +0100 Subject: [PATCH 10/11] Formatting fixes --- PWGUD/Core/UDHelpers.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PWGUD/Core/UDHelpers.h b/PWGUD/Core/UDHelpers.h index aefcd7739af..6dafd495e39 100644 --- a/PWGUD/Core/UDHelpers.h +++ b/PWGUD/Core/UDHelpers.h @@ -31,9 +31,9 @@ #include "TLorentzVector.h" +#include #include #include -#include // namespace with helpers for UD framework namespace udhelpers From ea76f12bc642d8ee94576d77694e3612ae8c0d08 Mon Sep 17 00:00:00 2001 From: sandor-lokos Date: Tue, 24 Mar 2026 22:10:37 +0100 Subject: [PATCH 11/11] PWGUD/Tasks/CMakeLists.txt COMPONENT_NAME warning -- fixed --- PWGUD/Tasks/CMakeLists.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PWGUD/Tasks/CMakeLists.txt b/PWGUD/Tasks/CMakeLists.txt index bd33c051da1..d8aa575b906 100644 --- a/PWGUD/Tasks/CMakeLists.txt +++ b/PWGUD/Tasks/CMakeLists.txt @@ -272,4 +272,4 @@ o2physics_add_dpl_workflow(sg-exclusive-jpsi-midrapidity o2physics_add_dpl_workflow(fitbit-mapping SOURCES upcTestFITBitMapping.cxx PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore - COMPONENT_NAME PWGUD) \ No newline at end of file + COMPONENT_NAME Analysis) \ No newline at end of file