Traccc core (rebased)

Core/include/Acts/Geometry/GeometryHierarchyMap.hpp

@@ -93,6 +93,16 @@ class GeometryHierarchyMap {
   /// Return the number of stored elements.
   Size size() const { return m_values.size(); }
 
+  // Return the geometry identifier - value pairs (i.e., it reconstructs the
+  // input).
+  const std::vector<InputElement> getElements() const {
+    std::vector<InputElement> res;
+    for (std::size_t i = 0; i < size(); i++) {
+      res.push_back({m_ids[i], m_values[i]});
+    }
+    return res;
+  }
+
   /// Access the geometry identifier for the i-th element with bounds check.
   ///
   /// @throws std::out_of_range for invalid indices

Examples/Algorithms/CMakeLists.txt

@@ -18,3 +18,4 @@ add_subdirectory(Vertexing)
 add_subdirectory_if(Alignment ACTS_BUILD_ALIGNMENT)
 add_subdirectory(Utilities)
 add_subdirectory(AmbiguityResolution)
+add_subdirectory_if(Traccc ACTS_BUILD_PLUGIN_TRACCC)

Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationAlgorithm.hpp

@@ -88,21 +88,25 @@ class DigitizationAlgorithm final : public IAlgorithm {
   using Digitizer = std::variant<CombinedDigitizer<0>, CombinedDigitizer<1>,
                                  CombinedDigitizer<2>, CombinedDigitizer<3>,
                                  CombinedDigitizer<4>>;
-
+                                 
   /// Configuration of the Algorithm
   DigitizationConfig m_cfg;
   /// Digitizers within geometry hierarchy
   Acts::GeometryHierarchyMap<Digitizer> m_digitizers;
   /// Geometric digtizer
   ActsFatras::Channelizer m_channelizer;
 
+  using CellsMap =
+      std::map<Acts::GeometryIdentifier, std::vector<Cluster::Cell>>;
+
   ReadDataHandle<SimHitContainer> m_simContainerReadHandle{this,
                                                            "SimHitContainer"};
 
   WriteDataHandle<IndexSourceLinkContainer> m_sourceLinkWriteHandle{
       this, "SourceLinks"};
   WriteDataHandle<MeasurementContainer> m_measurementWriteHandle{
       this, "Measurements"};
+  WriteDataHandle<CellsMap> m_cellsWriteHandle{this, "Cells"};
   WriteDataHandle<ClusterContainer> m_clusterWriteHandle{this, "Clusters"};
   WriteDataHandle<IndexMultimap<ActsFatras::Barcode>>
       m_measurementParticlesMapWriteHandle{this, "MeasurementParticlesMap"};

Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationConfig.hpp

@@ -139,6 +139,8 @@ class DigitizationConfig {
   std::string outputSourceLinks = "sourcelinks";
   /// Output measurements collection.
   std::string outputMeasurements = "measurements";
+  /// Output cells map (geoID -> collection of cells).
+  std::string outputCells = "cells";
   /// Output cluster collection.
   std::string outputClusters = "clusters";
   /// Output collection to map measured hits to contributing particles.
@@ -155,7 +157,7 @@ class DigitizationConfig {
   /// How close do parameters have to be to consider merged
   const double mergeNsigma;
   /// Consider clusters that share a corner as merged (8-cell connectivity)
-  const bool mergeCommonCorner;
+  bool mergeCommonCorner;
   /// Energy deposit threshold for accepting a hit
   /// For a generic readout frontend we assume 1000 e/h pairs, in Si each
   /// e/h-pair requiers on average an energy of 3.65 eV (PDG  review 2023,

Examples/Algorithms/Digitization/src/DigitizationAlgorithm.cpp

@@ -72,6 +72,7 @@ ActsExamples::DigitizationAlgorithm::DigitizationAlgorithm(
   m_simContainerReadHandle.initialize(m_cfg.inputSimHits);
   m_sourceLinkWriteHandle.initialize(m_cfg.outputSourceLinks);
   m_measurementWriteHandle.initialize(m_cfg.outputMeasurements);
+  m_cellsWriteHandle.initialize(m_cfg.outputCells);
   m_clusterWriteHandle.initialize(m_cfg.outputClusters);
   m_measurementParticlesMapWriteHandle.initialize(
       m_cfg.outputMeasurementParticlesMap);
@@ -154,6 +155,10 @@ ActsExamples::ProcessCode ActsExamples::DigitizationAlgorithm::execute(
   // Some statistics
   std::size_t skippedHits = 0;
 
+  // Some algorithms do the clusterization themselves such as the traccc chain.
+  // Thus we need to store the cell data from the simulation.
+  CellsMap cellsMap;
+
   ACTS_DEBUG("Starting loop over modules ...");
   for (const auto& simHitsGroup : groupByModule(simHits)) {
     // Manual pair unpacking instead of using
@@ -255,8 +260,18 @@ ActsExamples::ProcessCode ActsExamples::DigitizationAlgorithm::execute(
             moduleClusters.add(std::move(dParameters), simHitIdx);
           }
 
-          for (auto& [dParameters, simhits] :
-               moduleClusters.digitizedParameters()) {
+          auto digitizeParametersResult = moduleClusters.digitizedParameters();
+
+          // Store the data of the cells from the simulation.
+          std::vector<Cluster::Cell> cells;
+          for (auto& [dParameters, simhits] : digitizeParametersResult) {
+            for (auto cell : dParameters.cluster.channels) {
+              cells.push_back(std::move(cell));
+            }
+          }
+          cellsMap.insert({moduleGeoId, cells});
+
+          for (auto& [dParameters, simhits] : digitizeParametersResult) {
             // The measurement container is unordered and the index under which
             // the measurement will be stored is known before adding it.
             Index measurementIdx = measurements.size();
@@ -293,6 +308,7 @@ ActsExamples::ProcessCode ActsExamples::DigitizationAlgorithm::execute(
 
   m_sourceLinkWriteHandle(ctx, std::move(sourceLinks));
   m_measurementWriteHandle(ctx, std::move(measurements));
+  m_cellsWriteHandle(ctx, std::move(cellsMap));
   m_clusterWriteHandle(ctx, std::move(clusters));
   m_measurementParticlesMapWriteHandle(ctx, std::move(measurementParticlesMap));
   m_measurementSimHitsMapWriteHandle(ctx, std::move(measurementSimHitsMap));

Examples/Algorithms/Traccc/CMakeLists.txt

@@ -0,0 +1,2 @@
+add_subdirectory(Common)
+add_subdirectory(Host)

Examples/Algorithms/Traccc/Common/CMakeLists.txt

@@ -0,0 +1,26 @@
+add_library(
+  ActsExamplesTracccCommon SHARED
+  src/TracccChainAlgorithmBase.cpp
+  src/Conversion/CellMapConversion.cpp
+  src/Conversion/DigitizationConversion.cpp
+  src/Conversion/MeasurementConversion.cpp
+  src/Debug/Debug.cpp
+)
+
+target_include_directories(
+  ActsExamplesTracccCommon
+  PUBLIC $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/include>
+)
+
+target_link_libraries(
+  ActsExamplesTracccCommon
+  PUBLIC
+    ActsPluginTraccc
+    ActsExamplesFramework
+    ActsExamplesDigitization
+)
+
+install(
+  TARGETS ActsExamplesTracccCommon
+  LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR}
+)

Examples/Algorithms/Traccc/Common/include/ActsExamples/Traccc/Common/Conversion/CellMapConversion.hpp

@@ -0,0 +1,34 @@
+// This file is part of the Acts project.
+//
+// Copyright (C) 2023 CERN for the benefit of the Acts project
+//
+// This Source Code Form is subject to the terms of the Mozilla Public
+// License, v. 2.0. If a copy of the MPL was not distributed with this
+// file, You can obtain one at http://mozilla.org/MPL/2.0/.
+
+#pragma once
+
+// Acts include(s)
+#include "Acts/Geometry/GeometryIdentifier.hpp"
+
+// Acts Examples include(s)
+#include "ActsExamples/EventData/Cluster.hpp"
+
+// Traccc include(s)
+#include "traccc/edm/cell.hpp"
+
+// System include(s).
+#include <cstdint>
+#include <cstdlib>
+#include <map>
+#include <vector>
+
+namespace ActsExamples::Traccc::Common::Conversion {
+
+/// @brief Converts a "geometry ID -> generic cell collection type" map to a "geometry ID -> traccc cell collection" map.
+/// @note The function sets the module link of the cells in the output to 0.
+/// @return Map from geometry ID to its cell data (as a vector of traccc cell data)
+std::map<std::uint64_t, std::vector<traccc::cell>> tracccCellsMap(
+    const std::map<Acts::GeometryIdentifier, std::vector<Cluster::Cell>>& map);
+
+}  // namespace ActsExamples::Traccc::Common::Conversion

Examples/Algorithms/Traccc/Common/include/ActsExamples/Traccc/Common/Conversion/DigitizationConversion.hpp

@@ -0,0 +1,29 @@
+// This file is part of the Acts project.
+//
+// Copyright (C) 2023 CERN for the benefit of the Acts project
+//
+// This Source Code Form is subject to the terms of the Mozilla Public
+// License, v. 2.0. If a copy of the MPL was not distributed with this
+// file, You can obtain one at http://mozilla.org/MPL/2.0/.
+
+#pragma once
+
+// Acts include(s)
+#include "Acts/Geometry/GeometryHierarchyMap.hpp"
+
+// Acts Examples include(s)
+#include "ActsExamples/Digitization/DigitizationConfig.hpp"
+
+// Traccc include(s)
+#include "traccc/io/digitization_config.hpp"
+
+namespace ActsExamples::Traccc::Common::Conversion {
+
+/// @brief Creates a traccc digitalization config from an Acts geometry hierarchy map
+/// that contains the digitization configuration.
+/// @param config the Acts geometry hierarchy map that contains the digitization configuration.
+/// @return a traccc digitization config.
+traccc::digitization_config tracccConfig(
+    const Acts::GeometryHierarchyMap<DigiComponentsConfig>& config);
+
+}  // namespace ActsExamples::Traccc::Common::Conversion

Examples/Algorithms/Traccc/Common/include/ActsExamples/Traccc/Common/Conversion/MeasurementConversion.hpp

@@ -0,0 +1,175 @@
+// This file is part of the Acts project.
+//
+// Copyright (C) 2023 CERN for the benefit of the Acts project
+//
+// This Source Code Form is subject to the terms of the Mozilla Public
+// License, v. 2.0. If a copy of the MPL was not distributed with this
+// file, You can obtain one at http://mozilla.org/MPL/2.0/.
+
+#pragma once
+
+// Plugin include(s)
+#include "Acts/Plugins/Traccc/Detail/AlgebraConversion.hpp"
+
+// Acts include(s)
+#include "Acts/Definitions/Algebra.hpp"
+#include "Acts/Definitions/TrackParametrization.hpp"
+#include "Acts/EventData/Measurement.hpp"
+#include "Acts/EventData/SourceLink.hpp"
+#include "Acts/Geometry/GeometryIdentifier.hpp"
+
+// Acts Examples include(s)
+#include "ActsExamples/EventData/IndexSourceLink.hpp"
+
+// Detray include(s)
+#include "detray/core/detector.hpp"
+#include "detray/tracks/bound_track_parameters.hpp"
+
+// Traccc include(s)
+#include "traccc/definitions/qualifiers.hpp"
+#include "traccc/definitions/track_parametrization.hpp"
+#include "traccc/edm/measurement.hpp"
+#include "traccc/edm/track_state.hpp"
+
+// System include(s)
+#include <cstdint>
+#include <cstdlib>
+#include <memory>
+#include <variant>
+#include <vector>
+
+namespace ActsExamples::Traccc::Common::Conversion {
+
+/// @brief Converts a traccc bound index to an Acts bound index.
+/// @param tracccBoundIndex the traccc bound index.
+/// @returns an Acts bound index.
+Acts::BoundIndices boundIndex(const traccc::bound_indices tracccBoundIndex);
+
+/// @brief Creates an Acts measurement from a traccc measurement.
+/// @tparam the dimension of the Acts measurement (subspace size).
+/// @param m the traccc measurement.
+/// @param sl the Acts source link to use for the Acts measurement.
+/// @returns an Acts measurement with data copied from the traccc measurement
+/// and with its source link set to the one provided to the function.
+template <std::size_t dim>
+inline Acts::Measurement<Acts::BoundIndices, dim> measurement(
+    const traccc::measurement& m, const Acts::SourceLink sl) {
+  auto params = Acts::TracccPlugin::detail::toActsVector<dim>(m.local);
+  std::array<Acts::BoundIndices, dim> indices;
+  for (unsigned int i = 0; i < dim; i++) {
+    indices[i] = boundIndex(traccc::bound_indices(m.subs.get_indices()[i]));
+  }
+  auto cov = Eigen::DiagonalMatrix<Acts::ActsScalar, static_cast<int>(dim)>(
+                 Acts::TracccPlugin::detail::toActsVector<dim>(m.variance))
+                 .toDenseMatrix();
+  return Acts::Measurement<Acts::BoundIndices, dim>(std::move(sl), indices,
+                                                    params, cov);
+}
+
+/// @brief Creates an Acts bound variant measurement from a traccc measurement.
+/// Using recursion, the functions determines the dimension of the traccc
+/// measurement which is used for the Acts measurement that the bound variant
+/// measurement holds. The dimension must lie between [0; max_dim].
+/// @tparam max_dim the largest possible dimension of any measurement type in the variant (default = 4).
+/// @param m the traccc measurement.
+/// @param sl the Acts source link to use for the Acts measurement.
+/// @returns an Acts bound variant measurement with data copied from the traccc measurement
+/// and with its source link set to the one provided to the function.
+template <std::size_t max_dim = 4UL>
+inline Acts::BoundVariantMeasurement boundVariantMeasurement(
+    const traccc::measurement& m, const Acts::SourceLink sl) {
+  if constexpr (max_dim == 0UL) {
+    std::string errorMsg = "Invalid/mismatching measurement dimension: " +
+                           std::to_string(m.meas_dim);
+    throw std::runtime_error(errorMsg.c_str());
+  } else {
+    if (m.meas_dim == max_dim) {
+      return measurement<max_dim>(m, sl);
+    }
+    return boundVariantMeasurement<max_dim - 1>(m, sl);
+  }
+}
+
+/// @brief Gets the the local position of the measurement.
+/// @param measurement the Acts measurement.
+/// @returns A two-dimensional vector containing the local position.
+/// The first item in the vector is local position on axis 0 and
+/// I.e., [local position (axis 0), local position (axis 1)].
+/// @note if the dimension is less than 2 then the remaining values are set to 0.
+template <std::size_t dim>
+inline Acts::ActsVector<2> getLocal(
+    const Acts::Measurement<Acts::BoundIndices, dim>& measurement) {
+  traccc::scalar loc0 = 0;
+  traccc::scalar loc1 = 0;
+  if constexpr (dim > Acts::BoundIndices::eBoundLoc0) {
+    loc0 = measurement.parameters()(Acts::BoundIndices::eBoundLoc0);
+  }
+  if constexpr (dim > Acts::BoundIndices::eBoundLoc1) {
+    loc1 = measurement.parameters()(Acts::BoundIndices::eBoundLoc1);
+  }
+  return Acts::ActsVector<2>(loc0, loc1);
+}
+
+/// @brief Get the the local position of the measurement.
+/// @param measurement the Acts bound variant measurement.
+/// @return A two-dimensional vector containing the local position.
+/// I.e., [local position (axis 0), local position (axis 1)].
+/// @note if the dimension is less than 2 then the remaining values are set to 0.
+inline Acts::ActsVector<2> getLocal(
+    const Acts::BoundVariantMeasurement& measurement) {
+  return std::visit([](auto& m) { return getLocal(m); }, measurement);
+}
+
+/// @brief Get the the variance of the measurement.
+/// @param measurement the Acts measurement.
+/// @return A two-dimensional vector containing the variance.
+/// I.e., [variance (axis 0), variance (axis 1)].
+/// @note if the dimension is less than 2 then the remaining values are set to 0.
+template <std::size_t dim>
+inline Acts::ActsVector<2> getVariance(
+    const Acts::Measurement<Acts::BoundIndices, dim>& measurement) {
+  traccc::scalar var0 = 0;
+  traccc::scalar var1 = 0;
+  if constexpr (dim >= Acts::BoundIndices::eBoundLoc0) {
+    var0 = measurement.covariance()(Acts::BoundIndices::eBoundLoc0,
+                                    Acts::BoundIndices::eBoundLoc0);
+  }
+  if constexpr (dim > Acts::BoundIndices::eBoundLoc1) {
+    var1 = measurement.covariance()(Acts::BoundIndices::eBoundLoc1,
+                                    Acts::BoundIndices::eBoundLoc1);
+  }
+  return Acts::ActsVector<2>(var0, var1);
+}
+
+/// @brief Get the the variance of the measurement.
+/// @param measurement the Acts bound variant measurement.
+/// @return A two-dimensional vector containing the variance.
+/// I.e., [variance (axis 0), variance (axis 1)].
+/// @note if the dimension is less than 2 then the remaining values are set to 0.
+inline Acts::ActsVector<2> getVariance(
+    const Acts::BoundVariantMeasurement& measurement) {
+  return std::visit([](auto& m) { return getVariance(m); }, measurement);
+}
+
+/// @brief Converts traccc measurements to acts measurements.
+/// @param detector The detray detector,
+/// @param measurements The traccc measurements,
+/// @return A vector of Acts bound variant measurements.
+/// @note The type IndexSourceLink is used for the measurements' source links.
+template <typename detector_t, typename allocator_t>
+inline auto createActsMeasurements(
+    const detector_t& detector,
+    const std::vector<traccc::measurement, allocator_t>& measurements) {
+  std::vector<Acts::BoundVariantMeasurement> measurementContainer;
+  for (const traccc::measurement& m : measurements) {
+    Acts::GeometryIdentifier moduleGeoId(
+        detector.surface(m.surface_link).source);
+    Index measurementIdx = measurementContainer.size();
+    IndexSourceLink idxSourceLink{moduleGeoId, measurementIdx};
+    measurementContainer.push_back(
+        boundVariantMeasurement(m, Acts::SourceLink{idxSourceLink}));
+  }
+  return measurementContainer;
+}
+
+}  // namespace ActsExamples::Traccc::Common::Conversion