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