BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSFieldEMInterpolated1D.cc
1/*
2Beam Delivery Simulation (BDSIM) Copyright (C) Royal Holloway,
3University of London 2001 - 2023.
4
5This file is part of BDSIM.
6
7BDSIM is free software: you can redistribute it and/or modify
8it under the terms of the GNU General Public License as published
9by the Free Software Foundation version 3 of the License.
10
11BDSIM is distributed in the hope that it will be useful, but
12WITHOUT ANY WARRANTY; without even the implied warranty of
13MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14GNU General Public License for more details.
15
16You should have received a copy of the GNU General Public License
17along with BDSIM. If not, see <http://www.gnu.org/licenses/>.
18*/
19#include "BDSDimensionType.hh"
20#include "BDSFieldEMInterpolated1D.hh"
21#include "BDSInterpolator1D.hh"
22
23#include "G4ThreeVector.hh"
24#include "G4Types.hh"
25
26#include <utility>
27
28BDSFieldEMInterpolated1D::BDSFieldEMInterpolated1D(BDSInterpolator1D* eInterpolatorIn,
29 BDSInterpolator1D* bInterpolatorIn,
30 const G4Transform3D& offset,
31 G4double eScalingIn,
32 G4double bScalingIn):
33 BDSFieldEMInterpolated(eInterpolatorIn, bInterpolatorIn, offset, eScalingIn, bScalingIn),
34 eInterpolator(eInterpolatorIn),
35 bInterpolator(bInterpolatorIn),
36 eDimensionIndex((eInterpolatorIn->FirstDimension()).underlying()),
37 eTime((eInterpolatorIn->FirstDimension()).underlying() > 2),
38 bDimensionIndex((bInterpolatorIn->FirstDimension()).underlying()),
39 bTime((bInterpolatorIn->FirstDimension()).underlying() > 2)
40{;}
41
42BDSFieldEMInterpolated1D::~BDSFieldEMInterpolated1D()
43{
44 delete eInterpolator;
45 delete bInterpolator;
46}
47
48std::pair<G4ThreeVector,G4ThreeVector> BDSFieldEMInterpolated1D::GetField(const G4ThreeVector& position,
49 const G4double t) const
50{
51 G4double eCoordinate = 0;
52 if (eTime)
53 {eCoordinate = t;}
54 else
55 {eCoordinate = position[eDimensionIndex];}
56 G4double bCoordinate = 0;
57 if (bTime)
58 {bCoordinate = t;}
59 else
60 {bCoordinate = position[bDimensionIndex];}
61 G4ThreeVector e = eInterpolator->GetInterpolatedValue(eCoordinate) * EScaling();
62 G4ThreeVector b = bInterpolator->GetInterpolatedValue(bCoordinate) * BScaling();
63 return std::make_pair(b,e);
64}
BDSInterpolator1D * eInterpolator
E Inteprolator the field is based on.
virtual std::pair< G4ThreeVector, G4ThreeVector > GetField(const G4ThreeVector &position, const G4double t=0) const
Return the interpolated field value at a given point.
const G4bool eTime
E Cache of whether to use time coordinate.
const G4bool bTime
B Cache of whether to use time coordinate.
const G4int bDimensionIndex
B Integer index to dimension to use.
BDSInterpolator1D * bInterpolator
B Interpolator the field is based on.
const G4int eDimensionIndex
E Integer index to dimension to use.
Class to provide scaling and a base class pointer for interpolator fields.
G4double EScaling() const
Accessor.
G4double BScaling() const
Accessor.
Interface for all 1D interpolators.
G4ThreeVector GetInterpolatedValue(G4double x) const
Public interface to a 1D interpolator. Returns Geant4 type as that's what will be needed.