Json/src/MaterialMapJsonConverter.cpp

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