View on GitHub

GEARS | 齿轮组

Geant4 Example Application with Rich features and Small footprints

Download this project as a .zip file Download this project as a tar.gz file

YouTube Python Step point Detected energy Data analysis Detector hits

Output

Generally speaking, the visualization of detector geometry and the screen dump of a Geant4 application can be all regarded as output of a Geant4 simulation. Strictly speaking, the output of a Geant4 simulation includes histograms and/or ntuples of data generated during the simulation, which can be used to reveal statistical distributions of, for example, positions and energy depositions of interactions.

GEARS utilizes Geant4 analysis managers to provide four output formats: ROOT, HDF5, CSV, and AIDA XML.

The output file name can be chosen using the macro command:

/analysis/setFileName gears.root

One of the following suffix is needed to specify the output file format: .root, .hdf5, .csv, .xml. Note that the output is disabled by default. It will be enabled if the output file name is not empty. So this macro command also works as a switch. Without it, no output file will be created.

Screen dump

Geant4 can print out on screen detailed information of a simulation if you increase the verbose level of tracking using the macro command /tracking/verbose, for example,

 # turn on detailed information about tracking
 /tracking/verbose 2
 # run a few events for debugging
 /run/beamOn 10
 # turn off screen dump for fast simulation
 /tracking/verbose 0
 # run a lot of events
 /run/beamOn 100000

You can re-direct the screen dump to a file for detailed examination:

$ gears macro.mac | tee log
$ less log

ROOT

The ROOT version of ntuples is TTree, which is commonly called a tree. Simply put, a tree is a table. Each row is called an entry (or an event) and each column is called a branch (of course, we are dealing with a tree after all). If you simulate 1000 events using GEARS, your ROOT tree would have 1000 entries. If you save the first interaction position x, y, and z in the tree, your tree would have 3 branches, each holds 1000 values of x, y, or z. If an event has more than just one hit, you will have a few x, y, and z’s for each event (row), and your table would have one more dimension, labeled as Iteration$ in TTree, the depth of which may change event by event since each event may have a different number of hits.

The ROOT TTree offers the following features that are desired for analyzing data generated by GEARS:

Here are some example codes that can be run in a ROOT interactive session to generate histograms:

$ root gears.root # open the output file in root format using ROOT
root [] .ls
TFile**         gears.root
 TFile*         gears.root
   KEY: TTree    t;1     Geant4 step points
root [] t->GetEntries()
(long long) 5000
root [] t->Show(0)
======> EVENT:0
 n   = 235
 m   = 2
 trk = (vector<int>*)0x351a520
 stp = (vector<int>*)0x398c7b0
 vlm = (vector<int>*)0x2d4dba0
 pro = (vector<int>*)0x3524f60
 pdg = (vector<int>*)0x34302c0
 pid = (vector<int>*)0x3972e80
 xx  = (vector<double>*)0x3548ce0
 yy  = (vector<double>*)0x2d4ef10
 zz  = (vector<double>*)0x354e600
 dt  = (vector<double>*)0x3549db0
 de  = (vector<double>*)0x34b2d90
 dl  = (vector<double>*)0x2f63630
 l   = (vector<double>*)0x2f636f0
 x   = (vector<double>*)0x2d50760
 y   = (vector<double>*)0x351a450
 z   = (vector<double>*)0x3543f10
 t   = (vector<double>*)0x351eed0
 k   = (vector<double>*)0x3a907c0
 p   = (vector<double>*)0x2f307c0
 q   = (vector<double>*)0x2ae30c0
 et  = (vector<double>*)0x35293f0
root [] t->Draw("x","k*(pdg==22)")

Python

ROOT is less known than Python outside of the high energy physics community. The good news for people who are not familiar with ROOT is that GEARS does not depend on ROOT to compile or run, even though its output can be saved in ROOT TTree format, and that the analysis of GEARS output can be done in Python instead of ROOT thanks to the uproot Python package. The best way to get started with analyzing GEARS output in Python would be to follow the uproot tutorial.

If you are familiar with ROOT and would like to migrate to Python for analyzing GEARS output, you can find here a brief list of Python equivalence of ROOT commands:

Julia

Utilizing https://github.com/JuliaHEP/UpROOT.jl, one can load simulation results in ROOT format in Julia as well:

julia> import Pkg
julia> Pkg.add("UpROOT")
julia> using UpROOT
julia> file = TFile("gears.root")
julia> tree = file["t"]
julia> tree.x

Step point

Geant4 tracks a particle step by step as it passes through the simulated world until it goes out of it, gets absorbed in a material, or changes to other particles, as shown in the following figure. A step point is a point in a particle track where the particle is generated or changed.

tracks.png

Geant4 artificially adds a step point on the boundary of two volumes when the track goes from one volume into the other, even if there is no change of the particle.

A step point in GEARS contains the following information:

They are saved in separated C++ vectors (arrays with various sizes). Such a flat data structure and very short variable names are chosen on purpose to make plotting of those variables in a ROOT interactive session easy.

Notice that the variable n is the total number of step points recorded in each event. m is the maximal copy number of a sensitive volume

Process id

The physics process generating each step point is saved in a variable pro[i], where i is the index of the step point. It equals to (process type) * 1000 + (sub type). The Process types are defined in G4ProcessType.hh, sub types are defined in G4HadronicProcessType.hh, G4DecayProcessType.hh, G4EmProcessSubType.hh, G4TransportationProcessType.hh, G4FastSimulationProcessType.hh, G4OpProcessSubType.hh, etc. They can be found in http://www-geant4.kek.jp/lxr/find?string=Type.hh.

Particle id

The type of particle related to a step point is saved in a variable pdg. It is the same as the PDG encoding of the particle. A Google search will give more information about it. The name pid is used for the parent particle’s track id. A complete list of nuclei PDG encoding can be found here.

Record information of step 0

One cannot get access to step 0 (initStep printed on screen if /tracking/verbose is set to be more than 0) through G4UserSteppingAction, which only provides access to step 1 and afterwards. However, step 0 contains some very important information that is constantly requested by the user. For example, the energy of a gamma ray from a radioactive decay is only available from step 0. Such information can be easily displayed using the following ROOT command with the Output ROOT tree, t:

  // draw kinetic energy, "k", of a gamma (pdg==22)
  // created by radioactive decay process (pro==6210)
  t->Draw("k","pro==6210 && pdg==22")

This is achieved by using G4SteppingVerbose instead of G4UserSteppingAction for data recording. The former provides a function called TrackingStarted() to print verbose information about step 0 on screen if /tracking/verbose is set to be more than 0. It also provides a function called StepInfo() to print verbose information about steps after step 0 on screen if /tracking/verbose is more than 0. G4SteppingVerbose::StepInfo() is called right before G4UserSteppingAction::UserSteppingAction(G4Step) is called in G4SteppingManager::Stepping(), it hence can be used to fully replace the functionality of G4UserSteppingAction::UserSteppingAction(G4Step).

In fact, G4UserSteppingAction::UserSteppingAction(G4Step*) is not used at all in GEARS. The Output class inherits TrackingStarted() and StepInfo() from G4SteppingVerbose to record data from all steps. There is another advantage of using G4SteppingVerbose instead of G4UserSteppingAction for recording, that is, G4SteppingVerbose is provided as a globally available singleton, which can be easily accessed at different places in the codes using:

  G4VSteppingVerbose::GetInstance()

This is used in G4UserRunAction to open and close a TFile, in G4UserEventAction to fill a TTree.

The catch is that functions in G4SteppingVerbose will not be called in G4SteppingManager unless /tracking/verbose is set, which will print too much information on screen for a long run. This is solved in EventAction::BeginOfEventAction by turning on tracking verbose all the time so that all functions in G4SteppingVerbose will be called, while at the same time, turning on G4SteppingVerbose local verbose flag Silent to run them in silent mode.

Total energy

In addition to arrays of parameters of individual step points, a GEARS output file also contains an array of the total energies deposited in the sensitive volumes of the simulated detector. The index of the array is the same as the copy numbers of those volumes. The name of the array is called et. Since the copy number of a sensitive volume has to start from 1, et[0] is used to store the total energy deposited in all sensitive volumes. You can use the following command to draw the total energy deposited in a sensitive volume with a copy number 1:

root [] t->Draw("et[1]")

Data analysis

One can use the following command to generate gears.root in GEARS/tutorials/output/:

$ gears radiate.mac

radiate.mac demonstrates how to use Geant4 macro commands to save step points and total energies in sensitive volumes. It uses the detector geometry defined in detector.tg.

Here are some sample ROOT commands that one can use to generate plots from gears.root:

// draw tracks of the primary particle on x-y plane in event 1
root[] t->Draw("x:y", "trk==1","l", 1, 1)
// show physics processes creating secondary particles
root[] t->Draw("pro","trk>1 && stp==0")
// display hits distribution in volume with copy number 1
root[] t->Draw("x:y:z", "vlm==1")
// show how many types of particles are involved in the simulation
root[] t->Draw("pdg");
// show physics process related to gamma-rays
root[] t->Draw("pro", "pdg==22 && stp!=0");
// show spatial distributions of secondary particles
root[] t->Draw("x:y", "trk>1")
// plot dE/dx versus momentum
root[] t->Draw("de/dl:p")
// draw energy spectrum recorded by volume (detector) 1
root[] t->Draw("et[1]")

Combine step points to hits in detector

Many of the step points in a Geant4 simulation are very close to each other, especially those of charged particles, as they quickly lose energy through multiple scattering or ionizing surrounding atoms in a very small range. The space resolution of a real-life detector is normally not enough to resolve these details. As seen by such a detector, all nearby step points act as a single hit, the energy of which is a sum of deposited energies from all these step points, the position of which is an energy-weighted average of all these step point positions. The modeling of detector response normally start with combined hits instead of the original step points directly from Geant4 to save computational power.

A ROOT script combineStepPointsToHits.C, is provided to demonstrate the combining procedure. It takes the output from GEARS and save the combined hits to another ROOT file hits.root. Another ROOT script, drawHits.C, is used to demonstrate the final result shown in the following plot.

combined hits

These ROOT scripts can be directly executed without compilation:

$ root -q combineStepPointsToHits.C
$ root drawHits.C

Note that they have to be run after you execute gears radiate.mac, since they need the Geant4 output as input.