Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4HadronXSDataTable.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: G4HadronXSDataTable.cc 96794 2016-05-09 10:09:30Z gcosmo $
27 //
28 //---------------------------------------------------------------------------
29 //
30 // GEANT4 Class file
31 //
32 // Description: Data structure for cross sections per materials
33 //
34 // Author: V.Ivanchenko 31.05.2018
35 //
36 // Modifications:
37 //
38 //----------------------------------------------------------------------------
39 //
40 
41 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
42 
43 #include "G4HadronXSDataTable.hh"
44 #include "G4PhysicsLogVector.hh"
45 #include "G4Material.hh"
46 #include "G4MaterialTable.hh"
47 #include "G4DynamicParticle.hh"
49 
50 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
51 
54  const G4Material* mat,
55  G4int bins, G4double emin,
56  G4double emax, G4bool spline)
57 {
58  G4int n = mat->GetNumberOfElements();
59  nElmMinusOne = n - 1;
61  if(nElmMinusOne > 0) {
62  G4PhysicsVector* first = nullptr;
63  xSections.resize(n, first);
64  first = new G4PhysicsLogVector(emin,emax,bins);
65  first->SetSpline(spline);
66  xSections[0] = first;
67  for(G4int i=1; i<n; ++i) {
68  xSections[i] = new G4PhysicsVector(*first);
69  }
70  std::vector<double> temp;
71  temp.resize(n, 0.0);
72  for(G4int j=0; j<=bins; ++j) {
73  G4double cross = 0.0;
74  G4double e = first->Energy(j);
75  dp->SetKineticEnergy(e);
76  for(G4int i=0; i<n; ++i) {
77  cross += xs->GetCrossSection(dp, (*theElementVector)[i], mat);
78  temp[i] = cross;
79  }
80  G4double fact = (cross > 0.0) ? 1.0/cross : 0.0;
81  for(G4int i=0; i<n; ++i) {
82  G4double y = (i<n-1) ? temp[i]*fact : 1.0;
83  xSections[i]->PutValue(j, y);
84  }
85  }
86  }
87 }
88 
89 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
90 
92 {
93  if(nElmMinusOne > 0) {
94  for(G4int i=0; i<=nElmMinusOne; ++i) { delete xSections[i]; }
95  }
96 }
97 
98 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
99 
101 {}
102 
103 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
104 
106 {}
107 
108 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
109 
112  G4int bins, G4double emin, G4double emax,
113  G4bool spline)
114 {
116  if(nn > nMaterials) {
117  G4PhysicsLogVector* first = nullptr;
118  G4int sbins = std::max(10, bins/5);
120  for(size_t i=nMaterials; i<nn; ++i) {
121  const G4Material* mat = (*mtable)[i];
122  G4PhysicsVector* v = nullptr;
123  G4HadElementSelector* es = nullptr;
124  // create real vector only for complex materials
125  if(mat->GetNumberOfElements() > 1) {
126  if(!first) {
127  first = new G4PhysicsLogVector(emin, emax, bins);
128  first->SetSpline(spline);
129  v = first;
130  } else {
131  v = new G4PhysicsVector(*first);
132  }
133  for(G4int j=0; j<=bins; ++j) {
134  G4double e = first->Energy(j);
135  dp->SetKineticEnergy(e);
136  G4double cros = xs->ComputeCrossSection(dp, mat);
137  v->PutValue(j, cros);
138  }
139  elmSelectors[i] = new G4HadElementSelector(dp, xs, mat, sbins, emin, emax, spline);
140  }
141  xsData.push_back(v);
142  elmSelectors.push_back(es);
143  }
144  nMaterials = nn;
145  }
146 }
147 
148 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
149 
151 {
152  for(size_t i=0; i<nMaterials; ++i) {
153  delete xsData[i];
154  delete elmSelectors[i];
155  }
156 }
157 
158 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
159 
161 {}
162 
163 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
G4double Energy(size_t index) const
T max(const T t1, const T t2)
brief Return the largest of the two arguments
static size_t GetNumberOfMaterials()
Definition: G4Material.cc:600
static G4MaterialTable * GetMaterialTable()
Definition: G4Material.cc:593
void SetKineticEnergy(G4double aEnergy)
std::vector< G4PhysicsVector * > xSections
Float_t y
Definition: compare.C:6
G4double ComputeCrossSection(const G4DynamicParticle *, const G4Material *)
G4HadElementSelector(G4DynamicParticle *, G4CrossSectionDataStore *, const G4Material *, G4int bins, G4double emin, G4double emax, G4bool spline)
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:191
static const G4double emax
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
std::vector< G4HadElementSelector * > elmSelectors
void SetSpline(G4bool)
Float_t mat
void Initialise(G4DynamicParticle *, G4CrossSectionDataStore *, G4int bins, G4double emin, G4double emax, G4bool spline)
std::vector< G4Material * > G4MaterialTable
int G4int
Definition: G4Types.hh:78
Char_t n[5]
G4double GetCrossSection(const G4DynamicParticle *, const G4Material *)
std::vector< G4PhysicsVector * > xsData
const G4ElementVector * theElementVector
void PutValue(size_t index, G4double theValue)
size_t GetNumberOfElements() const
Definition: G4Material.hh:187