BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSBeamPipeFactoryPoints.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 "BDSBeamPipeFactoryPoints.hh"
20#include "BDSBeamPipe.hh"
21#include "BDSExtent.hh"
22#include "BDSUtilities.hh"
23
24#include "globals.hh" // geant4 globals / types
25#include "G4CutTubs.hh"
26#include "G4ExtrudedSolid.hh"
27#include "G4IntersectionSolid.hh"
28#include "G4SubtractionSolid.hh"
29#include "G4ThreeVector.hh"
30#include "G4TwoVector.hh"
31#include "G4VSolid.hh"
32
33#include <cmath>
34#include <vector>
35
36BDSBeamPipeFactoryPoints::BDSBeamPipeFactoryPoints()
37{
39}
40
41BDSBeamPipeFactoryPoints::~BDSBeamPipeFactoryPoints()
42{;}
43
45{
47 BDSBeamPipeFactoryBase::CleanUp();
48}
49
51{
53 extentX = 0;
54 extentY = 0;
55
56 vacuumEdge.clear();
57 beamPipeInnerEdge.clear();
58 beamPipeOuterEdge.clear();
59 containerEdge.clear();
61
62 beamPipeInnerSolid = nullptr;
63 beamPipeOuterSolid = nullptr;
64
65 pointsFile = "";
66 pointsUnit = "";
67}
68
69void BDSBeamPipeFactoryPoints::AppendPoint(std::vector<G4TwoVector>& vec,
70 G4double x,
71 G4double y)
72{
73 vec.emplace_back(G4TwoVector(x,y));
74}
75
76void BDSBeamPipeFactoryPoints::AppendAngle(std::vector<G4TwoVector>& vec,
77 G4double startAngle,
78 G4double finishAngle,
79 G4double radius,
80 G4int nPoints,
81 G4double xOffset,
82 G4double yOffset)
83{
84 AppendAngleEllipse(vec, startAngle, finishAngle, radius, radius,nPoints, xOffset, yOffset);
85}
86
87void BDSBeamPipeFactoryPoints::AppendAngleEllipse(std::vector<G4TwoVector>& vec,
88 G4double startAngle,
89 G4double finishAngle,
90 G4double radiusA,
91 G4double radiusB,
92 G4int nPoints,
93 G4double xOffset,
94 G4double yOffset)
95{
96 G4double diff = finishAngle - startAngle;
97 G4double delta = diff / (G4double)nPoints;
98 G4double ang = startAngle;
99 for (G4int i = 0; i < nPoints; i++)
100 { // l for local
101 G4double xl = xOffset + radiusA*std::sin(ang);
102 G4double yl = yOffset + radiusB*std::cos(ang);
103 AppendPoint(vec, xl, yl);
104 ang += delta;
105 }
106}
107
109 G4double length,
110 G4bool buildLongForIntersection)
111{
113 G4double zHalfLength = 0.5*length;
114 if (buildLongForIntersection)
115 {zHalfLength *= 1.5;}
116
117 G4TwoVector zOffsets(0,0); // the transverse offset of each plane from 0,0
118 G4double zScale = 1; // the scale at each end of the points = 1
119 vacuumSolid = new G4ExtrudedSolid(name + "_vacuum_solid", // name
120 vacuumEdge, // vector of TwoVector points
121 zHalfLength - lengthSafety, // half length for +- planes
122 zOffsets, zScale, // dx,dy offset for each face, scaling
123 zOffsets, zScale); // dx,dy offset for each face, scaling
124
125 beamPipeInnerSolid = new G4ExtrudedSolid(name + "_bp_inner_solid",
127 zHalfLength * 4,
128 zOffsets, zScale,
129 zOffsets, zScale);
130
131 beamPipeOuterSolid = new G4ExtrudedSolid(name + "_bp_outer_solid",
133 zHalfLength - 4*lengthSafety,
134 zOffsets, zScale,
135 zOffsets, zScale);
136
137 allSolids.insert(beamPipeInnerSolid);
138 allSolids.insert(beamPipeOuterSolid);
139
140 beamPipeSolid = new G4SubtractionSolid(name + "_pipe_solid", // name
141 beamPipeOuterSolid, // this
142 beamPipeInnerSolid); // minus this
143
144 containerSolid = new G4ExtrudedSolid(name + "_container_solid",
146 zHalfLength,
147 zOffsets, zScale,
148 zOffsets, zScale);
149
150 containerSubtractionSolid = new G4ExtrudedSolid(name + "_container_subtraction_solid",
152 zHalfLength*4,
153 zOffsets, zScale,
154 zOffsets, zScale);
155}
156
158 G4double length,
159 const G4ThreeVector& inputFace,
160 const G4ThreeVector& outputFace)
161{
162 // long length for unambiguous boolean - ensure no gaps in beam pipe geometry
163 // extra factor 2 to be safe
164 G4double angledVolumeLength = BDS::CalculateSafeAngledVolumeLength(inputFace, outputFace, length*2, intersectionRadius);
165
166 // create straight solids that are a bit long
167 CreateSolids(name + "_straight", angledVolumeLength, true);
168
169 // now intersect them with one G4CutTubs to get the angled faces
170 G4double zHalfLength = length*0.5 - lengthSafety;
171 G4double zHalfLengthContainer = length*0.5;
172
173 // check faces of angled volume don't intersect - if it can be built, remaining angled volumes can be built
174 CheckAngledVolumeCanBeBuilt(length, inputFace, outputFace, intersectionRadius, name);
175
176 G4VSolid* faceSolid = new G4CutTubs(name + "_face_solid", // name
177 0, // inner radius
178 intersectionRadius, // outer radius
179 zHalfLength, // z half length
180 0, // start angle
181 CLHEP::twopi, // sweep angle
182 inputFace, // input face normal
183 outputFace); // output face normal
184
185 // different only by container length > length of inner section
186 G4VSolid* faceSolidContainer = new G4CutTubs(name + "_cont_face_solid", // name
187 0, // inner radius
188 intersectionRadius, // outer radius
189 zHalfLengthContainer, // z half length
190 0, // start angle
191 CLHEP::twopi, // sweep angle
192 inputFace, // input face normal
193 outputFace); // output face normal
194
195 allSolids.insert(faceSolid);
196 allSolids.insert(faceSolidContainer);
197
198 // copy pointers to do intersection with then reassign member pointer
199 // to point to new solid
200 G4VSolid* vacuumTemp = vacuumSolid;
201 vacuumSolid = new G4IntersectionSolid(name + "_vacuum_solid",
202 vacuumTemp,
203 faceSolid);
204
205 G4VSolid* beamPipeTemp = beamPipeSolid;
206 beamPipeSolid = new G4IntersectionSolid(name + "_pipe_solid",
207 beamPipeTemp,
208 faceSolid);
209
210 G4VSolid* containerTemp = containerSolid;
211 containerSolid = new G4IntersectionSolid(name + "_container_solid",
212 containerTemp,
213 faceSolidContainer);
214
215 // container subtraction solid can just be long with flat edges as
216 // only used transversely
217}
218
220 G4double lengthIn,
221 G4double aper1In,
222 G4double aper2In,
223 G4double aper3In,
224 G4double aper4In,
225 G4Material* vacuumMaterialIn,
226 G4double beamPipeThicknessIn,
227 G4Material* beamPipeMaterialIn,
228 const G4String& pointsFileIn,
229 const G4String& pointsUnitIn)
230{
231 // clean up after last usage
232 CleanUp();
233
234 pointsFile = pointsFileIn;
235 pointsUnit = pointsUnitIn;
236
237 // generate extruded solid edges - provided by derived class
238 GeneratePoints(aper1In, aper2In, aper3In, aper4In, beamPipeThicknessIn);
239
240 // calculate and set the intersection solid radius
241 intersectionRadius = CalculateIntersectionRadius(aper1In, aper2In, aper3In, aper4In, beamPipeThicknessIn);
242
243 // create solids based on the member vectors of points
244 CreateSolids(nameIn, lengthIn);
245
246 return CommonFinalConstruction(nameIn, vacuumMaterialIn, beamPipeMaterialIn, lengthIn);
247}
248
250 G4double lengthIn,
251 const G4ThreeVector& inputFaceNormalIn,
252 const G4ThreeVector& outputFaceNormalIn,
253 G4double aper1In,
254 G4double aper2In,
255 G4double aper3In,
256 G4double aper4In,
257 G4Material* vacuumMaterialIn,
258 G4double beamPipeThicknessIn,
259 G4Material* beamPipeMaterialIn,
260 const G4String& pointsFileIn,
261 const G4String& pointsUnitIn)
262{
263 // clean up after last usage
264 CleanUp();
265
266 pointsFile = pointsFileIn;
267 pointsUnit = pointsUnitIn;
268
269 // generate extruded solid edges - provided by derived class
270 GeneratePoints(aper1In, aper2In, aper3In, aper4In, beamPipeThicknessIn);
271
272 inputFaceNormal = inputFaceNormalIn;
273 outputFaceNormal = outputFaceNormalIn;
274
275 // calculate and set the intersection solid radius
276 intersectionRadius = CalculateIntersectionRadius(aper1In, aper2In, aper3In, aper4In, beamPipeThicknessIn);
277
278 // create solids based on the member vectors of points
280
281 return CommonFinalConstruction(nameIn, vacuumMaterialIn, beamPipeMaterialIn, lengthIn);
282}
283
285 G4Material* vacuumMaterialIn,
286 G4Material* beamPipeMaterialIn,
287 G4double lengthIn)
288{
289 BDSBeamPipeFactoryBase::CommonConstruction(nameIn, vacuumMaterialIn,
290 beamPipeMaterialIn, lengthIn);
291
292 // record extents
293 BDSExtent ext = BDSExtent(extentX, extentY, lengthIn*0.5);
294
295 // calculate radius if a tube were to be place around it
296 G4double containerRadius = intersectionRadius;
297
298 BDSBeamPipe* aPipe = BuildBeamPipeAndRegisterVolumes(ext,containerRadius);
299
300 return aPipe;
301}
G4ThreeVector outputFaceNormal
For recording the face normals in the finished pipe component.
G4ThreeVector inputFaceNormal
For recording the face normals in the finished pipe component.
void CommonConstruction(const G4String &nameIn, G4Material *vacuumMaterialIn, G4Material *beamPipeMaterialIn, G4double length)
finalise beampipe construction
BDSBeamPipe * BuildBeamPipeAndRegisterVolumes(BDSExtent extent, G4double containerRadius, G4bool containerIsCircular=false)
build beampipe and register logical volumes
static void CheckAngledVolumeCanBeBuilt(G4double length, const G4ThreeVector &inputfaceAngle, const G4ThreeVector &outputfaceAngle, G4double horizontalWidth, const G4String &name)
check if a beam pipe volume with angled faces can be constructed
G4VSolid * containerSubtractionSolid
Longer (in length) version of container solid for unambiguous subtraction.
G4VSolid * beamPipeInnerSolid
Solid for inner surface of beam pipe.
std::vector< G4TwoVector > containerSubtractionEdge
Vector of x,y coordinates for the container subtraction volume.
std::vector< G4TwoVector > beamPipeInnerEdge
Vector of x,y coordinates for beam pipe inner edge.
virtual void CleanUp()
Clear member vectors and run base class clean up to clear pointers between runs.
G4VSolid * beamPipeOuterSolid
Solid for outer surface of beam pipe.
virtual void GeneratePoints(G4double aper1, G4double aper2, G4double aper3, G4double aper4, G4double beamPipeThickness, G4int pointsPerTwoPi=40)=0
G4double extentY
Extents prepared by GeneratePoints function as only it knows the extents.
virtual G4double CalculateIntersectionRadius(G4double aper1, G4double aper2, G4double aper3, G4double aper4, G4double beamPipeThickness)=0
void CleanUpPoints()
Actual clean up here for this class only.
virtual BDSBeamPipe * CommonFinalConstruction(const G4String &nameIn, G4Material *vacuumMaterialIn, G4Material *beamPipeMaterialIn, G4double lengthIn)
G4double extentX
Extents prepared by GeneratePoints function as only it knows the extents.
std::vector< G4TwoVector > vacuumEdge
Vector of x,y coordinates for vacuum extruded solid edge.
virtual BDSBeamPipe * CreateBeamPipe(const G4String &nameIn, G4double lengthIn, G4double aper1In, G4double aper2In, G4double aper3In, G4double aper4In, G4Material *vacuumMaterialIn, G4double beamPipeThicknessIn, G4Material *beamPipeMaterialIn, const G4String &pointsFileIn="", const G4String &pointsUnitIn="")
Required overloaded method from BDSBeamPipeFactoryBase to build straight piece of beam pipe.
void CreateSolidsAngled(const G4String &name, G4double length, const G4ThreeVector &inputFace, const G4ThreeVector &outputFace)
Create angled solids. Uses CreateSolids() method.
void AppendAngle(std::vector< G4TwoVector > &vec, G4double startAngle, G4double finishAngle, G4double radius, G4int nPoints=10, G4double xOffset=0, G4double yOffset=0)
void CreateSolids(const G4String &name, G4double length, G4bool buildLongForIntersection=false)
G4String pointsFile
Temporary copy of input variables.
std::vector< G4TwoVector > containerEdge
Vector of x,y coordinates for the container volume.
std::vector< G4TwoVector > beamPipeOuterEdge
Vector of x,y coordinates for beam pipe outer edge.
G4String pointsUnit
Temporary copy of input variables.
void AppendAngleEllipse(std::vector< G4TwoVector > &vec, G4double startAngle, G4double finishAngle, G4double radiusA, G4double radiusB, G4int nPoints=10, G4double xOffset=0, G4double yOffset=0)
Generate 2-vector points (and append them) about an ellipse.
A holder class for a piece of beam pipe geometry.
Definition: BDSBeamPipe.hh:45
Holder for +- extents in 3 dimensions.
Definition: BDSExtent.hh:39
G4double lengthSafety
Cache of global constants variable.
G4double CalculateSafeAngledVolumeLength(G4double angleIn, G4double angleOut, G4double length, G4double containerWidth, G4double containerHeight=0)
Calculate safe length of an angled volume so it fills the length of its container.