BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSBLMFactory.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 "BDSBLM.hh"
20#include "BDSBLMFactory.hh"
21#include "BDSBLMType.hh"
22#include "BDSDebug.hh"
23#include "BDSException.hh"
24#include "BDSExtent.hh"
25#include "BDSGeometryExternal.hh"
26#include "BDSGeometryFactory.hh"
27#include "BDSGlobalConstants.hh"
28#include "BDSMaterials.hh"
29#include "BDSUtilities.hh"
30
31#include "parser/blmplacement.h"
32
33#include "globals.hh"
34#include "G4Box.hh"
35#include "G4LogicalVolume.hh"
36#include "G4Orb.hh"
37#include "G4Tubs.hh"
38
39#include <map>
40
41BDSBLMFactory::BDSBLMFactory()
42{;}
43
44BDSBLMFactory::~BDSBLMFactory()
45{;}
46
47void BDSBLMFactory::PositiveFinite(G4double value,
48 const G4String& parameterName,
49 const G4String& blmName) const
50{
51 if (!BDS::IsFinite(value) || value < 0)
52 {throw BDSException(__METHOD_NAME__, "\"" + parameterName + "\" must be +ve finite for blm named \"" + blmName + "\"");}
53}
54
55BDSBLM* BDSBLMFactory::CreateBLM(const GMAD::BLMPlacement& bp,
56 G4VSensitiveDetector* sd)
57{
58 return CreateBLM(G4String(bp.name),
59 G4String(bp.geometryFile),
60 G4String(bp.geometryType),
61 G4String(bp.blmMaterial),
62 bp.blm1 * CLHEP::m,
63 bp.blm2 * CLHEP::m,
64 bp.blm3 * CLHEP::m,
65 bp.blm4 * CLHEP::m,
66 sd,
67 bp.bias);
68}
69
70BDSBLM* BDSBLMFactory::CreateBLM(const G4String& name,
71 const G4String& geometryFile,
72 const G4String& geometryType,
73 const G4String& material,
74 G4double blm1,
75 G4double blm2,
76 G4double blm3,
77 G4double /*blm4*/,
78 G4VSensitiveDetector* sd,
79 const G4String& bias)
80{
81 // if geometry file is specified then we load the external file.
82 BDSBLM* result = nullptr;
83 if (!geometryFile.empty())
84 {
85 BDSGeometryExternal* blmGeom = BDSGeometryFactory::Instance()->BuildGeometry(name, geometryFile);
86 result = new BDSBLM(blmGeom);
87 }
88 else
89 {
90 BDSBLMType shape = BDS::DetermineBLMType(geometryType);
91 switch (shape.underlying())
92 {
93 case BDSBLMType::cylinder:
94 {result = BuildBLMCylinder(name, material, blm1, blm2); break;}
95 case BDSBLMType::cube:
96 {result = BuildBLMCube(name, material, blm1, blm2, blm3); break;}
97 case BDSBLMType::sphere:
98 {result = BuildBLMSphere(name, material, blm1); break;}
99 default:
100 {break;}
101 }
102 }
103 if (!result)
104 {return result;}
105
106 // attach sensitivity
107 for (auto lv : result->GetAllLogicalVolumes())
108 {lv->SetSensitiveDetector(sd);}
109 result->GetContainerLogicalVolume()->SetSensitiveDetector(sd); // attach separately to the container
110
111 // biasing - name of objects attached here and built later
112 result->SetBias(bias);
113 return result;
114}
115
117 const G4String& material,
118 G4double halfLength,
119 G4double radius)
120{
121 PositiveFinite(halfLength, "blm1 (half length)", name);
122 PositiveFinite(radius, "blm2 (radius)", name);
123 G4Tubs* shape = new G4Tubs(name + "_solid",
124 0,
125 radius,
126 halfLength,
127 0,
128 CLHEP::twopi);
129 BDSExtent ext = BDSExtent(radius, radius, halfLength);
130 return CommonConstruction(name, material, shape, ext);
131}
132
134 const G4String& material,
135 G4double halfLengthX,
136 G4double halfLengthY,
137 G4double halfLengthZ)
138{
139 PositiveFinite(halfLengthX, "blm1 (half length in x)", name);
140 PositiveFinite(halfLengthY, "blm2 (half length in y)", name);
141 PositiveFinite(halfLengthZ, "blm3 (half length in z)", name);
142 G4Box* shape = new G4Box(name + "_solid", halfLengthX, halfLengthY, halfLengthZ);
143 BDSExtent ext = BDSExtent(halfLengthX, halfLengthY, halfLengthZ);
144 return CommonConstruction(name, material, shape, ext);
145}
146
148 const G4String& material,
149 G4double radius)
150{
151 PositiveFinite(radius, "blm1 (radius)", name);
152 G4Orb* shape = new G4Orb(name + "_solid", radius);
153 BDSExtent ext = BDSExtent(radius, radius, radius);
154 return CommonConstruction(name, material, shape, ext);
155}
156
158 const G4String& material,
159 G4VSolid* shape,
160 const BDSExtent& extent)
161{
162 G4Material* mat = BDSMaterials::Instance()->GetMaterial(material);
163 G4LogicalVolume* lv = new G4LogicalVolume(shape, mat, name + "_lv");
164 lv->SetUserLimits(BDSGlobalConstants::Instance()->DefaultUserLimits());
165 BDSBLM* result = new BDSBLM(shape, lv, extent);
166 return result;
167}
BDSBLM * BuildBLMCylinder(const G4String &name, const G4String &material, G4double halfLength, G4double radius)
Build the geometry for a cylinder.
BDSBLM * BuildBLMCube(const G4String &name, const G4String &material, G4double halfLengthX, G4double halfLengthY, G4double halfLengthZ)
Build the geometry for a cube.
BDSBLM * CommonConstruction(const G4String &name, const G4String &material, G4VSolid *shape, const BDSExtent &extent)
BDSBLM * BuildBLMSphere(const G4String &name, const G4String &material, G4double radius)
Build the geometry for a sphere.
A beam loss monitor.
Definition: BDSBLM.hh:42
void SetBias(const G4String &biasIn)
Definition: BDSBLM.hh:62
General exception with possible name of object and message.
Definition: BDSException.hh:35
Holder for +- extents in 3 dimensions.
Definition: BDSExtent.hh:39
G4LogicalVolume * GetContainerLogicalVolume() const
Accessor - see member for more info.
virtual std::set< G4LogicalVolume * > GetAllLogicalVolumes() const
Access all logical volumes belonging to this component.
A loaded piece of externally provided geometry.
BDSGeometryExternal * BuildGeometry(G4String componentName, const G4String &formatAndFilePath, std::map< G4String, G4Colour * > *colourMapping=nullptr, G4bool autoColour=true, G4double suggestedLength=0, G4double suggestedHorizontalWidth=0, std::vector< G4String > *namedVacuumVolumes=nullptr, G4bool makeSensitive=true, BDSSDType sensitivityType=BDSSDType::energydep, BDSSDType vacuumSensitivityType=BDSSDType::energydepvacuum, G4bool stripOuterVolumeAndMakeAssembly=false, G4UserLimits *userLimitsToAttachToAllLVs=nullptr, G4bool dontReloadGeometry=false)
static BDSGeometryFactory * Instance()
Singleton accessor.
static BDSGlobalConstants * Instance()
Access method.
static BDSMaterials * Instance()
Singleton pattern access.
Definition: BDSMaterials.cc:39
G4Material * GetMaterial(G4String material) const
Get material by name.
Improve type-safety of native enum data type in C++.
type underlying() const
return underlying value (can be used in switch statement)
blm for parser.
Definition: blmplacement.h:39
std::string name
Name of this samplerplacement.
Definition: blmplacement.h:41
BDSBLMType DetermineBLMType(G4String blmType)
function that gives corresponding enum value for string (case-insensitive)
Definition: BDSBLMType.cc:39
G4bool IsFinite(G4double value, G4double tolerance=std::numeric_limits< double >::epsilon())