File indexing completed on 2026-07-08 07:50:54
0001
0002
0003
0004
0005
0006
0007
0008
0009 #pragma once
0010
0011 #include "Acts/Definitions/Algebra.hpp"
0012 #include "Acts/Geometry/GeometryContext.hpp"
0013 #include "Acts/Surfaces/BoundaryTolerance.hpp"
0014 #include "Acts/Surfaces/CylinderBounds.hpp"
0015 #include "Acts/Surfaces/RegularSurface.hpp"
0016 #include "Acts/Surfaces/Surface.hpp"
0017 #include "Acts/Utilities/AnyGridView.hpp"
0018 #include "Acts/Utilities/AxisDefinitions.hpp"
0019 #include "Acts/Utilities/Grid.hpp"
0020 #include "Acts/Utilities/IAxis.hpp"
0021 #include "Acts/Utilities/Intersection.hpp"
0022
0023 #include <iostream>
0024 #include <limits>
0025 #include <vector>
0026
0027 namespace Acts {
0028
0029 using SurfaceVector = std::vector<const Surface*>;
0030
0031
0032
0033
0034
0035
0036
0037 class SurfaceArray {
0038 public:
0039
0040 struct ISurfaceGridLookup {
0041
0042
0043
0044 virtual void fill(const GeometryContext& gctx,
0045 const SurfaceVector& surfaces) = 0;
0046
0047
0048
0049
0050
0051
0052 virtual const SurfaceVector& lookup(const Vector3& position,
0053 const Vector3& direction) const = 0;
0054
0055
0056
0057
0058
0059 virtual SurfaceVector& lookup(std::size_t bin) = 0;
0060
0061
0062
0063
0064
0065 virtual const SurfaceVector& lookup(std::size_t bin) const = 0;
0066
0067
0068
0069
0070
0071
0072 virtual const SurfaceVector& neighbors(const Vector3& position,
0073 const Vector3& direction) const = 0;
0074
0075
0076
0077
0078 virtual std::size_t size() const = 0;
0079
0080
0081
0082
0083 virtual Vector3 getBinCenter(std::size_t bin) const = 0;
0084
0085
0086
0087
0088 virtual std::vector<const IAxis*> getAxes() const = 0;
0089
0090
0091
0092 virtual std::optional<AnyGridConstView<SurfaceVector>> getGridView()
0093 const = 0;
0094
0095
0096
0097 virtual const Surface* surfaceRepresentation() const = 0;
0098
0099
0100
0101
0102
0103
0104 virtual bool isValidBin(std::size_t bin) const = 0;
0105
0106
0107
0108
0109 virtual std::vector<AxisDirection> binningValues() const { return {}; };
0110
0111
0112 virtual ~ISurfaceGridLookup() = 0;
0113 };
0114
0115
0116
0117 template <class Axis1, class Axis2>
0118 struct SurfaceGridLookup : ISurfaceGridLookup {
0119
0120 using Grid_t = Grid<SurfaceVector, Axis1, Axis2>;
0121
0122
0123
0124
0125
0126
0127 SurfaceGridLookup(std::shared_ptr<RegularSurface> representative,
0128 double tolerance, std::tuple<Axis1, Axis2> axes,
0129 std::vector<AxisDirection> bValues = {})
0130 : m_representative(std::move(representative)),
0131 m_tolerance(tolerance),
0132 m_grid(std::move(axes)),
0133 m_binValues(std::move(bValues)) {
0134 m_neighborMap.resize(m_grid.size());
0135 }
0136
0137
0138
0139
0140
0141
0142
0143
0144
0145
0146 void fill(const GeometryContext& gctx,
0147 const SurfaceVector& surfaces) override {
0148 for (const Surface* surface : surfaces) {
0149 const std::size_t globalBin = fillSurfaceToBinMapping(gctx, *surface);
0150 if (globalBin == 0) {
0151 continue;
0152 }
0153
0154 fillBinToSurfaceMapping(gctx, *surface, globalBin);
0155 }
0156
0157 populateNeighborCache();
0158 }
0159
0160 const SurfaceVector& lookup(const Vector3& position,
0161 const Vector3& direction) const override {
0162 return m_grid.at(findGlobalBin(position, direction,
0163 std::numeric_limits<double>::infinity()));
0164 }
0165
0166
0167
0168
0169
0170 SurfaceVector& lookup(std::size_t bin) override { return m_grid.at(bin); }
0171
0172
0173
0174
0175
0176 const SurfaceVector& lookup(std::size_t bin) const override {
0177 return m_grid.at(bin);
0178 }
0179
0180
0181
0182
0183
0184
0185 const SurfaceVector& neighbors(const Vector3& position,
0186 const Vector3& direction) const override {
0187 return m_neighborMap.at(findGlobalBin(
0188 position, direction, std::numeric_limits<double>::infinity()));
0189 }
0190
0191
0192
0193
0194 std::size_t size() const override { return m_grid.size(); }
0195
0196
0197
0198
0199 std::vector<AxisDirection> binningValues() const override {
0200 return m_binValues;
0201 }
0202
0203
0204
0205
0206 Vector3 getBinCenter(std::size_t bin) const override {
0207 auto gctx = GeometryContext::dangerouslyDefaultConstruct();
0208 return getBinCenterImpl(gctx, bin);
0209 }
0210
0211
0212
0213
0214 std::vector<const IAxis*> getAxes() const override {
0215 auto arr = m_grid.axes();
0216 return std::vector<const IAxis*>(arr.begin(), arr.end());
0217 }
0218
0219 std::optional<AnyGridConstView<SurfaceVector>> getGridView()
0220 const override {
0221 return AnyGridConstView<SurfaceVector>{m_grid};
0222 }
0223
0224 const Surface* surfaceRepresentation() const override {
0225 return m_representative.get();
0226 }
0227
0228
0229
0230
0231
0232
0233 bool isValidBin(std::size_t bin) const override {
0234 std::array<std::size_t, 2> indices = m_grid.localBinsFromGlobalBin(bin);
0235 std::array<std::size_t, 2> nBins = m_grid.numLocalBins();
0236 for (std::size_t i = 0; i < indices.size(); ++i) {
0237 std::size_t idx = indices.at(i);
0238 if (idx <= 0 || idx >= nBins.at(i) + 1) {
0239 return false;
0240 }
0241 }
0242 return true;
0243 }
0244
0245 private:
0246
0247 std::size_t fillSurfaceToBinMapping(const GeometryContext& gctx,
0248 const Surface& surface) {
0249 const Vector3 pos = surface.referencePosition(gctx, AxisDirection::AxisR);
0250 const Vector3 normal = m_representative->normal(gctx, pos);
0251 const std::size_t globalBin = findGlobalBin(pos, normal, m_tolerance);
0252 if (globalBin != 0) {
0253 m_grid.at(globalBin).push_back(&surface);
0254 }
0255 return globalBin;
0256 };
0257
0258
0259 void fillBinToSurfaceMapping(const GeometryContext& gctx,
0260 const Surface& surface, std::size_t startBin) {
0261 const std::array<std::size_t, 2> startIndices =
0262 m_grid.localBinsFromGlobalBin(startBin);
0263 const auto startNeighborIndices =
0264 m_grid.neighborHoodIndices(startIndices, 1u);
0265
0266 std::set<std::size_t> visited({startBin});
0267 std::vector<std::size_t> queue(startNeighborIndices.begin(),
0268 startNeighborIndices.end());
0269
0270 while (!queue.empty()) {
0271 const std::size_t current = queue.back();
0272 queue.pop_back();
0273 if (visited.contains(current)) {
0274 continue;
0275 }
0276
0277 const std::array<std::size_t, 2> currentIndices =
0278 m_grid.localBinsFromGlobalBin(current);
0279 visited.insert(current);
0280
0281 const std::array<double, 2> gridLocal =
0282 m_grid.binCenter(currentIndices);
0283 const Vector2 surfaceLocal = gridToSurfaceLocal(gridLocal);
0284 const Vector3 normal = m_representative->normal(gctx, surfaceLocal);
0285 const Vector3 global =
0286 m_representative->localToGlobal(gctx, surfaceLocal, normal);
0287
0288 const Intersection3D intersection =
0289 surface.intersect(gctx, global, normal, BoundaryTolerance::None())
0290 .closest();
0291 if (!intersection.isValid() ||
0292 std::abs(intersection.pathLength()) > m_tolerance) {
0293 continue;
0294 }
0295 m_grid.at(current).push_back(&surface);
0296
0297 const auto neighborIndices =
0298 m_grid.neighborHoodIndices(currentIndices, 1u);
0299 queue.insert(queue.end(), neighborIndices.begin(),
0300 neighborIndices.end());
0301 }
0302 };
0303
0304 void populateNeighborCache() {
0305
0306 for (std::size_t i = 0; i < m_grid.size(); i++) {
0307 if (!isValidBin(i)) {
0308 continue;
0309 }
0310 const std::array<std::size_t, 2> indices =
0311 m_grid.localBinsFromGlobalBin(i);
0312 std::vector<const Surface*>& neighbors = m_neighborMap.at(i);
0313 neighbors.clear();
0314
0315 for (std::size_t idx : m_grid.neighborHoodIndices(indices, 1u)) {
0316 const std::vector<const Surface*>& binContent = m_grid.at(idx);
0317 std::copy(binContent.begin(), binContent.end(),
0318 std::back_inserter(neighbors));
0319 }
0320
0321 std::ranges::sort(neighbors);
0322 auto last = std::ranges::unique(neighbors);
0323 neighbors.erase(last.begin(), last.end());
0324 neighbors.shrink_to_fit();
0325 }
0326 }
0327
0328 Vector3 getBinCenterImpl(const GeometryContext& gctx,
0329 std::size_t bin) const {
0330 const std::array<double, 2> gridLocal =
0331 m_grid.binCenter(m_grid.localBinsFromGlobalBin(bin));
0332 const Vector2 surfaceLocal = gridToSurfaceLocal(gridLocal);
0333 return m_representative->localToGlobal(gctx, surfaceLocal);
0334 }
0335
0336 const CylinderBounds* getCylinderBounds() const {
0337 return dynamic_cast<const CylinderBounds*>(&m_representative->bounds());
0338 }
0339
0340 Vector2 gridToSurfaceLocal(std::array<double, 2> gridLocal) const {
0341 Vector2 surfaceLocal = Eigen::Map<Vector2>(gridLocal.data());
0342 if (const CylinderBounds* bounds = getCylinderBounds();
0343 bounds != nullptr) {
0344 surfaceLocal[0] *= bounds->get(CylinderBounds::eR);
0345 }
0346 return surfaceLocal;
0347 }
0348 std::array<double, 2> surfaceToGridLocal(Vector2 local) const {
0349 std::array<double, 2> gridLocal = {local[0], local[1]};
0350 if (const CylinderBounds* bounds = getCylinderBounds();
0351 bounds != nullptr) {
0352 gridLocal[0] /= bounds->get(CylinderBounds::eR);
0353 }
0354 return gridLocal;
0355 }
0356
0357 std::size_t findGlobalBin(const Vector3& position, const Vector3& direction,
0358 double tolerance) const {
0359 auto gctx = GeometryContext::dangerouslyDefaultConstruct();
0360
0361 const Intersection3D intersection =
0362 m_representative
0363 ->intersect(gctx, position, direction,
0364 BoundaryTolerance::Infinite())
0365 .closest();
0366 if (!intersection.isValid() ||
0367 std::abs(intersection.pathLength()) > tolerance) {
0368 return 0;
0369 }
0370 const Vector2 surfaceLocal =
0371 m_representative
0372 ->globalToLocal(gctx, intersection.position(), direction)
0373 .value();
0374 const std::array<double, 2> gridLocal = surfaceToGridLocal(surfaceLocal);
0375 return m_grid.globalBinFromPosition(gridLocal);
0376 }
0377
0378 std::shared_ptr<RegularSurface> m_representative;
0379 double m_tolerance{};
0380 Grid_t m_grid;
0381 std::vector<AxisDirection> m_binValues;
0382 std::vector<SurfaceVector> m_neighborMap;
0383 };
0384
0385
0386
0387 struct SingleElementLookup : ISurfaceGridLookup {
0388
0389
0390 explicit SingleElementLookup(SurfaceVector::value_type element)
0391 : m_element({element}) {}
0392
0393
0394
0395 explicit SingleElementLookup(const SurfaceVector& elements)
0396 : m_element(elements) {}
0397
0398
0399
0400 const SurfaceVector& lookup(const Vector3& ,
0401 const Vector3& ) const override {
0402 return m_element;
0403 }
0404
0405
0406
0407 SurfaceVector& lookup(std::size_t ) override { return m_element; }
0408
0409
0410
0411 const SurfaceVector& lookup(std::size_t ) const override {
0412 return m_element;
0413 }
0414
0415
0416
0417 const SurfaceVector& neighbors(
0418 const Vector3& ,
0419 const Vector3& ) const override {
0420 return m_element;
0421 }
0422
0423
0424
0425 std::size_t size() const override { return 1; }
0426
0427
0428
0429 Vector3 getBinCenter(std::size_t ) const override {
0430 return Vector3(0, 0, 0);
0431 }
0432
0433
0434
0435 std::vector<const IAxis*> getAxes() const override { return {}; }
0436
0437 std::optional<AnyGridConstView<SurfaceVector>> getGridView()
0438 const override {
0439 return std::nullopt;
0440 }
0441
0442 const Surface* surfaceRepresentation() const override { return nullptr; }
0443
0444
0445
0446 void fill(const GeometryContext& ,
0447 const SurfaceVector& ) override {}
0448
0449
0450
0451 bool isValidBin(std::size_t ) const override { return true; }
0452
0453 private:
0454 SurfaceVector m_element;
0455 };
0456
0457
0458
0459
0460
0461
0462
0463
0464 explicit SurfaceArray(std::unique_ptr<ISurfaceGridLookup> gridLookup,
0465 std::vector<std::shared_ptr<const Surface>> surfaces,
0466 const Transform3& transform = Transform3::Identity());
0467
0468
0469
0470 explicit SurfaceArray(std::shared_ptr<const Surface> srf);
0471
0472
0473
0474
0475
0476
0477 const SurfaceVector& at(const Vector3& position,
0478 const Vector3& direction) const {
0479 return p_gridLookup->lookup(position, direction);
0480 }
0481
0482
0483
0484
0485 SurfaceVector& at(std::size_t bin) { return p_gridLookup->lookup(bin); }
0486
0487
0488
0489
0490 const SurfaceVector& at(std::size_t bin) const {
0491 return p_gridLookup->lookup(bin);
0492 }
0493
0494
0495
0496
0497
0498
0499
0500
0501 const SurfaceVector& neighbors(const Vector3& position,
0502 const Vector3& direction) const {
0503 return p_gridLookup->neighbors(position, direction);
0504 }
0505
0506
0507
0508
0509 std::size_t size() const { return p_gridLookup->size(); }
0510
0511
0512
0513
0514 Vector3 getBinCenter(std::size_t bin) const {
0515 return p_gridLookup->getBinCenter(bin);
0516 }
0517
0518
0519
0520
0521
0522
0523 const SurfaceVector& surfaces() const { return m_surfacesRawPointers; }
0524
0525
0526
0527
0528
0529 std::vector<const IAxis*> getAxes() const { return p_gridLookup->getAxes(); }
0530
0531
0532
0533
0534
0535
0536 bool isValidBin(std::size_t bin) const {
0537 return p_gridLookup->isValidBin(bin);
0538 }
0539
0540
0541
0542 const Transform3& transform() const { return m_transform; }
0543
0544
0545
0546
0547 std::vector<AxisDirection> binningValues() const {
0548 return p_gridLookup->binningValues();
0549 };
0550
0551
0552
0553
0554
0555 std::ostream& toStream(const GeometryContext& gctx, std::ostream& sl) const;
0556
0557
0558
0559 const ISurfaceGridLookup& gridLookup() const { return *p_gridLookup; }
0560
0561 private:
0562 std::unique_ptr<ISurfaceGridLookup> p_gridLookup;
0563
0564 std::vector<std::shared_ptr<const Surface>> m_surfaces;
0565
0566
0567 SurfaceVector m_surfacesRawPointers;
0568
0569
0570 Transform3 m_transform;
0571 };
0572
0573 }