Optimal exploitation of hadronic final states played a key role in successes of all recent collider experiment in HEP, and the ability to use hadronic final states will continue to be one of the decisive issues during the analysis phase of the LHC experiments. Monte Carlo programs like Geant4 facilitate the use of hadronic final states, and have been developed for many years.
We give an overview of the Object Oriented frameworks for hadronic generators in Geant4, and illustrate the physics models underlying hadronic shower simulation on examples, including the three basic types of modeling; data-driven, parametrisation-driven, and theory-driven modeling, and their possible realisations in the Object Oriented component system of Geant4. We put particular focus on the level of extendibility that can and has been achieved by our Russian dolls approach to Object Oriented design, and the role and importance of the frameworks in a component system.
The purpose of this section is to explain the implementation frameworks used in and provided by Geant4 for hadronic shower simulation as in the 1.1 release of the program. The implementation frameworks follow the Russian dolls approach to implementation framework design. A top-level, very abstracting implementation framework provides the basic interface to the other Geant4 categories, and fulfils the most general use-case for hadronic shower simulation. It is refined for more specific use-cases by implementing a hierarchy of implementation frameworks, each level implementing the common logic of a particular use-case package in a concrete implementation of the interface specification of one framework level above, this way refining the granularity of abstraction and delegation. This defines the Russian dolls architectural pattern. Abstract classes are used as the delegation mechanism [1] .
All framework functional requirements were obtained through use-case analysis. In the following we present for each framework level the compressed use-cases, requirements, designs including the flexibility provided, and illustrate the framework functionality with examples. All design patterns cited are to be read as defined in [ Gamma1995 ].
There are two principal use-cases of the level 1 framework. A user will want to choose the processes used for his particular simulation run, and a physicist will want to write code for processes of his own and use these together with the rest of the system in a seamless manner.
Both requirements are implemented in a sub-set of the tracking-physics interface in Geant4}. The class diagram is shown in Figure 3.3.
All processes have a common base-class G4VProcess
,
from which a set of specialised classes are derived. Two of them are
used as base classes for hadronic processes at rest and in flight
(G4VDiscreteProcess
), and for processes like
radioactive decay where the same implementation can represent both these
extreme cases (G4VRestDiscreteProcess
).
Each of these classes declares two types of methods; one for calculating the time to interaction or the physical interaction length, allowing tracking to request the information necessary to decide on the process responsible for final state production, and one to compute the final state. These are pure virtual methods, and have to be implemented in each individual derived class, as enforced by the compiler.
Note on at-rest processes: starting with Geant4 version 9.6
-
when the Bertini and Fritiof final-state models have been extended
down to zero kinetic energy and used also for simulating the nuclear capture
at-rest - the at-rest processes derive from G4HadronicProcess
,
hence from G4VDiscreteProcess
, instead than from
G4VRestProcess
as in the initial design of at-rest processes.
This requires some adaptation a discrete process to handle
an at-rest one using top level interface G4VProcess
.
A different solution, under consideration but not yet
implemented, would be instead to have G4HadronicProcess
inheriting from G4VRestDiscreteProcess
: in this way,
G4HadronicProcess
, and therefore any theory-driven
final-state model, could be deployed for any kind of hadronic process,
including capture-at-rest processes and radioactive decays.
The functionality provided is through the use of process base-class pointers
in the tracking-physics interface, and the
G4Process-Manager
.
All functionality is implemented in abstract, and registration of derived process
classes with the G4Process-Manager
of an individual particle
allows for arbitrary combination of both Geant4 provided processes, and
user-implemented processes. This registration mechanism is a modification on a
Chain of Responsibility. It is outside the scope of the current paper, and its
description is available from
G4Manual.
At the next level of abstraction, only processes that occur for particles in flight are considered. For these, it is easily observed that the sources of cross sections and final state production are rarely the same. Also, different sources will come with different restrictions. The principal use-cases of the framework are addressing these commonalities. A user might want to combine different cross sections and final state or isotope production models as provided by Geant4, and a physicist might want to implement his own model for particular situation, and add cross-section data sets that are relevant for his particular analysis to the system in a seamless manner.
The above requirements are implemented in three framework components, one for cross-sections, final state production, and isotope production each. The class diagrams are shown in Figure 3.4 for the cross-section aspects, Figure 3.5 for the final state production aspects, and figure Figure 3.6 for the isotope production aspects.
Figure 3.4. Level 2 implementation framework of the hadronic category of Geant4; cross-section aspect.
Figure 3.5. Level 2 implementation framework of the hadronic category of Geant4; final state production aspect.
Figure 3.6. Level 2 implementation framework of the hadronic category of Geant4; isotope production aspect
The three parts are integrated in the G4Hadronic-Process
class, that serves as base-class for all hadronic processes of particles in flight.
Each hadronic process is derived from G4Hadronic-Process}
,
which holds a list of ``cross section data sets''.
The term ``data set'' is
representing an object that encapsulates methods and data for calculating
total cross sections for a given process in a certain range of validity.
The implementations may take any form. It can be a simple equation as well as
sophisticated parameterisations, or evaluated data.
All cross section data set classes are derived from the abstract class
G4VCrossSection-DataSet}
, which declares methods that allow
the process inquire, about the applicability of an individual data-set through
IsApplicable(const G4DynamicParticle*, const G4Element*)
,
and to delegate the calculation of the actual cross-section value through
GetCrossSection(const G4DynamicParticle*, const G4Element*)
.
The G4HadronicInteraction
base class is provided for final state
generation. It declares a minimal interface of only one pure virtual method:
G4VParticleChange* ApplyYourself(const G4Track &, G4Nucleus &)}.
G4HadronicProcess
provides a registry for final state
production models inheriting from G4Hadronic-Interaction
.
Again, final state production model is meant in very general terms. This can
be an implementation of a quark gluon string model
[QGSM],
a sampling code for ENDF/B data formats
[
ENDFForm
],
or a parametrisation describing only
neutron elastic scattering off Silicon up to 300~MeV.
For isotope production, a base class (G4VIsotope-Production
)
is provided. It declares a method
G4IsoResult * GetIsotope(const G4Track &, const G4Nucleus &)
that calculates and returns the isotope production information. Any concrete
isotope production model will inherit from this class, and implement the
method. Again, the modeling possibilities are not limited, and the
implementation of concrete production models is not restricted in any way.
By convention, the GetIsotope
method returns NULL, if the model
is not applicable for the current projectile target combination.
G4HadronicProcess
provides registering possibilities
for data sets. A default is provided covering all possible
conditions to some approximation. The process stores and retrieves
the data sets through a data store that acts like a FILO stack
(a Chain of Responsibility with a
First In Last Out decision strategy). This allows the user to map out
the entire parameter space by overlaying cross section data sets to optimise
the overall result. Examples are the cross sections for low energy neutron
transport. If these are registered last by the user, they will be used
whenever low energy neutrons are encountered. In all other conditions the
system falls back on the default, or data sets with earlier registration dates.
The fact that the registration is done through abstract base classes with no
side-effects allows the user to implement and use his own cross sections.
Examples are special reaction cross sections of
K0-nuclear interactions
that might be used for ∊/∊' analysis at LHC to control the
systematic error.
The G4HadronicProcess
class provides a
registration service for classes deriving from
G4Hadronic-Interaction
, and delegates final state
production to the applicable model.
G4Hadronic-Interaction
provides the functionality needed
to define and enforce the applicability of a particular model.
Models inheriting from G4Hadronic-Interaction
can be restricted in applicability in projectile
type and energy, and can be activated/deactivated for individual materials and
elements. This allows a user to use final state production models in
arbitrary combinations, and to write his own models for
final state production. The design is a variant of a Chain of Responsibility.
An example would be the likely CMS scenario - the combination of low energy
neutron transport with a quantum molecular dynamics
[QMD],
invariant phase space decay
[CHIPS],
and fast parametrised models for calorimeter materials, with user defined
modeling of interactions of spallation nucleons with the most abundant
tracker and calorimeter materials.
The G4HadronicProcess
by default calculates the isotope production
information from the final state given by the transport model. In addition, it
provides a registering mechanism for isotope production models that run in
parasitic mode to the transport models and inherit from
G4VIsotope-Production
. The registering mechanism behaves like a FILO
stack, again based on Chain of Responsibility. The models will be asked for
isotope production information in inverse order of registration. The first
model that returns a non-NULL value will be applied. In addition, the
G4Hadronic-Process
provides the basic infrastructure for accessing and
steering of isotope-production information. It allows to enable and disable
the calculation of isotope production information globally, or for individual
processes, and to retrieve the isotope production information through the
G4IsoParticleChange * GetIsotopeProductionInfo()}
method at the end of each step. The G4HadronicProcess
is a finite state
machine that will ensure the GetIsotope-ProductionInfo
returns a
non-zero value only at the first call after isotope production occurred. An
example of the use of this functionality is the study of activation of a
Germanium detector in a high precision, low background experiment.
Figure 3.7. Level 3 implementation framework of the hadronic category of Geant4; theoretical model aspect.
Geant4 provides at present one implementation framework for theory
driven models. The main use-case is that of a user wishing to use theoretical
models in general, and to use various intra-nuclear transport or pre-compound
models together with models simulating the initial interactions at very high
energies.
To provide the above flexibility, the following abstract base classes have been implemented:
G4VHighEnergyGenerator
G4VIntranuclearTransportModel
G4VPreCompoundModel
G4VExcitationHandler
In addition, the class G4TheoFS-Generator
is provided to orchestrate
interactions between these classes. The class diagram is shown in
Figure 3.7.
G4VHighEnergy-Generator
serves as base class for parton transport or
parton string models, and for Adapters to event generators. This class
declares two methods, Scatter
, and
GetWoundedNucleus
.
The base class for INT inherits from G4Hadronic-Interaction
,
making any concrete implementation usable as a stand-alone model. In doing so, it
re-declares the ApplyYourself
interface of
G4Hadronic-Interaction
,
and adds a second interface, Propagate
, for further propagation
after high energy interactions. Propagate
takes as arguments a
three-dimensional model of a wounded nucleus, and a set of secondaries with
energies and positions.
The base class for pre-equilibrium decay models, G4VPre-CompoundModel
,
inherits from G4Hadronic-Interaction
, again making any concrete
implementation usable as stand-alone model. It adds a pure virtual
DeExcite
method for further evolution of the system when
intra-nuclear transport assumptions break down.
DeExcite
takes a G4Fragment
,
the Geant4 representation of an excited nucleus, as argument.
The base class for evaporation phases, G4VExcitation-Handler
,
declares an abstract method, BreakItUP()
, for compound decay.
The G4TheoFSGenerator
class inherits from
G4Hadronic-Interaction
,
and hence can be registered as a model for final state production with a
hadronic process. It allows a concrete implementation of
G4VIntranuclear-TransportModel
and
G4VHighEnergy-Generator
to be
registered, and delegates initial interactions, and intra-nuclear transport
of the corresponding secondaries to the respective classes. The design is a
complex variant of a Strategy. The most spectacular application of this
pattern is the use of parton-string models for string excitation, quark
molecular dynamics for correlated string decay, and quantum molecular dynamics
for transport, a combination which promises to result in a coherent
description of quark gluon plasma in high energy nucleus-nucleus interactions.
The class G4VIntranuclearTransportModel
provides
registering mechanisms for concrete implementations of
G4VPreCompound-Model
, and provides
concrete intra-nuclear transports with the possibility of delegating
pre-compound decay to these models.
G4VPreCompoundModel
provides a registering mechanism
for compound decay through the
G4VExcitation-Handler
interface, and provides concrete
implementations with the possibility of delegating the decay of a compound
nucleus.
The concrete scenario of G4TheoFS-Generator
using a
dual parton model
and a classical cascade, which in turn uses an exciton pre-compound model that
delegates to an evaporation phase, would be the following:
G4TheoFS-Generator
receives the conditions of the interaction;
a primary particle and a nucleus. It asks the dual parton model to perform the
initial scatterings, and return the final state, along with the by then
damaged nucleus. The nucleus records all information on the damage sustained.
G4TheoFS-Generator
forwards all information to the classical cascade,
that propagates the particles in the already damaged nucleus, keeping track of
interactions, further damage to the nucleus, etc.. Once the cascade assumptions
break down, the cascade will collect the information of the current state of
the hadronic system, like excitation energy and number of excited particles,
and interpret it as a pre-compound system. It delegates the decay of this to
the exciton model. The exciton model will take the information provided, and
calculate
transitions and emissions, until the number of excitons in the system equals
the mean number of excitons expected in equilibrium for the current excitation
energy. Then it hands over to the
evaporation phase. The evaporation phase decays the residual nucleus, and the Chain of
Command rolls back to G4TheoFS-Generator
, accumulating the
information produced in the various levels of delegation.
Figure 3.8. Level 4 implementation framework of the hadronic category of Geant4; parton string aspect.
Figure 3.9. Level 4 implementation framework of the hadronic category of Geant4; intra-nuclear transport aspect.
The use-cases of this level are related to commonalities and detailed choices in string-parton models and cascade models. They are use-cases of an expert user wishing to alter details of a model, or a theoretical physicist, wishing to study details of a particular model.
To fulfil the requirements on string models, two abstract classes are provided,
the G4VParton-StringModel
and the
G4VString-Fragmentation
.
The base class for parton string models, G4VParton-StringModel
,
declares the GetStrings()
pure virtual method.
G4VString-Fragmentation
, the pure abstract base class for string
fragmentation, declares the interface for string fragmentation.
To fulfill the requirements on intra-nuclear transport, two abstract classes
are provided, G4V3DNucleus
, and G4VScatterer
.
At this point in time, the usage of these intra-nuclear transport related classes
by concrete codes is not enforced by designs, as the details of the
cascade loop are still model dependent, and more experience has to be gathered to achieve
standardisation. It is within the responsibility of the implementers of concrete
intra-nuclear transport
codes to use the abstract interfaces as defined in these classes.
The class diagram is shown in Figure 3.8 for the string parton model aspects, and in Figure 3.9 for the intra-nuclear transport.
Again variants of Strategy, Bridge and Chain of Responsibility are used.
G4VParton-StringModel
implements the initialisation of a
three-dimensional model of a nucleus, and the logic of scattering. It
delegates secondary production to string fragmentation through a
G4VString-Fragmentation
pointer. It provides a registering service for
the concrete string fragmentation, and delegates the string excitation to
derived classes. Selection of string excitation is through selection of
derived class. Selection of string fragmentation is through registration.
Figure 3.10. Level 5 implementation framework of the hadronic category of Geant4; string fragmentation aspect.
The use-case of this level is that of a user or theoretical physicist wishing to understand the systematic effects involved in combining various fragmentation functions with individual types of string fragmentation. Note that this framework level is meeting the current state of the art, making extensions and changes of interfaces in subsequent releases likely.
A base class for fragmentation functions,
G4VFragmentation-Function}
, is
provided. It declares the GetLightConeZ()
interface.
The design is a basic Strategy. The class diagram is shown in
Figure 3.10.
At this point in time, the usage of the
G4VFragmentation-Function
is not enforced by design,
but made available from the
G4VString-Fragmentation
to an implementer of a concrete string
decay. G4VString-Fragmentation
provides a registering
mechanism for the
concrete fragmentation function. It delegates the calculation of
zf of the
hadron to split of the string to the concrete implementation. Standardisation
in this area is expected.
For some applications Geant4 might not provide the most appropriate physics implementantion or, in fact, any physics implementation at all. In such cases, it is up to the user to develop the necessary processes, models and cross sections and integrate them into his version of the Geant4 toolkit. The user's process, model or cross section may then be used in a physics list to replace some of those already provided by the toolkit. This modularity requires that user classes be derived from a set of base classes which have been provided to aid integration with the toolkit and to spare the user from consideration of many details not related to the physics at hand.
Processes communicate with Geant4 tracking, telling it where or when an interaction is supposed to occur, and what is supposed to happen at that point. A hadronic process may be implemented directly, or through the use af a framework of classes that modularize physics functionality, make available several utilities and reduce unneccesary code duplication. In the latter, recommended, approach the user must in general develop as many as three classes: a process, a cross section and a model. Instances of the cross section and model classes must then be assigned to the process. In practice it is usually necessary to develop only a model class or a cross section class, since a number of processes are already provided by Geant4. Before writing any code, users should check that Geant4 has not already povided the necessary models, cross sections or processes.
A hadronic model is repsonsible for the generation of a set of final state four-vectors, given an initial projectile and target. New models should derive from the G4HadronicInteraction base class and at least two methods in this class must be implemented:
virtual G4HadFinalState* ApplyYourself(const G4HadProjectile& aProjectile, G4Nucleus& targetNucleus)
which is responsible for generating the final state of the interaction including the specification of all particle types and four-momenta, and
virtual G4bool IsApplicable(const G4HadProjectile& aProjectile, G4Nucleus& targetNucleus)
which is responsible for checking that the incident particle type and energy, and the Z and A of the target nucleus, can be adequately handled by the model. G4HadronicInteraction provides a number of utilities to aid in the implementation of these methods.
When implementing ApplyYourself(), the Get() methods of the G4HadProjectile and G4Nucleus classes provide all the initial state information necessary for the generation of the final state. For G4HadProjectile, GetDefinition() provides the particle type, and Get4Momentum() provides the total energy and momentum. For G4Nucleus, GetZ_asInt(), GetN_asInt() and GetA_asInt() provide Z, N and A, while AtomicMass() provides the mass. Additional utility methods are available for both G4HadProjectile and G4Nucleus .
The inputs to the model assume that the incident particle (G4HadProjectile) travels along the z axis and interacts with the target (G4Nucleus) which is at rest in the lab frame. Before invoking the ApplyYourself() method, the process rotates the direction of the projectile to be along the z axis and then performs the inverse rotation on the final state particles produced by the model.
The model must perform two additional transformations: into the CM frame, and back out of it after the interaction is complete.
It is thus the model developer's responsibility to:
Step 4) is accomplished by using the various Get() and Set() methods provided by the class G4HadFinalState. The developer must also decide whether the original projectile survives the interaction, disappears or is suspended. This is done with the SetStatusChange() method. If the projectile survives, the change in its energy and momentum must be set with the provided methods and it must be flagged as "isAlive". If the particle disappears it must be flagged as "stopAndKill".
Geant4 provides a large number of Lorentz transformation tools which may be used to complete steps 1) and 3).
How step 2) is accomplished is entirely up to the developer. This could be as simple as a look-up table which assigns a final state to an initial state, or as complex as a theoretical high energy generator. Typically, the user will have to provide methods of sampling final state multiplicities, energies and angles using random number generators provided by Geant4. For example, the cosine of the polar angle of an isotropic angular distribution could be sampled as follows:
G4double cosTheta = 2.*G4UniformRandom() - 1.;
The developer must also see to it that the model conserves energy and momentum. Currently the hadronic framework checks that final states do not exceed reasonably small limits of non-conservation.
For complex models it is recommended that the user become familiar with the Geant4 hadronic framework. This is covered in detail in the chapter on Extended Functionality in this manual. The framework uses the object-oriented principles of abstraction and re-use to provide a number of services to the developer. The part of the framework used will depend on the type of model. For example, high energy models can take advantage of already-developed string excitation and decay functions and medium energy models can use the intra-nuclear propagation base class and the nuclear de-excitation handler.
The is a straightforward, but important, method. Most models are quite specific in their range of use and the developer must codify this. It is recommended that this method test for ranges of projectile energy, particle type and target atomic number and weight, and return false when these ranges are exceeded.
New cross section sets should derive from G4VCrossSectionDataSet . This class serves as a container of cross section data and provides a number of access methods that must be implemented by the developer. The essential methods are:
G4double GetElementCrossSection(const G4DynamicParticle*, G4int Z)
which retrieves element-based cross sections,
G4double GetIsoCrossSection(const G4DynamicParticle*, G4int Z, G4int A)
which retrieves isotope-based cross sections,
G4bool IsElementApplicable(const G4DynamicParticle*, G4int Z)
which sets the Z range of the element-based data set,
G4bool IsIsoApplicable(const G4DynamicParticle*, G4int, G4int A)
which sets the Z and A range of the isotope-based data set, and
SetMinKinEnergy(G4double) SetMaxKinEnergy(G4double) GetMinKinEnergy() and GetMaxKinEnergy()
for defining the applicable energy range of the data set.
As mentioned above, it is preferable to add new physics in terms of a model, and assign the model to an existing process, rather than develop a new, specific process. Under certain circumstances though, a directly implemented process may be necessary. In that case it must derive from G4HadronicProcess and three methods of that class must be implemented:
virtual G4VParticleChange* PostStepDoIt(const G4Track&, const G4Step&) , virtual G4bool IsApplicable(const G4ParticleDefinition&) , and G4double GetMeanFreePath(const G4Track& aTrack, G4double, G4ForceCondition*).
PostStepDoIt() is responsible for generating the final state of an interaction given the track and step information. It must update the state of the track, flagging it as "Alive", "StopButAlive", "StopAndKill", "KillTrackAndSecondaries", "Suspend", or "PostponeToNextEvent". It is roughly analagous to the ApplyYourself() method in models.
IsApplicable() serves the same purpose in processes as it does in models.
GetMeanFreePath() gets the cross section as a function of particle type, energy, and target material, and converts it to a mean free path, which is in turn passed on to the tracking. This method can be quite simple:
G4double particle = aTrack.GetDynamicParticle(); G4double material = aTrack.GetMaterial(); return factor/theCrossSectionDataStore->GetCrossSection(particle, material);
provided an appropriate cross section data set is already available in the data store.
The above discussion refers to in-flight processes. At-rest hadronic processes do not currently derive from G4HadronicProcess, but from G4VRestProcess or G4VRestDiscreteProcess. As such they do not employ the full hadronic framework and must be implemented directly without models. The methods to be implemented are similar to those in the in-flight case:
virtual G4VParticleChange* AtRestDoIt(const G4Track&, const G4Step&) , and G4bool IsApplicable(const G4ParticleDefinition&) .
[1] The same can be achieved with template specialisations with slightly improved CPU performance but at the cost of significantly more complex designs and, with present compilers, significantly reduced portability.