Menu

Subscribe

[Tutorial] GEANT4 Introduction: Basics & Geometry

on research

Weihe Zeng
Tsinghua University
Nov. 27, 2017

Introduction to Monte Carlo

Monte Carlo methods (or Monte Carlo experiments) are a broad class of computational algorithms that rely on repeated random sampling to obtain numerical results. Their essential idea is using randomness to solve problems that might be deterministic in principle. (wikipedia)

casino

A Small Discussion

150%

Picture from wikipedia

Applications

  • High Energy / Nuclear Physics
  • Nuclear Reactor Design
  • Radiation Protection
  • Medical imaging and radiotherapy ...

60%

atlas@LHC

Building Blocks of Physics Simulation

  • Virtual Reality (World): (3 dimensional) geometry, material, environment;
  • Incident particles (Player): radioactive sources, high energy collision;
  • Physics Model (Rules): Physical laws, cross section;
// examples/basic/B1/exampleB1.cc
// Detector construction
runManager->SetUserInitialization(new B1DetectorConstruction());

// Physics list
G4VModularPhysicsList* physicsList = new QBBC;  
physicsList->SetVerboseLevel(1);  
runManager->SetUserInitialization(physicsList);

// User action initialization
runManager->SetUserInitialization(new B1ActionInitialization());  
/* in B1ActionInitialization 
  SetUserAction(new B1PrimaryGeneratorAction); */

A Simple Simulation

60%

World: Infinite space; Particles with random directions and velocity; The only interaction: annihilation of constant cross-section (or mean free path).

Intensity of particles after traveling a distance of $x$:

$$ I = I_{0} e^{-x / l} $$

Probability of annihilate at the distance $x$:

$$f(x) = \dfrac{1}{l} e^{-x/l}$$

Getting Started with GEANT4

Textbook: Page 2-48

int main(){  
  // construct the default run manager
  G4RunManager* runManager = new G4RunManager;
  // set mandatory initialization classes
  runManager->SetUserInitialization(new ExG4DetectorConstruction01);
  runManager->SetUserInitialization(new ExG4PhysicsList00);
  runManager->SetUserInitialization(new ExG4ActionInitialization01);
  // initialize G4 kernel
  runManager->Initialize();
  // get the pointer to the UI manager and set verbosities
  G4UImanager* UI = G4UImanager::GetUIpointer();
  UI->ApplyCommand("/run/verbose 1");
  UI->ApplyCommand("/event/verbose 1");
  UI->ApplyCommand("/tracking/verbose 1");
  // start a run
  int numberOfEvent = 3;
  runManager->BeamOn(numberOfEvent);
  // job termination
  delete runManager;
  return 0;
}

Geometry in GEANT4

Textbook: Chapter 4, page 90-167

  • Basic components
    • System of units
    • Material
    • Geometry
    • Transformation
  • How to define Materials
  • How to define Solids: G4Box, G4Tube as examples.
  • How to define Logical Volumes
  • How to define Physical Volumes and repeated structure
  • others: electromagnetic field, sensitive detector ...

Class Reference: http://www-geant4.kek.jp/Reference/10.03.p03

Basic components

System of Units Textbook: Chapter 3.3, Page 56-58

| Name | Basic Unit | | ----------------------- | ---------- | | millimeter | mm | | nanosecond | ns | | Mega electron Volt | MeV | | positron charge | eplus | | degree Kelvin | kelvin | | the amount of substance | mole | | radian | radian | | steradian | steradian |

#include "SystemOfUnits.h"
G4double Size=15*km, KineticEnergy=90.3*GeV, density=11*mg/cm3;

Material

Textbook: Chapter 4.2, page 140-145

  • Isotopes: the properties of atoms, $^{3}$H, $^{235}$U, $^{76}$Ge
    • atomic number,
    • number of nucleons.
  • Elements: the properties of elements, Uranium, Germanium, Lead
    • effective atomic number,
    • effective mass of a mole,
    • the number of isotopes ...
  • Materials: the macroscopic properties of matter, Air, Water
    • density,
    • physical state,
    • temperature,
    • pressure ...

Geometry

70%

  • Logical volumes (P110-112): Solid + Material
  • Physical volumes (P112-121): Logical volume + Spatial positioning

Geometry Hierarchy

Mother and daughter volumes:

  • A volume is placed in its mother volume (except World volume)
  • One or more volumes can be placed in a mother volume
  • Position and rotation of the daughter volume is described with respect to the local coordinate system of the mother volume

Overlap:

  • Daughter volumes cannot protrude from the mother volume
  • Daughter volumes cannot overlap

motherdaughter

Transformation

  • Position / Translation
G4ThreeVector ZTrans(0, 0, 50*cm);
  • Rotation
G4RotationMatrix* yRot = new G4RotationMatrix;  
// Rotates X and Z axes only
// Rotates 45 degrees
yRot->rotateY(PI/4.*rad);
  • Transform 3D
G4Transform3D( G4RotationMatrix &pRot,  
          const G4ThreeVector &tlate)

Define Materials

Build materials from scratch:

// define Elements
a = 1.01*g/mole;  
G4Element* elH = new G4Element(name="Hydrogen",  
                               symbol="H" , z= 1., a);
// define an Element from isotopes, by relative abundance
G4Isotope* U5 = new G4Isotope(name="U235", iz=92, n=235, a=235.01*g/mole);  
G4Isotope* U8 = new G4Isotope(name="U238", iz=92, n=238, a=238.03*g/mole);  
G4Element* elU = new G4Element(name="enriched Uranium", symbol="U", ncomponents=2);  
elU->AddIsotope(U5, abundance= 90.*perCent);  
elU->AddIsotope(U8, abundance= 10.*perCent);

// define a material from elements. case 1: chemical molecule
density = 1.000*g/cm3;  
G4Material* H2O = new G4Material(name="Water", density, ncomponents=2);  
H2O->AddElement(elH, natoms=2);  
H2O->AddElement(elO, natoms=1);  
// define a material from elements. case 2: mixture by fractional mass
density = 1.290*mg/cm3;  
G4Material* Air = new G4Material(name="Air " , density, ncomponents=2);  
Air->AddElement(elN, fractionmass=0.7);  
Air->AddElement(elO, fractionmass=0.3);

Use the internal Geant4 database (Textbook: Appendix 6, page 368-381)

#include "G4NistManager.hh"
G4NistManager* man = G4NistManager::Instance();  
man->SetVerbose(1);  
// define elements
G4Element* C = man->FindOrBuildElement("C");  
G4Element* Pb = man->FindOrBuildMaterial("Pb");  
// define pure NIST materials
G4Material* Al = man->FindOrBuildMaterial("G4_Al");  
G4Material* Cu = man->FindOrBuildMaterial("G4_Cu");  
// define NIST materials
G4Material* H2O = man->FindOrBuildMaterial("G4_WATER");  
G4Material* Sci = man->FindOrBuildMaterial("G4_PLASTIC_SC_VINYLTOLUENE");  
G4Material* SiO2 = man->FindOrBuildMaterial("G4_SILICON_DIOXIDE");  
G4Material* Air = man->FindOrBuildMaterial("G4_AIR");  
// HEP materials
G4Material* PbWO4 = man->FindOrBuildMaterial("G4_PbWO4");  
G4Material* lAr = man->FindOrBuildMaterial("G4_lAr");  
G4Material* vac = man->FindOrBuildMaterial("G4_Galactic");  
// define gas material at non STP conditions (T = 120K, P=0.5atm)
G4Material* coldAr = man->ConstructNewGasdMaterial(  
                 "ColdAr","G4_Ar",120.*kelvin,0.5*atmosphere);

Define Solids

G4Box(const G4String& pName,  
            G4double pX,    // half length
            G4double pY,    // half width
            G4double pZ)    // half hight
G4Tubs(const G4String& pName,  
              G4double pRMin,
              G4double pRMax,
              G4double pDz,
              G4double pSPhi,
              G4double pDPhi)

tube

Boolean Operations

G4Box* box = new G4Box("Box",20*mm,30*mm,40*mm);  
G4Tubs* cyl = new G4Tubs("Cylinder",0,50*mm,50*mm,0,twopi);

G4UnionSolid* union =  
  new G4UnionSolid("Box+Cylinder", box, cyl);
G4IntersectionSolid* intersection =  
  new G4IntersectionSolid("Box*Cylinder", box, cyl);
G4SubtractionSolid* subtraction =  
  new G4SubtractionSolid("Box-Cylinder", box, cyl);

// transform cyl's coordinate system
G4UnionSolid* unionMoved =  
  new G4UnionSolid("Box+CylinderMoved", box, cyl, yRot, zTrans);

Note: The constituent solids of a Boolean operation should possibly avoid be composed by sharing all or part of their surfaces. Moreover, the final Boolean solid should represent a single 'closed' solid.

Define Logical Volumes

G4LogicalVolume(  G4VSolid*       pSolid,  
                  G4Material*     pMaterial,
                  const G4String& Name,
                  G4FieldManager* pFieldMgr=0,
                  G4VSensitiveDetector* pSDetector=0,
                  G4UserLimits*         pULimits=0,
                  G4bool                Optimise=true )

For example. In B1DetectorConstruction.cc

G4LogicalVolume* logicWorld =  
    new G4LogicalVolume(solidWorld,          //its solid
                        world_mat,           //its material
                        "World");            //its name

Define Physical Volumes

G4PVPlacement(G4RotationMatrix*    pRot,  
              const G4ThreeVector& tlate,
              G4LogicalVolume* pCurrentLogical,
              const G4String&  pName,
              G4LogicalVolume* pMotherLogical, // 0 for world
              G4bool           pMany,
              G4int            pCopyNo,
              G4bool           pSurfChk=false )

G4VPhysicalVolume* physEnv = new G4PVPlacement(  
          0,                       //no rotation
        G4ThreeVector(10,10,5),  //at (10,10,5)
        logicEnv,                //its logical volume
        "Envelope",              //its name
        logicWorld,              //its mother  volume
        false,                   //For future use
        0,                       //copy number
        checkOverlaps);          //overlaps checking

Note: One logical volume can be placed more than once. Position and rotation is described with respect to mother volume

Repeated Volumes

  • Construct PV inside a loop;
  • Replicas: volumes differing only in their positioning, and completely filling the containing mother volume. It MUST be the mother's only daughter.
G4PVReplica( const G4String& pName,  
      G4LogicalVolume* pCurrentLogical,
      G4LogicalVolume* pMotherLogical,
      const EAxis pAxis,
      const G4int nReplicas,
      const G4double width,
      const G4double offset=0 )
  • Parameterised Volumes: check out example/basic/B2/B2c

Others

  • Electromagnetic Field: Chapter 4.3, page 145-152;
// create global field
G4ThreeVector fieldValue = G4ThreeVector();  
fMagFieldMessenger = new G4GlobalMagFieldMessenger(fieldValue);  
// for ppart of the volume
logicVolume->SetFieldManager(localFieldManager, true);

  • Sensitive Detectors: Chapter 4.4, page 152-161;

  • Color of volumes: Chapter 8.6, page 294-299;

// Instantiation of a set of visualization attributes
G4VisAttributes * calTubeVisAtt = new  
  G4VisAttributes(G4Colour(0.,1.,1.));
// Assignment visualization attributes to the logical volume
myTargetLog->SetVisAttributes(calTubeVisAtt);

Review of example B1

See text editor...

class B1DetectorConstruction:public G4VUserDetectorConstruction  
{
  public:
    B1DetectorConstruction();
    virtual ~B1DetectorConstruction();
    virtual G4VPhysicalVolume* Construct();
};


G4VPhysicalVolume* B1DetectorConstruction::Construct()  
{  
  /* ... */
  return physWorld;
}

Visualization

Textbook: Chapter 8, page 252-311

The many graphics systems that Geant4 supports are complementary to each other.

  • ASCIITree
  • OpenGL and QT
  • OpenInventor
  • HepRep
  • DAWN
  • VRML
  • ...

Choose the driver that meets your needs: Chapter 8.1.3, page 254.

GEANT4 commands related to visualization: Chapter 8.4, 275-286.

The Simplest Geometry Printer: ASCIITree

/vis/open ASCIITree
/vis/drawVolume
/vis/ASCIITree/Verbose 4
/vis/viewer/flush

Output:

# example/basic/B1
"World":0 / "World" / "World"(G4Box), 20.736 L  , 1.20479 mg/cm3 (G4_AIR)
    "Envelope":0 / "Envelope" / "Envelope"(G4Box), 12 L  , 1 g/cm3  (G4_WATER)
      "Shape1":0 / "Shape1" / "Shape1"(G4Cons), 1.75929 dL , 1.127 g/cm3  (G4_A-150_TISSUE)
      "Shape2":0 / "Shape2" / "Shape2"(G4Trd), 9.36 dL , 1.85 g/cm3  (G4_BONE_COMPACT_ICRU)
Calculating mass(es)...  
Overall volume of "World":0, is 20.736 L   and the daughter-included mass to unlimited depth is 12.8285 kg

Real-time Visualization with GUI controls: QT

See example/basic/B1/vis.mac for details.

60%

File Based Visualization: VRML

/vis/open VRML2FILE
/vis/drawVolume
/vis/viewer/flush

Recommended viewers: freeWRL, BS contact

80%

Homework

Deadline: 19:00, Dec. 4, 2017

70%

Let's find out who is the GEANT4 artist. Please modify examples/basic/B1/src/B1DetectorConstruction.cc and model a human (or robot), which should has at least one head, one body, a pair of arms and legs.

Submission: the modified source file and a screenshot of your masterpiece.