BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSAcceleratorComponent.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 "BDSAcceleratorComponent.hh"
20#include "BDSBeamPipeInfo.hh"
21#include "BDSDebug.hh"
22#include "BDSException.hh"
23#include "BDSFieldBuilder.hh"
24#include "BDSFieldInfo.hh"
25#include "BDSGlobalConstants.hh"
26#include "BDSMaterials.hh"
27#include "BDSTunnelInfo.hh"
28#include "BDSUtilities.hh"
29
30#include "G4AssemblyVolume.hh"
31#include "G4LogicalVolume.hh"
32#include "G4Material.hh"
33#include "G4ThreeVector.hh"
34#include "G4UserLimits.hh"
35
36#include <algorithm>
37#include <cmath>
38#include <string>
39
40G4Material* BDSAcceleratorComponent::emptyMaterial = nullptr;
41G4Material* BDSAcceleratorComponent::worldMaterial = nullptr;
46G4VisAttributes* BDSAcceleratorComponent::containerVisAttr = nullptr;
48
50 G4double arcLengthIn,
51 G4double angleIn,
52 const G4String& typeIn,
53 BDSBeamPipeInfo* beamPipeInfoIn,
54 const G4ThreeVector& inputFaceNormalIn,
55 const G4ThreeVector& outputFaceNormalIn,
56 BDSFieldInfo* fieldInfoIn):
57 BDSGeometryComponent(nullptr,nullptr),
58 name(nameIn),
59 arcLength(arcLengthIn),
60 type(typeIn),
61 angle(angleIn),
62 beamPipeInfo(beamPipeInfoIn),
63 endPieceBefore(nullptr),
64 endPieceAfter(nullptr),
65 userLimits(nullptr),
66 fieldInfo(fieldInfoIn),
67 copyNumber(-1), // -1 initialisation since it will be incremented when placed
68 inputFaceNormal(inputFaceNormalIn),
69 outputFaceNormal(outputFaceNormalIn)
70{
71 // initialise static members
72 if (!emptyMaterial)
73 {
74 const auto globals = BDSGlobalConstants::Instance(); // shortcut
75 emptyMaterial = BDSMaterials::Instance()->GetMaterial(globals->EmptyMaterial());
76 worldMaterial = BDSMaterials::Instance()->GetMaterial(globals->WorldMaterial());
77 lengthSafety = globals->LengthSafety();
78 lengthSafetyLarge = globals->LengthSafetyLarge();
79 checkOverlaps = globals->CheckOverlaps();
80 sensitiveOuter = globals->SensitiveOuter();
81 sensitiveVacuum = globals->StoreELossVacuum();
82 containerVisAttr = BDSGlobalConstants::Instance()->ContainerVisAttr();
83 }
84
85 // Prevent negative length components.
86 if (arcLength < 0)
87 {
88 G4String message = "Negative length for component named \"" + name + "\" with length " + std::to_string(arcLength);
89 throw BDSException(__METHOD_NAME__, message);
90 }
91
92 // calculate the chord length if the angle is finite
93 if (BDS::IsFinite(angleIn))
94 {chordLength = 2.0 * arcLengthIn * std::sin(0.5*angleIn) / angleIn;}
95 else
96 {chordLength = arcLengthIn;}
97#ifdef BDSDEBUG
98 G4cout << __METHOD_NAME__ << "angle: " << angleIn << G4endl;
99 G4cout << __METHOD_NAME__ << "arc length: " << arcLengthIn << G4endl;
100 G4cout << __METHOD_NAME__ << "chord length: " << chordLength << G4endl;
101#endif
102
103 initialised = false;
104}
105
106BDSAcceleratorComponent::~BDSAcceleratorComponent()
107{
108 delete beamPipeInfo;
109 // Don't delete usersLimits as could be globals one. If not it'll be registered
110 // with BDSGeometryComponent
111}
112
114{
115 if (initialised)
116 {return;} // protect against duplicated initialisation and memory leaks
117
118 Build();
119
120 // field construction must be done after all the geometry is constructed if the
121 // field is to propagate to the daughter volumes correctly.
122 if (fieldInfo)
123 {
125 containerLogicalVolume,
126 true);
127 }
128 initialised = true; // record that this component has been initialised
129}
130
132{
133 BuildContainerLogicalVolume(); // pure virtual provided by derived class
136 if (containerLogicalVolume)
137 {containerLogicalVolume->SetVisAttributes(containerVisAttr);}
138}
139
141{
142 fieldInfo = fieldInfoIn;
143}
144
146{
147 if (!BDS::IsFinite(angle))
148 {return 0;}
149 G4double bendingRadius = arcLength / angle;
150 G4double sagitta = bendingRadius - std::sqrt(std::pow(bendingRadius,2) - std::pow(chordLength, 2));
151 return sagitta;
152}
153
155{
156 G4ThreeVector zeroAngle = G4ThreeVector(0,0,-1);
157 G4ThreeVector cross = inputFaceNormal.cross(zeroAngle);
158 G4double det = cross.mag2();
159 return BDS::IsFinite(det); // if finite, there is an angle and if not, no angle
160}
161
163{
164 G4ThreeVector zeroAngle = G4ThreeVector(0,0,1);
165 G4ThreeVector cross = outputFaceNormal.cross(zeroAngle);
166 G4double det = cross.mag2();
167 return BDS::IsFinite(det);
168}
169
171{
172 // user limits
173 auto defaultUL = BDSGlobalConstants::Instance()->DefaultUserLimits();
174 //copy the default and update with the length of the object rather than the default 1m
175 G4UserLimits* ul = BDS::CreateUserLimits(defaultUL, std::max(chordLength, arcLength));
176 if (ul != defaultUL) // if it's not the default register it
177 {RegisterUserLimits(ul);}
178 userLimits = ul; // assign to member
179}
180
182{
183 if (!userLimits)
184 {return;}
185 if (containerLogicalVolume || containerAssembly)
186 {
187 if (containerIsAssembly && containerAssembly)
188 {AttachUserLimitsToAssembly(containerAssembly, userLimits);}
189 else if (containerLogicalVolume)
190 {containerLogicalVolume->SetUserLimits(userLimits);}
191 }
192}
193
195{
196 // get full set minus ones marked to be excluded completely from biasing
197 std::set<G4LogicalVolume*> result = GetAllBiasingVolumes();
198
199 for (auto lv : acceleratorVacuumLV)
200 {result.erase(lv);}
201 return result;
202}
203
205{
206 if (fieldInfo)
207 {fieldInfo->SetUsePlacementWorldTransform(true);}
208}
BDSFieldInfo * fieldInfo
Recipe for field that could overlay this whole component.
G4UserLimits * userLimits
Cache of user limits.
G4bool AngledInputFace() const
Whether the face normal is angled at all w.r.t. the incoming / outgoing reference trajectory.
G4double Sagitta() const
Calculate the sagitta - ie the distance between the chord and the arc at the centre.
G4ThreeVector outputFaceNormal
Output face unit normal vector in outgoing reference coordinate frame.
G4ThreeVector inputFaceNormal
Input face unit normal vector in incoming reference coordinate frame.
const G4double arcLength
Const protected member variable that may not be changed by derived classes.
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 G4bool sensitiveOuter
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.
virtual std::set< G4LogicalVolume * > GetAcceleratorMaterialLogicalVolumes() const
Return a set of logical volumes excluding the ones in the 'vacuum' set.
static G4Material * worldMaterial
Useful variable often used in construction.
BDSAcceleratorComponent()=delete
virtual void BuildContainerLogicalVolume()=0
void SetField(BDSFieldInfo *fieldInfoIn)
Set the field definition for the whole component.
virtual void AttachUserLimits() const
Doesn't change member variables, but may change their contents.
G4double chordLength
Protected member variable that can be modified by derived classes.
G4bool AngledOutputFace() const
Whether the face normal is angled at all w.r.t. the incoming / outgoing reference trajectory.
std::set< G4LogicalVolume * > acceleratorVacuumLV
A set of logical volumes we classify as 'vacuum' for biasing purposes.
static G4VisAttributes * containerVisAttr
Useful variable often used in construction.
static G4bool sensitiveVacuum
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()
Holder class for all information required to describe a beam pipe model.
General exception with possible name of object and message.
Definition: BDSException.hh:35
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
A generic geometry component for a bdsim model.
virtual std::set< G4LogicalVolume * > GetAllBiasingVolumes() const
Return all logical volumes that should be used for biasing minus any that are in the excluded set.
static void AttachUserLimitsToAssembly(G4AssemblyVolume *av, G4UserLimits *ul)
Utility function to apply user limits to an assembly volume as there's not interface.
void RegisterUserLimits(G4UserLimits *userLimit)
G4bool containerIsAssembly
True if the 'container' is really an assembly; false if an LV.
static BDSGlobalConstants * Instance()
Access method.
static BDSMaterials * Instance()
Singleton pattern access.
Definition: BDSMaterials.cc:39
G4Material * GetMaterial(G4String material) const
Get material by name.
G4UserLimits * CreateUserLimits(G4UserLimits *defaultUL, G4double length, G4double fraction=1.6)
G4bool IsFinite(G4double value, G4double tolerance=std::numeric_limits< double >::epsilon())