Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-05 08:09:37

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2016-2020 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/LayerCreator.hpp"
0010 
0011 #include "Acts/Definitions/Algebra.hpp"
0012 #include "Acts/Geometry/CylinderLayer.hpp"
0013 #include "Acts/Geometry/DiscLayer.hpp"
0014 #include "Acts/Geometry/Extent.hpp"
0015 #include "Acts/Geometry/Layer.hpp"
0016 #include "Acts/Geometry/PlaneLayer.hpp"
0017 #include "Acts/Geometry/ProtoLayer.hpp"
0018 #include "Acts/Geometry/SurfaceArrayCreator.hpp"
0019 #include "Acts/Surfaces/CylinderBounds.hpp"
0020 #include "Acts/Surfaces/RadialBounds.hpp"
0021 #include "Acts/Surfaces/RectangleBounds.hpp"
0022 #include "Acts/Surfaces/Surface.hpp"
0023 
0024 #include <algorithm>
0025 #include <array>
0026 #include <iterator>
0027 #include <ostream>
0028 #include <set>
0029 #include <utility>
0030 
0031 namespace Acts {
0032 class PlanarBounds;
0033 }  // namespace Acts
0034 
0035 using Acts::VectorHelpers::perp;
0036 using Acts::VectorHelpers::phi;
0037 
0038 Acts::LayerCreator::LayerCreator(const Acts::LayerCreator::Config& lcConfig,
0039                                  std::unique_ptr<const Logger> logger)
0040     : m_cfg(lcConfig), m_logger(std::move(logger)) {}
0041 
0042 void Acts::LayerCreator::setConfiguration(
0043     const Acts::LayerCreator::Config& lcConfig) {
0044   // @todo check consistency
0045   // copy the configuration
0046   m_cfg = lcConfig;
0047 }
0048 
0049 void Acts::LayerCreator::setLogger(std::unique_ptr<const Logger> newLogger) {
0050   m_logger = std::move(newLogger);
0051 }
0052 
0053 Acts::MutableLayerPtr Acts::LayerCreator::cylinderLayer(
0054     const GeometryContext& gctx,
0055     std::vector<std::shared_ptr<const Surface>> surfaces, std::size_t binsPhi,
0056     std::size_t binsZ, std::optional<ProtoLayer> _protoLayer,
0057     const Transform3& transform, std::unique_ptr<ApproachDescriptor> ad) const {
0058   ProtoLayer protoLayer =
0059       _protoLayer ? *_protoLayer : ProtoLayer(gctx, surfaces);
0060 
0061   // Remaining layer parameters - they include the envelopes
0062   double layerR = protoLayer.medium(binR);
0063   double layerZ = protoLayer.medium(binZ);
0064   double layerHalfZ = 0.5 * protoLayer.range(binZ);
0065   double layerThickness = protoLayer.range(binR);
0066 
0067   ACTS_VERBOSE("Creating a cylindrical Layer:");
0068   ACTS_VERBOSE(" - with layer R     = " << layerR);
0069   ACTS_VERBOSE(" - from R min/max   = " << protoLayer.min(binR, false) << " / "
0070                                         << protoLayer.max(binR, false));
0071   ACTS_VERBOSE(" - with R thickness = " << layerThickness);
0072   ACTS_VERBOSE("   - incl envelope  = " << protoLayer.envelope[binR][0u]
0073                                         << " / "
0074                                         << protoLayer.envelope[binR][1u]);
0075 
0076   ACTS_VERBOSE(" - with z min/max   = "
0077                << protoLayer.min(binZ, false) << " (-"
0078                << protoLayer.envelope[binZ][0u] << ") / "
0079                << protoLayer.max(binZ, false) << " (+"
0080                << protoLayer.envelope[binZ][1u] << ")");
0081 
0082   ACTS_VERBOSE(" - z center         = " << layerZ);
0083   ACTS_VERBOSE(" - halflength z     = " << layerHalfZ);
0084 
0085   // create the layer transforms if not given
0086   // we need to transform in case layerZ != 0, so that the layer will be
0087   // correctly defined using the halflength
0088   Translation3 addTranslation(0., 0., 0.);
0089   if (transform.isApprox(Transform3::Identity())) {
0090     addTranslation = Translation3(0., 0., layerZ);
0091     ACTS_VERBOSE(" - layer z shift  = " << -layerZ);
0092   }
0093 
0094   ACTS_VERBOSE(" - with phi min/max = " << protoLayer.min(binPhi, false)
0095                                         << " / "
0096                                         << protoLayer.max(binPhi, false));
0097   ACTS_VERBOSE(" - # of modules     = " << surfaces.size() << " ordered in ( "
0098                                         << binsPhi << " x " << binsZ << ")");
0099   std::unique_ptr<SurfaceArray> sArray;
0100   if (!surfaces.empty()) {
0101     sArray = m_cfg.surfaceArrayCreator->surfaceArrayOnCylinder(
0102         gctx, std::move(surfaces), binsPhi, binsZ, protoLayer);
0103 
0104     checkBinning(gctx, *sArray);
0105   }
0106 
0107   // create the layer and push it back
0108   std::shared_ptr<const CylinderBounds> cBounds(
0109       new CylinderBounds(layerR, layerHalfZ));
0110 
0111   // create the layer
0112   MutableLayerPtr cLayer = CylinderLayer::create(
0113       addTranslation * transform, cBounds, std::move(sArray), layerThickness,
0114       std::move(ad), active);
0115 
0116   if (!cLayer) {
0117     ACTS_ERROR("Creation of cylinder layer did not succeed!");
0118   }
0119   associateSurfacesToLayer(*cLayer);
0120 
0121   // now return
0122   return cLayer;
0123 }
0124 
0125 Acts::MutableLayerPtr Acts::LayerCreator::cylinderLayer(
0126     const GeometryContext& gctx,
0127     std::vector<std::shared_ptr<const Surface>> surfaces, BinningType bTypePhi,
0128     BinningType bTypeZ, std::optional<ProtoLayer> _protoLayer,
0129     const Transform3& transform, std::unique_ptr<ApproachDescriptor> ad) const {
0130   ProtoLayer protoLayer =
0131       _protoLayer ? *_protoLayer : ProtoLayer(gctx, surfaces);
0132 
0133   // remaining layer parameters
0134   double layerR = protoLayer.medium(binR);
0135   double layerZ = protoLayer.medium(binZ);
0136   double layerHalfZ = 0.5 * protoLayer.range(binZ);
0137   double layerThickness = protoLayer.range(binR);
0138 
0139   // adjust the layer radius
0140   ACTS_VERBOSE("Creating a cylindrical Layer:");
0141   ACTS_VERBOSE(" - with layer R     = " << layerR);
0142   ACTS_VERBOSE(" - from R min/max   = " << protoLayer.min(binR, false) << " / "
0143                                         << protoLayer.max(binR, false));
0144   ACTS_VERBOSE(" - with R thickness = " << layerThickness);
0145   ACTS_VERBOSE("   - incl envelope  = " << protoLayer.envelope[binR][0u]
0146                                         << " / "
0147                                         << protoLayer.envelope[binR][1u]);
0148   ACTS_VERBOSE(" - with z min/max   = "
0149                << protoLayer.min(binZ, false) << " (-"
0150                << protoLayer.envelope[binZ][0u] << ") / "
0151                << protoLayer.max(binZ, false) << " (+"
0152                << protoLayer.envelope[binZ][1u] << ")");
0153   ACTS_VERBOSE(" - z center         = " << layerZ);
0154   ACTS_VERBOSE(" - halflength z     = " << layerHalfZ);
0155 
0156   // create the layer transforms if not given
0157   // we need to transform in case layerZ != 0, so that the layer will be
0158   // correctly defined using the halflength
0159   // create the layer transforms if not given
0160   Translation3 addTranslation(0., 0., 0.);
0161   if (transform.isApprox(Transform3::Identity()) && bTypeZ == equidistant) {
0162     addTranslation = Translation3(0., 0., layerZ);
0163     ACTS_VERBOSE(" - layer z shift    = " << -layerZ);
0164   }
0165 
0166   ACTS_VERBOSE(" - with phi min/max = " << protoLayer.min(binPhi, false)
0167                                         << " / "
0168                                         << protoLayer.max(binPhi, false));
0169   ACTS_VERBOSE(" - # of modules     = " << surfaces.size() << "");
0170 
0171   // create the surface array
0172   std::unique_ptr<SurfaceArray> sArray;
0173   if (!surfaces.empty()) {
0174     sArray = m_cfg.surfaceArrayCreator->surfaceArrayOnCylinder(
0175         gctx, std::move(surfaces), bTypePhi, bTypeZ, protoLayer);
0176 
0177     checkBinning(gctx, *sArray);
0178   }
0179 
0180   // create the layer and push it back
0181   std::shared_ptr<const CylinderBounds> cBounds(
0182       new CylinderBounds(layerR, layerHalfZ));
0183 
0184   // create the layer
0185   MutableLayerPtr cLayer = CylinderLayer::create(
0186       addTranslation * transform, cBounds, std::move(sArray), layerThickness,
0187       std::move(ad), active);
0188 
0189   if (!cLayer) {
0190     ACTS_ERROR("Creation of cylinder layer did not succeed!");
0191   }
0192   associateSurfacesToLayer(*cLayer);
0193 
0194   // now return
0195   return cLayer;
0196 }
0197 
0198 Acts::MutableLayerPtr Acts::LayerCreator::discLayer(
0199     const GeometryContext& gctx,
0200     std::vector<std::shared_ptr<const Surface>> surfaces, std::size_t binsR,
0201     std::size_t binsPhi, std::optional<ProtoLayer> _protoLayer,
0202     const Transform3& transform, std::unique_ptr<ApproachDescriptor> ad) const {
0203   ProtoLayer protoLayer =
0204       _protoLayer ? *_protoLayer : ProtoLayer(gctx, surfaces);
0205 
0206   double layerZ = protoLayer.medium(binZ);
0207   double layerThickness = protoLayer.range(binZ);
0208 
0209   // adjust the layer radius
0210   ACTS_VERBOSE("Creating a disk Layer:");
0211   ACTS_VERBOSE(" - at Z position    = " << layerZ);
0212   ACTS_VERBOSE(" - from Z min/max   = " << protoLayer.min(binZ, false) << " / "
0213                                         << protoLayer.max(binZ, false));
0214   ACTS_VERBOSE(" - with Z thickness = " << layerThickness);
0215   ACTS_VERBOSE("   - incl envelope  = " << protoLayer.envelope[binZ][0u]
0216                                         << " / "
0217                                         << protoLayer.envelope[binZ][1u]);
0218   ACTS_VERBOSE(" - with R min/max   = "
0219                << protoLayer.min(binR, false) << " (-"
0220                << protoLayer.envelope[binR][0u] << ") / "
0221                << protoLayer.max(binR, false) << " (+"
0222                << protoLayer.envelope[binR][1u] << ")");
0223   ACTS_VERBOSE(" - with phi min/max = " << protoLayer.min(binPhi, false)
0224                                         << " / "
0225                                         << protoLayer.max(binPhi, false));
0226   ACTS_VERBOSE(" - # of modules    = " << surfaces.size() << " ordered in ( "
0227                                        << binsR << " x " << binsPhi << ")");
0228 
0229   // create the layer transforms if not given
0230   Translation3 addTranslation(0., 0., 0.);
0231   if (transform.isApprox(Transform3::Identity())) {
0232     addTranslation = Translation3(0., 0., layerZ);
0233   }
0234   // create the surface array
0235   std::unique_ptr<SurfaceArray> sArray;
0236   if (!surfaces.empty()) {
0237     sArray = m_cfg.surfaceArrayCreator->surfaceArrayOnDisc(
0238         gctx, std::move(surfaces), binsR, binsPhi, protoLayer, transform);
0239 
0240     checkBinning(gctx, *sArray);
0241   }
0242 
0243   // create the share disc bounds
0244   auto dBounds = std::make_shared<const RadialBounds>(protoLayer.min(binR),
0245                                                       protoLayer.max(binR));
0246 
0247   // create the layers
0248   // we use the same transform here as for the layer itself
0249   // for disk this is fine since we don't bin in Z, so does not matter
0250   MutableLayerPtr dLayer =
0251       DiscLayer::create(addTranslation * transform, dBounds, std::move(sArray),
0252                         layerThickness, std::move(ad), active);
0253 
0254   if (!dLayer) {
0255     ACTS_ERROR("Creation of disc layer did not succeed!");
0256   }
0257   associateSurfacesToLayer(*dLayer);
0258   // return the layer
0259   return dLayer;
0260 }
0261 
0262 Acts::MutableLayerPtr Acts::LayerCreator::discLayer(
0263     const GeometryContext& gctx,
0264     std::vector<std::shared_ptr<const Surface>> surfaces, BinningType bTypeR,
0265     BinningType bTypePhi, std::optional<ProtoLayer> _protoLayer,
0266     const Transform3& transform, std::unique_ptr<ApproachDescriptor> ad) const {
0267   ProtoLayer protoLayer =
0268       _protoLayer ? *_protoLayer : ProtoLayer(gctx, surfaces);
0269 
0270   double layerZ = protoLayer.medium(binZ);
0271   double layerThickness = protoLayer.range(binZ);
0272 
0273   // adjust the layer radius
0274   ACTS_VERBOSE("Creating a disk Layer:");
0275   ACTS_VERBOSE(" - at Z position    = " << layerZ);
0276   ACTS_VERBOSE(" - from Z min/max   = " << protoLayer.min(binZ, false) << " / "
0277                                         << protoLayer.max(binZ, false));
0278   ACTS_VERBOSE(" - with Z thickness = " << layerThickness);
0279   ACTS_VERBOSE("   - incl envelope  = " << protoLayer.envelope[binZ][0u]
0280                                         << " / "
0281                                         << protoLayer.envelope[binZ][1u]);
0282   ACTS_VERBOSE(" - with R min/max   = "
0283                << protoLayer.min(binR, false) << " (-"
0284                << protoLayer.envelope[binR][0u] << ") / "
0285                << protoLayer.max(binR, false) << " (+"
0286                << protoLayer.envelope[binR][1u] << ")");
0287   ACTS_VERBOSE(" - with phi min/max = " << protoLayer.min(binPhi, false)
0288                                         << " / "
0289                                         << protoLayer.max(binPhi, false));
0290   ACTS_VERBOSE(" - # of modules     = " << surfaces.size());
0291 
0292   // create the layer transforms if not given
0293   Translation3 addTranslation(0., 0., 0.);
0294   if (transform.isApprox(Transform3::Identity())) {
0295     addTranslation = Translation3(0., 0., layerZ);
0296   }
0297 
0298   // create the surface array
0299   std::unique_ptr<SurfaceArray> sArray;
0300   if (!surfaces.empty()) {
0301     sArray = m_cfg.surfaceArrayCreator->surfaceArrayOnDisc(
0302         gctx, std::move(surfaces), bTypeR, bTypePhi, protoLayer, transform);
0303 
0304     checkBinning(gctx, *sArray);
0305   }
0306 
0307   // create the shared disc bounds
0308   auto dBounds = std::make_shared<const RadialBounds>(protoLayer.min(binR),
0309                                                       protoLayer.max(binR));
0310 
0311   // create the layers
0312   MutableLayerPtr dLayer =
0313       DiscLayer::create(addTranslation * transform, dBounds, std::move(sArray),
0314                         layerThickness, std::move(ad), active);
0315   if (!dLayer) {
0316     ACTS_ERROR("Creation of disc layer did not succeed!");
0317   }
0318   associateSurfacesToLayer(*dLayer);
0319   // return the layer
0320   return dLayer;
0321 }
0322 
0323 Acts::MutableLayerPtr Acts::LayerCreator::planeLayer(
0324     const GeometryContext& gctx,
0325     std::vector<std::shared_ptr<const Surface>> surfaces, std::size_t bins1,
0326     std::size_t bins2, BinningValue bValue,
0327     std::optional<ProtoLayer> _protoLayer, const Transform3& transform,
0328     std::unique_ptr<ApproachDescriptor> ad) const {
0329   ProtoLayer protoLayer =
0330       _protoLayer ? *_protoLayer : ProtoLayer(gctx, surfaces);
0331 
0332   // remaining layer parameters
0333   double layerHalf1 = 0, layerHalf2 = 0, layerThickness = 0;
0334   switch (bValue) {
0335     case BinningValue::binX: {
0336       layerHalf1 = 0.5 * (protoLayer.max(binY) - protoLayer.min(binY));
0337       layerHalf2 = 0.5 * (protoLayer.max(binZ) - protoLayer.min(binZ));
0338       layerThickness = (protoLayer.max(binX) - protoLayer.min(binX));
0339       break;
0340     }
0341     case BinningValue::binY: {
0342       layerHalf1 = 0.5 * (protoLayer.max(binX) - protoLayer.min(binX));
0343       layerHalf2 = 0.5 * (protoLayer.max(binZ) - protoLayer.min(binZ));
0344       layerThickness = (protoLayer.max(binY) - protoLayer.min(binY));
0345       break;
0346     }
0347     case BinningValue::binZ: {
0348       layerHalf1 = 0.5 * (protoLayer.max(binX) - protoLayer.min(binX));
0349       layerHalf2 = 0.5 * (protoLayer.max(binY) - protoLayer.min(binY));
0350       layerThickness = (protoLayer.max(binZ) - protoLayer.min(binZ));
0351       break;
0352     }
0353     default:
0354       throw std::invalid_argument("Invalid binning value");
0355   }
0356 
0357   double centerX = 0.5 * (protoLayer.max(binX) + protoLayer.min(binX));
0358   double centerY = 0.5 * (protoLayer.max(binY) + protoLayer.min(binY));
0359   double centerZ = 0.5 * (protoLayer.max(binZ) + protoLayer.min(binZ));
0360 
0361   ACTS_VERBOSE("Creating a plane Layer:");
0362   ACTS_VERBOSE(" - with layer center     = "
0363                << "(" << centerX << ", " << centerY << ", " << centerZ << ")");
0364   ACTS_VERBOSE(" - from X min/max   = " << protoLayer.min(binX) << " / "
0365                                         << protoLayer.max(binX));
0366   ACTS_VERBOSE(" - from Y min/max   = " << protoLayer.min(binY) << " / "
0367                                         << protoLayer.max(binY));
0368   ACTS_VERBOSE(" - with Z thickness = " << layerThickness);
0369   ACTS_VERBOSE("   - incl envelope  = " << protoLayer.envelope[bValue][0u]
0370                                         << " / "
0371                                         << protoLayer.envelope[bValue][1u]);
0372 
0373   // create the layer transforms if not given
0374   // we need to transform in case centerX/centerY/centerZ != 0, so that the
0375   // layer will be correctly defined
0376   Translation3 addTranslation(0., 0., 0.);
0377   if (transform.isApprox(Transform3::Identity())) {
0378     addTranslation = Translation3(centerX, centerY, centerZ);
0379     ACTS_VERBOSE(" - layer shift  = "
0380                  << "(" << centerX << ", " << centerY << ", " << centerZ
0381                  << ")");
0382   }
0383 
0384   std::unique_ptr<SurfaceArray> sArray;
0385   if (!surfaces.empty()) {
0386     sArray = m_cfg.surfaceArrayCreator->surfaceArrayOnPlane(
0387         gctx, std::move(surfaces), bins1, bins2, bValue, protoLayer, transform);
0388 
0389     checkBinning(gctx, *sArray);
0390   }
0391 
0392   // create the layer and push it back
0393   std::shared_ptr<const PlanarBounds> pBounds(
0394       new RectangleBounds(layerHalf1, layerHalf2));
0395 
0396   // create the layer
0397   MutableLayerPtr pLayer =
0398       PlaneLayer::create(addTranslation * transform, pBounds, std::move(sArray),
0399                          layerThickness, std::move(ad), active);
0400 
0401   if (!pLayer) {
0402     ACTS_ERROR("Creation of plane layer did not succeed!");
0403   }
0404   associateSurfacesToLayer(*pLayer);
0405 
0406   // now return
0407   return pLayer;
0408 }
0409 
0410 void Acts::LayerCreator::associateSurfacesToLayer(Layer& layer) const {
0411   if (layer.surfaceArray() != nullptr) {
0412     auto surfaces = layer.surfaceArray()->surfaces();
0413 
0414     for (auto& surface : surfaces) {
0415       auto mutableSurface = const_cast<Surface*>(surface);
0416       mutableSurface->associateLayer(layer);
0417     }
0418   }
0419 }
0420 
0421 bool Acts::LayerCreator::checkBinning(const GeometryContext& gctx,
0422                                       const SurfaceArray& sArray) const {
0423   // do consistency check: can we access all sensitive surfaces
0424   // through the binning? If not, surfaces get lost and the binning does not
0425   // work
0426 
0427   ACTS_VERBOSE("Performing consistency check")
0428 
0429   std::vector<const Surface*> surfaces = sArray.surfaces();
0430   std::set<const Surface*> sensitiveSurfaces(surfaces.begin(), surfaces.end());
0431   std::set<const Surface*> accessibleSurfaces;
0432   std::size_t nEmptyBins = 0;
0433   std::size_t nBinsChecked = 0;
0434 
0435   // iterate over all bins
0436   std::size_t size = sArray.size();
0437   for (std::size_t b = 0; b < size; ++b) {
0438     std::vector<const Surface*> binContent = sArray.at(b);
0439     // we don't check under/overflow bins
0440     if (!sArray.isValidBin(b)) {
0441       continue;
0442     }
0443     for (const auto& srf : binContent) {
0444       accessibleSurfaces.insert(srf);
0445     }
0446     if (binContent.empty()) {
0447       nEmptyBins++;
0448     }
0449     nBinsChecked++;
0450   }
0451 
0452   std::vector<const Acts::Surface*> diff;
0453   std::set_difference(sensitiveSurfaces.begin(), sensitiveSurfaces.end(),
0454                       accessibleSurfaces.begin(), accessibleSurfaces.end(),
0455                       std::inserter(diff, diff.begin()));
0456 
0457   ACTS_VERBOSE(" - Checked " << nBinsChecked << " valid bins");
0458 
0459   if (nEmptyBins > 0) {
0460     ACTS_ERROR(" -- Not all bins point to surface. " << nEmptyBins << " empty");
0461   } else {
0462     ACTS_VERBOSE(" -- All bins point to a surface");
0463   }
0464 
0465   if (!diff.empty()) {
0466     ACTS_ERROR(
0467         " -- Not all sensitive surfaces are accessible through binning. "
0468         "sensitive: "
0469         << sensitiveSurfaces.size()
0470         << "    accessible: " << accessibleSurfaces.size());
0471 
0472     // print all inaccessibles
0473     ACTS_ERROR(" -- Inaccessible surfaces: ");
0474     for (const auto& srf : diff) {
0475       // have to choose BinningValue here
0476       Vector3 ctr = srf->binningPosition(gctx, binR);
0477       ACTS_ERROR(" Surface(x=" << ctr.x() << ", y=" << ctr.y()
0478                                << ", z=" << ctr.z() << ", r=" << perp(ctr)
0479                                << ", phi=" << phi(ctr) << ")");
0480     }
0481 
0482   } else {
0483     ACTS_VERBOSE(" -- All sensitive surfaces are accessible through binning.");
0484   }
0485 
0486   return nEmptyBins == 0 && diff.empty();
0487 }