Low-level Celeritas integration

This subsection contains details of importing Geant4 data into Celeritas.

Geant4 geometry utilities

class GeantGdmlLoader

Load a GDML file into memory.

The pointer treatment gives three options:

  • ignore leaves names as they are imported by Geant4’s GDML reader, which strips them from material/region names but leaves solid/logical/physical pointers in place.

  • truncate lets the Geant4 GDML remove the pointers, which cuts everything after 0x including suffixes like _refl added during volume construction.

  • remove uses a regular expression to remove pointers from volume names.

The detectors option reads auxiliary tags in the structure that have auxtype=SensDet and returns a multimap of strings to volume pointers.

inline G4VPhysicalVolume *celeritas::load_gdml(std::string const &filename)

Load a Geant4 geometry, excising pointers.

This provides a good default for using GDML in Celeritas.

Returns:

Geant4-owned world volume

inline void celeritas::save_gdml(G4VPhysicalVolume const *world, std::string const &out_filename)

Write a GDML file to the given filename.

inline std::unordered_set<G4LogicalVolume const*> celeritas::find_geant_volumes(std::unordered_set<std::string> names)

Find Geant4 logical volumes corresponding to a list of names.

If logical volumes with duplicate names are present, they will all show up in the output and a warning will be emitted. If one is missing, a RuntimeError will be raised.

static std::string_view const labels[] = {"Vol1", "Vol2"};
auto vols = find_geant_volumes(make_span(labels));

Geant4 physics interfaces

class GeantImporter : public celeritas::ImporterInterface

Load problem data directly from Geant4.

This can be used to circumvent ROOT as a serialization tool, whether to simplify the toolchain or to integrate better with user frameworks. As much data as possible is imported (subject to the data selection); downstream Celeritas classes will validate the imported data as needed.

GeantImporter import(GeantSetup("blah.gdml"));
ImportData data = import();
or to import from an existing, initialized Geant4 state:
GeantImport import(world_volume);
ImportData data = import();

class GeantSetup

Construct a Geant 4 run manager and populate internal Geant4 physics.

This is usually passed directly into GeantImporter . It hides Geant4 implementation details (including header files) from the rest of the code. It is safe to include even when Geant4 is unavailable!

The setup is targeted specifically for physics that Celeritas supports.

Geant4 physics options

struct GeantPhysicsOptions

Construction options for Geant physics.

These options attempt to default to our closest match to G4StandardEmPhysics. They are passed to the EmPhysicsList and FtfpBert physics lists to provide an easy way to set up physics options.

Gamma physics

bool compton_scattering = {true}

Enable Compton scattering.

bool photoelectric = {true}

Enable the photoelectric effect.

bool rayleigh_scattering = {true}

Enable Rayleigh scattering.

bool gamma_conversion = {true}

Enable electron pair production.

bool gamma_general = {false}

Use G4GammaGeneral instead of individual gamma processes.

Electron and positron physics

bool coulomb_scattering = {false}

Enable discrete Coulomb.

bool ionization = {true}

Enable e- and e+ ionization.

bool annihilation = {true}

Enable positron annihilation.

BremsModelSelection brems = {BremsModelSelection::all}

Enable bremsstrahlung and select a model.

MevEnergy seltzer_berger_limit = {1e3}

Upper limit for the Seltzer-Berger bremsstrahlung model.

MscModelSelection msc = {MscModelSelection::urban}

Enable multiple coulomb scattering and select a model.

RelaxationSelection relaxation = {RelaxationSelection::none}

Enable atomic relaxation and select a model.

Physics options

int em_bins_per_decade = {7}

Number of log-spaced bins per factor of 10 in energy.

bool eloss_fluctuation = {true}

Enable universal energy fluctuations.

bool lpm = {true}

Apply relativistic corrections for select models.

bool integral_approach = {true}

See PhysicsParamsOptions::disable_integral_xs.

Cutoff options

MevEnergy min_energy = {0.1 * 1e-3}

Lowest energy of any EM physics process.

MevEnergy max_energy = {100 * 1e6}

Highest energy of any EM physics process.

double linear_loss_limit = {0.01}

See PhysicsParamsOptions::linear_loss_limit.

MevEnergy lowest_electron_energy = {0.001}

Tracking cutoff kinetic energy for e-/e+.

MevEnergy lowest_muhad_energy = {0.001}

Tracking cutoff kinetic energy for muons/hadrons.

bool apply_cuts = {false}

Kill secondaries below the production cut.

double default_cutoff = {0.1 * units::centimeter}

Set the default production cut for all particle types [len].

Multiple scattering configuration

double msc_range_factor = {0.04}

e-/e+ range factor for MSC models

double msc_muhad_range_factor = {0.2}

Muon/hadron range factor for MSC models.

double msc_safety_factor = {0.6}

Safety factor for MSC models.

double msc_lambda_limit = {0.1 * units::centimeter}

Lambda limit for MSC models [len].

double msc_theta_limit = {constants::pi}

Polar angle limit between single and multiple Coulomb scattering.

double angle_limit_factor = {1}

Factor for dynamic computation of angular limit between SS and MSC.

bool msc_displaced = {true}

Whether lateral displacement is enabled for e-/e+ MSC.

bool msc_muhad_displaced = {false}

Whether lateral displacement is enabled for muon/hadron MSC.

MscStepLimitAlgorithm msc_step_algorithm = {MscStepLimitAlgorithm::safety}

Step limit algorithm for e-/e+ MSC models.

MscStepLimitAlgorithm msc_muhad_step_algorithm{MscStepLimitAlgorithm::minimal}

Step limit algorithm for muon/hadron MSC models.

NuclearFormFactorType form_factor = {NuclearFormFactorType::exponential}

Nuclear form factor model for Coulomm scattering.

Public Members

GeantMuonPhysicsOptions muon

Muon EM physics.

bool verbose = {false}

Print detailed Geant4 output.

GeantOpticalPhysicsOptions optical = {GeantOpticalPhysicsOptions::deactivated()}

Optical physics options.