Physics processes describe how particles interact with a material. Seven major categories of processes are provided by Geant4:
The generalization and abstraction of physics processes is a key issue in the design of Geant4. All physics processes are treated in the same manner from the tracking point of view. The Geant4 approach enables anyone to create a process and assign it to a particle type. This openness should allow the creation of processes for novel, domain-specific or customised purposes by individuals or groups of users.
Each process has two groups of methods which play an important
role in tracking, GetPhysicalInteractionLength
(GPIL) and
DoIt
. The GPIL method gives the step length from the
current space-time point to the next space-time point. It does this
by calculating the probability of interaction based on the
process's cross section information. At the end of this step the
DoIt
method should be invoked. The DoIt
method
implements the details of the interaction, changing the particle's
energy, momentum, direction and position, and producing secondary
tracks if required. These changes are recorded as
G4VParticleChange
objects(see
Particle Change).
G4VProcess
is the base class for all physics processes.
Each physics process must implement virtual methods of
G4VProcess
which describe the interaction (DoIt) and
determine when an interaction should occur (GPIL). In order to
accommodate various types of interactions G4VProcess
provides three DoIt
methods:
G4VParticleChange* AlongStepDoIt( const G4Track& track,
const G4Step& stepData )
This method is invoked while G4SteppingManager
is
transporting a particle through one step. The corresponding
AlongStepDoIt
for each defined process is applied for
every step regardless of which process produces the minimum step
length. Each resulting change to the track information is recorded
and accumulated in G4Step
. After all processes have been
invoked, changes due to AlongStepDoIt
are applied to
G4Track
, including the particle relocation and the safety
update. Note that after the invocation of AlongStepDoIt
,
the endpoint of the G4Track
object is in a new volume if the
step was limited by a geometric boundary. In order to obtain
information about the old volume, G4Step
must be accessed,
since it contains information about both endpoints of a step.
G4VParticleChange* PostStepDoIt( const G4Track& track,
const G4Step& stepData )
This method is invoked at the end point of a step, only if its
process has produced the minimum step length, or if the process is
forced to occur. G4Track
will be updated after each
invocation of PostStepDoIt
, in contrast to the
AlongStepDoIt
method.
G4VParticleChange* AtRestDoIt( const G4Track& track,
const G4Step& stepData )
This method is invoked only for stopped particles, and only if its process produced the minimum step length or the process is forced to occur.
For each of the above DoIt
methods G4VProcess
provides a corresponding pure virtual GPIL method:
G4double PostStepGetPhysicalInteractionLength( const
G4Track& track, G4double previousStepSize, G4ForceCondition*
condition )
This method generates the step length allowed by its process. It also provides a flag to force the interaction to occur regardless of its step length.
G4double AlongStepGetPhysicalInteractionLength( const
G4Track& track, G4double previousStepSize, G4double
currentMinimumStep, G4double& proposedSafety, G4GPILSelection*
selection )
This method generates the step length allowed by its process.
G4double AtRestGetPhysicalInteractionLength( const
G4Track& track, G4ForceCondition* condition )
This method generates the step length in time allowed by its process. It also provides a flag to force the interaction to occur regardless of its step length.
Other pure virtual methods in G4VProcess
follow:
virtual G4bool IsApplicable(const
G4ParticleDefinition&)
returns true if this process object is applicable to the particle type.
virtual void PreparePhysicsTable(const
G4ParticleDefinition&)
and
virtual void BuildPhysicsTable(const
G4ParticleDefinition&)
is messaged by the process manager, whenever cross section tables should be prepared and rebuilt due to changing cut-off values. It is not mandatory if the process is not affected by cut-off values.
virtual void StartTracking()
and
virtual void EndTracking()
are messaged by the tracking manager at the beginning and end of tracking the current track.
Specialized processes may be derived from seven additional
virtual base classes which are themselves derived from
G4VProcess
. Three of these classes are used for simple
processes:
G4VRestProcess
Processes using only the AtRestDoIt
method.
example: neutron capture
G4VDiscreteProcess
Processes using only the PostStepDoIt
method.
example: compton scattering, hadron inelastic interaction
The other four classes are provided for rather complex processes:
G4VContinuousDiscreteProcess
Processes using both AlongStepDoIt
and
PostStepDoIt
methods.
example: transportation, ionisation(energy loss and delta ray)
G4VRestDiscreteProcess
Processes using both AtRestDoIt
and
PostStepDoIt
methods.
example: positron annihilation, decay (both in flight and at rest)
G4VRestContinuousProcess
Processes using both AtRestDoIt
and
AlongStepDoIt
methods.
G4VRestContinuousDiscreteProcess
Processes using AtRestDoIt
,
AlongStepDoIt and
PostStepDoIt methods.
G4VParticleChange
and its descendants are used to store
the final state information of the track, including secondary
tracks, which has been generated by the DoIt
methods. The
instance of G4VParticleChange
is the only object whose
information is updated by the physics processes, hence it is
responsible for updating the step. The stepping manager collects
secondary tracks and only sends requests via particle change to
update G4Step
.
G4VParticleChange
is introduced as an abstract class. It
has a minimal set of methods for updating G4Step
and
handling secondaries. A physics process can therefore define its
own particle change derived from G4VParticleChange
. Three
pure virtual methods are provided,
virtual G4Step* UpdateStepForAtRest( G4Step* step)
,
virtual G4Step* UpdateStepForAlongStep( G4Step* step )
and
virtual G4Step* UpdateStepForPostStep( G4Step* step)
,
which correspond to the three DoIt
methods of
G4VProcess
. Each derived class should implement these
methods.
This section summarizes the electromagnetic (EM) physics processes which are provided with Geant4. Extended information are avalable at EM web pages. For details on the implementation of these processes please refer to the Physics Reference Manual.
To use the electromagnetic physics data files are needed. The user should set the environment variable G4LEDATA to the directory with this files. These files are distributed together with Geant4 and can be obtained via Geant4 download web page. For Geant4 version 10.1 G4EMLOW6.41 data set is required.
The following is a summary of the electromagnetic processes available in Geant4.
Photon processes
G4GammaConversion
)G4PhotoElectricEffect
)
G4ComptonScattering
)
G4RayleighScattering
)
G4GammaConversionToMuons
)
Electron/positron processes
G4eIonisation
)G4eBremsstrahlung
)G4eMultipleScattering
)G4eplusAnnihilation
)G4AnnihiToMuPair
)G4eeToHadrons
)
Muon processes
G4MuIonisation
)G4MuBremsstrahlung
) G4MuPairProduction
)G4MuMultipleScattering
)
Hadron/ion processes
G4hIonisation
)
G4ionIonisation
)
G4hhIonisation
)
G4mplIonisation
)
G4hMultipleScattering
)
G4hBremsstrahlung
)
G4hPairProduction
)
Coulomb scattering processes
G4CoulombScattering
)
G4ScreenedNuclearRecoil
)
Processes for simulation of polarized electron and gamma beams
G4PolarizedCompton
)
G4PolarizedGammaConversion
)
G4PolarizedPhotoElectricEffect
)
G4ePolarizedBremsstrahlung
)
G4ePolarizedIonisation
)
G4eplusPolarizedAnnihilation
)
Processes for simulation of X-rays and optical protons production by charged particles
G4SynchrotronRadiation
)G4TransitionRadiation
)G4Cerenkov
)G4Scintillation
)
The processes described above use physics model classes, which may be combined according to particle energy. It is possible to change the energy range over which different models are valid, and to apply other models specific to particle type, energy range, and G4Region. The following alternative models are available in the standard EM sub-library:
G4PAIModel
)G4PAIPhotModel
)G4BraggIonGasModel
)G4BetheBlochIonGasModel
)G4UrbanMscModel
)G4WentzelVIModel
)G4LowEWentzelVIModel
)
It is recommended to use physics constructor classes provided
with reference physics lists (in subdirectory source/physics_lists/constructors/electromagnetic
of the Geant4 source distribution):
default EM physics (class name G4EmStandardPhysics
)
optional EM physics providing fast but less acurate electron transport
due to
"Simple" method of step limitation by multiple scattering, reduced
step limitation by ionisation process and enabled "ApplyCuts" option
(class name G4EmStandardPhysics_option1
)
optional EM physics providing fast but less acurate electron transport
due to
"Simple" method of step limitation by multiple scattering and reduced
step limitation by ionisation process
(class name G4EmStandardPhysics_option2
)
EM physics for simulation with high accuracy due to
"UseDistanceToBoundary" multiple
scattering step limitation and usage of
G4UrbanMscModel
for all charged particles,
reduced finalRange
parameter of stepping function
optimized per particle type, alternative model
G4KleinNishinaModel
for Compton scattering,
enabled Rayleigh scattering, enabled fluorescence,
enabled nuclear stopping,
G4IonParameterisedLossModel
for ion
ionisation, 20 bins per energy decade of physics tables,
and 10 eV low-energy limit for tables
(class name G4EmStandardPhysics_option3
)
Combination of best EM models for simulation with high accuracy
includes "UseDistanceToBoundary" multiple
scattering step limitation,
RangeFactor = 0.02
,
reduced finalRange
parameter of stepping function
optimized per particle type,
enabled Rayleigh scattering, enabled fluorescence,
enabled nuclear stopping, enable accurate angular generator
for ionisation models,
G4LivermorePhotoElectricModel
,
G4LowEPComptonModel
below 20 MeV,
G4PenelopeGammaConversionModel
below 1 GeV,
G4PenelopeIonisationModel
fro electrons
and positrons below 100 keV,
G4IonParameterisedLossModel
for ion ionisation, and 20 bins per energy decade of physics tables,
(class name G4EmStandardPhysics_option4
)
Models based on Livermore data bases for electrons and gamma,
enabled Rayleigh scattering, enabled fluorescence,
enabled nuclear stopping, enable accurate angular generator
for ionisation models,
G4IonParameterisedLossModel
for ion ionisation, and 20 bins per energy decade of physics tables,
(G4EmLivermorePhysics
);
Models for simulation of linear polarized gamma
based on Livermore data bases for electrons
and gamma (G4EmLivermorePolarizedPhysics
);
Models based on Livermore data bases and new model for
Compton scattering G4LowEPComptonModel
,
new low-energy model of multiple scatetring
G4LowEWenzelMscModel
(G4EmLowEPPhysics
);
Penelope2008 models for electrons, positrons and
gamma,
enabled Rayleigh scattering, enabled fluorescence,
enabled nuclear stopping, enable accurate angular generator
for ionisation models,
G4IonParameterisedLossModel
for ion ionisation, and 20 bins per energy decade of physics tables,
(G4EmPenelopePhysics
);
Low-energy Geant4-DNA physics (G4EmDNAPhysics
).
Alternative low-energy Geant4-DNA physics (G4EmDNAPhysics_option1
). Uses G4LowEWentzelVIModel elastic scattering model for electrons.
Examples of the registration of these physics constructor and
construction of alternative combinations of options are shown
in basic, extended and advanced examples, which can be found in
the subdirectories examples/basic
,
examples/extended/electromagnetic
and
examples/advanced
of the Geant4 source distribution.
Examples illustrating the use
of electromagnetic processes are available as part of the Geant4
release.
Options are available for steering of
electromagnetic processes. These options may be invoked either by
UI commands or by the new C++ interface class G4EmParameters
.
The interface
G4EmParameters::Instance()
is thread safe, EM parameters are shared
between threads, and parameters are shared between all EM processes.
This class has the following public methods:
The old interface class G4EmProcessOptions
is still available and
all old methods are kept unchanged but it is recommended to use only
limited number of complicated methods of this class, which should be applied per process and called
in each thread (both master and workers):
The corresponding UI command can be accessed in the UI subdirectories "/process/eLoss", "/process/em", and "/process/msc". The following types of step limitation by multiple scattering are available:
G4EmCalculator is a class which provides access to cross sections and stopping powers. This class can be used anywhere in the user code provided the physics list has already been initialised (G4State_Idle). G4EmCalculator has "Get" methods which can be applied to materials for which physics tables are already built, and "Compute" methods which can be applied to any material defined in the application or existing in the Geant4 internal database. The public methods of this class are:
For these interfaces, particles, materials, or processes may be pointers or strings with names.
A physical interaction is described by a process class which can handle physics models, described by model classes. The following is a summary of the Low Energy Electromagnetic physics models available in Geant4. Further information is available in the web pages of the Geant4 Low Energy Electromagnetic Physics Working Group, accessible from the Geant4 web site, “who we are” section, then “working groups”.
The physics content of these models is documented in the Geant4 Physics Reference Manual. They are based on the Livermore data library, on the ICRU73 data tables or on the Penelope2008 Monte Carlo code. They adopt the same software design as the "standard" Geant4 electromagnetic models.
Examples of the registration of physics constructor with low-energy
electromagnetic models are shown in Geant4 extended examples
(examples/extended/electromagnetic
in the Geant4 source distribution).
Advanced examples (examples/advanced
in the Geant4 source distribution)
illustrate alternative instantiation of these processes.
Both are available as part of the Geant4 release.
Remember that production cuts for secondaries can be specified as range cuts, which are converted at initialisation time into energy thresholds for secondary gamma, electron, positron and proton production. The cut for proton is applied by elastic scattering processes to aal recoil ions.
A range cut value is set by default to 0.7 mm in Geant4 reference physics lists. This value can be specified in the optional SetCuts() method of the user Physics list or via UI commands. For eg. to set a range cut of 10 micrometers, one can use:
/run/setCut 0.01 mm
or, for a given particle type (for eg. electron),
/run/setCutForAGivenParticle e- 0.01 mm
If a range cut equivalent to an energy lower than 990 eV is specified, the energy cut is still set to 990 eV. In order to decrease this value (for eg. down to 250 eV, in order to simulate low energy emission lines of the fluorescence spectrum), one may use the following UI command before the "/run/initialize" command:
/cuts/setLowEdge 250 eV
or alternatively directly in the user Physics list, in the optional SetCuts() method, using:
G4ProductionCutsTable::GetProductionCutsTable()->SetEnergyRange(250*eV, 1*GeV);
A command is also available in order to disable usage of production threshold for fluorescence and Auger electron production:
/process/em/deexcitationIgnoreCut true
For part of EM processes it is possible to factorise out
sampling of secondary energy and direction. Using an interface
G4VEmModel
base class
SetAngularDistribution(G4VEmAngularDistribution*)
it is possible to substitute default angular generator of a model.
Angular generators in standard and lowenergy sub-packages follow the same
abstract interface.
For photoelectric models several angular generators are available:
For bremsstrahlung models following angular generators are available:
G4DipBustGenerator (default);
G4ModifiedTsai;
G4Generator2BS;
G4Generator2BN;
G4PenelopeBremsstrahlungAngular.
For models of ionisation a new optional angular generator is available:
G4DeltaAngle.
It may be useful to create more than one secondary at an interaction. For example, electrons incident on a target in a medical linac produce photons through bremsstrahlung. The variance reduction technique of bremsstrahlung splitting involves choosing N photons from the expected distribution, and assigning each a weight of 1/N.
Similarly, if the secondaries are not important, one can kill them with a survival probability of 1/N. The weight of the survivors is increased by a factor N. This is known as Russian roulette.
Neither biasing technique is applied if the resulting daughter particles would have a weight below 1/N, in the case of brem splitting, or above 1, in the case of Russian roulette.
These techniques can be enabled in Geant4 electromagnetics with the macro commands
/process/em/setSecBiasing processName Region factor energyLimit energyUnit
where: processName is the name of the process to apply the biasing to; Region is the region in which to apply biasing; factor is the inverse of the brem splitting or Russian roulette factor (1/N); energyLimit energyUnit is the high energy limit. If the first secondary has energy above this limit, no biasing is applied.
For example,
/process/em/setSecBiasing eBrem target 10 100 MeV
will result in electrons undergoing bremsstrahlung in the target region being split 10 times (if the first photon sampled has an energy less than 100 MeV).
Note that the biasing needs to be specified for each process individually. To apply Russian Roulette to daughter electrons from interactions of photons, issue the macro command for the processes phot, compt, conv.
Reference: BEAMnrc Users Manual, D.W.O Rogers, B. Walters, I. Kawrakow. NRCC Report PIRS-0509(A)revL, available at http://www.irs.inms.nrc.ca/inms/irs/BEAM/beamhome.html
Photon models
G4LivermorePhotoElectricModel
)G4LivermorePolarizedPhotoElectricModel
)G4LivermoreComptonModel
)G4LowEPComptonModel
)G4LivermorePolarizedComptonModel
)G4LivermoreRayleighModel
)G4LivermorePolarizedRayleighModel
)G4LivermoreGammaConversionModel
)G4LivermoreNuclearGammaConversionModel
)G4LivermoreGammaConversionModelRC
)G4LivermorePolarizedGammaConversionModel
)
Electron models
G4LivermoreBremsstrahlungModel
)G4LivermoreIonisationModel
)
Ionisation and delta ray production (class G4IonParametrisedLossModel)
The ion model uses ICRU 73 stopping powers, if corresponding ion-material combinations are covered by the ICRU 73 report (up to 1 GeV/nucleon), and otherwise applies a Bethe-Bloch based formalism. For compounds, ICRU 73 stopping powers are employed if the material name coincides with the name of Geant4 NIST materials (e.g. G4_WATER). Elemental materials are matched to the corresponding ICRU 73 stopping powers by means of the atomic number of the material. The material name may be arbitrary in this case. For a list of applicable materials, the user is referred to the ICRU 73 report.
The model requires data files to be copied by the user to his/her code repository. These files are distributed together with the Geant4 release. The user should set the environment variable G4LEDATA to the directory where he/she has copied the files.
The model is dedicated to be used with the G4ionIonisation process and its applicability is restricted to G4GenericIon particles. The ion model is not used by default by this process and must be instantiated and registered by the user:
G4ionIonisation* ionIoni = new G4ionIonisation(); ionIoni -> SetEmModel(new G4IonParametrisedLossModel());
Photon models
G4PenelopeComptonModel
)G4PenelopeRayleighModel
)G4PenelopePhotoElectricModel
)
Electron models
G4PenelopeBremsstrahlungModel
)G4PenelopeIonisationModel
)
Positron models
G4PenelopeBremsstrahlungModel
)G4PenelopeIonisationModel
)
All Penelope models can be applied up to a maximum energy of 100 GeV, although it is advisable not to use them above a few hundreds of MeV.
Options are available in the all Penelope Models, allowing to set (and retrieve) the verbosity level of the model, namely the amount of information which is printed on the screen.
SetVerbosityLevel(G4int)
GetVerbosityLevel()
The default verbosity level is 0 (namely, no textual output on the screen). The default value should be used in general for normal runs. Higher verbosity levels are suggested only for testing and debugging purposes.
The verbosity scale defined for all Penelope processes is the following:
0 = no printout on the screen (default)
1 = issue warnings only in the case of energy non-conservation in the final state (should never happen)
2 = reports full details on the energy budget in the final state
3 = writes also informations on cross section calculation, data file opening and sampling of atoms
4 = issues messages when entering in methods
The Geant4 low energy electromagnetic Physics package has been extended down to energies of a few electronVolts suitable for the simulation of radiation effects in liquid water for applications in micro/nanodosimetry at the cellular and sub-cellular level. These developments take place in the framework of the on-going Geant4-DNA project (see more in the Geant4-DNA web pages or in the web pages of the Geant4 Low Energy Electromagnetic Physics Working Group).
The Geant4-DNA process and model classes apply to electrons, protons, hydrogen, alpha particles and their charge states.
Elastic scattering :
Excitation
Ionisation
Attachment
Vibrational excitation
Excitation
Ionisation
Charge decrease
Excitation
Ionisation
Charge increase
Excitation
Ionisation
Charge increase
Excitation
Ionisation
Charge increase
Charge decrease
Excitation
Ionisation
Charge decrease
Ionisation
An example of the registration of these processes in a physics list is
given in the G4EmDNAPhysics constructor (in source/physics_lists/constructors/electromagnetic
in the Geant4 source distribution).
An example of the usage of this constructor in a physics list is given in the "dnaphysics"
extended example, which explains how to extract basic information from Geant4-DNA Physics processes.
The "microdosimetry" extended example illustrates how to combine Geant4-DNA processes with Standard electromagnetic processes (combination of discrete and condensed history Geant4 electromagnetic processes at different scales).
Since Geant4 release 10.1, Geant4-DNA can also be used for the modelling of water radiolysis (physico-chemistry and chemistry stages). Three extended examples, "chem1", "chem2" and "chem3" illustrate this. More information is available from the Geant4-DNA website.
To run the Geant4-DNA extension, data files need to be copied by the user to his/her code repository. These files are distributed together with the Geant4 release. The user should set the environment variable G4LEDATA to the directory where he/she has copied the files.
A full list of publications regarding Geant4-DNA is directly available from the Geant4-DNA website or from the Geant4@IN2P3 web site).
A unique interface named G4VAtomicDeexcitation is available in Geant4 for the simulation of atomic deexcitation using Standard, Low Energy and Very Low Energy electromagnetic processes. Atomic deexcitation includes fluorescence and Auger electron emission induced by photons, electrons and ions (PIXE); see more details in:
A. Mantero et al., PIXE Simulation in Geant4, X-Ray Spec., 40, 135-140, 2011.
It can be activated for processes producing vacancies in atomic shells. Currently these processes are the photoelectric effect, ionization and Compton scattering.
The activation of atomic deexcitation in continuous processes in a user physics list can be done through the following G4EmProcessOptions class methods, respectively to activate deexcitation processes, Auger effect and deexcitation from ions (PIXE):
G4EmProcessOptions::SetFluo(G4bool); G4EmProcessOptions::SetAuger(G4bool); G4EmProcessOptions::SetPIXE(G4bool);
Please note that in order to have Auger emission it is mandatory to activate Auger electron production for the region in which it is needed (World region included):
G4EmProcessOptions::SetDeexcitationActiveRegion(const G4String& , G4bool, G4bool, G4bool);
The same methods can be accessed interactively via UI commands:
/run/initialize /process/em/deexcitation region true true true /process/em/fluo true /process/em/auger true /process/em/pixe true
Fluorescence from photons and electrons is activated by default in Livermore and Penelope physics constructors, while Auger production and PIXE are not.
The alternative set of data by Bearden et al. (1967) for the modelling of fluorescence lines had been added to the G4LEDATA archive. This set can be selected in any physics list with:
G4AtomicTransitionManager::Instance()->SetFluoDirectory("fluor_Bearden");
Note that in MT mode, reinitialisation of Physics can be done using the following command:
/run/physicsModified
This is done automatically for user applications when atomic de-excitation commands are used.
The user can also select which cross section model to use in order to calculate shell ionisation cross sections for generating PIXE. Again, it is possible to use methods of the G4EmProcessOptions class in the user Physics list:
G4EmProcessOptions::SetPIXECrossSectionModel(const G4String&);
where the string can be "Empirical", "ECPSSR_FormFactor" or "ECPSSR_Analytical", or alternatively with the following UI commands:
/process/em/pixeXSmodel value
where value is equal to Empirical or ECPSSR_FormFactor or ECPSSR_Analytical.
Different shell cross sections models are available : "ECPSSR_Analytical" models derive from an analytical calculation of the ECPSSR theory (see A. Mantero et al., X-Ray Spec.40 (2011) 135-140) and it reproduces K and L shell cross sections over a wide range of energies; "ECPSSR_FormFactor" models derive from A. Taborda et al. calculations (see A. Taborda et al., X-Ray Spec. 40 (2011) 127-134) of ECPSSR values directly form Form Factors and it covers K, L shells on the range 0.1-100 MeV and M shells in the range 0.1-10 MeV; the "empirical" models are from Paul "reference values" (for protons and alphas for K-Shell) and Orlic empirical model for L shells (only for protons and ions with Z>2). The later ones are the models used by default. Out of the energy boundaries, "ECPSSR_Analytical" model is used. We recommend to use default settings if not sure what to use.
The TestEm5 extended/electromagetic example shows how to simulate atomic deexcitation (see for eg. the pixe.mac macro).
(Previously named Geant4-MuElec)
The Geant4 low energy electromagnetic Physics package has been extended down to energies of a few electronVolts suitable for the simulation of radiation effects in highly integrated microelectronic components.
The Geant4-MicroElec process and model classes apply to electrons, protons and heavy ions in silicon.
Elastic scattering :
Ionization
Ionisation
Ionization
A full list of publications regarding Geant4-MicroElec is directly available from the Geant4-MicroElec website.
A new Compton scattering model for unpolarised photons has been developed in the relativistic impulse approximation. The model was developed as an alternative to low energy electromagnetic Compton scattering models developed from Ribberfors' Compton scattering framework (Livermore, Penelope Compton models). The model class is named named G4LowEPComptonModel.
G4LowEPComptonModel has been added to the physics constructor G4EmStandardPhysics_option4, containing the most accurate models from the Standard and Low Energy Electromagnetic physics working groups.
The G4hImpactIonisation process deals with ionisation by impact of hadrons and alpha particles, and the following generation of PIXE (Particle Induced X-ray Emission). This process and related classes can be found in source/processes/electromagnetic/pii .
Further documentation about PIXE simulation with this process is available here.
A detailed description of the related physics features can be found in:
M. G. Pia et al., PIXE Simulation with Geant4, IEEE Trans. Nucl. Sci., vol. 56, no. 6, pp. 3614-3649, 2009.
A brief summary of the related physics features can be found in the Geant4 Physics Reference Manual.
An example of how to use this process is shown below. A more extensive example is available in the eRosita Geant4 advanced example (see examples/advanced/eRosita in your Geant4 installation source).
#include "G4hImpactIonisation.hh" [...] void eRositaPhysicsList::ConstructProcess() { [...] theParticleIterator->reset(); while( (*theParticleIterator)() ) { G4ParticleDefinition* particle = theParticleIterator->value(); G4ProcessManager* processManager = particle->GetProcessManager(); G4String particleName = particle->GetParticleName(); if (particleName == "proton") { // Instantiate the G4hImpactIonisation process G4hImpactIonisation* hIonisation = new G4hImpactIonisation(); // Select the cross section models to be applied for K, L and M shell vacancy creation // (here the ECPSSR model is selected for K, L and M shell; one can mix and match // different models for each shell) hIonisation->SetPixeCrossSectionK("ecpssr"); hIonisation->SetPixeCrossSectionL("ecpssr"); hIonisation->SetPixeCrossSectionM("ecpssr"); // Register the process with the processManager associated with protons processManager -> AddProcess(hIonisation, -1, 2, 2); } } }
The following cross section model options are available:
ecpssr
(based on the ECPSSR theory) ecpssr_hs
(based on the ECPSSR theory, with Hartree-Slater correction) ecpssr_ua
(based on the ECPSSR theory, with United Atom Hartree-Slater correction) ecpssr_he
(based on the ECPSSR theory, with high energy correction) pwba
(plane wave Born approximation) paul
(based on the empirical model by Paul and Sacher) kahoul
(based on the empirical model by Kahoul et al.) ecpssr
ecpssr_ua
pwba
miyagawa
(based on the empirical model by Miyagawa et al.) orlic
(based on the empirical model by Orlic et al.) sow
(based on the empirical model by Sow et al.) ecpssr
pwba
ecpssr
ecpssr_hs
pwba
paul
(based on the empirical model by Paul and Bolik) ecpssr
pwba
ecpssr
pwba
The G4hImpactIonisation process uses a PIXE Data Library.
The PIXE Data Library is distributed in the Geant4 G4PII data set, which must be downloaded along with Geant4 source code.
The G4PIIDATA environment variable must be defined to refer to the location of the G4PII PIXE data library in your filesystem; for instance, if you use a c-like shell:
setenv G4PIIDATA path_to_where_G4PII_has_been_downloaded
Further documentation about the PIXE Data Library is available here.
This section briefly introduces the hadronic physics processes installed in Geant4. For details of the implementation of hadronic interactions available in Geant4, please refer to the Physics Reference Manual.
Each hadronic process object (derived from
G4HadronicProcess
) may have one or more cross section data
sets associated with it. The term "data set" is meant, in a broad
sense, to be an object that encapsulates methods and data for
calculating total cross sections for a given process. The methods
and data may take many forms, from a simple equation using a few
hard-wired numbers to a sophisticated parameterisation using large
data tables. Cross section data sets are derived from the abstract
class G4VCrossSectionDataSet
, and are required to implement
the following methods:
G4bool IsApplicable( const G4DynamicParticle*, const G4Element* )
This method must return True
if the data set is able to
calculate a total cross section for the given particle and
material, and False
otherwise.
G4double GetCrossSection( const G4DynamicParticle*, const G4Element* )
This method, which will be invoked only if True
was
returned by IsApplicable
, must return a cross section, in
Geant4 default units, for the given particle and material.
void BuildPhysicsTable( const G4ParticleDefinition& )
This method may be invoked to request the data set to recalculate its internal database or otherwise reset its state after a change in the cuts or other parameters of the given particle type.
void DumpPhysicsTable( const G4ParticleDefinition& ) = 0
This method may be invoked to request the data set to print its internal database and/or other state information, for the given particle type, to the standard output stream.
Cross section data sets are used by the process for the
calculation of the physical interaction length. A given cross
section data set may only apply to a certain energy range, or may
only be able to calculate cross sections for a particular type of
particle. The class G4CrossSectionDataStore
has been
provided to allow the user to specify, if desired, a series of data
sets for a process, and to arrange the priority of data sets so
that the appropriate one is used for a given energy range,
particle, and material. It implements the following public
methods:
G4CrossSectionDataStore() ~G4CrossSectionDataStore()
and
G4double GetCrossSection( const G4DynamicParticle*, const G4Element* )
For a given particle and material, this method returns a cross
section value provided by one of the collection of cross section
data sets listed in the data store object. If there are no known
data sets, a G4Exception
is thrown and DBL_MIN
is
returned. Otherwise, each data set in the list is queried, in
reverse list order, by invoking its IsApplicable
method
for the given particle and material. The first data set object that
responds positively will then be asked to return a cross section
value via its GetCrossSection
method. If no data set
responds positively, a G4Exception
is thrown and
DBL_MIN
is returned.
void AddDataSet( G4VCrossSectionDataSet* aDataSet )
This method adds the given cross section data set to the end of the
list of data sets in the data store. For the evaluation of cross
sections, the list has a LIFO (Last In First Out) priority, meaning
that data sets added later to the list will have priority over
those added earlier to the list. Another way of saying this, is
that the data store, when given a GetCrossSection
request,
does the IsApplicable
queries in the reverse list order,
starting with the last data set in the list and proceeding to the
first, and the first data set that responds positively is used to
calculate the cross section.
void BuildPhysicsTable( const G4ParticleDefinition& aParticleType )
This method may be invoked to indicate to the data store that there
has been a change in the cuts or other parameters of the given
particle type. In response, the data store will invoke the
BuildPhysicsTable
of each of its data sets.
void DumpPhysicsTable( const G4ParticleDefinition& )
This method may be used to request the data store to invoke the
DumpPhysicsTable
method of each of its data sets.
The defaults for total cross section data and calculations have
been encapsulated in the singleton class
G4HadronCrossSections
. Each hadronic process:
G4HadronInelasticProcess
,
G4HadronElasticProcess
,
G4HadronFissionProcess
,
and G4HadronCaptureProcess
,
comes already equipped with a cross section data store and a
default cross section data set. The data set objects are really
just shells that invoke the singleton G4HadronCrossSections
to do the real work of calculating cross sections.
The default cross sections can be overridden in whole or in part
by the user. To this end, the base class G4HadronicProcess
has a ``get'' method:
G4CrossSectionDataStore* GetCrossSectionDataStore()
which gives public access to the data store for each process. The user's cross section data sets can be added to the data store according to the following framework:
G4Hadron...Process aProcess(...) MyCrossSectionDataSet myDataSet(...) aProcess.GetCrossSectionDataStore()->AddDataSet( &MyDataSet )
The added data set will override the default cross section data
whenever so indicated by its IsApplicable
method.
In addition to the ``get'' method, G4HadronicProcess
also
has the method
void SetCrossSectionDataStore( G4CrossSectionDataStore* )
which allows the user to completely replace the default data store with a new data store.
It should be noted that a process does not send any information
about itself to its associated data store (and hence data set)
objects. Thus, each data set is assumed to be formulated to
calculate cross sections for one and only one type of process. Of
course, this does not prevent different data sets from sharing
common data and/or calculation methods, as in the case of the
G4HadronCrossSections
class mentioned above. Indeed,
G4VCrossSectionDataSet
specifies only the abstract interface
between physics processes and their data sets, and leaves the user
free to implement whatever sort of underlying structure is
appropriate.
The current implementation of the data set
G4HadronCrossSections
reuses the total cross-sections for
inelastic and elastic scattering, radiative capture and fission as
used with GHEISHA to provide cross-sections
for calculation
of the respective mean free paths of a given particle in a given
material.
The cross section data for low energy neutron transport are
organized in a set of files that are read in by the corresponding
data set classes at time zero. Hereby the file system is used, in
order to allow highly granular access to the data. The ``root''
directory of the cross-section directory structure is accessed
through an environment variable, G4NEUTRONHPDATA
,
which is to be set by the user. The classes accessing the total
cross-sections of the individual processes, i.e., the cross-section
data set classes for low energy neutron transport, are
G4NeutronHPElasticData
,
G4NeutronHPCaptureData
,
G4NeutronHPFissionData
,
and G4NeutronHPInelasticData
.
For detailed descriptions of the low energy neutron total cross-sections, they may be registered by the user as described above with the data stores of the corresponding processes for neutron interactions.
It should be noted that using these total cross section classes does not require that the neutron_hp models also be used. It is up to the user to decide whethee this is desirable or not for his particular problem.
A prototype of the compact version of neutron cross sections derived from HP database
are provided with new classes G4NeutronHPElasticData
,
G4NeutronCaptureXS
,
G4NeutronElasticXS
,
and G4NeutronInelasticXS
.
The following process classes have been implemented:
pi- absorption (class name G4PionMinusAbsorptionAtRest
or G4PiMinusAbsorptionAtRest
)
kaon- absorption (class name G4KaonMinusAbsorptionAtRest
or G4KaonMinusAbsorption
)
neutron capture (class name G4NeutronCaptureAtRest
)
anti-proton annihilation (class name
G4AntiProtonAnnihilationAtRest
)
anti-neutron annihilation (class name
G4AntiNeutronAnnihilationAtRest
)
mu- capture (class name G4MuonMinusCaptureAtRest
)
alternative CHIPS model for any negativly charged particle
(class name G4QCaptureAtRest
)
Obviously the last process does not, strictly speaking, deal with a ``hadron at rest''. It does, nonetheless, share common features with the others in the above list because of the implementation model chosen. The differences between the alternative implementation for kaon and pion absorption concern the fast part of the emitted particle spectrum. G4PiMinusAbsorptionAtRest, and G4KaonMinusAbsorptionAtRest focus especially on a good description of this part of the spectrum.
All of these classes are derived from the abstract class
G4VRestProcess
. In addition to the constructor and
destructor methods, the following public methods of the abstract
class have been implemented for each of the above six
processes:
AtRestGetPhysicalInteractionLength( const G4Track&,
G4ForceCondition* )
This method returns the time taken before the interaction actually occurs. In all processes listed above, except for muon capture, a value of zero is returned. For the muon capture process the muon capture lifetime is returned.
AtRestDoIt( const G4Track&, const G4Step&)
This method generates the secondary particles produced by the process.
IsApplicable( const G4ParticleDefinition& )
This method returns the result of a check to see if the process is possible for a given particle.
Including a ``hadron at rest'' process for a particle, a pi- for example, into the Geant4 system is straightforward and can be done in the following way:
create a process:
theProcess = new G4PionMinusAbsorptionAtRest();
register the process with the particle's process manager:
theParticleDef = G4PionMinus::PionMinus(); G4ProcessManager* pman = theParticleDef->GetProcessManager(); pman->AddRestProcess( theProcess );
For hadrons in motion, there are four physics process classes. Table 5.1 shows each process and the particles for which it is relevant.
G4HadronElasticProcess
| pi+, pi-, K+, K0S, K0L, K-, p, p-bar, n, n-bar, lambda, lambda-bar, Sigma+, Sigma-, Sigma+-bar, Sigma--bar, Xi0, Xi-, Xi0-bar, Xi--bar |
G4HadronInelasticProcess
| pi+, pi-, K+, K0S, K0L, K-, p, p-bar, n, n-bar, lambda, lambda-bar, Sigma+, Sigma-, Sigma+-bar, Sigma--bar, Xi0, Xi-, Xi0-bar, Xi--bar |
G4HadronFissionProcess
| all |
G4CaptureProcess
| n, n-bar |
Table 5.1. Hadronic processes and relevant particles.
To register an inelastic process model for a particle, a proton for example, first get the pointer to the particle's process manager:
G4ParticleDefinition *theProton = G4Proton::ProtonDefinition(); G4ProcessManager *theProtonProcMan = theProton->GetProcessManager();
Create an instance of the particle's inelastic process:
G4ProtonInelasticProcess *theProtonIEProc = new G4ProtonInelasticProcess();
Create an instance of the model which determines the secondaries produced in the interaction, and calculates the momenta of the particles:
G4LEProtonInelastic *theProtonIE = new G4LEProtonInelastic();
Register the model with the particle's inelastic process:
theProtonIEProc->RegisterMe( theProtonIE );
Finally, add the particle's inelastic process to the list of discrete processes:
theProtonProcMan->AddDiscreteProcess( theProtonIEProc );
The particle's inelastic process class,
G4ProtonInelasticProcess
in the example above, derives from
the G4HadronicInelasticProcess
class, and simply defines the
process name and calls the G4HadronicInelasticProcess
constructor. All of the specific particle inelastic processes
derive from the G4HadronicInelasticProcess
class, which
calls the PostStepDoIt
function, which returns the
particle change object from the G4HadronicProcess
function
GeneralPostStepDoIt
. This class also gets the mean free
path, builds the physics table, and gets the microscopic cross
section. The G4HadronicInelasticProcess
class derives from
the G4HadronicProcess
class, which is the top level hadronic
process class. The G4HadronicProcess
class derives from the
G4VDiscreteProcess
class. The inelastic, elastic, capture,
and fission processes derive from the G4HadronicProcess
class. This pure virtual class also provides the energy range
manager object and the RegisterMe
access function.
A sample case for the proton's inelastic interaction model class
is shown in Example 5.1, where
G4LEProtonInelastic.hh
is the name of the include
file:
Example 5.1. An example of a proton inelastic interaction model class.
----------------------------- include file ------------------------------------------ #include "G4InelasticInteraction.hh" class G4LEProtonInelastic : public G4InelasticInteraction { public: G4LEProtonInelastic() : G4InelasticInteraction() { SetMinEnergy( 0.0 ); SetMaxEnergy( 25.*GeV ); } ~G4LEProtonInelastic() { } G4ParticleChange *ApplyYourself( const G4Track &aTrack, G4Nucleus &targetNucleus ); private: void CascadeAndCalculateMomenta( required arguments ); }; ----------------------------- source file ------------------------------------------ #include "G4LEProtonInelastic.hh" G4ParticleChange * G4LEProton Inelastic::ApplyYourself( const G4Track &aTrack, G4Nucleus &targetNucleus ) { theParticleChange.Initialize( aTrack ); const G4DynamicParticle *incidentParticle = aTrack.GetDynamicParticle(); // create the target particle G4DynamicParticle *targetParticle = targetNucleus.ReturnTargetParticle(); CascadeAndCalculateMomenta( required arguments ) { ... } return &theParticleChange; }
The CascadeAndCalculateMomenta
function is the bulk of
the model and is to be provided by the model's creator. It should
determine what secondary particles are produced in the interaction,
calculate the momenta for all the particles, and put this
information into the ParticleChange object which is
returned.
The G4LEProtonInelastic
class derives from the
G4InelasticInteraction
class, which is an abstract base
class since the pure virtual function ApplyYourself
is not
defined there. G4InelasticInteraction
itself derives from
the G4HadronicInteraction
abstract base class. This class is
the base class for all the model classes. It sorts out the energy
range for the models and provides class utilities. The
G4HadronicInteraction
class provides the
Set/GetMinEnergy
and the Set/GetMaxEnergy
functions which determine the minimum and maximum energy range for
the model. An energy range can be set for a specific element, a
specific material, or for general applicability:
void SetMinEnergy( G4double anEnergy, G4Element *anElement ) void SetMinEnergy( G4double anEnergy, G4Material *aMaterial ) void SetMinEnergy( const G4double anEnergy ) void SetMaxEnergy( G4double anEnergy, G4Element *anElement ) void SetMaxEnergy( G4double anEnergy, G4Material *aMaterial ) void SetMaxEnergy( const G4double anEnergy )
In Geant4, any model can be run together with any other model without the need for the implementation of a special interface, or batch suite, and the ranges of applicability for the different models can be steered at initialisation time. This way, highly specialised models (valid only for one material and particle, and applicable only in a very restricted energy range) can be used in the same application, together with more general code, in a coherent fashion.
Each model has an intrinsic range of applicability, and the model chosen for a simulation depends very much on the use-case. Consequently, there are no ``defaults''. However, physics lists are provided which specify sets of models for various purposes.
Three types of hadronic shower models have been implemented: parametrisation driven models, data driven models, and theory driven models.
Parametrisation driven models are used for all processes
pertaining to particles coming to rest, and interacting with the
nucleus. For particles in flight, two sets of models exist for
inelastic scattering; low energy, and high energy models. Both sets
are based originally on the GHEISHA
package of Geant3.21,
and the original approaches to primary interaction, nuclear
excitation, intra-nuclear cascade and evaporation is kept. The
models are located in the sub-directories
hadronics/models/low_energy
and
hadronics/models/high_energy
. The low energy models are
targeted towards energies below 20 GeV; the high energy models
cover the energy range from 20 GeV to O(TeV). Fission, capture and
coherent elastic scattering are also modeled through parametrised
models.
Data driven models are available for the transport of low
energy neutrons in matter in sub-directory
hadronics/models/neutron_hp
. The modeling is based
on the data formats of ENDF/B-VI,
and all distributions of this standard data format are implemented.
The data sets used are selected from data libraries that conform to
these standard formats. The file system is used in order to allow granular
access to, and flexibility in, the use of the cross sections for different
isotopes, and channels. The energy coverage of these models is from
thermal energies to 20 MeV.
Theory driven models are available for inelastic scattering in
a first implementation, covering the full energy range of LHC
experiments. They are located in sub-directory
hadronics/models/generator
. The current philosophy
implies the usage of parton string models at high energies, of
intra-nuclear transport models at intermediate energies, and of
statistical break-up models for de-excitation.
Nuclear models fail (sometimes catastrophically) at predicting with reasonable accuracies the nuclear cross sections of neutrons (and other particles). For this reason, all physical quantities relevant for an accurate modeling of nuclear reactions in Monte Carlo simulations need to be provided as a database which includes, ideally:
For the case of neutron induced reactions, such databases are called “evaluated data”, in the sense that they contain recommended values for different quantities that rely on compilations of experimental nuclear data and usually completed with theoretical predictions, benchmarked against available experimental data (i.e. integral and differential experiments) when possible. It should be noticed that the information available varies from isotope to isotope and can be incomplete or totally missing.
The G4NeutronHP package in GEANT4 allows using evaluated nuclear data libraries in the G4NDL format. GEANT4 users should know that any simulation involving neutrons with energies below 20 MeV and not using the G4NeutronHP package can lead to unreliable results. GEANT4 users are therefore encouraged to use it, although they should be aware of the limitations of using evaluated nuclear data libraries.
An example about how to implement the G4NeutronHP package into physics
list in a GEANT4 application can be found in the example case (among others
distributed with
GEANT4) extended/radioactivedecay/rdecay02
. Three
different processes are included in that example: elastic, capture and
inelastic. The inelastic reactions in G4NeutronHP are all reactions except
elastic, capture and fission, so fission should also be included in the
physics list, if needed, and it is done in the same way as it is done for
the other three.
The G4NeutronHP package must be used together with evaluated nuclear data libraries. They are distributed by the GEANT4 collaboration (http://geant4.web.cern.ch/geant4/support/download.shtml) and from the IAEA nuclear data web site (http://www-nds.iaea.org/geant4/) where a larger set of different libraries, including isotopes with Z > 92, is available.
The evaluated nuclear data libraries do differ and thus the results of the Monte Carlo simulations will depend on the library used. It is a safe practice to perform simulations with (at least) two different libraries for estimating the uncertainties associated to the nuclear data.
Together with a good implementation of the physics list, users must be very careful with the definition of the materials performed in a Monte Carlo simulation when low energy neutron transport is relevant. In contrast to other kind of simulations, the isotopic composition of the elements which compose the different materials can strongly affect the obtained simulation results. Because of this, it is strongly recommended to define specifically the isotopic composition of each element used in the simulation, as it is described in the GEANT4 user’s manual. In principle, such a practice is not mandatory if natural isotopic compositions are used, since GEANT4 contains them in their databases. However, by defining them explicitly some unexpected problems may be avoided and a better control of the simulation will be achieved.
It is highly recommended or mandatory to set the following UNIX environment variables running a GEANT4 application:
G4NEUTRONHPDATA
G4NEUTRONHP_SKIP_MISSING_ISOTOPES=1
It sets to zero the cross section of the isotopes which are not present in the neutron library. If GEANT4 doesn’t find an isotope, then it looks for the natural composition data of that element. Only if the element is not found then the cross section is set to zero. On the contrary, if this variable is not defined, GEANT4 looks then for the neutron data of another isotope close in Z and A, which will have completely different nuclear properties and lead to incorrect results (highly recommended).
G4NEUTRONHP_DO_NOT_ADJUST_FINAL_STATE=1
If this variable is not defined, a GEANT4 model that attempts to satisfy the energy and momentum conservation in some nuclear reactions, by generating artificial gamma rays. By setting such a variable one avoids the correction and leads to the result obtained with the ENDF-6 libraries. Even though energy and momentum conservation are desirable, the ENDF-6 libraries do not provide the necessary correlations between secondary particles for satisfying them in all cases. On the contrary, ENDF-6 libraries intrinsically violate energy and momentum conservation for several processes and have been built for preserving the overall average quantities such as average energy releases, average number of secondaries… (highly recommended).
AllowForHeavyElements=1
Activates the physics for isotopes with Z>92 (recommended).
The G4NDL format libraries are based on the ENDF-6 format libraries, which contain evaluated (i.e. recommended) nuclear data prepared for their use in transport codes. These data are essentially nuclear reaction cross sections together with the distribution in energy and angle of the secondary reaction products. As a consequence of how the data is written in the ENDF files, there are some features that may be or may be not expected in the results of a Monte Carlo calculation.
The information concerning the creation of the reaction products can be incomplete and/or uncorrelated, in the sense that is described below:
This applies when there is no information about how to generate a secondary particle. As an example, it is possible to have only the cross section data of an (n,p) reaction, without any information concerning the energy and angle of the secondary proton. In this case GEANT4 will produce the proton considering that it is emitted isotropically in the center of mass frame, with an energy which is deduced from assuming that the residual nucleus is in its ground state.
This applies when:
The energy and angle distributions of a reaction product may be uncorrelated. As a consequence, the reaction products can be generated with an unphysical energy-angle relationship.
The energy-angle distributions of different reaction products of a certain reaction are always uncorrelated. As an example, consider that in a (n, 2p) reaction at a certain neutron energy both resulting protons can be emitted with energies ranging from 0 to 5MeV. In this case the energy and angle of each proton will be sampled independently of the energy and angle of the other proton, so there will be events in which both protons will be emitted with energies close to 5 MeV and there will also be events in which both protons will be emitted with energies close to 0 MeV. As a consequence, energy and angular momentum won’t be conserved event by event. However, energy will be conserved in average and the resulting proton energy spectrum will be correctly produced.
There are some cases where several nuclear reactions are put
together as if they were a single reaction (MT=5
reaction, in ENDF-6 format nomenclature). In those cases the information
consists in a cross section, which is the sum of all of them, plus a
reaction product yield and energy-angle distributions for each secondary
particle. In this case the amount of each secondary particle produced has to
be sampled every time the reaction occurs, and it is done independently of
the amount of the other secondary particles produced.
Thus, in this case neither the energy and angular momentum nor the number of nucleons is conserved event by event, but all the quantities should be conserved in average. As a consequence, it is also not possible to deduce which are the residual nuclei produced, since no information is available concerning what are the specific nuclear reactions which take place. It has to be said that sometimes ENDF libraries include the residual nuclei as an outgoing particle. However, GEANT4 does not manage that information, at present. This situation is quite uncommon in neutron data libraries up to 20 MeV. However, it is quite common to find it in charged particle libraries below 20 MeV or in neutron libraries above 20 MeV.
As a consequence of what has been presented above, some general features can be expected in the results of a Monte Carlo calculation performed with the G4NeutronHP package:
The neutron transport, which means how the neutron looses energy in the collisions, when and how it is absorbed…, is quite trustable, since the main purpose of the ENDF neutron libraries is to perform this neutron transport.
The production of neutrons due to neutron induced nuclear reactions is usually trustable, with the exception of the energy-angle correlations when several neutrons are produced in the same nuclear reaction.
The results concerning the production of charged particles have to be always questioned. A look into the ENDF format library used can indicate which results are trustable and which are not. This can be done, for example, in http://t2.lanl.gov/data/data.html, among other websites.
The results concerning the production of γ-rays have to be questioned always. For example, the information on the number and energies of γ-rays emitted in the neutron capture process is incomplete for almost all the nuclei and is frequently also uncorrelated. When the information is available, it will be used, but one can obtain results which are quite far from reality on an event by event basis: the total energy of the cascade won’t be correct in many cases and only some specific γ-rays which are stored in the neutron databases will be emitted. If there isn’t any information concerning these γ-rays, GEANT4 will use a simple a model instead which is generally missing the relevant spectroscopic information. The results concerning the generation of residual nuclei (for example, in activation calculations) are usually trustable, with the exception of libraries with MT=5 reactions, as described above (2).
As a general conclusion, users should always be critical with the results obtained with Monte Carlo simulation codes, and this also applies to GEANT4. They have to anticipate which results can be trusted and which results should be questioned. For the particular case of the a closer look into the underlying evaluated nuclear datain the ENDF format libraries will allow to check what is the information available in a certain library for some specific isotope and a certain reaction. There are several public nuclear data sites like http://t2.lanl.gov/data/data.html.
The transport of very low energy neutrons (below 5 eV) has to be performed using the thermal neutron data libraries. At these energies, the fact that the nuclei are in atoms which form part of a certain molecule inside a material (crystal lattice, liquid, plastic…) plays an important role, since there can be a transference of momentum between the neutron and the whole structure of the material, not only with the nucleus. This is of particular importance for material used as neutron moderators, i.e., materials with low A (mass number) used to decrease the incident neutron energy in only a few collisions. Since the property is related to the nucleus in the material, as an example, there is the need for having different thermal libraries for Hydrogen in polyethylene, Hydrogen in water and so on.
If neutron collisions at these energies are relevant for the problem to be simulated, thermal libraries should be used for the materials if they are available. If they are not, the results obtained from the simulation will not be trustable in the neutron energy range below 5 eV, especially when using low mass elements in the simulation.
To use the thermal libraries the following lines should be included in the physics list:
G4HadronElasticProcess* theNeutronElasticProcess = new G4HadronElasticProcess; // Cross Section Data set G4NeutronHPElasticData* theHPElasticData = new G4NeutronHPElasticData; theNeutronElasticProcess->AddDataSet(theHPElasticData); G4NeutronHPThermalScatteringData* theHPThermalScatteringData = new G4NeutronHPThermalScatteringData; theNeutronElasticProcess->AddDataSet(theHPThermalScatteringData); // Models G4NeutronHPElastic* theNeutronElasticModel = new G4NeutronHPElastic; theNeutronElasticModel->SetMinEnergy(4.0*eV); theNeutronElasticProcess->RegisterMe(theNeutronElasticModel); G4NeutronHPThermalScattering* theNeutronThermalElasticModel = new G4NeutronHPThermalScattering; theNeutronThermalElasticModel->SetMaxEnergy(4.0*eV); theNeutronElasticProcess->RegisterMe(theNeutronThermalElasticModel); // Apply Processes to Process Manager of Neutron G4ProcessManager* pmanager = G4Neutron::Neutron()->GetProcessManager(); pmanager->AddDiscreteProcess(theNeutronElasticProcess);
And the materials should be defined with a specific name. For example, to use the thermal library for Hydrogen in water, the water should be defined as:
G4Element* elTSHW = new G4Element("TS_H_of_Water", "H_WATER", 1.0, 1.0079*g/mole); G4Material* matH2O_TS = new G4Material("Water_TS", density=1.0*g/cm3, ncomponents=2); matH2O_TS->AddElement(elTSHW,natoms=2); matH2O_TS->AddElement(elO,natoms=1);
where the important thing is the
name "TS_H_of_Water", which is a specific
name used by G4NeutronHP. In order to see which thermal libraries are
available, they can be found in the G4NDL4.0/ThermalScattering folder (or
equivalent, for other neutron libraries). Then, one has to look into the
G4NeutronHPThermalScatteringNames.cc
source file, under
source/processes/hadronic/models/neutron_hp/src
. There
are some lines similar to:
names.insert(std::pair<G4String,G4String>("TS_H_of_Water", "h_water"));
where "TS_H_of_Water" means Hydrogen in water. Names similar to "TS_H_of_Water" like "TS_C_of_Graphite" or "TS_H_of_Polyethylene" can be found and used in the same way as described above.
Nuclear reactions at intermediate energies (from a few MeV to a few GeV) are typically modelled in two stages. The first, fast reaction stage is described by a dynamical model (quantum molecular dynamics, intranuclear cascade, pre-compound, etc.) and often results in the production of one or several excited nuclei. The second reaction stage describes the de-excitation of the excited nuclei and it is usually handled by statistical de-excitation models. The models for the two reaction stages can in principle be chosen independently, but the current design of the Geant4 hadronics framework makes it difficult to do this at the physics-list level. However, another solution exists.
Geant4 provides several nuclear de-excitation models. The default one is
G4ExcitationHandler
, which is described in
detail in the
Physics Reference Manual. The
Bertini-style G4CascadeInterface
uses an internal
de-excitation model. The ABLA V3 model is also
available.
It is possible to replace the default de-excitation model with ABLA
V3 for any intranuclear-cascade model in Geant4 except
G4CascadeInterface
. The easiest way to do this is to
call the SetDeExcitation()
method of the relevant
intranuclear-cascade-model interface. This can be done even if you are using
one of the reference physics lists. The technique is the following.
For clarity's sake, assume you are using the FTFP_INCLXX
physics list, which uses INCL++, the Liege Intranuclear
Cascade model (G4INCLXXInterface
) at intermediate
energies. You can couple INCL++ to ABLA
V3 by adding a run action (Section 6.2.1) and adding the following code snippet to
BeginOfRunAction()
.
Example 5.2. Coupling the INCL++ model to ABLA V3
#include "G4HadronicInteraction.hh" #include "G4HadronicInteractionRegistry.hh" #include "G4INCLXXInterface.hh" #include "G4AblaInterface.hh" void MyRunAction::BeginOfRunAction(const G4Run*) { // Get hold of pointers to the INCL++ model interfaces std::vector<G4HadronicInteraction *> interactions = G4HadronicInteractionRegistry::Instance() ->FindAllModels(G4INCLXXInterfaceStore::GetInstance()->getINCLXXVersionName()); for(std::vector<G4HadronicInteraction *>::const_iterator iInter=interactions.begin(), e=interactions.end(); iInter!=e; ++iInter) { G4INCLXXInterface *theINCLInterface = static_cast<G4INCLXXInterface*>(*iInter); if(theINCLInterface) { // Instantiate the ABLA model G4HadronicInteraction *interaction = G4HadronicInteractionRegistry::Instance()->FindModel("ABLA"); G4AblaInterface *theAblaInterface = static_cast<G4AblaInterface*>(interaction); if(!theAblaInterface) theAblaInterface = new G4AblaInterface; // Couple INCL++ to ABLA G4cout << "Coupling INCLXX to ABLA" << G4endl; theINCLInterface->SetDeExcitation(theAblaInterface); } } }
This technique may be applied to any intranuclear-cascade model (i.e. models
that inherit from G4VIntraNuclearTransportModel
), except
G4CascadeInterface
. For example, if your physics list
relies on the Binary-Cascade model (e.g.
FTF_BIC), you'll need to do
// Get hold of a pointer to the Binary-Cascade model interface std::vector<G4HadronicInteraction *> interactions = G4HadronicInteractionRegistry::Instance() ->FindAllModels("Binary Cascade"); for(std::vector<G4HadronicInteraction *>::const_iterator iInter=interactions.begin(), e=interactions.end(); iInter!=e; ++iInter) { G4BinaryCascade *theBICInterface = static_cast<G4BinaryCascade*>(*iInter); if(theBICInterface) { // Instantiate ABLA V3 as in the example above // [...] // Couple BIC to ABLA theBICInterface->SetDeExcitation(theAblaInterface); } }
This section briefly introduces decay processes installed in Geant4. For details of the implementation of particle decays, please refer to the Physics Reference Manual.
Geant4 provides a G4Decay
class for both ``at rest'' and
``in flight'' particle decays. G4Decay
can be applied to all
particles except:
G4ParticleDefinition::thePDGMass <= 0
G4ParticleDefinition::thePDGLifeTime < 0
G4ParticleDefinition::fShortLivedFlag = True
Decay for some particles may be switched on or off by using
G4ParticleDefinition::SetPDGStable()
as well as
ActivateProcess()
and InActivateProcess()
methods of G4ProcessManager
.
G4Decay
proposes the step length (or step time for
AtRest
) according to the lifetime of the particle unless
PreAssignedDecayProperTime
is defined in
G4DynamicParticle
.
The G4Decay
class itself does not define decay modes of
the particle. Geant4 provides two ways of doing this:
G4DecayChannel
in G4DecayTable
,
andthePreAssignedDecayProducts
of
G4DynamicParticle
The G4Decay
class calculates the
PhysicalInteractionLength
and boosts decay products
created by G4VDecayChannel
or event generators. See below
for information on the determination of the decay modes.
An object of G4Decay
can be shared by particles.
Registration of the decay process to particles in the
ConstructPhysics
method of PhysicsList
(see Section 2.5.3)
is shown in Example 5.3.
Example 5.3.
Registration of the decay process to particles in the
ConstructPhysics
method of PhysicsList.
#include "G4Decay.hh" void MyPhysicsList::ConstructGeneral() { // Add Decay Process G4Decay* theDecayProcess = new G4Decay(); theParticleIterator->reset(); while( (*theParticleIterator)() ){ G4ParticleDefinition* particle = theParticleIterator->value(); G4ProcessManager* pmanager = particle->GetProcessManager(); if (theDecayProcess->IsApplicable(*particle)) { pmanager ->AddProcess(theDecayProcess); // set ordering for PostStepDoIt and AtRestDoIt pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep); pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest); } } }
Each particle has its G4DecayTable
, which stores information
on the decay modes of the particle. Each decay mode, with its
branching ratio, corresponds to an object of various ``decay
channel'' classes derived from G4VDecayChannel
. Default
decay modes are created in the constructors of particle classes.
For example, the decay table of the neutral pion has
G4PhaseSpaceDecayChannel
and
G4DalitzDecayChannel
as follows:
// create a decay channel G4VDecayChannel* mode; // pi0 -> gamma + gamma mode = new G4PhaseSpaceDecayChannel("pi0",0.988,2,"gamma","gamma"); table->Insert(mode); // pi0 -> gamma + e+ + e- mode = new G4DalitzDecayChannel("pi0",0.012,"e-","e+"); table->Insert(mode);
Decay modes and branching ratios defined in Geant4 are listed in Section 5.3.2.
Branching ratios and life time can be set in tracking time.
// set lifetime G4Neutron::Neutron()->SetPDGLifeTime(885.7*second); // allow neutron decay G4Neutron::Neutron()->SetPDGStable(false);
Branching ratios and life time can be modified by using user commands, also.
Example: Set 100% br for dalitz decay of pi0
Idle> /particle/select pi0 Idle> /particle/property/decay/select 0 Idle> /particle/property/decay/br 0 Idle> /particle/property/decay/select 1 Idle> /particle/property/decay/br 1 Idle> /particle/property/decay/dump G4DecayTable: pi0 0: BR: 0 [Phase Space] : gamma gamma 1: BR: 1 [Dalitz Decay] : gamma e- e+
Decays of heavy flavor particles such as B mesons are very complex,
with many varieties of decay modes and decay mechanisms. There are
many models for heavy particle decay provided by various event
generators and it is impossible to define all the decay modes of
heavy particles by using G4VDecayChannel
. In other words,
decays of heavy particles cannot be defined by the Geant4 decay
process, but should be defined by event generators or other
external packages. Geant4 provides two ways to do this:
pre-assigned decay mode
and external decayer
.
In the latter approach, the class G4VExtDecayer
is used
for the interface to an external package which defines decay modes
for a particle. If an instance of G4VExtDecayer
is attached
to G4Decay
, daughter particles will be generated by the
external decay handler.
In the former case, decays of heavy particles are simulated by
an event generator and the primary event contains the decay
information. G4VPrimaryGenerator
automatically attaches any
daughter particles to the parent particle as the
PreAssignedDecayProducts member of G4DynamicParticle
.
G4Decay
adopts these pre-assigned daughter particles instead
of asking G4VDecayChannel
to generate decay products.
In addition, the user may assign a pre-assigned
decay
time for a specific track in its rest frame (i.e. decay time is
defined in the proper time) by using the
G4PrimaryParticle::SetProperTime()
method.
G4VPrimaryGenerator
sets the PreAssignedDecayProperTime
member of G4DynamicParticle
. G4Decay
uses this decay time instead of the life time of the particle type.
Gamma-nuclear and lepto-nuclear reactions are handled in Geant4 as hybrid processes which typically require both electromagnetic and hadronic models for their implementation. While neutrino-induced reactions are not currently provided, the Geant4 hadronic framework is general enough to include their future implementation as a hybrid of weak and hadronic models.
The general scheme followed is to factor the full interaction into an electromagnetic (or weak) vertex, in which a virtual particle is generated, and a hadronic vertex in which the virtual particle interacts with a target nucleus. In most cases the hadronic vertex is implemented by an existing Geant4 model which handles the intra-nuclear propagation.
The cross sections for these processes are parameterizations, either directly of data or of theoretical distributions determined from the integration of lepton-nucleon cross sections double differential in energy loss and momentum transfer.
For the most part gammas can be treated as hadrons and in fact they interact
that way with the nucleus when the Bertini-style cascade
G4CascadeInterface
and QGSP models
are used. These models may be assigned to
G4PhotoNuclearProcess
as shown in the following partial
code:
G4TheoFSGenerator* theHEModel = new G4TheoFSGenerator; G4QGSModel* theStringModel = new G4QGSModel<G4GammaParticipants>; G4ExcitedStringDecay* theStringDecay = new G4ExcitedStringDecay(theFragmentation=new G4QGSMFragmentation); theStringModel->SetFragmentationModel(theStringDecay); theHEModel->SetHighEnergyGenerator(theStringModel); theHEModel->SetTransport(new G4GeneratorPrecompoundInterface); theHEModel->SetMinEnergy(8*GeV); G4CascadeInterface* theLEModel = new G4CascadeInterface; theLEModel->SetMaxEnergy(10*GeV); G4PhotoNuclearProcess* thePhotoNuclearProcess = new G4PhotoNuclearProcess; thePhotoNuclearProcess->RegisterMe(theLEModel); thePhotoNuclearProcess->RegisterMe(theHEModel); G4ProcessManager* procMan = G4Gamma::Gamma()->GetProcessManager(); procMan->AddDiscreteProcess(thePhotoNuclearProcess);
Electro-nuclear reactions in Geant4 are handled by the classes
G4ElectronNuclearProcess
and
G4PositronNuclearProcess
, which are both implmented by
G4ElectroVDNuclearModel
. This model consists of three
sub-models: code which generates the virtual photon from the lepton-nucleus
vertex, the Bertini-style cascade to handle the low and medium energy photons,
and the FTFP model to handle the high energy photons.
Muon-nuclear reactions are handled similarly. The process
G4MuonNuclearProcess
can be assigned the
G4MuonVDNuclearModel
which in turn is implemented by
three sub-models: virtual gamma generation code, Bertini-style cascade and
the FTFP model.
A photon is considered to be optical when its wavelength is much greater than the typical atomic spacing. In GEANT4 optical photons are treated as a class of particle distinct from their higher energy gamma cousins. This implementation allows the wave-like properties of electromagnetic radiation to be incorporated into the optical photon process. Because this theoretical description breaks down at higher energies, there is no smooth transition as a function of energy between the optical photon and gamma particle classes.
For the simulation of optical photons to work correctly in GEANT4, they must be imputed a linear polarization. This is unlike most other particles in GEANT4 but is automatically and correctly done for optical photons that are generated as secondaries by existing processes in GEANT4. Not so, if the user wishes to start optical photons as primary particles. In this case, the user must set the linear polarization using particle gun methods, the General Particle Source, or his/her PrimaryGeneratorAction. For an unpolarized source, the linear polarization should be sampled randomly for each new primary photon.
The GEANT4 catalogue of processes at optical wavelengths includes refraction and reflection at medium boundaries, bulk absorption, Mie and Rayleigh scattering. Processes which produce optical photons include the Cerenkov effect and scintillation. Optical photons are generated in GEANT4 without energy conservation and their energy must therefore not be tallied as part of the energy balance of an event.
The optical properties of the medium which are key to the
implementation of these types of processes are stored as entries in
a G4MaterialPropertiesTable
which is linked to the
G4Material
in question. These properties may be constants
or they may be expressed as a function of the photon's energy.
This table is a private data member of the G4Material
class. The G4MaterialPropertiesTable
is implemented as a
hash directory, in which each entry consists of a value and
a key. The key is used to quickly and efficiently retrieve
the corresponding value. All values in the dictionary are either
instantiations of G4double
or the class
G4MaterialPropertyVector
, and all keys are of type
G4String
.
A G4MaterialPropertyVector
is a typedef of
G4PhysicsOrderedFreeVector. The entries are a pair of numbers, which in
the case of an optical property, are the photon energy and corresponding
property value. It is possible for the user to
add as many material (optical) properties to the material as he
wishes using the methods supplied by the
G4MaterialPropertiesTable
class. An example of this is
shown in Example 5.4. In this example the
interpolation of the G4MaterialPropertyVector is to be done by a spline fit.
The default is a linear interpolation.
Example 5.4.
Optical properties added to a G4MaterialPropertiesTable
and linked to a G4Material
const G4int NUMENTRIES = 32; G4double ppckov[NUMENTRIES] = {2.034*eV, ......, 4.136*eV}; G4double rindex[NUMENTRIES] = {1.3435, ......, 1.3608}; G4double absorption[NUMENTRIES] = {344.8*cm, ......, 1450.0*cm]; G4MaterialPropertiesTable *MPT = new G4MaterialPropertiesTable(); MPT -> AddConstProperty("SCINTILLATIONYIELD",100./MeV); MPT -> AddProperty("RINDEX",ppckov,rindex,NUMENTRIES}->SetSpline(true); MPT -> AddProperty("ABSLENGTH",ppckov,absorption,NUMENTRIES}->SetSpline(true); scintillator -> SetMaterialPropertiesTable(MPT);
The radiation of Cerenkov light occurs when a charged particle moves through a dispersive medium faster than the group velocity of light in that medium. Photons are emitted on the surface of a cone, whose opening angle with respect to the particle's instantaneous direction decreases as the particle slows down. At the same time, the frequency of the photons emitted increases, and the number produced decreases. When the particle velocity drops below the local speed of light, the radiation ceases and the emission cone angle collapses to zero. The photons produced by this process have an inherent polarization perpendicular to the cone's surface at production.
The flux, spectrum, polarization and emission of Cerenkov
radiation in the AlongStepDoIt
method of the class
G4Cerenkov
follow well-known formulae, with two inherent
computational limitations. The first arises from step-wise
simulation, and the second comes from the requirement that
numerical integration calculate the average number of Cerenkov
photons per step. The process makes use of a
G4PhysicsTable
which contains incremental integrals to
expedite this calculation.
The time and position of Cerenkov photon emission are calculated
from quantities known at the beginning of a charged particle's
step. The step is assumed to be rectilinear even in the presence of
a magnetic field. The user may limit the step size by specifying a
maximum (average) number of Cerenkov photons created during the
step, using the SetMaxNumPhotonsPerStep(const G4int
NumPhotons)
method. The actual number generated will
necessarily be different due to the Poissonian nature of the
production. In the present implementation, the production density
of photons is distributed evenly along the particle's track
segment, even if the particle has slowed significantly during the
step. The step can also be limited with
the SetMaxBetaChangePerStep
method, where
the argument is the allowed change in percent).
The frequently very large number of secondaries produced in a
single step (about 300/cm in water), compelled the idea in
GEANT3.21 of suspending the primary particle until all its progeny
have been tracked. Despite the fact that GEANT4 employs dynamic
memory allocation and thus does not suffer from the limitations of
GEANT3.21 with its fixed large initial ZEBRA store, GEANT4
nevertheless provides for an analogous functionality with the
public method SetTrackSecondariesFirst
. An example of the
registration of the Cerenkov process is given in
Example 5.5.
Example 5.5.
Registration of the Cerenkov process in PhysicsList
.
#include "G4Cerenkov.hh" void ExptPhysicsList::ConstructOp(){ G4Cerenkov* theCerenkovProcess = new G4Cerenkov("Cerenkov"); G4int MaxNumPhotons = 300; theCerenkovProcess->SetTrackSecondariesFirst(true); theCerenkovProcess->SetMaxBetaChangePerStep(10.0); theCerenkovProcess->SetMaxNumPhotonsPerStep(MaxNumPhotons); theParticleIterator->reset(); while( (*theParticleIterator)() ){ G4ParticleDefinition* particle = theParticleIterator->value(); G4ProcessManager* pmanager = particle->GetProcessManager(); G4String particleName = particle->GetParticleName(); if (theCerenkovProcess->IsApplicable(*particle)) { pmanager->AddProcess(theCerenkovProcess); pmanager->SetProcessOrdering(theCerenkovProcess,idxPostStep); } } }
Every scintillating material has a characteristic light yield,
SCINTILLATIONYIELD
, and an intrinsic resolution,
RESOLUTIONSCALE
, which generally broadens the statistical
distribution of generated photons. A wider intrinsic resolution is
due to impurities which are typical for doped crystals like NaI(Tl)
and CsI(Tl). On the other hand, the intrinsic resolution can also
be narrower when the Fano factor plays a role. The actual number of
emitted photons during a step fluctuates around the mean number of
photons with a width given by
ResolutionScale*sqrt(MeanNumberOfPhotons)
. The average
light yield, MeanNumberOfPhotons
, has a linear dependence
on the local energy deposition, but it may be different for minimum
ionizing and non-minimum ionizing particles.
A scintillator is also characterized by its photon emission
spectrum and by the exponential decay of its time spectrum. In
GEANT4 the scintillator can have a fast and a slow component. The
relative strength of the fast component as a fraction of total
scintillation yield is given by the YIELDRATIO
.
Scintillation may be simulated by specifying these empirical
parameters for each material. It is sufficient to specify in the
user's DetectorConstruction
class a relative spectral
distribution as a function of photon energy for the scintillating
material. An example of this is shown in
Example 5.6
Example 5.6.
Specification of scintillation properties in
DetectorConstruction
.
const G4int NUMENTRIES = 9; G4double Scnt_PP[NUMENTRIES] = { 6.6*eV, 6.7*eV, 6.8*eV, 6.9*eV, 7.0*eV, 7.1*eV, 7.2*eV, 7.3*eV, 7.4*eV }; G4double Scnt_FAST[NUMENTRIES] = { 0.000134, 0.004432, 0.053991, 0.241971, 0.398942, 0.000134, 0.004432, 0.053991, 0.241971 }; G4double Scnt_SLOW[NUMENTRIES] = { 0.000010, 0.000020, 0.000030, 0.004000, 0.008000, 0.005000, 0.020000, 0.001000, 0.000010 }; G4Material* Scnt; G4MaterialPropertiesTable* Scnt_MPT = new G4MaterialPropertiesTable(); Scnt_MPT->AddProperty("FASTCOMPONENT", Scnt_PP, Scnt_FAST, NUMENTRIES); Scnt_MPT->AddProperty("SLOWCOMPONENT", Scnt_PP, Scnt_SLOW, NUMENTRIES); Scnt_MPT->AddConstProperty("SCINTILLATIONYIELD", 5000./MeV); Scnt_MPT->AddConstProperty("RESOLUTIONSCALE", 2.0); Scnt_MPT->AddConstProperty("FASTTIMECONSTANT", 1.*ns); Scnt_MPT->AddConstProperty("SLOWTIMECONSTANT", 10.*ns); Scnt_MPT->AddConstProperty("YIELDRATIO", 0.8); Scnt->SetMaterialPropertiesTable(Scnt_MPT);
In cases where the scintillation yield of a scintillator depends
on the particle type, different scintillation processes may be
defined for them. How this yield scales to the one specified for
the material is expressed with the
ScintillationYieldFactor
in the user's
PhysicsList
as shown in
Example 5.7.
In those cases where the fast to slow excitation ratio changes with particle
type, the method SetScintillationExcitationRatio
can be
called for each scintillation process (see the advanced
underground_physics example). This overwrites the
YieldRatio
obtained from the
G4MaterialPropertiesTable
.
Example 5.7.
Implementation of the scintillation process in
PhysicsList
.
G4Scintillation* theMuonScintProcess = new G4Scintillation("Scintillation"); theMuonScintProcess->SetTrackSecondariesFirst(true); theMuonScintProcess->SetScintillationYieldFactor(0.8); theParticleIterator->reset(); while( (*theParticleIterator)() ){ G4ParticleDefinition* particle = theParticleIterator->value(); G4ProcessManager* pmanager = particle->GetProcessManager(); G4String particleName = particle->GetParticleName(); if (theMuonScintProcess->IsApplicable(*particle)) { if (particleName == "mu+") { pmanager->AddProcess(theMuonScintProcess); pmanager->SetProcessOrderingToLast(theMuonScintProcess, idxAtRest); pmanager->SetProcessOrderingToLast(theMuonScintProcess, idxPostStep); } } }
A Gaussian-distributed number of photons is generated according
to the energy lost during the step. A resolution scale of 1.0
produces a statistical fluctuation around the average yield set
with AddConstProperty("SCINTILLATIONYIELD")
, while values
> 1 broaden the fluctuation. A value of zero produces no
fluctuation. Each photon's frequency is sampled from the empirical
spectrum. The photons originate evenly along the track segment and
are emitted uniformly into 4π with a random linear polarization
and at times characteristic for the scintillation component.
When there are multiple scintillators in the simulation and/or when the scintillation yield is a non-linear function of the energy deposited, the user can also define an array of total scintillation light yields as a function of the energy deposited and particle type. The available particles are protons, electrons, deuterons, tritons, alphas, and carbon ions. These are the particles known to significantly effect the scintillation light yield, of for example, BC501A (NE213/EJ301) liquid organic scintillator and BC420 plastic scintillator as function of energy deposited.
The method works as follows:
In the user's physics lists, the user must set a G4bool flag that allows scintillation light emission to depend on the energy deposited by particle type:
theScintProcess->SetScintillationByParticleType(true);
The user must also specify and add, via the AddProperty method of the MPT, the scintillation light yield as function of incident particle energy with new keywords, for example: PROTONSCINTILLATIONYIELD etc. and pairs of protonEnergy and scintLightYield.
Wavelength Shifting (WLS) fibers are used in many high-energy particle physics experiments. They absorb light at one wavelength and re-emit light at a different wavelength and are used for several reasons. For one, they tend to decrease the self-absorption of the detector so that as much light reaches the PMTs as possible. WLS fibers are also used to match the emission spectrum of the detector with the input spectrum of the PMT.
A WLS material is characterized by its photon absorption and
photon emission spectrum and by a possible time delay between the
absorption and re-emission of the photon. Wavelength Shifting may
be simulated by specifying these empirical parameters for each WLS
material in the simulation. It is sufficient to specify in the
user's DetectorConstruction
class a relative spectral
distribution as a function of photon energy for the WLS material.
WLSABSLENGTH is the absorption length of the material as a function
of the photon's energy. WLSCOMPONENT is the relative emission
spectrum of the material as a function of the photon's energy,
and WLSTIMECONSTANT accounts for any time delay which may occur
between absorption and re-emission of the photon. An example is
shown in Example 5.8.
Example 5.8.
Specification of WLS properties in DetectorConstruction
.
const G4int nEntries = 9; G4double PhotonEnergy[nEntries] = { 6.6*eV, 6.7*eV, 6.8*eV, 6.9*eV, 7.0*eV, 7.1*eV, 7.2*eV, 7.3*eV, 7.4*eV }; G4double RIndexFiber[nEntries] = { 1.60, 1.60, 1.60, 1.60, 1.60, 1.60, 1.60, 1.60, 1.60 }; G4double AbsFiber[nEntries] = {0.1*mm,0.2*mm,0.3*mm,0.4*cm,1.0*cm,10*cm,1.0*m,10.0*m,10.0*m}; G4double EmissionFiber[nEntries] = {0.0, 0.0, 0.0, 0.1, 0.5, 1.0, 5.0, 10.0, 10.0 }; G4Material* WLSFiber; G4MaterialPropertiesTable* MPTFiber = new G4MaterialPropertiesTable(); MPTFiber->AddProperty("RINDEX",PhotonEnergy,RIndexFiber,nEntries); MPTFiber->AddProperty("WLSABSLENGTH",PhotonEnergy,AbsFiber,nEntries); MPTFiber->AddProperty("WLSCOMPONENT",PhotonEnergy,EmissionFiber,nEntries); MPTFiber->AddConstProperty("WLSTIMECONSTANT", 0.5*ns); WLSFiber->SetMaterialPropertiesTable(MPTFiber);
The process is defined in the PhysicsList in the usual way. The process class name is G4OpWLS. It should be instantiated with theWLSProcess = new G4OpWLS("OpWLS") and attached to the process manager of the optical photon as a DiscreteProcess. The way the WLSTIMECONSTANT is used depends on the time profile method chosen by the user. If in the PhysicsList theWLSProcess->UseTimeGenerator("exponential") option is set, the time delay between absorption and re-emission of the photon is sampled from an exponential distribution, with the decay term equal to WLSTIMECONSTANT. If, on the other hand, theWLSProcess->UseTimeGenerator("delta") is chosen, the time delay is a delta function and equal to WLSTIMECONSTANT. The default is "delta" in case the G4OpWLS::UseTimeGenerator(const G4String name) method is not used.
The implementation of optical photon bulk absorption,
G4OpAbsorption
, is trivial in that the process merely
kills the particle. The procedure requires the user to fill the
relevant G4MaterialPropertiesTable
with empirical data for
the absorption length, using ABSLENGTH
as the property key
in the public method AddProperty
. The absorption length is
the average distance traveled by a photon before being absorpted by
the medium; i.e. it is the mean free path returned by the
GetMeanFreePath
method.
The differential cross section in Rayleigh scattering,
dσ/dω, is proportional
to 1+cos2(θ),
where θ is the polar of the new polarization vector with
respect to the old polarization vector. The G4OpRayleigh
scattering process samples this angle accordingly and then
calculates the scattered photon's new direction by requiring that
it be perpendicular to the photon's new polarization in such a way
that the final direction, initial and final polarizations are all
in one plane. This process thus depends on the particle's
polarization (spin). The photon's polarization is a data member of
the G4DynamicParticle
class.
A photon which is not assigned a polarization at production,
either via the SetPolarization
method of the
G4PrimaryParticle
class, or indirectly with the
SetParticlePolarization
method of the
G4ParticleGun
class, may not be Rayleigh scattered.
Optical photons produced by the G4Cerenkov
process have
inherently a polarization perpendicular to the cone's surface at
production. Scintillation photons have a random linear polarization
perpendicular to their direction.
The process requires a G4MaterialPropertiesTable
to be
filled by the user with Rayleigh scattering length data. The
Rayleigh scattering attenuation length is the average distance
traveled by a photon before it is Rayleigh scattered in the medium
and it is the distance returned by the GetMeanFreePath
method. The G4OpRayleigh
class provides a
RayleighAttenuationLengthGenerator
method which calculates
the attenuation coefficient of a medium following the
Einstein-Smoluchowski formula whose derivation requires the use of
statistical mechanics, includes temperature, and depends on the
isothermal compressibility of the medium. This generator is
convenient when the Rayleigh attenuation length is not known from
measurement but may be calculated from first principles using the
above material constants. For a medium named Water and no
Rayleigh scattering attenutation length specified by the user, the
program automatically calls the
RayleighAttenuationLengthGenerator
which calculates it for 10 degrees Celsius liquid water.
Mie Scattering (or Mie solution) is an analytical solution of Maxwell's equations for scattering of optical photons by spherical particles. It is significant only when the radius of the scattering object is of order of the wave length. The analytical expressions for Mie Scattering are very complicated since they are a series sum of Bessel functions. One common approximation made is call Henyey-Greenstein (HG). The implementation in Geant4 follows the HG approximation (for details see the Physics Reference Manual) and the treatment of polarization and momentum are similar to that of Rayleigh scattering. We require the final polarization direction to be perpendicular to the momentum direction. We also require the final momentum, initial polarization and final polarization to be in the same plane.
The process requires a G4MaterialPropertiesTable to be filled by the user with Mie scattering length data (entered with the name: MIEHG) analogous to Rayleigh scattering. The Mie scattering attenuation length is the average distance traveled by a photon before it is Mie scattered in the medium and it is the distance returned by the GetMeanFreePath method. In practice, the user not only needs to provide the attenuation length of Mie scattering, but also needs to provide the constant parameters of the approximation: g_f, g_b, and r_f. (with AddConstProperty and with the names: MIEHG_FORWARD, MIEHG_BACKWARD, and MIEHG_FORWARD_RATIO, respectively; see Novice Example N06.)
Reference: E. Hecht and A. Zajac, Optics [ Hecht1974 ]
For the simple case of a perfectly smooth interface between two
dielectric materials, all the user needs to provide are the
refractive indices of the two materials stored in their respective
G4MaterialPropertiesTable
. In all other cases, the optical
boundary process design relies on the concept of surfaces.
The information is split into two classes. One class in the
material category keeps information about the physical properties
of the surface itself, and a second class in the geometry category
holds pointers to the relevant physical and logical volumes
involved and has an association to the physical class. Surface
objects of the second type are stored in a related table and can be
retrieved by either specifying the two ordered pairs of physical
volumes touching at the surface, or by the logical volume entirely
surrounded by this surface. The former is called a border
surface while the latter is referred to as the skin
surface. This second type of surface is useful in situations
where a volume is coded with a reflector and is placed into many
different mother volumes. A limitation is that the skin surface can
only have one and the same optical property for all of the enclosed
volume's sides. The border surface is an ordered pair of physical
volumes, so in principle, the user can choose different optical
properties for photons arriving from the reverse side of the same
interface. For the optical boundary process to use a border
surface, the two volumes must have been positioned with
G4PVPlacement
. The ordered combination can exist at many
places in the simulation. When the surface concept is not needed,
and a perfectly smooth surface exists beteen two dielectic
materials, the only relevant property is the index of refraction, a
quantity stored with the material, and no restriction exists on how
the volumes were positioned.
The physical surface object also specifies which model the
boundary process should use to simulate interactions with that
surface. In addition, the physical surface can have a material
property table all its own. The usage of this table allows all
specular constants to be wavelength dependent. In case the surface
is painted or wrapped (but not a cladding), the table may include
the thin layer's index of refraction. This allows the simulation of
boundary effects at the intersection between the medium and the
surface layer, as well as the Lambertian reflection at the far side
of the thin layer. This occurs within the process itself and does
not invoke the G4Navigator
. Combinations of surface finish
properties, such as polished or
ground and front
painted or back painted, enumerate the different
situations which can be simulated.
When a photon arrives at a medium boundary its behavior depends on the nature of the two materials that join at that boundary. Medium boundaries may be formed between two dielectric materials or a dielectric and a metal. In the case of two dielectric materials, the photon can undergo total internal reflection, refraction or reflection, depending on the photon's wavelength, angle of incidence, and the refractive indices on both sides of the boundary. Furthermore, reflection and transmission probabilites are sensitive to the state of linear polarization. In the case of an interface between a dielectric and a metal, the photon can be absorbed by the metal or reflected back into the dielectric. If the photon is absorbed it can be detected according to the photoelectron efficiency of the metal.
As expressed in Maxwell's equations, Fresnel reflection and refraction are intertwined through their relative probabilities of occurrence. Therefore neither of these processes, nor total internal reflection, are viewed as individual processes deserving separate class implementation. Nonetheless, an attempt was made to adhere to the abstraction of having independent processes by splitting the code into different methods where practicable.
One implementation of the G4OpBoundaryProcess
class
employs the
UNIFIED model
[A. Levin and C. Moisan, A More Physical Approach
to Model the Surface Treatment of Scintillation Counters and its
Implementation into DETECT, TRIUMF Preprint TRI-PP-96-64, Oct.
1996] of the DETECT program [G.F. Knoll, T.F. Knoll and T.M.
Henderson, Light Collection Scintillation Detector Composites for
Neutron Detection, IEEE Trans. Nucl. Sci., 35 (1988) 872.]. It
applies to dielectric-dielectric interfaces and tries to provide a
realistic simulation, which deals with all aspects of surface
finish and reflector coating. The surface may be assumed as smooth
and covered with a metallized coating representing a specular
reflector with given reflection coefficient, or painted with a
diffuse reflecting material where Lambertian reflection occurs. The
surfaces may or may not be in optical contact with another
component and most importantly, one may consider a surface to be
made up of micro-facets with normal vectors that follow given
distributions around the nominal normal for the volume at the
impact point. For very rough surfaces, it is possible for the
photon to inversely aim at the same surface again after reflection
of refraction and so multiple interactions with the boundary are
possible within the process itself and without the need for
relocation by G4Navigator
.
The UNIFIED model (Figure 5.1) provides for a range of different reflection mechanisms. The specular lobe constant represents the reflection probability about the normal of a micro facet. The specular spike constant, in turn, illustrates the probability of reflection about the average surface normal. The diffuse lobe constant is for the probability of internal Lambertian reflection, and finally the back-scatter spike constant is for the case of several reflections within a deep groove with the ultimate result of exact back-scattering. The four probabilities must add up to one, with the diffuse lobe constant being implicit. The reader may consult the reference for a thorough description of the model.
Example 5.9.
Dielectric-dielectric surface properties
defined via the G4OpticalSurface
.
G4VPhysicalVolume* volume1; G4VPhysicalVolume* volume2; G4OpticalSurface* OpSurface = new G4OpticalSurface("name"); G4LogicalBorderSurface* Surface = new G4LogicalBorderSurface("name",volume1,volume2,OpSurface); G4double sigma_alpha = 0.1; OpSurface -> SetType(dielectric_dielectric); OpSurface -> SetModel(unified); OpSurface -> SetFinish(groundbackpainted); OpSurface -> SetSigmaAlpha(sigma_alpha); const G4int NUM = 2; G4double pp[NUM] = {2.038*eV, 4.144*eV}; G4double specularlobe[NUM] = {0.3, 0.3}; G4double specularspike[NUM] = {0.2, 0.2}; G4double backscatter[NUM] = {0.1, 0.1}; G4double rindex[NUM] = {1.35, 1.40}; G4double reflectivity[NUM] = {0.3, 0.5}; G4double efficiency[NUM] = {0.8, 0.1}; G4MaterialPropertiesTable* SMPT = new G4MaterialPropertiesTable(); SMPT -> AddProperty("RINDEX",pp,rindex,NUM); SMPT -> AddProperty("SPECULARLOBECONSTANT",pp,specularlobe,NUM); SMPT -> AddProperty("SPECULARSPIKECONSTANT",pp,specularspike,NUM); SMPT -> AddProperty("BACKSCATTERCONSTANT",pp,backscatter,NUM); SMPT -> AddProperty("REFLECTIVITY",pp,reflectivity,NUM); SMPT -> AddProperty("EFFICIENCY",pp,efficiency,NUM); OpSurface -> SetMaterialPropertiesTable(SMPT);
The original GEANT3.21 implementation of this process is also available via the GLISUR methods flag. [GEANT Detector Description and Simulation Tool, Application Software Group, Computing and Networks Division, CERN, PHYS260-6 tp 260-7.].
Example 5.10.
Dielectric metal surface properties defined via the
G4OpticalSurface
.
G4LogicalVolume* volume_log; G4OpticalSurface* OpSurface = new G4OpticalSurface("name"); G4LogicalSkinSurface* Surface = new G4LogicalSkinSurface("name",volume_log,OpSurface); OpSurface -> SetType(dielectric_metal); OpSurface -> SetFinish(ground); OpSurface -> SetModel(glisur); G4double polish = 0.8; G4MaterialPropertiesTable *OpSurfaceProperty = new G4MaterialPropertiesTable(); OpSurfaceProperty -> AddProperty("REFLECTIVITY",pp,reflectivity,NUM); OpSurfaceProperty -> AddProperty("EFFICIENCY",pp,efficiency,NUM); OpSurface -> SetMaterialPropertiesTable(OpSurfaceProperty);
The reflectivity off a metal surface can also be calculated by way of a complex index of refraction. Instead of storing the REFLECTIVITY directly, the user stores the real part (REALRINDEX) and the imaginary part (IMAGINARYRINDEX) as a function of photon energy separately in the G4MaterialPropertyTable. Geant4 then calculates the reflectivity depending on the incident angle, photon energy, degree of TE and TM polarization, and this complex refractive index.
The program defaults to the GLISUR model and polished
surface finish when no specific model and surface finish is
specified by the user. In the case of a dielectric-metal interface,
or when the GLISUR model is specified, the only surface finish
options available are polished or ground. For
dielectric-metal surfaces, the G4OpBoundaryProcess
also
defaults to unit reflectivity and zero detection efficiency. In
cases where the user specifies the UNIFIED model
(Figure 5.1), but does not
otherwise specify the model reflection probability constants, the
default becomes Lambertian reflection.
Martin Janecek and Bill Moses (Lawrence Berkeley National Laboratory)
built an instrument for measuring the angular reflectivity distribution
inside of BGO crystals with common surface treatments and reflectors
applied. These results have been incorporate into the Geant4 code. A
third class of reflection type besides dielectric_metal and
dielectric_dielectric is added: dielectric_LUT. The distributions have
been converted to 21 look-up-tables (LUT); so far for 1 scintillator
material (BGO) x 3 surface treatments x 7 reflector materials. The
modified code allows the user to specify the surface treatment
(rough-cut, chemically etched, or mechanically polished), the attached
reflector (Lumirror, Teflon, ESR film, Tyvek, or TiO2 paint), and the
bonding type (air-coupled or glued). The glue used is MeltMount, and the
ESR film used is VM2000. Each LUT consists of measured angular
distributions with 4º by 5º resolution in theta and phi, respectively,
for incidence angles from 0º to 90º degrees, in 1º-steps. The code might
in the future be updated by adding more LUTs, for instance, for other
scintillating materials (such as LSO or NaI). To use these LUT the user
has to download them from
Geant4 Software Download and set an environment variable,
G4REALSURFACEDATA
, to the directory of
geant4/data/RealSurface1.0
. For details see:
M. Janecek, W. W. Moses, IEEE Trans. Nucl. Sci. 57 (3) (2010) 964-970.
The enumeration G4OpticalSurfaceFinish has been extended to include (what follows should be a 2 column table):
polishedlumirrorair, // mechanically polished surface, with lumirror polishedlumirrorglue, // mechanically polished surface, with lumirror & meltmount polishedair, // mechanically polished surface polishedteflonair, // mechanically polished surface, with teflon polishedtioair, // mechanically polished surface, with tio paint polishedtyvekair, // mechanically polished surface, with tyvek polishedvm2000air, // mechanically polished surface, with esr film polishedvm2000glue, // mechanically polished surface, with esr film & meltmount etchedlumirrorair, // chemically etched surface, with lumirror etchedlumirrorglue, // chemically etched surface, with lumirror & meltmount etchedair, // chemically etched surface etchedteflonair, // chemically etched surface, with teflon etchedtioair, // chemically etched surface, with tio paint etchedtyvekair, // chemically etched surface, with tyvek etchedvm2000air, // chemically etched surface, with esr film etchedvm2000glue, // chemically etched surface, with esr film & meltmount groundlumirrorair, // rough-cut surface, with lumirror groundlumirrorglue, // rough-cut surface, with lumirror & meltmount groundair, // rough-cut surface groundteflonair, // rough-cut surface, with teflon groundtioair, // rough-cut surface, with tio paint groundtyvekair, // rough-cut surface, with tyvek groundvm2000air, // rough-cut surface, with esr film groundvm2000glue // rough-cut surface, with esr film & meltmount
To use a look-up-table, all the user needs to specify for an
G4OpticalSurface
is:
SetType(dielectric_LUT), SetModel(LUT)
and for example,
SetFinish(polishedtyvekair)
.
In this section we describe how to use the parameterization or "fast simulation" facilities of GEANT4. Examples are provided in the examples/novice/N05 directory.
The Geant4 parameterization facilities allow you to shortcut the detailed tracking in a given volume and for given particle types in order for you to provide your own implementation of the physics and of the detector response.
Parameterisations are bound to a
G4Region
object, which, in the case of fast simulation is also called an
envelope. Prior to release 8.0,
parameterisations were bound
to a G4LogicalVolume
, the root of a volume hierarchy.
These root volumes are now attributes of the G4Region
.
Envelopes often correspond to the volumes of sub-detectors:
electromagnetic calorimeters, tracking chambers, etc. With GEANT4
it is also possible to define envelopes by overlaying a parallel or
"ghost" geometry as discussed in Section 5.2.6.7.
In GEANT4, parameterisations have three main features. You must specify:
GEANT4 will message your parameterisation code for each step
starting in any root G4LogicalVolume (including daughters.
sub-daughters, etc. of this volume) of the G4Region
.
It will proceed by first asking the available parameterisations for
the current particle type if one of them (and only one) wants to
issue a trigger. If so it will invoke its parameterisation. In this
case, the tracking
will not apply physics
to the particle in the step. Instead, the UserSteppingAction will be
invoked.
Parameterisations look like a "user stepping action" but are more advanced because:
G4Region
to which it is bound;G4Region
, that is, any volume in which the track is
located;G4Region
in which the track is travelling.
The GEANT4 components which allow the implementation and control of parameterisations are:
G4VFastSimulationModel
This is the abstract class for the implementation of parameterisations. You must inherit from it to implement your concrete parameterisation model.
G4FastSimulationManager
The G4VFastSimulationModel
objects are attached to the
G4Region
through a G4FastSimulationManager
.
This object will manage the list of models and will message them at
tracking time.
G4Region/Envelope
As mentioned before, an envelope in GEANT4 is a
G4Region
.
The parameterisation is bound to the G4Region
by
setting a G4FastSimulationManager
pointer to it.
The figure below shows how the G4VFastSimulationModel
and G4FastSimulationManager
objects are bound to the
G4Region
. Then for all root G4LogicalVolume's held by
the G4Region, the fast simulation code is active.
G4FastSimulationManagerProcess
This is a G4VProcess
. It provides the interface
between the tracking and the parameterisation. It must be set in the
process list of the particles you want to parameterise.
G4GlobalFastSimulationManager
This a singleton class which provides the management of the
G4FastSimulationManager
objects and some ghost
facilities.
The G4VFastSimulationModel
class has two constructors.
The second one allows you to get started quickly:
G4VFastSimulationModel(
const G4String& aName)
:
Here aName
identifies the parameterisation model.
G4VFastSimulationModel(const G4String&
aName, G4Region*, G4bool IsUnique=false):
In addition to the model name, this constructor accepts
a G4Region
pointer. The
needed G4FastSimulationManager
object is
constructed if necessary, passing to it the G4Region pointer and the
boolean value. If it already exists, the model is simply added to this
manager. Note that the G4VFastSimulationModel
object will not keep track of
the G4Region
passed in the constructor. The
boolean argument is there for optimization purposes: if you know that
the G4Region
has a unique
root G4LogicalVolume
, uniquely placed, you can set
the boolean value to "true".
The G4VFastSimulationModel
has three pure virtual methods which
must be overriden in your concrete class:
G4VFastSimulationModel(
const G4String& aName):
Here aName identifies the parameterisation model.
G4bool ModelTrigger(
const G4FastTrack&):
You must return "true" when the dynamic conditions to trigger your parameterisation are fulfilled. G4FastTrack provides access to the current G4Track, gives simple access to the current root G4LogicalVolume related features (its G4VSolid, and G4AffineTransform references between the global and the root G4LogicalVolume local coordinates systems) and simple access to the position and momentum expressed in the root G4LogicalVolume coordinate system. Using these quantities and the G4VSolid methods, you can for example easily check how far you are from the root G4LogicalVolume boundary.
G4bool IsApplicable(
const G4ParticleDefinition&):
In your implementation, you must return "true" when your model is applicable to the G4ParticleDefinition passed to this method. The G4ParticleDefinition provides all intrinsic particle information (mass, charge, spin, name ...).
If you want to implement a model which is valid only for certain particle types, it is recommended for efficiency that you use the static pointer of the corresponding particle classes.
As an example, in a model valid for gammas only, the IsApplicable() method should take the form:
#include "G4Gamma.hh" G4bool MyGammaModel::IsApplicable(const G4ParticleDefinition& partDef) { return &partDef == G4Gamma::GammaDefinition(); }
G4bool ModelTrigger(
const G4FastTrack&):
You must return "true" when the dynamic conditions to trigger your parameterisation are fulfilled. The G4FastTrack provides access to the current G4Track, gives simple access to envelope related features (G4LogicalVolume, G4VSolid, and G4AffineTransform references between the global and the envelope local coordinates systems) and simple access to the position and momentum expressed in the envelope coordinate system. Using these quantities and the G4VSolid methods, you can for example easily check how far you are from the envelope boundary.
void DoIt(
const G4FastTrack&, G4FastStep&):
The details of your parameterisation will be implemented in this method. The G4FastTrack reference provides the input information, and the final state of the particles after parameterisation must be returned through the G4FastStep reference. Tracking for the final state particles is requested after your parameterisation has been invoked.
G4FastSimulationManager functionnalities regarding the use of ghost volumes are explained in Section 5.2.6.7.
G4FastSimulationManager(
G4Region *anEnvelope, G4bool IsUnique=false
):
This is the only constructor. You specify the G4Region by providing its pointer. The G4FastSimulationManager object will bind itself to this G4Region. If you know that this G4Region has a single root G4LogicalVolume, placed only once, you can set the IsUnique boolean to "true" to allow some optimization.
Note that if you choose to use the G4VFastSimulationModel(const G4String&, G4Region*, G4bool) constructor for your model, the G4FastSimulationManager will be constructed using the given G4Region* and G4bool values of the model constructor.
The following two methods provide the usual management functions.
void AddFastSimulationModel(
G4VFastSimulationModel*)
RemoveFastSimulationModel(
G4VFastSimulationModel*)
This G4VProcess serves as an interface between the tracking and the parameterisation. At tracking time, it collaborates with the G4FastSimulationManager of the current volume, if any, to allow the models to trigger. If no manager exists or if no model issues a trigger, the tracking goes on normally.
In the present implementation, you must set this process in the G4ProcessManager of the particles you parameterise to enable your parameterisation.
The processes ordering is:
[n-3] ... [n-2] Multiple Scattering [n-1] G4FastSimulationManagerProcess [ n ] G4Transportation
This ordering is important if you use ghost geometries, since the G4FastSimulationManagerProcess will provide navigation in the ghost world to limit the step on ghost boundaries.
The G4FastSimulationManager must be added to the process list of a particle as a continuous and discrete process if you use ghost geometries for this particle. You can add it as a discrete process if you don't use ghosts.
The following code registers the G4FastSimulationManagerProcess with all the particles as a discrete and continuous process:
void MyPhysicsList::addParameterisation() { G4FastSimulationManagerProcess* theFastSimulationManagerProcess = new G4FastSimulationManagerProcess(); theParticleIterator->reset(); while( (*theParticleIterator)() ) { G4ParticleDefinition* particle = theParticleIterator->value(); G4ProcessManager* pmanager = particle->GetProcessManager(); pmanager->AddProcess(theFastSimulationManagerProcess, -1, 0, 0); } }
This class is a singleton which can be accessed as follows:
#include "G4GlobalFastSimulationManager.hh" ... ... G4GlobalFastSimulationManager* globalFSM; globalFSM = G4GlobalFastSimulationManager::getGlobalFastSimulationManager(); ... ...
Presently, you will mainly need to use the GlobalFastSimulationManager if you use ghost geometries.
In some cases, volumes of the tracking geometry do not allow envelopes to be defined. This may be the case with a geometry coming from a CAD system. Since such a geometry is flat, a parallel geometry must be used to define the envelopes.
Another interesting case involves defining an envelope which groups the electromagnetic and hadronic calorimeters of a detector into one volume. This may be useful when parameterizing the interaction of charged pions. You will very likely not want electrons to see this envelope, which means that ghost geometries have to be organized by particle flavours.
Using ghost geometries implies some more overhead in the parameterisation mechanism for the particles sensitive to ghosts, since navigation is provided in the ghost geometry by the G4FastSimulationManagerProcess. Usually, however, only a few volumes will be placed in this ghost world, so that the geometry computations will remain rather cheap.
In the existing implementation (temporary implementation with G4Region but before parallel geometry implementation), you may only consider ghost G4Regions with just one root G4LogicalVolume. The G4GlobalFastSimulationManager provides the construction of the ghost geometry by making first an empty "clone" of the world for tracking provided by the construct() method of your G4VUserDetectorConstruction concrete class. You provide the placement of the G4Region root G4LogicalVolume relative to the ghost world coordinates in the G4FastSimulationManager objects. A ghost G4Region is recognized by the fact that its associated G4FastSimulationManager retains a non-empty list of placements.
The G4GlobalFastSimulationManager will then use both those placements and the IsApplicable() methods of the models attached to the G4FastSimulationManager objects to build the flavour-dependant ghost geometries.
Then at the beginning of the tracking of a particle, the appropriate ghost world, if any, will be selected.
The steps required to build one ghost G4Region are:
built the ghost G4Region : myGhostRegion;
build the root G4LogicalVolume: myGhostLogical, set it to myGhostRegion;
build a G4FastSimulationManager object, myGhostFSManager, giving myGhostRegion as argument of the constructor;
give to the G4FastSimulationManager the placement of the myGhostLogical, by invoking for the G4FastSimulationManager method:
AddGhostPlacement(G4RotationMatrix*, const G4ThreeVector&);
or:
AddGhostPlacement(G4Transform3D*);
where the rotation matrix and translation vector of the 3-D transformation describe the placement relative to the ghost world coordinates.
build your G4VFastSimulationModel objects and add them to the myGhostFSManager. The IsApplicable() methods of your models will be used by the G4GlobalFastSimulationManager to build the ghost geometries corresponding to a given particle type.
Invoke the G4GlobalFastSimulationManager method:
G4GlobalFastSimulationManager::getGlobalFastSimulationManager()-> CloseFastSimulation();
This last call will cause the G4GlobalFastSimulationManager to build the flavour-dependent ghost geometries. This call must be done before the RunManager closes the geometry. (It is foreseen that the run manager in the future will invoke the CloseFastSimulation() to synchronize properly with the closing of the geometry).
Visualization facilities are provided for ghosts geometries. After the CloseFastSimulation() invocation, it is possible to ask for the drawing of ghosts in an interactive session. The basic commands are:
/vis/draw/Ghosts particle_name
which makes the drawing of the ghost geometry associated with the particle specified by name in the command line.
/vis/draw/Ghosts
which draws all the ghost geometries.
This section describes how to use the Gflash library. Gflash is a concrete parameterization which is based on the equations and parameters of the original Gflash package from H1(hep-ex/0001020, Grindhammer & Peters, see physics manual) and uses the "fast simulation" facilities of GEANT4 described above. Briefly, whenever a e-/e+ particle enters the calorimeter, it is parameterized if it has a minimum energy and the shower is expected to be contained in the calorimeter (or " parameterization envelope"). If this is fulfilled the particle is killed, as well as all secondaries, and the energy is deposited according to the Gflash equations. An example, provided in examples/extended/parametrisation/gflash/, shows how to interface Gflash to your application. The simulation time is measured, so the user can immediately see the speed increase resulting from the use of Gflash.
To use Gflash "out of the box" the following steps are necessary:
The user must add the fast simulation process to his process manager:
void MyPhysicsList::addParameterisation() { G4FastSimulationManagerProcess* theFastSimulationManagerProcess = new G4FastSimulationManagerProcess(); theParticleIterator->reset(); while( (*theParticleIterator)() ) { G4ParticleDefinition* particle = theParticleIterator->value(); G4ProcessManager* pmanager = particle->GetProcessManager(); pmanager->AddProcess(theFastSimulationManagerProcess, -1, 0, 0); } }
The envelope in which the parameterization should be performed must be specified (below: G4Region m_calo_region) and the GFlashShowerModel must be assigned to this region. Furthermore, the classes GFlashParticleBounds (which provides thresholds for the parameterization like minimal energy etc.), GflashHitMaker(a helper class to generate hits in the sensitive detector) and GFlashHomoShowerParamterisation (which does the computations) must be constructed (by the user at the moment) and assigned to the GFlashShowerModel. Please note that at the moment only homogeneous calorimeters are supported.
m_theFastShowerModel = new GFlashShowerModel("fastShowerModel",m_calo_region); m_theParametrisation = new GFlashHomoShowerParamterisation(matManager->getMaterial(mat)); m_theParticleBounds = new GFlashParticleBounds(); m_theHMaker = new GFlashHitMaker(); m_theFastShowerModel->SetParametrisation(*m_theParametrisation); m_theFastShowerModel->SetParticleBounds(*m_theParticleBounds) ; m_theFastShowerModel->SetHitMaker(*m_theHMaker);
The user must also set the material of the calorimeter, since the computation depends on the material.
It is mandatory to use G4VGFlashSensitiveDetector as (additional) base class for the sensitive detector.
class ExGflashSensitiveDetector: public G4VSensitiveDetector ,public G4VGFlashSensitiveDetector
Here it is necessary to implement a separate interface, where the GFlash spots are processed.
(ProcessHits(G4GFlashSpot*aSpot ,G4TouchableHistory* ROhist))
A separate interface is used, because the Gflash spots naturally contain less information than the full simulation.
Since the parameters in the Gflash package are taken from fits to full simulations with Geant3, some retuning might be necessary for good agreement with Geant4 showers. For experiment-specific geometries some retuning might be necessary anyway. The tuning is quite complicated since there are many parameters (some correlated) and cannot be described here (see again hep-ex/0001020). For brave users the Gflash framework already forsees the possibility of passing a class with the (users) parameters,GVFlashHomoShowerTuning, to the GFlashHomoShowerParamterisation constructor. The default parameters are the original Gflash parameters:
GFlashHomoShowerParameterisation(G4Material * aMat, GVFlashHomoShowerTuning * aPar = 0);
Now there is also a preliminary implemenation of a parameterization for sampling calorimeters.
The user must specify the active and passive material, as well as the thickness of the active and passive layer.
The sampling structure of the calorimeter is taken into account by using an "effective medium" to compute the shower shape.
All material properties needed are calculated automatically. If tuning is required, the user can pass his own parameter set in the class GFlashSamplingShowerTuning. Here the user can also set his calorimeter resolution.
All in all the constructor looks the following:
GFlashSamplingShowerParamterisation(G4Material * Mat1, G4Material * Mat2,G4double d1,G4double d2, GVFlashSamplingShowerTuning * aPar = 0);
An implementation of some tools that should help the user to tune the parameterization is forseen.
To be delivered by J. Apostolakis (<John.Apostolakis@cern.ch>
).