From a5d697b2975deb33a7164098f2d4a51849c15a2e Mon Sep 17 00:00:00 2001 From: Felix Schlepper Date: Mon, 2 Nov 2026 11:15:08 +0100 Subject: [PATCH 1/2] ALICE3: TRK and FT3 topology dictionaries per layer, their builder and macro Signed-off-by: Felix Schlepper --- .../ALICE3/TRKFT3/common/CMakeLists.txt | 1 + .../TRKFT3/common/macros/CMakeLists.txt | 17 + .../common/macros/CheckTopologyDictionary.C | 341 ++++++++++++++++++ .../common/reconstruction/CMakeLists.txt | 7 + .../BuildTopologyDictionary.h | 91 +++++ .../TRKReconstruction/TopologyDictionary.h | 123 +++++++ .../src/BuildTopologyDictionary.cxx | 273 ++++++++++++++ .../src/TRKReconstructionLinkDef.h | 27 ++ .../reconstruction/src/TopologyDictionary.cxx | 152 ++++++++ 9 files changed, 1032 insertions(+) create mode 100644 Detectors/Upgrades/ALICE3/TRKFT3/common/macros/CMakeLists.txt create mode 100644 Detectors/Upgrades/ALICE3/TRKFT3/common/macros/CheckTopologyDictionary.C create mode 100644 Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/include/TRKReconstruction/BuildTopologyDictionary.h create mode 100644 Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/include/TRKReconstruction/TopologyDictionary.h create mode 100644 Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/src/BuildTopologyDictionary.cxx create mode 100644 Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/src/TRKReconstructionLinkDef.h create mode 100644 Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/src/TopologyDictionary.cxx diff --git a/Detectors/Upgrades/ALICE3/TRKFT3/common/CMakeLists.txt b/Detectors/Upgrades/ALICE3/TRKFT3/common/CMakeLists.txt index 6cd471f6f74de..e6accc8b86b30 100644 --- a/Detectors/Upgrades/ALICE3/TRKFT3/common/CMakeLists.txt +++ b/Detectors/Upgrades/ALICE3/TRKFT3/common/CMakeLists.txt @@ -12,3 +12,4 @@ add_subdirectory(simulation) add_subdirectory(reconstruction) add_subdirectory(workflow) +add_subdirectory(macros) diff --git a/Detectors/Upgrades/ALICE3/TRKFT3/common/macros/CMakeLists.txt b/Detectors/Upgrades/ALICE3/TRKFT3/common/macros/CMakeLists.txt new file mode 100644 index 0000000000000..9e43453d3e138 --- /dev/null +++ b/Detectors/Upgrades/ALICE3/TRKFT3/common/macros/CMakeLists.txt @@ -0,0 +1,17 @@ +# 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. + + +o2_add_test_root_macro(CheckTopologyDictionary.C + PUBLIC_LINK_LIBRARIES O2::TRKReconstruction + O2::DetectorsBase + O2::SimulationDataFormat + LABELS alice3 COMPILE_ONLY) diff --git a/Detectors/Upgrades/ALICE3/TRKFT3/common/macros/CheckTopologyDictionary.C b/Detectors/Upgrades/ALICE3/TRKFT3/common/macros/CheckTopologyDictionary.C new file mode 100644 index 0000000000000..6d7a3fac9147b --- /dev/null +++ b/Detectors/Upgrades/ALICE3/TRKFT3/common/macros/CheckTopologyDictionary.C @@ -0,0 +1,341 @@ +// Copyright 2019-2026 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. + +/// \file CheckTopologyDictionary.C +/// \brief Builds the TRK or FT3 topology dictionary from MC timeframes and checks it on others +/// +/// Each timeframe directory holds the clusters with their MC labels +/// (o2clus_trk.root, o2clus_ft3.root), the signal hits (sgn_HitsTRK.root, +/// sgn_HitsFT3.root) and the geometry. A signal cluster is taken with its +/// particle's crossing of the chip, except at the matrix's edge (it may have +/// lost pixels) and where the particle left more clusters on the chip than it +/// crossed it (pieces of one crossing). +/// +/// With training timeframes the dictionary is built from them and written to +/// \a dictionaryFile, else read from it. On the check timeframes it prints per +/// layer the pulls residual / width without and with the offsets (median and +/// normalised MAD) and writes per entry the residuals' count, median and +/// width beside the entry's (tree "entries" in \a output). +/// +/// root -b -q 'CheckTopologyDictionary.C+("TRK", "TRKTopologyDictionary.root", {"tf7", "tf8"}, {"tf2", "tf3"})' + +#if !defined(__CLING__) || defined(__ROOTCLING__) +#include "DataFormatsTRKFT3/Cluster.h" +#include "DataFormatsTRKFT3/Hit.h" +#include "DetectorsBase/GeometryManager.h" +#include "DetectorsCommonDataFormats/DetID.h" +#include "FT3Base/GeometryTGeo.h" +#include "MathUtils/Utils.h" +#include "SimulationDataFormat/MCCompLabel.h" +#include "SimulationDataFormat/MCTruthContainer.h" +#include "TRKBase/GeometryTGeo.h" +#include "TRKBase/Specs.h" +#include "TRKReconstruction/BuildTopologyDictionary.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#endif + +namespace check_dictionary +{ + +using o2::math_utils::Point3D; + +template +auto& geometry() +{ + if constexpr (DetID == o2::detectors::DetID::TRK) { + return *o2::trk::GeometryTGeo::Instance(); + } else { + return *o2::ft3::GeometryTGeo::Instance(); + } +} + +/// The median and normalised MAD of \a values, which it reorders. +std::pair medianAndWidth(std::vector& values) +{ + const auto median = [](std::vector& v) { + const auto middle = v.begin() + static_cast(v.size() / 2); + std::nth_element(v.begin(), middle, v.end()); + return v.size() % 2 == 1 ? *middle : 0.5f * (*std::max_element(v.begin(), middle) + *middle); + }; + if (values.empty()) { + return {0.f, 0.f}; + } + const float m = median(values); + std::vector deviations(values.size()); + std::ranges::transform(values, deviations.begin(), [m](float v) { return std::abs(v - m); }); + return {m, 1.4826f * median(deviations)}; +} + +/// The o2sim tree of \a path, held open by \a file. +TTree& treeOf(std::unique_ptr& file, const std::string& path) +{ + file.reset(TFile::Open(path.c_str())); + auto* tree = file && !file->IsZombie() ? file->Get("o2sim") : nullptr; + if (tree == nullptr) { + throw std::runtime_error(fmt::format("{}: no o2sim tree", path)); + } + return *tree; +} + +/// Reads entry 0 of \a tree's branch \a name into \a object. +template +void readBranch(TTree& tree, const std::string& name, T& object) +{ + auto* branch = tree.GetBranch(name.c_str()); + if (branch == nullptr) { + throw std::runtime_error(fmt::format("no branch {}", name)); + } + auto* pointer = &object; + branch->SetAddress(&pointer); + branch->GetEntry(0); + branch->ResetAddress(); +} + +/// Calls \a f(cluster, pattern, start, end) for the selected signal clusters +/// of \a directory, with their particle's crossing in O2's local frame. +template +void forEachCluster(const std::string& directory, F&& f) +{ + using Dictionary = o2::trkft3::TopologyDictionary; + using ClusterType = o2::trkft3::Cluster; + const std::string name = o2::detectors::DetID(DetID).getName(); + struct Segment { + int track, chip; + Point3D start, end; ///< global frame + }; + const auto keyOf = [](const Segment& segment) { return std::pair{segment.track, segment.chip}; }; + + // The signal hits, by event, sorted by particle and chip. + std::vector> events; + { + std::unique_ptr file; + auto& tree = treeOf(file, fmt::format("{}/sgn_Hits{}.root", directory, name)); + std::vector hits, *hitsPointer = &hits; + tree.SetBranchAddress(fmt::format("{}Hit", name).c_str(), &hitsPointer); + events.resize(tree.GetEntries()); + for (std::size_t event = 0; event < events.size(); ++event) { + tree.GetEntry(static_cast(event)); + for (const auto& hit : hits) { + events[event].push_back({hit.GetTrackID(), hit.GetDetectorID(), {hit.GetStartX(), hit.GetStartY(), hit.GetStartZ()}, {hit.GetX(), hit.GetY(), hit.GetZ()}}); + } + std::ranges::sort(events[event], {}, keyOf); + } + tree.ResetBranchAddresses(); + } + + std::unique_ptr file; + auto& tree = treeOf(file, fmt::format("{}/o2clus_{}.root", directory, DetID == o2::detectors::DetID::TRK ? "trk" : "ft3")); + for (int stream = 0; tree.GetBranch(fmt::format("{}ClusterComp_{}", name, stream).c_str()) != nullptr; ++stream) { + std::vector clusters; + std::vector patterns; + o2::dataformats::MCTruthContainer labels; + readBranch(tree, fmt::format("{}ClusterComp_{}", name, stream), clusters); + readBranch(tree, fmt::format("{}ClusterPatt_{}", name, stream), patterns); + readBranch(tree, fmt::format("{}ClusterMCTruth_{}", name, stream), labels); + + // The signal clusters with their particle's hits on the chip; the clusters of each particle and chip. + struct Candidate { + const ClusterType* cluster; + const unsigned char* pattern; + std::span segments; + std::tuple key; ///< event, track, chip + }; + std::vector candidates; + std::map, int> nClusters; + const unsigned char* next = patterns.data(); + for (std::size_t i = 0; i < clusters.size(); ++i) { + const auto& cluster = clusters[i]; + const auto* pattern = next; + next += 2 + ((pattern[0] * pattern[1] + 7) / 8); + const auto ofCluster = labels.getLabels(i); + const auto signal = std::find_if(ofCluster.begin(), ofCluster.end(), [](const auto& label) { return label.isValid() && !label.isNoise() && label.getSourceID() == 0; }); + if (signal == ofCluster.end()) { + continue; + } + const auto& segments = events.at(signal->getEventID()); + const auto range = std::ranges::equal_range(segments, std::pair{signal->getTrackID(), static_cast(cluster.chipID)}, {}, keyOf); + if (range.empty()) { + continue; + } + candidates.push_back({&cluster, pattern, {range.begin(), range.end()}, {signal->getEventID(), signal->getTrackID(), cluster.chipID}}); + ++nClusters[candidates.back().key]; + } + + for (const auto& candidate : candidates) { + const auto& cluster = *candidate.cluster; + const auto* pattern = candidate.pattern; + // Left out: at the matrix's edge, or one of several pieces of a crossing. + const bool vd = cluster.subDetID == 0; + const int nRows = vd ? o2::trk::constants::VD::petal::layer::nRows[cluster.layer] : o2::trk::constants::moduleMLOT::chip::nRows; + const int nCols = vd ? o2::trk::constants::VD::petal::layer::nCols : o2::trk::constants::moduleMLOT::chip::nCols; + if (cluster.row == 0 || cluster.col == 0 || cluster.row + pattern[0] >= nRows || cluster.col + pattern[1] >= nCols || + nClusters[candidate.key] > static_cast(candidate.segments.size())) { + continue; + } + // A particle crossing the chip more than once: the crossing whose middle is nearest the centre of gravity. + const auto& l2g = geometry().getMatrixL2G(cluster.chipID); + const auto cog = Dictionary::centreOfGravity(cluster, pattern); + float nearest = std::numeric_limits::max(); + Point3D start, end; + for (const auto& segment : candidate.segments) { + const auto a = l2g ^ segment.start, b = l2g ^ segment.end; + const auto fa = Dictionary::localToFlat(cluster, a), fb = Dictionary::localToFlat(cluster, b); + const float distance = std::hypot(0.5f * (fa.X() + fb.X()) - cog.X(), 0.5f * (fa.Z() + fb.Z()) - cog.Z()); + if (distance < nearest) { + nearest = distance; + start = a; + end = b; + } + } + f(cluster, pattern, start, end); + } + } + LOGP(info, "{}: read", directory); +} + +template +void check(const std::string& dictionaryFile, const std::vector& timeframes, const std::vector& training, const std::string& output, + bool applyOffsets, int nThreads) +{ + using Dictionary = o2::trkft3::TopologyDictionary; + const auto& directories = training.empty() ? timeframes : training; + if (directories.empty()) { + throw std::invalid_argument("no timeframe directories"); + } + o2::base::GeometryManager::loadGeometry(directories.front() + "/o2sim_geometry-aligned.root"); + geometry().fillMatrixCache(o2::math_utils::bit2Mask(o2::math_utils::TransformType::L2G)); + + std::unique_ptr dictionary; + if (!training.empty()) { + int nLayers = 0; + if constexpr (DetID == o2::detectors::DetID::TRK) { + nLayers = o2::trk::constants::VD::petal::nLayers + geometry().getNumberOfLayersMLOT(); + } else { + nLayers = geometry().getNumberOfLayers(); + } + o2::trkft3::BuildTopologyDictionary builder(nLayers, {.nThreads = nThreads}); + for (const auto& directory : training) { + forEachCluster(directory, [&](const auto& cluster, const auto* pattern, const auto& start, const auto& end) { builder.accountCluster(cluster, pattern, start, end); }); + } + builder.build(); + builder.setApplyOffsets(applyOffsets); + builder.saveDictionaryRoot(dictionaryFile); + dictionary = std::make_unique(builder.getDictionary()); + LOGP(info, "written to {}", dictionaryFile); + } else { + dictionary.reset(Dictionary::loadFrom(dictionaryFile)); + } + dictionary->print(); + if (timeframes.empty()) { + return; + } + const bool applied = dictionary->getApplyOffsets(); + + // The residuals hit - centre of gravity at the layer's depth, by layer and entry, across and along. + const int nLayers = dictionary->getNumberOfLayers(); + dictionary->setApplyOffsets(false); + std::vector, 2>>> residuals(nLayers); + for (const auto& directory : timeframes) { + forEachCluster(directory, [&](const auto& cluster, const auto* pattern, const auto& start, const auto& end) { + const int layer = Dictionary::layerOf(cluster); + if (layer < 0 || layer >= nLayers) { + throw std::runtime_error(fmt::format("{}: on no layer of the dictionary's {}", cluster.asString(), nLayers)); + } + const int id = dictionary->findID(layer, pattern); + const auto position = dictionary->getFlatCoordinates(cluster, pattern, id); + const auto a = Dictionary::localToFlat(cluster, start), b = Dictionary::localToFlat(cluster, end); + const float t = b.Y() != a.Y() ? (position.Y() - a.Y()) / (b.Y() - a.Y()) : -1.f; + if (t >= 0.f && t <= 1.f) { + auto& entry = residuals[layer][id]; + entry[0].push_back(a.X() + (t * (b.X() - a.X())) - position.X()); + entry[1].push_back(a.Z() + (t * (b.Z() - a.Z())) - position.Z()); + } + }); + } + + const std::unique_ptr out{TFile::Open(output.c_str(), "RECREATE")}; + auto* entries = new TTree("entries", "per layer and entry, across (x) and along (z): the residuals against the entry [cm]"); + int tLayer = 0, tId = 0, tCount = 0; + bool tGroup = false; + std::array tMedian{}, tWidth{}, tEntryOffset{}, tEntryWidth{}; + entries->Branch("layer", &tLayer); + entries->Branch("id", &tId); + entries->Branch("isGroup", &tGroup); + entries->Branch("count", &tCount); + for (int axis = 0; axis < 2; ++axis) { + const char* a = axis == 0 ? "x" : "z"; + entries->Branch(fmt::format("{}Median", a).c_str(), &tMedian[axis]); + entries->Branch(fmt::format("{}Width", a).c_str(), &tWidth[axis]); + entries->Branch(fmt::format("{}EntryOffset", a).c_str(), &tEntryOffset[axis]); + entries->Branch(fmt::format("{}EntryWidth", a).c_str(), &tEntryWidth[axis]); + } + LOGP(info, "pulls residual / width per layer, median and normalised MAD, without | with the offsets (the dictionary's: {})", applied ? "applied" : "not applied"); + for (int layer = 0; layer < nLayers; ++layer) { + std::array, 2> pulls, offsetPulls; + for (auto& [id, values] : residuals[layer]) { + const auto& entry = dictionary->getEntry(layer, id); + const std::array width{entry.xWidth, entry.zWidth}, offset{entry.xOffset, entry.zOffset}; + for (int axis = 0; axis < 2; ++axis) { + for (const float r : values[axis]) { + pulls[axis].push_back(r / width[axis]); + offsetPulls[axis].push_back((r - offset[axis]) / width[axis]); + } + std::tie(tMedian[axis], tWidth[axis]) = medianAndWidth(values[axis]); + tEntryOffset[axis] = offset[axis]; + tEntryWidth[axis] = width[axis]; + } + tLayer = layer; + tId = id; + tGroup = entry.isGroup; + tCount = static_cast(values[0].size()); + entries->Fill(); + } + const auto [x, xWidth] = medianAndWidth(pulls[0]); + const auto [z, zWidth] = medianAndWidth(pulls[1]); + const auto [xo, xoWidth] = medianAndWidth(offsetPulls[0]); + const auto [zo, zoWidth] = medianAndWidth(offsetPulls[1]); + LOGP(info, " layer {:2}, {:9} residuals: across {:6.3f} {:5.3f} | {:6.3f} {:5.3f}, along {:6.3f} {:5.3f} | {:6.3f} {:5.3f}", layer, pulls[0].size(), x, xWidth, xo, + xoWidth, z, zWidth, zo, zoWidth); + } + out->Write(); + LOGP(info, "written to {}", output); +} + +} // namespace check_dictionary + +/// \a detector "TRK" or "FT3"; \a timeframes the directories to check on, \a training those to build from (none: \a dictionaryFile is read); +/// \a applyOffsets whether a built dictionary adds the entries' median offsets; \a nThreads the depth scan's. +void CheckTopologyDictionary(const std::string& detector, const std::string& dictionaryFile, const std::vector& timeframes, + const std::vector& training = {}, const std::string& output = "CheckTopologyDictionary.root", bool applyOffsets = false, + int nThreads = 8) +{ + if (detector == "TRK") { + check_dictionary::check(dictionaryFile, timeframes, training, output, applyOffsets, nThreads); + } else if (detector == "FT3") { + check_dictionary::check(dictionaryFile, timeframes, training, output, applyOffsets, nThreads); + } else { + LOGP(fatal, "unknown detector {}: TRK or FT3", detector); + } +} diff --git a/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/CMakeLists.txt b/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/CMakeLists.txt index b2f8a954888fc..0dcbe5589e23d 100644 --- a/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/CMakeLists.txt +++ b/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/CMakeLists.txt @@ -16,6 +16,8 @@ endif() o2_add_library(TRKReconstruction TARGETVARNAME targetName SOURCES src/Clusterer.cxx + src/TopologyDictionary.cxx + src/BuildTopologyDictionary.cxx $<$:src/ClustererACTS.cxx> PUBLIC_LINK_LIBRARIES Microsoft.GSL::GSL @@ -30,3 +32,8 @@ o2_add_library(TRKReconstruction if(Acts_FOUND) target_compile_definitions(${targetName} PUBLIC O2_WITH_ACTS) endif() + +o2_target_root_dictionary(TRKReconstruction + HEADERS include/TRKReconstruction/TopologyDictionary.h + include/TRKReconstruction/BuildTopologyDictionary.h + LINKDEF src/TRKReconstructionLinkDef.h) diff --git a/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/include/TRKReconstruction/BuildTopologyDictionary.h b/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/include/TRKReconstruction/BuildTopologyDictionary.h new file mode 100644 index 0000000000000..8bbdbee378f83 --- /dev/null +++ b/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/include/TRKReconstruction/BuildTopologyDictionary.h @@ -0,0 +1,91 @@ +// Copyright 2019-2026 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. + +/// \file BuildTopologyDictionary.h +/// \brief Builds the TRK and FT3 topology dictionaries from MC clusters and their particles' crossings +/// +/// Each cluster is accounted with its particle's crossing of the chip: the +/// ends of its path through the silicon. Per layer, build first scans the +/// depth: the local y, on a grid of depthStep, at which the residuals' +/// widths across and along, added in quadrature, are smallest. It then cuts +/// every crossing at that depth: a pattern with at least minCount residuals +/// gets its own entry; a group takes the widths of its rare patterns, else of +/// the rare patterns with the same span class on the axis, else of the layer. + +#ifndef ALICEO2_TRKFT3_BUILDTOPOLOGYDICTIONARY_H +#define ALICEO2_TRKFT3_BUILDTOPOLOGYDICTIONARY_H + +#include "TRKReconstruction/TopologyDictionary.h" +#include +#include +#include +#include + +namespace o2::trkft3 +{ + +template +class BuildTopologyDictionary +{ + public: + using Dictionary = TopologyDictionary; + using ClusterType = typename Dictionary::ClusterType; + + struct Settings { + uint32_t minCount = 137; ///< residuals for an entry of its own: the width to 10 % + float depthStep = 0.5e-4f; ///< [cm] + int nThreads = 1; ///< for the depth scan + }; + + BuildTopologyDictionary(int nLayers, const Settings& settings); + + /// A cluster and its particle's crossing of the chip, in O2's local frame [cm]. + void accountCluster(const ClusterType& cluster, const unsigned char* pattern, const math_utils::Point3D& start, const math_utils::Point3D& end); + + /// Builds every layer; throws std::runtime_error for a layer with too few residuals. + void build(); + + void setApplyOffsets(bool apply) noexcept { mDictionary.setApplyOffsets(apply); } + [[nodiscard]] const Dictionary& getDictionary() const noexcept { return mDictionary; } + /// Writes the dictionary as "ccdb_object". + void saveDictionaryRoot(const std::string& fileName) const; + + private: + /// A cluster's crossing relative to its centre of gravity in the layer's flat frame [cm]. + struct Crossing { + uint32_t pattern = 0; ///< index among the layer's patterns + float x0 = 0.f; ///< start - centre of gravity, across + float y0 = 0.f; ///< start, depth + float z0 = 0.f; ///< start - centre of gravity, along + float dx = 0.f; ///< end - start + float dy = 0.f; + float dz = 0.f; + }; + struct LayerInput { + std::unordered_map patternIndex; ///< the pattern's bytes -> index + std::vector patterns; + std::vector crossings; + }; + + [[nodiscard]] float measureDepth(int layer, const LayerInput& input) const; + void makeEntries(int layer, const LayerInput& input, LayerDictionary& data) const; + + Settings mSettings; //! + std::vector mInput; //! + Dictionary mDictionary; //! +}; + +using TRKBuildTopologyDictionary = BuildTopologyDictionary; +using FT3BuildTopologyDictionary = BuildTopologyDictionary; + +} // namespace o2::trkft3 + +#endif diff --git a/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/include/TRKReconstruction/TopologyDictionary.h b/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/include/TRKReconstruction/TopologyDictionary.h new file mode 100644 index 0000000000000..3bb2b9983203a --- /dev/null +++ b/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/include/TRKReconstruction/TopologyDictionary.h @@ -0,0 +1,123 @@ +// Copyright 2019-2026 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. + +/// \file TopologyDictionary.h +/// \brief Positions and errors of TRK and FT3 clusters per layer and pixel pattern +/// +/// As the ITS's topology dictionary, measured from MC, but per layer: the +/// layer's depth (the local y positions are taken at) and, per common pattern +/// and per group of rare ones, the width (normalised median absolute +/// deviation) and median of the residuals hit - centre of gravity across +/// (local x) and along (local z). Every pattern finds an entry: its own, else +/// its group's. On the VD the dictionary works in the flat frame of +/// SegmentationChip::curvedToFlat; getClusterCoordinates returns O2's curved +/// local coordinates. + +#ifndef ALICEO2_TRKFT3_TOPOLOGYDICTIONARY_H +#define ALICEO2_TRKFT3_TOPOLOGYDICTIONARY_H + +#include "DataFormatsTRKFT3/Cluster.h" +#include "DetectorsCommonDataFormats/DetID.h" +#include "MathUtils/Cartesian.h" +#include +#include +#include +#include +#include + +namespace o2::trkft3 +{ + +template +class BuildTopologyDictionary; + +/// A common pattern, or the group of the rare patterns with the same span +/// classes (spans 1 to 5, 6 or more) on both axes. +struct DictionaryEntry { + std::vector pattern; ///< row span, column span, pixel bits as the clusteriser streams them; empty for a group + float xCOG = 0.f; ///< centre of gravity from the reference pixel across [pitches] + float zCOG = 0.f; ///< along [pitches] + float xWidth = 0.f; ///< the residuals' width across [cm] + float zWidth = 0.f; ///< along [cm] + float xOffset = 0.f; ///< the residuals' median across [cm] + float zOffset = 0.f; ///< along [cm] + bool isGroup = false; + + ClassDefNV(DictionaryEntry, 1); +}; + +struct LayerDictionary { + float depth = 0.f; ///< the local y positions are taken at [cm] + std::vector entries; ///< the common patterns, most frequent first, then the groups + std::unordered_multimap hashMap; ///< itsmft::ClusterTopology hash -> common pattern's entry + + ClassDefNV(LayerDictionary, 1); +}; + +template +class TopologyDictionary +{ + static_assert(DetID == o2::detectors::DetID::TRK || DetID == o2::detectors::DetID::FT3, "only TRK and FT3 dictionaries are supported"); + + public: + using ClusterType = Cluster; + + static constexpr int NSpanClasses = 6; ///< spans 1 to 5, 6 or more + static constexpr int NGroups = NSpanClasses * NSpanClasses; ///< one group per pair of span classes, row class major + + static constexpr int spanClass(int span) noexcept { return std::min(span, NSpanClasses) - 1; } + static constexpr int groupOf(int rowSpan, int colSpan) noexcept { return spanClass(rowSpan) * NSpanClasses + spanClass(colSpan); } + + /// A cluster's layer: TRK's VD layers, then ML/OT (GeometryTGeo::getLayerTRK); FT3's discs. -1 for none (a VD disc). + static int layerOf(const ClusterType& cluster) noexcept; + /// A point from O2's local frame to the layer's flat frame (the VD's unrolled; elsewhere the same). + static math_utils::Point3D localToFlat(const ClusterType& cluster, const math_utils::Point3D& local) noexcept; + /// The cluster's centre of gravity in the flat frame at y = 0. \a pattern: row span, column span, pixel bits. + static math_utils::Point3D centreOfGravity(const ClusterType& cluster, const unsigned char* pattern) noexcept; + + [[nodiscard]] int getNumberOfLayers() const noexcept { return static_cast(mLayers.size()); } + [[nodiscard]] const LayerDictionary& getLayer(int layer) const { return mLayers.at(layer); } + [[nodiscard]] const DictionaryEntry& getEntry(int layer, int id) const { return mLayers.at(layer).entries.at(id); } + + /// The pattern's entry on the layer: its own if it is common, else its group's. + [[nodiscard]] int findID(int layer, const unsigned char* pattern) const; + + /// The position in the flat frame: the entry's centre of gravity, plus its offsets if applied, at the layer's depth. + [[nodiscard]] math_utils::Point3D getFlatCoordinates(const ClusterType& cluster, const unsigned char* pattern, int id) const; + /// The same in O2's local frame: curved on the VD. + [[nodiscard]] math_utils::Point3D getClusterCoordinates(const ClusterType& cluster, const unsigned char* pattern, int id) const; + [[nodiscard]] math_utils::Point3D getClusterCoordinates(const ClusterType& cluster, const unsigned char* pattern) const + { + return getClusterCoordinates(cluster, pattern, findID(layerOf(cluster), pattern)); + } + + [[nodiscard]] bool getApplyOffsets() const noexcept { return mApplyOffsets; } + void setApplyOffsets(bool apply) noexcept { mApplyOffsets = apply; } + + void print() const; + + static TopologyDictionary* loadFrom(const std::string& fileName, const std::string& objName = "ccdb_object"); + + private: + friend class BuildTopologyDictionary; + + std::vector mLayers; + bool mApplyOffsets = false; + + ClassDefNV(TopologyDictionary, 1); +}; + +using TRKTopologyDictionary = TopologyDictionary; +using FT3TopologyDictionary = TopologyDictionary; + +} // namespace o2::trkft3 + +#endif diff --git a/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/src/BuildTopologyDictionary.cxx b/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/src/BuildTopologyDictionary.cxx new file mode 100644 index 0000000000000..18910fa9adac6 --- /dev/null +++ b/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/src/BuildTopologyDictionary.cxx @@ -0,0 +1,273 @@ +// Copyright 2019-2026 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. + +/// \file BuildTopologyDictionary.cxx + +#include "TRKReconstruction/BuildTopologyDictionary.h" +#include "DataFormatsITSMFT/ClusterPattern.h" +#include "DataFormatsITSMFT/ClusterTopology.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +namespace o2::trkft3 +{ + +namespace +{ + +/// The median of \a values, which it reorders. +float medianOf(std::span values) +{ + const auto middle = values.begin() + static_cast(values.size() / 2); + std::nth_element(values.begin(), middle, values.end()); + if (values.size() % 2 == 1) { + return *middle; + } + return 0.5f * (*std::max_element(values.begin(), middle) + *middle); +} + +/// The median of \a values and their normalised median absolute deviation +/// (a Gaussian's sigma); \a deviations a buffer at least their size. +std::pair medianAndWidth(std::span values, std::span deviations) +{ + const float median = medianOf(values); + for (std::size_t i = 0; i < values.size(); ++i) { + deviations[i] = std::abs(values[i] - median); + } + return {median, 1.4826f * medianOf(deviations.first(values.size()))}; +} + +std::pair medianAndWidth(std::vector& values) +{ + std::vector deviations(values.size()); + return medianAndWidth(values, deviations); +} + +} // namespace + +template +BuildTopologyDictionary::BuildTopologyDictionary(int nLayers, const Settings& settings) : mSettings(settings), mInput(nLayers) +{ + mDictionary.mLayers.resize(nLayers); +} + +template +void BuildTopologyDictionary::accountCluster(const ClusterType& cluster, const unsigned char* pattern, const math_utils::Point3D& start, + const math_utils::Point3D& end) +{ + const int layer = Dictionary::layerOf(cluster); + if (layer < 0 || layer >= static_cast(mInput.size())) { + throw std::runtime_error(fmt::format("{}: on no layer of the dictionary's {}", cluster.asString(), mInput.size())); + } + auto& input = mInput[layer]; + const auto nBytes = static_cast(2 + ((pattern[0] * pattern[1] + 7) / 8)); + std::string bytes(reinterpret_cast(pattern), nBytes); + const auto [index, inserted] = input.patternIndex.try_emplace(bytes, static_cast(input.patterns.size())); + if (inserted) { + input.patterns.push_back(std::move(bytes)); + } + const auto cog = Dictionary::centreOfGravity(cluster, pattern); + const auto a = Dictionary::localToFlat(cluster, start); + const auto b = Dictionary::localToFlat(cluster, end); + input.crossings.push_back({.pattern = index->second, + .x0 = a.X() - cog.X(), + .y0 = a.Y(), + .z0 = a.Z() - cog.Z(), + .dx = b.X() - a.X(), + .dy = b.Y() - a.Y(), + .dz = b.Z() - a.Z()}); +} + +template +float BuildTopologyDictionary::measureDepth(int layer, const LayerInput& input) const +{ + const auto& crossings = input.crossings; + float low = std::numeric_limits::max(), high = std::numeric_limits::lowest(); + for (const auto& crossing : crossings) { + low = std::min({low, crossing.y0, crossing.y0 + crossing.dy}); + high = std::max({high, crossing.y0, crossing.y0 + crossing.dy}); + } + const float step = mSettings.depthStep; + const int nSteps = static_cast(std::floor((high - low) / step)); + if (nSteps < 1) { + throw std::runtime_error(fmt::format("layer {}: its hits span {} um of silicon, less than the depth step", layer, 1e4 * (high - low))); + } + // The grid's points sit half a step inside the cells; the crossings that reach all of them are scanned. + const auto depthAt = [low, step](int k) { return low + ((static_cast(k) + 0.5f) * step); }; + const float first = depthAt(0), last = depthAt(nSteps - 1); + std::vector scanned; + for (std::size_t i = 0; i < crossings.size(); ++i) { + const auto& crossing = crossings[i]; + if (crossing.dy != 0.f && std::min(crossing.y0, crossing.y0 + crossing.dy) <= first && std::max(crossing.y0, crossing.y0 + crossing.dy) >= last) { + scanned.push_back(static_cast(i)); + } + } + if (scanned.size() < mSettings.minCount) { + throw std::runtime_error(fmt::format("layer {}: {} crossings through the silicon, fewer than {}", layer, scanned.size(), mSettings.minCount)); + } + + std::vector objective(nSteps); + const auto scan = [&](int thread, int nThreads) { + std::vector residualsX(scanned.size()), residualsZ(scanned.size()), deviations(scanned.size()); + for (int k = thread; k < nSteps; k += nThreads) { + const float depth = depthAt(k); + for (std::size_t i = 0; i < scanned.size(); ++i) { + const auto& crossing = crossings[scanned[i]]; + const float t = (depth - crossing.y0) / crossing.dy; + residualsX[i] = crossing.x0 + (t * crossing.dx); + residualsZ[i] = crossing.z0 + (t * crossing.dz); + } + objective[k] = std::hypot(medianAndWidth(residualsX, deviations).second, medianAndWidth(residualsZ, deviations).second); + } + }; + const int nThreads = std::clamp(mSettings.nThreads, 1, nSteps); + std::vector threads; + for (int thread = 1; thread < nThreads; ++thread) { + threads.emplace_back(scan, thread, nThreads); + } + scan(0, nThreads); + for (auto& thread : threads) { + thread.join(); + } + return depthAt(static_cast(std::ranges::min_element(objective) - objective.begin())); +} + +template +void BuildTopologyDictionary::makeEntries(int layer, const LayerInput& input, LayerDictionary& data) const +{ + constexpr int NSpanClasses = Dictionary::NSpanClasses; + const auto nPatterns = input.patterns.size(); + const auto spansOf = [&input](uint32_t pattern) { return std::pair{static_cast(input.patterns[pattern][0]), static_cast(input.patterns[pattern][1])}; }; + + // The residuals at the depth, by pattern. + std::vector> residualsX(nPatterns), residualsZ(nPatterns); + std::vector layerX, layerZ; + for (const auto& crossing : input.crossings) { + const float t = crossing.dy != 0.f ? (data.depth - crossing.y0) / crossing.dy : -1.f; + if (!(t >= 0.f && t <= 1.f)) { + continue; + } + const float x = crossing.x0 + (t * crossing.dx), z = crossing.z0 + (t * crossing.dz); + residualsX[crossing.pattern].push_back(x); + residualsZ[crossing.pattern].push_back(z); + layerX.push_back(x); + layerZ.push_back(z); + } + if (layerX.size() < mSettings.minCount) { + throw std::runtime_error(fmt::format("layer {}: {} residuals at its depth, fewer than {}", layer, layerX.size(), mSettings.minCount)); + } + + // The common patterns, most frequent first (equal counts by their bytes). + std::vector common; + for (uint32_t pattern = 0; pattern < nPatterns; ++pattern) { + if (residualsX[pattern].size() >= mSettings.minCount) { + common.push_back(pattern); + } + } + std::ranges::sort(common, [&](uint32_t a, uint32_t b) { + return residualsX[a].size() != residualsX[b].size() ? residualsX[a].size() > residualsX[b].size() : input.patterns[a] < input.patterns[b]; + }); + std::vector isCommon(nPatterns, false); + for (const auto pattern : common) { + isCommon[pattern] = true; + const auto* bytes = reinterpret_cast(input.patterns[pattern].data()); + const auto [rowSpan, colSpan] = spansOf(pattern); + DictionaryEntry entry; + entry.pattern.assign(bytes, bytes + input.patterns[pattern].size()); + float cogRow = 0.f, cogCol = 0.f; + itsmft::ClusterPattern::getCOG(rowSpan, colSpan, bytes + 2, cogRow, cogCol); + entry.xCOG = -cogRow; + entry.zCOG = cogCol; + std::tie(entry.xOffset, entry.xWidth) = medianAndWidth(residualsX[pattern]); + std::tie(entry.zOffset, entry.zWidth) = medianAndWidth(residualsZ[pattern]); + data.hashMap.emplace(itsmft::ClusterTopology::getCompleteHash(rowSpan, colSpan, bytes + 2), static_cast(data.entries.size())); + data.entries.push_back(std::move(entry)); + } + + // The rare patterns pooled by group and by span class on each axis. + std::vector> groupX(Dictionary::NGroups), groupZ(Dictionary::NGroups), classX(NSpanClasses), classZ(NSpanClasses); + for (uint32_t pattern = 0; pattern < nPatterns; ++pattern) { + if (isCommon[pattern]) { + continue; + } + const auto [rowSpan, colSpan] = spansOf(pattern); + const int group = Dictionary::groupOf(rowSpan, colSpan); + groupX[group].insert(groupX[group].end(), residualsX[pattern].begin(), residualsX[pattern].end()); + groupZ[group].insert(groupZ[group].end(), residualsZ[pattern].begin(), residualsZ[pattern].end()); + auto& cx = classX[Dictionary::spanClass(rowSpan)]; + auto& cz = classZ[Dictionary::spanClass(colSpan)]; + cx.insert(cx.end(), residualsX[pattern].begin(), residualsX[pattern].end()); + cz.insert(cz.end(), residualsZ[pattern].begin(), residualsZ[pattern].end()); + } + const auto wholeX = medianAndWidth(layerX), wholeZ = medianAndWidth(layerZ); + const auto pooled = [this](std::vector& group, std::vector& spanClass, const std::pair& whole) { + if (group.size() >= mSettings.minCount) { + return medianAndWidth(group); + } + if (spanClass.size() >= mSettings.minCount) { + return medianAndWidth(spanClass); + } + return whole; + }; + for (int group = 0; group < Dictionary::NGroups; ++group) { + DictionaryEntry entry; + entry.isGroup = true; + std::tie(entry.xOffset, entry.xWidth) = pooled(groupX[group], classX[group / NSpanClasses], wholeX); + std::tie(entry.zOffset, entry.zWidth) = pooled(groupZ[group], classZ[group % NSpanClasses], wholeZ); + data.entries.push_back(std::move(entry)); + } + + for (std::size_t id = 0; id < data.entries.size(); ++id) { + const auto& entry = data.entries[id]; + if (!(entry.xWidth > 0.f && entry.zWidth > 0.f && std::isfinite(entry.xWidth) && std::isfinite(entry.zWidth))) { + throw std::runtime_error(fmt::format("layer {}: entry {} has widths {} and {} cm", layer, id, entry.xWidth, entry.zWidth)); + } + } +} + +template +void BuildTopologyDictionary::build() +{ + for (int layer = 0; layer < static_cast(mInput.size()); ++layer) { + auto& input = mInput[layer]; + if (input.crossings.empty()) { + throw std::runtime_error(fmt::format("layer {} has no clusters", layer)); + } + auto& data = mDictionary.mLayers[layer]; + data = LayerDictionary{}; + data.depth = measureDepth(layer, input); + makeEntries(layer, input, data); + input = LayerInput{}; + LOGP(info, "{} dictionary layer {}: depth {:.2f} um, {} common patterns", o2::detectors::DetID(DetID).getName(), layer, 1e4 * data.depth, + data.entries.size() - Dictionary::NGroups); + } +} + +template +void BuildTopologyDictionary::saveDictionaryRoot(const std::string& fileName) const +{ + TFile output(fileName.c_str(), "recreate"); + output.WriteObjectAny(&mDictionary, Dictionary::Class(), "ccdb_object"); + output.Close(); +} + +template class BuildTopologyDictionary; +template class BuildTopologyDictionary; + +} // namespace o2::trkft3 diff --git a/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/src/TRKReconstructionLinkDef.h b/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/src/TRKReconstructionLinkDef.h new file mode 100644 index 0000000000000..45b77dfee7454 --- /dev/null +++ b/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/src/TRKReconstructionLinkDef.h @@ -0,0 +1,27 @@ +// Copyright 2019-2026 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. + +#ifdef __CLING__ + +#pragma link off all globals; +#pragma link off all classes; +#pragma link off all functions; + +#pragma link C++ class o2::trkft3::DictionaryEntry + ; +#pragma link C++ class std::vector < o2::trkft3::DictionaryEntry> + ; +#pragma link C++ class o2::trkft3::LayerDictionary + ; +#pragma link C++ class std::vector < o2::trkft3::LayerDictionary> + ; +#pragma link C++ class o2::trkft3::TopologyDictionary < o2::detectors::DetID::TRK> + ; +#pragma link C++ class o2::trkft3::TopologyDictionary < o2::detectors::DetID::FT3> + ; +#pragma link C++ class o2::trkft3::BuildTopologyDictionary < o2::detectors::DetID::TRK> + ; +#pragma link C++ class o2::trkft3::BuildTopologyDictionary < o2::detectors::DetID::FT3> + ; + +#endif diff --git a/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/src/TopologyDictionary.cxx b/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/src/TopologyDictionary.cxx new file mode 100644 index 0000000000000..d035604178bd4 --- /dev/null +++ b/Detectors/Upgrades/ALICE3/TRKFT3/common/reconstruction/src/TopologyDictionary.cxx @@ -0,0 +1,152 @@ +// Copyright 2019-2026 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. + +/// \file TopologyDictionary.cxx + +#include "TRKReconstruction/TopologyDictionary.h" +#include "DataFormatsITSMFT/ClusterPattern.h" +#include "DataFormatsITSMFT/ClusterTopology.h" +#include "TRKBase/SegmentationChip.h" +#include +#include +#include +#include +#include +#include + +namespace o2::trkft3 +{ + +using o2::trk::SegmentationChip; + +template +int TopologyDictionary::layerOf(const ClusterType& cluster) noexcept +{ + if (cluster.layer < 0) { + return -1; + } + if constexpr (DetID == o2::detectors::DetID::TRK) { + if (cluster.subDetID == 0) { + return cluster.layer; + } + if (cluster.subDetID == 1) { + return o2::trk::constants::VD::petal::nLayers + cluster.layer; + } + return -1; + } else { + return cluster.layer; + } +} + +template +math_utils::Point3D TopologyDictionary::localToFlat(const ClusterType& cluster, const math_utils::Point3D& local) noexcept +{ + if (cluster.subDetID != 0) { + return local; + } + const auto flat = SegmentationChip::curvedToFlat(cluster.layer, local.X(), local.Y()); + return {flat.X(), flat.Y(), local.Z()}; +} + +template +math_utils::Point3D TopologyDictionary::centreOfGravity(const ClusterType& cluster, const unsigned char* pattern) noexcept +{ + float cogRow = 0.f, cogCol = 0.f; + itsmft::ClusterPattern::getCOG(pattern[0], pattern[1], pattern + 2, cogRow, cogCol); + float x = 0.f, z = 0.f; + SegmentationChip::detectorToLocalUnchecked(cluster.row, cluster.col, x, z, cluster.subDetID, cluster.layer, cluster.layer); + const bool vd = cluster.subDetID == 0; + x -= cogRow * (vd ? SegmentationChip::PitchRowVD : SegmentationChip::PitchRowMLOT); + z += cogCol * (vd ? SegmentationChip::PitchColVD : SegmentationChip::PitchColMLOT); + return {x, 0.f, z}; +} + +template +int TopologyDictionary::findID(int layer, const unsigned char* pattern) const +{ + const auto& data = mLayers.at(layer); + const int rowSpan = pattern[0], colSpan = pattern[1]; + // The hash leaves out a large pattern's last bytes: compare them all. + const auto nBytes = static_cast(2 + ((rowSpan * colSpan + 7) / 8)); + const auto [first, last] = data.hashMap.equal_range(itsmft::ClusterTopology::getCompleteHash(rowSpan, colSpan, pattern + 2)); + for (auto it = first; it != last; ++it) { + const auto& entryPattern = data.entries[it->second].pattern; + if (entryPattern.size() == nBytes && std::memcmp(entryPattern.data(), pattern, nBytes) == 0) { + return it->second; + } + } + return static_cast(data.entries.size()) - NGroups + groupOf(rowSpan, colSpan); +} + +template +math_utils::Point3D TopologyDictionary::getFlatCoordinates(const ClusterType& cluster, const unsigned char* pattern, int id) const +{ + const int layer = layerOf(cluster); + const auto& data = mLayers.at(layer); + const auto& entry = data.entries.at(id); + math_utils::Point3D flat; + if (entry.isGroup) { + flat = centreOfGravity(cluster, pattern); + } else { + float x = 0.f, z = 0.f; + SegmentationChip::detectorToLocalUnchecked(cluster.row, cluster.col, x, z, cluster.subDetID, cluster.layer, cluster.layer); + const bool vd = cluster.subDetID == 0; + flat.SetX(x + entry.xCOG * (vd ? SegmentationChip::PitchRowVD : SegmentationChip::PitchRowMLOT)); + flat.SetZ(z + entry.zCOG * (vd ? SegmentationChip::PitchColVD : SegmentationChip::PitchColMLOT)); + } + if (mApplyOffsets) { + flat.SetX(flat.X() + entry.xOffset); + flat.SetZ(flat.Z() + entry.zOffset); + } + flat.SetY(data.depth); + return flat; +} + +template +math_utils::Point3D TopologyDictionary::getClusterCoordinates(const ClusterType& cluster, const unsigned char* pattern, int id) const +{ + auto local = getFlatCoordinates(cluster, pattern, id); + if (cluster.subDetID == 0) { + const auto curved = SegmentationChip::flatToCurved(cluster.layer, local.X(), local.Y()); + local.SetX(curved.X()); + local.SetY(curved.Y()); + } + return local; +} + +template +void TopologyDictionary::print() const +{ + LOGP(info, "{} TopologyDictionary: {} layers, offsets {}", o2::detectors::DetID(DetID).getName(), mLayers.size(), mApplyOffsets ? "applied" : "not applied"); + for (std::size_t layer = 0; layer < mLayers.size(); ++layer) { + LOGP(info, " layer {:2}: depth {:6.2f} um, {} common patterns", layer, 1e4 * mLayers[layer].depth, mLayers[layer].entries.size() - NGroups); + } +} + +template +TopologyDictionary* TopologyDictionary::loadFrom(const std::string& fileName, const std::string& objName) +{ + LOGP(info, "Loading {} TopologyDictionary from {} with name {}", o2::detectors::DetID(DetID).getName(), fileName, objName); + std::unique_ptr file{TFile::Open(fileName.c_str(), "READ")}; + if (!file || file->IsZombie()) { + throw std::runtime_error(fmt::format("Failed to open {}", fileName)); + } + auto* dictionary = file->Get>(objName.c_str()); + if (dictionary == nullptr) { + throw std::runtime_error(fmt::format("Failed to load {} from {}", objName, fileName)); + } + return dictionary; +} + +template class TopologyDictionary; +template class TopologyDictionary; + +} // namespace o2::trkft3 From 7e01900e21fdf08a93d6c2b9fa4e448de9207e0a Mon Sep 17 00:00:00 2001 From: Felix Schlepper Date: Mon, 2 Nov 2026 17:50:44 +0100 Subject: [PATCH 2/2] Add ENABLE_UPGRADES macro to CheckTopologyDictionary.C --- .../ALICE3/TRKFT3/common/macros/CheckTopologyDictionary.C | 2 ++ 1 file changed, 2 insertions(+) diff --git a/Detectors/Upgrades/ALICE3/TRKFT3/common/macros/CheckTopologyDictionary.C b/Detectors/Upgrades/ALICE3/TRKFT3/common/macros/CheckTopologyDictionary.C index 6d7a3fac9147b..f97a586373308 100644 --- a/Detectors/Upgrades/ALICE3/TRKFT3/common/macros/CheckTopologyDictionary.C +++ b/Detectors/Upgrades/ALICE3/TRKFT3/common/macros/CheckTopologyDictionary.C @@ -27,6 +27,8 @@ /// /// root -b -q 'CheckTopologyDictionary.C+("TRK", "TRKTopologyDictionary.root", {"tf7", "tf8"}, {"tf2", "tf3"})' +#define ENABLE_UPGRADES // DetID::TRK and FT3; the macro is not compiled with O2's definitions + #if !defined(__CLING__) || defined(__ROOTCLING__) #include "DataFormatsTRKFT3/Cluster.h" #include "DataFormatsTRKFT3/Hit.h"