Hands-on 5 Storing Hits

Introduction

Installation

Calorimeter Geometry Implementation

Sensitive Detector and Hit Definition

User Defined Event Action Class

Hit Visualisation

Important:
Every time you open a new shell you need to set your environment up correctly.

Introduction

By the end of this Hands-on you should be able to:

Define a simple Cesium Iodide calorimeter

Define a sensitive detector and create and update hit collections

Implement a user defined event action class

Visualise hits

This exercise builds on the HandsOn3 exercise. It is based on an experimental setup developed for bremsstrahlung benchmarking at UCSF. The simplified experiment is made up of the following components:

Vacuum beam pipe

Titanium beam window

Iron evacuation chamber window

Silicon monitor

Beryllium target

Electrons are generated with a default energy of 15 MeV in the beam pipe. The electrons interact with the target and generate bremsstrahlung photons. A square Cesium Iodide calorimeter, segmented into 100 5cm*5cm*25cm cells is also implemented. The diagram below demonstrates the simulated experiment.

The following files are provided with the exercise:

beamTest.cc - main program

BeamTestDetectorConstruction – material, geometry and sensitive detector definitions

BeamTestPrimaryGeneratorAction - primary particle generator

BeamTestPhysicsList - user defined physics list

BeamTestEventAction - G4UserEventAction class

BeamTestCellParameterisation – calorimeter cell parameterisation

BeamTestEmCalorimeter – sensitive detector for the calorimeter

BeamTestEmCalorimeterHit – calorimeter hit

All the geometries except for the calorimeter are implemented. The calorimeter geometry, sensitive detector and hits are implemented over the course of the exercise.

Installation

If using windows, make sure G4UI_USE_TCSH is not set.

To start this exercise, unpack HandsOn5.tgz to your working directory.

Compile and link it using following commands:

$ cd HandsOn5

$ make

Once you have successfully compiled and linked the code, try it to see how it runs:

$ ./bin/$G4SYSTEM/beamTest run.mac

This will produce files named G4Data0.heprep and G4Data1.heprep, which can be viewed by invoking the Wired visualization tool.

$ wired G4Data0.heprep

You should be able to view the following:

You can see that all the components with the exception of the calorimeter have been implemented.

Calorimeter Geometry Implementation

The calorimeter mother volume is defined as an air filled box with dimensions 0.5m*0.5m* 0.25m. It is placed along the z axis at a distance of 0.5m from the centre of the world volume. One hundred cells are placed within the mother volume. The cells are constructed of CsI, with dimensions 5cm*5cm*25cm. G4PVPlacement is used to construct the mother volume, while G4PVParameterised, with the parameterisation provided by BeamTestCellParameterisation, is used to construct the CsI crystals.

The implementation, along with visualisation attributes, is shown below.

BeamTestDetectorConstruction.cc
----- snipped -----
// HandsOn5: Create CsI calorimeter
#include "BeamTestCellParameterisation.hh"
#include "BeamTestEmCalorimeter.hh"
#include "G4Box.hh"
#include "G4Colour.hh"
----- snipped -----
// Define elements
G4Element* N = new G4Element("Nitrogen", symbol="N", z=7., a=14.01*g/mole);
G4Element* O = new G4Element("Oxygen", symbol="O", z=8., a=16.00*g/mole);
//HandsOn5: Define CsI
G4Element* Cs = new G4Element("Cesium", symbol="Cs", z=55., a=132.9*g/mole);
G4Element* I = new G4Element("Iodine", symbol="I", z=53., a=126.9*g/mole);
// CsI
G4Material* CsI = new G4Material("CsI", density=4.51*g/cm3, ncomponents=2);
CsI->AddElement(I, .5);
CsI->AddElement(Cs, .5);
// Define air
G4Material* air = new G4Material("Air", density= 1.290*mg/cm3, ncomponents=2);
air->AddElement(N, fractionmass=0.7);
air->AddElement(O, fractionmass=0.3);
----- snipped -----
fpWorldPhysical = new G4PVPlacement(0, // Rotation matrix pointer
G4ThreeVector(), // Translation vector
fpWorldLogical, // Logical volume
"World_Physical", // Name
0, // Mother volume
false, // Unused boolean parameter
0); // Copy number
////////////////////////////////////////////////////////////////////////
// HandsOn5: Create CsI calorimeter
// Mother volume
G4VSolid* calorimeterSolid = new G4Box("Calorimeter_Solid", // Name
.5*m, // x half length
.5*m, // y half length
.25*m) ; // z half length
G4LogicalVolume* calorimeterLogical =
new G4LogicalVolume(calorimeterSolid, // Solid
air, // Material
"Calorimeter_Logical"); // Name
new G4PVPlacement(0, // Rotation matrix pointer
G4ThreeVector(0.,0.,0.5*m), // Translation vector
calorimeterLogical, // Logical volume
"Calorimeter_Physical", // Name
fpWorldLogical, // Mother volume
false, // Unused boolean
0); // Copy number
// 100 rectangular CsI cells
G4Material* csi = G4Material::GetMaterial("CsI");
G4VSolid* cellSolid = new G4Box("Cell_Solid", // Name
5*cm, // x half length
5*cm, // y half length
25.*cm); // z half length
G4LogicalVolume* cellLogical
= new G4LogicalVolume(cellSolid, // Solid
csi, // Material
"Cell_Logical"); // Name
G4VPVParameterisation* cellParam = new BeamTestCellParameterisation();
new G4PVParameterised("Cell_Physical", // Name
cellLogical, // Logical volume
calorimeterLogical, // Mother volume
kXAxis, // Axis
100, // Number of replicas
cellParam); // Parameterisation
----- snipped -----
// Invisible world volume.
fpWorldLogical->SetVisAttributes(G4VisAttributes::Invisible);
// HandsOn5: Calorimeter attributes
// Invisible calorimeter mother volume
calorimeterLogical->SetVisAttributes(G4VisAttributes::Invisible);
// Calorimeter cells - green with transparency
G4VisAttributes* calorimeterAttributes =
new G4VisAttributes(G4Colour(0.0, 1.0, 0.0, 0.1));
calorimeterAttributes->SetForceSolid(true);
cellLogical->SetVisAttributes(calorimeterAttributes);

Implement the above and compile, link and run the program. Use Wired to check that all the geometries are now in place. The geometry should be visible in G4Data0.heprep, while G4Data1.heprep should also show trajectories.

Sensitive Detector and Hit Definition

We are interested in finding the energy deposited in each CsI cell. To do this we need to:

Create a sensitive detector (BeamTestEmCalorimeter)

Generate an associated hit collection to record the energy deposited in each cell

Register the detector with the sensitive detector manager (G4SDManager)

Attach the sensitive detector to the parameterised calorimeter cell volume

A hit class, BeamTestEmCalorimeterHit, inheriting from G4VHit, is provided with the example. A partially implemented sensitive detector class, BeamTestEmCalorimeter is also provided. We will start off by completing the BeamTestEmCalorimeter implementation.

On initialisation, BeamTestEmCalorimeter is responsible for creating a new collection of 100 BeamTestEmCalorimeterHit objects. Each BeamTestEmCalorimeterHit object accumulates the energy deposited in one CsI cell. The collection of hits is then registered with the G4HCofThisEvent (hit collection of this event) object. The hit collection data is updated in the ProcessHits(G4Step* aStep) method of BeamTestEmCalorimeter, where the energy deposited for that step is accessed, and the appropriate BeamTestEmCalorimeterHit object updated.

The implementation is shown below.

BeamTestEmCalorimeter.cc
----- snipped -----
#include "G4HCofThisEvent.hh"
// HandsOn5: Hit collection
#include "G4SDManager.hh"
#include "G4Step.hh"
#include "G4TouchableHistory.hh"
#include "G4Track.hh"
----- snipped -----
void BeamTestEmCalorimeter::Initialize(G4HCofThisEvent* hitsCollectionOfThisEvent)
{
// HandsOn5: Creating hit collection
// Create a new collection
fHitsCollection =
new BeamTestEmCalorimeterHitsCollection(SensitiveDetectorName, collectionName[0]);
if(fHitsCollectionID < 0)
fHitsCollectionID = G4SDManager::GetSDMpointer()->GetCollectionID(fHitsCollection);
// Add collection to the event
hitsCollectionOfThisEvent->AddHitsCollection(fHitsCollectionID, fHitsCollection);
// Initialise hits
G4int i(0);
for (i=0; i<100; i++) {
BeamTestEmCalorimeterHit* aHit = new BeamTestEmCalorimeterHit(i);
fHitsCollection->insert(aHit);
}
}
----- snipped -----
G4bool BeamTestEmCalorimeter::ProcessHits(G4Step* aStep, G4TouchableHistory*)
{
// HandsOn5: Accumulating hit data
// Get energy deposited in this step
G4double depositedEnergy = aStep->GetTotalEnergyDeposit();
if (0 == depositedEnergy) return true;
// Get volume and copy number
G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
G4TouchableHistory* theTouchable = (G4TouchableHistory*)(preStepPoint->GetTouchable());
G4VPhysicalVolume* thePhysical = theTouchable->GetVolume();
G4int copyNo = thePhysical->GetCopyNo();
// Get corresponding hit
BeamTestEmCalorimeterHit* aHit = (*fHitsCollection)[copyNo];
// Check to see if this is the first time the hit has been updated
if (!(aHit->GetLogicalVolume())) {
// Set volume information
aHit->SetLogicalVolume(thePhysical->GetLogicalVolume());
G4AffineTransform aTrans = theTouchable->GetHistory()->GetTopTransform();
aTrans.Invert();
aHit->SetRotation(aTrans.NetRotation());
aHit->SetPosition(aTrans.NetTranslation());
}
// Accumulate energy deposition
aHit->AddDepositedEnergy(depositedEnergy);
return true;
}

After implementing the above, check that it compiles. The sensitive detector also needs to be registered with the sensitive detector manager and attached to the CsI cell volume, as shown below.

BeamTestDetectorConstruction.cc
----- snipped -----
#include "G4Tubs.hh"
#include "G4VisAttributes.hh"
// HandsOn5: Defining sensitive detector
#include "G4SDManager.hh"
----- snipped -----
new G4PVParameterised("Cell_Physical", // Name
cellLogical, // Logical volume
calorimeterLogical, // Mother volume
kXAxis, // Axis
100, // Number of replicas
cellParam); // Parameterisation
////////////////////////////////////////////////////////////////////////
// HandsOn5: Defining sensitive detector
// Create a new BetamTestEmCalorimeter sensitive detector
G4VSensitiveDetector* detector = new BeamTestEmCalorimeter("Calorimeter");
// Get pointer to detector manager
G4SDManager* SDman = G4SDManager::GetSDMpointer();
// Register detector with manager
SDman->AddNewDetector(detector);
// Attach detector to volume defining calorimeter cells
cellLogical->SetSensitiveDetector(detector);
////////////////////////////////////////////////////////////////////////

After implementing the above, check that the program compiles.

User Defined Event Action Class

We can use a user action event class to access the hit data accumulated over a single event. The example includes a fully implemented class called BeamTestEventAction. This class inherits from G4UserEventAction. In the EndOfEventAction(const G4Event* event) method, the hit collection created above is accessed. The total energy deposited over all cells is summed and printed out.

We still need to register BeamTestEventAction with the run manager. The implementation is shown below.

beamTest.cc
----- snipped -----
#include "BeamTestDetectorConstruction.hh"
// HandsOn5: User defined event action
#include "BeamTestEventAction.hh"
----- snipped -----
// User action classes
runManager->SetUserAction(new BeamTestPrimaryGeneratorAction());
// HandsOn5: User defined event action
runManager->SetUserAction(new BeamTestEventAction);

After implementing the above, compile, link and run the program. Check that the total energy deposited in the calorimeter is printed out for each event as shown below.

Printout showing total energy deposited in the calorimeter
Voxelisation: top memory users:
Percent Memory Heads Nodes Pointers Total CPU Volume
------
55.48 0k 1 19 21 0.00 Calorimeter_Logical
44.52 0k 3 9 32 0.00 World_Logical
Start Run processing.
Energy deposited in calorimeter 0.963627 MeV
HepRepFile writing to G4Data1.heprep
Energy deposited in calorimeter 0 MeV
Energy deposited in calorimeter 0 MeV
Energy deposited in calorimeter 0.19467729 MeV
Energy deposited in calorimeter 3.4543694 MeV
Energy deposited in calorimeter 0 MeV
Energy deposited in calorimeter 0.11993501 MeV
Energy deposited in calorimeter 1.5138951 MeV
Energy deposited in calorimeter 1.2236197 MeV
Energy deposited in calorimeter 0 MeV
Energy deposited in calorimeter 0.035957673 MeV
Energy deposited in calorimeter 0.16419804 MeV
Energy deposited in calorimeter 0.75760033 MeV
Energy deposited in calorimeter 0.8852363 MeV
Energy deposited in calorimeter 0 MeV
Energy deposited in calorimeter 0 MeV
Energy deposited in calorimeter 0.12932751 MeV
Energy deposited in calorimeter 1.0485329 MeV
Energy deposited in calorimeter 0 MeV
Energy deposited in calorimeter 0 MeV
Run terminated.

Hit Visualisation

Hit visualization can be turned on through interactive commands. However, it is necessary to implement the Draw method of your hit class if you want to get any actual output. The example shown below draws a 5cm*5cm, magenta coloured tower, with a depth proportional to the deposited energy, on the end face of each CsI cell in which energy has been deposited.

BeamTestEmCalorimeterHit.cc
----- snipped -----
#include "G4UIcommand.hh"
#include "G4UnitsTable.hh"
// HandsOn5: Draw box
#include "G4Box.hh"
#include "G4Colour.hh"
#include "G4ParticleGun.hh"
#include "G4VisAttributes.hh"
----- snipped -----
void BeamTestEmCalorimeterHit::Draw()
{
G4VVisManager* pVVisManager = G4VVisManager::GetConcreteInstance();
if(pVVisManager & (fDepositedEnergy>0.)) {
// HandsOn5: Draw a box with depth propotional to the energy deposition
G4double scale = BeamTestPrimaryGeneratorAction::Gun()->GetParticleEnergy();
G4double depth = (50.*cm)*(fDepositedEnergy*MeV)/(scale*MeV);
// Back face of box is flat against front face of calorimeter cell
double z = fPosition.z() + 25.*cm;
G4ThreeVector myPos(fPosition.x(), fPosition.y(), z+depth);
G4Transform3D trans(fRotation.inverse(), myPos);
G4VisAttributes attribs;
// Magenta with transparency
G4Colour colour(1., 0., 1., 0.6);
attribs.SetColour(colour);
attribs.SetForceSolid(true);
G4Box box("MyBox", 5.*cm, 5.*cm, depth);
pVVisManager->Draw(box, attribs, trans);
}
}

After implementing the above, make the following modifications to run.mac to activate hit drawing and increase the primary electron energy to 500 MeV (to get a more interesting shower in the calorimeter).

.

run.mac
----- snipped -----
# Add trajectories to the visualization.
/vis/scene/add/trajectories
# HandsOn5: add hits to scene
/vis/scene/add/hits
# Accumulate multiple events in one picture.
#/vis/scene/endOfEventAction accumulate
# Trajectory colouring scheme
/vis/modeling/trajectories/create/drawByCharge
/vis/modeling/trajectories/drawByCharge-0/set -1 blue
/vis/modeling/trajectories/drawByCharge-0/set 1 blue
/vis/modeling/trajectories/drawByCharge-0/set 0 red
#HandsOn5: change primary particle energy
/gun/energy 500 MeV

Compile, link and run the program for around 20 events. A .heprep file should be now be produced for each event. Using Wired to view the results, you should see something like the picture shown below.

You can also can also view output using the OGLIX, OGLSX or OGLXm drivers. The picture below was produced using OGLXm.

You can download the complete source of this exercise from HandsOn5_complete.tgz

Geant4 Tutorial Course

Gean4.v8.0p01

May 2006

Jane Tinslay – Heavily based on a tutorial given at a SLAC in March 2006 by Tsukasa Aso