BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSFieldESinusoid.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 "BDSFieldESinusoid.hh"
20#include "BDSMagnetStrength.hh"
21#include "BDSUtilities.hh"
22
23#include "G4ThreeVector.hh"
24#include "G4Types.hh"
25
26#include "CLHEP/Units/SystemOfUnits.h"
27
28#include <cmath>
29
30BDSFieldESinusoid::BDSFieldESinusoid(BDSMagnetStrength const* strength,
31 G4double brho):
32 BDSFieldESinusoid((*strength)["efield"],
33 G4ThreeVector( (*strength)["ex"], (*strength)["ey"], (*strength)["ez"]),
34 (*strength)["frequency"],
35 (*strength)["phase"])
36{
37 G4int sign = BDS::Sign(brho);
38 eField *= sign;
39}
40
41BDSFieldESinusoid::BDSFieldESinusoid(G4double eFieldAmplitudeIn,
42 const G4ThreeVector& unitDirectionIn,
43 G4double frequencyIn,
44 G4double phaseOffsetIn):
45 eField(eFieldAmplitudeIn),
46 unitDirection(unitDirectionIn.unit()),
47 angularFrequency(CLHEP::twopi*frequencyIn),
48 phase(phaseOffsetIn)
49{
51}
52
53G4ThreeVector BDSFieldESinusoid::GetField(const G4ThreeVector& /*position*/,
54 const G4double t) const
55{
56 G4double ef = eField * std::cos(angularFrequency*t + phase);
57 G4ThreeVector field = unitDirection * ef;
58 return field;
59}
A sinusoidal electric (only) field that doesn't vary with position. Uses cosine.
G4double angularFrequency
Angular frequency of field.
virtual G4ThreeVector GetField(const G4ThreeVector &position, const G4double t) const
Accessor for field value.
G4double eField
Amplitude of electric field in V/m.
const G4ThreeVector unitDirection
Unit vector for direction of field.
G4double phase
Phase in radians.
G4bool finiteStrength
Flag to cache whether finite nor not.
Definition: BDSFieldE.hh:83
Efficient storage of magnet strengths.
G4bool IsFinite(G4double value, G4double tolerance=std::numeric_limits< double >::epsilon())