BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSBunchHaloFlatSigma.hh
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#ifndef BDSBUNCHHALOFLATSIGMA_H
20#define BDSBUNCHHALOFLATSIGMA_H
21#include "BDSBunch.hh"
22#include "BDSBunchType.hh"
23#include "BDSException.hh"
24
25#include "G4Transform3D.hh"
26#include "G4Types.hh"
27
28#include <algorithm>
29#include <cmath>
30#include <iterator>
31#include <memory>
32#include <numeric>
33#include <vector>
34
35namespace GMAD
36{
37 class Beam;
38}
39
40namespace BDS
41{
50{
51 G4double action;
52 G4double angle;
53};
54
63{
64 G4double position;
65 G4double momentum;
66};
67
76{
77 G4double alpha;
78 G4double beta;
79};
80
81PhaseSpaceCoord PhaseSpaceCoordFromActionAngle(ActionAngleCoord aa, const TwissPair& tp);
82
83template <typename T>
84std::vector<G4double> CumulativeDistances(T pairs)
85{
86 if (pairs.size() < 2)
87 {throw BDSException("At least two points needed to calculated adjacent distances.");}
88 std::vector<G4double> distances;
89
90 auto it = std::begin(pairs);
91 auto next = std::next(it);
92 for (; next != std::end(pairs); ++it, ++next)
93 {
94 G4double x0 = it->position;
95 G4double xp0 = it->momentum;
96
97 G4double x1 = next->position;
98 G4double xp1 = next->momentum;
99
100 G4double xd = std::pow((x0 - x1), 2);
101 G4double xpd = std::pow((xp0 - xp1), 2);
102
103 G4double distance = std::sqrt(xd + xpd);
104
105 distances.push_back(distance);
106 }
107 std::vector<G4double> result;
108 std::partial_sum(std::begin(distances), std::end(distances), std::back_inserter(result));
109 return result;
110}
111
120{
121public:
122 EllipsePointGenerator(G4double actionIn,
123 const TwissPair& tp);
124 ~EllipsePointGenerator() = default;
125
126 PhaseSpaceCoord GetRandomPointOnEllipse() const;
127 inline double EllipsePerimeter() const {return pathLengths.back();};
128
129private:
130 G4double PathLengthToAngle(G4double pathLength) const;
131
132 G4double action;
133 TwissPair twisspair;
134 std::vector<G4double> angles;
135 std::vector<G4double> pathLengths;
136};
137
138} // namespace BDS
139
147{
148public:
150 virtual ~BDSBunchHaloFlatSigma(){};
151 virtual void SetOptions(const BDSParticleDefinition* beamParticle,
152 const GMAD::Beam& beam,
153 const BDSBunchType& distrType,
154 G4Transform3D beamlineTransformIn = G4Transform3D::Identity,
155 const G4double beamlineS = 0);
156 virtual void CheckParameters();
158
159private:
161 G4double alphaX;
162 G4double alphaY;
163 G4double betaX;
164 G4double betaY;
165 G4double emitX;
166 G4double emitY;
167 G4double sigmaX;
168 G4double sigmaY;
170
171 G4double haloNSigmaXInner;
172 G4double haloNSigmaXOuter;
173 G4double haloNSigmaYInner;
174 G4double haloNSigmaYOuter;
175};
176
177#endif
Bunch halo distribution where the PDF is uniformly distribution in sigma.
G4double alphaX
Twiss parameter.
G4double sigmaY
Twiss parameter.
G4double betaY
Twiss parameter.
virtual void SetOptions(const BDSParticleDefinition *beamParticle, const GMAD::Beam &beam, const BDSBunchType &distrType, G4Transform3D beamlineTransformIn=G4Transform3D::Identity, const G4double beamlineS=0)
G4double sigmaX
Twiss parameter.
G4double emitX
Twiss parameter.
G4double betaX
Twiss parameter.
G4double emitY
Twiss parameter.
virtual BDSParticleCoordsFull GetNextParticleLocal()
G4double alphaY
Twiss parameter.
The base class for bunch distribution generators.
Definition: BDSBunch.hh:47
G4double beamlineS
Beamline initial S position.
Definition: BDSBunch.hh:219
General exception with possible name of object and message.
Definition: BDSException.hh:35
A set of particle coordinates including energy and weight.
Wrapper for particle definition.
Improve type-safety of native enum data type in C++.
Class for generating points uniformly on ellipse perimeters via interpolation. Part of implementation...
Beam loader.
Definition: Beam.hh:37
Beam class.
Definition: beam.h:44
Return either G4Tubs or G4CutTubs depending on flat face.
Parser namespace for GMAD language. Combination of Geant4 and MAD.
Simple struct for storing action/angle pairs to aid readability. Implementation detail.
Simple struct for storing position/momentum pairs to aid readability. Implementation detail.
Simple struct for storing Twiss alpha/beta pairs to aid readability. Implementation detail.