BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
BDSBeamPipeFactoryRectEllipse.cc
1/*
2Beam Delivery Simulation (BDSIM) Copyright (C) Royal Holloway,
3University of London 2001 - 2022.
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 "BDSBeamPipeFactoryBase.hh"
20#include "BDSBeamPipeFactoryElliptical.hh"
21#include "BDSBeamPipeFactoryRectEllipse.hh"
22#include "BDSBeamPipe.hh"
23#include "BDSDebug.hh"
24#include "BDSException.hh"
25#include "BDSExtent.hh"
26#include "BDSUtilities.hh"
27
28#include "globals.hh"
29#include "G4Box.hh"
30#include "G4CutTubs.hh"
31#include "G4EllipticalTube.hh"
32#include "G4IntersectionSolid.hh"
33#include "G4LogicalVolume.hh"
34#include "G4SubtractionSolid.hh"
35#include "G4ThreeVector.hh"
36#include "G4VSolid.hh"
37
38#include "CLHEP/Units/SystemOfUnits.h"
39
40#include <cmath>
41#include <set>
42#include <string>
43#include <utility> // for std::pair
44
45BDSBeamPipeFactoryRectEllipse::BDSBeamPipeFactoryRectEllipse()
46{;}
47
49 G4double lengthIn,
50 G4double aper1In,
51 G4double aper2In,
52 G4double aper3In,
53 G4double aper4In,
54 G4Material* vacuumMaterialIn,
55 G4double beamPipeThicknessIn,
56 G4Material* beamPipeMaterialIn,
57 const G4String& /*pointsFileIn*/,
58 const G4String& /*pointsUnitIn*/)
59{
60 // clean up after last usage
61 CleanUp();
62
63 // CERN use rectellipse for everything because you know, why use different aperture types
64 // tolerate the case where it's a) equal rect and ellipse parameters resulting in an ellipse
65 // and b) the ellipse params are less than rect ones so also just an ellipse
66 if (((aper1In == aper3In) && (aper2In == aper4In)) ||
67 ( (aper3In < aper1In) && (aper4In < aper2In) ) )
68 {
70 BDSBeamPipe* result = ef->CreateBeamPipe(nameIn, lengthIn, aper1In, aper2In, aper3In, aper4In,
71 vacuumMaterialIn, beamPipeThicknessIn, beamPipeMaterialIn);
72 delete ef;
73 return result;
74 }
75
76 // build the solids
77 //vacuum cylindrical solid (circular cross-section)
78 G4VSolid* vacCylSolid = new G4EllipticalTube(nameIn + "_vacuum_ellipsoid", // name
79 aper3In, // horizontal semi-axis
80 aper4In, // vertical semi-axis
81 lengthIn*0.5-lengthSafety); // half length
82 //vacuum box solid (rectangular cross-section)
83 G4VSolid* vacRectSolid = new G4Box(nameIn + "_vacuum_box", // name
84 aper1In, // x half width
85 aper2In, // y half width
86 4*lengthIn); // z full width (long for unambiguous intersection)
87 allSolids.insert(vacCylSolid);
88 allSolids.insert(vacRectSolid);
89 //intersection of both of these gives the desired shape
90 vacuumSolid = new G4IntersectionSolid(nameIn + "_vacuum_solid", // name
91 vacCylSolid, // solid 1
92 vacRectSolid); // solid 2
93
94 //beampipe solid
95 //beampipe inner edge for subtraction (actually just like vacuum + lengthSafety)
96 G4VSolid* bpInnerCylSolid = new G4EllipticalTube(nameIn + "_pipe_inner_ellipsoid",// name
97 aper3In + lengthSafetyLarge, // horizontal semi-axis
98 aper4In + lengthSafetyLarge, // vertical semi-axis
99 1.5*lengthIn); // length big for unambiguous subtraction (but < outerlength)
100 //beampipe inner edge box solid (rectangular cross-section)
101 G4VSolid* bpInnerRectSolid = new G4Box(nameIn + "_pipe_inner_box", // name
102 aper1In + lengthSafetyLarge,// x half width
103 aper2In + lengthSafetyLarge,// y half width
104 1.7*lengthIn); // z long for unambiguous intersection
105 //beampipe inner intersection - 1.5*length long which is > half length for unambiguous subtraction later
106 G4VSolid* bpInnerSolid = new G4IntersectionSolid(nameIn + "_pipe_inner_solid", // name
107 bpInnerCylSolid, // solid 1
108 bpInnerRectSolid); // solid 2
109
110 //beampipe outer edge for subtraction (actually just like vacuum + lengthSafety)
111 G4double extraWidth = beamPipeThicknessIn + lengthSafetyLarge;
112 G4VSolid* bpOuterCylSolid = new G4EllipticalTube(nameIn + "_pipe_inner_ellipsoid",// name
113 aper3In + extraWidth, // horizontal semi-axis
114 aper4In + extraWidth, // hotizontal semi-axis
115 (lengthIn*0.5)-lengthSafety); // half length
116 //beampipe outer edge box solid (rectangular cross-section)
117 G4VSolid* bpOuterRectSolid = new G4Box(nameIn + "_pipe_inner_box", // name
118 aper1In + extraWidth, // x half width
119 aper2In + extraWidth, // y half width
120 lengthIn); // z full width (long for unambiguous intersection)
121 G4VSolid* bpOuterSolid = new G4IntersectionSolid(nameIn + "_pipe_inner_solid", // name
122 bpOuterCylSolid, // solid 1
123 bpOuterRectSolid); // solid 2
124
125 allSolids.insert(bpInnerCylSolid);
126 allSolids.insert(bpInnerRectSolid);
127 allSolids.insert(bpInnerSolid);
128 allSolids.insert(bpOuterCylSolid);
129 allSolids.insert(bpOuterRectSolid);
130 allSolids.insert(bpOuterSolid);
131
132 //beampipe final subtraction between outer and inner edge
133 beamPipeSolid = new G4SubtractionSolid(nameIn + "_pipe_solid", // name
134 bpOuterSolid, // this
135 bpInnerSolid); // minus this
136
137 //container cylindrical solid (circular cross-section)
138 G4VSolid* contCylSolid = new G4EllipticalTube(nameIn + "_vacuum_ellipsoid", // name
139 aper3In + extraWidth + lengthSafetyLarge, // horizontal semi-axis
140 aper4In + extraWidth + lengthSafetyLarge, // vertical semi-axis
141 lengthIn*0.5); // half length
142 //vacuum box solid (rectangular cross-section)
143 G4VSolid* contRectSolid = new G4Box(nameIn + "_vacuum_box", // name
144 aper1In + extraWidth + lengthSafetyLarge, // x half width
145 aper2In + extraWidth + lengthSafetyLarge, // y half width
146 lengthIn); // z full width (long for unambiguous intersection)
147
148 allSolids.insert(contCylSolid);
149 allSolids.insert(contRectSolid);
150
151 //intersection of both of these gives the desired shape
152 containerSolid = new G4IntersectionSolid(nameIn + "_vacuum_solid", // name
153 contCylSolid, // solid 1
154 contRectSolid); // solid 2
155
156 G4double width = aper3In + extraWidth + lengthSafetyLarge;
157 G4double height = aper2In + extraWidth + lengthSafetyLarge;
158
159 CreateContainerSubtractionSolid(nameIn, lengthIn, beamPipeThicknessIn, aper1In, aper2In, aper3In, aper4In);
160
161 return CommonFinalConstruction(nameIn, vacuumMaterialIn, beamPipeMaterialIn, lengthIn, width, height);
162}
163
165 G4double lengthIn,
166 const G4ThreeVector& inputFaceNormalIn,
167 const G4ThreeVector& outputFaceNormalIn,
168 G4double aper1In,
169 G4double aper2In,
170 G4double aper3In,
171 G4double aper4In,
172 G4Material* vacuumMaterialIn,
173 G4double beamPipeThicknessIn,
174 G4Material* beamPipeMaterialIn,
175 const G4String& /*pointsFileIn*/,
176 const G4String& /*pointsUnitIn*/)
177{
178 // clean up after last usage
179 CleanUp();
180
181 // CERN use rectellipse for everything because you know, why use different aperture types
182 // tolerate the case where it's a) equal rect and ellipse parameters resulting in an ellipse
183 // and b) the ellipse params are less than rect ones so also just an ellipse
184 if (((aper1In == aper3In) && (aper2In == aper4In)) ||
185 ( (aper3In < aper1In) && (aper4In < aper2In) ) )
186 {
188 BDSBeamPipe* result = ef->CreateBeamPipe(nameIn, lengthIn, inputFaceNormalIn, outputFaceNormalIn,
189 aper1In, aper2In, aper3In, aper4In,
190 vacuumMaterialIn, beamPipeThicknessIn, beamPipeMaterialIn);
191 delete ef;
192 return result;
193 }
194
195 inputFaceNormal = inputFaceNormalIn;
196 outputFaceNormal = outputFaceNormalIn;
197
198 G4double width = std::max(aper1In,aper3In) + beamPipeThicknessIn + lengthSafetyLarge;
199 G4double height = std::max(aper2In,aper4In) + beamPipeThicknessIn + lengthSafetyLarge;
200
201 CreateGeneralAngledSolids(nameIn, lengthIn, aper1In, aper2In, aper3In, aper4In,
202 beamPipeThicknessIn, inputFaceNormal, outputFaceNormal);
203 CreateContainerSubtractionSolid(nameIn, lengthIn, beamPipeThicknessIn, aper1In,
204 aper2In, aper3In, aper4In);
205
206 return CommonFinalConstruction(nameIn, vacuumMaterialIn, beamPipeMaterialIn,
207 lengthIn, width, height);
208}
209
211 G4Material* vacuumMaterialIn,
212 G4Material* beamPipeMaterialIn,
213 G4double lengthIn,
214 G4double containerWidthIn,
215 G4double containerHeightIn)
216{
217 BDSBeamPipeFactoryBase::CommonConstruction(nameIn, vacuumMaterialIn,
218 beamPipeMaterialIn, lengthIn);
219
220 // record extents
221 BDSExtent ext = BDSExtent(containerWidthIn, containerHeightIn, lengthIn*0.5);
222
223 // build the BDSBeamPipe instance and return it
224 return BuildBeamPipeAndRegisterVolumes(ext,containerWidthIn);
225}
226
230 G4double lengthIn,
231 G4double aper1In,
232 G4double aper2In,
233 G4double aper3In,
234 G4double aper4In,
235 G4double beamPipeThicknessIn,
236 const G4ThreeVector& inputfaceIn,
237 const G4ThreeVector& outputfaceIn)
238{
239 // long length for unambiguous boolean - ensure no gaps in beam pipe geometry
240 G4double angledVolumeLength = BDS::CalculateSafeAngledVolumeLength(inputfaceIn, outputfaceIn, lengthIn, aper1In);
241
242 // build the solids
243 // we can get the rectangular ellipse as in the straight case with the intersection of
244 // an elliptical tube (always solid) and then we can use another intersection solid
245 // with a larger (wider and taller) G4CutTubs to get the angled faces.
246 //vacuum cylindrical solid (circular cross-section)
247 G4VSolid* vacCylSolid = new G4EllipticalTube(nameIn + "_vacuum_ellipsoid",// name
248 aper3In, // horizontal semi-axis
249 aper4In, // vertical semi-axis
250 angledVolumeLength); // long length for unambiguous intersection
251 //vacuum box solid (rectangular cross-section)
252 G4VSolid* vacRectSolid = new G4Box(nameIn + "_vacuum_box", // name
253 aper1In, // x half width
254 aper2In, // y half width
255 angledVolumeLength); // long length for unambiguous intersection
256 //intersection of both of these gives the desired shape
257 G4VSolid* longVacuumSolid = new G4IntersectionSolid(nameIn + "_vacuum_long_solid", // name
258 vacCylSolid, // solid 1
259 vacRectSolid); // solid 2
260 //prepare angled face large cylinder for intersection get angled faces
261 //we can actually use this for the beampipe too later on - whew
262 G4double angledFaceRadius = (std::max(std::max(aper1In,aper2In),std::max(aper3In,aper4In)) + beamPipeThicknessIn) * 1.1;
263
264 // check faces of angled volume don't intersect - if it can be built, remaining angled volumes can be built
265 CheckAngledVolumeCanBeBuilt(lengthIn, inputfaceIn, outputfaceIn, angledFaceRadius, nameIn);
266
267 G4VSolid* vacuumAngledSolid = new G4CutTubs(nameIn + "_pipe_angled_faces", // name
268 0, // inner radius
269 angledFaceRadius, // outer radius
270 (lengthIn*0.5) - (lengthSafety), // accurate length
271 0, // rotation start angle
272 CLHEP::twopi, // rotation sweep angle
273 inputfaceIn, // input face normal
274 outputfaceIn); // output face normal
275
276 allSolids.insert(vacCylSolid);
277 allSolids.insert(vacRectSolid);
278 allSolids.insert(longVacuumSolid);
279 allSolids.insert(vacuumAngledSolid);
280
281 vacuumSolid = new G4IntersectionSolid(nameIn + "_vacuum_solid", // name
282 longVacuumSolid, // solid 1
283 vacuumAngledSolid); // solid 2
284
285 //beampipe cylindrical solid (circular cross-section)
286 //beampipe inner edge for subtraction (actually just like vacuum + lengthSafety)
287 if (beamPipeThicknessIn <= lengthSafetyLarge)
288 {
289 G4String message = "beam pipe thickness in element \"" + nameIn
290 + "\" is thinner than minimum " + std::to_string(lengthSafetyLarge / CLHEP::mm) + "mm.";
291 throw BDSException(__METHOD_NAME__, message);
292 }
293 G4VSolid* bpInnerCylSolid = new G4EllipticalTube(nameIn + "_pipe_inner_ellipsoid", // name
294 aper3In + lengthSafetyLarge, // horizontal semi-axis
295 aper4In + lengthSafetyLarge, // vertical semi-axis
296 angledVolumeLength); // length big for unambiguous subtraction (but < outerlength)
297 //beampipe inner edge box solid (rectangular cross-section)
298 G4VSolid* bpInnerRectSolid = new G4Box(nameIn + "_pipe_inner_box", // name
299 aper1In + lengthSafetyLarge, // x half width
300 aper2In + lengthSafetyLarge, // y half width
301 angledVolumeLength); // long length for unambiguous intersection
302 //beampipe inner intersection - 1.5*length long which is > half length for unambiguous subtraction later
303 G4VSolid* bpInnerSolid = new G4IntersectionSolid(nameIn + "_pipe_inner_solid", // name
304 bpInnerCylSolid, // solid 1
305 bpInnerRectSolid); // solid 2
306
307 //beampipe outer edge for subtraction (actually just like vacuum + lengthSafety)
308 //this length should be less than bpInnerSolid above but longer than the actual length for later intersection
309 G4double extraWidth = beamPipeThicknessIn + lengthSafetyLarge;
310 G4VSolid* bpOuterCylSolid = new G4EllipticalTube(nameIn + "_pipe_inner_ellipsoid", // name
311 aper3In + extraWidth, // horizontal semi-axis
312 aper4In + extraWidth, // vertical semi-axis
313 lengthIn); // length
314 //beampipe outer edge box solid (rectangular cross-section)
315 G4VSolid* bpOuterRectSolid = new G4Box(nameIn + "_pipe_inner_box", // name
316 aper1In + extraWidth, // x half width
317 aper2In + extraWidth, // y half width
318 angledVolumeLength); // long length for unambiguous intersection
319 G4VSolid* bpOuterSolid = new G4IntersectionSolid(nameIn + "_pipe_inner_solid", // name
320 bpOuterCylSolid, // solid 1
321 bpOuterRectSolid); // solid 2
322 //correct solid shape but too long
323 G4VSolid* longBeamPipeSolid = new G4SubtractionSolid(nameIn + "_long_pipe_solid", // name
324 bpOuterSolid, // this
325 bpInnerSolid); // minus this
326
327 allSolids.insert(bpInnerCylSolid);
328 allSolids.insert(bpInnerRectSolid);
329 allSolids.insert(bpInnerSolid);
330 allSolids.insert(bpOuterCylSolid);
331 allSolids.insert(bpOuterRectSolid);
332 allSolids.insert(bpOuterSolid);
333
334 //final beampipe solid with correct shape and angled faces
335 beamPipeSolid = new G4IntersectionSolid(nameIn + "_pipe_solid", // name
336 longBeamPipeSolid, // solid1
337 vacuumAngledSolid); // solid2
338
339 //container solid
340 //container cylindrical solid (circular cross-section)
341 G4VSolid* contCylSolid = new G4EllipticalTube(nameIn + "_vacuum_ellipsoid", // name
342 aper3In + extraWidth + lengthSafetyLarge, // horizontal semi-axis
343 aper4In + extraWidth + lengthSafetyLarge, // vertical semi-axis
344 angledVolumeLength); // length
345 // length both long and different from next solid for unamibiguous intersection
346
347 //vacuum box solid (rectangular cross-section)
348 G4VSolid* contRectSolid = new G4Box(nameIn + "_vacuum_box", // name
349 aper1In + extraWidth + lengthSafetyLarge, // x half width
350 aper2In + extraWidth + lengthSafetyLarge, // y half width
351 lengthIn); // z full width (long for unambiguous intersection)
352 //intersection of both of these gives the desired shape
353 G4VSolid* longContainerSolid = new G4IntersectionSolid(nameIn + "_long_container_solid", // name
354 contCylSolid, // solid 1
355 contRectSolid); // solid 2
356
357 allSolids.insert(contCylSolid);
358 allSolids.insert(contRectSolid);
359 allSolids.insert(longContainerSolid);
360
361 containerSolid = new G4IntersectionSolid(nameIn + "_container_solid", // name
362 longContainerSolid, // solid 1
363 vacuumAngledSolid); // solid 2
364}
365
366void BDSBeamPipeFactoryRectEllipse::CreateContainerSubtractionSolid(const G4String& nameIn,
367 G4double& lengthIn,
368 G4double& beamPipeThicknessIn,
369 G4double& aper1In,
370 G4double& aper2In,
371 G4double& aper3In,
372 G4double& aper4In)
373{
374 //container cylindrical solid (circular cross-section)
375 G4VSolid* contSubCylSolid = new G4EllipticalTube(nameIn + "_subtraction_ellipsoid", // name
376 aper3In + beamPipeThicknessIn + 3*lengthSafetyLarge, // horizontal semi-axis
377 aper4In + beamPipeThicknessIn + 3*lengthSafetyLarge, // vertical semi-axis
378 2*lengthIn); // long length for unambiguous subtraction
379 //vacuum box solid (rectangular cross-section)
380 G4VSolid* contSubRectSolid = new G4Box(nameIn + "_subtraction_box", // name
381 aper1In + beamPipeThicknessIn + 2*lengthSafetyLarge, // x half width
382 aper2In + beamPipeThicknessIn + 2*lengthSafetyLarge, // y half width
383 1.7*lengthIn); // z full width (long for unambiguous intersection)
384 allSolids.insert(contSubCylSolid);
385 allSolids.insert(contSubRectSolid);
386
387 //intersection of both of these gives the desired shape
388 containerSubtractionSolid = new G4IntersectionSolid(nameIn + "_subtraction_solid", // name
389 contSubCylSolid, // solid 1
390 contSubRectSolid); // solid 2
391}
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.
Factory for elliptical beam pipes.
virtual BDSBeamPipe * CreateBeamPipe(const G4String &nameIn, G4double lengthIn, G4double aper1=0, G4double aper2=0, G4double aper3=0, G4double aper4=0, G4Material *vacuumMaterialIn=nullptr, G4double beamPipeThicknessIn=0, G4Material *beamPipeMaterialIn=nullptr, const G4String &pointsFileIn="", const G4String &pointsUnitIn="")
create a flat ended beampipe
BDSBeamPipe * CommonFinalConstruction(const G4String &nameIn, G4Material *vacuumMaterialIn, G4Material *beamPipeMaterialIn, G4double lengthIn, G4double widthIn, G4double heightIn)
virtual BDSBeamPipe * CreateBeamPipe(const G4String &nameIn, G4double lengthIn, G4double aper1=0, G4double aper2=0, G4double aper3=0, G4double aper4=0, G4Material *vacuumMaterialIn=nullptr, G4double beamPipeThicknessIn=0, G4Material *beamPipeMaterialIn=nullptr, const G4String &pointsFileIn="", const G4String &pointsUnitIn="")
create a flat ended beampipe
void CreateGeneralAngledSolids(const G4String &nameIn, G4double lengthIn, G4double aper1In, G4double aper2In, G4double aper3In, G4double aper4In, G4double beamPipeThicknessIn, const G4ThreeVector &inputfaceIn, const G4ThreeVector &outputfaceIn)
A holder class for a piece of beam pipe geometry.
Definition: BDSBeamPipe.hh:45
General exception with possible name of object and message.
Definition: BDSException.hh:35
Holder for +- extents in 3 dimensions.
Definition: BDSExtent.hh:39
G4double lengthSafety
Cache of global constants variable.
G4double lengthSafetyLarge
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.