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