Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4ParticleHPThermalScatteringData.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 // Thermal Neutron Scattering
27 // Koi, Tatsumi (SCCS/SLAC)
28 //
29 // Class Description
30 // Cross Sections for a high precision (based on evaluated data
31 // libraries) description of themal neutron scattering below 4 eV;
32 // Based on Thermal neutron scattering files
33 // from the evaluated nuclear data files ENDF/B-VI, Release2
34 // To be used in your physics list in case you need this physics.
35 // In this case you want to register an object of this class with
36 // the corresponding process.
37 // Class Description - End
38 
39 // 15-Nov-06 First implementation is done by T. Koi (SLAC/SCCS)
40 // 070625 implement clearCurrentXSData to fix memory leaking by T. Koi
41 // P. Arce, June-2014 Conversion neutron_hp to particle_hp
42 //
43 
44 #include <list>
45 #include <algorithm>
46 
48 #include "G4ParticleHPManager.hh"
49 
50 #include "G4SystemOfUnits.hh"
51 #include "G4Neutron.hh"
52 #include "G4ElementTable.hh"
53 
54 #include "G4Threading.hh"
55 
57 :G4VCrossSectionDataSet("NeutronHPThermalScatteringData")
58 ,coherent(NULL)
59 ,incoherent(NULL)
60 ,inelastic(NULL)
61 {
62 // Upper limit of neutron energy
63  emax = 4*eV;
64  SetMinKinEnergy( 0*MeV );
65  SetMaxKinEnergy( emax );
66 
67  ke_cache = 0.0;
68  xs_cache = 0.0;
69  element_cache = NULL;
70  material_cache = NULL;
71 
72  indexOfThermalElement.clear();
73 
75 }
76 
78 {
79 
80  clearCurrentXSData();
81 
82  delete names;
83 }
84 
86  G4int /*Z*/ , G4int /*A*/ ,
87  const G4Element* element ,
88  const G4Material* material )
89 {
90  G4double eKin = dp->GetKineticEnergy();
91  if ( eKin > 4.0*eV //GetMaxKinEnergy()
92  || eKin < 0 //GetMinKinEnergy()
93  || dp->GetDefinition() != G4Neutron::Neutron() ) return false;
94 
95  if ( dic.find( std::pair < const G4Material* , const G4Element* > ( (G4Material*)NULL , element ) ) != dic.end()
96  || dic.find( std::pair < const G4Material* , const G4Element* > ( material , element ) ) != dic.end() ) return true;
97 
98  return false;
99 
100 // return IsApplicable( dp , element );
101 /*
102  G4double eKin = dp->GetKineticEnergy();
103  if ( eKin > 4.0*eV //GetMaxKinEnergy()
104  || eKin < 0 //GetMinKinEnergy()
105  || dp->GetDefinition() != G4Neutron::Neutron() ) return false;
106  return true;
107 */
108 }
109 
111  G4int /*Z*/ , G4int /*A*/ ,
112  const G4Isotope* /*iso*/ ,
113  const G4Element* element ,
114  const G4Material* material )
115 {
116  //if ( dp->GetKineticEnergy() == ke_cache && element == element_cache && material == material_cache ) return xs_cache;
117 
118  ke_cache = dp->GetKineticEnergy();
119  element_cache = element;
120  material_cache = material;
121  //G4double xs = GetCrossSection( dp , element , material->GetTemperature() );
122  G4double xs = GetCrossSection( dp , element , material );
123  xs_cache = xs;
124  return xs;
125  //return GetCrossSection( dp , element , material->GetTemperature() );
126 }
127 
128 void G4ParticleHPThermalScatteringData::clearCurrentXSData()
129 {
130  std::map< G4int , std::map< G4double , G4ParticleHPVector* >* >::iterator it;
131  std::map< G4double , G4ParticleHPVector* >::iterator itt;
132 
133  if ( coherent != NULL ) {
134  for ( it = coherent->begin() ; it != coherent->end() ; it++ )
135  {
136  if ( it->second != NULL )
137  {
138  for ( itt = it->second->begin() ; itt != it->second->end() ; itt++ )
139  {
140  delete itt->second;
141  }
142  }
143  delete it->second;
144  }
145  coherent->clear();
146  }
147 
148  if ( incoherent != NULL ) {
149  for ( it = incoherent->begin() ; it != incoherent->end() ; it++ )
150  {
151  if ( it->second != NULL )
152  {
153  for ( itt = it->second->begin() ; itt != it->second->end() ; itt++ )
154  {
155  delete itt->second;
156  }
157  }
158  delete it->second;
159  }
160  incoherent->clear();
161  }
162 
163  if ( inelastic != NULL ) {
164  for ( it = inelastic->begin() ; it != inelastic->end() ; it++ )
165  {
166  if ( it->second != NULL )
167  {
168  for ( itt = it->second->begin() ; itt != it->second->end() ; itt++ )
169  {
170  delete itt->second;
171  }
172  }
173  delete it->second;
174  }
175  inelastic->clear();
176  }
177 
178 }
179 
180 
181 
183 {
184  G4bool result = false;
185 
186  G4double eKin = aP->GetKineticEnergy();
187  // Check energy
188  if ( eKin < emax )
189  {
190  // Check Particle Species
191  if ( aP->GetDefinition() == G4Neutron::Neutron() )
192  {
193  // anEle is one of Thermal elements
194  G4int ie = (G4int) anEle->GetIndex();
195  std::vector < G4int >::iterator it;
196  for ( it = indexOfThermalElement.begin() ; it != indexOfThermalElement.end() ; it++ )
197  {
198  if ( ie == *it ) return true;
199  }
200  }
201  }
202 
203 /*
204  if ( names->IsThisThermalElement ( anEle->GetName() ) )
205  {
206  // Check energy and projectile species
207  G4double eKin = aP->GetKineticEnergy();
208  if ( eKin < emax && aP->GetDefinition() == G4Neutron::Neutron() ) result = true;
209  }
210 */
211  return result;
212 }
213 
214 
216 {
217 
218  if ( &aP != G4Neutron::Neutron() )
219  throw G4HadronicException(__FILE__, __LINE__, "Attempt to use NeutronHP data for particles other than neutrons!!!");
220 
221  //std::map < std::pair < G4Material* , const G4Element* > , G4int > dic;
222  //
223  dic.clear();
224  if ( G4Threading::IsMasterThread() ) clearCurrentXSData();
225 
226  std::map < G4String , G4int > co_dic;
227 
228  //Searching Nist Materials
229  static G4ThreadLocal G4MaterialTable* theMaterialTable = 0 ; if (!theMaterialTable) theMaterialTable= G4Material::GetMaterialTable();
230  size_t numberOfMaterials = G4Material::GetNumberOfMaterials();
231  for ( size_t i = 0 ; i < numberOfMaterials ; i++ )
232  {
233  G4Material* material = (*theMaterialTable)[i];
234  size_t numberOfElements = material->GetNumberOfElements();
235  for ( size_t j = 0 ; j < numberOfElements ; j++ )
236  {
237  const G4Element* element = material->GetElement(j);
238  if ( names->IsThisThermalElement ( material->GetName() , element->GetName() ) )
239  {
240  G4int ts_ID_of_this_geometry;
241  G4String ts_ndl_name = names->GetTS_NDL_Name( material->GetName() , element->GetName() );
242  if ( co_dic.find ( ts_ndl_name ) != co_dic.end() )
243  {
244  ts_ID_of_this_geometry = co_dic.find ( ts_ndl_name ) -> second;
245  }
246  else
247  {
248  ts_ID_of_this_geometry = co_dic.size();
249  co_dic.insert ( std::pair< G4String , G4int >( ts_ndl_name , ts_ID_of_this_geometry ) );
250  }
251 
252  //G4cout << "Neutron HP Thermal Scattering Data : Registering a material-element pair of "
253  // << material->GetName() << " " << element->GetName()
254  // << " as internal thermal scattering id of " << ts_ID_of_this_geometry << "." << G4endl;
255 
256  dic.insert( std::pair < std::pair < G4Material* , const G4Element* > , G4int > ( std::pair < G4Material* , const G4Element* > ( material , element ) , ts_ID_of_this_geometry ) );
257  }
258  }
259  }
260 
261  //Searching TS Elements
262  static G4ThreadLocal G4ElementTable* theElementTable = 0 ; if (!theElementTable) theElementTable= G4Element::GetElementTable();
263  size_t numberOfElements = G4Element::GetNumberOfElements();
264  //size_t numberOfThermalElements = 0;
265  for ( size_t i = 0 ; i < numberOfElements ; i++ )
266  {
267  const G4Element* element = (*theElementTable)[i];
268  if ( names->IsThisThermalElement ( element->GetName() ) )
269  {
270  if ( names->IsThisThermalElement ( element->GetName() ) )
271  {
272  G4int ts_ID_of_this_geometry;
273  G4String ts_ndl_name = names->GetTS_NDL_Name( element->GetName() );
274  if ( co_dic.find ( ts_ndl_name ) != co_dic.end() )
275  {
276  ts_ID_of_this_geometry = co_dic.find ( ts_ndl_name ) -> second;
277  }
278  else
279  {
280  ts_ID_of_this_geometry = co_dic.size();
281  co_dic.insert ( std::pair< G4String , G4int >( ts_ndl_name , ts_ID_of_this_geometry ) );
282  }
283 
284  //G4cout << "Neutron HP Thermal Scattering: Registering an element of "
285  // << material->GetName() << " " << element->GetName()
286  // << " as internal thermal scattering id of " << ts_ID_of_this_geometry << "." << G4endl;
287 
288  dic.insert( std::pair < std::pair < const G4Material* , const G4Element* > , G4int > ( std::pair < const G4Material* , const G4Element* > ( (G4Material*)NULL , element ) , ts_ID_of_this_geometry ) );
289  }
290  }
291  }
292 
293  G4cout << G4endl;
294  G4cout << "Neutron HP Thermal Scattering Data: Following material-element pairs and/or elements are registered." << G4endl;
295  for ( std::map < std::pair < const G4Material* , const G4Element* > , G4int >::iterator it = dic.begin() ; it != dic.end() ; it++ )
296  {
297  if ( it->first.first != NULL )
298  {
299  G4cout << "Material " << it->first.first->GetName() << " - Element " << it->first.second->GetName() << ", internal thermal scattering id " << it->second << G4endl;
300  }
301  else
302  {
303  G4cout << "Element " << it->first.second->GetName() << ", internal thermal scattering id " << it->second << G4endl;
304  }
305  }
306  G4cout << G4endl;
307 
308 
309  //G4cout << "Neutron HP Thermal Scattering Data: Following NDL thermal scattering files are assigned to the internal thermal scattering id." << G4endl;
310  //for ( std::map < G4String , G4int >::iterator it = co_dic.begin() ; it != co_dic.end() ; it++ )
311  //{
312  // G4cout << "NDL file name " << it->first << ", internal thermal scattering id " << it->second << G4endl;
313  //}
314 
316 
317  coherent = hpmanager->GetThermalScatteringCoherentCrossSections();
318  incoherent = hpmanager->GetThermalScatteringIncoherentCrossSections();
319  inelastic = hpmanager->GetThermalScatteringInelasticCrossSections();
320 
321  if ( G4Threading::IsMasterThread() ) {
322 
323  if ( coherent == NULL ) coherent = new std::map< G4int , std::map< G4double , G4ParticleHPVector* >* >;
324  if ( incoherent == NULL ) incoherent = new std::map< G4int , std::map< G4double , G4ParticleHPVector* >* >;
325  if ( inelastic == NULL ) inelastic = new std::map< G4int , std::map< G4double , G4ParticleHPVector* >* >;
326 
327 
328  // Read Cross Section Data files
329 
330  G4String dirName;
331  if ( !getenv( "G4NEUTRONHPDATA" ) )
332  throw G4HadronicException(__FILE__, __LINE__, "Please setenv G4NEUTRONHPDATA to point to the neutron cross-section files.");
333  G4String baseName = getenv( "G4NEUTRONHPDATA" );
334 
335  dirName = baseName + "/ThermalScattering";
336 
337  G4String ndl_filename;
338  G4String full_name;
339 
340  for ( std::map < G4String , G4int >::iterator it = co_dic.begin() ; it != co_dic.end() ; it++ )
341  {
342 
343  ndl_filename = it->first;
344  G4int ts_ID = it->second;
345 
346  // Coherent
347  full_name = dirName + "/Coherent/CrossSection/" + ndl_filename;
348  std::map< G4double , G4ParticleHPVector* >* coh_amapTemp_EnergyCross = readData( full_name );
349  coherent->insert ( std::pair < G4int , std::map< G4double , G4ParticleHPVector* >* > ( ts_ID , coh_amapTemp_EnergyCross ) );
350 
351  // Incoherent
352  full_name = dirName + "/Incoherent/CrossSection/" + ndl_filename;
353  std::map< G4double , G4ParticleHPVector* >* incoh_amapTemp_EnergyCross = readData( full_name );
354  incoherent->insert ( std::pair < G4int , std::map< G4double , G4ParticleHPVector* >* > ( ts_ID , incoh_amapTemp_EnergyCross ) );
355 
356  // Inelastic
357  full_name = dirName + "/Inelastic/CrossSection/" + ndl_filename;
358  std::map< G4double , G4ParticleHPVector* >* inela_amapTemp_EnergyCross = readData( full_name );
359  inelastic->insert ( std::pair < G4int , std::map< G4double , G4ParticleHPVector* >* > ( ts_ID , inela_amapTemp_EnergyCross ) );
360 
361  }
362  hpmanager->RegisterThermalScatteringCoherentCrossSections( coherent );
363  hpmanager->RegisterThermalScatteringIncoherentCrossSections( incoherent );
364  hpmanager->RegisterThermalScatteringInelasticCrossSections( inelastic );
365  }
366 }
367 
368 
369 
370 std::map< G4double , G4ParticleHPVector* >* G4ParticleHPThermalScatteringData::readData ( G4String full_name )
371 {
372 
373  std::map< G4double , G4ParticleHPVector* >* aData = new std::map< G4double , G4ParticleHPVector* >;
374 
375  //std::ifstream theChannel( full_name.c_str() );
376  std::istringstream theChannel;
377  G4ParticleHPManager::GetInstance()->GetDataStream(full_name,theChannel);
378 
379  //G4cout << "G4ParticleHPThermalScatteringData " << name << G4endl;
380 
381  G4int dummy;
382  while ( theChannel >> dummy ) // MF // Loop checking, 11.05.2015, T. Koi
383  {
384  theChannel >> dummy; // MT
385  G4double temp;
386  theChannel >> temp;
387  G4ParticleHPVector* anEnergyCross = new G4ParticleHPVector;
388  G4int nData;
389  theChannel >> nData;
390  anEnergyCross->Init ( theChannel , nData , eV , barn );
391  aData->insert ( std::pair < G4double , G4ParticleHPVector* > ( temp , anEnergyCross ) );
392  }
393  //theChannel.close();
394 
395  return aData;
396 
397 }
398 
399 
400 
402 {
403  if( &aP != G4Neutron::Neutron() )
404  throw G4HadronicException(__FILE__, __LINE__, "Attempt to use NeutronHP data for particles other than neutrons!!!");
405 // G4cout << "G4ParticleHPThermalScatteringData::DumpPhysicsTable still to be implemented"<<G4endl;
406 }
407 
408 //#include "G4Nucleus.hh"
409 //#include "G4NucleiPropertiesTable.hh"
410 //#include "G4Neutron.hh"
411 //#include "G4Electron.hh"
412 
413 
414 
415 /*
416 G4double G4ParticleHPThermalScatteringData::GetCrossSection( const G4DynamicParticle* aP , const G4Element*anE , G4double aT )
417 {
418 
419  G4double result = 0;
420  const G4Material* aM = NULL;
421 
422  G4int iele = anE->GetIndex();
423 
424  if ( dic.find( std::pair < const G4Material* , const G4Element* > ( (G4Material*)NULL , anE ) ) != dic.end() )
425  {
426  iele = dic.find( std::pair < const G4Material* , const G4Element* > ( (G4Material*)NULL , anE ) )->second;
427  }
428  else if ( dic.find( std::pair < const G4Material* , const G4Element* > ( aM , anE ) ) != dic.end() )
429  {
430  iele = dic.find( std::pair < const G4Material* , const G4Element* > ( aM , anE ) )->second;
431  }
432  else
433  {
434  return result;
435  }
436 
437  G4double Xcoh = GetX ( aP , aT , coherent.find(iele)->second );
438  G4double Xincoh = GetX ( aP , aT , incoherent.find(iele)->second );
439  G4double Xinela = GetX ( aP , aT , inelastic.find(iele)->second );
440 
441  result = Xcoh + Xincoh + Xinela;
442 
443  //G4cout << "G4ParticleHPThermalScatteringData::GetCrossSection Tot= " << result/barn << " Coherent= " << Xcoh/barn << " Incoherent= " << Xincoh/barn << " Inelastic= " << Xinela/barn << G4endl;
444 
445  return result;
446 
447 }
448 */
449 
451 {
452  G4double result = 0;
453 
454  G4int ts_id =getTS_ID( aM , anE );
455 
456  if ( ts_id == -1 ) return result;
457 
458  G4double aT = aM->GetTemperature();
459 
460  G4double Xcoh = GetX ( aP , aT , coherent->find(ts_id)->second );
461  G4double Xincoh = GetX ( aP , aT , incoherent->find(ts_id)->second );
462  G4double Xinela = GetX ( aP , aT , inelastic->find(ts_id)->second );
463 
464  result = Xcoh + Xincoh + Xinela;
465 
466  //G4cout << "G4ParticleHPThermalScatteringData::GetCrossSection Tot= " << result/barn << " Coherent= " << Xcoh/barn << " Incoherent= " << Xincoh/barn << " Inelastic= " << Xinela/barn << G4endl;
467 
468  return result;
469 }
470 
471 
473 {
474  G4double result = 0;
475  G4int ts_id = getTS_ID( aM , anE );
476  G4double aT = aM->GetTemperature();
477  result = GetX ( aP , aT , inelastic->find( ts_id )->second );
478  return result;
479 }
480 
482 {
483  G4double result = 0;
484  G4int ts_id = getTS_ID( aM , anE );
485  G4double aT = aM->GetTemperature();
486  result = GetX ( aP , aT , coherent->find( ts_id )->second );
487  return result;
488 }
489 
491 {
492  G4double result = 0;
493  G4int ts_id = getTS_ID( aM , anE );
494  G4double aT = aM->GetTemperature();
495  result = GetX ( aP , aT , incoherent->find( ts_id )->second );
496  return result;
497 }
498 
499 
500 
501 G4int G4ParticleHPThermalScatteringData::getTS_ID ( const G4Material* material , const G4Element* element )
502 {
503  G4int result = -1;
504  if ( dic.find( std::pair < const G4Material* , const G4Element* > ( (G4Material*)NULL , element ) ) != dic.end() )
505  return dic.find( std::pair < const G4Material* , const G4Element* > ( (G4Material*)NULL , element ) )->second;
506  if ( dic.find( std::pair < const G4Material* , const G4Element* > ( material , element ) ) != dic.end() )
507  return dic.find( std::pair < const G4Material* , const G4Element* > ( material , element ) )->second;
508  return result;
509 }
510 
511 
512 
513 
514 G4double G4ParticleHPThermalScatteringData::GetX ( const G4DynamicParticle* aP, G4double aT , std::map < G4double , G4ParticleHPVector* >* amapTemp_EnergyCross )
515 {
516 
517  G4double result = 0;
518  if ( amapTemp_EnergyCross->size() == 0 ) return result;
519 
520 
521  G4double eKinetic = aP->GetKineticEnergy();
522 
523  if ( amapTemp_EnergyCross->size() == 1 ) {
524  if ( std::fabs ( aT - amapTemp_EnergyCross->begin()->first ) / amapTemp_EnergyCross->begin()->first > 0.1 ) {
525  G4cout << "G4ParticleHPThermalScatteringData:: The temperature of material ("
526  << aT/kelvin << "K) is different more than 10% from temperature of thermal scattering file expected ("
527  << amapTemp_EnergyCross->begin()->first << "K). Result may not be reliable."
528  << G4endl;
529  }
530  result = amapTemp_EnergyCross->begin()->second->GetXsec ( eKinetic );
531  return result;
532  }
533 
534  std::map< G4double , G4ParticleHPVector* >::iterator it;
535  for ( it = amapTemp_EnergyCross->begin() ; it != amapTemp_EnergyCross->end() ; it++ ) {
536  if ( aT < it->first ) break;
537  }
538  //if ( it == amapTemp_EnergyCross->begin() && it != amapTemp_EnergyCross->end() ) it++; // lower than the first
539  //if ( it != amapTemp_EnergyCross->begin() && it == amapTemp_EnergyCross->end() ) it--; // upper than the last
540  if ( it == amapTemp_EnergyCross->begin() ) {
541  it++; // lower than the first
542  } else if ( it == amapTemp_EnergyCross->end() ) {
543  it--; // upper than the last
544  }
545 
546  G4double TH = it->first;
547  G4double XH = it->second->GetXsec ( eKinetic );
548 
549  //G4cout << "G4ParticleHPThermalScatteringData::GetX TH " << TH << " E " << eKinetic << " XH " << XH << G4endl;
550 
551  if ( it != amapTemp_EnergyCross->begin() ) it--;
552  G4double TL = it->first;
553  G4double XL = it->second->GetXsec ( eKinetic );
554 
555  //G4cout << "G4ParticleHPThermalScatteringData::GetX TL " << TL << " E " << eKinetic << " XL " << XL << G4endl;
556 
557  if ( TH == TL )
558  throw G4HadronicException(__FILE__, __LINE__, "Thermal Scattering Data Error!");
559 
560  G4double T = aT;
561  G4double X = ( XH - XL ) / ( TH - TL ) * ( T - TL ) + XL;
562  result = X;
563 
564  return result;
565 }
566 
567 
569 {
570  names->AddThermalElement( nameG4Element , filename );
571 }
573 {
574  outFile << "High Precision cross data based on thermal scattering data in evaluated nuclear data libraries for neutrons below 5eV on specific materials\n" ;
575 }
G4double G4ParticleHPJENDLHEData::G4double result
static G4ParticleHPManager * GetInstance()
void RegisterThermalScatteringIncoherentCrossSections(std::map< G4int, std::map< G4double, G4ParticleHPVector * > * > *val)
G4double GetKineticEnergy() const
G4int first(char) const
void DumpPhysicsTable(const G4ParticleDefinition &)
const G4String & GetName() const
Definition: G4Material.hh:178
static G4MaterialTable * GetMaterialTable()
Definition: G4Material.cc:587
std::vector< G4Material * > G4MaterialTable
void GetDataStream(G4String, std::istringstream &iss)
void RegisterThermalScatteringInelasticCrossSections(std::map< G4int, std::map< G4double, G4ParticleHPVector * > * > *val)
G4ParticleDefinition * GetDefinition() const
static constexpr double second
Definition: G4SIunits.hh:157
G4bool IsApplicable(const G4DynamicParticle *, const G4Element *)
#define G4ThreadLocal
Definition: tls.hh:89
const G4Element * GetElement(G4int iel) const
Definition: G4Material.hh:202
int G4int
Definition: G4Types.hh:78
G4double GetIsoCrossSection(const G4DynamicParticle *, G4int, G4int, const G4Isotope *, const G4Element *, const G4Material *)
G4double GetCoherentCrossSection(const G4DynamicParticle *, const G4Element *, const G4Material *)
G4double GetIncoherentCrossSection(const G4DynamicParticle *, const G4Element *, const G4Material *)
string material
Definition: eplot.py:19
void Init(std::istream &aDataFile, G4int total, G4double ux=1., G4double uy=1.)
G4GLOB_DLL std::ostream G4cout
static size_t GetNumberOfElements()
Definition: G4Element.cc:405
std::map< G4int, std::map< G4double, G4ParticleHPVector * > * > * GetThermalScatteringIncoherentCrossSections()
void SetMinKinEnergy(G4double value)
bool G4bool
Definition: G4Types.hh:79
size_t GetIndex() const
Definition: G4Element.hh:182
static constexpr double eV
Definition: G4SIunits.hh:215
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
static size_t GetNumberOfMaterials()
Definition: G4Material.cc:594
static constexpr double kelvin
Definition: G4SIunits.hh:281
G4bool IsIsoApplicable(const G4DynamicParticle *, G4int, G4int, const G4Element *, const G4Material *)
virtual void CrossSectionDescription(std::ostream &) const
void SetMaxKinEnergy(G4double value)
G4double GetTemperature() const
Definition: G4Material.hh:182
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
G4bool IsMasterThread()
Definition: G4Threading.cc:146
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
G4double GetCrossSection(const G4DynamicParticle *, const G4Element *, const G4Material *)
std::vector< G4Element * > G4ElementTable
double G4double
Definition: G4Types.hh:76
void RegisterThermalScatteringCoherentCrossSections(std::map< G4int, std::map< G4double, G4ParticleHPVector * > * > *val)
G4double GetInelasticCrossSection(const G4DynamicParticle *, const G4Element *, const G4Material *)
const G4String & GetName() const
Definition: G4Element.hh:127
static constexpr double barn
Definition: G4SIunits.hh:105
static G4ElementTable * GetElementTable()
Definition: G4Element.cc:398
void BuildPhysicsTable(const G4ParticleDefinition &)
std::map< G4int, std::map< G4double, G4ParticleHPVector * > * > * GetThermalScatteringInelasticCrossSections()
std::map< G4int, std::map< G4double, G4ParticleHPVector * > * > * GetThermalScatteringCoherentCrossSections()