Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4LENDCrossSection.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 // Class Description
27 // Cross Section for LEND (Low Energy Nuclear Data)
28 // LEND is Geant4 interface for GIDI (General Interaction Data Interface)
29 // which gives a discription of nuclear and atomic reactions, such as
30 // Binary collision cross sections
31 // Particle number multiplicity distributions of reaction products
32 // Energy and angular distributions of reaction products
33 // Derived calculational constants
34 // GIDI is developped at Lawrence Livermore National Laboratory
35 // Class Description - End
36 
37 // 071025 First implementation done by T. Koi (SLAC/SCCS)
38 // 101118 Name modifications for release T. Koi (SLAC/PPA)
39 
40 #include "G4LENDCrossSection.hh"
41 #include "G4Pow.hh"
42 #include "G4SystemOfUnits.hh"
43 #include "G4ElementTable.hh"
44 #include "G4HadronicException.hh"
45 
46 //TK110811
48  const G4Element* /*element*/ , const G4Material* /*material*/ )
49 {
50  G4double eKin = dp->GetKineticEnergy();
51  if ( dp->GetDefinition() != proj ) return false;
52  if ( eKin > GetMaxKinEnergy() || eKin < GetMinKinEnergy() ) return false;
53 
54  return true;
55 }
56 
58  const G4Isotope* isotope , const G4Element* /*elment*/ , const G4Material* material )
59 {
60 
61  G4double xs = 0.0;
62  G4double ke = dp->GetKineticEnergy();
63  G4double temp = material->GetTemperature();
64  G4int iM = isotope->Getm();
65 
66  G4GIDI_target* aTarget = usedTarget_map.find( lend_manager->GetNucleusEncoding( iZ , iA , iM ) )->second->GetTarget();
67 
68  xs = getLENDCrossSection ( aTarget , ke , temp );
69 
70  return xs;
71 }
72 
73 
74 /*
75 G4bool G4LENDCrossSection::IsApplicable(const G4DynamicParticle*aP, const G4Element*)
76 {
77  G4bool result = true;
78  G4double eKin = aP->GetKineticEnergy();
79  if( eKin > GetMaxKinEnergy() || aP->GetDefinition() != proj ) result = false;
80  return result;
81 }
82 */
83 
86 {
87 
88  proj = NULL; //will be set in an inherited class
89  //default_evaluation = "endl99";
90  //default_evaluation = "ENDF.B-VII.0";
91  default_evaluation = "ENDF/BVII.1";
92 
93  allow_nat = false;
94  allow_any = false;
95 
96  SetMinKinEnergy( 0*MeV );
97  SetMaxKinEnergy( 20*MeV );
98 
99  lend_manager = G4LENDManager::GetInstance();
100 
101 }
102 
104 {
105 
106  for ( std::map< G4int , G4LENDUsedTarget* >::iterator
107  it = usedTarget_map.begin() ; it != usedTarget_map.end() ; it ++ )
108  {
109  delete it->second;
110  }
111 
112 }
113 
115 {
117 }
118 
120 {
121 
122  if ( &aP != proj )
123  throw G4HadronicException(__FILE__, __LINE__, "Attempt to use LEND data for particles other than neutrons!!!");
124 
125  G4cout << G4endl;
126  G4cout << "Dump Cross Sections of " << GetName() << G4endl;
127  G4cout << "(Pointwise cross-section at 300 Kelvin.)" << G4endl;
128  G4cout << G4endl;
129 
130  G4cout << "Target informaiton " << G4endl;
131 
132  for ( std::map< G4int , G4LENDUsedTarget* >::iterator
133  it = usedTarget_map.begin() ; it != usedTarget_map.end() ; it ++ )
134  {
135  G4cout
136  << "Wanted " << it->second->GetWantedEvaluation()
137  << ", Z= " << it->second->GetWantedZ()
138  << ", A= " << it->second->GetWantedA()
139  << "; Actual " << it->second->GetActualEvaluation()
140  << ", Z= " << it->second->GetActualZ()
141  << ", A= " << it->second->GetActualA()
142  << ", " << it->second->GetTarget()
143  << G4endl;
144 
145  G4int ie = 0;
146 
147  G4GIDI_target* aTarget = it->second->GetTarget();
148  G4double aT = 300;
149  for ( ie = 0 ; ie < 130 ; ie++ )
150  {
151  G4double ke = 1.0e-5 * G4Pow::GetInstance()->powA ( 10.0 , ie/10.0 ) *eV;
152 
153  if ( ke < 20*MeV )
154  {
155  G4cout << " " << GetName() << ", cross section at " << ke/eV << " [eV] = " << getLENDCrossSection ( aTarget , ke , aT )/barn << " [barn] " << G4endl;
156  }
157  }
158  G4cout << G4endl;
159 
160  }
161 
162 }
163 
164 
165 /*
166 //110810
167 //G4double G4LENDCrossSection::GetCrossSection(const G4DynamicParticle* aP , const G4Element* anElement , G4double aT)
168 G4double G4LENDCrossSection::GetCrossSection(const G4DynamicParticle* aP , int iZ , const G4Material* aMat)
169 {
170 
171 //110810
172  G4double aT = aMat->GetTemperature();
173  G4Element* anElement = lend_manager->GetNistElementBuilder()->FindOrBuildElement( iZ );
174 
175  G4double ke = aP->GetKineticEnergy();
176  G4double XS = 0.0;
177 
178  G4int numberOfIsotope = anElement->GetNumberOfIsotopes();
179 
180  if ( numberOfIsotope > 0 )
181  {
182  // User Defined Abundances
183  for ( G4int i_iso = 0 ; i_iso < numberOfIsotope ; i_iso++ )
184  {
185 
186  G4int iZ = anElement->GetIsotope( i_iso )->GetZ();
187  G4int iA = anElement->GetIsotope( i_iso )->GetN();
188  G4double ratio = anElement->GetRelativeAbundanceVector()[i_iso];
189 
190  G4GIDI_target* aTarget = usedTarget_map.find( lend_manager->GetNucleusEncoding( iZ , iA ) )->second->GetTarget();
191  XS += ratio*getLENDCrossSection ( aTarget , ke , aT );
192 
193  }
194  }
195  else
196  {
197  // Natural Abundances
198  G4NistElementBuilder* nistElementBuild = lend_manager->GetNistElementBuilder();
199  G4int iZ = int ( anElement->GetZ() );
200  G4int numberOfNistIso = nistElementBuild->GetNumberOfNistIsotopes( int ( anElement->GetZ() ) );
201 
202  G4int Nfirst = nistElementBuild->GetNistFirstIsotopeN( iZ );
203  for ( G4int i = 0 ; i < numberOfNistIso ; i++ )
204  {
205  G4int iA = Nfirst + i;
206  G4double ratio = nistElementBuild->GetIsotopeAbundance( iZ , iA );
207  if ( ratio > 0.0 )
208  {
209  G4GIDI_target* aTarget = usedTarget_map.find( lend_manager->GetNucleusEncoding( iZ , iA ) )->second->GetTarget();
210  XS += ratio*getLENDCrossSection ( aTarget , ke , aT );
211  //G4cout << ke/eV << " " << iZ << " " << iMass << " " << aTarget << " " << getLENDCrossSection ( aTarget , ke , aT ) << G4endl;
212  }
213  }
214  }
215 
216  //G4cout << "XS= " << XS << G4endl;
217  return XS;
218 }
219 
220 
221 
222 //110810
223 //G4double G4LENDCrossSection::GetIsoCrossSection(const G4DynamicParticle* dp, const G4Isotope* isotope, G4double aT )
224 G4double G4LENDCrossSection::GetIsoCrossSection(const G4DynamicParticle* dp, const G4Isotope* isotope, const G4Material* aMat)
225 {
226 
227 //110810
228  G4double aT = aMat->GetTemperature();
229 
230  G4double ke = dp->GetKineticEnergy();
231 
232  G4int iZ = isotope->GetZ();
233  G4int iA = isotope->GetN();
234 
235  G4GIDI_target* aTarget = usedTarget_map.find( lend_manager->GetNucleusEncoding( iZ , iA ) )->second->GetTarget();
236 
237  return getLENDCrossSection ( aTarget , ke , aT );
238 
239 }
240 
241 
242 
243 //110810
244 //G4double G4LENDCrossSection::GetZandACrossSection(const G4DynamicParticle* dp, G4int iZ, G4int iA, G4double aT)
245 G4double G4LENDCrossSection::GetZandACrossSection(const G4DynamicParticle* dp, G4int iZ, G4int iA, const G4Material* aMat)
246 {
247 
248 //110810
249  G4double aT = aMat->GetTemperature();
250 
251  G4double ke = dp->GetKineticEnergy();
252 
253  G4GIDI_target* aTarget = usedTarget_map.find( lend_manager->GetNucleusEncoding( iZ , iA ) )->second->GetTarget();
254 
255  return getLENDCrossSection ( aTarget , ke , aT );
256 
257 }
258 */
259 
260 
261 
262 void G4LENDCrossSection::recreate_used_target_map()
263 {
264  for ( std::map< G4int , G4LENDUsedTarget* >::iterator
265  it = usedTarget_map.begin() ; it != usedTarget_map.end() ; it ++ )
266  {
267  delete it->second;
268  }
269  usedTarget_map.clear();
270 
272 }
273 
274 
275 
277 {
278 
279  lend_manager->RequestChangeOfVerboseLevel( verboseLevel );
280 
281  size_t numberOfElements = G4Element::GetNumberOfElements();
282  static const G4ElementTable* theElementTable = G4Element::GetElementTable();
283 
284  for ( size_t i = 0 ; i < numberOfElements ; ++i )
285  {
286 
287  const G4Element* anElement = (*theElementTable)[i];
288  G4int numberOfIsotope = anElement->GetNumberOfIsotopes();
289 
290  if ( numberOfIsotope > 0 )
291  {
292  // User Defined Abundances
293  for ( G4int i_iso = 0 ; i_iso < numberOfIsotope ; i_iso++ )
294  {
295  G4int iZ = anElement->GetIsotope( i_iso )->GetZ();
296  G4int iA = anElement->GetIsotope( i_iso )->GetN();
297  G4int iIsomer = anElement->GetIsotope( i_iso )->Getm();
298 
299  //G4LENDUsedTarget* aTarget = new G4LENDUsedTarget ( G4Neutron::Neutron() , default_evaluation , iZ , iA );
300  G4LENDUsedTarget* aTarget = new G4LENDUsedTarget ( proj , default_evaluation , iZ , iA , iIsomer );
301  if ( allow_nat == true ) aTarget->AllowNat();
302  if ( allow_any == true ) aTarget->AllowAny();
303  usedTarget_map.insert( std::pair< G4int , G4LENDUsedTarget* > ( lend_manager->GetNucleusEncoding( iZ , iA , iIsomer ) , aTarget ) );
304  }
305  }
306  else
307  {
308  // Natural Abundances
309  G4NistElementBuilder* nistElementBuild = lend_manager->GetNistElementBuilder();
310  G4int iZ = int ( anElement->GetZ() );
311  //G4cout << nistElementBuild->GetNumberOfNistIsotopes( int ( anElement->GetZ() ) ) << G4endl;
312  G4int numberOfNistIso = nistElementBuild->GetNumberOfNistIsotopes( int ( anElement->GetZ() ) );
313 
314  for ( G4int ii = 0 ; ii < numberOfNistIso ; ii++ )
315  {
316  //G4cout << nistElementBuild->GetIsotopeAbundance( iZ , nistElementBuild->GetNistFirstIsotopeN( iZ ) + i ) << G4endl;
317  if ( nistElementBuild->GetIsotopeAbundance( iZ , nistElementBuild->GetNistFirstIsotopeN( iZ ) + ii ) > 0 )
318  {
319  G4int iMass = nistElementBuild->GetNistFirstIsotopeN( iZ ) + ii;
320  //G4cout << iZ << " " << nistElementBuild->GetNistFirstIsotopeN( iZ ) + i << " " << nistElementBuild->GetIsotopeAbundance ( iZ , iMass ) << G4endl;
321  G4int iIsomer = 0;
322 
323  G4LENDUsedTarget* aTarget = new G4LENDUsedTarget ( proj , default_evaluation , iZ , iMass );
324  if ( allow_nat == true ) aTarget->AllowNat();
325  if ( allow_any == true ) aTarget->AllowAny();
326  usedTarget_map.insert( std::pair< G4int , G4LENDUsedTarget* > ( lend_manager->GetNucleusEncoding( iZ , iMass , iIsomer ) , aTarget ) );
327 
328  }
329 
330  }
331  }
332  }
333 
334  G4cout << "Dump UsedTarget for " << GetName() << G4endl;
335  //G4cout << "Requested Evaluation, Z , A -> Actual Evaluation, Z , A(0=Nat) , Pointer of Target" << G4endl;
336  G4cout << "Requested Evaluation, Z , A -> Actual Evaluation, Z , A(0=Nat) " << G4endl;
337  for ( std::map< G4int , G4LENDUsedTarget* >::iterator
338  it = usedTarget_map.begin() ; it != usedTarget_map.end() ; it ++ )
339  {
340  G4cout
341  << " " << it->second->GetWantedEvaluation()
342  << ", " << it->second->GetWantedZ()
343  << ", " << it->second->GetWantedA()
344  << " -> " << it->second->GetActualEvaluation()
345  << ", " << it->second->GetActualZ()
346  << ", " << it->second->GetActualA()
347  //<< ", " << it->second->GetTarget()
348  << G4endl;
349  }
350 
351 }
352 
353  // elow ehigh xs_elow xs_ehigh ke (ke < elow)
355 {
356  //XS propotinal to 1/v at low energy -> 1/root(E)
357  //XS = a * 1/root(E) + b
358  G4double a = ( y2 - y1 ) / ( 1/std::sqrt(x2) - 1/std::sqrt(x1) );
359  G4double b = y1 - a * 1/std::sqrt(x1);
360  G4double result = a * 1/std::sqrt(ke) + b;
361  return result;
362 }
G4double G4ParticleHPJENDLHEData::G4double result
size_t GetNumberOfIsotopes() const
Definition: G4Element.hh:159
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
virtual G4double getLENDCrossSection(G4GIDI_target *, G4double, G4double)
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
void DumpPhysicsTable(const G4ParticleDefinition &)
G4double GetUltraLowEnergyExtrapolatedXS(G4double, G4double, G4double, G4double, G4double)
G4double GetKineticEnergy() const
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
G4NistElementBuilder * GetNistElementBuilder()
G4double GetZ() const
Definition: G4Element.hh:131
G4double GetIsotopeAbundance(G4int Z, G4int N) const
G4int GetNumberOfNistIsotopes(G4int Z) const
G4ParticleDefinition * GetDefinition() const
G4ParticleDefinition * proj
static constexpr double second
Definition: G4SIunits.hh:157
const G4String & GetName() const
int G4int
Definition: G4Types.hh:78
G4LENDCrossSection(const G4String name="")
G4bool IsIsoApplicable(const G4DynamicParticle *, G4int, G4int, const G4Element *, const G4Material *)
G4int GetNucleusEncoding(G4int iZ, G4int iA, G4int iM)
string material
Definition: eplot.py:19
tuple b
Definition: test.py:12
G4int GetN() const
Definition: G4Isotope.hh:94
G4GLOB_DLL std::ostream G4cout
static size_t GetNumberOfElements()
Definition: G4Element.cc:405
G4int Getm() const
Definition: G4Isotope.hh:100
void SetMinKinEnergy(G4double value)
bool G4bool
Definition: G4Types.hh:79
G4int GetNistFirstIsotopeN(G4int Z) const
typedef int(XMLCALL *XML_NotStandaloneHandler)(void *userData)
static constexpr double eV
Definition: G4SIunits.hh:215
G4int GetZ() const
Definition: G4Isotope.hh:91
static G4LENDManager * GetInstance()
void SetMaxKinEnergy(G4double value)
const G4Isotope * GetIsotope(G4int iso) const
Definition: G4Element.hh:170
G4double GetTemperature() const
Definition: G4Material.hh:182
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
std::vector< G4Element * > G4ElementTable
double G4double
Definition: G4Types.hh:76
G4bool RequestChangeOfVerboseLevel(G4int)
static constexpr double barn
Definition: G4SIunits.hh:105
static G4ElementTable * GetElementTable()
Definition: G4Element.cc:398
G4double GetIsoCrossSection(const G4DynamicParticle *, G4int, G4int, const G4Isotope *, const G4Element *, const G4Material *)
void BuildPhysicsTable(const G4ParticleDefinition &)