Geant4  10.01.p03
G4HadronicInteraction.cc
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 // $Id: G4HadronicInteraction.cc 70714 2013-06-05 07:39:59Z gcosmo $
27 //
28 // Hadronic Interaction base class
29 // original by H.P. Wellisch
30 // modified by J.L. Chuma, TRIUMF, 21-Mar-1997
31 // Last modified: 04-Apr-1997
32 // reimplemented 1.11.2003 JPW.
33 // 23-Jan-2009 V.Ivanchenko move constructor and destructor to the body
34 
35 #include <iostream>
36 
37 #include "G4HadronicInteraction.hh"
38 #include "G4SystemOfUnits.hh"
40 #include "G4HadronicException.hh"
41 
43  verboseLevel(0), theMinEnergy(0.0), theMaxEnergy(25.0*GeV),
44  isBlocked(false), recoilEnergyThreshold(0.0), theModelName(modelName),
45  epCheckLevels(DBL_MAX, DBL_MAX)
46 {
48 }
49 
50 
52 {
54 }
55 
56 
57 G4double
60 {
61  return 0.0;
62 }
63 
65  const G4Material *aMaterial, const G4Element *anElement ) const
66 {
67  if( IsBlocked(aMaterial) ) { return 0.0; }
68  if( IsBlocked(anElement) ) { return 0.0; }
69  size_t length = theMinEnergyListElements.size();
70  if(0 < length) {
71  for(size_t i=0; i<length; ++i ) {
72  if( anElement == theMinEnergyListElements[i].second )
73  { return theMinEnergyListElements[i].first; }
74  }
75  }
76  length = theMinEnergyList.size();
77  if(0 < length) {
78  for(size_t i=0; i<length; ++i ) {
79  if( aMaterial == theMinEnergyList[i].second )
80  { return theMinEnergyList[i].first; }
81  }
82  }
83  if(IsBlocked()) { return 0.0; }
84  if( verboseLevel > 1 ) {
85  G4cout << "*** Warning from HadronicInteraction::GetMinEnergy" << G4endl
86  << " material " << aMaterial->GetName()
87  << " not found in min energy List" << G4endl;
88  }
89  return theMinEnergy;
90 }
91 
93  const G4Element *anElement )
94 {
95  if( IsBlocked(anElement) ) {
96  G4cout << "*** Warning from HadronicInteraction::SetMinEnergy" << G4endl
97  << " The model is not active for the Element "
98  << anElement->GetName() << "." << G4endl;
99  }
100  size_t length = theMinEnergyListElements.size();
101  if(0 < length) {
102  for(size_t i=0; i<length; ++i ) {
103  if( anElement == theMinEnergyListElements[i].second )
104  {
105  theMinEnergyListElements[i].first = anEnergy;
106  return;
107  }
108  }
109  }
110  theMinEnergyListElements.push_back(std::pair<G4double, const G4Element *>(anEnergy, anElement));
111 }
112 
114  const G4Material *aMaterial )
115 {
116  if( IsBlocked(aMaterial) ) {
117  G4cout << "*** Warning from HadronicInteraction::SetMinEnergy" << G4endl
118  << " The model is not active for the Material "
119  << aMaterial->GetName() << "." << G4endl;
120  }
121  size_t length = theMinEnergyList.size();
122  if(0 < length) {
123  for(size_t i=0; i<length; ++i ) {
124  if( aMaterial == theMinEnergyList[i].second )
125  {
126  theMinEnergyList[i].first = anEnergy;
127  return;
128  }
129  }
130  }
131  theMinEnergyList.push_back(std::pair<G4double, const G4Material *>(anEnergy, aMaterial));
132 }
133 
135  const G4Element *anElement ) const
136 {
137  if( IsBlocked(aMaterial) ) { return 0.0; }
138  if( IsBlocked(anElement) ) { return 0.0; }
139  size_t length = theMaxEnergyListElements.size();
140  if(0 < length) {
141  for(size_t i=0; i<length; ++i ) {
142  if( anElement == theMaxEnergyListElements[i].second )
143  { return theMaxEnergyListElements[i].first; }
144  }
145  }
146  length = theMaxEnergyList.size();
147  if(0 < length) {
148  for(size_t i=0; i<length; ++i ) {
149  if( aMaterial == theMaxEnergyList[i].second )
150  { return theMaxEnergyList[i].first; }
151  }
152  }
153  if(IsBlocked()) { return 0.0; }
154  if( verboseLevel > 1 ) {
155  G4cout << "*** Warning from HadronicInteraction::GetMaxEnergy" << G4endl
156  << " material " << aMaterial->GetName()
157  << " not found in min energy List" << G4endl;
158  }
159  return theMaxEnergy;
160 }
161 
163  const G4Element *anElement )
164 {
165  if( IsBlocked(anElement) ) {
166  G4cout << "*** Warning from HadronicInteraction::SetMaxEnergy" << G4endl
167  << "Warning: The model is not active for the Element "
168  << anElement->GetName() << "." << G4endl;
169  }
170  size_t length = theMaxEnergyListElements.size();
171  if(0 < length) {
172  for(size_t i=0; i<length; ++i ) {
173  if( anElement == theMaxEnergyListElements[i].second )
174  {
175  theMaxEnergyListElements[i].first = anEnergy;
176  return;
177  }
178  }
179  }
180  theMaxEnergyListElements.push_back(std::pair<G4double, const G4Element *>(anEnergy, anElement));
181 }
182 
184  const G4Material *aMaterial )
185 {
186  if( IsBlocked(aMaterial) ) {
187  G4cout << "*** Warning from HadronicInteraction::SetMaxEnergy" << G4endl
188  << "Warning: The model is not active for the Material "
189  << aMaterial->GetName() << "." << G4endl;
190  }
191  size_t length = theMaxEnergyList.size();
192  if(0 < length) {
193  for(size_t i=0; i<length; ++i ) {
194  if( aMaterial == theMaxEnergyList[i].second )
195  {
196  theMaxEnergyList[i].first = anEnergy;
197  return;
198  }
199  }
200  }
201  theMaxEnergyList.push_back(std::pair<G4double, const G4Material *>(anEnergy, aMaterial));
202 }
203 
205 {
206  theBlockedList.push_back(aMaterial);
207 }
208 
210 {
211  theBlockedListElements.push_back(anElement);
212 }
213 
214 
216 {
217  for (size_t i=0; i<theBlockedList.size(); ++i) {
218  if (aMaterial == theBlockedList[i]) return true;
219  }
220  return false;
221 }
222 
223 
225 {
226  for (size_t i=0; i<theBlockedListElements.size(); ++i) {
227  if (anElement == theBlockedListElements[i]) return true;
228  }
229  return false;
230 }
231 
232 const std::pair<G4double, G4double> G4HadronicInteraction::GetFatalEnergyCheckLevels() const
233 {
234  // default level of Check
235  return std::pair<G4double, G4double>(2.*perCent, 1. * GeV);
236 }
237 
238 std::pair<G4double, G4double>
240 {
241  return epCheckLevels;
242 }
243 
244 
245 void G4HadronicInteraction::ModelDescription(std::ostream& outFile) const
246 {
247  outFile << "The description for this model has not been written yet.\n";
248 }
249 
250 /*
251 G4HadronicInteraction::G4HadronicInteraction(const G4HadronicInteraction &right )
252 {
253  *this = right;
254 }
255 
256 const G4HadronicInteraction&
257 G4HadronicInteraction::operator=(const G4HadronicInteraction &right )
258 {
259  G4String text = "unintended use of G4HadronicInteraction::operator=";
260  throw G4HadronicException(__FILE__, __LINE__, text);
261  return right;
262 }
263  */
264 /* end of file */
265 
std::vector< std::pair< G4double, const G4Element * > > theMaxEnergyListElements
G4double GetMinEnergy() const
G4bool isBlocked(ParticleList const &p, Nucleus const *const n)
Check Pauli blocking.
std::vector< const G4Element * > theBlockedListElements
virtual G4double SampleInvariantT(const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
std::vector< std::pair< G4double, const G4Material * > > theMaxEnergyList
G4double GetMaxEnergy() const
virtual const std::pair< G4double, G4double > GetFatalEnergyCheckLevels() const
const G4String & GetName() const
Definition: G4Material.hh:178
virtual void ModelDescription(std::ostream &outFile) const
void RemoveMe(G4HadronicInteraction *aModel)
void RegisterMe(G4HadronicInteraction *aModel)
int G4int
Definition: G4Types.hh:78
std::vector< std::pair< G4double, const G4Material * > > theMinEnergyList
void SetMinEnergy(G4double anEnergy)
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
virtual std::pair< G4double, G4double > GetEnergyMomentumCheckLevels() const
std::vector< std::pair< G4double, const G4Element * > > theMinEnergyListElements
static const double second
Definition: G4SIunits.hh:138
static const double GeV
Definition: G4SIunits.hh:196
static const double perCent
Definition: G4SIunits.hh:296
G4HadronicInteraction(const G4String &modelName="HadronicModel")
std::pair< G4double, G4double > epCheckLevels
static G4HadronicInteractionRegistry * Instance()
void SetMaxEnergy(const G4double anEnergy)
std::vector< const G4Material * > theBlockedList
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
const G4String & GetName() const
Definition: G4Element.hh:127
#define DBL_MAX
Definition: templates.hh:83
void DeActivateFor(const G4Material *aMaterial)