Skip to content

Commit

Permalink
Fix loading of GDML files with reflection (#1626)
Browse files Browse the repository at this point in the history
* Instantiate UI pointers to maintain logging consistency
* Master UI pointer is only set after MT manager starts
* Add optional SD reading to GDML loader
* Update scoped mem name
* Add option to remove pointer names more carefully
* Move and rename load/save gdml methods
* Combine orange Geant helper loaders
* Add GDML loader test
* Change default GDML load to excise
* Attach reflected SDs
* Address feedback
  • Loading branch information
sethrj authored Feb 17, 2025
1 parent 83e6546 commit 3fe3999
Show file tree
Hide file tree
Showing 27 changed files with 798 additions and 284 deletions.
76 changes: 23 additions & 53 deletions app/celer-g4/DetectorConstruction.cc
Original file line number Diff line number Diff line change
Expand Up @@ -10,8 +10,6 @@
#include <G4ChordFinder.hh>
#include <G4Exception.hh>
#include <G4FieldManager.hh>
#include <G4GDMLAuxStructType.hh>
#include <G4GDMLParser.hh>
#include <G4LogicalVolume.hh>
#include <G4MagneticField.hh>
#include <G4SDManager.hh>
Expand All @@ -23,6 +21,7 @@
#include "corecel/io/Logger.hh"
#include "corecel/io/OutputRegistry.hh"
#include "corecel/math/ArrayUtils.hh"
#include "geocel/GeantGdmlLoader.hh"
#include "celeritas/Quantities.hh"
#include "celeritas/field/RZMapFieldInput.hh"
#include "celeritas/field/RZMapFieldParams.hh"
Expand Down Expand Up @@ -78,9 +77,27 @@ G4VPhysicalVolume* DetectorConstruction::Construct()
== (GlobalSetup::Instance()->input().sd_type
!= SensitiveDetectorType::none));

auto geo = this->construct_geo();
CELER_ASSERT(geo.world);
detectors_ = std::move(geo.detectors);
std::string const& filename = GlobalSetup::Instance()->GetGeometryFile();
CELER_VALIDATE(!filename.empty(),
<< "no GDML input file was specified (geometry_file)");

auto& sd = GlobalSetup::Instance()->GetSDSetupOptions();
auto loaded = [&] {
GeantGdmlLoader::Options opts;
opts.detectors = sd.enabled;
GeantGdmlLoader load(opts);
return load(filename);
}();

if (sd.enabled && loaded.detectors.empty())
{
CELER_LOG(warning)
<< R"(No sensitive detectors were found in the GDML file)";
sd.enabled = false;
}

CELER_ASSERT(loaded.world);
detectors_ = std::move(loaded.detectors);

if (!detectors_.empty()
&& GlobalSetup::Instance()->input().sd_type
Expand Down Expand Up @@ -110,54 +127,7 @@ G4VPhysicalVolume* DetectorConstruction::Construct()
CELER_ASSERT(field.along_step);
GlobalSetup::Instance()->SetAlongStepFactory(std::move(field.along_step));
mag_field_ = std::move(field.g4field);
return geo.world.release();
}

//---------------------------------------------------------------------------//
/*!
* Construct shared geometry information.
*/
auto DetectorConstruction::construct_geo() const -> GeoData
{
CELER_LOG_LOCAL(status) << "Loading detector geometry";

G4GDMLParser gdml_parser;
gdml_parser.SetStripFlag(true);

std::string const& filename = GlobalSetup::Instance()->GetGeometryFile();
CELER_VALIDATE(!filename.empty(),
<< "no GDML input file was specified (geometry_file)");
constexpr bool validate_gdml_schema = false;
gdml_parser.Read(filename, validate_gdml_schema);

auto& sd = GlobalSetup::Instance()->GetSDSetupOptions();

MapDetectors detectors;
if (sd.enabled)
{
// Find sensitive detectors
for (auto const& lv_vecaux : *gdml_parser.GetAuxMap())
{
for (G4GDMLAuxStructType const& aux : lv_vecaux.second)
{
if (aux.type == "SensDet")
{
detectors.insert({aux.value, lv_vecaux.first});
}
}
}

if (detectors.empty())
{
CELER_LOG(warning) << "No sensitive detectors were found in the "
"GDML file";
sd.enabled = false;
}
}

// Claim ownership of world volume and pass it to the caller
return {std::move(detectors),
UPPhysicalVolume{gdml_parser.GetWorldVolume()}};
return loaded.world;
}

//---------------------------------------------------------------------------//
Expand Down
7 changes: 0 additions & 7 deletions app/celer-g4/DetectorConstruction.hh
Original file line number Diff line number Diff line change
Expand Up @@ -46,17 +46,11 @@ class DetectorConstruction final : public G4VUserDetectorConstruction
private:
//// TYPES ////

using UPPhysicalVolume = std::unique_ptr<G4VPhysicalVolume>;
using MapDetectors = std::multimap<std::string, G4LogicalVolume*>;
using AlongStepFactory = SetupOptions::AlongStepFactory;
using SPMagneticField = std::shared_ptr<G4MagneticField>;
using SPSimpleCalo = std::shared_ptr<GeantSimpleCalo>;

struct GeoData
{
MapDetectors detectors;
UPPhysicalVolume world;
};
struct FieldData
{
AlongStepFactory along_step;
Expand All @@ -73,7 +67,6 @@ class DetectorConstruction final : public G4VUserDetectorConstruction

//// METHODS ////

GeoData construct_geo() const;
FieldData construct_field() const;

template<class F>
Expand Down
4 changes: 3 additions & 1 deletion doc/implementation/geant4-interface.rst
Original file line number Diff line number Diff line change
Expand Up @@ -168,7 +168,9 @@ This subsection contains details of importing Geant4 data into Celeritas.
Geant4 geometry utilities
^^^^^^^^^^^^^^^^^^^^^^^^^

.. doxygenfunction:: celeritas::load_geant_geometry_native
.. doxygenclass:: celeritas::GeantGdmlLoader
.. doxygenfunction:: celeritas::load_gdml
.. doxygenfunction:: celeritas::save_gdml
.. doxygenfunction:: celeritas::find_geant_volumes

Geant4 physics interfaces
Expand Down
5 changes: 3 additions & 2 deletions src/accel/SharedParams.cc
Original file line number Diff line number Diff line change
Expand Up @@ -43,6 +43,7 @@
#include "corecel/sys/ScopedMem.hh"
#include "corecel/sys/ScopedProfiling.hh"
#include "corecel/sys/ThreadId.hh"
#include "geocel/GeantGdmlLoader.hh"
#include "geocel/GeantUtils.hh"
#include "geocel/g4/GeantGeoParams.hh"
#include "celeritas/Types.hh"
Expand Down Expand Up @@ -519,8 +520,8 @@ void SharedParams::initialize_core(SetupOptions const& options)
"when manually loading a geometry (the "
"'geometry_file' option is also set)");

write_geant_geometry(GeantImporter::get_world_volume(),
options.geometry_output_file);
save_gdml(GeantImporter::get_world_volume(),
options.geometry_output_file);
}

CoreParams::Input params;
Expand Down
7 changes: 4 additions & 3 deletions src/celeritas/ext/GeantSetup.cc
Original file line number Diff line number Diff line change
Expand Up @@ -25,7 +25,7 @@
#include "corecel/io/ScopedTimeLog.hh"
#include "corecel/sys/ScopedMem.hh"
#include "corecel/sys/ScopedProfiling.hh"
#include "geocel/GeantGeoUtils.hh"
#include "geocel/GeantGdmlLoader.hh"
#include "geocel/GeantUtils.hh"
#include "geocel/ScopedGeantExceptionHandler.hh"
#include "geocel/ScopedGeantLogger.hh"
Expand Down Expand Up @@ -104,8 +104,9 @@ GeantSetup::GeantSetup(std::string const& gdml_filename, Options options)
{
CELER_LOG(status) << "Initializing Geant4 geometry and physics list";

// Load GDML and save a copy of the pointer
world_ = load_geant_geometry_native(gdml_filename);
// Load GDML and reference the world pointer
// TODO: pass GdmlLoader options through SetupOptions
world_ = load_gdml(gdml_filename);
CELER_ASSERT(world_);

// Construct the geometry
Expand Down
1 change: 1 addition & 0 deletions src/geocel/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -43,6 +43,7 @@ if(CELERITAS_USE_Geant4)
celeritas_add_object_library(geocel_geant4
GeantUtils.cc
GeantGeoUtils.cc
GeantGdmlLoader.cc
ScopedGeantExceptionHandler.cc
ScopedGeantLogger.cc
g4/GeantGeoParams.cc
Expand Down
177 changes: 177 additions & 0 deletions src/geocel/GeantGdmlLoader.cc
Original file line number Diff line number Diff line change
@@ -0,0 +1,177 @@
//------------------------------- -*- C++ -*- -------------------------------//
// Copyright Celeritas contributors: see top-level COPYRIGHT file for details
// SPDX-License-Identifier: (Apache-2.0 OR MIT)
//---------------------------------------------------------------------------//
//! \file geocel/GeantGdmlLoader.cc
//---------------------------------------------------------------------------//
#include "GeantGdmlLoader.hh"

#include <regex>
#include <G4GDMLAuxStructType.hh>
#include <G4GDMLParser.hh>
#include <G4LogicalVolumeStore.hh>
#include <G4PhysicalVolumeStore.hh>
#include <G4ReflectionFactory.hh>
#include <G4SolidStore.hh>
#include <G4Version.hh>

#include "corecel/io/ScopedTimeLog.hh"
#include "corecel/sys/ScopedMem.hh"

#include "ScopedGeantExceptionHandler.hh"
#include "ScopedGeantLogger.hh"

namespace celeritas
{
namespace
{
//---------------------------------------------------------------------------//
bool search_pointer(std::string const& s, std::smatch& ptr_match)
{
// Search either for th end of the expression, or an underscore that likely
// indicates a _refl or _PV suffix
static std::regex const ptr_regex{"0x[0-9a-f]{4,16}(?=_|$)"};
return std::regex_search(s, ptr_match, ptr_regex);
}

//---------------------------------------------------------------------------//
//! Remove pointer address from inside geometry names
template<class StoreT>
void remove_pointers(StoreT& obj_store)
{
std::smatch sm;
for (auto* obj : obj_store)
{
if (!obj)
continue;

std::string const& name = obj->GetName();
if (search_pointer(name, sm))
{
std::string new_name = name.substr(0, sm.position());
new_name += name.substr(sm.position() + sm.length());
obj->SetName(new_name);
}
}
#if G4VERSION_NUMBER >= 1100
obj_store.UpdateMap(); // Added by geommng-V10-07-00
#endif
}

//---------------------------------------------------------------------------//
} // namespace

//---------------------------------------------------------------------------//
/*!
* Load a gdml input file, creating a pointer owned by Geant4.
*
* Geant4's constructors for physical/logical volumes register \c this pointers
* in the "volume stores" which can be cleared with
* \c celeritas::reset_geant_geometry .
*
* Note that material and element names (at least as
* of [email protected]) are \em always stripped: only volumes and solids keep
* their extension.
*/
auto GeantGdmlLoader::operator()(std::string const& filename) const -> Result
{
CELER_EXPECT(!filename.empty());
CELER_LOG(info) << "Loading Geant4 geometry from GDML at " << filename;

if (!G4Threading::IsMasterThread())
{
// Always-on debug assertion (not a "runtime" error but a
// subtle programming logic error that always causes a crash)
CELER_DEBUG_FAIL(
"Geant4 geometry cannot be loaded from a worker thread", internal);
}

ScopedMem record_mem("GeantGdmlLoader.load");
ScopedTimeLog scoped_time;

ScopedGeantLogger scoped_logger;
ScopedGeantExceptionHandler scoped_exceptions;

G4GDMLParser gdml_parser;
gdml_parser.SetStripFlag(opts_.pointers == PointerTreatment::truncate);

gdml_parser.Read(filename, /* validate_gdml_schema = */ false);

Result result;
result.world = gdml_parser.GetWorldVolume();

if (opts_.detectors)
{
// Find sensitive detectors
auto const* refl_factory = G4ReflectionFactory::Instance();
CELER_ASSERT(refl_factory);

for (auto&& [lv, vecaux] : *gdml_parser.GetAuxMap())
{
for (G4GDMLAuxStructType const& aux : vecaux)
{
std::string const& sd_name = aux.value;
if (aux.type == "SensDet")
{
result.detectors.insert({sd_name, lv});
if (auto* refl_lv = refl_factory->GetReflectedLV(lv))
{
// Add the reflected volume as well
result.detectors.insert({sd_name, refl_lv});
}
}
}
}
}

if (opts_.pointers == PointerTreatment::remove)
{
remove_pointers(*G4SolidStore::GetInstance());
remove_pointers(*G4PhysicalVolumeStore::GetInstance());
remove_pointers(*G4LogicalVolumeStore::GetInstance());
}

CELER_ENSURE(result.world);
return result;
}

//---------------------------------------------------------------------------//
/*!
* Write a GDML file to the given filename.
*/
void save_gdml(G4VPhysicalVolume const* world, std::string const& out_filename)
{
CELER_EXPECT(world);

CELER_LOG(info) << "Writing Geant4 geometry to GDML at " << out_filename;
ScopedMem record_mem("save_gdml");
ScopedTimeLog scoped_time;

ScopedGeantLogger scoped_logger;
ScopedGeantExceptionHandler scoped_exceptions;

G4GDMLParser parser;
parser.SetOverlapCheck(false);

if (!world->GetLogicalVolume()->GetRegion())
{
CELER_LOG(warning) << "Geant4 regions have not been set up: skipping "
"export of energy cuts and regions";
}
else
{
parser.SetEnergyCutsExport(true);
parser.SetRegionExport(true);
}

parser.SetSDExport(true);
parser.SetStripFlag(false);
#if G4VERSION_NUMBER >= 1070
parser.SetOutputFileOverwrite(true);
#endif

parser.Write(out_filename, world, /* append_pointers = */ true);
}

//---------------------------------------------------------------------------//
} // namespace celeritas
Loading

0 comments on commit 3fe3999

Please sign in to comment.