Cell Sorting, main branch (2026.02.20.)

core/CMakeLists.txt

@@ -1,6 +1,6 @@
 # TRACCC library, part of the ACTS project (R&D line)
 #
-# (c) 2021-2025 CERN for the benefit of the ACTS project
+# (c) 2021-2026 CERN for the benefit of the ACTS project
 #
 # Mozilla Public License Version 2.0
 
@@ -78,6 +78,8 @@ traccc_add_library( traccc_core core TYPE SHARED
   "src/clusterization/measurement_sorting_algorithm.cpp"
   "include/traccc/clusterization/clusterization_algorithm.hpp"
   "src/clusterization/clusterization_algorithm.cpp"
+  "include/traccc/clusterization/silicon_cell_sorting_algorithm.hpp"
+  "src/clusterization/silicon_cell_sorting_algorithm.cpp"
   # Finding algorithmic code
   "include/traccc/finding/candidate_link.hpp"
   "include/traccc/finding/finding_config.hpp"

core/include/traccc/clusterization/silicon_cell_sorting_algorithm.hpp

@@ -0,0 +1,59 @@
+/** TRACCC library, part of the ACTS project (R&D line)
+ *
+ * (c) 2026 CERN for the benefit of the ACTS project
+ *
+ * Mozilla Public License Version 2.0
+ */
+
+#pragma once
+
+// Library include(s).
+#include "traccc/edm/silicon_cell_collection.hpp"
+#include "traccc/utils/algorithm.hpp"
+#include "traccc/utils/messaging.hpp"
+
+// VecMem include(s).
+#include <vecmem/memory/memory_resource.hpp>
+
+// System include(s).
+#include <functional>
+#include <memory>
+
+namespace traccc::host {
+
+/// Algorithm sorting the detector cell information
+///
+/// Clusterization algorithms may (and do) rely on cells being "strictly
+/// sorted". This algorithm can be used to make sure that this would be the
+/// case.
+///
+class silicon_cell_sorting_algorithm
+    : public algorithm<edm::silicon_cell_collection::host(
+          const edm::silicon_cell_collection::const_view&)>,
+      public messaging {
+
+    public:
+    /// Constructor
+    ///
+    /// @param mr The memory resource to use for the algorithm
+    /// @param logger The logger to use for the algorithm
+    ///
+    silicon_cell_sorting_algorithm(
+        vecmem::memory_resource& mr,
+        std::unique_ptr<const Logger> logger = getDummyLogger().clone());
+
+    /// Callable operator performing the sorting on a container
+    ///
+    /// @param cells The cells to sort
+    /// @return The sorted cells
+    ///
+    output_type operator()(
+        const edm::silicon_cell_collection::const_view& cells) const override;
+
+    private:
+    /// The memory resource to use
+    std::reference_wrapper<vecmem::memory_resource> m_mr;
+
+};  // class silicon_cell_sorting_algorithm
+
+}  // namespace traccc::host

core/src/clusterization/silicon_cell_sorting_algorithm.cpp

@@ -0,0 +1,48 @@
+/** TRACCC library, part of the ACTS project (R&D line)
+ *
+ * (c) 2026 CERN for the benefit of the ACTS project
+ *
+ * Mozilla Public License Version 2.0
+ */
+
+// Library include(s).
+#include "traccc/clusterization/silicon_cell_sorting_algorithm.hpp"
+
+// System include(s).
+#include <algorithm>
+#include <numeric>
+
+namespace traccc::host {
+
+silicon_cell_sorting_algorithm::silicon_cell_sorting_algorithm(
+    vecmem::memory_resource& mr, std::unique_ptr<const Logger> logger)
+    : messaging(std::move(logger)), m_mr(mr) {}
+
+silicon_cell_sorting_algorithm::output_type
+silicon_cell_sorting_algorithm::operator()(
+    const edm::silicon_cell_collection::const_view& cells_view) const {
+
+    // Create a device container on top of the view.
+    const edm::silicon_cell_collection::const_device cells{cells_view};
+
+    // Create a vector of cell indices, which would be sorted.
+    vecmem::vector<unsigned int> indices(cells.size(), &(m_mr.get()));
+    std::iota(indices.begin(), indices.end(), 0u);
+
+    // Sort the indices according to the cells.
+    std::sort(indices.begin(), indices.end(),
+              [&](unsigned int lhs, unsigned int rhs) {
+                  return cells.at(lhs) < cells.at(rhs);
+              });
+
+    // Fill an output container with the sorted cells.
+    edm::silicon_cell_collection::host result{m_mr.get()};
+    for (unsigned int i : indices) {
+        result.push_back(cells.at(i));
+    }
+
+    // Return the sorted cells.
+    return result;
+}
+
+}  // namespace traccc::host

device/alpaka/CMakeLists.txt

@@ -44,6 +44,8 @@ traccc_add_alpaka_library( traccc_alpaka alpaka TYPE SHARED
   "src/clusterization/clusterization_algorithm.cpp"
   "include/traccc/alpaka/clusterization/measurement_sorting_algorithm.hpp"
   "src/clusterization/measurement_sorting_algorithm.cpp"
+  "include/traccc/alpaka/clusterization/silicon_cell_sorting_algorithm.hpp"
+  "src/clusterization/silicon_cell_sorting_algorithm.cpp"
   # Seeding code
   "include/traccc/alpaka/seeding/silicon_pixel_spacepoint_formation_algorithm.hpp"
   "src/seeding/silicon_pixel_spacepoint_formation_algorithm.cpp"

device/alpaka/include/traccc/alpaka/clusterization/silicon_cell_sorting_algorithm.hpp

@@ -0,0 +1,58 @@
+/** TRACCC library, part of the ACTS project (R&D line)
+ *
+ * (c) 2026 CERN for the benefit of the ACTS project
+ *
+ * Mozilla Public License Version 2.0
+ */
+
+#pragma once
+
+// Local include(s).
+#include "traccc/alpaka/utils/algorithm_base.hpp"
+
+// Project include(s).
+#include "traccc/device/algorithm_base.hpp"
+#include "traccc/edm/silicon_cell_collection.hpp"
+#include "traccc/utils/algorithm.hpp"
+#include "traccc/utils/messaging.hpp"
+
+namespace traccc::alpaka {
+
+/// Algorithm sorting the detector cell information
+///
+/// Clusterization algorithms may (and do) rely on cells being "strictly
+/// sorted". This algorithm can be used to make sure that this would be the
+/// case.
+///
+class silicon_cell_sorting_algorithm
+    : public algorithm<edm::silicon_cell_collection::buffer(
+          const edm::silicon_cell_collection::const_view&)>,
+      device::algorithm_base,
+      alpaka::algorithm_base,
+      public messaging {
+
+    public:
+    /// Constructor for the algorithm
+    ///
+    /// @param mr The memory resource(s) to use in the algorithm
+    /// @param copy The copy object to use for copying data between device
+    ///             and host memory blocks
+    /// @param q The Alpaka queue to perform the operations in
+    /// @param config The clustering configuration
+    ///
+    silicon_cell_sorting_algorithm(
+        const traccc::memory_resource& mr, ::vecmem::copy& copy,
+        alpaka::queue& q,
+        std::unique_ptr<const Logger> logger = getDummyLogger().clone());
+
+    /// Callable operator performing the sorting on a container
+    ///
+    /// @param cells The cells to sort
+    /// @return The sorted cells
+    ///
+    output_type operator()(
+        const edm::silicon_cell_collection::const_view& cells) const override;
+
+};  // class silicon_cell_sorting_algorithm
+
+}  // namespace traccc::alpaka

device/alpaka/src/clusterization/silicon_cell_sorting_algorithm.cpp

@@ -0,0 +1,109 @@
+/** TRACCC library, part of the ACTS project (R&D line)
+ *
+ * (c) 2026 CERN for the benefit of the ACTS project
+ *
+ * Mozilla Public License Version 2.0
+ */
+
+// Local include(s).
+#include "traccc/alpaka/clusterization/silicon_cell_sorting_algorithm.hpp"
+
+#include "../utils/get_queue.hpp"
+#include "../utils/parallel_algorithms.hpp"
+#include "../utils/thread_id.hpp"
+
+// Project include(s).
+#include "traccc/clusterization/device/silicon_cell_sorter.hpp"
+#include "traccc/clusterization/device/sorting_index_filler.hpp"
+
+namespace traccc::alpaka {
+namespace kernels {
+
+/// Kernel filling the output buffer with sorted cells.
+struct fill_sorted_silicon_cells {
+    /// @param[in] acc Alpaka accelerator object
+    /// @param[in] input_view View of the input cells
+    /// @param[out] output_view View of the output cells
+    /// @param[in] sorted_indices_view View of the sorted cell indices
+    ///
+    template <typename TAcc>
+    ALPAKA_FN_ACC void operator()(
+        TAcc const& acc,
+        const edm::silicon_cell_collection::const_view input_view,
+        edm::silicon_cell_collection::view output_view,
+        const vecmem::data::vector_view<const unsigned int> sorted_indices_view)
+        const {
+
+        // Create the device objects.
+        const edm::silicon_cell_collection::const_device input{input_view};
+        edm::silicon_cell_collection::device output{output_view};
+        const vecmem::device_vector<const unsigned int> sorted_indices{
+            sorted_indices_view};
+
+        // Stop early if we can.
+        const unsigned int index = details::thread_id1{acc}.getGlobalThreadId();
+        if (index >= input.size()) {
+            return;
+        }
+
+        // Copy one measurement into the correct position.
+        output.at(index) = input.at(sorted_indices.at(index));
+    }
+};  // struct fill_sorted_silicon_cells
+
+}  // namespace kernels
+
+silicon_cell_sorting_algorithm::silicon_cell_sorting_algorithm(
+    const traccc::memory_resource& mr, vecmem::copy& copy, alpaka::queue& q,
+    std::unique_ptr<const Logger> logger)
+    : device::algorithm_base(mr, copy),
+      alpaka::algorithm_base(q),
+      messaging(std::move(logger)) {}
+
+auto silicon_cell_sorting_algorithm::operator()(
+    const edm::silicon_cell_collection::const_view& cells_view) const
+    -> output_type {
+
+    // Exit early if there are no cells.
+    if (cells_view.capacity() == 0) {
+        return {};
+    }
+
+    // Get a convenience variable for the queue that we'll be using.
+    auto& aqueue = details::get_queue(queue());
+
+    // Create a vector of cell indices, which would be sorted.
+    vecmem::data::vector_buffer<unsigned int> indices(cells_view.capacity(),
+                                                      mr().main);
+    copy().setup(indices)->wait();
+    details::for_each(aqueue, mr(), indices.ptr(),
+                      indices.ptr() + indices.capacity(),
+                      device::sorting_index_filler{indices});
+
+    // Sort the indices according to the (correct) order of the cells.
+    details::sort(aqueue, mr(), indices.ptr(),
+                  indices.ptr() + indices.capacity(),
+                  device::silicon_cell_sorter{cells_view});
+
+    // Create the output buffer.
+    output_type result{cells_view.capacity(), mr().main,
+                       cells_view.size().capacity()
+                           ? vecmem::data::buffer_type::resizable
+                           : vecmem::data::buffer_type::fixed_size};
+    copy().setup(result)->ignore();
+    copy()(cells_view.size(), result.size())->ignore();
+
+    // Fill it with the sorted cells.
+    const unsigned int BLOCK_SIZE = warp_size() * 8;
+    const unsigned int n_blocks =
+        (cells_view.capacity() + BLOCK_SIZE - 1) / BLOCK_SIZE;
+    auto workDiv = makeWorkDiv<Acc>(BLOCK_SIZE, n_blocks);
+    ::alpaka::exec<Acc>(aqueue, workDiv, kernels::fill_sorted_silicon_cells{},
+                        cells_view, vecmem::get_data(result),
+                        vecmem::get_data(indices));
+
+    // Return the sorted buffer.
+    return result;
+}
+
+}  // namespace traccc::alpaka

device/common/CMakeLists.txt

@@ -34,6 +34,7 @@ traccc_add_library( traccc_device_common device_common
    "include/traccc/clusterization/device/barcode_based_sorter.hpp"
    "include/traccc/clusterization/device/ccl_kernel_definitions.hpp"
    "include/traccc/clusterization/device/ccl_kernel.hpp"
+   "include/traccc/clusterization/device/silicon_cell_sorter.hpp"
    "include/traccc/clusterization/device/sorting_index_filler.hpp"
    "include/traccc/clusterization/device/clusterization_algorithm.hpp"
    "src/clusterization/clusterization_algorithm.cpp"

device/common/include/traccc/clusterization/device/silicon_cell_sorter.hpp

@@ -0,0 +1,51 @@
+/** TRACCC library, part of the ACTS project (R&D line)
+ *
+ * (c) 2026 CERN for the benefit of the ACTS project
+ *
+ * Mozilla Public License Version 2.0
+ */
+
+#pragma once
+
+// Project include(s).
+#include "traccc/definitions/qualifiers.hpp"
+#include "traccc/edm/silicon_cell_collection.hpp"
+
+namespace traccc::device {
+
+/// Functor to help with sorting cells on a device
+class silicon_cell_sorter {
+
+    public:
+    /// Constructor, capturing the silicon cells to sort indices by
+    silicon_cell_sorter(const edm::silicon_cell_collection::const_view& cells)
+        : m_cells(cells) {}
+
+    /// Index comparison operator
+    ///
+    /// The logic is a bit convoluted here. :-( The index sequence/vector that
+    /// we sort may be larger than the cell collection. (Whenever we are dealing
+    /// with resizable cell collections. Which may happen sometimes.) In this
+    /// case the behaviour should be that for the non-existent cells the indices
+    /// should be left unchanged. Which the following logic should do...
+    ///
+    TRACCC_HOST_DEVICE bool operator()(unsigned int lhs,
+                                       unsigned int rhs) const {
+
+        const edm::silicon_cell_collection::const_device cells{m_cells};
+        if (lhs >= cells.size()) {
+            return false;
+        }
+        if (rhs >= cells.size()) {
+            return true;
+        }
+        return cells.at(lhs) < cells.at(rhs);
+    }
+
+    private:
+    /// The cells to sort indices by
+    edm::silicon_cell_collection::const_view m_cells;
+
+};  // class silicon_cell_sorter
+
+}  // namespace traccc::device

device/cuda/CMakeLists.txt

@@ -51,6 +51,8 @@ traccc_add_library( traccc_cuda cuda TYPE SHARED
   "src/clusterization/clusterization_algorithm.cu"
   "include/traccc/cuda/clusterization/measurement_sorting_algorithm.hpp"
   "src/clusterization/measurement_sorting_algorithm.cu"
+  "include/traccc/cuda/clusterization/silicon_cell_sorting_algorithm.hpp"
+  "src/clusterization/silicon_cell_sorting_algorithm.cu"
   "src/clusterization/kernels/ccl_kernel.cu"
   "src/clusterization/kernels/ccl_kernel.cuh"
   "src/clusterization/kernels/reify_cluster_data.cu"