src/Geometry/KDTreeTrackingGeometryBuilder.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/Geometry/KDTreeTrackingGeometryBuilder.hpp"
0010
0011 #include "Acts/Geometry/CylinderLayer.hpp"
0012 #include "Acts/Geometry/DiscLayer.hpp"
0013 #include "Acts/Geometry/Extent.hpp"
0014 #include "Acts/Geometry/ITrackingVolumeHelper.hpp"
0015 #include "Acts/Geometry/Layer.hpp"
0016 #include "Acts/Geometry/LayerCreator.hpp"
0017 #include "Acts/Geometry/Polyhedron.hpp"
0018 #include "Acts/Geometry/ProtoLayer.hpp"
0019 #include "Acts/Geometry/TrackingGeometry.hpp"
0020 #include "Acts/Geometry/TrackingVolume.hpp"
0021 #include "Acts/Geometry/Volume.hpp"
0022 #include "Acts/Geometry/VolumeBounds.hpp"
0023 #include "Acts/Surfaces/CylinderBounds.hpp"
0024 #include "Acts/Surfaces/RadialBounds.hpp"
0025 #include "Acts/Surfaces/Surface.hpp"
0026 #include "Acts/Surfaces/SurfaceArray.hpp"
0027 #include "Acts/Surfaces/SurfaceBounds.hpp"
0028 #include "Acts/Utilities/BinningType.hpp"
0029 #include "Acts/Utilities/RangeXD.hpp"
0030
0031 #include <cstddef>
0032 #include <optional>
0033 #include <ostream>
0034 #include <stdexcept>
0035 #include <utility>
0036
0037 Acts::KDTreeTrackingGeometryBuilder::KDTreeTrackingGeometryBuilder(
0038     const Acts::KDTreeTrackingGeometryBuilder::Config& cfg,
0039     std::unique_ptr<const Logger> logger)
0040     : m_cfg(cfg), m_logger(std::move(logger)) {
0041   m_cfg.protoDetector.harmonize();
0042 }
0043
0044 std::unique_ptr<const Acts::TrackingGeometry>
0045 Acts::KDTreeTrackingGeometryBuilder::trackingGeometry(
0046     const GeometryContext& gctx) const {
0047   using MeasuredSurface =
0048       std::pair<std::array<double, 2u>, std::shared_ptr<Surface>>;
0049   // Prepare all the surfaces
0050   std::vector<MeasuredSurface> surfacesMeasured;
0051   surfacesMeasured.reserve(m_cfg.surfaces.size());
0052   for (auto& s : m_cfg.surfaces) {
0053     auto ext = s->polyhedronRepresentation(gctx, 1u).extent();
0054     surfacesMeasured.push_back(MeasuredSurface{
0055         std::array<double, 2u>{ext.medium(AxisDirection::AxisZ),
0056                                ext.medium(AxisDirection::AxisR)},
0057         s});
0058   }
0059
0060   // Create the KDTree
0061   ACTS_INFO("Full KDTree has " << surfacesMeasured.size() << " surfaces.");
0062   SurfaceKDT surfaceKDT(std::move(surfacesMeasured));
0063
0064   // Walk through the proto builder
0065   auto protoWorld = m_cfg.protoDetector.worldVolume;
0066
0067   KDTreeTrackingGeometryBuilder::Cache cCache;
0068
0069   auto worldTrackingVolume =
0070       translateVolume(cCache, gctx, surfaceKDT, protoWorld);
0071
0072   ACTS_INFO("Retrieved " << cCache.surfaceCounter
0073                          << " surfaces from the KDTree.");
0074
0075   // return the geometry to the service
0076   return std::make_unique<const Acts::TrackingGeometry>(worldTrackingVolume);
0077 }
0078
0079 std::shared_ptr<Acts::TrackingVolume>
0080 Acts::KDTreeTrackingGeometryBuilder::translateVolume(
0081     Cache& cCache, const GeometryContext& gctx, const SurfaceKDT& kdt,
0082     const ProtoVolume& ptVolume, const std::string& indent) const {
0083   ACTS_DEBUG(indent << "Processing ProtoVolume: " << ptVolume.name);
0084   std::vector<std::shared_ptr<const TrackingVolume>> translatedVolumes = {};
0085
0086   // Volume extent
0087   auto rangeR = ptVolume.extent.range(Acts::AxisDirection::AxisR);
0088   auto rangeZ = ptVolume.extent.range(Acts::AxisDirection::AxisZ);
0089
0090   // Simple gap volume
0091   if (!ptVolume.container.has_value()) {
0092     ACTS_VERBOSE(indent << "> empty volume to be built");
0093     MutableTrackingVolumeVector mtv = {};
0094     auto tVolume = m_cfg.trackingVolumeHelper->createGapTrackingVolume(
0095         gctx, mtv, nullptr, rangeR.min(), rangeR.max(), rangeZ.min(),
0096         rangeZ.max(), 0, false, ptVolume.name);
0097     ACTS_DEBUG(indent << "> translated into gap volume bounds: "
0098                       << tVolume->volumeBounds());
0099     return tVolume;
0100   } else {
0101     // Get the container information
0102     auto& cts = ptVolume.container.value();
0103
0104     // Container information must be present
0105     if (cts.constituentVolumes.empty()) {
0106       throw std::invalid_argument(
0107           "KDTreeTrackingGeometryBuilder: no constituents given.");
0108     }
0109
0110     // This volume is a volume container
0111     if (!cts.layerContainer) {
0112       ACTS_VERBOSE(indent << "> volume container with "
0113                           << cts.constituentVolumes.size() << " constituents.");
0114       for (auto& cVolume : cts.constituentVolumes) {
0115         auto dtVolume = translateVolume(cCache, gctx, kdt, cVolume,
0116                                         indent + m_cfg.hierarchyIndent);
0117         translatedVolumes.push_back(dtVolume);
0118       }
0119       auto tVolume = m_cfg.trackingVolumeHelper->createContainerTrackingVolume(
0120           gctx, translatedVolumes);
0121       ACTS_DEBUG(indent << "> translated into container volume bounds: "
0122                         << tVolume->volumeBounds());
0123       return tVolume;
0124     } else {
0125       // This volume is a layer container
0126       std::vector<std::shared_ptr<const Layer>> layers = {};
0127       ACTS_VERBOSE(indent << "> layer container with "
0128                           << cts.constituentVolumes.size() << " layers.");
0129       for (auto& plVolume : cts.constituentVolumes) {
0130         if (plVolume.internal.has_value()) {
0131           layers.push_back(translateLayer(cCache, gctx, kdt, plVolume,
0132                                           indent + m_cfg.hierarchyIndent));
0133         } else {
0134           ACTS_WARNING(indent << "> layer type volume has no internal "
0135                                  "description, layer not built.");
0136         }
0137       }
0138       // Create a new tracking volume with those layers
0139       auto tVolume = m_cfg.trackingVolumeHelper->createTrackingVolume(
0140           gctx, layers, {}, nullptr, rangeR.min(), rangeR.max(), rangeZ.min(),
0141           rangeZ.max(), ptVolume.name);
0142       ACTS_DEBUG(indent << "> translated into bounds: "
0143                         << tVolume->volumeBounds());
0144       return tVolume;
0145     }
0146   }
0147 }
0148
0149 /// @return a new tracking volume
0150 std::shared_ptr<const Acts::Layer>
0151 Acts::KDTreeTrackingGeometryBuilder::translateLayer(
0152     Cache& cCache, const GeometryContext& gctx, const SurfaceKDT& kdt,
0153     const ProtoVolume& plVolume, const std::string& indent) const {
0154   ACTS_DEBUG(indent + "Processing ProtoVolume: " << plVolume.name);
0155
0156   // This is only called if the volume has internal structure
0157   auto& its = plVolume.internal.value();
0158
0159   // Try to pull from the kd tree
0160   RangeXD<2u, double> zrRange;
0161   zrRange[0u] = plVolume.extent.range(Acts::AxisDirection::AxisZ);
0162   zrRange[1u] = plVolume.extent.range(Acts::AxisDirection::AxisR);
0163
0164   auto layerSurfaces = kdt.rangeSearchWithKey(zrRange);
0165   ACTS_VERBOSE(indent + ">> looking z/r range = " << zrRange.toString());
0166   ACTS_VERBOSE(indent + ">> found " << layerSurfaces.size()
0167                                     << " surfaces in the KDTree.");
0168   cCache.surfaceCounter += layerSurfaces.size();
0169
0170   std::shared_ptr<const Acts::Layer> tLayer = nullptr;
0171
0172   if (layerSurfaces.size() == 1u) {
0173     auto surface = layerSurfaces[0u].second;
0174     const auto& transform = surface->transform(gctx);
0175     if (its.layerType == Acts::Surface::SurfaceType::Cylinder) {
0176       ACTS_VERBOSE(indent +
0177                    ">> creating cylinder layer from a single surface.");
0178       // Get the bounds
0179       auto cylinderBounds =
0180           dynamic_cast<const CylinderBounds*>(&(surface->bounds()));
0181       auto cylinderBoundsClone =
0182           std::make_shared<const CylinderBounds>(*cylinderBounds);
0183       auto cylinderLayer =
0184           CylinderLayer::create(transform, cylinderBoundsClone, nullptr, 1.);
0185       cylinderLayer->assignSurfaceMaterial(surface->surfaceMaterialSharedPtr());
0186       tLayer = cylinderLayer;
0187     } else if (its.layerType == Acts::Surface::SurfaceType::Disc) {
0188       ACTS_VERBOSE(indent +
0189                    ">> creating cylinder layer from a single surface.");
0190       // Get the bounds
0191       auto radialBounds =
0192           dynamic_cast<const RadialBounds*>(&(surface->bounds()));
0193       auto radialBoundsClone =
0194           std::make_shared<const RadialBounds>(*radialBounds);
0195       auto discLayer =
0196           DiscLayer::create(transform, radialBoundsClone, nullptr, 1.);
0197       discLayer->assignSurfaceMaterial(surface->surfaceMaterialSharedPtr());
0198       tLayer = discLayer;
0199     } else {
0200       throw std::invalid_argument(
0201           "KDTreeTrackingGeometryBuilder: layer type is neither cylinder nor "
0202           "disk.");
0203     }
0204
0205   } else if (layerSurfaces.size() > 1u) {
0206     // Make a const collection out of the surfaces
0207     std::vector<std::shared_ptr<const Surface>> cLayerSurfaces;
0208     cLayerSurfaces.reserve(layerSurfaces.size());
0209     for (const auto& s : layerSurfaces) {
0210       cLayerSurfaces.push_back(s.second);
0211     }
0212
0213     Acts::BinningType bType0 = Acts::equidistant;
0214     Acts::BinningType bType1 = Acts::equidistant;
0215     std::size_t bins0 = 0;
0216     std::size_t bins1 = 0;
0217     // In case explicit binning is given
0218     if (its.surfaceBinning.size() == 2u) {
0219       bType0 = its.surfaceBinning[0u].type;
0220       bType1 = its.surfaceBinning[1u].type;
0221       // In case explicit bin numbers are given in addition
0222       if (bType0 == Acts::equidistant && bType1 == Acts::equidistant &&
0223           its.surfaceBinning[0u].bins() > 1u &&
0224           its.surfaceBinning[1u].bins() > 1u) {
0225         bins0 = its.surfaceBinning[0u].bins();
0226         bins1 = its.surfaceBinning[1u].bins();
0227         ACTS_VERBOSE(indent + ">> binning provided externally to be "
0228                      << bins0 << " x " << bins1 << ".");
0229       }
0230     }
0231
0232     Acts::ProtoLayer pLayer(gctx, cLayerSurfaces);
0233     pLayer.envelope = plVolume.extent.envelope();
0234     if (its.layerType == Acts::Surface::SurfaceType::Cylinder) {
0235       ACTS_VERBOSE(indent + ">> creating cylinder layer with "
0236                    << cLayerSurfaces.size() << " surfaces.");
0237       // Forced equidistant or auto-binned
0238       tLayer = (bins0 * bins1 > 0)
0239                    ? m_cfg.layerCreator->cylinderLayer(gctx, cLayerSurfaces,
0240                                                        bins0, bins1, pLayer)
0241                    : m_cfg.layerCreator->cylinderLayer(gctx, cLayerSurfaces,
0242                                                        bType0, bType1, pLayer);
0243
0244     } else if (its.layerType == Acts::Surface::SurfaceType::Disc) {
0245       ACTS_VERBOSE(indent + ">> creating disc layer with "
0246                    << cLayerSurfaces.size() << " surfaces.");
0247       // Forced equidistant or auto-binned
0248       tLayer = (bins0 * bins1 > 0)
0249                    ? m_cfg.layerCreator->discLayer(gctx, cLayerSurfaces, bins0,
0250                                                    bins1, pLayer)
0251
0252                    : m_cfg.layerCreator->discLayer(gctx, cLayerSurfaces, bType0,
0253                                                    bType1, pLayer);
0254     } else {
0255       throw std::invalid_argument(
0256           "KDTreeTrackingGeometryBuilder: layer type is neither cylinder nor "
0257           "disk.");
0258     }
0259   }
0260   if (tLayer != nullptr && tLayer->representingVolume() != nullptr) {
0261     ACTS_DEBUG(indent << "> translated into layer bounds: "
0262                       << tLayer->representingVolume()->volumeBounds());
0263   } else {
0264     throw std::runtime_error(
0265         "KDTreeTrackingGeometryBuilder: layer was not built.");
0266   }
0267
0268   return tLayer;
0269 }