Creating Geometry

Precepts

Units may be specified but default to Geant4 ones (e.g. mm, rad).
Rotations are made using Tait-Bryan angles (rotation about reference x,y,z axes).
A Registry object should be used to hold all things related in a model and passed into the constructors of most objects.
GDML-like full lengths are used instead of typically half lengths

Units

The default units are mm and rad for length and angle.
Most constructors will take units as an optional key word argument (‘kwarg’).
The kwargs are typically lunit for length unit and aunit for angle unit.
Units are defined in pyg4ometry.gdml.Units.py.

The following units (as strings) are accepted:

Unit	Value
nm	1e-6
um	1e-3
mm	1
cm	10
m	1e3
km	1e6
deg	\(\pi/180\)
degree	\(\pi/180\)
rad	1
radian	1
mrad	1e-3
urad	1e-6
eV	1e-3
keV	1
MeV	1e+3
none	1
ns	1e-9
us	1e-6
ms	1e-3
s	1

Examples:

reg = pyg4ometry.geant4.Registry()
boxSolid = pyg4ometry.genat4.solid.Box("aBox", 10, 20, 30, reg)

This defines a box with the default units (none specified), so mm. We can specify them:

boxSolid = pyg4ometry.genat4.solid.Box("aBox", 10, 20, 30, reg, "cm")

Geant4 Python Scripting

Making use of pyg4ometry requires the following modules

import pyg4ometry

To make a simple geometry of a box located at the origin

# load pyg4ometry
import pyg4ometry

# registry to store gdml data
reg  = pyg4ometry.geant4.Registry()

# world solid and logical
ws   = pyg4ometry.geant4.solid.Box("ws",50,50,50,reg)
wl   = pyg4ometry.geant4.LogicalVolume(ws,"G4_Galactic","wl",reg)
reg.setWorld(wl.name)

# box placed at origin
b1   = pyg4ometry.geant4.solid.Box("b1",10,10,10,reg)
b1_l = pyg4ometry.geant4.LogicalVolume(b1,"G4_Fe","b1_l",reg)
b1_p = pyg4ometry.geant4.PhysicalVolume([0,0,0],[0,0,0],b1_l,"b1_p",wl,reg)

# visualise geometry
v = pyg4ometry.visualisation.VtkViewer()
v.addLogicalVolume(wl)
v.addAxes(20)
v.view()

Here is the vtk visualiser output of the above example

GDML Defines

In GDML there are multiple define objects that can be used parameterise geometry, materials etc. These can be used as variables or definitions and mean that any equations used will be retained in GDML output. For example a GDML constant can be created in the following way

# registry to store gdml data
reg = pyg4ometry.geant4.Registry()

# constant called x
x = pyg4ometry.gdml.Constant("x",10,reg)

The normal set of mathematical operations in python can be performed and evaluated

y = 2*x + 10
y.eval()

>> 30

The constant x can of course be changed and y re-evaluated

x.setExpression(5)
y.eval()

>> 20

Note

Standard mathematical functions can be used with GDML defines (Constant, Variable, etc). So sin, cos, tan, exp and so on, but pyg4ometry functions have to be used

x  = pyg4ometry.gdml.Constant("x",10,reg)
cx = pyg4ometry.gdml.cos(x)

So the box example above can be rewritten using constants

# load pyg4ometry
import pyg4ometry

# registry to store gdml data
reg  = pyg4ometry.geant4.Registry()

bx = pyg4ometry.gdml.Constant("bx","10",reg,True)
by = pyg4ometry.gdml.Constant("by",2*bx,reg,True)
bz = pyg4ometry.gdml.Constant("bz",2*by,reg,True)

# world solid and logical
ws   = pyg4ometry.geant4.solid.Box("ws",50,50,50,reg)
wl   = pyg4ometry.geant4.LogicalVolume(ws,"G4_Galactic","wl",reg)

# box placed at origin
b1   = pyg4ometry.geant4.solid.Box("b1",bx,by,bz,reg)
b1_l = pyg4ometry.geant4.LogicalVolume(b1,"G4_Fe","b1_l",reg)
b1_p = pyg4ometry.geant4.PhysicalVolume([0,0,0],[0,0,0],b1_l,"b1_p",wl,reg)

# visualise geometry
v = pyg4ometry.visualisation.VtkViewer()
v.addLogicalVolume(wl)
v.addAxes(20)
v.view()

Note

All GDML defines (Constant, Variable, etc) can be used in the construction of other pyg4ometry classes interchangeably instead of floats or strings (where strings are either numbers or a GDML expression)

Warning

Avoid reassigning variables used as defines, this can have unexpected consequences so for example

b1   = pyg4ometry.geant4.solid.Box("b1",bx,by,bz,reg)
b1.pX = 20              # do not do this
b1.pX.setExpression(20) # rather do this

Solids

The python geant4 solids match the Geant4 constructors as much possible (different constructor signatures are not supported in python). For example looking at the G4Box class

pyg4ometry.geant4.solid.Box(name, pX, pY, pZ, registry, lunit)

G4Box(const G4String& pName, G4double  pX, G4double  pY, G4double pZ)

A full list of solids can be found in Geant4 solids.

Warning

The parameters stick to the GDML convention of full lengths opposed to half lengths.

Materials

As with solids materials are defined in a similar way to Geant4 C++. Python does not have overloaded constructors, so unique signatures are needed, in contrast to Geant4.

To define a material from the Geant4 predefined (e.g. NIST) materials

import pyg4ometry.geant4 as _g4
wm = _g4.MaterialPredefined("G4_Galactic")
bm = _g4.MaterialPredefined("G4_Fe")

To define a single element in terms of atomic number, atomic mass and density.

import pyg4ometry.geant4 as _g4
wm = _g4.MaterialSingleElement("galactic",1,1.008,1e-25,reg)   # low density hydrogen
bm = _g4.MaterialSingleElement("iron",26,55.8452,7.874,reg)    # iron at near room temp

To define a compound two elements using the mass fraction

import pyg4ometry.geant4 as _g4
wm = _g4.MaterialCompound("air",1.290e-3,2,reg)
ne = _g4.ElementSimple("nitrogen","N",7,14.01)
oe = _g4.ElementSimple("oxygen","O",8,16.0)
wm.add_element_massfraction(ne,0.7)
wm.add_element_massfraction(oe,0.3)
bm = _g4.MaterialSingleElement("iron",26,55.8452,7.874,reg)    # iron at near room temp

To define a compound using number of atoms

import pyg4ometry.geant4 as _g4
bm = _g4.MaterialCompound("plastic",1.38,3,reg)    # Generic PET C_10 H_8 O_4
he = _g4.ElementSimple("hydrogen","H",1,1.008)
ce = _g4.ElementSimple("carbon","C",6,12.0096)
oe = _g4.ElementSimple("oxygen","O",8,16.0)
bm.add_element_natoms(he,8)
bm.add_element_natoms(ce,10)
bm.add_element_natoms(oe,4)

Material as a mixture of materials

import pyg4ometry.geant4 as _g4
bm     = _g4.MaterialCompound("YellowBrass_C26800", 8.14, 2, reg)
copper = _g4.MaterialPredefined("G4_Cu")
zinc   = _g4.MaterialPredefined("G4_Zn")
bm.add_material(copper, 0.67)
bm.add_material(zinc, 0.33)

Example of elements formed by isotopes

import pyg4ometry.geant4 as _g4
u235 = _g4.Isotope("U235", 92, 235, 235.044)
u238 = _g4.Isotope("U238", 92, 238, 238.051)
uranium = _g4.ElementIsotopeMixture("uranium", "U", 2)
uranium.add_isotope(u235, 0.00716)
uranium.add_isotope(u238, 0.99284)
bm = _g4.MaterialCompound("natural_uranium", 19.1, 1, reg)
bm.add_element_massfraction(uranium, 1)

NIST Materials

Geant4 has many predefined materials according to the NIST database. Their name typically starts with G4_. These typically can be used with MaterialPredefined and we do not need to specify the full composition - Geant4 will find them at run time.

However, in the case of conversion to FLUKA, these are fully expanded according to their definition in Geant4 based on a cache in pyg4ometry of the material compositions generated using BDSIM from Geant4 (10.7.p01 as of writing). Should the user wish to use these, they can be accessed from the functions in the geant4 module.

import pyg4ometry
nistHydrogenElement = pyg4ometry.geant4.nist_element_2geant4Element('G4_H')

Note, an ‘element’ cannot be used as a ‘material’ in a logical volume. We must upgrade it to a material for that. The NIST elements contain the appropriate mixture of natural isotopes and can be used in MaterialCompound as demonstrated above.

Alternatively, we can access the NIST materials and materials of elements.

import pyg4ometry
nistHydrogenMaterial = pyg4ometry.geant4.nist_material_2geant4Material('G4_H')
nistConcreteMaterial = pyg4ometry.geant4.nist_material_2geant4Material('G4_CONCRETE')

Detector Construction

This largely proceeds in exactly the same way as in G4 or GDML. Hierarchy of solids, booleans, logical, physical (replica, division, param) volumes.

Create registry to hold everything
Create solids
Create logical volumes
Place logical volumes (construct physical volumes)
Visualise
Check
Export

Transformations & Physical Volumes

Transformations in 3D are essential for the easy placement of solids in a CSG tree or LV placement. There is not a specific transformation class like in Geant4. The matrices and vectors used for placements are here typically Numpy arrays or matrices.

Geant4 has two possible constructors for a physical volume. These provide active and passive transformations. In pyg4ometry, only one is provided.

The transform in a physical volume first translates the placed logical volume with respect to the mother logical, then rotates it.

The physical volume class is documented here: Geant4 module, but an example is shown here.

import pyg4ometry
r = pyg4ometry.geant4.Registry()
vacuum = _g4.MaterialPredefined("G4_Galactic")
water = _g4.MaterialPredefined("G4_WATER")
worldSolid = pyg4ometry.geant4.solid.Box("world_solid", 100, 100, 100, reg)
boxSolid = pyg4ometry.geant4.solid.Box("box_solid", 10, 20, 40, reg)
worldLV = pyg4ometry.geant4.LogicalVolume(worldSolid, vacuum, "world_lv", reg)
boxLV = pyg4ometry.geant4.LogicalVolume(boxSolid, water, "box_lv", reg)

pyg4ometry.geant4.PhysicalVolume([0,0,0],
                                 [0,0,0],
                                 boxLV,
                                 "box_pv",
                                 worldLV,
                                 reg)

This creates a box of water inside a box of vacuum. The box of water is 10 x 20 x 50 mm long (note mm are the default length units), and it is placed with no offset and no rotation (i.e. at the centre) of the world volume. Alternatively:

import numpy as np
pyg4ometry.geant4.PhysicalVolume([0,np.pi/3.0,0],
                                 [0,0,0],
                                 boxLV,
                                 "box_pv",
                                 worldLV,
                                 reg)

In this case, the box is placed with no offset but with a rotation of \(\pi/3\) radians about the y axis of the world box.

Note

The rotations are Tait-Bryan angles, which are rotations about the reference x,y,z axes. i.e. if there is a rotation about both x and y, these are independent and it is not a compound frame that is rotated. These are commonly thought of like an aircraft and called pitch, yaw and tilt.

There are utility functions for translation between different transformations in pyg4ometry.transformation. See Transformation.

Optical Surfaces

Optical surfaces can be created in a similar way as in Geant4 C++. A pyg4ometry.geant4.solid.OpticalSurface instance holds all the needed properties of the surface (including extra properties, e.g. for optical processes). This is then assigned to the surface between either

two physical volumes: pyg4ometry.geant4.BorderSurface, or
a logical volume and all its neighbouring volumes: pyg4ometry.geant4.SkinSurface.

opa = _g4.solid.OpticalSurface("AirSurface", finish="polished", model="glisur", surf_type="dielectric_dielectric", value="1", registry=reg)
opw = _g4.solid.OpticalSurface("WaterSurface", finish="ground", model="unified", surf_type="dielectric_dielectric", value="0", registry=reg)

_g4.SkinSurface("AirSurface", air_lv, opa, reg)
_g4.BorderSurface("WaterSurface", water_phys, world_phys, opw, reg)

Properties of Materials and Optical Surfaces

Materials and optical surfaces support adding properties that can be used by Geant4 to influence processes, e.g. for scintillation, refraction or other optical processes.

In the GDML, a matrix is used to hold the value(s) of the property.

addProperty(name, matrix) - Add a property based on an existing pyg4ometry.gdml.Matrix object.
addVecProperty(name, e, v, eunit='eV', vunit='') - Add a property based on a energy vector and a value vector.
addConstProperty(name, value, vunit='')- Add a property that has only one constant value.

Units can be specified by setting the parameters eunit for the energy vector and vunit for the values. The given vectors are expected to be homogeneous in their units.

Note

Optical properties can only use units (or combinations of units) that are also defined in pyg4ometry. If needed, additional units can be added: pyg4ometry.gdml.Units.units['ps'] = 1e-12.

scint = _g4.Material(...)
scint.addConstProperty('SCINTILLATIONTIMECONSTANT1', 2.5, vunit='ns')
scint.addConstProperty('SCINTILLATIONYIELD', 8000, vunit='/MeV')
scint.addVecProperty('RINDEX', [1, 10], [1.3, 1.05])

FLUKA Geometry Creation

In a very similar way to geant4 geometry authoring it is possible to use pyg4ometry to create fluka output. To create a simple region consisting of a single body

import pyg4ometry.convert as convert
import pyg4ometry.visualisation as vi
from pyg4ometry.fluka import RPP, Region, Zone, FlukaRegistry

freg = FlukaRegistry()

rpp = RPP("RPP_BODY", 0, 10, 0, 10, 0, 10, flukaregistry=freg)
z = Zone()
z.addIntersection(rpp)
region = Region("RPP_REG", material="COPPER")
region.addZone(z)
freg.addRegion(region)

greg = convert.fluka2Geant4(freg)
greg.getWorldVolume().clipSolid()

v = vi.VtkViewer()
v.addAxes(length=20)
v.addLogicalVolume(greg.getWorldVolume())
v.view()