| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
>> 26 // $Id: G4VCrossSectionDataSet.cc 66241 2012-12-13 18:34:42Z gunter $
>> 27 //
26 // ------------------------------------------- 28 // -------------------------------------------------------------------
27 // 29 //
28 // GEANT4 Class file 30 // GEANT4 Class file
29 // 31 //
>> 32 //
30 // File name: G4VCrossSectionDataSet 33 // File name: G4VCrossSectionDataSet
31 // 34 //
32 // Author F.W. Jones, TRIUMF, 20-JAN-97 35 // Author F.W. Jones, TRIUMF, 20-JAN-97
33 // 36 //
34 // Modifications: 37 // Modifications:
35 // 23.01.2009 V.Ivanchenko move constructor an 38 // 23.01.2009 V.Ivanchenko move constructor and destructor to source
36 // 12.08.2011 G.Folger, V.Ivanchenko, T.Koi, D 39 // 12.08.2011 G.Folger, V.Ivanchenko, T.Koi, D.Wright redesign the class
37 // 40 //
38 41
39 #include "G4VCrossSectionDataSet.hh" 42 #include "G4VCrossSectionDataSet.hh"
40 #include "G4SystemOfUnits.hh" 43 #include "G4SystemOfUnits.hh"
41 #include "G4CrossSectionDataSetRegistry.hh" 44 #include "G4CrossSectionDataSetRegistry.hh"
>> 45 #include "G4DynamicParticle.hh"
42 #include "G4Material.hh" 46 #include "G4Material.hh"
43 #include "G4Element.hh" 47 #include "G4Element.hh"
44 #include "G4Isotope.hh" 48 #include "G4Isotope.hh"
45 #include "G4NistManager.hh" 49 #include "G4NistManager.hh"
>> 50 #include "G4HadronicException.hh"
>> 51 #include "G4HadTmpUtil.hh"
46 #include "Randomize.hh" 52 #include "Randomize.hh"
47 #include "G4HadronicParameters.hh" <<
48 53
49 G4VCrossSectionDataSet::G4VCrossSectionDataSet 54 G4VCrossSectionDataSet::G4VCrossSectionDataSet(const G4String& nam) :
50 name(nam), << 55 verboseLevel(0),minKinEnergy(0.0),maxKinEnergy(100*TeV),name(nam)
51 maxKinEnergy(G4HadronicParameters::Instance( <<
52 { 56 {
53 registry = G4CrossSectionDataSetRegistry::In << 57 G4CrossSectionDataSetRegistry::Instance()->Register(this);
54 registry->Register(this); <<
55 } 58 }
56 59
57 G4VCrossSectionDataSet::~G4VCrossSectionDataSe 60 G4VCrossSectionDataSet::~G4VCrossSectionDataSet()
58 { 61 {
59 registry->DeRegister(this); << 62 G4CrossSectionDataSetRegistry::Instance()->DeRegister(this);
60 } 63 }
61 64
62 G4bool 65 G4bool
63 G4VCrossSectionDataSet::IsElementApplicable(co 66 G4VCrossSectionDataSet::IsElementApplicable(const G4DynamicParticle*,
64 G4int, 67 G4int,
65 const G4Material*) 68 const G4Material*)
66 { 69 {
67 return false; 70 return false;
68 } 71 }
69 72
70 G4bool 73 G4bool
71 G4VCrossSectionDataSet::IsIsoApplicable(const 74 G4VCrossSectionDataSet::IsIsoApplicable(const G4DynamicParticle*,
72 G4int, 75 G4int, G4int,
73 const 76 const G4Element*,
74 const 77 const G4Material*)
75 { 78 {
76 return false; 79 return false;
77 } 80 }
78 81
79 G4double 82 G4double
80 G4VCrossSectionDataSet::ComputeCrossSection(co 83 G4VCrossSectionDataSet::ComputeCrossSection(const G4DynamicParticle* part,
81 co << 84 const G4Element* elm,
82 co << 85 const G4Material* mat)
83 { 86 {
84 G4int Z = elm->GetZasInt(); << 87 G4int Z = G4lrint(elm->GetZ());
85 88
86 if (IsElementApplicable(part, Z, mat)) { 89 if (IsElementApplicable(part, Z, mat)) {
87 return GetElementCrossSection(part, Z, mat 90 return GetElementCrossSection(part, Z, mat);
88 } 91 }
89 92
90 // isotope-wise cross section making sum ove 93 // isotope-wise cross section making sum over available
91 // isotope cross sections, which may be inco 94 // isotope cross sections, which may be incomplete, so
92 // the result is corrected 95 // the result is corrected
93 std::size_t nIso = elm->GetNumberOfIsotopes( << 96 G4int nIso = elm->GetNumberOfIsotopes();
94 G4double fact = 0.0; 97 G4double fact = 0.0;
95 G4double xsec = 0.0; 98 G4double xsec = 0.0;
>> 99 G4Isotope* iso = 0;
96 100
97 // user-defined isotope abundances << 101 if (0 < nIso) {
98 const G4IsotopeVector* isoVector = elm->GetI << 102
99 const G4double* abundVector = elm->GetRelati << 103 // user-defined isotope abundances
100 for (std::size_t j=0; j<nIso; ++j) { << 104 G4IsotopeVector* isoVector = elm->GetIsotopeVector();
101 const G4Isotope* iso = (*isoVector)[j]; << 105 G4double* abundVector = elm->GetRelativeAbundanceVector();
102 G4int A = iso->GetN(); << 106
103 if(abundVector[j] > 0.0 && IsIsoApplicable << 107 for (G4int j = 0; j<nIso; ++j) {
104 fact += abundVector[j]; << 108 iso = (*isoVector)[j];
105 xsec += abundVector[j]*GetIsoCrossSectio << 109 G4int A = iso->GetN();
>> 110 if(abundVector[j] > 0.0 && IsIsoApplicable(part, Z, A, elm, mat)) {
>> 111 fact += abundVector[j];
>> 112 xsec += abundVector[j]*GetIsoCrossSection(part, Z, A, iso, elm, mat);
>> 113 }
106 } 114 }
107 } <<
108 return (fact > 0.0) ? xsec/fact : 0.0; <<
109 } <<
110 115
111 G4double << 116 } else {
112 G4VCrossSectionDataSet::ComputeCrossSectionPer << 117
113 G4double kinEnergy, G4 << 118 // natural isotope abundances
114 const G4ParticleDefini << 119 G4NistManager* nist = G4NistManager::Instance();
115 const G4Element* elm, << 120 G4int n0 = nist->GetNistFirstIsotopeN(Z);
116 { << 121 G4int nn = nist->GetNumberOfNistIsotopes(Z);
117 G4int Z = elm->GetZasInt(); << 122 for (G4int A = n0; A < n0+nn; ++A) {
118 std::size_t nIso = elm->GetNumberOfIsotopes( << 123 G4double abund = nist->GetIsotopeAbundance(Z, A);
119 G4double xsec = 0.0; << 124 if(abund > 0.0 && IsIsoApplicable(part, Z, A, elm, mat)) {
120 const G4IsotopeVector* isoVector = elm->GetI << 125 fact += abund;
121 const G4double* abundVector = elm->GetRelati << 126 xsec += abund*GetIsoCrossSection(part, Z, A, iso, elm, mat);
122 for (std::size_t j=0; j<nIso; ++j) { << 127 }
123 const G4Isotope* iso = (*isoVector)[j]; << 128 }
124 G4int A = iso->GetN(); <<
125 xsec += abundVector[j]* <<
126 ComputeIsoCrossSection(kinEnergy, loge, <<
127 } 129 }
>> 130 if(fact > 0.0) { xsec /= fact; }
128 return xsec; 131 return xsec;
129 } 132 }
130 133
131 G4double 134 G4double
132 G4VCrossSectionDataSet::GetElementCrossSection 135 G4VCrossSectionDataSet::GetElementCrossSection(const G4DynamicParticle* dynPart,
133 G4int Z, 136 G4int Z,
134 const G4Material* mat) 137 const G4Material* mat)
135 { 138 {
136 G4ExceptionDescription ed; << 139 G4cout << "G4VCrossSectionDataSet::GetCrossSection per element ERROR: "
137 ed << "GetElementCrossSection is not impleme << 140 << " there is no cross section for "
138 << "Particle: " << dynPart->GetDefinition << 141 << dynPart->GetDefinition()->GetParticleName()
139 << " Ekin(MeV)= " << dynPart->GetKineti << 142 << " E(MeV)= " << dynPart->GetKineticEnergy()/MeV;
140 if(nullptr != mat) { ed << " material: " << << 143 if(mat) { G4cout << " inside " << mat->GetName(); }
141 ed << " target Z= " << Z << G4endl; << 144 G4cout << " for Z= " << Z << G4endl;
142 G4Exception("G4VCrossSectionDataSet::GetElem << 145 throw G4HadronicException(__FILE__, __LINE__,
143 FatalException, ed); << 146 "G4VCrossSectionDataSet::GetElementCrossSection is absent");
144 return 0.0; 147 return 0.0;
145 } 148 }
146 149
147 G4double 150 G4double
148 G4VCrossSectionDataSet::GetIsoCrossSection(con 151 G4VCrossSectionDataSet::GetIsoCrossSection(const G4DynamicParticle* dynPart,
149 G4int Z, G4int A, 152 G4int Z, G4int A,
150 const G4Isotope*, 153 const G4Isotope*,
151 const G4Element* elm, 154 const G4Element* elm,
152 const G4Material* mat) 155 const G4Material* mat)
153 { 156 {
154 G4ExceptionDescription ed; << 157 G4cout << "G4VCrossSectionDataSet::GetCrossSection per isotope ERROR: "
155 ed << "GetIsoCrossSection is not implemented << 158 << " there is no cross section for "
156 << "Particle: " << dynPart->GetDefinition << 159 << dynPart->GetDefinition()->GetParticleName()
157 << " Ekin(MeV)= " << dynPart->GetKineti << 160 << " E(MeV)= " << dynPart->GetKineticEnergy()/MeV;
158 if(nullptr != mat) { ed << " material: " << << 161 if(mat) { G4cout << " inside " << mat->GetName(); }
159 if(nullptr != elm) { ed << " element: " << e << 162 if(elm) { G4cout << " for " << elm->GetName(); }
160 ed << " target Z= " << Z << " A= " << A << G << 163 G4cout << " Z= " << Z << " A= " << A << G4endl;
161 G4Exception("G4VCrossSectionDataSet::GetIsoC << 164 throw G4HadronicException(__FILE__, __LINE__,
162 FatalException, ed); << 165 "G4VCrossSectionDataSet::GetIsoCrossSection is absent");
163 return 0.0; <<
164 } <<
165 <<
166 G4double <<
167 G4VCrossSectionDataSet::ComputeIsoCrossSection <<
168 const G4ParticleDefinition* p <<
169 G4int Z, G4int A, <<
170 const G4Isotope*, <<
171 const G4Element* elm, <<
172 const G4Material* mat) <<
173 { <<
174 G4ExceptionDescription ed; <<
175 ed << "GetIsoCrossSection is not implemented <<
176 << "Particle: " << pd->GetParticleName() <<
177 << " Ekin(MeV)= " << kinEnergy/CLHEP::M <<
178 if(nullptr != mat) { ed << " material: " << <<
179 if(nullptr != elm) { ed << " element: " << e <<
180 ed << " target Z= " << Z << " A= " << A << G <<
181 G4Exception("G4VCrossSectionDataSet::GetIsoC <<
182 FatalException, ed); <<
183 return 0.0; 166 return 0.0;
184 } 167 }
185 168
186 const G4Isotope* << 169 G4Isotope*
187 G4VCrossSectionDataSet::SelectIsotope(const G4 << 170 G4VCrossSectionDataSet::SelectIsotope(const G4Element* anElement, G4double)
188 G4double <<
189 { 171 {
190 G4int nIso = (G4int)anElement->GetNumberOfIs << 172 G4int nIso = anElement->GetNumberOfIsotopes();
191 const G4Isotope* iso = anElement->GetIsotope << 173 G4IsotopeVector* isoVector = anElement->GetIsotopeVector();
>> 174 G4Isotope* iso = (*isoVector)[0];
192 175
193 // more than 1 isotope 176 // more than 1 isotope
194 if(1 < nIso) { 177 if(1 < nIso) {
195 const G4double* abundVector = anElement->G << 178 G4double* abundVector = anElement->GetRelativeAbundanceVector();
196 G4double sum = 0.0; 179 G4double sum = 0.0;
197 G4double q = G4UniformRand(); 180 G4double q = G4UniformRand();
198 for (G4int j=0; j<nIso; ++j) { << 181 for (G4int j = 0; j<nIso; ++j) {
199 sum += abundVector[j]; 182 sum += abundVector[j];
200 if(q <= sum) { 183 if(q <= sum) {
201 iso = anElement->GetIsotope(j); << 184 iso = (*isoVector)[j];
202 break; 185 break;
203 } 186 }
204 } 187 }
205 } 188 }
206 return iso; 189 return iso;
207 } 190 }
208 191
209 void G4VCrossSectionDataSet::BuildPhysicsTable 192 void G4VCrossSectionDataSet::BuildPhysicsTable(const G4ParticleDefinition&)
210 {} 193 {}
211 194
212 void G4VCrossSectionDataSet::DumpPhysicsTable( 195 void G4VCrossSectionDataSet::DumpPhysicsTable(const G4ParticleDefinition&)
213 {} 196 {}
214 197
215 void G4VCrossSectionDataSet::CrossSectionDescr 198 void G4VCrossSectionDataSet::CrossSectionDescription(std::ostream& outFile) const
216 { 199 {
217 outFile << "The description for this cross s 200 outFile << "The description for this cross section data set has not been written yet.\n";
218 } 201 }
219 202