Geant4 11.2.2
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4HadronicProcessStore.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// -------------------------------------------------------------------
28//
29// GEANT4 Class header file
30//
31//
32// File name: G4HadronicProcessStore
33//
34// Author: Vladimir Ivanchenko
35//
36// Creation date: 09.05.2008
37//
38// Modifications:
39//
40//
41// Class Description:
42//
43
44// -------------------------------------------------------------------
45//
46
47#ifndef G4HadronicProcessStore_h
48#define G4HadronicProcessStore_h 1
49
50
51#include "globals.hh"
52#include "G4DynamicParticle.hh"
53#include "G4ThreeVector.hh"
54#include "G4HadronicProcess.hh"
59#include <map>
60#include <vector>
61#include <iostream>
62
63class G4Element;
66
68{
69
71
72public:
73
75
77
78 void Clean();
80 const G4ParticleDefinition* particle,
81 G4double kineticEnergy,
82 const G4VProcess* process,
83 const G4Element* element,
84 const G4Material* material=nullptr);
85
87 const G4ParticleDefinition* particle,
88 G4double kineticEnergy,
89 const G4VProcess* process,
90 const G4Material* material);
91
93 const G4ParticleDefinition *aParticle,
94 G4double kineticEnergy,
95 const G4Material *material);
96
98 const G4ParticleDefinition *aParticle,
99 G4double kineticEnergy,
100 const G4Element *anElement,
101 const G4Material* mat=nullptr);
102
104 const G4ParticleDefinition *aParticle,
105 G4double kineticEnergy,
106 G4int Z, G4int A);
107
109 const G4ParticleDefinition *aParticle,
110 G4double kineticEnergy,
111 const G4Material *material);
112
114 const G4ParticleDefinition *aParticle,
115 G4double kineticEnergy,
116 const G4Element *anElement, const G4Material* mat=0);
117
119 const G4ParticleDefinition *aParticle,
120 G4double kineticEnergy,
121 G4int Z, G4int A);
122
124 const G4ParticleDefinition *aParticle,
125 G4double kineticEnergy,
126 const G4Material *material);
127
129 const G4ParticleDefinition *aParticle,
130 G4double kineticEnergy,
131 const G4Element *anElement,
132 const G4Material* mat=nullptr);
133
135 const G4ParticleDefinition *aParticle,
136 G4double kineticEnergy,
137 G4int Z, G4int A);
138
140 const G4ParticleDefinition *aParticle,
141 G4double kineticEnergy,
142 const G4Material *material);
143
145 const G4ParticleDefinition *aParticle,
146 G4double kineticEnergy,
147 const G4Element *anElement,
148 const G4Material* mat=nullptr);
149
151 const G4ParticleDefinition *aParticle,
152 G4double kineticEnergy,
153 G4int Z, G4int A);
154
156 const G4ParticleDefinition *aParticle,
157 G4double kineticEnergy,
158 const G4Material *material);
159
161 const G4ParticleDefinition *aParticle,
162 G4double kineticEnergy,
163 const G4Element *anElement,
164 const G4Material* mat=nullptr);
165
167 const G4ParticleDefinition *aParticle,
168 G4double kineticEnergy,
169 G4int Z, G4int A);
170
171 // register/deregister processes following G4HadronicProcess interface
173
175 const G4ParticleDefinition*);
176
179
181
182 // register/deregister processes following only G4VProcess interface
184
186 const G4ParticleDefinition*);
187
189
190 void SetBuildXSTable(G4bool val);
191
192 G4bool GetBuildXSTable() const;
193
194 void PrintInfo(const G4ParticleDefinition*);
195
196 void Dump(G4int level);
197 void DumpHtml();
198 void PrintHtml(const G4ParticleDefinition*, std::ofstream&);
199 void PrintModelHtml(const G4HadronicInteraction * model) const;
200
201 void SetVerbose(G4int val);
203 // these methods are obsolete and will be removed
204
206 G4HadronicProcessType subType);
207
208 // Energy-momentum non-conservation limits and reporting
209 void SetEpReportLevel(G4int level);
210
211 void SetProcessAbsLevel(G4double absoluteLevel);
212
213 void SetProcessRelLevel(G4double relativeLevel);
214
215private:
216
217 // constructor
219
220 // print process info
221 void Print(G4int idxProcess, G4int idxParticle);
222
223 G4String HtmlFileName(const G4String &) const;
224
225 static G4ThreadLocal G4HadronicProcessStore* instance;
226
227 typedef const G4ParticleDefinition* PD;
228 typedef G4HadronicProcess* HP;
229 typedef G4HadronicInteraction* HI;
230
231 // hadronic processes following G4HadronicProcess interface
232 std::vector<G4HadronicProcess*> process;
233 std::vector<G4HadronicInteraction*> model;
234 std::vector<G4String> modelName;
235 std::vector<PD> particle;
236 std::vector<G4int> wasPrinted;
237
238 std::multimap<PD,HP> p_map;
239 std::multimap<HP,HI> m_map;
240
241 // hadronic processes following only G4VProcess interface
242 std::vector<G4VProcess*> extraProcess;
243 std::multimap<PD,G4VProcess*> ep_map;
244
246
247 // counters and options
248 G4int n_proc = 0;
249 G4int n_model = 0;
250 G4int n_part = 0;
251 G4int n_extra = 0;
252
253 G4bool buildTableStart = true;
254 G4bool buildXSTable = false;
255
256 // cache
257 HP currentProcess = nullptr;
258 PD currentParticle = nullptr;
259 PD theGenericIon;
260
261 G4DynamicParticle localDP;
262
263 G4HadronicEPTestMessenger* theEPTestMessenger;
264};
265
266#endif
G4HadronicProcessType
double G4double
Definition G4Types.hh:83
bool G4bool
Definition G4Types.hh:86
int G4int
Definition G4Types.hh:85
const G4double A[17]
void DeRegister(G4HadronicProcess *)
G4double GetCaptureCrossSectionPerAtom(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4Element *anElement, const G4Material *mat=nullptr)
G4double GetCaptureCrossSectionPerIsotope(const G4ParticleDefinition *aParticle, G4double kineticEnergy, G4int Z, G4int A)
G4double GetCrossSectionPerVolume(const G4ParticleDefinition *particle, G4double kineticEnergy, const G4VProcess *process, const G4Material *material)
G4double GetCaptureCrossSectionPerVolume(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4Material *material)
G4HadronicProcess * FindProcess(const G4ParticleDefinition *, G4HadronicProcessType subType)
void RegisterParticle(G4HadronicProcess *, const G4ParticleDefinition *)
G4double GetChargeExchangeCrossSectionPerVolume(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4Material *material)
void PrintHtml(const G4ParticleDefinition *, std::ofstream &)
void SetProcessAbsLevel(G4double absoluteLevel)
G4double GetChargeExchangeCrossSectionPerIsotope(const G4ParticleDefinition *aParticle, G4double kineticEnergy, G4int Z, G4int A)
G4double GetFissionCrossSectionPerIsotope(const G4ParticleDefinition *aParticle, G4double kineticEnergy, G4int Z, G4int A)
G4double GetInelasticCrossSectionPerAtom(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4Element *anElement, const G4Material *mat=nullptr)
G4double GetInelasticCrossSectionPerIsotope(const G4ParticleDefinition *aParticle, G4double kineticEnergy, G4int Z, G4int A)
G4double GetInelasticCrossSectionPerVolume(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4Material *material)
void SetProcessRelLevel(G4double relativeLevel)
void DeRegisterExtraProcess(G4VProcess *)
G4double GetFissionCrossSectionPerAtom(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4Element *anElement, const G4Material *mat=nullptr)
G4double GetChargeExchangeCrossSectionPerAtom(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4Element *anElement, const G4Material *mat=nullptr)
void RegisterExtraProcess(G4VProcess *)
G4double GetElasticCrossSectionPerVolume(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4Material *material)
void RegisterParticleForExtraProcess(G4VProcess *, const G4ParticleDefinition *)
static G4HadronicProcessStore * Instance()
G4double GetElasticCrossSectionPerAtom(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4Element *anElement, const G4Material *mat=0)
void RegisterInteraction(G4HadronicProcess *, G4HadronicInteraction *)
G4double GetCrossSectionPerAtom(const G4ParticleDefinition *particle, G4double kineticEnergy, const G4VProcess *process, const G4Element *element, const G4Material *material=nullptr)
void Register(G4HadronicProcess *)
void PrintModelHtml(const G4HadronicInteraction *model) const
G4double GetElasticCrossSectionPerIsotope(const G4ParticleDefinition *aParticle, G4double kineticEnergy, G4int Z, G4int A)
void PrintInfo(const G4ParticleDefinition *)
G4double GetFissionCrossSectionPerVolume(const G4ParticleDefinition *aParticle, G4double kineticEnergy, const G4Material *material)
#define G4ThreadLocal
Definition tls.hh:77