BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSDetectorConstruction.hh
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#ifndef BDSDETECTORCONSTRUCTION_H
20#define BDSDETECTORCONSTRUCTION_H
21
22#include "BDSExtent.hh"
23#include "BDSExtentGlobal.hh"
24
25#include "globals.hh" // geant4 types / globals
26#include "G4Transform3D.hh"
27#include "G4Version.hh"
28#include "G4VUserDetectorConstruction.hh"
29
30#include <list>
31#include <map>
32#include <set>
33#include <string>
34#include <vector>
35
36class G4LogicalVolume;
37class G4Region;
38class G4VPhysicalVolume;
39
40namespace GMAD {
41 class BLMPlacement;
42 struct Element;
43 template<typename T> class FastList;
44 class Placement;
45 class Query;
46 class SamplerPlacement;
47 class ScorerMesh;
48}
49
51class BDSBeamline;
52class BDSBeamlineSet;
54class BDSFieldObjects;
57class BDSSamplerInfo;
58
59#if G4VERSION_NUMBER > 1009
61#endif
62
74class BDSDetectorConstruction: public G4VUserDetectorConstruction
75{
76public:
77 explicit BDSDetectorConstruction(BDSComponentFactoryUser* userComponentFactoryIn = nullptr);
79
83
86 virtual G4VPhysicalVolume* Construct();
87
89 virtual void ConstructSDandField();
90
94 void BuildPhysicsBias();
95
97 inline void SetDesignParticle(const BDSParticleDefinition* defIn) {designParticle = defIn;}
98
101
103 const std::vector<BDSFieldQueryInfo*>& FieldQueries() const {return fieldQueries;}
104
113 static void PlaceBeamlineInWorld(BDSBeamline* beamline,
114 G4VPhysicalVolume* containerPV,
115 G4bool checkOverlaps = false,
116 G4bool setRegions = false,
117 G4bool registerInfo = false,
118 G4bool useCLPlacementTransform = false,
119 G4bool useIncrementalCopyNumbers = false,
120 G4bool registerPlacementNamesForOutput = false);
121
126 static G4Transform3D CreatePlacementTransform(const GMAD::Placement& placement,
127 const BDSBeamline* beamLine,
128 G4double* S = nullptr,
129 BDSExtent* placementExtent = nullptr,
130 const G4String& objectTypeForErrorMsg = "placement");
131
132 // Create a scorermesh placement transform. Turns the scorermesh into a
134 static G4Transform3D CreatePlacementTransform(const GMAD::ScorerMesh& scorerMesh,
135 const BDSBeamline* beamLine,
136 G4double* S = nullptr);
137
140 static G4Transform3D CreatePlacementTransform(const GMAD::SamplerPlacement& samplerPlacement,
141 const BDSBeamline* beamLine,
142 G4double* S = nullptr);
143
145 static G4Transform3D CreatePlacementTransform(const GMAD::BLMPlacement& blmPlacement,
146 const BDSBeamline* beamLine,
147 G4double* S = nullptr,
148 BDSExtent* blmExtent = nullptr);
149
151 static G4Transform3D CreatePlacementTransform(const GMAD::Query& queryPlacement,
152 const BDSBeamline* beamLine,
153 G4double* S = nullptr);
154
156 static G4ThreeVector SideToLocalOffset(const GMAD::Placement& placement,
157 const BDSBeamline* beamLine,
158 const BDSExtent& placementExtent);
159
162 static BDSSamplerInfo* BuildSamplerInfo(const GMAD::Element* element);
163
166 G4bool BuildSamplerWorld() const {return nSamplers > 0;}
167
168 G4bool BuildPlacementFieldsWorld() const {return buildPlacementFieldsWorld;}
169
171 static std::vector<BDSFieldQueryInfo*> PrepareFieldQueries(const BDSBeamline* mainBeamline);
172
173private:
177
181
184
186 void InitialiseRegions();
187
189 void InitialiseApertures();
190
192 void BuildBeamlines();
193
198 const G4String& name,
199 const G4Transform3D& initialTransform = G4Transform3D(),
200 G4double initialS = 0.0,
201 G4bool beamlineIsCircular = false,
202 G4bool isPlacementBeamline= false);
203
205 void BuildTunnel();
206
209 G4VPhysicalVolume* BuildWorld();
210
212 void ComponentPlacement(G4VPhysicalVolume* worldPV);
213
216 G4bool UnsuitableFirstElement(std::list<GMAD::Element>::const_iterator element);
217
220
224
227
231
232#if G4VERSION_NUMBER > 1009
234 BDSBOptrMultiParticleChangeCrossSection* BuildCrossSectionBias(const std::list<std::string>& biasList,
235 const std::list<std::string>& defaultBias,
236 const G4String& elementName);
237
240
242 void VerboseSensitivity() const;
244 void PrintSensitiveDetectorsOfLV(const G4LogicalVolume* lv, G4int currentDepth) const;
245
247 std::vector<BDSBOptrMultiParticleChangeCrossSection*> biasObjects;
248 std::map<G4String, BDSBOptrMultiParticleChangeCrossSection*> biasSetObjects;
249#endif
250
251#ifdef BDSDEBUG
252 bool debug = true;
253#else
254 bool debug = false;
255#endif
256
258 G4bool verbose;
261
264
266 std::vector<BDSFieldObjects*> fields;
267
268 G4bool circular;
274
275 BDSComponentFactoryUser* userComponentFactory;
276
277 G4int nSamplers;
278 G4bool buildPlacementFieldsWorld;
279
281 std::set<G4LogicalVolume*> worldContentsLogicalVolumes;
282 std::set<G4LogicalVolume*> worldVacuumLogicalVolumes;
283 G4LogicalVolume* worldLogicalVolume;
284
285 std::vector<BDSFieldQueryInfo*> fieldQueries;
286
287 // for developer checks only
288#ifdef BDSCHECKUSERLIMITS
289 void PrintUserLimitsSummary(const G4VPhysicalVolume* world) const;
290 void PrintUserLimitsPV(const G4VPhysicalVolume* aPV, G4double globalMinEK) const;
291#endif
292};
293
294#endif
295
A holder class for all representations of the accelerator model created in BDSIM.
Simple struct to return a beamline plus associated beam lines.
A vector of BDSBeamlineElement instances - a beamline.
Definition: BDSBeamline.hh:61
Factory for user specified accelerator components.
Class that constructs a Geant4 model of an accelerator.
G4bool circular
Whether or not we're building a circular machine.
BDSBeamlineSet BuildBeamline(const GMAD::FastList< GMAD::Element > &beamLine, const G4String &name, const G4Transform3D &initialTransform=G4Transform3D(), G4double initialS=0.0, G4bool beamlineIsCircular=false, G4bool isPlacementBeamline=false)
void SetDesignParticle(const BDSParticleDefinition *defIn)
Set the design particle definition.
BDSDetectorConstruction & operator=(const BDSDetectorConstruction &)=delete
assignment and copy constructor not implemented nor used
BDSBOptrMultiParticleChangeCrossSection * BuildCrossSectionBias(const std::list< std::string > &biasList, const std::list< std::string > &defaultBias, const G4String &elementName)
Function that creates physics biasing cross section.
void InitialiseRegions()
Create and set parameters for various G4Regions.
std::vector< BDSBOptrMultiParticleChangeCrossSection * > biasObjects
List of bias objects - for memory management.
BDSExtent CalculateExtentOfScorerMesh(const GMAD::ScorerMesh &sm) const
Calculate local extent of scorer mesh in 3D.
BDSAcceleratorModel * acceleratorModel
Accelerator model pointer.
static G4ThreeVector SideToLocalOffset(const GMAD::Placement &placement, const BDSBeamline *beamLine, const BDSExtent &placementExtent)
Attach component with extent2 to component with extent1 with placement.
G4bool verbose
Variable copied from global constants.
static std::vector< BDSFieldQueryInfo * > PrepareFieldQueries(const BDSBeamline *mainBeamline)
Prepare field queries from parser information.
G4VPhysicalVolume * BuildWorld()
void VerboseSensitivity() const
Print out the sensitivity of every single volume so far constructed in the world.
void InitialiseApertures()
Create all aperture definitions from parser and store in BDSAcceleratorModel.
void CountPlacementFields()
Count number of fields required for placements.
G4bool canSampleAngledFaces
Whether the integrator set permits sampling elements with angled faces.
BDSExtent worldExtent
Record of the world extent.
std::set< G4LogicalVolume * > worldContentsLogicalVolumes
Cache of possibly loaded logical volumes from a world geometry file - used for biasing.
static BDSSamplerInfo * BuildSamplerInfo(const GMAD::Element *element)
BDSExtent CalculateExtentOfSamplerPlacement(const GMAD::SamplerPlacement &sp) const
Calculate local extent of custom user sampler.
void ConstructScoringMeshes()
Construct scoring meshes.
BDSExtentGlobal CalculateExtentOfScorerMeshes(const BDSBeamline *bl) const
G4bool checkOverlaps
Variable copied from global constants.
G4int nSamplers
Count of number of samplers to be built.
void ComponentPlacement(G4VPhysicalVolume *worldPV)
Place beam line, tunnel beam line, end pieces and placements in world.
const std::vector< BDSFieldQueryInfo * > & FieldQueries() const
Access vector of query objects.
void BuildBeamlines()
Build the main beam line and then any other required beam lines.
const BDSParticleDefinition * designParticle
Particle definition all components are built w.r.t. Includes rigidity etc.
void BuildTunnel()
Build the tunnel around the already constructed flat beam line.
BDSExtent WorldExtent() const
Public access to the world extent.
virtual void ConstructSDandField()
Construct sensitive detectors and fields.
std::vector< BDSFieldObjects * > fields
All fields.
void PrintSensitiveDetectorsOfLV(const G4LogicalVolume *lv, G4int currentDepth) const
Recursive function to print out each sensitive detector name.
static void PlaceBeamlineInWorld(BDSBeamline *beamline, G4VPhysicalVolume *containerPV, G4bool checkOverlaps=false, G4bool setRegions=false, G4bool registerInfo=false, G4bool useCLPlacementTransform=false, G4bool useIncrementalCopyNumbers=false, G4bool registerPlacementNamesForOutput=false)
virtual G4VPhysicalVolume * Construct()
static G4Transform3D CreatePlacementTransform(const GMAD::Placement &placement, const BDSBeamline *beamLine, G4double *S=nullptr, BDSExtent *placementExtent=nullptr, const G4String &objectTypeForErrorMsg="placement")
G4bool UnsuitableFirstElement(std::list< GMAD::Element >::const_iterator element)
BDSExtentGlobal CalculateExtentOfSamplerPlacements(const BDSBeamline *beamLine) const
Holder for +- extents in 3 dimensions with a rotation and translation.
Holder for +- extents in 3 dimensions.
Definition: BDSExtent.hh:39
A holder for all the Geant4 field related objects.
Holder class for all information required for a field query.
Wrapper for particle definition.
All info required to build a sampler but not place it.
blm for parser.
Definition: blmplacement.h:39
List with Efficient Lookup.
Definition: fastlist.h:42
Placement class for parser.
Definition: placement.h:41
Query structure class for parser.
Definition: query.h:37
Sampler placement class for parser.
ScorerMesh class for parser.
Definition: scorermesh.h:38
Parser namespace for GMAD language. Combination of Geant4 and MAD.
Element class.
Definition: element.h:43