Back to home page

EIC code displayed by LXR
 
 

    

Warning, file /acts/Core/src/Geometry/SurfaceArrayCreator.cpp was not indexed or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).

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/SurfaceArrayCreator.hpp"
0010 
0011 #include "Acts/Definitions/Algebra.hpp"
0012 #include "Acts/Surfaces/CylinderSurface.hpp"
0013 #include "Acts/Surfaces/DiscSurface.hpp"
0014 #include "Acts/Surfaces/PlanarBounds.hpp"
0015 #include "Acts/Surfaces/PlaneSurface.hpp"
0016 #include "Acts/Surfaces/RectangleBounds.hpp"
0017 #include "Acts/Surfaces/Surface.hpp"
0018 #include "Acts/Surfaces/SurfaceArray.hpp"
0019 #include "Acts/Utilities/AxisDefinitions.hpp"
0020 #include "Acts/Utilities/BinningType.hpp"
0021 #include "Acts/Utilities/Helpers.hpp"
0022 
0023 #include <algorithm>
0024 #include <numbers>
0025 #include <stdexcept>
0026 
0027 namespace Acts {
0028 
0029 using VectorHelpers::perp;
0030 using VectorHelpers::phi;
0031 
0032 SurfaceArray SurfaceArrayCreator::surfaceArrayOnCylinder(
0033     const GeometryContext& gctx,
0034     std::vector<std::shared_ptr<const Surface>> surfaces, std::size_t binsPhi,
0035     std::size_t binsZ, std::optional<ProtoLayer> protoLayerOpt,
0036     const Transform3& transform, std::uint8_t maxNeighborDistance) const {
0037   using enum AxisDirection;
0038 
0039   const std::vector<const Surface*> surfacesRaw = unpackSmartPointers(surfaces);
0040   // Check if we have proto layer, else build it
0041   const ProtoLayer protoLayer =
0042       protoLayerOpt ? *protoLayerOpt : ProtoLayer(gctx, surfacesRaw);
0043 
0044   ACTS_VERBOSE("Creating a SurfaceArray on a cylinder");
0045   ACTS_VERBOSE(" -- with " << surfaces.size() << " surfaces.");
0046   ACTS_VERBOSE(" -- with phi x z  = " << binsPhi << " x " << binsZ << " = "
0047                                       << binsPhi * binsZ << " bins.");
0048 
0049   Transform3 fullTransform = transform;
0050   const auto pAxisPhi =
0051       createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Closed,
0052                             AxisPhi, protoLayer, fullTransform, binsPhi);
0053   const auto pAxisZ =
0054       createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Bound, AxisZ,
0055                             protoLayer, fullTransform, binsZ);
0056 
0057   const double R = protoLayer.medium(AxisR, true);
0058   const double halfZ = protoLayer.range(AxisZ, true) * 0.5;
0059   const double layerTolerance = protoLayer.range(AxisR) * 0.5;
0060 
0061   auto surface = Surface::makeShared<CylinderSurface>(fullTransform, R, halfZ);
0062   ACTS_VERBOSE("- projection surface is: " << surface->toString(gctx));
0063 
0064   return SurfaceArray(gctx, std::move(surfaces), std::move(surface),
0065                       layerTolerance, {*pAxisPhi, *pAxisZ},
0066                       maxNeighborDistance);
0067 }
0068 
0069 SurfaceArray SurfaceArrayCreator::surfaceArrayOnCylinder(
0070     const GeometryContext& gctx,
0071     std::vector<std::shared_ptr<const Surface>> surfaces, BinningType bTypePhi,
0072     BinningType bTypeZ, std::optional<ProtoLayer> protoLayerOpt,
0073     const Transform3& transform, std::uint8_t maxNeighborDistance) const {
0074   using enum AxisDirection;
0075 
0076   const std::vector<const Surface*> surfacesRaw = unpackSmartPointers(surfaces);
0077   // check if we have proto layer, else build it
0078   const ProtoLayer protoLayer =
0079       protoLayerOpt ? *protoLayerOpt : ProtoLayer(gctx, surfacesRaw);
0080 
0081   const double R = protoLayer.medium(AxisR, true);
0082   const double halfZ = protoLayer.range(AxisZ, true) * 0.5;
0083   const double layerTolerance = protoLayer.range(AxisR) * 0.5;
0084 
0085   std::unique_ptr<const IAxis> pAxisPhi;
0086   std::unique_ptr<const IAxis> pAxisZ;
0087 
0088   Transform3 fullTransform = transform;
0089 
0090   if (bTypePhi == equidistant) {
0091     pAxisPhi =
0092         createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Closed,
0093                               AxisPhi, protoLayer, fullTransform, 0);
0094   } else {
0095     pAxisPhi = createVariableAxis(gctx, surfacesRaw, AxisBoundaryType::Closed,
0096                                   AxisPhi, protoLayer, fullTransform);
0097   }
0098 
0099   if (bTypeZ == equidistant) {
0100     pAxisZ = createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Bound,
0101                                    AxisZ, protoLayer, fullTransform);
0102   } else {
0103     pAxisZ = createVariableAxis(gctx, surfacesRaw, AxisBoundaryType::Bound,
0104                                 AxisZ, protoLayer, fullTransform);
0105   }
0106 
0107   auto surface = Surface::makeShared<CylinderSurface>(fullTransform, R, halfZ);
0108 
0109   const std::size_t bins0 = pAxisPhi->getNBins();
0110   const std::size_t bins1 = pAxisZ->getNBins();
0111   ACTS_VERBOSE("Creating a SurfaceArray on a cylinder");
0112   ACTS_VERBOSE(" -- with " << surfaces.size() << " surfaces.");
0113   ACTS_VERBOSE(" -- with phi x z  = " << bins0 << " x " << bins1 << " = "
0114                                       << bins0 * bins1 << " bins.");
0115 
0116   return SurfaceArray(gctx, std::move(surfaces), std::move(surface),
0117                       layerTolerance, {*pAxisPhi, *pAxisZ},
0118                       maxNeighborDistance);
0119 }
0120 
0121 SurfaceArray SurfaceArrayCreator::surfaceArrayOnDisc(
0122     const GeometryContext& gctx,
0123     std::vector<std::shared_ptr<const Surface>> surfaces, std::size_t binsR,
0124     std::size_t binsPhi, std::optional<ProtoLayer> protoLayerOpt,
0125     const Transform3& transform, std::uint8_t maxNeighborDistance) const {
0126   using enum AxisDirection;
0127 
0128   const std::vector<const Surface*> surfacesRaw = unpackSmartPointers(surfaces);
0129   // check if we have proto layer, else build it
0130   const ProtoLayer protoLayer =
0131       protoLayerOpt ? *protoLayerOpt : ProtoLayer(gctx, surfacesRaw);
0132 
0133   ACTS_VERBOSE("Creating a SurfaceArray on a disc");
0134 
0135   Transform3 fullTransform = transform;
0136   const auto pAxisR =
0137       createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Bound, AxisR,
0138                             protoLayer, fullTransform, binsR);
0139   const auto pAxisPhi =
0140       createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Closed,
0141                             AxisPhi, protoLayer, fullTransform, binsPhi);
0142 
0143   const double Z = protoLayer.medium(AxisZ, true);
0144   const double Rmin = protoLayer.min(AxisR, true);
0145   const double Rmax = protoLayer.max(AxisR, true);
0146   const double layerThickness = protoLayer.range(AxisZ) * 0.5;
0147   ACTS_VERBOSE("- z-position of disc estimated as " << Z);
0148   ACTS_VERBOSE("- full transform is \n" << fullTransform.matrix());
0149 
0150   if (fullTransform.translation().norm() < s_transformEquivalentTolerance) {
0151     ACTS_VERBOSE(
0152         "input transform does not have translation: putting projection surface "
0153         "at center of gravity in z");
0154     fullTransform.translate(Vector3::UnitZ() * Z);
0155   }
0156 
0157   auto surface = Surface::makeShared<DiscSurface>(fullTransform, Rmin, Rmax);
0158   ACTS_VERBOSE("- projection surface is: " << surface->toString(gctx));
0159 
0160   const std::size_t bins0 = pAxisR->getNBins();
0161   const std::size_t bins1 = pAxisPhi->getNBins();
0162   ACTS_VERBOSE(" -- with " << surfaces.size() << " surfaces.");
0163   ACTS_VERBOSE(" -- with r x phi  = " << bins0 << " x " << bins1 << " = "
0164                                       << bins0 * bins1 << " bins.");
0165 
0166   return SurfaceArray(gctx, std::move(surfaces), std::move(surface),
0167                       layerThickness, {*pAxisR, *pAxisPhi},
0168                       maxNeighborDistance);
0169 }
0170 
0171 SurfaceArray SurfaceArrayCreator::surfaceArrayOnDisc(
0172     const GeometryContext& gctx,
0173     std::vector<std::shared_ptr<const Surface>> surfaces, BinningType bTypeR,
0174     BinningType bTypePhi, std::optional<ProtoLayer> protoLayerOpt,
0175     const Transform3& transform, std::uint8_t maxNeighborDistance) const {
0176   using enum AxisDirection;
0177 
0178   const std::vector<const Surface*> surfacesRaw = unpackSmartPointers(surfaces);
0179   // check if we have proto layer, else build it
0180   const ProtoLayer protoLayer =
0181       protoLayerOpt ? *protoLayerOpt : ProtoLayer(gctx, surfacesRaw);
0182 
0183   ACTS_VERBOSE("Creating a SurfaceArray on a disc");
0184 
0185   std::unique_ptr<const IAxis> pAxisPhi;
0186   std::unique_ptr<const IAxis> pAxisR;
0187 
0188   Transform3 fullTransform = transform;
0189   Transform3 inverseTransform = transform.inverse();
0190 
0191   if (bTypeR == equidistant) {
0192     pAxisR = createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Bound,
0193                                    AxisR, protoLayer, fullTransform);
0194   } else {
0195     pAxisR = createVariableAxis(gctx, surfacesRaw, AxisBoundaryType::Bound,
0196                                 AxisR, protoLayer, fullTransform);
0197   }
0198 
0199   // if we have more than one R ring, we need to figure out
0200   // the number of phi bins.
0201   if (pAxisR->getNBins() > 1) {
0202     // more than one R-Ring, we need to adjust
0203     // this FORCES equidistant binning
0204     std::vector<std::vector<const Surface*>> phiModules(pAxisR->getNBins());
0205     for (const auto& srf : surfacesRaw) {
0206       const Vector3 bpos =
0207           inverseTransform * srf->referencePosition(gctx, AxisR);
0208       const std::size_t bin =
0209           pAxisR->getBin(perp(bpos)) - 1;  // subtract underflow bin
0210       phiModules.at(bin).push_back(srf);
0211     }
0212 
0213     std::vector<std::size_t> nPhiModules;
0214     const auto& matcher = m_cfg.surfaceMatcher;
0215     const auto equal = [&gctx, &matcher](const Surface* a, const Surface* b) {
0216       return matcher(gctx, AxisPhi, a, b);
0217     };
0218 
0219     std::transform(
0220         phiModules.begin(), phiModules.end(), std::back_inserter(nPhiModules),
0221         [&equal,
0222          this](const std::vector<const Surface*>& surfaces_) -> std::size_t {
0223           return this->findKeySurfaces(surfaces_, equal).size();
0224         });
0225 
0226     // @FIXME: Problem: phi binning runs rotation to optimize
0227     // for bin edges. This FAILS after this modification, since
0228     // the bin count is the one from the lowest module-count bin,
0229     // but the rotation is done considering all bins.
0230     // This might be resolved through bin completion, but not sure.
0231     // @TODO: check in extrapolation
0232     const std::size_t nBinsPhi =
0233         *std::min_element(nPhiModules.begin(), nPhiModules.end());
0234     pAxisPhi =
0235         createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Closed,
0236                               AxisPhi, protoLayer, fullTransform, nBinsPhi);
0237 
0238   } else {
0239     // use regular determination
0240     if (bTypePhi == equidistant) {
0241       pAxisPhi =
0242           createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Closed,
0243                                 AxisPhi, protoLayer, fullTransform, 0);
0244     } else {
0245       pAxisPhi = createVariableAxis(gctx, surfacesRaw, AxisBoundaryType::Closed,
0246                                     AxisPhi, protoLayer, fullTransform);
0247     }
0248   }
0249 
0250   const double Z = protoLayer.medium(AxisZ, true);
0251   const double Rmin = protoLayer.min(AxisR, true);
0252   const double Rmax = protoLayer.max(AxisR, true);
0253   const double layerThickness = protoLayer.range(AxisZ) * 0.5;
0254   ACTS_VERBOSE("- z-position of disc estimated as " << Z);
0255 
0256   if (fullTransform.translation().norm() < s_transformEquivalentTolerance) {
0257     ACTS_VERBOSE(
0258         "input transform does not have translation: putting projection surface "
0259         "at center of gravity in z");
0260     fullTransform.translate(Vector3::UnitZ() * Z);
0261   }
0262 
0263   auto surface = Surface::makeShared<DiscSurface>(fullTransform, Rmin, Rmax);
0264 
0265   // get the number of bins
0266   const std::size_t bins0 = pAxisR->getNBins();
0267   const std::size_t bins1 = pAxisPhi->getNBins();
0268   ACTS_VERBOSE(" -- with " << surfaces.size() << " surfaces.");
0269   ACTS_VERBOSE(" -- with r x phi  = " << bins0 << " x " << bins1 << " = "
0270                                       << bins0 * bins1 << " bins.");
0271 
0272   return SurfaceArray(gctx, std::move(surfaces), std::move(surface),
0273                       layerThickness, {*pAxisR, *pAxisPhi},
0274                       maxNeighborDistance);
0275 }
0276 
0277 /// SurfaceArrayCreator interface method - create an array on a plane
0278 SurfaceArray SurfaceArrayCreator::surfaceArrayOnPlane(
0279     const GeometryContext& gctx,
0280     std::vector<std::shared_ptr<const Surface>> surfaces, std::size_t bins1,
0281     std::size_t bins2, AxisDirection aDir,
0282     std::optional<ProtoLayer> protoLayerOpt, const Transform3& transform,
0283     std::uint8_t maxNeighborDistance) const {
0284   using enum AxisDirection;
0285 
0286   const std::vector<const Surface*> surfacesRaw = unpackSmartPointers(surfaces);
0287   // check if we have proto layer, else build it
0288   const ProtoLayer protoLayer =
0289       protoLayerOpt ? *protoLayerOpt : ProtoLayer(gctx, surfacesRaw);
0290 
0291   ACTS_VERBOSE("Creating a SurfaceArray on a plance");
0292   ACTS_VERBOSE(" -- with " << surfaces.size() << " surfaces.");
0293   ACTS_VERBOSE(" -- with " << bins1 << " x " << bins2 << " = " << bins1 * bins2
0294                            << " bins.");
0295   Transform3 fullTransform = transform;
0296 
0297   const double layerTolerance = protoLayer.range(aDir) * 0.5;
0298 
0299   // Axis along the binning
0300   switch (aDir) {
0301     case AxisX: {
0302       const auto pAxis1 =
0303           createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Bound,
0304                                 AxisY, protoLayer, fullTransform, bins1);
0305       const auto pAxis2 =
0306           createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Bound,
0307                                 AxisZ, protoLayer, fullTransform, bins2);
0308       auto surface = Surface::makeShared<PlaneSurface>(
0309           fullTransform,
0310           std::make_shared<RectangleBounds>(
0311               Vector2(protoLayer.min(AxisY), protoLayer.min(AxisZ)),
0312               Vector2(protoLayer.max(AxisY), protoLayer.max(AxisZ))));
0313       return SurfaceArray(gctx, std::move(surfaces), std::move(surface),
0314                           layerTolerance, {*pAxis1, *pAxis2},
0315                           maxNeighborDistance);
0316     }
0317     case AxisY: {
0318       const auto pAxis1 =
0319           createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Bound,
0320                                 AxisX, protoLayer, fullTransform, bins1);
0321       const auto pAxis2 =
0322           createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Bound,
0323                                 AxisZ, protoLayer, fullTransform, bins2);
0324       auto surface = Surface::makeShared<PlaneSurface>(
0325           fullTransform,
0326           std::make_shared<RectangleBounds>(
0327               Vector2(protoLayer.min(AxisX), protoLayer.min(AxisY)),
0328               Vector2(protoLayer.max(AxisX), protoLayer.max(AxisY))));
0329       return SurfaceArray(gctx, std::move(surfaces), std::move(surface),
0330                           layerTolerance, {*pAxis1, *pAxis2},
0331                           maxNeighborDistance);
0332     }
0333     case AxisZ: {
0334       const auto pAxis1 =
0335           createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Bound,
0336                                 AxisX, protoLayer, fullTransform, bins1);
0337       const auto pAxis2 =
0338           createEquidistantAxis(gctx, surfacesRaw, AxisBoundaryType::Bound,
0339                                 AxisY, protoLayer, fullTransform, bins2);
0340       auto surface = Surface::makeShared<PlaneSurface>(
0341           fullTransform,
0342           std::make_shared<RectangleBounds>(
0343               Vector2(protoLayer.min(AxisX), protoLayer.min(AxisY)),
0344               Vector2(protoLayer.max(AxisX), protoLayer.max(AxisY))));
0345       return SurfaceArray(gctx, std::move(surfaces), std::move(surface),
0346                           layerTolerance, {*pAxis1, *pAxis2},
0347                           maxNeighborDistance);
0348     }
0349     default:
0350       break;
0351   }
0352 
0353   throw std::invalid_argument(
0354       "SurfaceArrayCreator::surfaceArrayOnPlane: Invalid binning direction");
0355 }
0356 
0357 std::vector<const Surface*> SurfaceArrayCreator::findKeySurfaces(
0358     const std::vector<const Surface*>& surfaces,
0359     const std::function<bool(const Surface*, const Surface*)>& equal) const {
0360   std::vector<const Surface*> keys;
0361   for (const auto& srfA : surfaces) {
0362     bool exists = false;
0363     for (const auto& srfB : keys) {
0364       if (equal(srfA, srfB)) {
0365         exists = true;
0366         break;
0367       }
0368     }
0369     if (!exists) {
0370       keys.push_back(srfA);
0371     }
0372   }
0373 
0374   return keys;
0375 }
0376 
0377 std::size_t SurfaceArrayCreator::determineBinCount(
0378     const GeometryContext& gctx, const std::vector<const Surface*>& surfaces,
0379     AxisDirection aDir) const {
0380   const auto& matcher = m_cfg.surfaceMatcher;
0381   const auto equal = [&gctx, &aDir, &matcher](const Surface* a,
0382                                               const Surface* b) {
0383     return matcher(gctx, aDir, a, b);
0384   };
0385   const std::vector<const Surface*> keys = findKeySurfaces(surfaces, equal);
0386 
0387   return keys.size();
0388 }
0389 
0390 std::unique_ptr<const IAxis> SurfaceArrayCreator::createVariableAxis(
0391     const GeometryContext& gctx, const std::vector<const Surface*>& surfaces,
0392     AxisBoundaryType aBoundaryType, AxisDirection aDir,
0393     const ProtoLayer& protoLayer, Transform3& transform) const {
0394   using enum AxisDirection;
0395 
0396   if (surfaces.empty()) {
0397     throw std::logic_error(
0398         "No surfaces handed over for creating arbitrary bin utility!");
0399   }
0400   // BinningOption is open for z and r, in case of phi binning reset later
0401   // the vector with the binning Values (boundaries for each bin)
0402 
0403   // bind matcher with binning type
0404   const auto& matcher = m_cfg.surfaceMatcher;
0405   // find the key surfaces
0406   const auto equal = [&gctx, &aDir, &matcher](const Surface* a,
0407                                               const Surface* b) {
0408     return matcher(gctx, aDir, a, b);
0409   };
0410   std::vector<const Surface*> keys = findKeySurfaces(surfaces, equal);
0411 
0412   std::vector<double> binEdges;
0413   if (aDir == AxisPhi) {
0414     std::stable_sort(keys.begin(), keys.end(),
0415                      [&gctx](const Surface* a, const Surface* b) {
0416                        return (phi(a->referencePosition(gctx, AxisPhi)) <
0417                                phi(b->referencePosition(gctx, AxisPhi)));
0418                      });
0419 
0420     const double maxPhi =
0421         0.5 * (phi(keys.at(0)->referencePosition(gctx, AxisPhi)) +
0422                phi(keys.at(1)->referencePosition(gctx, AxisPhi)));
0423 
0424     // create rotation, so that maxPhi is +pi
0425     const double angle = -(std::numbers::pi + maxPhi);
0426     transform = transform * AngleAxis3(angle, Vector3::UnitZ());
0427 
0428     // iterate over all key surfaces, and use their mean position as aDirs,
0429     // but
0430     // rotate using transform from before
0431     double previous = phi(keys.at(0)->referencePosition(gctx, AxisPhi));
0432     // go through key surfaces
0433     for (std::size_t i = 1; i < keys.size(); i++) {
0434       const Surface* surface = keys.at(i);
0435       // create central binning values which is the mean of the center
0436       // positions in the binning direction of the current and previous
0437       // surface
0438       const double edge =
0439           0.5 * (previous + phi(surface->referencePosition(gctx, AxisPhi))) +
0440           angle;
0441       binEdges.push_back(edge);
0442       previous = phi(surface->referencePosition(gctx, AxisPhi));
0443     }
0444 
0445     // segments
0446     constexpr unsigned int segments = 72;
0447 
0448     // get the bounds of the last surfaces
0449     const Surface* backSurface = keys.back();
0450     const PlanarBounds* backBounds =
0451         dynamic_cast<const PlanarBounds*>(&(backSurface->bounds()));
0452     if (backBounds == nullptr) {
0453       ACTS_ERROR(
0454           "Given SurfaceBounds are not planar - not implemented for "
0455           "other bounds yet! ");
0456     }
0457     // get the global vertices
0458     const std::vector<Vector3> backVertices =
0459         makeGlobalVertices(gctx, *backSurface, backBounds->vertices(segments));
0460     const double maxBValue = phi(*std::ranges::max_element(
0461         backVertices, {}, [](const Vector3& v) { return phi(v); }));
0462 
0463     binEdges.push_back(maxBValue);
0464 
0465     binEdges.push_back(std::numbers::pi);
0466 
0467   } else if (aDir == AxisZ) {
0468     std::ranges::stable_sort(keys, {}, [&gctx](const Surface* s) {
0469       return s->referencePosition(gctx, AxisZ).z();
0470     });
0471 
0472     binEdges.push_back(protoLayer.min(AxisZ));
0473     binEdges.push_back(protoLayer.max(AxisZ));
0474 
0475     // the z-center position of the previous surface
0476     double previous = keys.front()->referencePosition(gctx, AxisZ).z();
0477     // go through key surfaces
0478     for (auto surface = keys.begin() + 1; surface != keys.end(); surface++) {
0479       // create central binning values which is the mean of the center
0480       // positions in the binning direction of the current and previous
0481       // surface
0482       binEdges.push_back(
0483           0.5 * (previous + (*surface)->referencePosition(gctx, AxisZ).z()));
0484       previous = (*surface)->referencePosition(gctx, AxisZ).z();
0485     }
0486   } else {  // AxisR
0487     std::ranges::stable_sort(keys, {}, [&gctx](const Surface* s) {
0488       return perp(s->referencePosition(gctx, AxisR));
0489     });
0490 
0491     binEdges.push_back(protoLayer.min(AxisR));
0492     binEdges.push_back(protoLayer.max(AxisR));
0493 
0494     // the r-center position of the previous surface
0495     double previous = perp(keys.front()->referencePosition(gctx, AxisR));
0496 
0497     // go through key surfaces
0498     for (auto surface = keys.begin() + 1; surface != keys.end(); surface++) {
0499       // create central binning values which is the mean of the center
0500       // positions in the binning direction of the current and previous
0501       // surface
0502       binEdges.push_back(
0503           0.5 * (previous + perp((*surface)->referencePosition(gctx, AxisR))));
0504       previous = perp((*surface)->referencePosition(gctx, AxisR));
0505     }
0506   }
0507   std::ranges::sort(binEdges);
0508   ACTS_VERBOSE("Create variable binning Axis for binned SurfaceArray");
0509   ACTS_VERBOSE("    AxisDirection: " << aDir);
0510   ACTS_VERBOSE("    Number of bins: " << (binEdges.size() - 1));
0511   ACTS_VERBOSE("    (Min/Max) = (" << binEdges.front() << "/"
0512                                        << binEdges.back() << ")");
0513 
0514   return IAxis::createVariable(aBoundaryType, binEdges, aDir);
0515 }
0516 
0517 std::unique_ptr<const IAxis> SurfaceArrayCreator::createEquidistantAxis(
0518     const GeometryContext& gctx, const std::vector<const Surface*>& surfaces,
0519     AxisBoundaryType aBoundaryType, AxisDirection aDir,
0520     const ProtoLayer& protoLayer, Transform3& transform,
0521     std::size_t nBins) const {
0522   using enum AxisDirection;
0523 
0524   if (surfaces.empty()) {
0525     throw std::logic_error(
0526         "No surfaces handed over for creating equidistant axis!");
0527   }
0528   // check the binning type first
0529 
0530   double minimum = protoLayer.min(aDir, false);
0531   double maximum = protoLayer.max(aDir, false);
0532 
0533   std::size_t binNumber = 0;
0534   if (nBins == 0) {
0535     // determine bin count
0536     binNumber = determineBinCount(gctx, surfaces, aDir);
0537   } else {
0538     // use bin count
0539     binNumber = nBins;
0540   }
0541 
0542   // bind matcher & context with binning type
0543   auto matcher = m_cfg.surfaceMatcher;
0544 
0545   // now check the binning value
0546   if (aDir == AxisPhi) {
0547     minimum = protoLayer.min(AxisPhi, true);
0548     maximum = protoLayer.max(AxisPhi, true);
0549 
0550     if (m_cfg.doPhiBinningOptimization) {
0551       minimum = -std::numbers::pi;
0552       maximum = std::numbers::pi;
0553 
0554       // Phi binning
0555       // set the binning option for phi
0556       // sort first in phi
0557       const Surface* maxElem =
0558           *std::ranges::max_element(surfaces, {}, [&gctx](const Surface* s) {
0559             return phi(s->referencePosition(gctx, AxisR));
0560           });
0561 
0562       // rotate to max phi module plus one half step
0563       // this should make sure that phi wrapping at +- pi
0564       // never falls on a module center
0565       const double surfaceMax = phi(maxElem->referencePosition(gctx, AxisR));
0566       const double gridStep =
0567           2 * std::numbers::pi / static_cast<double>(binNumber);
0568       const double gridMax = std::numbers::pi - 0.5 * gridStep;
0569       const double angle = gridMax - surfaceMax;
0570 
0571       // replace given transform ref
0572       transform = transform * AngleAxis3(angle, Vector3::UnitZ());
0573     }
0574   }
0575 
0576   // assign the bin size
0577   ACTS_VERBOSE("Create equidistant binning Axis for binned SurfaceArray");
0578   ACTS_VERBOSE("    AxisDirection: " << aDir);
0579   ACTS_VERBOSE("    Number of bins: " << binNumber);
0580   ACTS_VERBOSE("    (Min/Max) = (" << minimum << "/" << maximum << ")");
0581 
0582   return IAxis::createEquidistant(aBoundaryType, minimum, maximum, binNumber,
0583                                   aDir);
0584 }
0585 
0586 std::vector<Vector3> SurfaceArrayCreator::makeGlobalVertices(
0587     const GeometryContext& gctx, const Surface& surface,
0588     const std::vector<Vector2>& locVertices) const {
0589   std::vector<Vector3> globVertices;
0590   for (auto& vertex : locVertices) {
0591     Vector3 globVertex = surface.localToGlobal(gctx, vertex, Vector3());
0592     globVertices.push_back(globVertex);
0593   }
0594   return globVertices;
0595 }
0596 
0597 }  // namespace Acts
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!