BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
BDSGeometryFactoryBase.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 "BDSColours.hh"
20#include "BDSColourFromMaterial.hh"
21#include "BDSGeometryExternal.hh"
22#include "BDSGeometryFactoryBase.hh"
23#include "BDSGlobalConstants.hh"
24#include "BDSUtilities.hh"
25
26#include "globals.hh"
27#include "G4AssemblyVolume.hh"
28#include "G4Colour.hh"
29#include "G4LogicalVolume.hh"
30#include "G4RotationMatrix.hh"
31#include "G4String.hh"
32#include "G4VisAttributes.hh"
33#include "G4VPhysicalVolume.hh"
34
35#include <algorithm>
36#include <map>
37#include <set>
38#include <vector>
39
40BDSGeometryFactoryBase::BDSGeometryFactoryBase()
41{
43}
44
45BDSGeometryFactoryBase::~BDSGeometryFactoryBase()
46{;}
47
48std::set<G4VisAttributes*> BDSGeometryFactoryBase::ApplyColourMapping(std::set<G4LogicalVolume*>& lvsIn,
49 std::map<G4String, G4Colour*>* mapping,
50 G4bool autoColour,
51 const G4String& prefixToStripFromName)
52{
53 std::set<G4VisAttributes*> visAttributes; // empty set
54
55 // no mapping, just return.
56 if (!mapping && !autoColour)
57 {return visAttributes;}
58
59 if (mapping)
60 {
61 if (mapping->size() == 1)
62 {// only one colour for all - simpler
63 G4VisAttributes* vis = new G4VisAttributes(*BDSColours::Instance()->GetColour("gdml"));
64 vis->SetVisibility(true);
65 visAttributes.insert(vis);
66 for (auto lv : lvsIn)
67 {lv->SetVisAttributes(*vis);}
68 return visAttributes;
69 }
70
71 // else iterate over all lvs and required vis attributes
72 // prepare required vis attributes
73 std::map<G4String, G4VisAttributes*> attMap;
74 for (const auto& it : *mapping)
75 {
76 G4VisAttributes* vis = new G4VisAttributes(*(it.second));
77 vis->SetVisibility(true);
78 visAttributes.insert(vis);
79 attMap[it.first] = vis;
80 }
81
82 for (auto lv : lvsIn)
83 {// iterate over all volumes
84 const G4String& name = lv->GetName();
85 for (const auto& it : attMap)
86 {// iterate over all mappings to find first one that matches substring
87 if (BDS::StrContains(name, it.first))
88 {
89 lv->SetVisAttributes(it.second);
90 break;
91 }
92 }
93 }
94 }
95
96 if (autoColour)
97 {
98 for (auto lv : lvsIn)
99 {
100 const G4VisAttributes* existingVis = lv->GetVisAttributes();
101 if (!existingVis)
102 {
103 G4Colour* c = BDSColourFromMaterial::Instance()->GetColour(lv->GetMaterial(), prefixToStripFromName);
104 G4VisAttributes* vis = new G4VisAttributes(*c);
105 vis->SetVisibility(true);
106 visAttributes.insert(vis);
107 lv->SetVisAttributes(vis);
108 }
109 }
110 }
111
112 return visAttributes;
113}
114
115void BDSGeometryFactoryBase::ApplyUserLimits(const std::set<G4LogicalVolume*>& lvsIn,
116 G4UserLimits* userLimits)
117{
118 for (auto& lv : lvsIn)
119 {lv->SetUserLimits(userLimits);}
120}
121
123{
124 CleanUpBase();
125}
126
127G4String BDSGeometryFactoryBase:: PreprocessedName(const G4String& objectName,
128 const G4String& /*acceleratorComponentName*/) const
129{return objectName;}
130
131std::set<G4LogicalVolume*> BDSGeometryFactoryBase::GetVolumes(const std::set<G4LogicalVolume*>& allLVs,
132 std::vector<G4String>* volumeNames,
133 G4bool preprocessFile,
134 const G4String& componentName) const
135{
136 if (!volumeNames)
137 {return std::set<G4LogicalVolume*>();}
138
139 std::vector<G4String> expectedVolumeNames;
140 if (preprocessFile)
141 {
142 // transform the names to those the preprocessor would produce
143 expectedVolumeNames.resize(volumeNames->size());
144 std::transform(volumeNames->begin(),
145 volumeNames->end(),
146 expectedVolumeNames.begin(),
147 [&](const G4String& n){return PreprocessedName(n, componentName);});
148 }
149 else
150 {expectedVolumeNames = *volumeNames;}
151
152 std::set<G4LogicalVolume*> volsMatched;
153 for (const auto& en : expectedVolumeNames)
154 {
155 for (auto lv : allLVs)
156 {
157 if (lv->GetName() == en)
158 {volsMatched.insert(lv);}
159 }
160 }
161 return volsMatched;
162}
163
165{
167 xmin = 0;
168 xmax = 0;
169 ymin = 0;
170 ymax = 0;
171 zmin = 0;
172 zmax = 0;
173}
174
176{
177 xmin = std::min(ext.XNeg(), xmin);
178 xmax = std::max(ext.XPos(), xmax);
179 ymin = std::min(ext.YNeg(), ymin);
180 ymax = std::max(ext.YPos(), ymax);
181 zmin = std::min(ext.ZNeg(), zmin);
182 zmax = std::max(ext.ZPos(), zmax);
183}
184
185void BDSGeometryFactoryBase::ExpandExtent(G4double x0, G4double rx,
186 G4double y0, G4double ry,
187 G4double z0, G4double rz)
188{
189 xmin = std::min(x0-rx, xmin);
190 xmax = std::max(x0+rx, xmax);
191 ymin = std::min(y0-ry, ymin);
192 ymax = std::max(y0+ry, ymax);
193 zmin = std::min(z0-rz, zmin);
194 zmax = std::max(z0+rz, zmax);
195}
196
197void BDSGeometryFactoryBase::ExpandExtent(G4double x0, G4double xLow, G4double xHigh,
198 G4double y0, G4double yLow, G4double yHigh,
199 G4double z0, G4double zLow, G4double zHigh)
200{
201 xmin = std::min(x0+xLow, xmin);
202 xmax = std::max(x0+xHigh, xmax);
203 ymin = std::min(y0+yLow, ymin);
204 ymax = std::max(y0+yHigh, ymax);
205 zmin = std::min(z0+zLow, zmin);
206 zmax = std::max(z0+zHigh, zmax);
207}
static BDSColourFromMaterial * Instance()
Singleton pattern.
G4Colour * GetColour(const G4Material *material, const G4String &prefixToStripFromName="")
Get colour from name.
static BDSColours * Instance()
singleton pattern
Definition: BDSColours.cc:33
Holder for +- extents in 3 dimensions.
Definition: BDSExtent.hh:39
G4double XPos() const
Accessor.
Definition: BDSExtent.hh:66
G4double ZPos() const
Accessor.
Definition: BDSExtent.hh:70
G4double XNeg() const
Accessor.
Definition: BDSExtent.hh:65
G4double YNeg() const
Accessor.
Definition: BDSExtent.hh:67
G4double YPos() const
Accessor.
Definition: BDSExtent.hh:68
G4double ZNeg() const
Accessor.
Definition: BDSExtent.hh:69
virtual void FactoryBaseCleanUp()
Empty containers for next use - factories are never deleted so can't rely on scope.
G4double ymax
Extent in one dimension.
void ExpandExtent(const BDSExtent &extent)
Expand the acuumulated extents using a (possibly asymmetric) extent instance.
virtual G4String PreprocessedName(const G4String &objectName, const G4String &acceleratorComponentName) const
std::set< G4LogicalVolume * > GetVolumes(const std::set< G4LogicalVolume * > &allLVs, std::vector< G4String > *volumeNames, G4bool preprocessGDML, const G4String &componentName) const
Get the volumes that match the name. Volume names are matched exactly and are case sensitive.
virtual void CleanUp()
Virtual clean up that derived classes can override that calls CleanUpBase().
G4double xmin
Extent in one dimension.
G4double ymin
Extent in one dimension.
G4double zmin
Extent in one dimension.
G4double zmax
Extent in one dimension.
G4double xmax
Extent in one dimension.
virtual void ApplyUserLimits(const std::set< G4LogicalVolume * > &lvsIn, G4UserLimits *userLimits)
Attach a set of user limits to every logical volume supplied.
virtual std::set< G4VisAttributes * > ApplyColourMapping(std::set< G4LogicalVolume * > &lvs, std::map< G4String, G4Colour * > *mapping, G4bool autoColour, const G4String &preprocessPrefix="")
G4bool StrContains(const G4String &str, const G4String &test)
Utility function to simplify lots of syntax changes for pedantic g4 changes.
Definition: BDSUtilities.cc:66