io/csv/intersection2D.hpp

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008
0009 #pragma once
0010
0011 // Project include(s)
0012 #include "detray/io/csv/dfe.hpp"
0013 #include "detray/io/utils/create_path.hpp"
0014 #include "detray/navigation/intersection/intersection.hpp"
0015 #include "detray/utils/ranges.hpp"
0016
0017 // System include(s).
0018 #include <cstdint>
0019 #include <filesystem>
0020
0021 namespace detray::io::csv {
0022
0023 /// Type to read the data of a line-surface intersection
0024 struct intersection2D {
0025   unsigned int track_id = 0;
0026   std::uint64_t identifier = 0ul;
0027   unsigned int type = 0u;
0028   unsigned int transform_index = 0u;
0029   unsigned int mask_id = 0u;
0030   unsigned int mask_index = 0u;
0031   unsigned int n_masks = 1u;
0032   unsigned int material_id = 0u;
0033   unsigned int material_index = 0u;
0034   double loc_0 = 0.;
0035   double loc_1 = 0.;
0036   double path = 0.;
0037   unsigned int volume_link = 0u;
0038   int direction = 0;
0039   int status = 0;
0040
0041   DFE_NAMEDTUPLE(intersection2D, track_id, identifier, type, transform_index,
0042                  mask_id, mask_index, material_id, material_index, loc_0, loc_1,
0043                  path, volume_link, direction, status);
0044 };
0045
0046 /// Read intersections from csv file
0047 /// @returns vector of intersections
0048 template <typename detector_t>
0049 inline auto read_intersection2D(const std::string &file_name) {
0050   using algebra_t = typename detector_t::algebra_type;
0051   using scalar_t = dscalar<algebra_t>;
0052   using surface_t = typename detector_t::surface_type;
0053   using nav_link_t = typename surface_t::navigation_link;
0054   using mask_link_t = typename surface_t::mask_link;
0055   using material_link_t = typename surface_t::material_link;
0056   using mask_id_t = typename detector_t::masks::id;
0057   using material_id_t = typename detector_t::material::id;
0058
0059   using intersection_t =
0060       detray::intersection2D<surface_t, algebra_t, intersection::contains_pos>;
0061
0062   dfe::NamedTupleCsvReader<io::csv::intersection2D> inters_reader(file_name);
0063
0064   io::csv::intersection2D inters_data{};
0065   std::vector<std::vector<intersection_t>> intersections_per_track;
0066
0067   // Read the data
0068   while (inters_reader.read(inters_data)) {
0069     // Add new intersection to correct track
0070     dindex trk_index{inters_data.track_id};
0071     while (intersections_per_track.size() <= trk_index) {
0072       intersections_per_track.push_back({});
0073     }
0074
0075     // Read the intersection
0076     intersection_t inters{};
0077
0078     mask_link_t mask_link{};
0079     mask_link.set_id(static_cast<mask_id_t>(inters_data.mask_id));
0080     if constexpr (detray::concepts::interval<
0081                       typename mask_link_t::index_type>) {
0082       mask_link.set_index({inters_data.mask_index, inters_data.n_masks});
0083     } else {
0084       mask_link.set_index(inters_data.mask_index);
0085     }
0086
0087     material_link_t material_link{
0088         static_cast<material_id_t>(inters_data.material_id),
0089         inters_data.material_index};
0090     inters.set_surface({inters_data.transform_index, mask_link, material_link,
0091                         dindex_invalid, surface_id::e_unknown});
0092     inters.surface().set_identifier(
0093         geometry::identifier{inters_data.identifier});
0094     inters.set_local({static_cast<scalar_t>(inters_data.loc_0),
0095                       static_cast<scalar_t>(inters_data.loc_1), 0.f});
0096     inters.set_path(static_cast<scalar_t>(inters_data.path));
0097     inters.set_volume_link(static_cast<nav_link_t>(inters_data.volume_link));
0098     inters.set_direction(static_cast<bool>(inters_data.direction));
0099     inters.set_status(static_cast<intersection::status>(inters_data.status));
0100
0101     // Add to collection
0102     intersections_per_track[trk_index].push_back(inters);
0103   }
0104
0105   // Check the result
0106   if (intersections_per_track.empty()) {
0107     throw std::invalid_argument(
0108         "ERROR: csv reader: Failed to read intersection data");
0109   }
0110
0111   return intersections_per_track;
0112 }
0113
0114 /// Write intersections to csv file
0115 template <typename intersection_t>
0116 inline void write_intersection2D(
0117     const std::string &file_name,
0118     const std::vector<std::vector<intersection_t>> &intersections_per_track,
0119     const bool replace = true) {
0120   using sf_desc_t = decltype(intersections_per_track.at(0).at(0).surface());
0121   using mask_link_t = typename sf_desc_t::mask_link;
0122
0123   // Don't write over existing data
0124   std::string inters_file_name{file_name};
0125   if (!replace && io::file_exists(file_name)) {
0126     inters_file_name = io::alt_file_name(file_name);
0127   } else {
0128     // Make sure the output directories exit
0129     io::create_path(std::filesystem::path{inters_file_name}.parent_path());
0130   }
0131
0132   dfe::NamedTupleCsvWriter<io::csv::intersection2D> inters_writer(
0133       inters_file_name);
0134
0135   for (const auto &[track_idx, intersections] :
0136        detray::views::enumerate(intersections_per_track)) {
0137     // Skip empty traces
0138     if (intersections.empty()) {
0139       continue;
0140     }
0141
0142     for (const auto &inters : intersections) {
0143       io::csv::intersection2D inters_data{};
0144
0145       inters_data.track_id = track_idx;
0146       inters_data.identifier = inters.surface().identifier().value();
0147       inters_data.type =
0148           static_cast<unsigned int>(inters.surface().identifier().id());
0149       inters_data.transform_index = inters.surface().transform();
0150       inters_data.mask_id =
0151           static_cast<unsigned int>(inters.surface().mask().id());
0152       const auto mask_index = inters.surface().mask().index();
0153       if constexpr (detray::concepts::interval<
0154                         typename mask_link_t::index_type>) {
0155         inters_data.mask_index = static_cast<unsigned int>(mask_index.lower());
0156         inters_data.n_masks = static_cast<unsigned int>(mask_index.size());
0157       } else {
0158         inters_data.mask_index = static_cast<unsigned int>(mask_index);
0159         inters_data.n_masks = 1u;
0160       }
0161       inters_data.material_id =
0162           static_cast<unsigned int>(inters.surface().material().id());
0163       inters_data.material_index = inters.surface().material().index();
0164       inters_data.loc_0 = inters.local()[0];
0165       inters_data.loc_1 = inters.local()[1];
0166       inters_data.path = inters.path();
0167       inters_data.volume_link = inters.volume_link();
0168       inters_data.direction = static_cast<int>(inters.is_along());
0169       inters_data.status = static_cast<int>(inters.status());
0170
0171       inters_writer.append(inters_data);
0172     }
0173   }
0174 }
0175
0176 }  // namespace detray::io::csv