Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
MaterialMapJsonConverter.cpp
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file MaterialMapJsonConverter.cpp
1 // This file is part of the Acts project.
2 //
3 // Copyright (C) 2017-2023 CERN for the benefit of the Acts project
4 //
5 // This Source Code Form is subject to the terms of the Mozilla Public
6 // License, v. 2.0. If a copy of the MPL was not distributed with this
7 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
8 
9 #include "Acts/Plugins/Json/MaterialMapJsonConverter.hpp"
10 
11 #include "Acts/Definitions/Tolerance.hpp"
12 #include "Acts/Geometry/ApproachDescriptor.hpp"
13 #include "Acts/Geometry/BoundarySurfaceT.hpp"
14 #include "Acts/Geometry/CuboidVolumeBounds.hpp"
15 #include "Acts/Geometry/CutoutCylinderVolumeBounds.hpp"
16 #include "Acts/Geometry/CylinderVolumeBounds.hpp"
17 #include "Acts/Geometry/GeometryHierarchyMap.hpp"
18 #include "Acts/Geometry/GeometryIdentifier.hpp"
19 #include "Acts/Geometry/Layer.hpp"
20 #include "Acts/Geometry/TrackingGeometry.hpp"
21 #include "Acts/Geometry/TrackingVolume.hpp"
22 #include "Acts/Geometry/VolumeBounds.hpp"
23 #include "Acts/Material/ISurfaceMaterial.hpp"
24 #include "Acts/Material/IVolumeMaterial.hpp"
25 #include "Acts/Material/ProtoSurfaceMaterial.hpp"
26 #include "Acts/Material/ProtoVolumeMaterial.hpp"
27 #include "Acts/Plugins/Json/ITrackingGeometryJsonDecorator.hpp"
28 #include "Acts/Plugins/Json/IVolumeMaterialJsonDecorator.hpp"
29 #include "Acts/Plugins/Json/MaterialJsonConverter.hpp"
30 #include "Acts/Plugins/Json/SurfaceJsonConverter.hpp"
31 #include "Acts/Plugins/Json/VolumeJsonConverter.hpp"
32 #include "Acts/Surfaces/RectangleBounds.hpp"
33 #include "Acts/Surfaces/Surface.hpp"
34 #include "Acts/Surfaces/SurfaceArray.hpp"
35 #include "Acts/Utilities/BinUtility.hpp"
36 #include "Acts/Utilities/BinnedArray.hpp"
37 #include "Acts/Utilities/BinningType.hpp"
38 #include <Acts/Surfaces/AnnulusBounds.hpp>
39 #include <Acts/Surfaces/CylinderBounds.hpp>
40 #include <Acts/Surfaces/RadialBounds.hpp>
41 #include <Acts/Surfaces/SurfaceBounds.hpp>
42 #include <Acts/Surfaces/TrapezoidBounds.hpp>
43 
44 #include <algorithm>
45 #include <cmath>
46 #include <cstddef>
47 #include <map>
48 #include <stdexcept>
49 
50 namespace Acts {
51 // specialisations of decoration helper function
52 // to pick correct objects from the container object
53 template <>
54 inline void decorateJson<Acts::SurfaceAndMaterialWithContext>(
55  const ITrackingGeometryJsonDecorator* decorator,
56  const Acts::SurfaceAndMaterialWithContext& src, nlohmann::json& dest) {
57  if (decorator != nullptr && std::get<0>(src) != nullptr) {
58  decorator->decorate(*std::get<0>(src), dest);
59  }
60 }
61 template <>
62 inline void decorateJson<Acts::TrackingVolumeAndMaterial>(
63  const ITrackingGeometryJsonDecorator* decorator,
64  const Acts::TrackingVolumeAndMaterial& src, nlohmann::json& dest) {
65  if (decorator != nullptr && src.first != nullptr) {
66  decorator->decorate(*src.first, dest);
67  }
68 }
69 
70 template <>
71 inline void decorateJson<Acts::IVolumeMaterial>(
72  const IVolumeMaterialJsonDecorator* decorator,
73  const Acts::IVolumeMaterial* src, nlohmann::json& dest) {
74  if (decorator != nullptr && src != nullptr) {
75  decorator->decorate(*src, dest);
76  }
77 }
78 template <>
79 inline void decorateJson<Acts::ISurfaceMaterial>(
80  const IVolumeMaterialJsonDecorator* decorator,
81  const Acts::ISurfaceMaterial* src, nlohmann::json& dest) {
82  if (decorator != nullptr && src != nullptr) {
83  decorator->decorate(*src, dest);
84  }
85 }
86 } // namespace Acts
87 
88 namespace {
89 
90 Acts::SurfaceAndMaterialWithContext defaultSurfaceMaterial(
91  const std::shared_ptr<const Acts::Surface>& surface,
92  const Acts::GeometryContext& context) {
93  if (surface->surfaceMaterialSharedPtr() != nullptr) {
94  return {surface, surface->surfaceMaterialSharedPtr(), context};
95  }
96  Acts::BinUtility bUtility;
97  // Check which type of bounds is associated to the surface
98  const Acts::SurfaceBounds& surfaceBounds = surface->bounds();
99  const Acts::RadialBounds* radialBounds =
100  dynamic_cast<const Acts::RadialBounds*>(&surfaceBounds);
101  const Acts::CylinderBounds* cylinderBounds =
102  dynamic_cast<const Acts::CylinderBounds*>(&surfaceBounds);
103  const Acts::AnnulusBounds* annulusBounds =
104  dynamic_cast<const Acts::AnnulusBounds*>(&surfaceBounds);
105  const Acts::RectangleBounds* rectangleBounds =
106  dynamic_cast<const Acts::RectangleBounds*>(&surfaceBounds);
107  const Acts::TrapezoidBounds* trapezoidBounds =
108  dynamic_cast<const Acts::TrapezoidBounds*>(&surfaceBounds);
109 
110  if (radialBounds != nullptr) {
111  bUtility += Acts::BinUtility(
112  1,
113  radialBounds->get(Acts::RadialBounds::eAveragePhi) -
114  radialBounds->get(Acts::RadialBounds::eHalfPhiSector),
115  radialBounds->get(Acts::RadialBounds::eAveragePhi) +
116  radialBounds->get(Acts::RadialBounds::eHalfPhiSector),
117  (radialBounds->get(Acts::RadialBounds::eHalfPhiSector) - M_PI) <
118  Acts::s_epsilon
119  ? Acts::closed
120  : Acts::open,
121  Acts::binPhi);
122  bUtility += Acts::BinUtility(1, radialBounds->rMin(), radialBounds->rMax(),
123  Acts::open, Acts::binR);
124  }
125  if (cylinderBounds != nullptr) {
126  bUtility += Acts::BinUtility(
127  1,
128  cylinderBounds->get(Acts::CylinderBounds::eAveragePhi) -
129  cylinderBounds->get(Acts::CylinderBounds::eHalfPhiSector),
130  cylinderBounds->get(Acts::CylinderBounds::eAveragePhi) +
131  cylinderBounds->get(Acts::CylinderBounds::eHalfPhiSector),
132  (cylinderBounds->get(Acts::CylinderBounds::eHalfPhiSector) - M_PI) <
133  Acts::s_epsilon
134  ? Acts::closed
135  : Acts::open,
136  Acts::binPhi);
137  bUtility += Acts::BinUtility(
138  1, -1 * cylinderBounds->get(Acts::CylinderBounds::eHalfLengthZ),
139  cylinderBounds->get(Acts::CylinderBounds::eHalfLengthZ), Acts::open,
140  Acts::binZ);
141  }
142  if (annulusBounds != nullptr) {
143  bUtility +=
144  Acts::BinUtility(1, annulusBounds->get(Acts::AnnulusBounds::eMinPhiRel),
145  annulusBounds->get(Acts::AnnulusBounds::eMaxPhiRel),
146  Acts::open, Acts::binPhi);
147  bUtility += Acts::BinUtility(1, annulusBounds->rMin(),
148  annulusBounds->rMax(), Acts::open, Acts::binR);
149  }
150  if (rectangleBounds != nullptr) {
151  bUtility +=
152  Acts::BinUtility(1, rectangleBounds->get(Acts::RectangleBounds::eMinX),
153  rectangleBounds->get(Acts::RectangleBounds::eMaxX),
154  Acts::open, Acts::binX);
155  bUtility +=
156  Acts::BinUtility(1, rectangleBounds->get(Acts::RectangleBounds::eMinY),
157  rectangleBounds->get(Acts::RectangleBounds::eMaxY),
158  Acts::open, Acts::binY);
159  }
160  if (trapezoidBounds != nullptr) {
161  double halfLengthX =
162  std::max(trapezoidBounds->get(Acts::TrapezoidBounds::eHalfLengthXnegY),
163  trapezoidBounds->get(Acts::TrapezoidBounds::eHalfLengthXposY));
164  bUtility += Acts::BinUtility(1, -1 * halfLengthX, halfLengthX, Acts::open,
165  Acts::binX);
166  bUtility += Acts::BinUtility(
167  1, -1 * trapezoidBounds->get(Acts::TrapezoidBounds::eHalfLengthY),
168  trapezoidBounds->get(Acts::TrapezoidBounds::eHalfLengthY), Acts::open,
169  Acts::binY);
170  }
171  return {surface, std::make_shared<Acts::ProtoSurfaceMaterial>(bUtility),
172  context};
173 }
174 
175 Acts::TrackingVolumeAndMaterial defaultVolumeMaterial(
176  const Acts::TrackingVolume* volume) {
177  Acts::BinUtility bUtility;
178  if (volume->volumeMaterialSharedPtr() != nullptr) {
179  return {volume, volume->volumeMaterialSharedPtr()};
180  }
181  // Check which type of bound is associated to the volume
182  auto cyBounds = dynamic_cast<const Acts::CylinderVolumeBounds*>(
183  &(volume->volumeBounds()));
184  auto cutcylBounds = dynamic_cast<const Acts::CutoutCylinderVolumeBounds*>(
185  &(volume->volumeBounds()));
186  auto cuBounds =
187  dynamic_cast<const Acts::CuboidVolumeBounds*>(&(volume->volumeBounds()));
188 
189  if (cyBounds != nullptr) {
190  bUtility +=
191  Acts::BinUtility(1, cyBounds->get(Acts::CylinderVolumeBounds::eMinR),
192  cyBounds->get(Acts::CylinderVolumeBounds::eMaxR),
193  Acts::open, Acts::binR);
194  bUtility += Acts::BinUtility(
195  1, -cyBounds->get(Acts::CylinderVolumeBounds::eHalfPhiSector),
196  cyBounds->get(Acts::CylinderVolumeBounds::eHalfPhiSector),
197  (cyBounds->get(Acts::CylinderVolumeBounds::eHalfPhiSector) - M_PI) <
198  Acts::s_epsilon
199  ? Acts::closed
200  : Acts::open,
201  Acts::binPhi);
202  bUtility += Acts::BinUtility(
203  1, -cyBounds->get(Acts::CylinderVolumeBounds::eHalfLengthZ),
204  cyBounds->get(Acts::CylinderVolumeBounds::eHalfLengthZ), Acts::open,
205  Acts::binZ);
206  }
207  if (cutcylBounds != nullptr) {
208  bUtility += Acts::BinUtility(
209  1, cutcylBounds->get(Acts::CutoutCylinderVolumeBounds::eMinR),
210  cutcylBounds->get(Acts::CutoutCylinderVolumeBounds::eMaxR), Acts::open,
211  Acts::binR);
212  bUtility += Acts::BinUtility(1, -M_PI, M_PI, Acts::closed, Acts::binPhi);
213  bUtility += Acts::BinUtility(
214  1, -cutcylBounds->get(Acts::CutoutCylinderVolumeBounds::eHalfLengthZ),
215  cutcylBounds->get(Acts::CutoutCylinderVolumeBounds::eHalfLengthZ),
216  Acts::open, Acts::binZ);
217  } else if (cuBounds != nullptr) {
218  bUtility += Acts::BinUtility(
219  1, -cuBounds->get(Acts::CuboidVolumeBounds::eHalfLengthX),
220  cuBounds->get(Acts::CuboidVolumeBounds::eHalfLengthX), Acts::open,
221  Acts::binX);
222  bUtility += Acts::BinUtility(
223  1, -cuBounds->get(Acts::CuboidVolumeBounds::eHalfLengthY),
224  cuBounds->get(Acts::CuboidVolumeBounds::eHalfLengthY), Acts::open,
225  Acts::binY);
226  bUtility += Acts::BinUtility(
227  1, -cuBounds->get(Acts::CuboidVolumeBounds::eHalfLengthZ),
228  cuBounds->get(Acts::CuboidVolumeBounds::eHalfLengthZ), Acts::open,
229  Acts::binZ);
230  }
231  return {volume, std::make_shared<Acts::ProtoVolumeMaterial>(bUtility)};
232 }
233 } // namespace
234 
235 Acts::MaterialMapJsonConverter::MaterialMapJsonConverter(
236  const Acts::MaterialMapJsonConverter::Config& config,
237  Acts::Logging::Level level)
238  : m_cfg(config),
239  m_logger{getDefaultLogger("MaterialMapJsonConverter", level)},
240  m_volumeMaterialConverter(m_volumeName),
241  m_volumeConverter(m_volumeName),
242  m_surfaceMaterialConverter(m_surfaceName),
243  m_surfaceConverter(m_surfaceName) {}
244 
247 nlohmann::json Acts::MaterialMapJsonConverter::materialMapsToJson(
248  const DetectorMaterialMaps& maps,
249  const IVolumeMaterialJsonDecorator* decorator) {
250  VolumeMaterialMap volumeMap = maps.second;
251  std::vector<std::pair<GeometryIdentifier, const IVolumeMaterial*>>
252  mapVolumeInit;
253  for (auto it = volumeMap.begin(); it != volumeMap.end(); it++) {
254  mapVolumeInit.push_back({it->first, it->second.get()});
255  }
256  GeometryHierarchyMap<const IVolumeMaterial*> hierarchyVolumeMap(
257  mapVolumeInit);
258  nlohmann::json materialVolume =
259  m_volumeMaterialConverter.toJson(hierarchyVolumeMap, decorator);
260  SurfaceMaterialMap surfaceMap = maps.first;
261  std::vector<std::pair<GeometryIdentifier, const ISurfaceMaterial*>>
262  mapSurfaceInit;
263  for (auto it = surfaceMap.begin(); it != surfaceMap.end(); it++) {
264  mapSurfaceInit.push_back({it->first, it->second.get()});
265  }
266  GeometryHierarchyMap<const ISurfaceMaterial*> hierarchySurfaceMap(
267  mapSurfaceInit);
268  nlohmann::json materialSurface =
269  m_surfaceMaterialConverter.toJson(hierarchySurfaceMap, decorator);
270  nlohmann::json materialMap;
271  materialMap["Volumes"] = materialVolume;
272  materialMap["Surfaces"] = materialSurface;
273  return materialMap;
274 }
275 
276 Acts::MaterialMapJsonConverter::DetectorMaterialMaps
277 Acts::MaterialMapJsonConverter::jsonToMaterialMaps(
278  const nlohmann::json& materialmap) {
279  nlohmann::json materialVolume = materialmap["Volumes"];
280  GeometryHierarchyMap<const IVolumeMaterial*> hierarchyVolumeMap =
281  m_volumeMaterialConverter.fromJson(materialVolume);
282  VolumeMaterialMap volumeMap;
283  for (size_t i = 0; i < hierarchyVolumeMap.size(); i++) {
284  std::shared_ptr<const IVolumeMaterial> volumePointer(
285  hierarchyVolumeMap.valueAt(i));
286  volumeMap.insert({hierarchyVolumeMap.idAt(i), std::move(volumePointer)});
287  }
288  nlohmann::json materialSurface = materialmap["Surfaces"];
289  GeometryHierarchyMap<const ISurfaceMaterial*> hierarchySurfaceMap =
290  m_surfaceMaterialConverter.fromJson(materialSurface);
291  SurfaceMaterialMap surfaceMap;
292  for (size_t i = 0; i < hierarchySurfaceMap.size(); i++) {
293  std::shared_ptr<const ISurfaceMaterial> surfacePointer(
294  hierarchySurfaceMap.valueAt(i));
295  surfaceMap.insert({hierarchySurfaceMap.idAt(i), std::move(surfacePointer)});
296  }
297 
298  Acts::MaterialMapJsonConverter::DetectorMaterialMaps maps = {surfaceMap,
299  volumeMap};
300 
301  // Return the filled maps
302  return maps;
303 }
304 
305 nlohmann::json Acts::MaterialMapJsonConverter::trackingGeometryToJson(
306  const Acts::TrackingGeometry& tGeometry,
307  const ITrackingGeometryJsonDecorator* decorator) {
308  std::vector<std::pair<GeometryIdentifier, Acts::TrackingVolumeAndMaterial>>
309  volumeHierarchy;
310  std::vector<
311  std::pair<GeometryIdentifier, Acts::SurfaceAndMaterialWithContext>>
312  surfaceHierarchy;
313  convertToHierarchy(volumeHierarchy, surfaceHierarchy,
314  tGeometry.highestTrackingVolume());
315  GeometryHierarchyMap<Acts::TrackingVolumeAndMaterial> hierarchyVolumeMap(
316  volumeHierarchy);
317  nlohmann::json jsonVolumes =
318  m_volumeConverter.toJson(hierarchyVolumeMap, decorator);
319  GeometryHierarchyMap<Acts::SurfaceAndMaterialWithContext> hierarchySurfaceMap(
320  surfaceHierarchy);
321  nlohmann::json jsonSurfaces =
322  m_surfaceConverter.toJson(hierarchySurfaceMap, decorator);
323  nlohmann::json hierarchyMap;
324  hierarchyMap["Volumes"] = jsonVolumes;
325  hierarchyMap["Surfaces"] = jsonSurfaces;
326  return hierarchyMap;
327 }
328 
329 void Acts::MaterialMapJsonConverter::convertToHierarchy(
330  std::vector<std::pair<GeometryIdentifier, Acts::TrackingVolumeAndMaterial>>&
331  volumeHierarchy,
332  std::vector<
333  std::pair<GeometryIdentifier, Acts::SurfaceAndMaterialWithContext>>&
334  surfaceHierarchy,
335  const Acts::TrackingVolume* tVolume) {
336  auto sameId =
337  [tVolume](
338  const std::pair<GeometryIdentifier, Acts::TrackingVolumeAndMaterial>&
339  pair) { return (tVolume->geometryId() == pair.first); };
340  if (std::find_if(volumeHierarchy.begin(), volumeHierarchy.end(), sameId) !=
341  volumeHierarchy.end()) {
342  // this volume was already visited
343  return;
344  }
345  if ((tVolume->volumeMaterial() != nullptr ||
346  m_cfg.processNonMaterial == true) &&
347  m_cfg.processVolumes == true) {
348  volumeHierarchy.push_back(
349  {tVolume->geometryId(), defaultVolumeMaterial(tVolume)});
350  }
351  // there are confined volumes
352  if (tVolume->confinedVolumes() != nullptr) {
353  // get through the volumes
354  auto& volumes = tVolume->confinedVolumes()->arrayObjects();
355  // loop over the volumes
356  for (auto& vol : volumes) {
357  // recursive call
358  convertToHierarchy(volumeHierarchy, surfaceHierarchy, vol.get());
359  }
360  }
361  // there are dense volumes
362  if (m_cfg.processDenseVolumes && !tVolume->denseVolumes().empty()) {
363  // loop over the volumes
364  for (auto& vol : tVolume->denseVolumes()) {
365  // recursive call
366  convertToHierarchy(volumeHierarchy, surfaceHierarchy, vol.get());
367  }
368  }
369  if (tVolume->confinedLayers() != nullptr) {
370  // get the layers
371  auto& layers = tVolume->confinedLayers()->arrayObjects();
372  // loop over the layers
373  for (auto& lay : layers) {
374  if (m_cfg.processRepresenting == true) {
375  auto& layRep = lay->surfaceRepresentation();
376  if ((layRep.surfaceMaterial() != nullptr ||
377  m_cfg.processNonMaterial == true) &&
378  layRep.geometryId() != GeometryIdentifier()) {
379  surfaceHierarchy.push_back(
380  {layRep.geometryId(),
381  defaultSurfaceMaterial(layRep.getSharedPtr(), m_cfg.context)});
382  }
383  }
384  if (lay->approachDescriptor() != nullptr &&
385  m_cfg.processApproaches == true) {
386  for (auto& asf : lay->approachDescriptor()->containedSurfaces()) {
387  if (asf->surfaceMaterial() != nullptr ||
388  m_cfg.processNonMaterial == true) {
389  surfaceHierarchy.push_back(
390  {asf->geometryId(),
391  defaultSurfaceMaterial(asf->getSharedPtr(), m_cfg.context)});
392  }
393  }
394  }
395  if (lay->surfaceArray() != nullptr && m_cfg.processSensitives == true) {
396  for (auto& ssf : lay->surfaceArray()->surfaces()) {
397  if (ssf->surfaceMaterial() != nullptr ||
398  m_cfg.processNonMaterial == true) {
399  auto sp = ssf->getSharedPtr();
400  auto sm = defaultSurfaceMaterial(sp, m_cfg.context);
401  auto id = ssf->geometryId();
402 
403  std::pair p{id, sm};
404  surfaceHierarchy.push_back(p);
405  }
406  }
407  }
408  }
409  }
410  // Let's finally check the boundaries
411  for (auto& bsurf : tVolume->boundarySurfaces()) {
412  // the surface representation
413  auto& bssfRep = bsurf->surfaceRepresentation();
414  if (bssfRep.geometryId().volume() == tVolume->geometryId().volume() &&
415  m_cfg.processBoundaries == true) {
416  if (bssfRep.surfaceMaterial() != nullptr ||
417  m_cfg.processNonMaterial == true) {
418  surfaceHierarchy.push_back(
419  {bssfRep.geometryId(),
420  defaultSurfaceMaterial(bssfRep.getSharedPtr(), m_cfg.context)});
421  }
422  }
423  }
424 }