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/Utilities/BinningType.hpp"
0027 #include "Acts/Utilities/RangeXD.hpp"
0028
0029 #include <cstddef>
0030 #include <optional>
0031 #include <ostream>
0032 #include <stdexcept>
0033 #include <utility>
0034
0035 namespace Acts {
0036
0037 KDTreeTrackingGeometryBuilder::KDTreeTrackingGeometryBuilder(
0038     const 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 TrackingGeometry>
0045 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 TrackingGeometry>(worldTrackingVolume);
0077 }
0078
0079 std::shared_ptr<TrackingVolume> KDTreeTrackingGeometryBuilder::translateVolume(
0080     Cache& cCache, const GeometryContext& gctx, const SurfaceKDT& kdt,
0081     const ProtoVolume& ptVolume, const std::string& indent) const {
0082   ACTS_DEBUG(indent << "Processing ProtoVolume: " << ptVolume.name);
0083   std::vector<std::shared_ptr<const TrackingVolume>> translatedVolumes = {};
0084
0085   // Volume extent
0086   auto rangeR = ptVolume.extent.range(AxisDirection::AxisR);
0087   auto rangeZ = ptVolume.extent.range(AxisDirection::AxisZ);
0088
0089   // Simple gap volume
0090   if (!ptVolume.container.has_value()) {
0091     ACTS_VERBOSE(indent << "> empty volume to be built");
0092     MutableTrackingVolumeVector mtv = {};
0093     auto tVolume = m_cfg.trackingVolumeHelper->createGapTrackingVolume(
0094         gctx, mtv, nullptr, rangeR.min(), rangeR.max(), rangeZ.min(),
0095         rangeZ.max(), 0, false, ptVolume.name);
0096     ACTS_DEBUG(indent << "> translated into gap volume bounds: "
0097                       << tVolume->volumeBounds());
0098     return tVolume;
0099   } else {
0100     // Get the container information
0101     auto& cts = ptVolume.container.value();
0102
0103     // Container information must be present
0104     if (cts.constituentVolumes.empty()) {
0105       throw std::invalid_argument(
0106           "KDTreeTrackingGeometryBuilder: no constituents given.");
0107     }
0108
0109     // This volume is a volume container
0110     if (!cts.layerContainer) {
0111       ACTS_VERBOSE(indent << "> volume container with "
0112                           << cts.constituentVolumes.size() << " constituents.");
0113       for (auto& cVolume : cts.constituentVolumes) {
0114         auto dtVolume = translateVolume(cCache, gctx, kdt, cVolume,
0115                                         indent + m_cfg.hierarchyIndent);
0116         translatedVolumes.push_back(dtVolume);
0117       }
0118       auto tVolume = m_cfg.trackingVolumeHelper->createContainerTrackingVolume(
0119           gctx, translatedVolumes);
0120       ACTS_DEBUG(indent << "> translated into container volume bounds: "
0121                         << tVolume->volumeBounds());
0122       return tVolume;
0123     } else {
0124       // This volume is a layer container
0125       std::vector<std::shared_ptr<const Layer>> layers = {};
0126       ACTS_VERBOSE(indent << "> layer container with "
0127                           << cts.constituentVolumes.size() << " layers.");
0128       for (auto& plVolume : cts.constituentVolumes) {
0129         if (plVolume.internal.has_value()) {
0130           layers.push_back(translateLayer(cCache, gctx, kdt, plVolume,
0131                                           indent + m_cfg.hierarchyIndent));
0132         } else {
0133           ACTS_WARNING(indent << "> layer type volume has no internal "
0134                                  "description, layer not built.");
0135         }
0136       }
0137       // Create a new tracking volume with those layers
0138       auto tVolume = m_cfg.trackingVolumeHelper->createTrackingVolume(
0139           gctx, layers, {}, nullptr, rangeR.min(), rangeR.max(), rangeZ.min(),
0140           rangeZ.max(), ptVolume.name);
0141       ACTS_DEBUG(indent << "> translated into bounds: "
0142                         << tVolume->volumeBounds());
0143       return tVolume;
0144     }
0145   }
0146 }
0147
0148 /// @return a new tracking volume
0149 std::shared_ptr<const Layer> 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, double> zrRange;
0159   zrRange[0u] = plVolume.extent.range(AxisDirection::AxisZ);
0160   zrRange[1u] = plVolume.extent.range(AxisDirection::AxisR);
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 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 == 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 == 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     BinningType bType0 = equidistant;
0212     BinningType bType1 = 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 == equidistant && bType1 == 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     ProtoLayer pLayer(gctx, cLayerSurfaces);
0231     pLayer.envelope = plVolume.extent.envelope();
0232     if (its.layerType == 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 == 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 }
0268
0269 }  // namespace Acts