src/Geometry/KDTreeTrackingGeometryBuilder.cpp

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