Hands-on 4   Scoring

 

l    Introduction

l    Installation

l    Scoring geometry implementation

l    Scorer definition

l    User defined run class

l    Applying filters

 

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:

 

l       Implement a spherical parameterised geometry to be used for scoring

l       Construct a sensitive detector (G4MultiFunctionalDetector) and attach scorers

l       Implement a user defined run class

l       Construct scorer filters based on particle type and energy

 

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:

 

l       Vacuum beam pipe

l       Titanium beam window

l       Iron evacuation chamber window

l       Silicon monitor

l       Beryllium target

 

Electrons are generated with an energy of 15 MeV in the beam pipe. The electrons interact with the target and generate bremsstrahlung photons. A spherical, parameterised, air filled scoring geometry is also implemented. The number of photons crossing this scoring volume is counted in ten theta bins and ten energy bins. The diagram below demonstrates the simulated experiment.

 

 

 

The following files are provided with the exercise:

 

l       beamTest.cc - main program

l       BeamTestDetectorConstruction – material, geometry , sensitive detector and scorer definitions

l       BeamTestPrimaryGeneratorAction - primary particle generator

l       BeamTestPhysicsList - user defined physics list

l       BeamTestRun - user defined run class

l       BeamTestRunAction - G4UserRunAction class

l       BeamTestScoreParameterisation – scoring volume parameterisation

 

All the geometries except the spherical scoring geometry are implemented. The scoring geometry, sensitive detector, scorers and filters 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 HandsOn4.tgz to your working directory.

 

Compile and link it using following commands:

 

     $ cd HandsOn4

     $ 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 scoring volume have been implemented.

 

 

Scoring Geometry Implementation

 

The mother scoring volume is defined as a thin air filled sphere. It has an inner radius of 1 m, an outer radius of 1.001 m, and is centered in the world volume. Ten scoring rings are placed within the mother scoring volume. Each ring covers 6 degrees in theta and is constructed of air with an inner radius of 1m and an outer radius of 1.001 m. G4PVPlacement is used to construct the mother sphere, while G4PVParameterised is used to construct the rings. BeamTestScoreParameterisation  is used for the ring parameterisation.

 

The implementation, along with visualisation attributes is shown below.

 

BeamTestDetectorConstruction.cc

----- snipped -----

 

#include "G4Tubs.hh"

#include "G4VisAttributes.hh"

 

// HandsOn4: Scoring volumes

#include "BeamTestScoreParameterisation.hh"

#include "G4PVParameterised.hh"

#include "G4Sphere.hh"

----- snipped -----

 

  ////////////////////////////////////////////////////////////////////////

  // HandsOn4: Scoring volumes

  // Mother sphere

  G4double innerRadius = 1.0*m;

  G4double outerRadius = 1.001*m;

 

  G4VSolid* sphereSolid = new G4Sphere("Sphere_Solid",   // Name

                                         innerRadius,    // Inner radius

                                         outerRadius,    // Outer radius

                                         0.*deg,         // Starting phi

                                         360.*deg,       // Delta phi

                                         0.*deg,         // Starting theta

                                         180.*deg);      // Delta theta

 

  G4LogicalVolume* sphereLogical =

    new G4LogicalVolume(sphereSolid, air, "Sphere_Logical");

 

 

  new G4PVPlacement(0, G4ThreeVector(), sphereLogical, "Sphere_Physical",

                     fpWorldLogical, false, 0);

 

  // 10 scoring rings

  G4double deltaTheta = 6.*deg;

 

  G4VSolid* scoreSolid = new G4Sphere("Score_Solid",   // Name

                                        innerRadius,   // Inner radius

                                        outerRadius,   // Outer radius

                                        0.*deg,        // Starting phi

                                        360.*deg,      // Delta phi    

                                        0.*deg,        // Starting theta

                                        deltaTheta);   // Delta theta

 

  G4LogicalVolume* fScoreLogical = 

    new G4LogicalVolume(scoreSolid, air, "scoreLog");

 

  BeamTestScoreParameterisation* param =

    new BeamTestScoreParameterisation(innerRadius, outerRadius, deltaTheta);

 

  new G4PVParameterised("scorePhys",   // Name

                         fScoreLogical, // Logical volume

                         sphereLogical, // Mother volume

                         kZAxis,        // Axis

                         10,            // Number of replicas

                         param);        // Parameterisation

----- snipped -----

 

  // HandsOn4: Scoring volume attributes

  // Mother sphere

  G4VisAttributes* sphereAttributes =

        new G4VisAttributes(G4Colour(1.0,1.0,1.0,0.2));

  sphereAttributes->SetVisibility(false);

  sphereLogical->SetVisAttributes(sphereAttributes);

 

  // Ring scoring volumes. Green with transparancy.

  G4VisAttributes* scoreAttributes =

        new G4VisAttributes(G4Colour(0.0,1.0,0.0,0.2));

  scoreAttributes->SetVisibility(true);

  fScoreLogical->SetVisAttributes(scoreAttributes);

 

 

Implement the above and compile, link and run the program. Use Wired to check that all the scoring geometries are now in place. You should be able to see something like:

 

 

Note that the mother scoring sphere has been made invisible by default. You can change this by checking the Sphere_Physical[0] box.

 

 

Scorer Definition

 

We are interested in counting the number of photons per unit area which cross the inner surface of each scoring ring. To do this we need to:

 

l       Create a G4MultiFunctionalDetector detector

l       Register the detector with the sensitive detector manager (G4SDManager)

l       Attach the detector to the parameterised scoring volume

l       Create a primitive scorer for scoring surface current

l       Register the primitive scorer with the detector

 

 We will use the G4PSSphereSurfaceCurrent scorer to count the number of photons crossing each scoring ring. Since the ten scoring rings are replicas, we only need to create one scorer. This scorer will automatically record the current as a function of replica copy number.

 

The implementation is shown below.

 

BeamTestDetectorConstruction.cc

----- snipped -----

 

// HandsOn4: Scoring components

#include "G4MultiFunctionalDetector.hh"

#include "G4PSSphereSurfaceCurrent.hh"

#include "G4SDManager.hh"

----- snipped -----

 

  new G4PVParameterised("scorePhys",    // Name

                         fScoreLogical, // Logical volume

                         sphereLogical, // Mother volume

                         kZAxis,        // Axis

                         10,            // Number of replicas

                         param);        // Parameterisation

 

  // HandsOn4: Construct scoring components

  SetupScoring(fScoreLogical);

 

  ////////////////////////////////////////////////////////////////////////

  // Target  (Default parameters for Be )

  G4Material* beryllium = G4Material::GetMaterial("Beryllium");

----- snipped -----

 

void BeamTestDetectorConstruction::SetupScoring(G4LogicalVolume* scoringVolume)

{

  // HandsOn4: Construct scoring components

 

  // Create a new sensitive detector named "MyDetector"

  G4MultiFunctionalDetector* detector =

    new G4MultiFunctionalDetector("MyDetector");

 

  // Get pointer to detector manager

  G4SDManager* manager = G4SDManager::GetSDMpointer();  

 

  // Register detector with manager

  manager->AddNewDetector(detector);

 

  // Attach detector to scoring volume

  scoringVolume->SetSensitiveDetector(detector);

 

  // Create a primitive Scorer named myScorer

  G4PSSphereSurfaceCurrent* scorer = new G4PSSphereSurfaceCurrent("myScorer",

                                                                   fCurrent_In);

 

  // Register scorer with detector

  detector->RegisterPrimitive(scorer);

 

  G4cout<<"Created G4MultiFunctionalDetector named "

        <<detector->GetName()<<", and a G4PSSphereSurfaceCurrent scorer "

        <<"named "<<scorer->GetName()<<G4endl;

}

 

 

 

After implementing the above, compile, link and run the program. Check that the sensitive detector and scorer names are printed out, as shown below.

 

Printout showing detector and scorer details

----- snipped -----

Material:   Vacuum     density:  0.010 mg/cm3  temperature: 273.15 K  pressure:   0.02 atm  RadLength:  36.786 km

   --->  Element: Nitrogen (N)   Z =  7.0   N =  14.0   A =  14.01 g/mole  fractionMass:  70.00 %  Abundance  72.71 %

   --->  Element: Oxygen (O)   Z =  8.0   N =  16.0   A =  16.00 g/mole  fractionMass:  30.00 %  Abundance  27.29 %

 

 

Created G4MultiFunctionalDetector named MyDetector, and a G4PSSphereSurfaceCurrent scorer named myScorer

BeamTestPhysicsList::SetCuts:CutLength : 1 mm

 

----- snipped -----

 

 

User Defined Run Class

 

The scorer implemented in the previous section will accumulate data over a single event, but we are more interested in data accumulated over an entire run. Therefore, we need to create our own run class, inheriting from G4Run. This class will be responsible for merging the scorer data for all events.

 

A partially implemented user run class called BeamTestRun is provided with the example. Also included is a BeamTestRunAction class which inherits from G4UserRunAction. BeamTestRunAction is responsible for creating an instance of BeamTestRun, and is an example of one of the optional user action classes. BeamTestRunAction needs to be registered with the run manager in the main program, as shown below.

 

beamTest.cc

----- snipped -----

 

#include "BeamTestPrimaryGeneratorAction.hh"

// HandsOn4: User defined run action

#include "BeamTestRunAction.hh"

#include "G4RunManager.hh"

----- snipped -----

 

  // User action classes

  runManager->SetUserAction(new BeamTestPrimaryGeneratorAction());

  // HandsOn4: User defined run action

  runManager->SetUserAction(new BeamTestRunAction);

 

  // Initialize G4 kernel

 

 

After implementing the above, compile, link and run the program. Check that you get printout shown below.

 

Printout showing user defined BeamTestRun class creation

Voxelisation: top memory users:

    Percent     Memory      Heads    Nodes   Pointers    Total CPU    Volume

    -------   --------     ------   ------   --------   ----------    ----------

      67.84          0k         3       14         36         0.00    World_Logical

      32.16          0k         1        7         20         0.00    Sphere_Logical

Creating user define run class BeamTestRun

Start Run processing.

HepRepFile writing to G4Data1.heprep

 

The event by event scorer data merging still needs to be implemented in BeamTestRun. This class has a data member, fMap, which is used to hold the accumulated data. fMap is just a map between hit collection idÕs and scorer data (G4THitsMap<G4double>). Merging of scorer data is done in the BeamTestRun::RecordEvent(const G4Event* anEvent) method, and involves:

 

l       Accessing the hit collection for the event (G4HCofThisEvent object)

l       Using the G4HCofThisEvent object to access the G4THitsMap<double> data for each registered scorer

l       Accumulating the hit map data

 

 The implementation is shown below.

 

BeamTestRun.cc

----- snipped -----

 

void BeamTestRun::RecordEvent(const G4Event* anEvent)

{

   numberOfEvent++;

 

  // HandsOn4: Merging event data

  // Get the hits collection

  G4HCofThisEvent* eventHitCollection = anEvent->GetHCofThisEvent();

 

  if (!eventHitCollection) return;

 

  // Update our private fMap

  std::map< G4int, G4THitsMap<G4double>* >::iterator iter = fMap.begin();

 

  while (iter != fMap.end()) {

    G4int id = iter->first;

   

    // Get the hit collection corresponding to "id"

    G4THitsMap<G4double>* eventHitsMap

      = dynamic_cast< G4THitsMap<G4double>* >(eventHitCollection->GetHC(id));

 

    // Expect this to exist

    assert (0 != eventHitsMap);

 

    // Accumulate event data into our G4THitsMap<G4double> map

    *(iter->second) += *eventHitsMap;

 

    iter++;

  }

}

 

After implementing the above, modify the run.mac macro to run for 1000 events. If you donÕt want to visualise the results, disable the visualisation system by adding the command /vis/disable to run.mac. This will also speed up the processing.

 

Compile, link and run the program. Check that you get data printed out as shown below. Theta bins corresponding to the ten scoring rings are listed in the first column. The corresponding number of photons per unit area (mm-2) crossing the scoring ring are listed the second column.

 

Printout showing current printout in theta bin

Number of Events Processed:1000 events.

   Theta myScorer

       0  0.00227

       1  0.00115

       2 0.000476

       3 0.000327

       4 0.000244

       5 0.000162

       6 0.000111

       7 9.89e-05

       8 9.78e-05

       9 8.34e-05

Graphics systems deleted.

Visualization Manager deleting...

 

 

 

Applying Filters

 

So far we have only used one scorer to count all tracks crossing each of the scoring rings. It is possible to attach a filter to a scorer, so that only hits passing the filter will be processed.  For example, if we are only interested in scoring photons within a given energy range, we would attach a suitably configured G4SDParticleWithEnergyFilter filter to the scorer.  If we are interested in scoring a number of different energy ranges we can create multiple scorers, each with its own filter. This is demonstrated below, where ten scorers filtering on photons  in 1.2 MeV energy bins are generated and registered with the G4MultiFunctionalDetector.

 

BeamTestDetectorConstruction.cc

----- snipped -----

 

// HandsOn4: Filter

#include "G4SDParticleWithEnergyFilter.hh"

----- snipped -----

 

void BeamTestDetectorConstruction::SetupScoring(G4LogicalVolume* scoringVolume)

{

  // HandsOn4: Construct scoring components

 

  // Create a new sensitive detector named "MyDetector"

  G4MultiFunctionalDetector* detector =

    new G4MultiFunctionalDetector("MyDetector");

 

  // Get pointer to detector manager

  G4SDManager* manager = G4SDManager::GetSDMpointer();  

 

  // Register detector with manager

  manager->AddNewDetector(detector);

 

  // Attach detector to scoring volume

  scoringVolume->SetSensitiveDetector(detector);

 

  // HandsOn4: Applying filters to scorers

  G4int nEnergyBins = 10;

  G4int i(0);

   

  for (i=0; i<nEnergyBins; i++) {

   

    char name[10];

    std::sprintf(name,"Scorer%i", i);

 

    // Create a primitive Scorer

    G4PSSphereSurfaceCurrent* scorer =

      new G4PSSphereSurfaceCurrent(name, fCurrent_In);

  

    // Create a filter

    G4SDParticleWithEnergyFilter* filter =

      new G4SDParticleWithEnergyFilter("Gamma Filter");

 

    // Filter will pass only gammas

    filter->add("gamma");

 

    // Set energy range of gammas filter will accept

    G4double minEnergy = i*1.2*MeV;

    G4double maxEnergy = (i+1)*1.2*MeV;

    filter->SetKineticEnergy(minEnergy, maxEnergy);

 

    // Attach filter to scorer. Scorer will only process hits which pass

    // the filter

    scorer->SetFilter(filter);

   

    G4cout<<"Created scorer with name "<<name;

    G4cout<<", covering energy range "<<minEnergy*MeV;

    G4cout<<" : "<<maxEnergy*MeV<<" MeV"<<G4endl;

   

    // Register scorer with detector

    detector->RegisterPrimitive(scorer);

  }

}

 

 

After implementing the above, compile, link and run the program for around 10000 events. Check that you see the new scorers registered, and that you see the current printed out for each theta and energy bin, as shown below.

 

Printout showing current in each theta and energy bin

----- snipped -----

 

   --->  Element: Nitrogen (N)   Z =  7.0   N =  14.0   A =  14.01 g/mole  fractionMass:  70.00 %  Abundance  72.71 %

   --->  Element: Oxygen (O)   Z =  8.0   N =  16.0   A =  16.00 g/mole  fractionMass:  30.00 %  Abundance  27.29 %

 

 

Created scorer with name Scorer0, covering energy range 0 : 1.2 MeV

Created scorer with name Scorer1, covering energy range 1.2 : 2.4 MeV

Created scorer with name Scorer2, covering energy range 2.4 : 3.6 MeV

Created scorer with name Scorer3, covering energy range 3.6 : 4.8 MeV

Created scorer with name Scorer4, covering energy range 4.8 : 6 MeV

Created scorer with name Scorer5, covering energy range 6 : 7.2 MeV

Created scorer with name Scorer6, covering energy range 7.2 : 8.4 MeV

Created scorer with name Scorer7, covering energy range 8.4 : 9.6 MeV

Created scorer with name Scorer8, covering energy range 9.6 : 10.8 MeV

Created scorer with name Scorer9, covering energy range 10.8 : 12 MeV

BeamTestPhysicsList::SetCuts:CutLength : 1 mm

----- snipped -----

 

Number of Events Processed:10000 events.

   Theta  Scorer0  Scorer1  Scorer2  Scorer3  Scorer4  Scorer5  Scorer6  Scorer7  Scorer8  Scorer9

       0   0.0132  0.00378  0.00215  0.00093 0.000581 0.000436 0.000407 0.000378 0.000349 0.000174

       1  0.00724  0.00145 0.000573 0.000457 0.000185 0.000214 0.000126 0.000156 5.83e-05 1.94e-05

       2   0.0041 0.000693  0.00027   0.0002 0.000106 0.000112 4.11e-05 2.35e-05 1.76e-05 1.76e-05

       3  0.00272  0.00031  0.00014 6.36e-05 2.97e-05 2.55e-05 2.55e-05 4.24e-06 8.49e-06 4.24e-06

       4  0.00201 0.000191 6.03e-05 2.34e-05 1.34e-05 2.01e-05        0    1e-05 3.35e-06        0

       5  0.00138 0.000131 4.75e-05 1.68e-05 2.79e-06 5.58e-06        0        0        0        0

       6  0.00113 8.94e-05 3.14e-05 9.66e-06 4.83e-06        0        0 2.42e-06        0        0

       7  0.00101 4.52e-05 6.45e-06  4.3e-06        0        0        0        0        0 2.15e-06

       8 0.000927 4.11e-05 3.91e-06 3.91e-06 1.96e-06 1.96e-06        0        0        0        0

       9 0.000722 2.72e-05 7.25e-06 1.81e-06        0        0        0        0        0        0

Graphics systems deleted.


An output file called output.csv (Comma Separated Variables) should also be produced. This can be loaded into for example, OpenOffice or Excel, to make a plot of the data. The plot below shows the data for 104 events plotted for a selection of energy bins.

 

 

You can download the complete source of this exercise from HandsOn4_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