accel/detail/MagneticFieldUtils.hh

0001 //------------------------------- -*- C++ -*- -------------------------------//
0002 // Copyright Celeritas contributors: see top-level COPYRIGHT file for details
0003 // SPDX-License-Identifier: (Apache-2.0 OR MIT)
0004 //---------------------------------------------------------------------------//
0005 //! \file accel/detail/MagneticFieldUtils.hh
0006 //! \brief Utilities for magnetic field map creation
0007 //---------------------------------------------------------------------------//
0008 #pragma once
0009
0010 #include <G4Field.hh>
0011 #include <corecel/Assert.hh>
0012
0013 #include "corecel/Types.hh"
0014 #include "corecel/cont/Array.hh"
0015 #include "corecel/data/HyperslabIndexer.hh"
0016 #include "geocel/GeantGeoUtils.hh"
0017
0018 namespace celeritas
0019 {
0020 namespace detail
0021 {
0022 //---------------------------------------------------------------------------//
0023 /*!
0024  * Common field sampling setup and execution.
0025  *
0026  * This function encapsulates the shared pattern of:
0027  * 1. Getting the G4 field
0028  * 2. Setting up the HyperslabIndexer
0029  * 3. Sampling the field on a grid
0030  *
0031  * \param field_data Output parameter array to store field values (must be
0032  *                   pre-allocated with size equal to the product of all dims)
0033  * \param dims Grid dimensions
0034  * \param calc_position Callable that computes position given
0035  *                           grid indices. Must have signature:
0036  *                           Array<G4double, 4>(size_type, size_type,
0037  *                           size_type) returning [x, y, z, 0] coordinates
0038  * \param convert_field Callable that converts field from G4 to
0039  *                        native units in the correct coordinate space. Must
0040  *                        have signature: void(Array<G4double, 3> const&,
0041  *                        real_type*) taking G4 field [Bx, By, Bz] and writing
0042  *                        converted values to output pointer
0043  */
0044 template<typename PositionCalc, typename FieldConverter>
0045 inline void setup_and_sample_field(real_type* field_data,
0046                                    Array<size_type, 4> const& dims,
0047                                    PositionCalc const& calc_position,
0048                                    FieldConverter const& convert_field)
0049 {
0050     HyperslabIndexer const flat_index{dims};
0051     G4Field const* g4field = celeritas::geant_field();
0052     CELER_VALIDATE(g4field,
0053                    << "no Geant4 global field has been set: cannot build "
0054                       "magnetic field map");
0055
0056     Array<G4double, 3> bfield;
0057
0058     for (size_type i = 0; i < dims[0]; ++i)
0059     {
0060         for (size_type j = 0; j < dims[1]; ++j)
0061         {
0062             for (size_type k = 0; k < dims[2]; ++k)
0063             {
0064                 // Calculate position for this grid point
0065                 Array<G4double, 4> pos = calc_position(i, j, k);
0066
0067                 // Sample field at this position
0068                 g4field->GetFieldValue(pos.data(), bfield.data());
0069
0070                 // Convert and store field values
0071                 auto* cur_bfield = field_data + flat_index(i, j, k, 0);
0072                 convert_field(bfield, cur_bfield);
0073             }
0074         }
0075     }
0076 }
0077
0078 //---------------------------------------------------------------------------//
0079 }  // namespace detail
0080 }  // namespace celeritas