From 13abbc541a7502a964805235fd318e05a0d1f750 Mon Sep 17 00:00:00 2001
From: Joana Niermann <joana.niermann@cern.ch>
Date: Tue, 14 Nov 2023 15:14:23 +0100
Subject: [PATCH] Initial draft

---
 core/include/detray/tools/grid_builder.hpp    |  4 +-
 .../detray/tools/material_map_builder.hpp     | 82 +++++++++++++++++++
 core/include/detray/utils/type_list.hpp       |  1 +
 .../io/common/{ => detail}/grid_reader.hpp    | 53 +++++-------
 .../detray/io/common/detector_reader.hpp      |  6 +-
 .../io/common/homogeneous_material_reader.hpp |  5 ++
 .../detray/io/common/material_map_reader.hpp  | 76 +++++++++++++++++
 .../detray/io/common/surface_grid_reader.hpp  | 53 ++++++++++++
 io/include/detray/io/json/json_reader.hpp     | 16 ++--
 .../detectors/create_telescope_detector.hpp   | 18 ++--
 10 files changed, 265 insertions(+), 49 deletions(-)
 create mode 100644 core/include/detray/tools/material_map_builder.hpp
 rename io/include/detray/io/common/{ => detail}/grid_reader.hpp (90%)
 create mode 100644 io/include/detray/io/common/material_map_reader.hpp
 create mode 100644 io/include/detray/io/common/surface_grid_reader.hpp

diff --git a/core/include/detray/tools/grid_builder.hpp b/core/include/detray/tools/grid_builder.hpp
index 0da37dca89..a96dba15d0 100644
--- a/core/include/detray/tools/grid_builder.hpp
+++ b/core/include/detray/tools/grid_builder.hpp
@@ -32,10 +32,12 @@ namespace detray {
 template <typename detector_t, typename grid_t,
           typename bin_filler_t = detail::fill_by_pos,
           typename grid_factory_t = grid_factory_type<grid_t>>
-class grid_builder final : public volume_decorator<detector_t> {
+class grid_builder : public volume_decorator<detector_t> {
 
     public:
     using scalar_type = typename detector_t::scalar_type;
+    using detector_type = detector_t;
+    using value_type = typename detector_type::surface_type;
 
     // TODO: nullptr can lead to exceptions, remove in the full implementation
     DETRAY_HOST
diff --git a/core/include/detray/tools/material_map_builder.hpp b/core/include/detray/tools/material_map_builder.hpp
new file mode 100644
index 0000000000..773ade1386
--- /dev/null
+++ b/core/include/detray/tools/material_map_builder.hpp
@@ -0,0 +1,82 @@
+/** Detray library, part of the ACTS project (R&D line)
+ *
+ * (c) 2022-2023 CERN for the benefit of the ACTS project
+ *
+ * Mozilla Public License Version 2.0
+ */
+
+#pragma once
+
+// Project include(s).
+#include "detray/geometry/detector_volume.hpp"
+#include "detray/geometry/surface.hpp"
+#include "detray/tools/bin_association.hpp"
+#include "detray/tools/bin_fillers.hpp"
+#include "detray/tools/grid_factory.hpp"
+#include "detray/tools/surface_factory_interface.hpp"
+#include "detray/tools/volume_builder.hpp"
+#include "detray/tools/volume_builder_interface.hpp"
+
+// System include(s)
+#include <array>
+#include <cassert>
+#include <memory>
+#include <vector>
+
+namespace detray {
+
+/// @brief Build a grid of a certain shape.
+///
+/// Decorator class to a volume builder that adds a grid as the volumes
+/// geometry accelerator structure.
+template <typename detector_t,
+          typename bin_filler_t = detail::fill_by_pos>
+class material_map_builder final {
+
+    public:
+    using scalar_type = typename detector_t::scalar_type;
+    using detector_type = detector_t;
+    using value_type = typename detector_type::surface_type;
+
+    struct bin_data {
+        n_axis::multi_bin<2u> local_bin;
+        material_slab<scalar_type> mat;
+    };
+
+    // TODO: nullptr can lead to exceptions, remove in the full implementation
+    DETRAY_HOST
+    material_map_builder(
+        std::unique_ptr<volume_builder_interface<detector_t>> vol_builder)
+        : grid_builder_t(std::move(vol_builder)) {}
+
+    /// Add the volume and the grid to the detector @param det
+    DETRAY_HOST
+    auto build(detector_t &det, typename detector_t::geometry_context ctx = {})
+        -> typename detector_t::volume_type * override {
+
+        /// Build the material grid for every surface that has material
+        for (auto& [idx, sf_desc] : detray::views::enumerate(this->surfaces())) {
+            // Only map material to surfaces
+            if (bin_data_per_surface.at(idx).empty()) {
+                continue;
+            }
+        }
+        // Add the grid to the detector and link it to its surface
+        constexpr auto gid{detector_t::materials::template get_id<grid_t>()};
+        det.material_store().template push_back<gid>(m_grid);
+        sf_desc.material() = 
+            {gid, det.material_store().template size<gid>() - 1};
+
+        // Call the underlying volume builder(s) and give the volume to the
+        // next decorator
+        return volume_decorator<detector_t>::build(det, ctx);
+    }
+
+    private:
+    /// The surface this material map belongs to (index is volume local)
+    std::map<dindex, std::vector<bin_data>> bin_data_per_surface;
+    grid_factory<material_slab<scalar_type>, simple_serializer, replacer> m_factory{};
+    bin_filler_t m_bin_filler{};
+};
+
+}  // namespace detray
diff --git a/core/include/detray/utils/type_list.hpp b/core/include/detray/utils/type_list.hpp
index e5a6cefd2c..90bfa03c37 100644
--- a/core/include/detray/utils/type_list.hpp
+++ b/core/include/detray/utils/type_list.hpp
@@ -19,6 +19,7 @@
 namespace detray::types {
 
 /// @brief type list implementation
+/// @see https://www.codingwiththomas.com/blog/getting-started-with-typelists
 template <typename... Ts>
 struct type_list {};
 
diff --git a/io/include/detray/io/common/grid_reader.hpp b/io/include/detray/io/common/detail/grid_reader.hpp
similarity index 90%
rename from io/include/detray/io/common/grid_reader.hpp
rename to io/include/detray/io/common/detail/grid_reader.hpp
index 83838b5e04..c3436c0b66 100644
--- a/io/include/detray/io/common/grid_reader.hpp
+++ b/io/include/detray/io/common/detail/grid_reader.hpp
@@ -13,7 +13,6 @@
 #include "detray/io/common/io_interface.hpp"
 #include "detray/io/common/payloads.hpp"
 #include "detray/tools/detector_builder.hpp"
-#include "detray/tools/grid_builder.hpp"
 #include "detray/tools/grid_factory.hpp"
 #include "detray/utils/ranges.hpp"
 #include "detray/utils/type_list.hpp"
@@ -25,11 +24,11 @@
 #include <type_traits>
 #include <vector>
 
-namespace detray {
+namespace detray::detail {
 
 /// @brief Abstract base class for surface grid readers
-template <class detector_t,
-          typename value_t = typename detector_t::surface_type,
+template <class detector_t, typename value_t,
+          template<typename, typename, typename> class grid_builder_t,
           typename CAP = std::integral_constant<std::size_t, 9>,
           typename DIM = std::integral_constant<std::size_t, 2>>
 class grid_reader : public reader_interface<detector_t> {
@@ -41,30 +40,19 @@ class grid_reader : public reader_interface<detector_t> {
     using algebra_t = typename detector_t::transform3;
     using scalar_t = typename algebra_t::scalar_type;
 
-    // Bin content type to be written into the grid bin payloads
-    // For detector surface descriptors, write only the surface index
-    using content_t = std::conditional_t<
-        std::is_same_v<value_t, typename detector_t::surface_type>, std::size_t,
-        value_t>;
-
     static constexpr std::size_t Dim{DIM()};
     static constexpr std::size_t bin_capacity{CAP()};
 
-    protected:
-    /// Tag the reader as "surface_grids"
-    inline static const std::string tag = "surface_grids";
-
     public:
     /// Same constructors for this class as for base_type
     using base_type::base_type;
 
     protected:
-    /// Deserialize the detector grids @param grids_data from their IO
-    /// payload
+    /// Deserialize the detector grids @param grids_data from their IO payload
+    template <typename content_t>
     static void deserialize(
         detector_builder<typename detector_t::metadata, volume_builder>
             &det_builder,
-        typename detector_t::name_map &,
         const detector_grids_payload<content_t> &grids_data) {
 
         // Deserialize the grids volume by volume
@@ -200,7 +188,7 @@ class grid_reader : public reader_interface<detector_t> {
     ///
     /// @param grid_data grid IO payload (read from file)
     /// @param det_builder gather the grid data and build the final volume
-    template <typename bounds_ts, typename binning_ts,
+    template <typename bounds_ts, typename binning_ts, typename content_t,
               std::enable_if_t<types::size<bounds_ts> == Dim and
                                    types::size<binning_ts> == Dim,
                                bool> = true>
@@ -215,9 +203,8 @@ class grid_reader : public reader_interface<detector_t> {
                     "Unknown accelerator id in geometry file!");
             } else {
                 throw std::invalid_argument(
-                    "Given accelerator id could not be matched to a grid "
-                    "type: " +
-                    std::to_string(static_cast<std::int64_t>(acc_link)));
+                    "Given accelerator id could not be matched to a grid type: "
+                    + std::to_string(static_cast<std::int64_t>(acc_link)));
             }
         };
 
@@ -227,7 +214,7 @@ class grid_reader : public reader_interface<detector_t> {
         constexpr auto binnings = binning_ts{};
 
         // Check only 2-dimensional grid types
-        if constexpr (Dim >= 2) {
+        if constexpr (Dim == 2) {
             switch (grid_data.acc_link.type) {
                 // rectangle, trapezoid, (triangle) grids
                 case io::detail::acc_type::cartesian2_grid: {
@@ -272,7 +259,7 @@ class grid_reader : public reader_interface<detector_t> {
     }
 
     /// @brief End of recursion: build the grid from the @param grid_data
-    template <typename local_frame_t, typename... bounds_ts,
+    template <typename local_frame_t, typename content_t, typename... bounds_ts,
               typename... binning_ts,
               std::enable_if_t<sizeof...(bounds_ts) == Dim and
                                    sizeof...(binning_ts) == Dim,
@@ -295,17 +282,18 @@ class grid_reader : public reader_interface<detector_t> {
 
         // The compiler will instantiate this function for all possible types of
         // grids: Only proceed, if the grid type is known by the detector
-        if constexpr (detector_t::accel::template is_defined<grid_t>()) {
+        if constexpr (detector_t::accel::template is_defined<grid_t>() ||
+                      detector_t::materials::template is_defined<grid_t>()) {
 
             // Decorate the current volume builder with the grid
-            using grid_builder_t =
-                grid_builder<detector_t, grid_t, detail::fill_by_pos>;
+            using builder_t =
+                grid_builder_t<detector_t, grid_t, detail::fill_by_pos>;
 
             const auto volume_idx{
                 base_type::deserialize(grid_data.volume_link)};
             auto v_builder =
-                det_builder.template decorate<grid_builder_t>(volume_idx);
-            auto vgr_builder = dynamic_cast<grid_builder_t *>(v_builder);
+                det_builder.template decorate<builder_t>(volume_idx);
+            auto vgr_builder = dynamic_cast<builder_t *>(v_builder);
 
             // Initialize the grid axes
             std::vector<std::size_t> n_bins{};
@@ -324,7 +312,7 @@ class grid_reader : public reader_interface<detector_t> {
 
             vgr_builder->init_grid(spans, n_bins, ax_bin_edges);
 
-            value_t empty_sf{};
+            value_t entry{};
             n_axis::multi_bin<Dim> mbin;
             for (const auto &bin_data : grid_data.bins) {
 
@@ -340,9 +328,10 @@ class grid_reader : public reader_interface<detector_t> {
 
                 // For now assume surfaces ids as the only grid input
                 for (const auto c : bin_data.content) {
-                    empty_sf.set_volume(volume_idx);
-                    empty_sf.set_index(static_cast<dindex>(c));
-                    vgr_builder->get().populate(mbin, empty_sf);
+
+                    /*entry.set_volume(volume_idx);
+                    entry.set_index(static_cast<dindex>(c));*/
+                    vgr_builder->get().populate(mbin, serializer(entry));
                 }
             }
         }
diff --git a/io/include/detray/io/common/detector_reader.hpp b/io/include/detray/io/common/detector_reader.hpp
index 56c4171fe8..da6ffb4eff 100644
--- a/io/include/detray/io/common/detector_reader.hpp
+++ b/io/include/detray/io/common/detector_reader.hpp
@@ -97,9 +97,11 @@ auto assemble_reader(const io::detector_reader_config& cfg) noexcept(false) {
                 readers.template add<json_homogeneous_material_reader>(
                     file_name);
 
+            } else if (header.tag == "material_maps") {
+                readers.template add<json_material_maps_reader>(file_name);
+
             } else if (header.tag == "surface_grids") {
-                using surface_t = typename detector_t::surface_type;
-                readers.template add<json_grid_reader, surface_t,
+                readers.template add<json_surface_grid_reader,
                                      std::integral_constant<std::size_t, CAP>,
                                      std::integral_constant<std::size_t, DIM>>(
                     file_name);
diff --git a/io/include/detray/io/common/homogeneous_material_reader.hpp b/io/include/detray/io/common/homogeneous_material_reader.hpp
index ec782dedf2..6d9e6f91e9 100644
--- a/io/include/detray/io/common/homogeneous_material_reader.hpp
+++ b/io/include/detray/io/common/homogeneous_material_reader.hpp
@@ -24,6 +24,9 @@
 
 namespace detray {
 
+template <typename detector_t>
+class material_map_reader;
+
 /// @brief Abstract base class for a homogeneous material reader.
 template <class detector_t>
 class homogeneous_material_reader : public reader_interface<detector_t> {
@@ -34,6 +37,8 @@ class homogeneous_material_reader : public reader_interface<detector_t> {
     using material_type = io::detail::material_type;
     using scalar_type = typename detector_t::scalar_type;
 
+    friend class material_map_writer<detector_t>;
+
     protected:
     /// Tag the reader as "homogeneous material"
     inline static const std::string tag = "homogeneous_material";
diff --git a/io/include/detray/io/common/material_map_reader.hpp b/io/include/detray/io/common/material_map_reader.hpp
new file mode 100644
index 0000000000..438aa4a920
--- /dev/null
+++ b/io/include/detray/io/common/material_map_reader.hpp
@@ -0,0 +1,76 @@
+/** Detray library, part of the ACTS project (R&D line)
+ *
+ * (c) 2023 CERN for the benefit of the ACTS project
+ *
+ * Mozilla Public License Version 2.0
+ */
+
+#pragma once
+
+// Project include(s)
+#include "detray/io/common/detail/grid_reader.hpp"
+#include "detray/io/common/io_interface.hpp"
+#include "detray/io/common/payloads.hpp"
+#include "detray/materials/material_slab.hpp"
+#include "detray/tools/detector_builder.hpp"
+#include "detray/tools/material_map_builder.hpp"
+
+// System include(s)
+#include <stdexcept>
+#include <string>
+
+namespace detray {
+
+/// @brief Abstract base class for surface grid readers
+template <class detector_t,
+          typename DIM = std::integral_constant<std::size_t, 2>>
+class material_map_reader : public detail::grid_reader<detector_t, material_slab<typename detector_t::surface_type>, material_map_builder, 1u, DIM> {
+
+    using scalar_t = typename detector_t::surface_type;
+    using grid_reader_t = detail::grid_reader<detector_t, material_slab<scalar_t>, material_map_builder, 1u, DIM>;
+    using base_type = grid_reader_t;
+
+    using material_reader_t = homogeneous_material_reader<detector_t>;
+
+    protected:
+    /// Tag the reader as "surface_grids"
+    inline static const std::string tag = "material_maps";
+
+    public:
+    /// Same constructors for this class as for base_type
+    using base_type::base_type;
+
+    protected:
+    /// Deserialize the detector grids @param grids_data from their IO
+    /// payload
+    static void deserialize(
+        detector_builder<typename detector_t::metadata, volume_builder>
+            &det_builder,
+        typename detector_t::name_map &,
+        const detector_grids_payload<material_slab_payload> &grids_data) {
+
+        // Go through the full grid deserialization once volume material is 
+        // added, until then, giving the explicit bounds and binnings should be
+        // enough
+        using regular_t = regular<host_container_types, scalar_t>;
+        using bounds2D_ts = types::type_list<closed<0u>, closed<1u>>;
+        using binning2D_ts = types::type_list<regular_t, regular_t>;
+
+        for (const grid_payload<material_slab_payload>& g_data : grids_data) {
+            // Start form finding the grid local frame and then build the grid
+            if constexpr (DIM == 2) {
+                grid_reader_t::deserialize<bounds2D_ts, binning2D_ts>(g_data, det_builder);
+            } else if constexpr (DIM == 3) {
+                using bounds3D_ts = types::push_back<bounds2D_ts, closed<2u>>;
+                using binning3D_ts = types::push_back<binning2D_ts, regular_t>;
+
+                grid_reader_t::deserialize<bounds3D_ts, binning3D_ts>(g_data, det_builder);
+            } else {
+                throw std::invalid_argument("No 1D material grid type defined in detray");
+            }
+        }
+    }
+
+};
+
+}  // namespace detray
diff --git a/io/include/detray/io/common/surface_grid_reader.hpp b/io/include/detray/io/common/surface_grid_reader.hpp
new file mode 100644
index 0000000000..0d4db91497
--- /dev/null
+++ b/io/include/detray/io/common/surface_grid_reader.hpp
@@ -0,0 +1,53 @@
+/** Detray library, part of the ACTS project (R&D line)
+ *
+ * (c) 2023 CERN for the benefit of the ACTS project
+ *
+ * Mozilla Public License Version 2.0
+ */
+
+#pragma once
+
+// Project include(s)
+#include "detray/io/common/detail/grid_reader.hpp"
+#include "detray/io/common/io_interface.hpp"
+#include "detray/io/common/payloads.hpp"
+#include "detray/tools/grid_builder.hpp"
+#include "detray/tools/detector_builder.hpp"
+
+// System include(s)
+#include <string>
+
+namespace detray {
+
+/// @brief Abstract base class for surface grid readers
+template <class detector_t,
+          typename CAP = std::integral_constant<std::size_t, 9>,
+          typename DIM = std::integral_constant<std::size_t, 2>>
+class surface_grid_reader : public detail::grid_reader<detector_t, typename detector_t::surface_type, grid_builder, CAP, DIM> {
+
+    using grid_reader_t = detail::grid_reader<detector_t, typename detector_t::surface_type, grid_builder, CAP, DIM>;
+    using base_type = grid_reader_t;
+
+    protected:
+    /// Tag the reader as "surface_grids"
+    inline static const std::string tag = "surface_grids";
+
+    public:
+    /// Same constructors for this class as for base_type
+    using base_type::base_type;
+
+    protected:
+    /// Deserialize the detector grids @param grids_data from their IO
+    /// payload
+    static void deserialize(
+        detector_builder<typename detector_t::metadata, volume_builder>
+            &det_builder,
+        typename detector_t::name_map &,
+        const detector_grids_payload<std::size_t> &grids_data) {
+
+        grid_reader_t::deserialize(det_builder, grids_data);
+    }
+
+};
+
+}  // namespace detray
diff --git a/io/include/detray/io/json/json_reader.hpp b/io/include/detray/io/json/json_reader.hpp
index 6c62a2b078..066375dd6a 100644
--- a/io/include/detray/io/json/json_reader.hpp
+++ b/io/include/detray/io/json/json_reader.hpp
@@ -11,8 +11,9 @@
 #include "detray/io/common/detail/file_handle.hpp"
 #include "detray/io/common/detail/utils.hpp"
 #include "detray/io/common/geometry_reader.hpp"
-#include "detray/io/common/grid_reader.hpp"
+#include "detray/io/common/surface_grid_reader.hpp"
 #include "detray/io/common/homogeneous_material_reader.hpp"
+#include "detray/io/common/material_maps_reader.hpp"
 #include "detray/io/json/json.hpp"
 #include "detray/io/json/json_serializers.hpp"
 #include "detray/tools/detector_builder.hpp"
@@ -95,12 +96,17 @@ template <typename detector_t>
 using json_homogeneous_material_reader =
     json_reader<detector_t, homogeneous_material_reader>;
 
-/// Reads a homogeneous material descritption from file in json format
+/// Reads the material maps from file in json format
+template <typename detector_t,
+          typename DIM = std::integral_constant<std::size_t, 2>>
+using json_material_map_reader =
+    json_reader<detector_t, material_map_reader, DIM>;
+
+/// Reads the surface grids from file in json format
 template <typename detector_t,
-          typename value_t = typename detector_t::surface_type,
           typename CAP = std::integral_constant<std::size_t, 2>,
           typename DIM = std::integral_constant<std::size_t, 2>>
-using json_grid_reader =
-    json_reader<detector_t, grid_reader, value_t, CAP, DIM>;
+using json_surface_grid_reader =
+    json_reader<detector_t, surface_grid_reader, CAP, DIM>;
 
 }  // namespace detray
diff --git a/utils/include/detray/detectors/create_telescope_detector.hpp b/utils/include/detray/detectors/create_telescope_detector.hpp
index 9603269344..0400e64640 100644
--- a/utils/include/detray/detectors/create_telescope_detector.hpp
+++ b/utils/include/detray/detectors/create_telescope_detector.hpp
@@ -138,13 +138,12 @@ struct tel_det_config {
 }  // namespace
 
 /// Builds a detray geometry that contains only one volume with one type of
-/// surfaces, where the last surface is the portal that leaves the telescope.
-/// The detector is auto-constructed by following a trajectory state through
-/// space. The trajectory and the given distances determine the positions of
-/// the plane surfaces. The dis
+/// surfaces. The detector is auto-constructed by following a trajectory
+/// through space, along which the surfaces are placed. The portals are built
+/// from the bounding box around the sensors.
 ///
-/// @tparam mask_t the type of mask for the telescope surfaces
-/// @tparam trajectory_t the type of the pilot trajectory
+/// @tparam mask_shape_t the type of mask for the telescope surfaces
+/// @tparam trajectory_t the type of the pilot trajectory (ray/helix)
 ///
 /// @param resource the memory resource for the detector containers
 /// @param cfg configuration struct of the telescope detector
@@ -157,10 +156,11 @@ inline auto create_telescope_detector(
     const tel_det_config<mask_shape_t, trajectory_t> &cfg = {
         20.f * unit<scalar>::mm, 20.f * unit<scalar>::mm}) {
 
-    detector_builder<telescope_metadata<mask_shape_t>, volume_builder>
-        det_builder;
+    using builder_t = 
+        detector_builder<telescope_metadata<mask_shape_t>, volume_builder>;
+    using detector_t = typename builder_t::detector_type;
 
-    using detector_t = typename decltype(det_builder)::detector_type;
+    builder_t det_builder;
 
     // Detector and volume names
     typename detector_t::name_map name_map = {{0u, "telescope_detector"},