Table of Contents
- Setup
- Visualize and Check Geometry
- Running Geant4 simulation
- Using detector description in analysis
- References
Setup
Note: this tutorial component assume you are in an eic-shell
singularity session, and that you
load the pre-compiled detector plugins from the container with
source /opt/detector/setup.sh
Start from the top level directory we created the Quick Start, which is ~/eic
in our example.
cd ~/eic
First start by cloning part 1 of the tutorial.
git clone https://eicweb.phy.anl.gov/EIC/tutorials/ip6_tutorial_1.git part1cd part1
Inspect files
Look around see what's what.
The files
compact/materials.xml
,
compact/elements.xml
,
and
gem_tracker.xml
are the compact detector description files. The first two files are included in
the main gem_tracker.xml
file. Note these files are often found in
directories named compact (but not always).
Compile the detector library
cmake -B build -S . -DCMAKE_INSTALL_PREFIX=$ATHENA_PREFIX -DCMAKE_CXX_STANDARD=17cmake --build build -j4 -- install
Visualize and Check Geometry
Look at disk tracker
dd_web_display --export gem_tracker.xml
Copy the generated file detector_geometry.root
to your local machine.
Then open the web viewer, click on the ...
and open the local root file.
Note this detector is not a good example to copy. It is just for instructional purposes. Examine the compact description, specifically this part:
<detector id="2" name="GEMTracker" vis="RedVis" type="my_GEMTracker" readout="GEMTrackerHits" ><layer id="1" z="-100 *cm" inner_r="40*cm" outer_r="120*cm" phi0_offset="0.0*deg" /><layer id="2" z="-80 *cm" inner_r="30*cm" outer_r="90*cm" phi0_offset="0.0*deg" /><layer id="3" z="-60 *cm" inner_r="20*cm" outer_r="70*cm" phi0_offset="0.0*deg" /><layer id="4" z="-40 *cm" inner_r="10*cm" outer_r="20.0*cm" phi0_offset="0.0*deg" /><layer id="5" z=" 40 *cm" inner_r="10*cm" outer_r="20.0*cm" phi0_offset="0.0*deg" /><layer id="6" z=" 60 *cm" inner_r="25*cm" outer_r="70.0*cm" phi0_offset="0.0*deg" /><layer id="7" z=" 80 *cm" inner_r="30*cm" outer_r="90.0*cm" phi0_offset="0.0*deg" /><layer id="8" z="100 *cm" inner_r="40*cm" outer_r="100.0*cm" phi0_offset="0.0*deg" /></detector>
Check for overlaps
checkOverlaps -t 0.0001 -c gem_tracker.xml
...Info in <TGeoNodeMatrix::CheckOverlaps>: Checking overlaps for world_volume and daughters within 0.0001Check overlaps: [==========] 11 [100.00 %] 00:00Info in <TGeoNodeMatrix::CheckOverlaps>: Number of illegal overlaps/extrusions : 0
Look at all the constants
Look at the constants as defined in the compact description (xml files).
npdet_info dump gem_tracker.xml
This should output something like this:
CrossingAngle = 0.020 = 0.020*radForwardTrackerPlane_z0 = 400.000 = 400*cmPlace_Center = 0.000 = 0*cmcompact_checksum = 1924675351.000 = 1924675351tracker_region_rmax = 200.000 = 2.0*mtracker_region_zmax = 400.000 = 4.0*mworld_side = 1000.000 = 10*mworld_x = 1000.000 = world_sideworld_y = 1000.000 = world_sideworld_z = 10000.000 = 10*world_side
Or if there are too many constants to dump, try searching:
npdet_info search world --all gem_tracker.xml
world_side = 1000.000 = 10*mworld_x = 1000.000 = world_sideworld_y = 1000.000 = world_sideworld_z = 10000.000 = 10*world_side
This is very useful with a good subsystem naming convention.
Running Geant4 simulation
npsim
is a clone of DD4hep's ddsim
but has a different output plugin.
npsim --runType run --enableG4GPS \--macroFile gps.mac \--compactFile ./gem_tracker.xml \--outputFile gem_tracker_sim.root
This will run geant4 using the General Particle Source (GPS) tool. The output will be a root file containing the generated events and GEM tracker hits. Look at gps.mac.
Refer to the GPS Documentation for more setup possibilities.
Using detector description in analysis
Simulation output data model
The geant4 output data model used by npsim
is described in a single yaml file.
Note that this is purposefully factorized from the larger EIC data model, eicd
, which is used
at every step post geant4.
Cell ID to position
root -b -q scripts/tutorial1_hit_position.cxx+
Look at the generated plot in results
:
This is the local position in the segmentation cell (i.e. pixel, strip, readout pad, etc.).
Detector handle
At the top of the script you can see these lines:
dd4hep::Detector& detector = dd4hep::Detector::getInstance();detector.fromCompact("gem_tracker.xml");dd4hep::rec::CellIDPositionConverter cellid_converter(detector);
This accomplishes the following:
- gets the main DD4hep Detector instance
- loads the compact detector file
- initializes the position converter tool (which provides thread safe access) Later on you can see its use while looping over hits:
auto pos0 = h.position;auto pos1 = cellid_converter.position(h.cellID);
Here pos0
is a dd4pod tracker hit
which has position
, a
dd4pod four vector.
This is the position (and time) of the Geant4 step through the sensitive tracking volume1.
pos1
is a three vector of the smallest sensitive detector segment.
The converter tool is used along with the hit's cellID
(i.e. unique channel number) to look up this position2.
- 1 Technically it is the average between the pre and post step points.
- 2 This is the first example showing how to use a single source of geometry information (i.e. the compact description file).
Accessing segmentation cell size
root -b -q scripts/tutorial2_cell_size.cxx+
The output should contain a lot of statements like this:
...Segmentation-Cell Position : 39,102,-100dim 1, 3,...
The following snippet shows how to get cell size information from within the hit loop.
auto context = cellid_converter.findContext( h.cellID ) ;dd4hep::Readout r = cellid_converter.findReadout( context->element ) ;dd4hep::Segmentation seg = r.segmentation() ;auto cell_dim = seg.cellDimensions(h.cellID);
Here cell_dim
is a std::vector<double>
passed by the specific segmentation implementation.
Where is the segmentation defined?
Look at the readouts
tag of gem_tracker.xml
and you will see this:
<readout name="GEMTrackerHits"><segmentation type="CartesianGridXY" grid_size_x="1*cm" grid_size_y="3*cm" /><id>system:8,barrel:2,layer:4,module:12,sensor:2,x:32:-16,y:-16</id><!--<segmentation type="PolarGridRPhi" grid_size_phi="3.0*degree" grid_size_r="5.0*cm"/><id>system:5,barrel:3,layer:4,module:5,r:32:-16,phi:-16</id>--></readout>
Note the grid_size_x
and grid_size_y
attributes with values corresponding to the sizes above.
In the CartesianGridXY
case, the cell sizes are always the same.
Change the segmentation
- Swap in the comment on segmentation and id to create the
PolarGridRPhi
segmentation. Should look like this:
<readout name="GEMTrackerHits"><segmentation type="PolarGridRPhi" grid_size_phi="3.0*degree" grid_size_r="5.0*cm"/><id>system:5,barrel:3,layer:4,module:5,r:32:-16,phi:-16</id></readout>
- You can now re-run the npsim command above, as we changed the detector pixel layout.
npsim --runType run --enableG4GPS \--macroFile gps.mac \--compactFile ./gem_tracker.xml \--outputFile gem_tracker_sim.root
- Re-run the root script:
root -b -q scripts/tutorial2_cell_size.cxx+
The output should look similar
...Segmentation-Cell Position : -96.44711531372378,-62.6334890267281,-100dim 5, 6.02139,Hit Position : 34.53880143830617,-11.266496282927323,-60Segmentation-Cell Position : 33.28697807033037,-10.815594803123158,-60dim 5, 1.8326,Hit Position : 46.356572126864805,-13.944556949561113,-80Segmentation-Cell Position : 47.552825814757675,-15.450849718747369,-80dim 5, 2.61799,...
Note that the dimensions of the cell in one direction is changing the rho coordinate! This is expected because the cell has dimensions rho =x and r*dphi = y.
Here we point out that access to this information is very useful. For example, a tracking covariance matrix could be computed using this method, again coming (dynamically) from a single source geometry description.
Id Specification
The cell ID is a 64bit integer with field names assigned to a few bits. Looking at the id specification we see the following:
<segmentation type="PolarGridRPhi" grid_size_phi="3.0*degree" grid_size_r="5.0*cm"/><id>system:5,barrel:3,layer:4,module:5,r:32:-16,phi:-16</id>
This means there are 5 bits associated with system (this is the <detector>
's id
). The field is mandatory but can differ from 5 bits.
The subsequent fields are arbitrary but should result in a unique 64 bit integer.
Here 3 bits for the "barrel", 4 bits for a "layer", 5 bits for the "module".
More bit fields could be allocated if desired (in part 2 it will become clear why we might want to do this).
The last two bit fields start at the 32nd bit, leaving 15 unused bits. These "r" and "phi" fields are not associated with physical volumes but rather the 2D segmentation above. The segmentation is a virtual geometry applied the smallest geometry volumes (nodes). So the last 32 bits are split between these two fields and are signed.
Note the concept of a VolumeID
is essentially equal to CellID
but where the bits in CellID
associated with a segmentation are set to zero. VolumeID
s are associated with actual constructed and placed volumes.
Run the example:
root -b -q scripts/tutorial3_id_spec.cxx+
You should see this (with the R-phi segmentation)
ID specification:system:0:5,barrel:5:3,layer:8:4,module:12:5,r:32:-16,phi:48:-16"layer" index is 2.
or this with the XY grid segmentation:
ID specification:system:0:8,barrel:8:2,layer:10:4,module:14:12,sensor:26:2,x:32:-16,y:48:-16"layer" index is 2.
This comes from the snippet at the initialization stage:
fmt::print("ID specification:\n");auto decoder = detector.readout("GEMTrackerHits").idSpec().decoder();fmt::print("{}\n", decoder->fieldDescription());auto layer_index = decoder->index("layer");fmt::print(" \"layer\" index is {}.\n", layer_index);
This gets an index for the named field "layer" and initializes a decoder. See BitFieldCoder documentation documentation for more details.
Later in the loop over event hits.
auto detector_layer = decoder->get(h.cellID, layer_index);if ((detector_layer !=4 ) && (detector_layer !=5 )) {continue;}
Note the layer_index
is computed once in the beginning when the decoder is also initialized.
This means it will be fast because it only parses the ID specification string once.
The statement above selects for the interior layers number 4 and 5:
Look at the generated plot in results
:
Homework
Can you adjust the compact description so that layers 4 and 5 have similar number of hits as the other layers?
References
- NPdet
- BitFieldCoder documentation documentation
- List of DD4hep Segmentations
- List of DD4hep detector constructors
- Compact description examples, and more examples