BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSFieldEInterpolated3D.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 "BDSFieldEInterpolated3D.hh"
21#include "BDSInterpolator3D.hh"
22
23#include "G4ThreeVector.hh"
24#include "G4Types.hh"
25
26BDSFieldEInterpolated3D::BDSFieldEInterpolated3D(BDSInterpolator3D* interpolatorIn,
27 const G4Transform3D& offset,
28 G4double eScalingIn):
29 BDSFieldEInterpolated(interpolatorIn, offset, eScalingIn),
30 interpolator(interpolatorIn),
31 firstDimensionIndex((interpolatorIn->FirstDimension()).underlying()),
32 firstTime((interpolatorIn->FirstDimension()).underlying() > 2),
33 secondDimensionIndex((interpolatorIn->SecondDimension()).underlying()),
34 secondTime((interpolatorIn->SecondDimension()).underlying() > 2),
35 thirdDimensionIndex((interpolatorIn->ThirdDimension()).underlying()),
36 thirdTime((interpolatorIn->ThirdDimension()).underlying() > 2)
37{;}
38
39BDSFieldEInterpolated3D::~BDSFieldEInterpolated3D()
40{
41 delete interpolator;
42}
43
44G4ThreeVector BDSFieldEInterpolated3D::GetField(const G4ThreeVector& position,
45 const G4double t) const
46{
47 G4double fCoordinate = 0;
48 if (firstTime)
49 {fCoordinate = t;}
50 else
51 {fCoordinate = position[firstDimensionIndex];}
52 G4double sCoordinate = 0;
53 if (secondTime)
54 {sCoordinate = t;}
55 else
56 {sCoordinate = position[secondDimensionIndex];}
57 G4double tCoordinate = 0; // 't' for third
58 if (thirdTime)
59 {tCoordinate = t;}
60 else
61 {tCoordinate = position[thirdDimensionIndex];}
62 return interpolator->GetInterpolatedValue(fCoordinate, sCoordinate, tCoordinate) * EScaling();
63}
virtual G4ThreeVector GetField(const G4ThreeVector &position, const G4double t=0) const
Return the interpolated field value at a given point.
const G4bool secondTime
Cache of whether to use time coordinate.
const G4int thirdDimensionIndex
Integer index to dimension to use.
const G4int secondDimensionIndex
Integer index to dimension to use.
const G4bool firstTime
Cache of whether to use time coordinate.
const G4bool thirdTime
Cache of whether to use time coordinate.
BDSInterpolator3D * interpolator
Interpolator the field is based on.
const G4int firstDimensionIndex
Integer index to dimension to use.
Class to provide scaling and a base class pointer for interpolator fields.
Interface for all 3D interpolators.
G4ThreeVector GetInterpolatedValue(G4double x, G4double y, G4double z) const
Public interface to a 3D interpolator. Returns Geant4 type as that's what will be needed.