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/Utilities/BinningType.hpp"
0028 #include "Acts/Utilities/RangeXD.hpp"
0029
0030 #include <cstddef>
0031 #include <optional>
0032 #include <ostream>
0033 #include <stdexcept>
0034 #include <utility>
0035
0036 namespace Acts {
0037
0038 KDTreeTrackingGeometryBuilder::KDTreeTrackingGeometryBuilder(
0039     const KDTreeTrackingGeometryBuilder::Config& cfg,
0040     std::unique_ptr<const Logger> logger)
0041     : m_cfg(cfg), m_logger(std::move(logger)) {
0042   m_cfg.protoDetector.harmonize();
0043 }
0044
0045 std::unique_ptr<const TrackingGeometry>
0046 KDTreeTrackingGeometryBuilder::trackingGeometry(
0047     const GeometryContext& gctx) const {
0048   using MeasuredSurface =
0049       std::pair<std::array<double, 2u>, std::shared_ptr<Surface>>;
0050   // Prepare all the surfaces
0051   std::vector<MeasuredSurface> surfacesMeasured;
0052   surfacesMeasured.reserve(m_cfg.surfaces.size());
0053   for (auto& s : m_cfg.surfaces) {
0054     auto ext = s->polyhedronRepresentation(gctx, 1u).extent();
0055     surfacesMeasured.push_back(MeasuredSurface{
0056         std::array<double, 2u>{ext.medium(AxisDirection::AxisZ),
0057                                ext.medium(AxisDirection::AxisR)},
0058         s});
0059   }
0060
0061   // Create the KDTree
0062   ACTS_INFO("Full KDTree has " << surfacesMeasured.size() << " surfaces.");
0063   SurfaceKDT surfaceKDT(std::move(surfacesMeasured));
0064
0065   // Walk through the proto builder
0066   auto protoWorld = m_cfg.protoDetector.worldVolume;
0067
0068   KDTreeTrackingGeometryBuilder::Cache cCache;
0069
0070   auto worldTrackingVolume =
0071       translateVolume(cCache, gctx, surfaceKDT, protoWorld);
0072
0073   ACTS_INFO("Retrieved " << cCache.surfaceCounter
0074                          << " surfaces from the KDTree.");
0075
0076   // return the geometry to the service
0077   return std::make_unique<const TrackingGeometry>(worldTrackingVolume);
0078 }
0079
0080 std::shared_ptr<TrackingVolume> 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(AxisDirection::AxisR);
0088   auto rangeZ = ptVolume.extent.range(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 Layer> KDTreeTrackingGeometryBuilder::translateLayer(
0151     Cache& cCache, const GeometryContext& gctx, const SurfaceKDT& kdt,
0152     const ProtoVolume& plVolume, const std::string& indent) const {
0153   ACTS_DEBUG(indent + "Processing ProtoVolume: " << plVolume.name);
0154
0155   // This is only called if the volume has internal structure
0156   auto& its = plVolume.internal.value();
0157
0158   // Try to pull from the kd tree
0159   RangeXD<2u, double> zrRange;
0160   zrRange[0u] = plVolume.extent.range(AxisDirection::AxisZ);
0161   zrRange[1u] = plVolume.extent.range(AxisDirection::AxisR);
0162
0163   auto layerSurfaces = kdt.rangeSearchWithKey(zrRange);
0164   ACTS_VERBOSE(indent + ">> looking z/r range = " << zrRange.toString());
0165   ACTS_VERBOSE(indent + ">> found " << layerSurfaces.size()
0166                                     << " surfaces in the KDTree.");
0167   cCache.surfaceCounter += layerSurfaces.size();
0168
0169   std::shared_ptr<const Layer> tLayer = nullptr;
0170
0171   if (layerSurfaces.size() == 1u) {
0172     auto surface = layerSurfaces[0u].second;
0173     const auto& transform = surface->transform(gctx);
0174     if (its.layerType == Surface::SurfaceType::Cylinder) {
0175       ACTS_VERBOSE(indent +
0176                    ">> creating cylinder layer from a single surface.");
0177       // Get the bounds
0178       auto cylinderBounds =
0179           dynamic_cast<const CylinderBounds*>(&(surface->bounds()));
0180       auto cylinderBoundsClone =
0181           std::make_shared<const CylinderBounds>(*cylinderBounds);
0182       auto cylinderLayer =
0183           CylinderLayer::create(transform, cylinderBoundsClone, nullptr, 1.);
0184       cylinderLayer->assignSurfaceMaterial(surface->surfaceMaterialSharedPtr());
0185       tLayer = cylinderLayer;
0186     } else if (its.layerType == Surface::SurfaceType::Disc) {
0187       ACTS_VERBOSE(indent +
0188                    ">> creating cylinder layer from a single surface.");
0189       // Get the bounds
0190       auto radialBounds =
0191           dynamic_cast<const RadialBounds*>(&(surface->bounds()));
0192       auto radialBoundsClone =
0193           std::make_shared<const RadialBounds>(*radialBounds);
0194       auto discLayer =
0195           DiscLayer::create(transform, radialBoundsClone, nullptr, 1.);
0196       discLayer->assignSurfaceMaterial(surface->surfaceMaterialSharedPtr());
0197       tLayer = discLayer;
0198     } else {
0199       throw std::invalid_argument(
0200           "KDTreeTrackingGeometryBuilder: layer type is neither cylinder nor "
0201           "disk.");
0202     }
0203
0204   } else if (layerSurfaces.size() > 1u) {
0205     // Make a const collection out of the surfaces
0206     std::vector<std::shared_ptr<const Surface>> cLayerSurfaces;
0207     cLayerSurfaces.reserve(layerSurfaces.size());
0208     for (const auto& s : layerSurfaces) {
0209       cLayerSurfaces.push_back(s.second);
0210     }
0211
0212     BinningType bType0 = equidistant;
0213     BinningType bType1 = equidistant;
0214     std::size_t bins0 = 0;
0215     std::size_t bins1 = 0;
0216     // In case explicit binning is given
0217     if (its.surfaceBinning.size() == 2u) {
0218       bType0 = its.surfaceBinning[0u].type;
0219       bType1 = its.surfaceBinning[1u].type;
0220       // In case explicit bin numbers are given in addition
0221       if (bType0 == equidistant && bType1 == equidistant &&
0222           its.surfaceBinning[0u].bins() > 1u &&
0223           its.surfaceBinning[1u].bins() > 1u) {
0224         bins0 = its.surfaceBinning[0u].bins();
0225         bins1 = its.surfaceBinning[1u].bins();
0226         ACTS_VERBOSE(indent + ">> binning provided externally to be "
0227                      << bins0 << " x " << bins1 << ".");
0228       }
0229     }
0230
0231     ProtoLayer pLayer(gctx, cLayerSurfaces);
0232     pLayer.envelope = plVolume.extent.envelope();
0233     if (its.layerType == Surface::SurfaceType::Cylinder) {
0234       ACTS_VERBOSE(indent + ">> creating cylinder layer with "
0235                    << cLayerSurfaces.size() << " surfaces.");
0236       // Forced equidistant or auto-binned
0237       tLayer = (bins0 * bins1 > 0)
0238                    ? m_cfg.layerCreator->cylinderLayer(gctx, cLayerSurfaces,
0239                                                        bins0, bins1, pLayer)
0240                    : m_cfg.layerCreator->cylinderLayer(gctx, cLayerSurfaces,
0241                                                        bType0, bType1, pLayer);
0242
0243     } else if (its.layerType == Surface::SurfaceType::Disc) {
0244       ACTS_VERBOSE(indent + ">> creating disc layer with "
0245                    << cLayerSurfaces.size() << " surfaces.");
0246       // Forced equidistant or auto-binned
0247       tLayer = (bins0 * bins1 > 0)
0248                    ? m_cfg.layerCreator->discLayer(gctx, cLayerSurfaces, bins0,
0249                                                    bins1, pLayer)
0250
0251                    : m_cfg.layerCreator->discLayer(gctx, cLayerSurfaces, bType0,
0252                                                    bType1, pLayer);
0253     } else {
0254       throw std::invalid_argument(
0255           "KDTreeTrackingGeometryBuilder: layer type is neither cylinder nor "
0256           "disk.");
0257     }
0258   }
0259   if (tLayer != nullptr && tLayer->representingVolume() != nullptr) {
0260     ACTS_DEBUG(indent << "> translated into layer bounds: "
0261                       << tLayer->representingVolume()->volumeBounds());
0262   } else {
0263     throw std::runtime_error(
0264         "KDTreeTrackingGeometryBuilder: layer was not built.");
0265   }
0266
0267   return tLayer;
0268 }
0269
0270 }  // namespace Acts