BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
BDSBeamPipeInfo.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 "BDSBeamPipeInfo.hh"
20#include "BDSBeamPipeType.hh"
21#include "BDSDebug.hh"
22#include "BDSException.hh"
23#include "BDSExtent.hh"
24#include "BDSMaterials.hh"
25#include "BDSUtilities.hh"
26
27#include "globals.hh" // geant4 types / globals
28#include "G4Material.hh"
29
30
32 G4double aper1In,
33 G4double aper2In,
34 G4double aper3In,
35 G4double aper4In,
36 G4Material* vacuumMaterialIn,
37 G4double beamPipeThicknessIn,
38 G4Material* beamPipeMaterialIn,
39 const G4ThreeVector& inputFaceNormalIn,
40 const G4ThreeVector& outputFaceNormalIn,
41 const G4String& pointsFileNameIn,
42 const G4String& pointsUnitIn):
43 beamPipeType(beamPipeTypeIn),
44 aper1(aper1In), aper2(aper2In), aper3(aper3In), aper4(aper4In),
45 aperOffsetX(0), aperOffsetY(0),
46 vacuumMaterial(vacuumMaterialIn),
47 beamPipeThickness(beamPipeThicknessIn),
48 beamPipeMaterial(beamPipeMaterialIn),
49 inputFaceNormal(inputFaceNormalIn),
50 outputFaceNormal(outputFaceNormalIn),
51 pointsFileName(pointsFileNameIn),
52 pointsUnit(pointsUnitIn)
53{
55}
56
57BDSBeamPipeInfo::BDSBeamPipeInfo(const G4String& beamPipeTypeIn,
58 G4double aper1In,
59 G4double aper2In,
60 G4double aper3In,
61 G4double aper4In,
62 const G4String& vacuumMaterialIn,
63 G4double beamPipeThicknessIn,
64 const G4String& beamPipeMaterialIn,
65 const G4ThreeVector& inputFaceNormalIn,
66 const G4ThreeVector& outputFaceNormalIn):
67 aper1(aper1In), aper2(aper2In), aper3(aper3In), aper4(aper4In),
68 aperOffsetX(0), aperOffsetY(0),
69 beamPipeThickness(beamPipeThicknessIn),
70 inputFaceNormal(inputFaceNormalIn),
71 outputFaceNormal(outputFaceNormalIn),
72 pointsFileName(""),
73 pointsUnit("mm")
74{
78
79 if (beamPipeType == BDSBeamPipeType::pointsfile)
80 {CheckAndSetPointsInfo(beamPipeTypeIn);}
82}
83
85 const G4String& beamPipeTypeIn,
86 G4double aper1In,
87 G4double aper2In,
88 G4double aper3In,
89 G4double aper4In,
90 const G4String& vacuumMaterialIn,
91 G4double beamPipeThicknessIn,
92 const G4String& beamPipeMaterialIn,
93 const G4ThreeVector& inputFaceNormalIn,
94 const G4ThreeVector& outputFaceNormalIn):
95 aperOffsetX(0), aperOffsetY(0),
96 inputFaceNormal(inputFaceNormalIn),
97 outputFaceNormal(outputFaceNormalIn)
98{
99 if (beamPipeTypeIn.empty())
100 {
101 beamPipeType = defaultInfo->beamPipeType;
102 pointsFileName = defaultInfo->pointsFileName; // copy even if empty
103 pointsUnit = defaultInfo->pointsUnit;
104 }
105 else
106 {
108 if (beamPipeType == BDSBeamPipeType::pointsfile)
109 {CheckAndSetPointsInfo(beamPipeTypeIn);}
110 }
111
112 if (!BDS::IsFinite(aper1In))
113 {aper1 = defaultInfo->aper1;}
114 else
115 {aper1 = aper1In;}
116 if (!BDS::IsFinite(aper2In))
117 {aper2 = defaultInfo->aper2;}
118 else
119 {aper2 = aper2In;}
120 if (!BDS::IsFinite(aper3In))
121 {aper3 = defaultInfo->aper3;}
122 else
123 {aper3 = aper3In;}
124 if (!BDS::IsFinite(aper4In))
125 {aper4 = defaultInfo->aper4;}
126 else
127 {aper4 = aper4In;}
128 if (!BDS::IsFinite(beamPipeThicknessIn))
129 {beamPipeThickness = defaultInfo->beamPipeThickness;}
130 else
131 {beamPipeThickness = beamPipeThicknessIn;}
132
133 if (vacuumMaterialIn == "")
134 {vacuumMaterial = defaultInfo->vacuumMaterial;}
135 else
136 {vacuumMaterial = BDSMaterials::Instance()->GetMaterial(vacuumMaterialIn);}
137 if (beamPipeMaterialIn == "")
138 {beamPipeMaterial = defaultInfo->beamPipeMaterial;}
139 else
140 {beamPipeMaterial = BDSMaterials::Instance()->GetMaterial(beamPipeMaterialIn);}
141
143}
144
145void BDSBeamPipeInfo::CheckAndSetPointsInfo(const G4String& beamPipeTypeIn)
146{
147 auto typeAndFileName = BDS::SplitOnColon(beamPipeTypeIn); // find first colon
148 G4String fname = typeAndFileName.second;
149 if (BDS::StrContains(fname, ":"))
150 {// optional second colon with units after it
151 auto fileNameAndUnit = BDS::SplitOnColon(fname);
152 pointsFileName = fileNameAndUnit.first;
153 pointsUnit = fileNameAndUnit.second;
154 }
155 else
156 {
157 pointsFileName = fname;
158 pointsUnit = "mm";
159 }
160}
161
163{
164 switch (beamPipeType.underlying())
165 {
166 case BDSBeamPipeType::circular:
167 {InfoOKForCircular(); break;}
168 case BDSBeamPipeType::elliptical:
169 {InfoOKForElliptical(); break;}
170 case BDSBeamPipeType::rectangular:
171 {InfoOKForRectangular(); break;}
172 case BDSBeamPipeType::lhc:
173 {InfoOKForLHC(); break;}
174 case BDSBeamPipeType::lhcdetailed:
175 {InfoOKForLHCDetailed(); break;}
176 case BDSBeamPipeType::rectellipse:
177 {InfoOKForRectEllipse(); break;}
178 case BDSBeamPipeType::racetrack:
179 {InfoOKForRaceTrack(); break;}
180 case BDSBeamPipeType::octagonal:
181 {InfoOKForOctagonal(); break;}
182 case BDSBeamPipeType::clicpcl:
183 {InfoOKForClicPCL(); break;}
184 default:
186 }
187}
188
190{
191 G4double extX = 0;
192 G4double extY = 0;
193 switch (beamPipeType.underlying())
194 {
195 case BDSBeamPipeType::circular:
196 case BDSBeamPipeType::circularvacuum:
197 {
198 extX = aper1;
199 extY = aper1;
200 break;
201 }
202 case BDSBeamPipeType::elliptical:
203 case BDSBeamPipeType::rectangular:
204 case BDSBeamPipeType::octagonal:
205 {
206 extX = aper1;
207 extY = aper2;
208 break;
209 }
210 case BDSBeamPipeType::lhc:
211 case BDSBeamPipeType::lhcdetailed:
212 case BDSBeamPipeType::rectellipse:
213 {
214 extX = std::min(aper1, aper3);
215 extY = std::min(aper2, aper3);
216 break;
217 }
218 case BDSBeamPipeType::racetrack:
219 {
220 extX = aper1 + aper3;
221 extY = aper2 + aper3;
222 break;
223 }
224 case BDSBeamPipeType::clicpcl:
225 {// this one is asymmetric so done separately
226 G4double extentX = aper1 + beamPipeThickness;
227 G4double extentYLow = -(std::abs(aper3) + beamPipeThickness);
228 G4double extentYHigh = aper2 + aper4 + beamPipeThickness;
229 BDSExtent ext = BDSExtent(-extentX, extentX,
230 extentYLow, extentYHigh,
231 0,0);
232 return ext;
233 break;
234 }
235 default:break;
236 }
237 BDSExtent ext = BDSExtent(extX, extY, 0);
238 return ext;
239}
240
242{
243 BDSExtent extentInner = ExtentInner();
244 G4double extX = extentInner.XPos(); // +ve values
245 G4double extY = extentInner.YPos();
246 extX += beamPipeThickness;
247 extY += beamPipeThickness;
248 return BDSExtent(extX, extY, 0);
249}
250
252{
253 BDSExtent ext = Extent();
254 return ext.MaximumAbsTransverse();
255}
256
258{
259 BDSExtent ext = ExtentInner();
260 return ext.MinimumAbsTransverse();
261}
262
264 G4bool setAper2,
265 G4bool setAper3,
266 G4bool setAper4)
267{
268#ifdef BDSDEBUG
269 G4cout << __METHOD_NAME__ << G4endl;
270 G4cout << "aper1: " << aper1 << " check it? " << setAper1 << G4endl;
271 G4cout << "aper2: " << aper2 << " check it? " << setAper2 << G4endl;
272 G4cout << "aper3: " << aper3 << " check it? " << setAper3 << G4endl;
273 G4cout << "aper4: " << aper4 << " check it? " << setAper4 << G4endl;
274#endif
275 G4bool shouldExit = false;
276 if (setAper1)
277 {
278 if (!BDS::IsFinite(aper1))
279 {G4cerr << "\"aper1\" not set, but required to be" << G4endl; shouldExit = true;}
280 }
281
282 if (setAper2)
283 {
284 if (!BDS::IsFinite(aper2))
285 {G4cerr << "\"aper2\" not set, but required to be" << G4endl; shouldExit = true;}
286 }
287
288 if (setAper3)
289 {
290 if (!BDS::IsFinite(aper3))
291 {G4cerr << "\"aper3\" not set, but required to be" << G4endl; shouldExit = true;}
292 }
293
294 if (setAper4)
295 {
296 if (!BDS::IsFinite(aper4))
297 {G4cerr << "\"aper4\" not set, but required to be" << G4endl; shouldExit = true;}
298 }
299
300 if (shouldExit)
301 {throw BDSException(__METHOD_NAME__, "aperture parameter missing");}
302}
303
305{
306 CheckRequiredParametersSet(true, false, false, false);
307}
308
310{
311 CheckRequiredParametersSet(true, true, false, false);
312}
313
315{
316 CheckRequiredParametersSet(true, true, false, false);
317}
318
320{
321 CheckRequiredParametersSet(true, true, true, false);
322
323 if ((aper3 > aper1) && (aper2 < aper3))
324 {throw BDSException(__METHOD_NAME__, "\"aper3\" > \"aper1\" (or \"beamPipeRadius\") for lhc aperture model - will not produce desired shape");}
325
326 if ((aper3 > aper2) && (aper1 < aper3))
327 {throw BDSException(__METHOD_NAME__, "\"aper3\" > \"aper2\" (or \"beamPipeRadius\") for lhc aperture model - will not produce desired shape");}
328}
329
331{
332 InfoOKForLHC();
333}
334
336{
337 CheckRequiredParametersSet(true, true, true, true);
338
339 // TODO
340}
341
343{
344 CheckRequiredParametersSet(true, true, true, false);
345}
346
348{
349 CheckRequiredParametersSet(true, true, true, true);
350
351 if (aper3 >= aper1)
352 {throw BDSException(__METHOD_NAME__, "aper3 is >= aper1 - invalid for an octagonal aperture");}
353 if (aper4 >= aper2)
354 {throw BDSException(__METHOD_NAME__, "aper4 is >= aper2 - invalid for an octagonal aperture");}
355}
356
358{
359 CheckRequiredParametersSet(true, true, true, false);
360}
Holder class for all information required to describe a beam pipe model.
G4double IndicativeRadiusInner() const
Return an indicative inner extent for the beam pipe vacuum.
BDSExtent Extent() const
G4Material * beamPipeMaterial
Public member for direct access.
G4double aper3
Public member for direct access.
G4String pointsUnit
Public member for direct access.
BDSBeamPipeInfo()=delete
Deleted default constructor to ensure one of supplied constructors is used.
void InfoOKForLHC()
Aperture info check for lhc aperture.
G4double aper1
Public member for direct access.
void InfoOKForOctagonal()
Aperture info check for octagonal aperture.
G4Material * vacuumMaterial
Public member for direct access.
void CheckAndSetPointsInfo(const G4String &beamPipeTypeIn)
Parse the type string to extract the file name and the optional units and assign to member variables.
void InfoOKForClicPCL()
Aperture info check for CLIC PCL aperture.
G4double beamPipeThickness
Public member for direct access.
void InfoOKForRaceTrack()
Aperture info check for racetrack aperture.
void InfoOKForElliptical()
Aperture info check for elliptical aperture.
BDSExtent ExtentInner() const
Return an extent for just the raw aperture.
void InfoOKForRectEllipse()
Aperture info check for rectellipse aperture.
G4String pointsFileName
Public member for direct access.
G4double aper4
Public member for direct access.
void InfoOKForLHCDetailed()
Aperture info check for lhc detailed aperture.
void InfoOKForCircular()
Aperture info check for circular aperture.
void CheckRequiredParametersSet(G4bool setAper1, G4bool setAper2, G4bool setAper3, G4bool setAper4)
void InfoOKForRectangular()
Aperture info check for rectangular aperture.
G4double IndicativeRadius() const
Return an indicative extent of the beam pipe - typically the maximum of x or y extent.
G4double aper2
Public member for direct access.
BDSBeamPipeType beamPipeType
Public member for direct access.
General exception with possible name of object and message.
Definition: BDSException.hh:35
Holder for +- extents in 3 dimensions.
Definition: BDSExtent.hh:39
G4double XPos() const
Accessor.
Definition: BDSExtent.hh:66
G4double MaximumAbsTransverse() const
Return the maximum absolute value considering only x,y.
Definition: BDSExtent.cc:171
G4double MinimumAbsTransverse() const
Return the minimum absolute value considering only x,y.
Definition: BDSExtent.cc:177
G4double YPos() const
Accessor.
Definition: BDSExtent.hh:68
static BDSMaterials * Instance()
Singleton pattern access.
Definition: BDSMaterials.cc:38
G4Material * GetMaterial(G4String material) const
Get material by name.
type underlying() const
return underlying value (can be used in switch statement)
std::pair< G4String, G4String > SplitOnColon(const G4String &formatAndPath)
G4bool StrContains(const G4String &str, const G4String &test)
Utility function to simplify lots of syntax changes for pedantic g4 changes.
Definition: BDSUtilities.cc:66
BDSBeamPipeType DetermineBeamPipeType(G4String apertureType)
Function that gives corresponding enum value for string (case-insensitive).
G4bool IsFinite(G4double value, G4double tolerance=std::numeric_limits< double >::epsilon())