BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSMagnet.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 "BDSBeamPipe.hh"
20#include "BDSBeamPipeFactory.hh"
21#include "BDSBeamPipeInfo.hh"
22#include "BDSBeamPipeType.hh"
23#include "BDSFieldBuilder.hh"
24#include "BDSFieldInfo.hh"
25#include "BDSIntegratorType.hh"
26#include "BDSMagnetGeometryType.hh"
27#include "BDSMagnetOuter.hh"
28#include "BDSMagnetOuterInfo.hh"
29#include "BDSMagnetOuterFactory.hh"
30#include "BDSMagnetOuterFactoryLHC.hh"
31#include "BDSMagnetStrength.hh"
32#include "BDSMagnetType.hh"
33#include "BDSMagnet.hh"
34#include "BDSUtilities.hh"
35
36#include "G4Box.hh"
37#include "G4LogicalVolume.hh"
38#include "G4Material.hh"
39#include "G4PVPlacement.hh"
40#include "G4String.hh"
41#include "G4ThreeVector.hh"
42#include "G4Transform3D.hh"
43#include "G4Types.hh"
44#include "G4VPhysicalVolume.hh"
45
46#include "CLHEP/Units/SystemOfUnits.h"
47
48class G4Userlimits;
49
51 const G4String& nameIn,
52 G4double lengthIn,
53 BDSBeamPipeInfo* beamPipeInfoIn,
54 BDSMagnetOuterInfo* magnetOuterInfoIn,
55 BDSFieldInfo* vacuumFieldInfoIn,
56 G4double angleIn,
57 BDSFieldInfo* outerFieldInfoIn,
58 G4bool isThinIn):
59 BDSAcceleratorComponent(nameIn, lengthIn, angleIn, typeIn.ToString(), beamPipeInfoIn),
60 magnetType(typeIn),
61 magnetOuterInfo(magnetOuterInfoIn),
62 vacuumFieldInfo(vacuumFieldInfoIn),
63 outerFieldInfo(outerFieldInfoIn),
64 beampipe(nullptr),
65 placeBeamPipe(false),
66 magnetOuterOffset(G4ThreeVector(0,0,0)),
67 outer(nullptr),
68 beamPipePlacementTransform(G4Transform3D()),
69 isThin(isThinIn)
70{
71 magnetOuterInfo->name += "_outer";
72 horizontalWidth = magnetOuterInfoIn->horizontalWidth;
74
75 beampipe = nullptr;
76 outer = nullptr;
77
78 placeBeamPipe = false;
79
80 // It's not possible to build advanced outer geometry for a very thin magnet.
81 if (lengthIn < 1e-4*CLHEP::m) // 100um minimum length for geometry
82 {magnetOuterInfo->geometryType = BDSMagnetGeometryType::none;}
83 // No beam pipe geometry for really short 'magnets'
84 if (lengthIn < 1e-6*CLHEP::m)
85 {
86 GetBeamPipeInfo()->beamPipeType = BDSBeamPipeType::circularvacuum;
87 isThin = true;
88 }
89}
90
92{
93 G4String result = "none";
94 switch (typeIn.underlying())
95 {
96 case BDSMagnetType::hkicker:
97 case BDSMagnetType::vkicker:
98 case BDSMagnetType::muonspoiler:
99 case BDSMagnetType::sectorbend:
100 case BDSMagnetType::rectangularbend:
101 {result = "field"; break;}
102 case BDSMagnetType::quadrupole:
103 {result = "k1"; break;}
104 case BDSMagnetType::sextupole:
105 {result = "k2"; break;}
106 case BDSMagnetType::octupole:
107 {result = "k3"; break;}
108 case BDSMagnetType::decapole:
109 {result = "k4"; break;}
110 default:
111 {break;} // leave as none without complaint
112 };
113
114 return result;
115}
116
117void BDSMagnet::SetInputFaceNormal(const G4ThreeVector& input)
118{
119 if (outer)
120 {outer->SetInputFaceNormal(input);}
121}
122
123void BDSMagnet::SetOutputFaceNormal(const G4ThreeVector& output)
124{
125 if (outer)
126 {outer->SetOutputFaceNormal(output);}
127}
128
129G4String BDSMagnet::Material() const
130{
131 if (magnetOuterInfo)
132 {return magnetOuterInfo->outerMaterial->GetName();}
133 else
135}
136
138{
142 BuildOuter();
143 // Instead of BDSAcceleratorComponent::Build just call BuildContainerLogicalVolume
144 // to control user limits in the case where there is no container and we just inherit
145 // the beam pipe container
147 BuildOuterField(); // must be done when the containerLV exists
148 PlaceComponents(); // place things (if needed) in container
149}
150
152{
153 beampipe = BDSBeamPipeFactory::Instance()->CreateBeamPipe(name+"_bp",
156
158
161
163
167}
168
170{
171 if (vacuumFieldInfo)
172 {
173 if (beamPipePlacementTransform != G4Transform3D::Identity)
174 {
175 G4Transform3D newFieldTransform = vacuumFieldInfo->Transform() * beamPipePlacementTransform;
176 vacuumFieldInfo->SetTransform(newFieldTransform);
177 }
178 // can use containerLV for field as we don't construct any geometry with thin elements.
179 if (isThin)
180 {
183 true);
184 }
185 else
186 {
189 true);
190 }
191 }
192}
193
195{
196 G4double outerLength = chordLength - 2*lengthSafety;
197 // The outerFieldInfo is required only so we don't reuse geometry when we need a unique field on it.
200 outerLength,
202 beampipe);
203
204 if (outer)
205 {
206 // copy necessary bits out of BDSGeometryComponent that holds
207 // container information for whole magnet object provided by
208 // magnet outer factory.
210 containerSolid = container->GetContainerSolid()->Clone();
211 G4ThreeVector contOffset = container->GetPlacementOffset();
212 // set the main offset of the whole magnet which is placed w.r.t. the
213 // zero coordinate of the container solid
214 SetPlacementOffset(contOffset);
215
217 InheritExtents(container, contOffset); // update extents
218
219 // Only clear after extents etc have been used
221
222 endPieceBefore = outer->EndPieceBefore();
223 endPieceAfter = outer->EndPieceAfter();
224
228 }
229}
230
232{
233 if (outer && outerFieldInfo)
234 {
235 // determine key for this specific magnet instance
236 G4String scalingKey = DetermineScalingKey(magnetType);
237
238 BDSMagnetStrength* scalingStrength = vacuumFieldInfo ? vacuumFieldInfo->MagnetStrength() : nullptr;
239 G4LogicalVolume* vol = outer->GetContainerLogicalVolume();
241 vol,
242 true,
243 scalingStrength,
244 scalingKey);
245 // Attach to the container but don't propagate to daughter volumes. This ensures
246 // any gap between the beam pipe and the outer also has a field.
248 containerLogicalVolume,
249 false,
250 scalingStrength,
251 scalingKey);
252
253 // in the case of LHC-style geometry, override the second beam pipe (which is a daughter of the outer)
254 // field to be the opposite sign of the main vacuum field (same strength)
255 auto mgt = magnetOuterInfo->geometryType;
256 if (mgt == BDSMagnetGeometryType::lhcleft || mgt == BDSMagnetGeometryType::lhcright)
257 {
258 std::set<BDSGeometryComponent*> daughtersSet = outer->GetAllDaughters();
259 std::vector<BDSGeometryComponent*> daughters(daughtersSet.begin(), daughtersSet.end());
260 if (daughters.size() == 1 && vacuumFieldInfo)
261 {
262 BDSFieldInfo* secondBPField = new BDSFieldInfo(*vacuumFieldInfo);
263 G4double sign = mgt == BDSMagnetGeometryType::lhcleft ? 1.0 : -1.0;
264 secondBPField->Translate(G4ThreeVector(sign * BDSMagnetOuterFactoryLHC::beamSeparation, 0, 0));
265 (*(secondBPField->MagnetStrength()))["field"] *= -1; // flip the sign
266 if (BDS::IsFinite((*(secondBPField->MagnetStrength()))["k1"]))
267 {(*(secondBPField->MagnetStrength()))["k1"] *= -1;}
268 secondBPField->SetIntegratorType(BDSIntegratorType::g4classicalrk4);
270 daughters[0]->GetContainerLogicalVolume(),
271 true);
272 }
273 }
274 }
275}
276
278{
279 // note beam pipe is not optional!
280 if (outer)
281 {//build around that
282 // container solid will have been updated in BuildOuter if the outer exists
283 containerLogicalVolume = new G4LogicalVolume(containerSolid,
285 name + "_container_lv");
286
287 // user limits - provided by BDSAcceleratorComponent
288 containerLogicalVolume->SetUserLimits(userLimits);
289 containerLogicalVolume->SetVisAttributes(containerVisAttr);
290
291 placeBeamPipe = true;
292 }
293 else
294 {
295 // use beam pipe container volume as ours and no need to place beam pipe
296 containerSolid = beampipe->GetContainerSolid();
297 containerLogicalVolume = beampipe->GetContainerLogicalVolume();
299 placeBeamPipe = false;
300 }
301}
302
304{
305 if (placeBeamPipe)
306 {
307 G4ThreeVector beamPipeOffset = -1*GetPlacementOffset();
308 // place beampipe
309 G4PVPlacement* beamPipePV = new G4PVPlacement(nullptr, // rotation
310 beamPipeOffset, // position in container
311 beampipe->GetContainerLogicalVolume(), // its logical volume
312 name + "_beampipe_pv", // its name
313 containerLogicalVolume, // its mother volume
314 false, // no boolean operation
315 0, // copy number
317
318 RegisterPhysicalVolume(beamPipePV);
319 }
320
321 if (outer)
322 {
323 //G4ThreeVector placementOffset = magnetOuterOffset + outer->GetPlacementOffset();
324 G4ThreeVector outerOffset = outer->GetPlacementOffset();
325
326 // place outer volume
327 G4PVPlacement* magnetOuterPV = new G4PVPlacement(nullptr, // rotation
328 outerOffset, // at normally (0,0,0)
329 outer->GetContainerLogicalVolume(), // its logical volume
330 name+"_outer_pv", // its name
331 containerLogicalVolume, // its mother volume
332 false, // no boolean operation
333 0, // copy number
335
336 RegisterPhysicalVolume(magnetOuterPV);
337 }
338}
339
341{
342 delete outerFieldInfo;
343 outerFieldInfo = outerFieldInfoIn;
344}
345
346void BDSMagnet::SetVacuumField(BDSFieldInfo* vacuumFieldInfoIn)
347{
348 delete vacuumFieldInfo;
349 vacuumFieldInfo = vacuumFieldInfoIn;
350}
351
352BDSMagnet::~BDSMagnet()
353{
354 delete magnetOuterInfo;
355 delete vacuumFieldInfo;
356 delete outerFieldInfo;
357}
358
360{
362 if (vacuumFieldInfo)
363 {vacuumFieldInfo->SetUsePlacementWorldTransform(true);}
364 if (outerFieldInfo)
365 {outerFieldInfo->SetUsePlacementWorldTransform(true);}
366}
Abstract class that represents a component of an accelerator.
G4UserLimits * userLimits
Cache of user limits.
virtual G4String Material() const
Return the name of a material associated with the component - ie the primary material.
void SetAcceleratorVacuumLogicalVolume(G4LogicalVolume *accVacLVIn)
const G4String name
Const protected member variable that may not be changed by derived classes.
static G4Material * emptyMaterial
Useful variable often used in construction.
static G4double lengthSafety
Useful variable often used in construction.
static G4bool checkOverlaps
Useful variable often used in construction.
G4double angle
Protected member variable that can be modified by derived classes.
BDSBeamPipeInfo * GetBeamPipeInfo() const
G4double chordLength
Protected member variable that can be modified by derived classes.
static G4VisAttributes * containerVisAttr
Useful variable often used in construction.
BDSBeamPipeInfo * beamPipeInfo
Optional beam pipe recipe that is written out to the survey if it exists.
virtual void SetFieldUsePlacementWorldTransform()
static BDSBeamPipeFactory * Instance()
Singleton accessor.
Holder class for all information required to describe a beam pipe model.
BDSBeamPipeType beamPipeType
Public member for direct access.
G4LogicalVolume * GetVacuumLogicalVolume() const
Access the vacuum volume to set fields and limits.
Definition: BDSBeamPipe.hh:65
G4ThreeVector InputFaceNormal() const
Accessor.
Definition: BDSBeamPipe.hh:73
G4ThreeVector OutputFaceNormal() const
Accessor.
Definition: BDSBeamPipe.hh:74
std::set< G4LogicalVolume * > GetVolumesForField() const
Definition: BDSBeamPipe.cc:50
void RegisterFieldForConstruction(const BDSFieldInfo *info, G4LogicalVolume *logicalVolume, const G4bool propagateToDaughters=false, const BDSMagnetStrength *magnetStrengthForScaling=nullptr, const G4String &scalingKey="none")
static BDSFieldBuilder * Instance()
Singleton pattern accessor.
All info required to build complete field of any type.
Definition: BDSFieldInfo.hh:66
void SetTransform(const G4Transform3D &transformIn)
Set the field definition transform.
G4Transform3D Transform() const
Transform for the field definition only.
void Translate(const G4ThreeVector &translationIn)
BDSMagnetStrength * MagnetStrength() const
Accessor.
A generic geometry component for a bdsim model.
virtual std::set< BDSGeometryComponent * > GetAllDaughters() const
Accessor - see member for more info.
G4LogicalVolume * GetContainerLogicalVolume() const
Accessor - see member for more info.
void InheritExtents(BDSGeometryComponent const *const anotherComponent)
Update the extents of this object with those of another object.
G4Transform3D GetPlacementTransform() const
Accessor - see member for more info.
void SetPlacementOffset(const G4ThreeVector &offsetIn)
Set the offset from 0,0,0 that the object should ideally be placed in its parent.
void RegisterDaughter(BDSGeometryComponent *anotherComponent)
void RegisterPhysicalVolume(G4VPhysicalVolume *physicalVolume)
G4VSolid * GetContainerSolid() const
Accessor - see member for more info.
G4ThreeVector GetPlacementOffset() const
Accessor - see member for more info.
static const G4double beamSeparation
Used in many places - make it a constant in the code and put here as most relevant.
static BDSMagnetOuterFactory * Instance()
Singleton accessor.
BDSMagnetOuter * CreateMagnetOuter(BDSMagnetType magnetType, BDSMagnetOuterInfo *outerInfo, G4double outerLength, G4double chordLength, BDSBeamPipe *beampipe)
Holder struct of all information required to create the outer geometry of a magnet.
BDSGeometryComponent * GetMagnetContainer() const
void ClearMagnetContainer()
Clear the memory of the now unneeded magnet container object.
void SetInputFaceNormal(const G4ThreeVector &input)
Setter for face normals.
G4ThreeVector InputFaceNormal() const
Accessor.
BDSSimpleComponent * EndPieceAfter() const
Access the end piece.
G4ThreeVector OutputFaceNormal() const
Accessor.
void SetOutputFaceNormal(const G4ThreeVector &output)
Setter for face normals.
BDSSimpleComponent * EndPieceBefore() const
Access the end piece.
Efficient storage of magnet strengths.
virtual void SetOutputFaceNormal(const G4ThreeVector &output)
Update face normal and also to beam pipe and magnet outer.
Definition: BDSMagnet.cc:123
G4Transform3D beamPipePlacementTransform
Definition: BDSMagnet.hh:162
BDSFieldInfo * vacuumFieldInfo
Field information for vacuum field.
Definition: BDSMagnet.hh:132
G4bool isThin
Boolean to store if the element is thin - will have no geometry constructed.
Definition: BDSMagnet.hh:165
G4double horizontalWidth
For outer volume construction.
Definition: BDSMagnet.hh:146
BDSMagnetType magnetType
Magnet type.
Definition: BDSMagnet.hh:126
BDSMagnetOuterInfo * magnetOuterInfo
Model information for the outer volume construction.
Definition: BDSMagnet.hh:129
virtual void SetFieldUsePlacementWorldTransform()
Override function as we have different field recipe objects.
Definition: BDSMagnet.cc:359
BDSBeamPipe * beampipe
The constructed beampipe.
Definition: BDSMagnet.hh:138
virtual void BuildBeampipe()
Definition: BDSMagnet.cc:151
virtual void BuildVacuumField()
Construct the field for the vacuum and attach it.
Definition: BDSMagnet.cc:169
G4bool placeBeamPipe
Definition: BDSMagnet.hh:143
virtual void BuildContainerLogicalVolume()
Definition: BDSMagnet.cc:277
BDSFieldInfo * outerFieldInfo
Field information for outer magnetic field (optional)
Definition: BDSMagnet.hh:135
BDSMagnet()=delete
Private default constructor to force the use of the supplied one.
virtual void Build()
Definition: BDSMagnet.cc:137
BDSMagnetOuter * outer
The assembled outer magnet geometry.
Definition: BDSMagnet.hh:158
virtual G4String Material() const
Accessor to outer material if it exists.
Definition: BDSMagnet.cc:129
virtual void BuildOuter()
Definition: BDSMagnet.cc:194
static G4String DetermineScalingKey(BDSMagnetType typeIn)
Definition: BDSMagnet.cc:91
void SetOuterField(BDSFieldInfo *outerFieldInfoIn)
@ { Delete existing field info and replace.
Definition: BDSMagnet.cc:340
virtual void SetInputFaceNormal(const G4ThreeVector &input)
Update face normal and also to beam pipe and magnet outer.
Definition: BDSMagnet.cc:117
virtual void BuildOuterField()
Construct the magnetic field for the outer magnet geometry.
Definition: BDSMagnet.cc:231
G4double containerRadius
Definition: BDSMagnet.hh:150
virtual void PlaceComponents()
Definition: BDSMagnet.cc:303
type underlying() const
return underlying value (can be used in switch statement)
G4ThreeVector RotateToReferenceFrame(G4ThreeVector faceNormal, G4double fullAngle)
G4bool IsFinite(G4double value, G4double tolerance=std::numeric_limits< double >::epsilon())