Json/src/MaterialMapJsonConverter.cpp

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 #include "Acts/Plugins/Json/MaterialMapJsonConverter.hpp"
0010
0011 #include "Acts/Definitions/Tolerance.hpp"
0012 #include "Acts/Geometry/ApproachDescriptor.hpp"
0013 #include "Acts/Geometry/BoundarySurfaceT.hpp"
0014 #include "Acts/Geometry/CuboidVolumeBounds.hpp"
0015 #include "Acts/Geometry/CutoutCylinderVolumeBounds.hpp"
0016 #include "Acts/Geometry/CylinderVolumeBounds.hpp"
0017 #include "Acts/Geometry/GeometryHierarchyMap.hpp"
0018 #include "Acts/Geometry/GeometryIdentifier.hpp"
0019 #include "Acts/Geometry/Layer.hpp"
0020 #include "Acts/Geometry/TrackingGeometry.hpp"
0021 #include "Acts/Geometry/TrackingVolume.hpp"
0022 #include "Acts/Geometry/VolumeBounds.hpp"
0023 #include "Acts/Material/ISurfaceMaterial.hpp"
0024 #include "Acts/Material/IVolumeMaterial.hpp"
0025 #include "Acts/Material/ProtoSurfaceMaterial.hpp"
0026 #include "Acts/Material/ProtoVolumeMaterial.hpp"
0027 #include "Acts/Plugins/Json/ITrackingGeometryJsonDecorator.hpp"
0028 #include "Acts/Plugins/Json/IVolumeMaterialJsonDecorator.hpp"
0029 #include "Acts/Plugins/Json/MaterialJsonConverter.hpp"
0030 #include "Acts/Plugins/Json/SurfaceJsonConverter.hpp"
0031 #include "Acts/Plugins/Json/VolumeJsonConverter.hpp"
0032 #include "Acts/Surfaces/RectangleBounds.hpp"
0033 #include "Acts/Surfaces/Surface.hpp"
0034 #include "Acts/Surfaces/SurfaceArray.hpp"
0035 #include "Acts/Utilities/BinUtility.hpp"
0036 #include "Acts/Utilities/BinnedArray.hpp"
0037 #include "Acts/Utilities/BinningType.hpp"
0038 #include <Acts/Surfaces/AnnulusBounds.hpp>
0039 #include <Acts/Surfaces/CylinderBounds.hpp>
0040 #include <Acts/Surfaces/RadialBounds.hpp>
0041 #include <Acts/Surfaces/SurfaceBounds.hpp>
0042 #include <Acts/Surfaces/TrapezoidBounds.hpp>
0043
0044 #include <algorithm>
0045 #include <cmath>
0046 #include <cstddef>
0047 #include <map>
0048 #include <numbers>
0049 #include <stdexcept>
0050
0051 namespace Acts {
0052 // specialisations of decoration helper function
0053 // to pick correct objects from the container object
0054 template <>
0055 inline void decorateJson<Acts::SurfaceAndMaterialWithContext>(
0056     const ITrackingGeometryJsonDecorator* decorator,
0057     const Acts::SurfaceAndMaterialWithContext& src, nlohmann::json& dest) {
0058   if (decorator != nullptr && std::get<0>(src) != nullptr) {
0059     decorator->decorate(*std::get<0>(src), dest);
0060   }
0061 }
0062 template <>
0063 inline void decorateJson<Acts::TrackingVolumeAndMaterial>(
0064     const ITrackingGeometryJsonDecorator* decorator,
0065     const Acts::TrackingVolumeAndMaterial& src, nlohmann::json& dest) {
0066   if (decorator != nullptr && src.first != nullptr) {
0067     decorator->decorate(*src.first, dest);
0068   }
0069 }
0070
0071 template <>
0072 inline void decorateJson<Acts::IVolumeMaterial>(
0073     const IVolumeMaterialJsonDecorator* decorator,
0074     const Acts::IVolumeMaterial* src, nlohmann::json& dest) {
0075   if (decorator != nullptr && src != nullptr) {
0076     decorator->decorate(*src, dest);
0077   }
0078 }
0079 template <>
0080 inline void decorateJson<Acts::ISurfaceMaterial>(
0081     const IVolumeMaterialJsonDecorator* decorator,
0082     const Acts::ISurfaceMaterial* src, nlohmann::json& dest) {
0083   if (decorator != nullptr && src != nullptr) {
0084     decorator->decorate(*src, dest);
0085   }
0086 }
0087 }  // namespace Acts
0088
0089 namespace {
0090
0091 Acts::SurfaceAndMaterialWithContext defaultSurfaceMaterial(
0092     const std::shared_ptr<const Acts::Surface>& surface,
0093     const Acts::GeometryContext& context) {
0094   if (surface->surfaceMaterialSharedPtr() != nullptr) {
0095     return {surface, surface->surfaceMaterialSharedPtr(), context};
0096   }
0097   Acts::BinUtility bUtility;
0098   // Check which type of bounds is associated to the surface
0099   const Acts::SurfaceBounds& surfaceBounds = surface->bounds();
0100   const Acts::RadialBounds* radialBounds =
0101       dynamic_cast<const Acts::RadialBounds*>(&surfaceBounds);
0102   const Acts::CylinderBounds* cylinderBounds =
0103       dynamic_cast<const Acts::CylinderBounds*>(&surfaceBounds);
0104   const Acts::AnnulusBounds* annulusBounds =
0105       dynamic_cast<const Acts::AnnulusBounds*>(&surfaceBounds);
0106   const Acts::RectangleBounds* rectangleBounds =
0107       dynamic_cast<const Acts::RectangleBounds*>(&surfaceBounds);
0108   const Acts::TrapezoidBounds* trapezoidBounds =
0109       dynamic_cast<const Acts::TrapezoidBounds*>(&surfaceBounds);
0110
0111   if (radialBounds != nullptr) {
0112     bUtility += Acts::BinUtility(
0113         1,
0114         radialBounds->get(Acts::RadialBounds::eAveragePhi) -
0115             radialBounds->get(Acts::RadialBounds::eHalfPhiSector),
0116         radialBounds->get(Acts::RadialBounds::eAveragePhi) +
0117             radialBounds->get(Acts::RadialBounds::eHalfPhiSector),
0118         (radialBounds->get(Acts::RadialBounds::eHalfPhiSector) -
0119          std::numbers::pi) < Acts::s_epsilon
0120             ? Acts::closed
0121             : Acts::open,
0122         Acts::AxisDirection::AxisPhi);
0123     bUtility += Acts::BinUtility(1, radialBounds->rMin(), radialBounds->rMax(),
0124                                  Acts::open, Acts::AxisDirection::AxisR);
0125   }
0126   if (cylinderBounds != nullptr) {
0127     bUtility += Acts::BinUtility(
0128         1,
0129         cylinderBounds->get(Acts::CylinderBounds::eAveragePhi) -
0130             cylinderBounds->get(Acts::CylinderBounds::eHalfPhiSector),
0131         cylinderBounds->get(Acts::CylinderBounds::eAveragePhi) +
0132             cylinderBounds->get(Acts::CylinderBounds::eHalfPhiSector),
0133         (cylinderBounds->get(Acts::CylinderBounds::eHalfPhiSector) -
0134          std::numbers::pi) < Acts::s_epsilon
0135             ? Acts::closed
0136             : Acts::open,
0137         Acts::AxisDirection::AxisPhi);
0138     bUtility += Acts::BinUtility(
0139         1, -1 * cylinderBounds->get(Acts::CylinderBounds::eHalfLengthZ),
0140         cylinderBounds->get(Acts::CylinderBounds::eHalfLengthZ), Acts::open,
0141         Acts::AxisDirection::AxisZ);
0142   }
0143   if (annulusBounds != nullptr) {
0144     bUtility +=
0145         Acts::BinUtility(1, annulusBounds->get(Acts::AnnulusBounds::eMinPhiRel),
0146                          annulusBounds->get(Acts::AnnulusBounds::eMaxPhiRel),
0147                          Acts::open, Acts::AxisDirection::AxisPhi);
0148     bUtility += Acts::BinUtility(1, static_cast<float>(annulusBounds->rMin()),
0149                                  static_cast<float>(annulusBounds->rMax()),
0150                                  Acts::open, Acts::AxisDirection::AxisR);
0151   }
0152   if (rectangleBounds != nullptr) {
0153     bUtility +=
0154         Acts::BinUtility(1, rectangleBounds->get(Acts::RectangleBounds::eMinX),
0155                          rectangleBounds->get(Acts::RectangleBounds::eMaxX),
0156                          Acts::open, Acts::AxisDirection::AxisX);
0157     bUtility +=
0158         Acts::BinUtility(1, rectangleBounds->get(Acts::RectangleBounds::eMinY),
0159                          rectangleBounds->get(Acts::RectangleBounds::eMaxY),
0160                          Acts::open, Acts::AxisDirection::AxisY);
0161   }
0162   if (trapezoidBounds != nullptr) {
0163     double halfLengthX =
0164         std::max(trapezoidBounds->get(Acts::TrapezoidBounds::eHalfLengthXnegY),
0165                  trapezoidBounds->get(Acts::TrapezoidBounds::eHalfLengthXposY));
0166     bUtility += Acts::BinUtility(1, -1 * halfLengthX, halfLengthX, Acts::open,
0167                                  Acts::AxisDirection::AxisX);
0168     bUtility += Acts::BinUtility(
0169         1, -1 * trapezoidBounds->get(Acts::TrapezoidBounds::eHalfLengthY),
0170         trapezoidBounds->get(Acts::TrapezoidBounds::eHalfLengthY), Acts::open,
0171         Acts::AxisDirection::AxisY);
0172   }
0173   return {surface, std::make_shared<Acts::ProtoSurfaceMaterial>(bUtility),
0174           context};
0175 }
0176
0177 Acts::TrackingVolumeAndMaterial defaultVolumeMaterial(
0178     const Acts::TrackingVolume* volume) {
0179   Acts::BinUtility bUtility;
0180   if (volume->volumeMaterialPtr() != nullptr) {
0181     return {volume, volume->volumeMaterialPtr()};
0182   }
0183   // Check which type of bound is associated to the volume
0184   auto cyBounds = dynamic_cast<const Acts::CylinderVolumeBounds*>(
0185       &(volume->volumeBounds()));
0186   auto cutcylBounds = dynamic_cast<const Acts::CutoutCylinderVolumeBounds*>(
0187       &(volume->volumeBounds()));
0188   auto cuBounds =
0189       dynamic_cast<const Acts::CuboidVolumeBounds*>(&(volume->volumeBounds()));
0190
0191   if (cyBounds != nullptr) {
0192     bUtility +=
0193         Acts::BinUtility(1, cyBounds->get(Acts::CylinderVolumeBounds::eMinR),
0194                          cyBounds->get(Acts::CylinderVolumeBounds::eMaxR),
0195                          Acts::open, Acts::AxisDirection::AxisR);
0196     bUtility += Acts::BinUtility(
0197         1, -cyBounds->get(Acts::CylinderVolumeBounds::eHalfPhiSector),
0198         cyBounds->get(Acts::CylinderVolumeBounds::eHalfPhiSector),
0199         (cyBounds->get(Acts::CylinderVolumeBounds::eHalfPhiSector) -
0200          std::numbers::pi) < Acts::s_epsilon
0201             ? Acts::closed
0202             : Acts::open,
0203         Acts::AxisDirection::AxisPhi);
0204     bUtility += Acts::BinUtility(
0205         1, -cyBounds->get(Acts::CylinderVolumeBounds::eHalfLengthZ),
0206         cyBounds->get(Acts::CylinderVolumeBounds::eHalfLengthZ), Acts::open,
0207         Acts::AxisDirection::AxisZ);
0208   }
0209   if (cutcylBounds != nullptr) {
0210     bUtility += Acts::BinUtility(
0211         1, cutcylBounds->get(Acts::CutoutCylinderVolumeBounds::eMinR),
0212         cutcylBounds->get(Acts::CutoutCylinderVolumeBounds::eMaxR), Acts::open,
0213         Acts::AxisDirection::AxisR);
0214     bUtility += Acts::BinUtility(1, -std::numbers::pi_v<float>,
0215                                  std::numbers::pi_v<float>, Acts::closed,
0216                                  Acts::AxisDirection::AxisPhi);
0217     bUtility += Acts::BinUtility(
0218         1, -cutcylBounds->get(Acts::CutoutCylinderVolumeBounds::eHalfLengthZ),
0219         cutcylBounds->get(Acts::CutoutCylinderVolumeBounds::eHalfLengthZ),
0220         Acts::open, Acts::AxisDirection::AxisZ);
0221   } else if (cuBounds != nullptr) {
0222     bUtility += Acts::BinUtility(
0223         1, -cuBounds->get(Acts::CuboidVolumeBounds::eHalfLengthX),
0224         cuBounds->get(Acts::CuboidVolumeBounds::eHalfLengthX), Acts::open,
0225         Acts::AxisDirection::AxisX);
0226     bUtility += Acts::BinUtility(
0227         1, -cuBounds->get(Acts::CuboidVolumeBounds::eHalfLengthY),
0228         cuBounds->get(Acts::CuboidVolumeBounds::eHalfLengthY), Acts::open,
0229         Acts::AxisDirection::AxisY);
0230     bUtility += Acts::BinUtility(
0231         1, -cuBounds->get(Acts::CuboidVolumeBounds::eHalfLengthZ),
0232         cuBounds->get(Acts::CuboidVolumeBounds::eHalfLengthZ), Acts::open,
0233         Acts::AxisDirection::AxisZ);
0234   }
0235   return {volume, std::make_shared<Acts::ProtoVolumeMaterial>(bUtility)};
0236 }
0237 }  // namespace
0238
0239 Acts::MaterialMapJsonConverter::MaterialMapJsonConverter(
0240     const Acts::MaterialMapJsonConverter::Config& config,
0241     Acts::Logging::Level level)
0242     : m_cfg(config),
0243       m_logger{getDefaultLogger("MaterialMapJsonConverter", level)},
0244       m_volumeMaterialConverter(m_volumeName),
0245       m_volumeConverter(m_volumeName),
0246       m_surfaceMaterialConverter(m_surfaceName),
0247       m_surfaceConverter(m_surfaceName) {}
0248
0249 /// Convert method
0250 ///
0251 nlohmann::json Acts::MaterialMapJsonConverter::materialMapsToJson(
0252     const DetectorMaterialMaps& maps,
0253     const IVolumeMaterialJsonDecorator* decorator) {
0254   VolumeMaterialMap volumeMap = maps.second;
0255   std::vector<std::pair<GeometryIdentifier, const IVolumeMaterial*>>
0256       mapVolumeInit;
0257   for (const auto& [key, value] : volumeMap) {
0258     mapVolumeInit.push_back({key, value.get()});
0259   }
0260   GeometryHierarchyMap<const IVolumeMaterial*> hierarchyVolumeMap(
0261       mapVolumeInit);
0262   nlohmann::json materialVolume =
0263       m_volumeMaterialConverter.toJson(hierarchyVolumeMap, decorator);
0264   SurfaceMaterialMap surfaceMap = maps.first;
0265   std::vector<std::pair<GeometryIdentifier, const ISurfaceMaterial*>>
0266       mapSurfaceInit;
0267   for (const auto& [key, value] : surfaceMap) {
0268     mapSurfaceInit.push_back({key, value.get()});
0269   }
0270   GeometryHierarchyMap<const ISurfaceMaterial*> hierarchySurfaceMap(
0271       mapSurfaceInit);
0272   nlohmann::json materialSurface =
0273       m_surfaceMaterialConverter.toJson(hierarchySurfaceMap, decorator);
0274   nlohmann::json materialMap;
0275   materialMap["Volumes"] = materialVolume;
0276   materialMap["Surfaces"] = materialSurface;
0277   return materialMap;
0278 }
0279
0280 Acts::MaterialMapJsonConverter::DetectorMaterialMaps
0281 Acts::MaterialMapJsonConverter::jsonToMaterialMaps(
0282     const nlohmann::json& materialmap) {
0283   nlohmann::json materialVolume = materialmap["Volumes"];
0284   GeometryHierarchyMap<const IVolumeMaterial*> hierarchyVolumeMap =
0285       m_volumeMaterialConverter.fromJson(materialVolume);
0286   VolumeMaterialMap volumeMap;
0287   for (std::size_t i = 0; i < hierarchyVolumeMap.size(); i++) {
0288     std::shared_ptr<const IVolumeMaterial> volumePointer(
0289         hierarchyVolumeMap.valueAt(i));
0290     volumeMap.insert({hierarchyVolumeMap.idAt(i), std::move(volumePointer)});
0291   }
0292   nlohmann::json materialSurface = materialmap["Surfaces"];
0293   GeometryHierarchyMap<const ISurfaceMaterial*> hierarchySurfaceMap =
0294       m_surfaceMaterialConverter.fromJson(materialSurface);
0295   SurfaceMaterialMap surfaceMap;
0296   for (std::size_t i = 0; i < hierarchySurfaceMap.size(); i++) {
0297     std::shared_ptr<const ISurfaceMaterial> surfacePointer(
0298         hierarchySurfaceMap.valueAt(i));
0299     surfaceMap.insert({hierarchySurfaceMap.idAt(i), std::move(surfacePointer)});
0300   }
0301
0302   Acts::MaterialMapJsonConverter::DetectorMaterialMaps maps = {surfaceMap,
0303                                                                volumeMap};
0304
0305   // Return the filled maps
0306   return maps;
0307 }
0308
0309 nlohmann::json Acts::MaterialMapJsonConverter::trackingGeometryToJson(
0310     const Acts::TrackingGeometry& tGeometry,
0311     const ITrackingGeometryJsonDecorator* decorator) {
0312   std::vector<std::pair<GeometryIdentifier, Acts::TrackingVolumeAndMaterial>>
0313       volumeHierarchy;
0314   std::vector<
0315       std::pair<GeometryIdentifier, Acts::SurfaceAndMaterialWithContext>>
0316       surfaceHierarchy;
0317   convertToHierarchy(volumeHierarchy, surfaceHierarchy,
0318                      tGeometry.highestTrackingVolume());
0319   GeometryHierarchyMap<Acts::TrackingVolumeAndMaterial> hierarchyVolumeMap(
0320       volumeHierarchy);
0321   nlohmann::json jsonVolumes =
0322       m_volumeConverter.toJson(hierarchyVolumeMap, decorator);
0323   GeometryHierarchyMap<Acts::SurfaceAndMaterialWithContext> hierarchySurfaceMap(
0324       surfaceHierarchy);
0325   nlohmann::json jsonSurfaces =
0326       m_surfaceConverter.toJson(hierarchySurfaceMap, decorator);
0327   nlohmann::json hierarchyMap;
0328   hierarchyMap["Volumes"] = jsonVolumes;
0329   hierarchyMap["Surfaces"] = jsonSurfaces;
0330   return hierarchyMap;
0331 }
0332
0333 void Acts::MaterialMapJsonConverter::convertToHierarchy(
0334     std::vector<std::pair<GeometryIdentifier, Acts::TrackingVolumeAndMaterial>>&
0335         volumeHierarchy,
0336     std::vector<
0337         std::pair<GeometryIdentifier, Acts::SurfaceAndMaterialWithContext>>&
0338         surfaceHierarchy,
0339     const Acts::TrackingVolume* tVolume) {
0340   auto sameId = [tVolume](const auto& pair) {
0341     return (tVolume->geometryId() == pair.first);
0342   };
0343   if (std::ranges::any_of(volumeHierarchy, sameId)) {
0344     // this volume was already visited
0345     return;
0346   }
0347   if ((tVolume->volumeMaterial() != nullptr ||
0348        m_cfg.processNonMaterial == true) &&
0349       m_cfg.processVolumes == true) {
0350     volumeHierarchy.push_back(
0351         {tVolume->geometryId(), defaultVolumeMaterial(tVolume)});
0352   }
0353   // there are confined volumes
0354   if (tVolume->confinedVolumes() != nullptr) {
0355     // get through the volumes
0356     auto& volumes = tVolume->confinedVolumes()->arrayObjects();
0357     // loop over the volumes
0358     for (auto& vol : volumes) {
0359       // recursive call
0360       convertToHierarchy(volumeHierarchy, surfaceHierarchy, vol.get());
0361     }
0362   }
0363   // there are dense volumes
0364   if (m_cfg.processDenseVolumes && !tVolume->denseVolumes().empty()) {
0365     // loop over the volumes
0366     for (auto& vol : tVolume->denseVolumes()) {
0367       // recursive call
0368       convertToHierarchy(volumeHierarchy, surfaceHierarchy, vol.get());
0369     }
0370   }
0371   if (tVolume->confinedLayers() != nullptr) {
0372     // get the layers
0373     auto& layers = tVolume->confinedLayers()->arrayObjects();
0374     // loop over the layers
0375     for (auto& lay : layers) {
0376       if (m_cfg.processRepresenting == true) {
0377         auto& layRep = lay->surfaceRepresentation();
0378         if ((layRep.surfaceMaterial() != nullptr ||
0379              m_cfg.processNonMaterial == true) &&
0380             layRep.geometryId() != GeometryIdentifier()) {
0381           surfaceHierarchy.push_back(
0382               {layRep.geometryId(),
0383                defaultSurfaceMaterial(layRep.getSharedPtr(), m_cfg.context)});
0384         }
0385       }
0386       if (lay->approachDescriptor() != nullptr &&
0387           m_cfg.processApproaches == true) {
0388         for (auto& asf : lay->approachDescriptor()->containedSurfaces()) {
0389           if (asf->surfaceMaterial() != nullptr ||
0390               m_cfg.processNonMaterial == true) {
0391             surfaceHierarchy.push_back(
0392                 {asf->geometryId(),
0393                  defaultSurfaceMaterial(asf->getSharedPtr(), m_cfg.context)});
0394           }
0395         }
0396       }
0397       if (lay->surfaceArray() != nullptr && m_cfg.processSensitives == true) {
0398         for (auto& ssf : lay->surfaceArray()->surfaces()) {
0399           if (ssf->surfaceMaterial() != nullptr ||
0400               m_cfg.processNonMaterial == true) {
0401             auto sp = ssf->getSharedPtr();
0402             auto sm = defaultSurfaceMaterial(sp, m_cfg.context);
0403             auto id = ssf->geometryId();
0404
0405             std::pair p{id, sm};
0406             surfaceHierarchy.push_back(p);
0407           }
0408         }
0409       }
0410     }
0411   }
0412   // Let's finally check the boundaries
0413   for (auto& bsurf : tVolume->boundarySurfaces()) {
0414     // the surface representation
0415     auto& bssfRep = bsurf->surfaceRepresentation();
0416     if (bssfRep.geometryId().volume() == tVolume->geometryId().volume() &&
0417         m_cfg.processBoundaries == true) {
0418       if (bssfRep.surfaceMaterial() != nullptr ||
0419           m_cfg.processNonMaterial == true) {
0420         surfaceHierarchy.push_back(
0421             {bssfRep.geometryId(),
0422              defaultSurfaceMaterial(bssfRep.getSharedPtr(), m_cfg.context)});
0423       }
0424     }
0425   }
0426 }