src/Geometry/GridPortalLinkMerging.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/Definitions/Tolerance.hpp"
0010 #include "Acts/Geometry/CompositePortalLink.hpp"
0011 #include "Acts/Geometry/GridPortalLink.hpp"
0012 #include "Acts/Surfaces/CylinderSurface.hpp"
0013 #include "Acts/Surfaces/DiscSurface.hpp"
0014 #include "Acts/Surfaces/PlaneSurface.hpp"
0015 #include "Acts/Surfaces/RadialBounds.hpp"
0016
0017 #include <algorithm>
0018 #include <cassert>
0019 #include <limits>
0020 #include <memory>
0021 #include <stdexcept>
0022
0023 namespace Acts {
0024
0025 namespace {
0026
0027 template <typename... Args>
0028 std::unique_ptr<GridPortalLink> makeGrid(
0029     const std::shared_ptr<RegularSurface>& surface, AxisDirection direction,
0030     const Logger& logger, std::tuple<Args...> args, const IAxis* otherAxis,
0031     bool prepend) {
0032   // @TODO: PlaneSurface support
0033
0034   ACTS_VERBOSE("Make resulting merged grid");
0035
0036   // This is to make it possible to construct Axis with the tuple arguments
0037   auto axisFactory = []<typename... Ts>(Ts&&... axisArgs) {
0038     return Axis{std::forward<Ts>(axisArgs)...};
0039   };
0040
0041   // Avoid copy-pasting identical code twice below
0042   auto makeGrid = [&](auto boundaryType) -> std::unique_ptr<GridPortalLink> {
0043     auto axisArgs = std::tuple_cat(std::tuple{boundaryType}, std::move(args));
0044     auto merged = std::apply(axisFactory, std::move(axisArgs));
0045
0046     if (otherAxis == nullptr) {
0047       // No other axis
0048       ACTS_VERBOSE("    ~> merged axis: " << merged);
0049       return GridPortalLink::make(surface, direction, std::move(merged));
0050     } else {
0051       return otherAxis->visit(
0052           [&](const auto& axis) -> std::unique_ptr<GridPortalLink> {
0053             if (prepend) {
0054               ACTS_VERBOSE("    ~> other axis (prepend): " << axis);
0055               ACTS_VERBOSE("    ~> merged axis: " << merged);
0056               return GridPortalLink::make(surface, axis, std::move(merged));
0057             } else {
0058               ACTS_VERBOSE("    ~> merged axis: " << merged);
0059               ACTS_VERBOSE("    ~> other axis (append): " << axis);
0060               return GridPortalLink::make(surface, std::move(merged), axis);
0061             }
0062           });
0063     }
0064   };
0065
0066   // Check if we're in the cylinder or disc case, and the resulting bounds wrap
0067   // around and should have closed binning
0068   if (direction == AxisDirection::AxisPhi ||
0069       direction == AxisDirection::AxisRPhi) {
0070     if (const auto* cylinder =
0071             dynamic_cast<const CylinderSurface*>(surface.get());
0072         cylinder != nullptr) {
0073       if (cylinder->bounds().coversFullAzimuth()) {
0074         return makeGrid(AxisClosed);
0075       }
0076     } else if (const auto* disc =
0077                    dynamic_cast<const DiscSurface*>(surface.get());
0078                disc != nullptr) {
0079       const auto& radialBounds =
0080           dynamic_cast<const RadialBounds&>(disc->bounds());
0081       if (radialBounds.coversFullAzimuth()) {
0082         return makeGrid(AxisClosed);
0083       }
0084     }
0085   }
0086
0087   return makeGrid(AxisBound);
0088 }
0089
0090 std::unique_ptr<GridPortalLink> mergeVariable(
0091     const std::shared_ptr<RegularSurface>& mergedSurface, const IAxis& axisA,
0092     const IAxis& axisB, double /*tolerance*/, AxisDirection direction,
0093     const Logger& logger, const IAxis* otherAxis, bool prepend) {
0094   ACTS_VERBOSE("Variable merge: direction is " << direction);
0095
0096   ACTS_VERBOSE("~> axis a: " << axisA);
0097   ACTS_VERBOSE("~> axis b: " << axisB);
0098
0099   std::vector<double> binEdges;
0100
0101   binEdges.reserve(axisA.getNBins() + axisB.getNBins() + 1);
0102
0103   auto edgesA = axisA.getBinEdges();
0104
0105   if (direction == AxisDirection::AxisR) {
0106     ACTS_VERBOSE("Performing asymmetric merge");
0107     std::ranges::copy(edgesA, std::back_inserter(binEdges));
0108
0109   } else {
0110     ACTS_VERBOSE("Performing symmetrized merge");
0111     double halfWidth =
0112         (axisA.getMax() - axisA.getMin() + axisB.getMax() - axisB.getMin()) /
0113         2.0;
0114     ACTS_VERBOSE("    ~> half width: " << halfWidth);
0115
0116     double shift = axisA.getMax() - halfWidth;
0117     ACTS_VERBOSE("    ~> shift: " << shift);
0118
0119     std::ranges::transform(edgesA, std::back_inserter(binEdges),
0120                            [&](double edge) { return edge + shift; });
0121   }
0122
0123   double stitchPoint = binEdges.back();
0124   auto edgesB = axisB.getBinEdges();
0125   std::transform(
0126       std::next(edgesB.begin()), edgesB.end(), std::back_inserter(binEdges),
0127       [&](double edge) { return edge - axisB.getMin() + stitchPoint; });
0128
0129   return makeGrid(mergedSurface, direction, logger,
0130                   std::tuple{std::move(binEdges)}, otherAxis, prepend);
0131 }
0132
0133 std::unique_ptr<GridPortalLink> mergeEquidistant(
0134     const std::shared_ptr<RegularSurface>& mergedSurface, const IAxis& axisA,
0135     const IAxis& axisB, double tolerance, AxisDirection direction,
0136     const Logger& logger, const IAxis* otherAxis, bool prepend) {
0137   ACTS_VERBOSE("===> potentially equidistant merge: checking bin widths");
0138
0139   ACTS_VERBOSE("~> axis a: " << axisA);
0140   ACTS_VERBOSE("~> axis b: " << axisB);
0141
0142   double binsWidthA =
0143       (axisA.getMax() - axisA.getMin()) / static_cast<double>(axisA.getNBins());
0144   double binsWidthB =
0145       (axisB.getMax() - axisB.getMin()) / static_cast<double>(axisB.getNBins());
0146
0147   ACTS_VERBOSE("  ~> binWidths: " << binsWidthA << " vs " << binsWidthB);
0148
0149   if (std::abs(binsWidthA - binsWidthB) < tolerance) {
0150     ACTS_VERBOSE("==> binWidths same: " << binsWidthA);
0151
0152     double min = std::numeric_limits<double>::signaling_NaN();
0153     double max = std::numeric_limits<double>::signaling_NaN();
0154
0155     if (direction == AxisDirection::AxisR) {
0156       ACTS_VERBOSE("Performing asymmetric merge");
0157       min = axisA.getMin();
0158       max = axisB.getMax();
0159     } else {
0160       ACTS_VERBOSE("Performing symmetrized merge");
0161
0162       double halfWidth =
0163           (axisA.getMax() - axisA.getMin() + axisB.getMax() - axisB.getMin()) /
0164           2.0;
0165
0166       min = -halfWidth;
0167       max = halfWidth;
0168     }
0169
0170     return makeGrid(mergedSurface, direction, logger,
0171                     std::tuple{min, max, axisA.getNBins() + axisB.getNBins()},
0172                     otherAxis, prepend);
0173
0174   } else {
0175     ACTS_VERBOSE("==> binWidths differ: " << binsWidthA << " vs " << binsWidthB
0176                                           << " ~> variable merge");
0177
0178     std::unique_ptr<GridPortalLink> mergedPortalLink =
0179         mergeVariable(mergedSurface, axisA, axisB, tolerance, direction, logger,
0180                       otherAxis, prepend);
0181     return mergedPortalLink;
0182   }
0183 }
0184
0185 std::unique_ptr<GridPortalLink> colinearMerge(
0186     const std::shared_ptr<RegularSurface>& mergedSurface, const IAxis& axisA,
0187     const IAxis& axisB, double tolerance, AxisDirection direction,
0188     const Logger& logger, const IAxis* otherAxis, bool prepend) {
0189   AxisType aType = axisA.getType();
0190   AxisType bType = axisB.getType();
0191   if (axisA.getBoundaryType() != axisB.getBoundaryType()) {
0192     ACTS_WARNING("AxisBoundaryTypes are different");
0193     return nullptr;
0194   }
0195
0196   if (axisA.getBoundaryType() != AxisBoundaryType::Bound) {
0197     ACTS_WARNING("AxisBoundaryType is not Bound, cannot do colinear merge");
0198     return nullptr;
0199   }
0200
0201   // convenience only
0202   auto mergeVariableLocal = [&] {
0203     return mergeVariable(mergedSurface, axisA, axisB, tolerance, direction,
0204                          logger, otherAxis, prepend);
0205   };
0206
0207   if (aType == AxisType::Equidistant && bType == AxisType::Equidistant) {
0208     auto mergedPortalLink =
0209         mergeEquidistant(mergedSurface, axisA, axisB, tolerance, direction,
0210                          logger, otherAxis, prepend);
0211     return mergedPortalLink;
0212   } else if (aType == AxisType::Variable && bType == AxisType::Variable) {
0213     ACTS_VERBOSE("===> variable merge");
0214     auto mergedPortalLink = mergeVariableLocal();
0215     return mergedPortalLink;
0216   } else if (aType == AxisType::Equidistant && bType == AxisType::Variable) {
0217     ACTS_VERBOSE("===> mixed merged");
0218     auto mergedPortalLink = mergeVariableLocal();
0219     return mergedPortalLink;
0220   } else {
0221     ACTS_VERBOSE("===> mixed merged");
0222     auto mergedPortalLink = mergeVariableLocal();
0223     return mergedPortalLink;
0224   }
0225 }
0226
0227 std::unique_ptr<PortalLinkBase> mergeGridPortals(
0228     const GridPortalLink* a, const GridPortalLink* b,
0229     const RegularSurface* surfaceA, const RegularSurface* surfaceB,
0230     AxisDirection loc0, AxisDirection loc1, AxisDirection direction,
0231     const Logger& logger) {
0232   assert(surfaceA != nullptr);
0233   assert(surfaceB != nullptr);
0234   assert(surfaceA->type() == surfaceB->type());
0235   assert(a->dim() == 2 || a->dim() == 1);
0236   assert(a->dim() == b->dim());
0237
0238   ACTS_VERBOSE(" - a: " << a->surface().bounds());
0239   ACTS_VERBOSE(" - b: " << b->surface().bounds());
0240
0241   // This tolerance is used checking bin equivalence. It's not intended to be
0242   // user configurable, as we don't foresee cases where the tolerance would have
0243   // to be adjusteed to the input geometry.
0244   constexpr auto tolerance = s_onSurfaceTolerance;
0245
0246   std::shared_ptr<RegularSurface> mergedSurface = nullptr;
0247   bool reversed = false;
0248
0249   if (const auto* cylinderA = dynamic_cast<const CylinderSurface*>(surfaceA);
0250       cylinderA != nullptr) {
0251     std::tie(mergedSurface, reversed) =
0252         cylinderA->mergedWith(dynamic_cast<const CylinderSurface&>(*surfaceB),
0253                               direction, true, logger);
0254   } else if (const auto* discA = dynamic_cast<const DiscSurface*>(surfaceA);
0255              discA != nullptr) {
0256     std::tie(mergedSurface, reversed) = discA->mergedWith(
0257         dynamic_cast<const DiscSurface&>(*surfaceB), direction, true, logger);
0258   } else if (const auto* planeA = dynamic_cast<const PlaneSurface*>(surfaceA);
0259              planeA != nullptr) {
0260     throw std::logic_error{"Plane surfaces not implemented yet"};
0261   } else {
0262     throw std::invalid_argument{"Unsupported surface type"};
0263   }
0264
0265   ACTS_VERBOSE("Merged surface: " << mergedSurface->bounds());
0266   ACTS_VERBOSE("~> reversed? " << std::boolalpha << reversed);
0267
0268   // Normalize ordering of grid portals and surfaces: a is always at lower
0269   // range than b
0270   if (reversed) {
0271     ACTS_VERBOSE("Swapping grid portals and surfaces after merging");
0272     std::swap(surfaceA, surfaceB);
0273     std::swap(a, b);
0274   }
0275
0276   // We do this here, because this is the only place where we've normalized
0277   // the ordering of grids just above: a is always "before" b, so we can
0278   // iterate over the bins in a unified way to fill the merged grid.
0279   auto fillGrid = [&](std::unique_ptr<GridPortalLink> merged) {
0280     if (merged != nullptr) {
0281       ACTS_VERBOSE("Post processing merged grid: " << merged->grid());
0282       GridPortalLink::fillMergedGrid(*a, *b, *merged, direction, logger);
0283     }
0284     return merged;
0285   };
0286
0287   if (a->dim() == 1) {
0288     ACTS_VERBOSE("Merge two 1D GridPortalLinks on " << surfaceA->name()
0289                                                     << "s in " << direction);
0290     if (a->direction() != b->direction()) {
0291       ACTS_VERBOSE("GridPortalLinks have different directions");
0292       ACTS_VERBOSE("=> cross merge");
0293
0294       // Find the one which is already binned in the merging direction
0295       const auto* aligned = a->direction() == direction ? a : b;
0296       const auto* alignedSurface =
0297           a->direction() == direction ? surfaceA : surfaceB;
0298       const auto* other = a->direction() == direction ? b : a;
0299       const auto* otherSurface =
0300           a->direction() == direction ? surfaceB : surfaceA;
0301
0302       // Extend the aligned one by the other one's axis
0303       auto aligned2D = aligned->extendTo2d(other->grid().axes().front());
0304       // Edtend the other one with a single bin
0305       auto other2D = other->extendTo2d(nullptr);
0306
0307       assert(aligned2D != nullptr);
0308       assert(other2D != nullptr);
0309
0310       ACTS_VERBOSE("Expanded grids:");
0311       ACTS_VERBOSE(" - aligned: " << aligned2D->grid());
0312       ACTS_VERBOSE(" - other: " << other2D->grid());
0313
0314       // @Fix ordering (a,b) should already be in good order, to save us one additional roundtrip
0315       if (a->direction() == direction) {
0316         return mergeGridPortals(aligned2D.get(), other2D.get(), alignedSurface,
0317                                 otherSurface, loc0, loc1, direction, logger);
0318       } else {
0319         return mergeGridPortals(other2D.get(), aligned2D.get(), otherSurface,
0320                                 alignedSurface, loc0, loc1, direction, logger);
0321       }
0322     }
0323
0324     const auto& axisA = *a->grid().axes().front();
0325     const auto& axisB = *b->grid().axes().front();
0326
0327     if (direction == loc0) {
0328       ACTS_VERBOSE("Grids are binned along " << a->direction());
0329       if (a->direction() == loc0) {
0330         ACTS_VERBOSE("=> colinear merge");
0331
0332         return fillGrid(colinearMerge(mergedSurface, axisA, axisB, tolerance,
0333                                       loc0, logger, nullptr, false));
0334
0335       } else {
0336         ACTS_VERBOSE("=> parallel merge");
0337
0338         auto a2D = a->extendTo2d(nullptr);
0339         auto b2D = b->extendTo2d(nullptr);
0340         assert(a2D != nullptr);
0341         assert(b2D != nullptr);
0342
0343         ACTS_VERBOSE("Expanded grids:");
0344         ACTS_VERBOSE(" - a: " << a2D->grid());
0345         ACTS_VERBOSE(" - b: " << b2D->grid());
0346
0347         return mergeGridPortals(a2D.get(), b2D.get(), surfaceA, surfaceB, loc0,
0348                                 loc1, direction, logger);
0349       }
0350     } else if (direction == loc1) {
0351       ACTS_VERBOSE("Grids are binned along " << a->direction());
0352       if (a->direction() == loc1) {
0353         ACTS_VERBOSE("=> colinear merge");
0354
0355         return fillGrid(colinearMerge(mergedSurface, axisA, axisB, tolerance,
0356                                       loc1, logger, nullptr, false));
0357
0358       } else {
0359         ACTS_VERBOSE("=> parallel merge");
0360         auto a2D = a->extendTo2d(nullptr);
0361         auto b2D = b->extendTo2d(nullptr);
0362         assert(a2D != nullptr);
0363         assert(b2D != nullptr);
0364
0365         ACTS_VERBOSE("Expanded grids:");
0366         ACTS_VERBOSE(" - a: " << a2D->grid());
0367         ACTS_VERBOSE(" - b: " << b2D->grid());
0368         return mergeGridPortals(a2D.get(), b2D.get(), surfaceA, surfaceB, loc0,
0369                                 loc1, direction, logger);
0370       }
0371
0372     } else {
0373       ACTS_ERROR("Invalid binning direction: " << direction);
0374       throw std::invalid_argument{"Invalid binning direction"};
0375     }
0376   } else {
0377     ACTS_VERBOSE("Merging two 2D GridPortalLinks on "
0378                  << surfaceA->name() << "s in " << direction << " direction");
0379
0380     const auto& loc0AxisA = *a->grid().axes().front();
0381     const auto& loc1AxisA = *a->grid().axes().back();
0382     const auto& loc0AxisB = *b->grid().axes().front();
0383     const auto& loc1AxisB = *b->grid().axes().back();
0384
0385     if (direction == loc0) {
0386       ACTS_VERBOSE("=> colinear merge along " << loc0);
0387       ACTS_VERBOSE("--> Checking if " << loc1 << " axes are identical");
0388
0389       if (loc1AxisA != loc1AxisB) {
0390         ACTS_WARNING("    ~> "
0391                      << loc1
0392                      << " axes are not identical, falling back to composite "
0393                         "merging");
0394         return nullptr;
0395       }
0396       ACTS_VERBOSE("    ~> they are!");
0397
0398       ACTS_VERBOSE("    ~> " << loc1 << " axis: " << loc1AxisA);
0399       return fillGrid(colinearMerge(mergedSurface, loc0AxisA, loc0AxisB,
0400                                     tolerance, loc0, logger, &loc1AxisA,
0401                                     false));
0402
0403     } else if (direction == loc1) {
0404       ACTS_VERBOSE("=> colinear merge along " << loc1);
0405       ACTS_VERBOSE("--> Checking if " << loc0 << " axes are identical");
0406
0407       if (loc0AxisA != loc0AxisB) {
0408         ACTS_WARNING("    ~> "
0409                      << loc0
0410                      << " axes are not identical, falling back to composite "
0411                         "merging");
0412         return nullptr;
0413       }
0414       ACTS_VERBOSE("    ~> they are!");
0415
0416       ACTS_VERBOSE("    ~> " << loc0 << " axis: " << loc0AxisA);
0417       return fillGrid(colinearMerge(mergedSurface, loc1AxisA, loc1AxisB,
0418                                     tolerance, loc1, logger, &loc0AxisA, true));
0419
0420     } else {
0421       ACTS_ERROR("Invalid binning direction: " << a->direction());
0422       throw std::invalid_argument{"Invalid binning direction"};
0423     }
0424   }
0425 }
0426
0427 std::unique_ptr<PortalLinkBase> mergeGridPortals(const GridPortalLink* a,
0428                                                  const GridPortalLink* b,
0429                                                  AxisDirection direction,
0430                                                  const Logger& logger) {
0431   using enum AxisDirection;
0432   assert(a->dim() == 2 || a->dim() == 1);
0433   assert(b->dim() == 2 || b->dim() == 1);
0434
0435   if (a->dim() < b->dim()) {
0436     return mergeGridPortals(b, a, direction, logger);
0437   }
0438
0439   ACTS_VERBOSE("Merging GridPortalLinks along " << direction << ":");
0440   ACTS_VERBOSE(" - a: " << a->grid() << " along: " << a->direction());
0441   ACTS_VERBOSE(" - b: " << b->grid() << " along: " << b->direction());
0442
0443   const auto* cylinder = dynamic_cast<const CylinderSurface*>(&a->surface());
0444   const auto* disc = dynamic_cast<const DiscSurface*>(&a->surface());
0445
0446   if (a->dim() == b->dim()) {
0447     ACTS_VERBOSE("Grid both have same dimension: " << a->dim());
0448
0449     if (cylinder != nullptr) {
0450       return mergeGridPortals(
0451           a, b, cylinder, &dynamic_cast<const CylinderSurface&>(b->surface()),
0452           AxisRPhi, AxisZ, direction, logger);
0453     } else if (disc != nullptr) {
0454       return mergeGridPortals(a, b, disc,
0455                               &dynamic_cast<const DiscSurface&>(b->surface()),
0456                               AxisR, AxisPhi, direction, logger);
0457     } else {
0458       // @TODO: Support PlaneSurface
0459       ACTS_VERBOSE("Surface type is not supported here, falling back");
0460       return nullptr;
0461     }
0462   } else {
0463     ACTS_VERBOSE("Grids have different dimension, extending rhs to 2D");
0464     const IAxis* otherAxis = nullptr;
0465
0466     if (b->direction() == direction) {
0467       ACTS_VERBOSE("1D grid is binned in merging direction " << direction);
0468       ACTS_VERBOSE("~> Adding complementary axis");
0469
0470       if (cylinder != nullptr) {
0471         otherAxis = direction == AxisRPhi ? a->grid().axes().back()
0472                                           : a->grid().axes().front();
0473       } else if (disc != nullptr) {
0474         otherAxis = direction == AxisR ? a->grid().axes().back()
0475                                        : a->grid().axes().front();
0476       } else {
0477         ACTS_VERBOSE("Surface type is not supported here, falling back");
0478         return nullptr;
0479       }
0480     } else {
0481       ACTS_VERBOSE("1D grid is binned in complementary direction");
0482     }
0483
0484     auto b2D = b->extendTo2d(otherAxis);
0485     ACTS_VERBOSE("-> new grid: " << b2D->grid());
0486     return mergeGridPortals(a, b2D.get(), direction, logger);
0487   }
0488 }
0489
0490 }  // namespace
0491
0492 void GridPortalLink::fillMergedGrid(const GridPortalLink& a,
0493                                     const GridPortalLink& b,
0494                                     GridPortalLink& merged,
0495                                     AxisDirection direction,
0496                                     const Logger& logger) {
0497   ACTS_VERBOSE("Filling merged grid");
0498   assert(a.dim() == b.dim());
0499   assert(a.direction() == b.direction());
0500   const auto locBinsA = a.numLocalBins();
0501   const auto locBinsB = b.numLocalBins();
0502
0503   ACTS_VERBOSE("a: " << a.grid());
0504   ACTS_VERBOSE("b: " << b.grid());
0505   ACTS_VERBOSE("merged: " << merged.grid());
0506
0507   if (a.dim() == 1) {
0508     std::size_t nBinsA = locBinsA.at(0);
0509     std::size_t nBinsB = locBinsB.at(0);
0510
0511     for (std::size_t i = 1; i <= nBinsA; ++i) {
0512       merged.atLocalBins({i}) = a.atLocalBins({i});
0513     }
0514     for (std::size_t i = 1; i <= nBinsB; ++i) {
0515       merged.atLocalBins({nBinsA + i}) = b.atLocalBins({i});
0516     }
0517   } else {
0518     if (a.direction() == direction) {
0519       std::size_t nBinsA = locBinsA.at(0);
0520       std::size_t nBinsB = locBinsB.at(0);
0521       std::size_t nBinsCommon = locBinsB.at(1);
0522
0523       for (std::size_t i = 1; i <= nBinsA; ++i) {
0524         for (std::size_t j = 1; j <= nBinsCommon; ++j) {
0525           merged.atLocalBins({i, j}) = a.atLocalBins({i, j});
0526         }
0527       }
0528
0529       for (std::size_t i = 1; i <= nBinsB; ++i) {
0530         for (std::size_t j = 1; j <= nBinsCommon; ++j) {
0531           std::size_t ti = i + nBinsA;
0532           merged.atLocalBins({ti, j}) = b.atLocalBins({i, j});
0533         }
0534       }
0535     } else {
0536       std::size_t nBinsA = locBinsA.at(1);
0537       std::size_t nBinsB = locBinsB.at(1);
0538       std::size_t nBinsCommon = locBinsB.at(0);
0539
0540       for (std::size_t i = 1; i <= nBinsCommon; ++i) {
0541         for (std::size_t j = 1; j <= nBinsA; ++j) {
0542           merged.atLocalBins({i, j}) = a.atLocalBins({i, j});
0543         }
0544       }
0545
0546       for (std::size_t i = 1; i <= nBinsCommon; ++i) {
0547         for (std::size_t j = 1; j <= nBinsB; ++j) {
0548           std::size_t tj = j + nBinsA;
0549           merged.atLocalBins({i, tj}) = b.atLocalBins({i, j});
0550         }
0551       }
0552     }
0553   }
0554 }
0555
0556 std::unique_ptr<PortalLinkBase> GridPortalLink::merge(const GridPortalLink& a,
0557                                                       const GridPortalLink& b,
0558                                                       AxisDirection direction,
0559                                                       const Logger& logger) {
0560   ACTS_VERBOSE("Merging two GridPortalLinks");
0561
0562   checkMergePreconditions(a, b, direction);
0563
0564   auto merged = mergeGridPortals(&a, &b, direction, logger);
0565   if (merged == nullptr) {
0566     ACTS_WARNING("Grid merging failed returning nullptr");
0567     return nullptr;
0568   }
0569
0570   return merged;
0571 }
0572
0573 }  // namespace Acts