Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4HadronicInteraction.hh
Go to the documentation of this file.
1//
2// ********************************************************************
3// * License and Disclaimer *
4// * *
5// * The Geant4 software is copyright of the Copyright Holders of *
6// * the Geant4 Collaboration. It is provided under the terms and *
7// * conditions of the Geant4 Software License, included in the file *
8// * LICENSE and available at http://cern.ch/geant4/license . These *
9// * include a list of copyright holders. *
10// * *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work make any representation or warranty, express or implied, *
14// * regarding this software system or assume any liability for its *
15// * use. Please see the license in the file LICENSE and URL above *
16// * for the full disclaimer and the limitation of liability. *
17// * *
18// * This code implementation is the result of the scientific and *
19// * technical work of the GEANT4 collaboration. *
20// * By using, copying, modifying or distributing the software (or *
21// * any work based on the software) you agree to acknowledge its *
22// * use in resulting scientific publications, and indicate your *
23// * acceptance of all terms of the Geant4 Software license. *
24// ********************************************************************
25//
26//
27// $Id$
28//
29// Hadronic Interaction abstract base class
30// This class is the base class for the model classes.
31// It sorts out the energy-range for the models and provides
32// class utilities.
33// original by H.P. Wellisch
34// Modified by J.L.Chuma, TRIUMF, 21-Mar-1997
35// Last modified: 3-Apr-1997
36// Added units to energy initialization: J.L. Chuma 04-Apr-97
37// Modified by J.L.Chuma, 05-May-97 to Initialize theBlockedCounter
38// Modified by J.L.Chuma, 08-Jul-97 to implement the Nucleus changes
39// Adding a registry for memory management of hadronic models, HPW 22-Mar-99
40// 23-Jan-2009 V.Ivanchenko move constructor and destructor to the body
41// and reorder methods in the header
42// 29-Jun-2009 V.Ivanchenko add SampleInvariantT method
43// 29-Aug-2009 V.Ivanchenko moved G4ReactionDynamics to G4InelasticInteraction,
44// add const pointers, and added recoilEnergyThreshold
45// member and accesors
46
47// Class Description
48// This is the base class for all hadronic interaction models in geant4.
49// If you want to implement a new way of producing a final state, please,
50// inherit from here.
51// Class Description - End
52
53#ifndef G4HadronicInteraction_h
54#define G4HadronicInteraction_h 1
55
56#include "G4HadFinalState.hh"
57#include "G4Material.hh"
58#include "G4Nucleus.hh"
59#include "G4Track.hh"
60#include "G4HadProjectile.hh"
61#include "G4ReactionDynamics.hh"
62
64{
65public: // With description
66
67 G4HadronicInteraction(const G4String& modelName = "HadronicModel");
68
69 virtual ~G4HadronicInteraction();
70
72 G4Nucleus & targetNucleus ) = 0;
73 // The interface to implement for final state production code.
74
76 G4double plab,
77 G4int Z, G4int A);
78 // The interface to implement sampling of scattering or change exchange
79
80 virtual G4bool IsApplicable(const G4HadProjectile &/*aTrack*/,
81 G4Nucleus & /*targetNucleus*/)
82 { return true;}
83
84 inline G4double GetMinEnergy() const
85 { return theMinEnergy; }
86
87 G4double GetMinEnergy( const G4Material *aMaterial,
88 const G4Element *anElement ) const;
89
90 inline void SetMinEnergy( G4double anEnergy )
91 { theMinEnergy = anEnergy; }
92
93 void SetMinEnergy( G4double anEnergy, const G4Element *anElement );
94
95 void SetMinEnergy( G4double anEnergy, const G4Material *aMaterial );
96
97 inline G4double GetMaxEnergy() const
98 { return theMaxEnergy; }
99
100 G4double GetMaxEnergy( const G4Material *aMaterial,
101 const G4Element *anElement ) const;
102
103 inline void SetMaxEnergy( const G4double anEnergy )
104 { theMaxEnergy = anEnergy; }
105
106 void SetMaxEnergy( G4double anEnergy, const G4Element *anElement );
107
108 void SetMaxEnergy( G4double anEnergy, const G4Material *aMaterial );
109
111 { return this; }
112
113 inline G4int GetVerboseLevel() const
114 { return verboseLevel; }
115
116 inline void SetVerboseLevel( G4int value )
117 { verboseLevel = value; }
118
119 inline const G4String& GetModelName() const
120 { return theModelName; }
121
122 void DeActivateFor(const G4Material* aMaterial);
123
124 inline void ActivateFor( const G4Material *aMaterial )
125 {
126 Block();
127 SetMaxEnergy(GetMaxEnergy(), aMaterial);
128 SetMinEnergy(GetMinEnergy(), aMaterial);
129 }
130
131 void DeActivateFor( const G4Element *anElement );
132 inline void ActivateFor( const G4Element *anElement )
133 {
134 Block();
135 SetMaxEnergy(GetMaxEnergy(), anElement);
136 SetMinEnergy(GetMinEnergy(), anElement);
137 }
138
139 G4bool IsBlocked( const G4Material *aMaterial ) const;
140 G4bool IsBlocked( const G4Element *anElement) const;
141
143 { recoilEnergyThreshold = val; }
144
146 { return recoilEnergyThreshold;}
147
148 inline G4bool operator==(const G4HadronicInteraction &right ) const
149 { return ( this == (G4HadronicInteraction *) &right ); }
150
151 inline G4bool operator!=(const G4HadronicInteraction &right ) const
152 { return ( this != (G4HadronicInteraction *) &right ); }
153
154 virtual const std::pair<G4double, G4double> GetFatalEnergyCheckLevels() const;
155
156 virtual std::pair<G4double, G4double> GetEnergyMomentumCheckLevels() const;
157
158 inline void SetEnergyMomentumCheckLevels(G4double relativeLevel, G4double absoluteLevel)
159 { epCheckLevels.first = relativeLevel;
160 epCheckLevels.second = absoluteLevel; }
161
162 virtual void ModelDescription(std::ostream& outFile) const ; //=0;
163
164private:
165
167 const G4HadronicInteraction& operator=(const G4HadronicInteraction &right);
168
169protected:
170
171 inline void SetModelName(const G4String& nam)
172 { theModelName = nam; }
173
174 inline G4bool IsBlocked() const { return isBlocked;}
175 inline void Block() { isBlocked = true; }
176
178 // the G4HadFinalState object which is modified and returned
179 // by address by the ApplyYourself method,
180 // (instead of aParticleChange as found in G4VProcess)
181
183 // control flag for output messages
184 // 0: silent
185 // 1: warning messages
186 // 2: more
187 // (instead of verboseLevel as found in G4VProcess)
188
189 // these two have global validity energy range
192
194
195private:
196
197 G4double recoilEnergyThreshold;
198
199 G4String theModelName;
200
201 std::pair<G4double, G4double> epCheckLevels;
202
203 std::vector<std::pair<G4double, const G4Material *> > theMinEnergyList;
204 std::vector<std::pair<G4double, const G4Material *> > theMaxEnergyList;
205 std::vector<std::pair<G4double, const G4Element *> > theMinEnergyListElements;
206 std::vector<std::pair<G4double, const G4Element *> > theMaxEnergyListElements;
207 std::vector<const G4Material *> theBlockedList;
208 std::vector<const G4Element *> theBlockedListElements;
209};
210
211#endif
double G4double
Definition: G4Types.hh:64
int G4int
Definition: G4Types.hh:66
bool G4bool
Definition: G4Types.hh:67
void DeActivateFor(const G4Material *aMaterial)
const G4HadronicInteraction * GetMyPointer() const
virtual const std::pair< G4double, G4double > GetFatalEnergyCheckLevels() const
void SetMinEnergy(G4double anEnergy)
void SetModelName(const G4String &nam)
virtual std::pair< G4double, G4double > GetEnergyMomentumCheckLevels() const
void SetVerboseLevel(G4int value)
virtual void ModelDescription(std::ostream &outFile) const
virtual G4bool IsApplicable(const G4HadProjectile &, G4Nucleus &)
void ActivateFor(const G4Element *anElement)
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)=0
void ActivateFor(const G4Material *aMaterial)
G4bool operator==(const G4HadronicInteraction &right) const
G4double GetRecoilEnergyThreshold() const
const G4String & GetModelName() const
void SetRecoilEnergyThreshold(G4double val)
G4bool operator!=(const G4HadronicInteraction &right) const
void SetEnergyMomentumCheckLevels(G4double relativeLevel, G4double absoluteLevel)
void SetMaxEnergy(const G4double anEnergy)
virtual G4double SampleInvariantT(const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)