ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
InterpolatedBFieldMap.hpp
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file InterpolatedBFieldMap.hpp
1 // This file is part of the Acts project.
2 //
3 // Copyright (C) 2016-2018 CERN for the benefit of the Acts project
4 //
5 // This Source Code Form is subject to the terms of the Mozilla Public
6 // License, v. 2.0. If a copy of the MPL was not distributed with this
7 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
8 
9 #pragma once
10 
11 #include <functional>
12 #include <optional>
13 #include <vector>
18 
19 namespace Acts {
20 
29 template <typename G>
31  public:
32  using Grid_t = G;
33  using FieldType = typename Grid_t::value_type;
34  static constexpr size_t DIM_POS = Grid_t::DIM;
35 
42  struct FieldCell {
44  static constexpr unsigned int N = 1 << DIM_POS;
45 
46  public:
59  std::array<double, DIM_POS> lowerLeft,
60  std::array<double, DIM_POS> upperRight,
61  std::array<Vector3D, N> fieldValues)
62  : m_transformPos(std::move(transformPos)),
63  m_lowerLeft(std::move(lowerLeft)),
64  m_upperRight(std::move(upperRight)),
65  m_fieldValues(std::move(fieldValues)) {}
66 
74  // defined in Interpolation.hpp
77  }
78 
84  bool isInside(const Vector3D& position) const {
85  const auto& gridCoordinates = m_transformPos(position);
86  for (unsigned int i = 0; i < DIM_POS; ++i) {
87  if (gridCoordinates[i] < m_lowerLeft.at(i) ||
88  gridCoordinates[i] >= m_upperRight.at(i)) {
89  return false;
90  }
91  }
92  return true;
93  }
94 
95  private:
97  std::function<ActsVectorD<DIM_POS>(const Vector3D&)> m_transformPos;
98 
100  std::array<double, DIM_POS> m_lowerLeft;
101 
103  std::array<double, DIM_POS> m_upperRight;
104 
109  std::array<Vector3D, N> m_fieldValues;
110  };
111 
121  std::function<ActsVectorD<DIM_POS>(const Vector3D&)> transformPos,
122  std::function<Vector3D(const FieldType&, const Vector3D&)>
123  transformBField,
124  Grid_t grid)
125  : m_transformPos(std::move(transformPos)),
126  m_transformBField(std::move(transformBField)),
127  m_grid(std::move(grid)) {}
128 
137  return m_transformBField(m_grid.interpolate(m_transformPos(position)),
138  position);
139  }
140 
149  const auto& gridPosition = m_transformPos(position);
150  const auto& indices = m_grid.localBinsFromPosition(gridPosition);
151  const auto& lowerLeft = m_grid.lowerLeftBinEdge(indices);
152  const auto& upperRight = m_grid.upperRightBinEdge(indices);
153 
154  // loop through all corner points
155  constexpr size_t nCorners = 1 << DIM_POS;
156  std::array<Vector3D, nCorners> neighbors;
157  const auto& cornerIndices = m_grid.closestPointsIndices(gridPosition);
158 
159  size_t i = 0;
160  for (size_t index : cornerIndices) {
161  neighbors.at(i++) = m_transformBField(m_grid.at(index), position);
162  }
163 
164  return FieldCell(m_transformPos, lowerLeft, upperRight,
165  std::move(neighbors));
166  }
167 
171  std::vector<size_t> getNBins() const {
172  auto nBinsArray = m_grid.numLocalBins();
173  return std::vector<size_t>(nBinsArray.begin(), nBinsArray.end());
174  }
175 
179  std::vector<double> getMin() const {
180  auto minArray = m_grid.minPosition();
181  return std::vector<double>(minArray.begin(), minArray.end());
182  }
183 
187  std::vector<double> getMax() const {
188  auto maxArray = m_grid.maxPosition();
189  return std::vector<double>(maxArray.begin(), maxArray.end());
190  }
191 
197  bool isInside(const Vector3D& position) const {
198  return m_grid.isInside(m_transformPos(position));
199  }
200 
204  const Grid_t& getGrid() const { return m_grid; }
205 
206  private:
208  std::function<ActsVectorD<DIM_POS>(const Vector3D&)> m_transformPos;
212  std::function<Vector3D(const FieldType&, const Vector3D&)> m_transformBField;
215 };
216 
232 template <typename Mapper_t>
234  public:
236  struct Config {
237  Config(Mapper_t m) : mapper(m) {}
238 
243  double scale = 1.;
244 
250  Mapper_t mapper;
251  };
252 
253  struct Cache {
257  Cache(std::reference_wrapper<const MagneticFieldContext> /*mcfg*/) {}
258 
259  std::optional<typename Mapper_t::FieldCell> fieldCell;
260  bool initialized = false;
261  };
262 
266  InterpolatedBFieldMap(Config config) : m_config(std::move(config)) {}
267 
271  Config getConfiguration() const { return m_config; }
272 
279  return m_config.mapper.getField(position);
280  }
281 
289  Vector3D getField(const Vector3D& position, Cache& cache) const {
290  if (!cache.fieldCell || !(*cache.fieldCell).isInside(position)) {
291  cache.fieldCell = getFieldCell(position);
292  }
293  return (*cache.fieldCell).getField(position);
294  }
295 
305  ActsMatrixD<3, 3>& /*derivative*/) const {
306  return m_config.mapper.getField(position);
307  }
308 
320  ActsMatrixD<3, 3>& /*derivative*/,
321  Cache& /*cache*/) const {
322  return m_config.mapper.getField(position);
323  }
324 
328  double getScale() const { return m_config.scale; }
329 
333  Mapper_t getMapper() const { return m_config.mapper; }
334 
340  bool isInside(const Vector3D& position) const {
341  return m_config.mapper.isInside(position);
342  }
343 
348 
349  private:
357  typename Mapper_t::FieldCell getFieldCell(const Vector3D& position) const {
358  return m_config.mapper.getFieldCell(position);
359  }
360 
363 };
364 
365 } // namespace Acts