Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4EnergyRangeManager.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: G4EnergyRangeManager.cc 98067 2016-07-01 16:33:54Z gcosmo $
27 //
28  // Hadronic Process: Energy Range Manager
29  // original by H.P. Wellisch
30  // modified by J.L. Chuma, TRIUMF, 22-Nov-1996
31  // Last modified: 24-Mar-1997
32  // fix in the counter-hndling: H.P. Wellisch 04-Apr-97
33  // throw an exception if no model found: J.L. Chuma 04-Apr-97
34 
35 #include "G4EnergyRangeManager.hh"
36 #include "Randomize.hh"
37 #include "G4HadronicException.hh"
38 
40  : theHadronicInteractionCounter(0)
41 {}
42 
44 {}
45 
47 {
48  theHadronicInteractionCounter = right.theHadronicInteractionCounter;
49  theHadronicInteraction = right.theHadronicInteraction;
50 }
51 
54 {
55  if (this != &right) {
56  theHadronicInteractionCounter = right.theHadronicInteractionCounter;
57  theHadronicInteraction = right.theHadronicInteraction;
58  }
59  return *this;
60 }
61 
63 {
64  if(!a) { return; }
65  if(0 < theHadronicInteractionCounter) {
66  for(G4int i=0; i<theHadronicInteractionCounter; ++i) {
67  if(a == theHadronicInteraction[i]) { return; }
68  }
69  }
70  theHadronicInteraction.push_back(a);
71  ++theHadronicInteractionCounter;
72 }
73 
74 
77  G4Nucleus & aTargetNucleus,
78  const G4Material* aMaterial,
79  const G4Element* anElement) const
80 {
81  if(0 == theHadronicInteractionCounter) {
82  throw G4HadronicException(__FILE__, __LINE__,
83  "GetHadronicInteraction: NO MODELS STORED");
84  }
85 
86  G4double kineticEnergy = aHadProjectile.GetKineticEnergy();
87  // For ions, get kinetic energy per nucleon
88  if ( aHadProjectile.GetDefinition()->GetBaryonNumber() > 1.5 ) {
89  kineticEnergy /= aHadProjectile.GetDefinition()->GetBaryonNumber();
90  }
91 
92  G4int cou = 0, memory = 0, memor2 = 0;
93  G4double emi1 = 0.0, ema1 = 0.0, emi2 = 0.0, ema2 = 0.0;
94 
95  for (G4int i = 0; i<theHadronicInteractionCounter; ++i) {
96  if ( theHadronicInteraction[i]->IsApplicable( aHadProjectile, aTargetNucleus ) ) {
97  G4double low = theHadronicInteraction[i]->GetMinEnergy( aMaterial, anElement );
98  // Work-around for particles with 0 kinetic energy, which still
99  // require a model to return a ParticleChange
100  //if (low == 0.) low = -DBL_MIN;
101  G4double high = theHadronicInteraction[i]->GetMaxEnergy( aMaterial, anElement );
102  if (low <= kineticEnergy && high > kineticEnergy) {
103  ++cou;
104  emi2 = emi1;
105  ema2 = ema1;
106  emi1 = low;
107  ema1 = high;
108  memor2 = memory;
109  memory = i;
110  }
111  }
112  }
113 
114  G4int mem = -1;
115  G4double rand;
116  switch (cou) {
117  case 0:
118  G4cout<<"G4EnergyRangeManager:GetHadronicInteraction: counter="
119  <<theHadronicInteractionCounter<<", Ek="
120  <<kineticEnergy<<", Material = "<<aMaterial->GetName()
121  <<", Element = "
122  <<anElement->GetName()<<G4endl;
123  for( G4int j=0; j<theHadronicInteractionCounter; ++j)
124  {
125  G4HadronicInteraction* HInt=theHadronicInteraction[j];
126  G4cout<<"*"<<j<<"* low=" <<HInt->GetMinEnergy(aMaterial,anElement)
127  <<", high="<<HInt->GetMaxEnergy(aMaterial,anElement)<<G4endl;
128  }
129  throw G4HadronicException(__FILE__, __LINE__,
130  "GetHadronicInteraction: No Model found");
131  return 0;
132  case 1:
133  mem = memory;
134  break;
135  case 2:
136  if( (emi2<=emi1 && ema2>=ema1) || (emi2>=emi1 && ema2<=ema1) )
137  {
138  G4cout<<"G4EnergyRangeManager:GetHadronicInteraction: counter="
139  <<theHadronicInteractionCounter<<", Ek="
140  <<kineticEnergy<<", Material = "<<aMaterial->GetName()
141  <<", Element = "
142  <<anElement->GetName()<<G4endl;
143  for( G4int j=0; j<theHadronicInteractionCounter; ++j)
144  {
145  G4HadronicInteraction* HInt=theHadronicInteraction[j];
146  G4cout<<"*"<<j<<"* low=" <<HInt->GetMinEnergy(aMaterial,anElement)
147  <<", high="<<HInt->GetMaxEnergy(aMaterial,anElement)<<G4endl;
148  }
149  throw G4HadronicException(__FILE__, __LINE__,
150  "GetHadronicInteraction: Energy ranges of two models fully overlapping");
151  }
152  rand = G4UniformRand();
153  if( emi1 < emi2 )
154  {
155  if( (ema1-kineticEnergy) < rand*(ema1-emi2) ) {
156  mem = memor2;
157  } else {
158  mem = memory;
159  }
160  } else {
161  if( (ema2-kineticEnergy) < rand*(ema2-emi1) ) {
162  mem = memory;
163  } else {
164  mem = memor2;
165  }
166  }
167  break;
168  default:
169  throw G4HadronicException(__FILE__, __LINE__,
170  "GetHadronicInteraction: More than two competing models in this energy range");
171  }
172 
173  return theHadronicInteraction[mem];
174 }
175 
176 
179  const G4Material* aMaterial,
180  const G4Element* anElement) const
181 {
182  if(0 == theHadronicInteractionCounter) {
183  throw G4HadronicException(__FILE__, __LINE__,
184  "GetHadronicInteraction: NO MODELS STORED");
185  }
186  G4int cou = 0, memory = 0, memor2 = 0;
187  G4double emi1 = 0.0, ema1 = 0.0, emi2 = 0.0, ema2 = 0.0;
188 
189  for (G4int i = 0; i<theHadronicInteractionCounter; ++i) {
190  G4double low = theHadronicInteraction[i]->GetMinEnergy( aMaterial, anElement );
191  // Work-around for particles with 0 kinetic energy, which still
192  // require a model to return a ParticleChange
193  //if (low == 0.) low = -DBL_MIN;
194  G4double high = theHadronicInteraction[i]->GetMaxEnergy( aMaterial, anElement );
195  if (low <= kineticEnergy && high > kineticEnergy) {
196  ++cou;
197  emi2 = emi1;
198  ema2 = ema1;
199  emi1 = low;
200  ema1 = high;
201  memor2 = memory;
202  memory = i;
203  }
204  }
205 
206  G4int mem = -1;
207  G4double rand;
208  switch (cou) {
209  case 0:
210  G4cout<<"G4EnergyRangeManager:GetHadronicInteraction: counter="
211  <<theHadronicInteractionCounter<<", Ek="
212  <<kineticEnergy<<", Material = "<<aMaterial->GetName()
213  <<", Element = "
214  <<anElement->GetName()<<G4endl;
215  for( G4int j=0; j<theHadronicInteractionCounter; ++j)
216  {
217  G4HadronicInteraction* HInt=theHadronicInteraction[j];
218  G4cout<<"*"<<j<<"* low=" <<HInt->GetMinEnergy(aMaterial,anElement)
219  <<", high="<<HInt->GetMaxEnergy(aMaterial,anElement)<<G4endl;
220  }
221  throw G4HadronicException(__FILE__, __LINE__,
222  "GetHadronicInteraction: No Model found");
223  return 0;
224  case 1:
225  mem = memory;
226  break;
227  case 2:
228  if( (emi2<=emi1 && ema2>=ema1) || (emi2>=emi1 && ema2<=ema1) )
229  {
230  G4cout<<"G4EnergyRangeManager:GetHadronicInteraction: counter="
231  <<theHadronicInteractionCounter<<", Ek="
232  <<kineticEnergy<<", Material = "<<aMaterial->GetName()
233  <<", Element = "
234  <<anElement->GetName()<<G4endl;
235  for( G4int j=0; j<theHadronicInteractionCounter; ++j)
236  {
237  G4HadronicInteraction* HInt=theHadronicInteraction[j];
238  G4cout<<"*"<<j<<"* low=" <<HInt->GetMinEnergy(aMaterial,anElement)
239  <<", high="<<HInt->GetMaxEnergy(aMaterial,anElement)<<G4endl;
240  }
241  throw G4HadronicException(__FILE__, __LINE__,
242  "GetHadronicInteraction: Energy ranges of two models fully overlapping");
243  }
244  rand = G4UniformRand();
245  if( emi1 < emi2 )
246  {
247  if( (ema1-kineticEnergy) < rand*(ema1-emi2) ) {
248  mem = memor2;
249  } else {
250  mem = memory;
251  }
252  } else {
253  if( (ema2-kineticEnergy) < rand*(ema2-emi1) ) {
254  mem = memory;
255  } else {
256  mem = memor2;
257  }
258  }
259  break;
260  default:
261  throw G4HadronicException(__FILE__, __LINE__,
262  "GetHadronicInteraction: More than two competing models in this energy range");
263  }
264 
265  return theHadronicInteraction[mem];
266 }
267 
268 std::vector<G4HadronicInteraction*>&
270 {
271  return theHadronicInteraction;
272 }
273 
274 #include "G4SystemOfUnits.hh"
276 {
277  G4cout << "G4EnergyRangeManager " << this << G4endl;
278  for (G4int i = 0 ; i < theHadronicInteractionCounter; i++) {
279  G4cout << " HadronicModel " << i <<":"
280  << theHadronicInteraction[i]->GetModelName() << G4endl;
281  if (verbose > 0) {
282  G4cout << " Minimum Energy "
283  << theHadronicInteraction[i]->GetMinEnergy()/GeV << " [GeV], "
284  << "Maximum Energy "
285  << theHadronicInteraction[i]->GetMaxEnergy()/GeV << " [GeV]"
286  << G4endl;
287  }
288  }
289 }
290 
291 void
293 {
294  for ( std::vector<G4HadronicInteraction*>::iterator
295  it = theHadronicInteraction.begin() ; it != theHadronicInteraction.end() ; it++ ) {
296  (*it)->BuildPhysicsTable( aParticleType );
297  }
298 }
299  /* end of file */
300 
G4double GetMinEnergy() const
void RegisterMe(G4HadronicInteraction *a)
G4double GetMaxEnergy() const
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
void Dump(G4int verbose=0)
G4HadronicInteraction * GetHadronicInteraction(const G4HadProjectile &aHadProjectile, G4Nucleus &aTargetNucleus, const G4Material *aMaterial, const G4Element *anElement) const
const G4String & GetName() const
Definition: G4Material.hh:178
int G4int
Definition: G4Types.hh:78
std::vector< G4HadronicInteraction * > & GetHadronicInteractionList()
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
const G4ParticleDefinition * GetDefinition() const
G4EnergyRangeManager & operator=(const G4EnergyRangeManager &right)
G4double GetKineticEnergy() const
void BuildPhysicsTable(const G4ParticleDefinition &)
static constexpr double GeV
Definition: G4SIunits.hh:217
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
const G4String & GetName() const
Definition: G4Element.hh:127