11 #include <framework/logging/Logger.h>
12 #include <framework/gearbox/GearDir.h>
13 #include <framework/database/IntervalOfValidity.h>
14 #include <geometry/GeometryManager.h>
15 #include <geometry/Materials.h>
16 #include <geometry/CreatorManager.h>
17 #include <geometry/CreatorBase.h>
18 #include <geometry/utilities.h>
19 #include <geometry/dbobjects/GeoConfiguration.h>
20 #include <geometry/bfieldmap/BFieldMap.h>
21 #include <geometry/bfieldmap/BFieldFrameworkInterface.h>
22 #include <framework/dbobjects/MagneticField.h>
23 #include <framework/database/DBStore.h>
26 #include "G4ThreeVector.hh"
27 #include "G4LogicalVolume.hh"
28 #include "G4PVPlacement.hh"
30 #include "G4RunManager.hh"
31 #include "G4GeometryManager.hh"
32 #include "G4PhysicalVolumeStore.hh"
33 #include "G4LogicalVolumeStore.hh"
34 #include "G4SolidStore.hh"
35 #include "G4RegionStore.hh"
36 #include "G4SurfaceProperty.hh"
37 #include "G4LogicalBorderSurface.hh"
38 #include "G4LogicalSkinSurface.hh"
39 #include "G4VisAttributes.hh"
40 #include "G4VoxelLimits.hh"
43 #include "Geant4GM/volumes/Factory.h"
44 #include "RootGM/volumes/Factory.h"
45 #include "VGM/volumes/IPlacement.h"
46 #include "TGeoManager.h"
57 using namespace gearbox;
69 double getTopMinSize(
const std::string& name, EAxis axis, G4LogicalVolume* volume,
double current)
71 G4VoxelLimits dummyLimits;
73 for (
int i = 0; i < volume->GetNoDaughters(); ++i) {
74 G4VPhysicalVolume* daughter = volume->GetDaughter(i);
75 G4AffineTransform trans(daughter->GetRotation(), daughter->GetTranslation());
76 double vmin{0}, vmax{0};
77 daughter->GetLogicalVolume()->GetSolid()->CalculateExtent((EAxis) axis, dummyLimits, trans, vmin, vmax);
78 extent = std::max({extent, std::fabs(vmin), std::fabs(vmax)});
80 B2DEBUG(100,
"Current global volume size in " << name <<
": " << std::fixed
81 << std::setprecision(std::numeric_limits<double>::max_digits10)
82 << current <<
", needed: " << extent);
83 if (current < extent && (extent - current) > Unit::um) {
85 B2WARNING(
"Global volume not large enough in " << name <<
" direction, enlarging from "
86 << current <<
" mm to " << extent <<
" mm");
88 B2DEBUG(10,
"Setting global volume size to " << name <<
"= +-" << extent <<
" mm");
104 void GeometryManager::clear()
106 B2DEBUG(50,
"Cleaning up Geometry");
107 for (G4VisAttributes* visAttr : m_VisAttributes)
delete visAttr;
108 m_VisAttributes.clear();
109 for (
CreatorBase* creator : m_creators)
delete creator;
111 m_topVolume =
nullptr;
113 G4GeometryManager::GetInstance()->OpenGeometry();
114 G4PhysicalVolumeStore::Clean();
115 G4LogicalVolumeStore::Clean();
116 G4SolidStore::Clean();
117 G4LogicalBorderSurface::CleanSurfaceTable();
118 G4LogicalSkinSurface::CleanSurfaceTable();
119 G4SurfaceProperty::CleanSurfacePropertyTable();
121 Materials::getInstance().clear();
126 std::string detectorName;
128 detectorName = detectorDir.
getString(
"Name");
129 B2DEBUG(50,
"Creating geometry for detector: " << detectorName);
130 }
catch (gearbox::PathEmptyError& e) {
131 B2FATAL(
"Could not read detector name, make sure gearbox is connected and "
132 << detectorDir.
getPath() <<
" points to the geometry description");
135 const double globalWidth = detectorDir.
getLength(
"Global/width", 0);
136 const double globalHeight = detectorDir.
getLength(
"Global/height", 0);
137 const double globalLength = detectorDir.
getLength(
"Global/length", 0);
138 const std::string globalMaterial = detectorDir.
getString(
"Global/material",
"Air");
140 GeoConfiguration config(detectorName, globalWidth, globalHeight, globalLength, globalMaterial);
143 Materials& materials = Materials::getInstance();
145 for (
const GearDir& mat : matlist.getNodes(
"Material")) {
147 config.addMaterial(material);
151 std::set<std::string> componentNames = m_components;
152 std::set<std::string> excludedNames = m_excluded;
153 std::set<std::string> additionalNames = m_additional;
154 if (!m_components.empty() && !m_additional.empty()) {
155 B2WARNING(
"Additional components are ignored when a list of components is provided.");
159 for (
const GearDir& component : detectorDir.
getNodes(
"DetectorComponent")) {
163 name = component.getString(
"@name");
164 creatorName = component.getString(
"Creator");
165 }
catch (gearbox::PathEmptyError& e) {
166 B2ERROR(
"Could not find required element Name or Creator for " << component.getPath());
169 const bool isDefault = component.getBool(
"@isDefault",
true);
172 componentNames.erase(name);
173 excludedNames.erase(name);
174 additionalNames.erase(name);
175 if (!m_components.empty() && m_components.count(name) == 0) {
176 B2DEBUG(50,
"DetectorComponent " << name <<
" not in list of components, skipping");
179 if (m_components.empty() && !isDefault) {
180 if (m_additional.count(name) == 0) {
181 B2DEBUG(50,
"DectorComponent " << name <<
" is not enabled by default, skipping");
184 B2DEBUG(50,
"DectorComponent " << name <<
" is enabled in addition to the default components");
187 if (!m_excluded.empty() && m_excluded.count(name) > 0) {
188 B2DEBUG(10,
"DetectorComponent " << name <<
" in list of excluded components, skipping");
192 string libraryName = component.getString(
"Creator/@library",
"");
194 CreatorBase* creator = CreatorManager::getCreator(creatorName, libraryName);
198 B2ERROR(
"Could not load creator " << creatorName <<
" from " << libraryName);
201 config.addComponent({name, creatorName, libraryName});
207 auto checkRemaining = [](
const std::string & type,
const std::set<std::string>& remainingNames) {
208 for (
const std::string& name : remainingNames) {
209 B2ERROR(
"Geometry '" << name <<
"' is specified in list of "
210 << type <<
" but could not be found");
214 checkRemaining(
"components", componentNames);
215 checkRemaining(
"excluded components", excludedNames);
216 checkRemaining(
"additional components", additionalNames);
224 createGeometry(config, type,
false);
234 BFieldMap::Instance().clear();
238 DBStore::Instance().addConstantOverride(
"MagneticField", fieldmap,
false);
242 G4RunManager* runManager = G4RunManager::GetRunManager();
243 if (runManager) runManager->ReinitializeGeometry(
true,
true);
246 B2DEBUG(10,
"Creating geometry for detector: " << config.getName());
249 Materials& materials = Materials::getInstance();
250 for (
const GeoMaterial& mat : config.getMaterials()) {
256 G4Material* top_mat = Materials::get(config.getGlobalMaterial());
257 double xHalfLength = config.getGlobalWidth() / Unit::cm * CLHEP::cm;
258 double yHalfLength = config.getGlobalHeight() / Unit::cm * CLHEP::cm;
259 double zHalfLength = config.getGlobalLength() / Unit::cm * CLHEP::cm;
263 G4Box* top_box =
new G4Box(
"Top", xHalfLength > 0 ? xHalfLength : 1,
264 yHalfLength > 0 ? yHalfLength : 1,
265 zHalfLength > 0 ? zHalfLength : 1);
266 G4LogicalVolume* top_log =
new G4LogicalVolume(top_box, top_mat,
"Top",
nullptr,
nullptr,
nullptr);
268 m_topVolume =
new G4PVPlacement(
nullptr, G4ThreeVector(), top_log,
"Top",
nullptr,
false, 0);
269 B2DEBUG(10,
"Created top volume with x= +-" << config.getGlobalWidth() <<
" cm, y= +-"
270 << config.getGlobalHeight() <<
" cm, z= +-" << config.getGlobalLength() <<
" cm");
273 auto getDensityScale = [
this](
const std::string & name) {
274 std::optional<double> scale;
276 if (
auto it = m_densityScaling.find(name); it != m_densityScaling.end()) {
281 auto globalScale = getDensityScale(
"*");
283 for (
const GeoComponent& component : config.getComponents()) {
284 CreatorBase* creator = CreatorManager::getCreator(component.getCreator(), component.getLibrary());
287 auto componentScale = getDensityScale(component.getName());
288 if (componentScale or globalScale) {
289 double scale = globalScale ? *globalScale : 1.0;
290 if (componentScale) scale *= *componentScale;
291 Materials::getInstance().setDensityScale(scale);
294 int oldSolids = G4SolidStore::GetInstance()->size();
295 int oldLogical = G4LogicalVolumeStore::GetInstance()->size();
296 int oldPhysical = G4PhysicalVolumeStore::GetInstance()->size();
298 if (!useDB)
throw CreatorBase::DBNotImplemented();
299 creator->
createFromDB(component.getName(), *top_log, type);
300 B2DEBUG(50,
"called creator " << component.getCreator() <<
" to create component " << component.getName() <<
" from DB");
302 }
catch (CreatorBase::DBNotImplemented& e) {
303 GearDir parameters = Gearbox::getInstance().getDetectorComponent(component.getName());
304 creator->
create(parameters, *top_log, type);
305 B2DEBUG(50,
"called creator " << component.getCreator() <<
" to create component " << component.getName() <<
" from Gearbox");
307 int newSolids = G4SolidStore::GetInstance()->size() - oldSolids;
308 int newLogical = G4LogicalVolumeStore::GetInstance()->size() - oldLogical;
309 int newPhysical = G4PhysicalVolumeStore::GetInstance()->size() - oldPhysical;
310 B2DEBUG(50,
"DetectorComponent " << component.getName() <<
" created " << newSolids
311 <<
" solids, " << newLogical <<
" logical volumes and "
312 << newPhysical <<
" physical volumes");
313 m_creators.push_back(creator);
315 Materials::getInstance().resetDensityScale();
316 if (m_assignRegions) {
318 G4Region* region {
nullptr};
321 for (
int i = 0; i < top_log->GetNoDaughters(); ++i) {
322 G4LogicalVolume* vol = top_log->GetDaughter(i)->GetLogicalVolume();
324 if (!vol->GetRegion()) {
326 if (!region) region = G4RegionStore::GetInstance()->FindOrCreateRegion(component.getName());
327 vol->SetRegion(region);
329 region->AddRootLogicalVolume(vol);
336 int newSolids = G4SolidStore::GetInstance()->size();
337 int newLogical = G4LogicalVolumeStore::GetInstance()->size();
338 int newPhysical = G4PhysicalVolumeStore::GetInstance()->size();
339 B2DEBUG(50,
"Created a total of " << newSolids <<
" solids, " << newLogical
340 <<
" logical volumes and " << newPhysical <<
" physical volumes");
343 top_box->SetXHalfLength(getTopMinSize(
"x", kXAxis, top_log, xHalfLength));
344 top_box->SetYHalfLength(getTopMinSize(
"y", kYAxis, top_log, yHalfLength));
345 top_box->SetZHalfLength(getTopMinSize(
"z", kZAxis, top_log, zHalfLength));
347 B2DEBUG(50,
"Optimizing geometry and creating lookup tables ...");
348 G4GeometryManager::GetInstance()->CloseGeometry(
true, LogSystem::Instance().isLevelEnabled(LogConfig::c_Debug, 200, PACKAGENAME()));
351 void GeometryManager::createTGeoRepresentation()
354 B2ERROR(
"No Geometry found, please create a geometry before converting it to ROOT::TGeo");
358 Geant4GM::Factory g4Factory;
359 if (LogSystem::Instance().isLevelEnabled(LogConfig::c_Debug, 200, PACKAGENAME())) {
360 g4Factory.SetDebug(1);
362 g4Factory.Import(m_topVolume);
363 RootGM::Factory rtFactory;
364 rtFactory.SetIgnore(1);
365 if (LogSystem::Instance().isLevelEnabled(LogConfig::c_Debug, 200, PACKAGENAME())) {
366 rtFactory.SetDebug(1);
368 g4Factory.Export(&rtFactory);
369 gGeoManager->CloseGeometry();
370 delete g4Factory.Top();
371 delete rtFactory.Top();
374 G4VisAttributes* GeometryManager::newVisAttributes()
376 m_VisAttributes.push_back(
new G4VisAttributes());
377 return m_VisAttributes.back();