Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4NeutronElasticXS.cc
이 파일의 문서화 페이지로 가기
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: G4NeutronElasticXS.cc 110744 2018-06-12 06:35:40Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class file
31 //
32 //
33 // File name: G4NeutronElasticXS
34 //
35 // Author Ivantchenko, Geant4, 3-Aug-09
36 //
37 // Modifications:
38 //
39 
40 #include "G4NeutronElasticXS.hh"
41 #include "G4Neutron.hh"
42 #include "G4DynamicParticle.hh"
43 #include "G4Element.hh"
44 #include "G4ElementTable.hh"
45 #include "G4PhysicsLogVector.hh"
46 #include "G4PhysicsVector.hh"
48 #include "G4HadronNucleonXsc.hh"
49 #include "G4NistManager.hh"
50 #include "G4Proton.hh"
51 
52 #include <iostream>
53 #include <fstream>
54 #include <sstream>
55 
56 // factory
57 #include "G4CrossSectionFactory.hh"
58 //
60 
61 using namespace std;
62 
65 
66 #ifdef G4MULTITHREADED
67  G4Mutex G4NeutronElasticXS::neutronElasticXSMutex = G4MUTEX_INITIALIZER;
68 #endif
69 
71  : G4VCrossSectionDataSet(Default_Name()),
73 {
74  // verboseLevel = 0;
75  if(verboseLevel > 0){
76  G4cout << "G4NeutronElasticXS::G4NeutronElasticXS Initialise for Z < "
77  << MAXZEL << G4endl;
78  }
81  isMaster = false;
82 }
83 
85 {
86  delete fNucleon;
87  delete ggXsection;
88  if(isMaster) {
89  for(G4int i=0; i<MAXZEL; ++i) {
90  delete data[i];
91  data[i] = nullptr;
92  }
93  }
94 }
95 
96 void G4NeutronElasticXS::CrossSectionDescription(std::ostream& outFile) const
97 {
98  outFile << "G4NeutronElasticXS calculates the neutron elastic scattering\n"
99  << "cross section on nuclei using data from the high precision\n"
100  << "neutron database. These data are simplified and smoothed over\n"
101  << "the resonance region in order to reduce CPU time.\n"
102  << "For high energies Glauber-Gribiv cross section is used.\n";
103 }
104 
105 G4bool
107  G4int, const G4Material*)
108 {
109  return true;
110 }
111 
112 G4double
114  G4int ZZ, const G4Material*)
115 {
116  G4double xs = 0.0;
117  G4double ekin = aParticle->GetKineticEnergy();
118 
119  G4int Z = (ZZ >= MAXZEL) ? MAXZEL - 1 : ZZ;
120 
121  G4PhysicsVector* pv = data[Z];
122  // G4cout << "G4NeutronElasticXS::GetCrossSection e= " << ekin
123  // << " Z= " << Z << G4endl;
124 
125  // element was not initialised
126  if(!pv) { return xs; }
127 
128  if(ekin <= pv->Energy(0)) { return (*pv)[0]; }
129 
130  if(ekin <= pv->GetMaxEnergy()) {
131  xs = pv->Value(ekin);
132  } else if(1 == Z) {
135  } else {
136  G4int Amean =
137  G4lrint(G4NistManager::Instance()->GetAtomicMassAmu(Z));
138  ggXsection->GetIsoCrossSection(aParticle, Z, Amean);
140  }
141 
142  if(verboseLevel > 0){
143  G4cout << "Z= " << Z << " Ekin(MeV)= " << ekin/CLHEP::MeV
144  << ", nElmXSel(bn)= " << xs/CLHEP::barn
145  << G4endl;
146  }
147  return xs;
148 }
149 
150 void
152 {
153  if(verboseLevel > 0){
154  G4cout << "G4NeutronElasticXS::BuildPhysicsTable for "
155  << p.GetParticleName() << G4endl;
156  }
157  if(p.GetParticleName() != "neutron") {
159  ed << p.GetParticleName() << " is a wrong particle type -"
160  << " only neutron is allowed";
161  G4Exception("G4NeutronElasticXS::BuildPhysicsTable(..)","had012",
162  FatalException, ed, "");
163  return;
164  }
165 
166  if(!data[1]) {
167 #ifdef G4MULTITHREADED
168  G4MUTEXLOCK(&neutronElasticXSMutex);
169  if(!data[1]) {
170 #endif
171  isMaster = true;
172 #ifdef G4MULTITHREADED
173  }
174  G4MUTEXUNLOCK(&neutronElasticXSMutex);
175 #endif
176  }
177 
178  // it is possible re-initialisation for the second run
179  if(isMaster) {
180 
181  // check environment variable
182  // Build the complete string identifying the file with the data set
183  char* path = getenv("G4NEUTRONXSDATA");
184 
185  G4DynamicParticle* dynParticle =
187 
188  // Access to elements
189  const G4ElementTable* theElmTable = G4Element::GetElementTable();
190  size_t numOfElm = G4Element::GetNumberOfElements();
191  for(size_t i=0; i<numOfElm; ++i) {
192  G4int Z = ((*theElmTable)[i])->GetZasInt();
193  if(Z >= MAXZEL) { Z = MAXZEL-1; }
194  //G4cout << "Z= " << Z << G4endl;
195  // Initialisation
196  if(!data[Z]) { Initialise(Z, dynParticle, path); }
197  }
198  delete dynParticle;
199  }
200 }
201 
202 void
204  const char* p)
205 {
206  if(data[Z]) { return; }
207  const char* path = p;
208  if(!p) {
209  // check environment variable
210  // Build the complete string identifying the file with the data set
211  path = getenv("G4NEUTRONXSDATA");
212  if (!path) {
213  G4Exception("G4NeutronElasticXS::Initialise(..)","had013",
215  "Environment variable G4NEUTRONXSDATA is not defined");
216  return;
217  }
218  }
219 
220  // upload data from file
221  data[Z] = new G4PhysicsLogVector();
222 
223  std::ostringstream ost;
224  ost << path << "/neutron/el" << Z ;
225  std::ifstream filein(ost.str().c_str());
226  if (!(filein)) {
228  ed << "Data file <" << ost.str().c_str()
229  << "> is not opened!";
230  G4Exception("G4NeutronElasticXS::Initialise(..)","had014",
231  FatalException, ed, "Check G4NEUTRONXSDATA");
232  return;
233  }else{
234  if(verboseLevel > 1) {
235  G4cout << "file " << ost.str()
236  << " is opened by G4NeutronElasticXS" << G4endl;
237  }
238 
239  // retrieve data from DB
240  if(!data[Z]->Retrieve(filein, true)) {
242  ed << "Data file <" << ost.str().c_str()
243  << "> is not retrieved!";
244  G4Exception("G4NeutronElasticXS::Initialise(..)","had015",
245  FatalException, ed, "Check G4NEUTRONXSDATA");
246  return;
247  }
248 
249  // smooth transition
250  G4double sig1 = (*(data[Z]))[data[Z]->GetVectorLength()-1];
251  dp->SetKineticEnergy(data[Z]->GetMaxEnergy());
252  G4double sig2 = 0.0;
253  if(1 == Z) {
256  } else {
257  G4int Amean =
258  G4lrint(G4NistManager::Instance()->GetAtomicMassAmu(Z));
259  ggXsection->GetIsoCrossSection(dp, Z, Amean);
261  }
262  if(sig2 > 0.) { coeff[Z] = sig1/sig2; }
263  }
264 }
G4double GetIsoCrossSection(const G4DynamicParticle *, G4int Z, G4int A, const G4Isotope *iso=0, const G4Element *elm=0, const G4Material *mat=0)
CLHEP::Hep3Vector G4ThreeVector
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:45
void SetKineticEnergy(G4double aEnergy)
virtual G4bool IsElementApplicable(const G4DynamicParticle *, G4int Z, const G4Material *)
#define G4endl
Definition: G4ios.hh:61
const char * p
Definition: xmltok.h:285
static G4PhysicsVector * data[MAXZEL]
G4HadronNucleonXsc * fNucleon
const G4String & GetParticleName() const
const G4int MAXZEL
G4double Value(G4double theEnergy, size_t &lastidx) const
G4double GetElasticHadronNucleonXsc()
const G4ParticleDefinition * proton
const XML_Char const XML_Char * data
Definition: expat.h:268
Float_t Z
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
#define G4MUTEX_INITIALIZER
Definition: G4Threading.hh:88
G4double GetHadronNucleonXscPDG(const G4DynamicParticle *, const G4ParticleDefinition *)
G4ComponentGGHadronNucleusXsc * ggXsection
std::vector< G4Element * > G4ElementTable
static constexpr double MeV
virtual void CrossSectionDescription(std::ostream &) const
#define G4_DECLARE_XS_FACTORY(cross_section)
virtual void BuildPhysicsTable(const G4ParticleDefinition &)
#define G4MUTEXUNLOCK(mutex)
Definition: G4Threading.hh:234
void Initialise(G4int Z, G4DynamicParticle *dp, const char *)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.hh:65
int G4lrint(double ad)
Definition: templates.hh:151
int G4int
Definition: G4Types.hh:78
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
static G4double coeff[MAXZEL]
virtual G4double GetElementCrossSection(const G4DynamicParticle *, G4int Z, const G4Material *mat=nullptr)
#define G4MUTEXLOCK(mutex)
Definition: G4Threading.hh:233
G4double GetKineticEnergy() const
G4GLOB_DLL std::ostream G4cout
static size_t GetNumberOfElements()
Definition: G4Element.cc:405
static constexpr double barn
Definition: SystemOfUnits.h:85
static G4NistManager * Instance()
static G4ElementTable * GetElementTable()
Definition: G4Element.cc:398
std::mutex G4Mutex
Definition: G4Threading.hh:84