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 | // |
---|
27 | // $Id: G4QIsotope.cc,v 1.10 2007/11/28 14:14:36 mkossov Exp $ |
---|
28 | // GEANT4 tag $Name: geant4-09-02-ref-02 $ |
---|
29 | // |
---|
30 | // ---------------- G4QIsotope class ---------------- |
---|
31 | // by Mikhail Kossov, December 2003. |
---|
32 | // class G4QIsotope gives Isotopes for Elements (CHIPS branch of G4) |
---|
33 | // ------------------------------------------------------------------ |
---|
34 | // **************************************************************************************** |
---|
35 | // ********** This CLASS is temporary moved from the photolepton_hadron directory ********* |
---|
36 | // ******* DO NOT MAKE ANY CHANGE! With time it'll move back to photolepton...(M.K.) ****** |
---|
37 | // **************************************************************************************** |
---|
38 | // |
---|
39 | // 1 2 3 4 5 6 7 8 9 |
---|
40 | //34567890123456789012345678901234567890123456789012345678901234567890123456789012345678901 |
---|
41 | |
---|
42 | //#define debug |
---|
43 | //#define cdebug |
---|
44 | //#define pdebug |
---|
45 | //#define ppdebug |
---|
46 | //#define sdebug |
---|
47 | |
---|
48 | #include "G4QIsotope.hh" |
---|
49 | #include <cmath> |
---|
50 | using namespace std; |
---|
51 | |
---|
52 | vector<vector<pair<G4int,G4double>*>*>G4QIsotope::natElements;//NaturalElems |
---|
53 | vector<vector<pair<G4int,G4double>*>*>G4QIsotope::natSumAbund;//NaturalSumAb |
---|
54 | vector<vector<pair<G4int,G4double>*>*>G4QIsotope::natIsoCrosS;//CSOfNatElems |
---|
55 | vector<pair<G4int,vector<pair<G4int,G4double>*>*>*>G4QIsotope::newElems; |
---|
56 | vector<pair<G4int,vector<pair<G4int,G4double>*>*>*>G4QIsotope::newSumAb; |
---|
57 | vector<pair<G4int,vector<pair<G4int,G4double>*>*>*>G4QIsotope::newIsoCS; |
---|
58 | |
---|
59 | G4QIsotope::G4QIsotope() |
---|
60 | { |
---|
61 | #ifdef cdebug |
---|
62 | G4cout<<"G4QIsotope::Constructor is called"<<G4endl; |
---|
63 | #endif |
---|
64 | vector<vector<pair<G4int,G4double> >*> natEl; |
---|
65 | #ifdef cdebug |
---|
66 | G4cout<<"G4QIsotope::Constructor natEl is booked"<<G4endl; |
---|
67 | #endif |
---|
68 | vector<pair<G4int,G4double> >*a0=new vector<pair<G4int,G4double> >; |
---|
69 | #ifdef cdebug |
---|
70 | G4cout<<"G4QIsotope::Constructor a0 is booked"<<G4endl; |
---|
71 | #endif |
---|
72 | a0->push_back(make_pair(1,1.)); |
---|
73 | #ifdef cdebug |
---|
74 | G4cout<<"G4QIsotope::Constructor a0 is filled by a pair"<<G4endl; |
---|
75 | #endif |
---|
76 | natEl.push_back(a0); |
---|
77 | #ifdef cdebug |
---|
78 | G4cout<<"G4QIsotope::Constructor a0 is filled in natEl"<<G4endl; |
---|
79 | #endif |
---|
80 | // If an error is found in this initialization, please, correct the simple tree below |
---|
81 | vector<pair<G4int,G4double> >*a1=new vector<pair<G4int,G4double> >; |
---|
82 | a1->push_back(make_pair(0,.99985)); |
---|
83 | a1->push_back(make_pair(1,1.)); |
---|
84 | natEl.push_back(a1); |
---|
85 | vector<pair<G4int,G4double> >*a2=new vector<pair<G4int,G4double> >; |
---|
86 | a2->push_back(make_pair(2,.999999863)); |
---|
87 | a2->push_back(make_pair(1,1.)); |
---|
88 | natEl.push_back(a2); |
---|
89 | vector<pair<G4int,G4double> >*a3=new vector<pair<G4int,G4double> >; |
---|
90 | a3->push_back(make_pair(4,.925)); |
---|
91 | a3->push_back(make_pair(3,1.)); |
---|
92 | natEl.push_back(a3); |
---|
93 | vector<pair<G4int,G4double> >*a4=new vector<pair<G4int,G4double> >; |
---|
94 | a4->push_back(make_pair(5,1.)); |
---|
95 | natEl.push_back(a4); |
---|
96 | vector<pair<G4int,G4double> >*a5=new vector<pair<G4int,G4double> >; |
---|
97 | a5->push_back(make_pair(6,.801)); |
---|
98 | a5->push_back(make_pair(5,1.)); |
---|
99 | natEl.push_back(a5); |
---|
100 | vector<pair<G4int,G4double> >*a6=new vector<pair<G4int,G4double> >; |
---|
101 | a6->push_back(make_pair(6,.989)); |
---|
102 | a6->push_back(make_pair(7,1.)); |
---|
103 | natEl.push_back(a6); |
---|
104 | vector<pair<G4int,G4double> >*a7=new vector<pair<G4int,G4double> >; |
---|
105 | a7->push_back(make_pair(7,.9963)); |
---|
106 | a7->push_back(make_pair(8,1.)); |
---|
107 | natEl.push_back(a7); |
---|
108 | vector<pair<G4int,G4double> >*a8=new vector<pair<G4int,G4double> >; |
---|
109 | a8->push_back(make_pair(8,.9976)); |
---|
110 | a8->push_back(make_pair(10,.9996)); |
---|
111 | a8->push_back(make_pair(9,1.)); |
---|
112 | natEl.push_back(a8); |
---|
113 | vector<pair<G4int,G4double> >*a9=new vector<pair<G4int,G4double> >; |
---|
114 | a9->push_back(make_pair(10,1.)); |
---|
115 | natEl.push_back(a9); |
---|
116 | vector<pair<G4int,G4double> >*b0=new vector<pair<G4int,G4double> >; |
---|
117 | b0->push_back(make_pair(10,.9948)); |
---|
118 | b0->push_back(make_pair(11,.9975)); |
---|
119 | b0->push_back(make_pair(12,1.)); |
---|
120 | natEl.push_back(b0); |
---|
121 | vector<pair<G4int,G4double> >*b1=new vector<pair<G4int,G4double> >; |
---|
122 | b1->push_back(make_pair(12,1.)); |
---|
123 | natEl.push_back(b1); |
---|
124 | vector<pair<G4int,G4double> >*b2=new vector<pair<G4int,G4double> >; |
---|
125 | b2->push_back(make_pair(12,.7899)); |
---|
126 | b2->push_back(make_pair(13,.8899)); |
---|
127 | b2->push_back(make_pair(14,1.)); |
---|
128 | natEl.push_back(b2); |
---|
129 | vector<pair<G4int,G4double> >*b3=new vector<pair<G4int,G4double> >; |
---|
130 | b3->push_back(make_pair(14,1.)); |
---|
131 | natEl.push_back(b3); |
---|
132 | vector<pair<G4int,G4double> >*b4=new vector<pair<G4int,G4double> >; |
---|
133 | b4->push_back(make_pair(14,.9223)); |
---|
134 | b4->push_back(make_pair(15,.969)); |
---|
135 | b4->push_back(make_pair(16,1.)); |
---|
136 | natEl.push_back(b4); |
---|
137 | vector<pair<G4int,G4double> >*b5=new vector<pair<G4int,G4double> >; |
---|
138 | b5->push_back(make_pair(16,1.)); |
---|
139 | natEl.push_back(b5); |
---|
140 | vector<pair<G4int,G4double> >*b6=new vector<pair<G4int,G4double> >; |
---|
141 | b6->push_back(make_pair(16,.9502)); |
---|
142 | b6->push_back(make_pair(18,.9923)); |
---|
143 | b6->push_back(make_pair(17,.9998)); |
---|
144 | b6->push_back(make_pair(20,1.)); |
---|
145 | natEl.push_back(b6); |
---|
146 | vector<pair<G4int,G4double> >*b7=new vector<pair<G4int,G4double> >; |
---|
147 | b7->push_back(make_pair(18,.7577)); |
---|
148 | b7->push_back(make_pair(20,1.)); |
---|
149 | natEl.push_back(b7); |
---|
150 | vector<pair<G4int,G4double> >*b8=new vector<pair<G4int,G4double> >; |
---|
151 | b8->push_back(make_pair(22,.996)); |
---|
152 | b8->push_back(make_pair(18,.99937)); |
---|
153 | b8->push_back(make_pair(20,1.)); |
---|
154 | natEl.push_back(b8); |
---|
155 | vector<pair<G4int,G4double> >*b9=new vector<pair<G4int,G4double> >; |
---|
156 | b9->push_back(make_pair(20,.932581)); |
---|
157 | b9->push_back(make_pair(22,.999883)); |
---|
158 | b9->push_back(make_pair(21,1.)); |
---|
159 | natEl.push_back(b9); |
---|
160 | vector<pair<G4int,G4double> >*c0=new vector<pair<G4int,G4double> >; |
---|
161 | c0->push_back(make_pair(20,.96941)); |
---|
162 | c0->push_back(make_pair(24,.99027)); |
---|
163 | c0->push_back(make_pair(22,.99674)); |
---|
164 | c0->push_back(make_pair(28,.99861)); |
---|
165 | c0->push_back(make_pair(23,.99996)); |
---|
166 | c0->push_back(make_pair(26,1.)); |
---|
167 | natEl.push_back(c0); |
---|
168 | vector<pair<G4int,G4double> >*c1=new vector<pair<G4int,G4double> >; |
---|
169 | c1->push_back(make_pair(24,1.)); |
---|
170 | natEl.push_back(c1); |
---|
171 | vector<pair<G4int,G4double> >*c2=new vector<pair<G4int,G4double> >; |
---|
172 | c2->push_back(make_pair(1,.738)); |
---|
173 | c2->push_back(make_pair(1,.818)); |
---|
174 | c2->push_back(make_pair(1,.891)); |
---|
175 | c2->push_back(make_pair(1,.946)); |
---|
176 | c2->push_back(make_pair(1,1.)); |
---|
177 | natEl.push_back(c2); |
---|
178 | vector<pair<G4int,G4double> >*c3=new vector<pair<G4int,G4double> >; |
---|
179 | c3->push_back(make_pair(28,.9975)); |
---|
180 | c3->push_back(make_pair(27,1.)); |
---|
181 | natEl.push_back(c3); |
---|
182 | vector<pair<G4int,G4double> >*c4=new vector<pair<G4int,G4double> >; |
---|
183 | c4->push_back(make_pair(28,.8379)); |
---|
184 | c4->push_back(make_pair(29,.9329)); |
---|
185 | c4->push_back(make_pair(26,.97635)); |
---|
186 | c4->push_back(make_pair(30,1.)); |
---|
187 | natEl.push_back(c4); |
---|
188 | vector<pair<G4int,G4double> >*c5=new vector<pair<G4int,G4double> >; |
---|
189 | c5->push_back(make_pair(30,1.)); |
---|
190 | natEl.push_back(c5); |
---|
191 | vector<pair<G4int,G4double> >*c6=new vector<pair<G4int,G4double> >; |
---|
192 | c6->push_back(make_pair(30,.9172)); |
---|
193 | c6->push_back(make_pair(28,.9762)); |
---|
194 | c6->push_back(make_pair(31,.9972)); |
---|
195 | c6->push_back(make_pair(32,1.)); |
---|
196 | natEl.push_back(c6); |
---|
197 | vector<pair<G4int,G4double> >*c7=new vector<pair<G4int,G4double> >; |
---|
198 | c7->push_back(make_pair(32,1.)); |
---|
199 | natEl.push_back(c7); |
---|
200 | vector<pair<G4int,G4double> >*c8=new vector<pair<G4int,G4double> >; |
---|
201 | c8->push_back(make_pair(30,.68077)); |
---|
202 | c8->push_back(make_pair(32,.943)); |
---|
203 | c8->push_back(make_pair(34,.97934)); |
---|
204 | c8->push_back(make_pair(33,.99074)); |
---|
205 | c8->push_back(make_pair(36,1.)); |
---|
206 | natEl.push_back(c8); |
---|
207 | vector<pair<G4int,G4double> >*c9=new vector<pair<G4int,G4double> >; |
---|
208 | c9->push_back(make_pair(34,.6917)); |
---|
209 | c9->push_back(make_pair(36,1.)); |
---|
210 | natEl.push_back(c9); |
---|
211 | vector<pair<G4int,G4double> >*d0=new vector<pair<G4int,G4double> >; |
---|
212 | d0->push_back(make_pair(34,.486)); |
---|
213 | d0->push_back(make_pair(36,.765)); |
---|
214 | d0->push_back(make_pair(38,.953)); |
---|
215 | d0->push_back(make_pair(37,.994)); |
---|
216 | d0->push_back(make_pair(40,1.)); |
---|
217 | natEl.push_back(d0); |
---|
218 | vector<pair<G4int,G4double> >*d1=new vector<pair<G4int,G4double> >; |
---|
219 | d1->push_back(make_pair(38,.60108)); |
---|
220 | d1->push_back(make_pair(40,1.)); |
---|
221 | natEl.push_back(d1); |
---|
222 | vector<pair<G4int,G4double> >*d2=new vector<pair<G4int,G4double> >; |
---|
223 | d2->push_back(make_pair(42,.3594)); |
---|
224 | d2->push_back(make_pair(40,.6360)); |
---|
225 | d2->push_back(make_pair(38,.8484)); |
---|
226 | d2->push_back(make_pair(41,.9256)); |
---|
227 | d2->push_back(make_pair(44,1.)); |
---|
228 | natEl.push_back(d2); |
---|
229 | vector<pair<G4int,G4double> >*d3=new vector<pair<G4int,G4double> >; |
---|
230 | d3->push_back(make_pair(42,1.)); |
---|
231 | natEl.push_back(d3); |
---|
232 | vector<pair<G4int,G4double> >*d4=new vector<pair<G4int,G4double> >; |
---|
233 | d4->push_back(make_pair(46,.4961)); |
---|
234 | d4->push_back(make_pair(44,.7378)); |
---|
235 | d4->push_back(make_pair(42,.8274)); |
---|
236 | d4->push_back(make_pair(48,.9148)); |
---|
237 | d4->push_back(make_pair(43,.9911)); |
---|
238 | d4->push_back(make_pair(40,1.)); |
---|
239 | natEl.push_back(d4); |
---|
240 | vector<pair<G4int,G4double> >*d5=new vector<pair<G4int,G4double> >; |
---|
241 | d5->push_back(make_pair(44,.5069)); |
---|
242 | d5->push_back(make_pair(46,1.)); |
---|
243 | natEl.push_back(d5); |
---|
244 | vector<pair<G4int,G4double> >*d6=new vector<pair<G4int,G4double> >; |
---|
245 | d6->push_back(make_pair(48,.57)); |
---|
246 | d6->push_back(make_pair(50,.743)); |
---|
247 | d6->push_back(make_pair(46,.859)); |
---|
248 | d6->push_back(make_pair(47,.974)); |
---|
249 | d6->push_back(make_pair(44,.9965)); |
---|
250 | d6->push_back(make_pair(42,1.)); |
---|
251 | natEl.push_back(d6); |
---|
252 | vector<pair<G4int,G4double> >*d7=new vector<pair<G4int,G4double> >; |
---|
253 | d7->push_back(make_pair(48,.7217)); |
---|
254 | d7->push_back(make_pair(50,1.)); |
---|
255 | natEl.push_back(d7); |
---|
256 | vector<pair<G4int,G4double> >*d8=new vector<pair<G4int,G4double> >; |
---|
257 | d8->push_back(make_pair(50,.8258)); |
---|
258 | d8->push_back(make_pair(48,.9244)); |
---|
259 | d8->push_back(make_pair(49,.9944)); |
---|
260 | d8->push_back(make_pair(46,1.)); |
---|
261 | natEl.push_back(d8); |
---|
262 | vector<pair<G4int,G4double> >*d9=new vector<pair<G4int,G4double> >; |
---|
263 | d9->push_back(make_pair(50,1.)); |
---|
264 | natEl.push_back(d9); |
---|
265 | vector<pair<G4int,G4double> >*e0=new vector<pair<G4int,G4double> >; |
---|
266 | e0->push_back(make_pair(50,.5145)); |
---|
267 | e0->push_back(make_pair(54,.6883)); |
---|
268 | e0->push_back(make_pair(53,.8598)); |
---|
269 | e0->push_back(make_pair(51,.972)); |
---|
270 | e0->push_back(make_pair(56,1.)); |
---|
271 | natEl.push_back(e0); |
---|
272 | vector<pair<G4int,G4double> >*e1=new vector<pair<G4int,G4double> >; |
---|
273 | e1->push_back(make_pair(52,1.)); |
---|
274 | natEl.push_back(e1); |
---|
275 | vector<pair<G4int,G4double> >*e2=new vector<pair<G4int,G4double> >; |
---|
276 | e2->push_back(make_pair(56,.2413)); |
---|
277 | e2->push_back(make_pair(54,.4081)); |
---|
278 | e2->push_back(make_pair(53,.5673)); |
---|
279 | e2->push_back(make_pair(50,.7157)); |
---|
280 | e2->push_back(make_pair(58,.8120)); |
---|
281 | e2->push_back(make_pair(55,.9075)); |
---|
282 | e2->push_back(make_pair(52,1.)); |
---|
283 | natEl.push_back(e2); |
---|
284 | vector<pair<G4int,G4double> >*e3=new vector<pair<G4int,G4double> >; |
---|
285 | e3->push_back(make_pair(55,1.)); |
---|
286 | natEl.push_back(e3); |
---|
287 | vector<pair<G4int,G4double> >*e4=new vector<pair<G4int,G4double> >; |
---|
288 | e4->push_back(make_pair(58,.316)); |
---|
289 | e4->push_back(make_pair(60,.502)); |
---|
290 | e4->push_back(make_pair(57,.673)); |
---|
291 | e4->push_back(make_pair(55,.8)); |
---|
292 | e4->push_back(make_pair(56,.926)); |
---|
293 | e4->push_back(make_pair(52,.9814)); |
---|
294 | e4->push_back(make_pair(54,1.)); |
---|
295 | natEl.push_back(e4); |
---|
296 | vector<pair<G4int,G4double> >*e5=new vector<pair<G4int,G4double> >; |
---|
297 | e5->push_back(make_pair(58,1.)); |
---|
298 | natEl.push_back(e5); |
---|
299 | vector<pair<G4int,G4double> >*e6=new vector<pair<G4int,G4double> >; |
---|
300 | e6->push_back(make_pair(60,.2733)); |
---|
301 | e6->push_back(make_pair(62,.5379)); |
---|
302 | e6->push_back(make_pair(59,.7612)); |
---|
303 | e6->push_back(make_pair(55,.8784)); |
---|
304 | e6->push_back(make_pair(58,.9898)); |
---|
305 | e6->push_back(make_pair(56,1.)); |
---|
306 | natEl.push_back(e6); |
---|
307 | vector<pair<G4int,G4double> >*e7=new vector<pair<G4int,G4double> >; |
---|
308 | e7->push_back(make_pair(60,.51839)); |
---|
309 | e7->push_back(make_pair(62,1.)); |
---|
310 | natEl.push_back(e7); |
---|
311 | vector<pair<G4int,G4double> >*e8=new vector<pair<G4int,G4double> >; |
---|
312 | e8->push_back(make_pair(66,.2873)); |
---|
313 | e8->push_back(make_pair(64,.5286)); |
---|
314 | e8->push_back(make_pair(59,.6566)); |
---|
315 | e8->push_back(make_pair(62,.7815)); |
---|
316 | e8->push_back(make_pair(65,.9037)); |
---|
317 | e8->push_back(make_pair(68,.9786)); |
---|
318 | e8->push_back(make_pair(58,.9911)); |
---|
319 | e8->push_back(make_pair(60,1.)); |
---|
320 | natEl.push_back(e8); |
---|
321 | vector<pair<G4int,G4double> >*e9=new vector<pair<G4int,G4double> >; |
---|
322 | e9->push_back(make_pair(66,.9577)); |
---|
323 | e9->push_back(make_pair(64,1.)); |
---|
324 | natEl.push_back(e9); |
---|
325 | vector<pair<G4int,G4double> >*f0=new vector<pair<G4int,G4double> >; |
---|
326 | f0->push_back(make_pair(70,.3259)); |
---|
327 | f0->push_back(make_pair(68,.5681)); |
---|
328 | f0->push_back(make_pair(66,.7134)); |
---|
329 | f0->push_back(make_pair(69,.7992)); |
---|
330 | f0->push_back(make_pair(67,.8760)); |
---|
331 | f0->push_back(make_pair(74,.9339)); |
---|
332 | f0->push_back(make_pair(72,.9802)); |
---|
333 | f0->push_back(make_pair(62,.9899)); |
---|
334 | f0->push_back(make_pair(64,1.)); |
---|
335 | //f0->push_back(make_pair(64,.9964)); |
---|
336 | //f0->push_back(make_pair(65,1.)); // Nine isotopes is the maximum, so Sn155 is out |
---|
337 | natEl.push_back(f0); |
---|
338 | vector<pair<G4int,G4double> >*f1=new vector<pair<G4int,G4double> >; |
---|
339 | f1->push_back(make_pair(70,.5736)); |
---|
340 | f1->push_back(make_pair(72,1.)); |
---|
341 | natEl.push_back(f1); |
---|
342 | vector<pair<G4int,G4double> >*f2=new vector<pair<G4int,G4double> >; |
---|
343 | f2->push_back(make_pair(78,.3387)); |
---|
344 | f2->push_back(make_pair(76,.6557)); |
---|
345 | f2->push_back(make_pair(74,.8450)); |
---|
346 | f2->push_back(make_pair(73,.9162)); |
---|
347 | f2->push_back(make_pair(72,.9641)); |
---|
348 | f2->push_back(make_pair(70,.9900)); |
---|
349 | f2->push_back(make_pair(71,.99905)); |
---|
350 | f2->push_back(make_pair(68,1.)); |
---|
351 | natEl.push_back(f2); |
---|
352 | vector<pair<G4int,G4double> >*f3=new vector<pair<G4int,G4double> >; |
---|
353 | f3->push_back(make_pair(74,1.)); |
---|
354 | natEl.push_back(f3); |
---|
355 | vector<pair<G4int,G4double> >*f4=new vector<pair<G4int,G4double> >; |
---|
356 | f4->push_back(make_pair(78,.269)); |
---|
357 | f4->push_back(make_pair(75,.533)); |
---|
358 | f4->push_back(make_pair(77,.745)); |
---|
359 | f4->push_back(make_pair(80,.849)); |
---|
360 | f4->push_back(make_pair(82,.938)); |
---|
361 | f4->push_back(make_pair(76,.979)); |
---|
362 | f4->push_back(make_pair(74,.9981)); |
---|
363 | f4->push_back(make_pair(70,.9991)); |
---|
364 | f4->push_back(make_pair(72,1.)); |
---|
365 | natEl.push_back(f4); |
---|
366 | vector<pair<G4int,G4double> >*f5=new vector<pair<G4int,G4double> >; |
---|
367 | f5->push_back(make_pair(78,1.)); |
---|
368 | natEl.push_back(f5); |
---|
369 | vector<pair<G4int,G4double> >*f6=new vector<pair<G4int,G4double> >; |
---|
370 | f6->push_back(make_pair(82,.717)); |
---|
371 | f6->push_back(make_pair(81,.8293)); |
---|
372 | f6->push_back(make_pair(80,.9078)); |
---|
373 | f6->push_back(make_pair(79,.97373)); |
---|
374 | f6->push_back(make_pair(78,.99793)); |
---|
375 | f6->push_back(make_pair(74,.99899)); |
---|
376 | f6->push_back(make_pair(76,1.)); |
---|
377 | natEl.push_back(f6); |
---|
378 | vector<pair<G4int,G4double> >*f7=new vector<pair<G4int,G4double> >; |
---|
379 | f7->push_back(make_pair(82,.999098)); |
---|
380 | f7->push_back(make_pair(81,1.)); |
---|
381 | natEl.push_back(f7); |
---|
382 | vector<pair<G4int,G4double> >*f8=new vector<pair<G4int,G4double> >; |
---|
383 | f8->push_back(make_pair(82,.8843)); |
---|
384 | f8->push_back(make_pair(84,.9956)); |
---|
385 | f8->push_back(make_pair(80,.9981)); |
---|
386 | f8->push_back(make_pair(78,1.)); |
---|
387 | natEl.push_back(f8); |
---|
388 | vector<pair<G4int,G4double> >*f9=new vector<pair<G4int,G4double> >; |
---|
389 | f9->push_back(make_pair(82,1.)); |
---|
390 | natEl.push_back(f9); |
---|
391 | vector<pair<G4int,G4double> >*g0=new vector<pair<G4int,G4double> >; |
---|
392 | g0->push_back(make_pair(82,.2713)); |
---|
393 | g0->push_back(make_pair(84,.5093)); |
---|
394 | g0->push_back(make_pair(86,.6812)); |
---|
395 | g0->push_back(make_pair(83,.8030)); |
---|
396 | g0->push_back(make_pair(85,.8860)); |
---|
397 | g0->push_back(make_pair(88,.9436)); |
---|
398 | g0->push_back(make_pair(90,1.)); |
---|
399 | natEl.push_back(g0); |
---|
400 | vector<pair<G4int,G4double> >*g1=new vector<pair<G4int,G4double> >; |
---|
401 | g1->push_back(make_pair(85,1.)); |
---|
402 | natEl.push_back(g1); |
---|
403 | vector<pair<G4int,G4double> >*g2=new vector<pair<G4int,G4double> >; |
---|
404 | g2->push_back(make_pair(90,.267)); |
---|
405 | g2->push_back(make_pair(92,.494)); |
---|
406 | g2->push_back(make_pair(85,.644)); |
---|
407 | g2->push_back(make_pair(87,.782)); |
---|
408 | g2->push_back(make_pair(86,.895)); |
---|
409 | g2->push_back(make_pair(88,.969)); |
---|
410 | g2->push_back(make_pair(82,1.)); |
---|
411 | natEl.push_back(g2); |
---|
412 | vector<pair<G4int,G4double> >*g3=new vector<pair<G4int,G4double> >; |
---|
413 | g3->push_back(make_pair(90,.522)); |
---|
414 | g3->push_back(make_pair(89,1.)); |
---|
415 | natEl.push_back(g3); |
---|
416 | vector<pair<G4int,G4double> >*g4=new vector<pair<G4int,G4double> >; |
---|
417 | g4->push_back(make_pair(94,.2484)); |
---|
418 | g4->push_back(make_pair(96,.4670)); |
---|
419 | g4->push_back(make_pair(92,.6717)); |
---|
420 | g4->push_back(make_pair(93,.8282)); |
---|
421 | g4->push_back(make_pair(91,.9762)); |
---|
422 | g4->push_back(make_pair(90,.9980)); |
---|
423 | g4->push_back(make_pair(88,1.)); |
---|
424 | natEl.push_back(g4); |
---|
425 | vector<pair<G4int,G4double> >*g5=new vector<pair<G4int,G4double> >; |
---|
426 | g5->push_back(make_pair(94,1.)); |
---|
427 | natEl.push_back(g5); |
---|
428 | vector<pair<G4int,G4double> >*g6=new vector<pair<G4int,G4double> >; |
---|
429 | g6->push_back(make_pair(98,.282)); |
---|
430 | g6->push_back(make_pair(96,.537)); |
---|
431 | g6->push_back(make_pair(97,.786)); |
---|
432 | g6->push_back(make_pair(95,.975)); |
---|
433 | g6->push_back(make_pair(94,.9984)); |
---|
434 | g6->push_back(make_pair(92,.9994)); |
---|
435 | g6->push_back(make_pair(90,1.)); |
---|
436 | natEl.push_back(g6); |
---|
437 | vector<pair<G4int,G4double> >*g7=new vector<pair<G4int,G4double> >; |
---|
438 | g7->push_back(make_pair(98,1.)); |
---|
439 | natEl.push_back(g7); |
---|
440 | vector<pair<G4int,G4double> >*g8=new vector<pair<G4int,G4double> >; |
---|
441 | g8->push_back(make_pair( 98,.3360)); |
---|
442 | g8->push_back(make_pair(100,.6040)); |
---|
443 | g8->push_back(make_pair( 99,.8335)); |
---|
444 | g8->push_back(make_pair(102,.9825)); |
---|
445 | g8->push_back(make_pair( 96,.9986)); |
---|
446 | g8->push_back(make_pair( 94,1.)); |
---|
447 | natEl.push_back(g8); |
---|
448 | vector<pair<G4int,G4double> >*g9=new vector<pair<G4int,G4double> >; |
---|
449 | g9->push_back(make_pair(100,1.)); |
---|
450 | natEl.push_back(g9); |
---|
451 | vector<pair<G4int,G4double> >*h0=new vector<pair<G4int,G4double> >; |
---|
452 | h0->push_back(make_pair(104,.3180)); |
---|
453 | h0->push_back(make_pair(102,.5370)); |
---|
454 | h0->push_back(make_pair(103,.6982)); |
---|
455 | h0->push_back(make_pair(101,.8412)); |
---|
456 | h0->push_back(make_pair(106,.9682)); |
---|
457 | h0->push_back(make_pair(100,.9987)); |
---|
458 | h0->push_back(make_pair( 98,1.)); |
---|
459 | natEl.push_back(h0); |
---|
460 | vector<pair<G4int,G4double> >*h1=new vector<pair<G4int,G4double> >; |
---|
461 | h1->push_back(make_pair(104,.9741)); |
---|
462 | h1->push_back(make_pair(105,1.)); |
---|
463 | natEl.push_back(h1); |
---|
464 | vector<pair<G4int,G4double> >*h2=new vector<pair<G4int,G4double> >; |
---|
465 | h2->push_back(make_pair(108,.35100)); |
---|
466 | h2->push_back(make_pair(106,.62397)); |
---|
467 | h2->push_back(make_pair(105,.81003)); |
---|
468 | h2->push_back(make_pair(107,.94632)); |
---|
469 | h2->push_back(make_pair(104,.99838)); |
---|
470 | h2->push_back(make_pair(102,1.)); |
---|
471 | natEl.push_back(h2); |
---|
472 | vector<pair<G4int,G4double> >*h3=new vector<pair<G4int,G4double> >; |
---|
473 | h3->push_back(make_pair(108,.99988)); |
---|
474 | h3->push_back(make_pair(107,1.)); |
---|
475 | natEl.push_back(h3); |
---|
476 | vector<pair<G4int,G4double> >*h4=new vector<pair<G4int,G4double> >; |
---|
477 | h4->push_back(make_pair(110,.307)); |
---|
478 | h4->push_back(make_pair(112,.593)); |
---|
479 | h4->push_back(make_pair(108,.856)); |
---|
480 | h4->push_back(make_pair(109,.9988)); |
---|
481 | h4->push_back(make_pair(106,1.)); |
---|
482 | natEl.push_back(h4); |
---|
483 | vector<pair<G4int,G4double> >*h5=new vector<pair<G4int,G4double> >; |
---|
484 | h5->push_back(make_pair(112,.626)); |
---|
485 | h5->push_back(make_pair(110,1.)); |
---|
486 | natEl.push_back(h5); |
---|
487 | vector<pair<G4int,G4double> >*h6=new vector<pair<G4int,G4double> >; |
---|
488 | h6->push_back(make_pair(116,.410)); |
---|
489 | h6->push_back(make_pair(114,.674)); |
---|
490 | h6->push_back(make_pair(113,.835)); |
---|
491 | h6->push_back(make_pair(112,.968)); |
---|
492 | h6->push_back(make_pair(111,.984)); |
---|
493 | h6->push_back(make_pair(110,.9998)); |
---|
494 | h6->push_back(make_pair(108,1.)); |
---|
495 | natEl.push_back(h6); |
---|
496 | vector<pair<G4int,G4double> >*h7=new vector<pair<G4int,G4double> >; |
---|
497 | h7->push_back(make_pair(116,.627)); |
---|
498 | h7->push_back(make_pair(114,1.)); |
---|
499 | natEl.push_back(h7); |
---|
500 | vector<pair<G4int,G4double> >*h8=new vector<pair<G4int,G4double> >; |
---|
501 | h8->push_back(make_pair(117,.338)); |
---|
502 | h8->push_back(make_pair(116,.667)); |
---|
503 | h8->push_back(make_pair(118,.920)); |
---|
504 | h8->push_back(make_pair(120,.992)); |
---|
505 | h8->push_back(make_pair(114,.9999)); |
---|
506 | h8->push_back(make_pair(112,1.)); |
---|
507 | natEl.push_back(h8); |
---|
508 | vector<pair<G4int,G4double> >*h9=new vector<pair<G4int,G4double> >; |
---|
509 | h9->push_back(make_pair(118,1.)); |
---|
510 | natEl.push_back(h9); |
---|
511 | vector<pair<G4int,G4double> >*i0=new vector<pair<G4int,G4double> >; |
---|
512 | i0->push_back(make_pair(122,.2986)); |
---|
513 | i0->push_back(make_pair(120,.5296)); |
---|
514 | i0->push_back(make_pair(119,.6983)); |
---|
515 | i0->push_back(make_pair(121,.8301)); |
---|
516 | i0->push_back(make_pair(118,.9298)); |
---|
517 | i0->push_back(make_pair(124,.9985)); |
---|
518 | i0->push_back(make_pair(116,1.)); |
---|
519 | natEl.push_back(i0); |
---|
520 | vector<pair<G4int,G4double> >*i1=new vector<pair<G4int,G4double> >; |
---|
521 | i1->push_back(make_pair(124,.70476)); |
---|
522 | i1->push_back(make_pair(122,1.)); |
---|
523 | natEl.push_back(i1); |
---|
524 | vector<pair<G4int,G4double> >*i2=new vector<pair<G4int,G4double> >; |
---|
525 | i2->push_back(make_pair(126,.524)); |
---|
526 | i2->push_back(make_pair(124,.765)); |
---|
527 | i2->push_back(make_pair(125,.986)); |
---|
528 | i2->push_back(make_pair(122,1.)); |
---|
529 | natEl.push_back(i2); |
---|
530 | vector<pair<G4int,G4double> >*i3=new vector<pair<G4int,G4double> >; |
---|
531 | i3->push_back(make_pair(126,1.)); |
---|
532 | natEl.push_back(i3); |
---|
533 | vector<pair<G4int,G4double> >*i4=new vector<pair<G4int,G4double> >; |
---|
534 | i4->push_back(make_pair(125,1.)); |
---|
535 | natEl.push_back(i4); |
---|
536 | vector<pair<G4int,G4double> >*i5=new vector<pair<G4int,G4double> >; |
---|
537 | i5->push_back(make_pair(136,1.)); |
---|
538 | natEl.push_back(i5); |
---|
539 | vector<pair<G4int,G4double> >*i6=new vector<pair<G4int,G4double> >; |
---|
540 | i6->push_back(make_pair(136,1.)); |
---|
541 | natEl.push_back(i6); |
---|
542 | vector<pair<G4int,G4double> >*i7=new vector<pair<G4int,G4double> >; |
---|
543 | i7->push_back(make_pair(138,1.)); |
---|
544 | natEl.push_back(i7); |
---|
545 | vector<pair<G4int,G4double> >*i8=new vector<pair<G4int,G4double> >; |
---|
546 | i8->push_back(make_pair(138,1.)); |
---|
547 | natEl.push_back(i8); |
---|
548 | vector<pair<G4int,G4double> >*i9=new vector<pair<G4int,G4double> >; |
---|
549 | i9->push_back(make_pair(142,1.)); |
---|
550 | natEl.push_back(i9); |
---|
551 | vector<pair<G4int,G4double> >*j0=new vector<pair<G4int,G4double> >; |
---|
552 | j0->push_back(make_pair(142,1.)); |
---|
553 | natEl.push_back(j0); |
---|
554 | vector<pair<G4int,G4double> >*j1=new vector<pair<G4int,G4double> >; |
---|
555 | j1->push_back(make_pair(140,1.)); |
---|
556 | natEl.push_back(j1); |
---|
557 | vector<pair<G4int,G4double> >*j2=new vector<pair<G4int,G4double> >; |
---|
558 | j2->push_back(make_pair(146,.992745)); |
---|
559 | j2->push_back(make_pair(143,.999945)); |
---|
560 | j2->push_back(make_pair(142,1.)); |
---|
561 | natEl.push_back(j2); |
---|
562 | vector<pair<G4int,G4double> >*j3=new vector<pair<G4int,G4double> >; |
---|
563 | j3->push_back(make_pair(144,1.)); |
---|
564 | natEl.push_back(j3); |
---|
565 | vector<pair<G4int,G4double> >*j4=new vector<pair<G4int,G4double> >; |
---|
566 | j4->push_back(make_pair(150,1.)); |
---|
567 | natEl.push_back(j4); |
---|
568 | vector<pair<G4int,G4double> >*j5=new vector<pair<G4int,G4double> >; |
---|
569 | j5->push_back(make_pair(148,1.)); |
---|
570 | natEl.push_back(j5); |
---|
571 | vector<pair<G4int,G4double> >*j6=new vector<pair<G4int,G4double> >; |
---|
572 | j6->push_back(make_pair(151,1.)); |
---|
573 | natEl.push_back(j6); |
---|
574 | vector<pair<G4int,G4double> >*j7=new vector<pair<G4int,G4double> >; |
---|
575 | j7->push_back(make_pair(150,1.)); |
---|
576 | natEl.push_back(j7); |
---|
577 | vector<pair<G4int,G4double> >*j8=new vector<pair<G4int,G4double> >; |
---|
578 | j8->push_back(make_pair(153,1.)); |
---|
579 | natEl.push_back(j8); |
---|
580 | vector<pair<G4int,G4double> >*j9=new vector<pair<G4int,G4double> >; |
---|
581 | j9->push_back(make_pair(157,1.)); |
---|
582 | natEl.push_back(j9); |
---|
583 | vector<pair<G4int,G4double> >*k0=new vector<pair<G4int,G4double> >; |
---|
584 | k0->push_back(make_pair(157,1.)); |
---|
585 | natEl.push_back(k0); |
---|
586 | vector<pair<G4int,G4double> >*k1=new vector<pair<G4int,G4double> >; |
---|
587 | k1->push_back(make_pair(157,1.)); |
---|
588 | natEl.push_back(k1); |
---|
589 | vector<pair<G4int,G4double> >*k2=new vector<pair<G4int,G4double> >; |
---|
590 | k2->push_back(make_pair(157,1.)); |
---|
591 | natEl.push_back(k2); |
---|
592 | vector<pair<G4int,G4double> >*k3=new vector<pair<G4int,G4double> >; |
---|
593 | k3->push_back(make_pair(157,1.)); |
---|
594 | natEl.push_back(k3); |
---|
595 | vector<pair<G4int,G4double> >*k4=new vector<pair<G4int,G4double> >; |
---|
596 | k4->push_back(make_pair(157,1.)); |
---|
597 | natEl.push_back(k4); |
---|
598 | vector<pair<G4int,G4double> >*k5=new vector<pair<G4int,G4double> >; |
---|
599 | k5->push_back(make_pair(157,1.)); |
---|
600 | natEl.push_back(k5); |
---|
601 | vector<pair<G4int,G4double> >*k6=new vector<pair<G4int,G4double> >; |
---|
602 | k6->push_back(make_pair(157,1.)); |
---|
603 | natEl.push_back(k6); |
---|
604 | vector<pair<G4int,G4double> >*k7=new vector<pair<G4int,G4double> >; |
---|
605 | k7->push_back(make_pair(155,1.)); |
---|
606 | natEl.push_back(k7); |
---|
607 | vector<pair<G4int,G4double> >*k8=new vector<pair<G4int,G4double> >; |
---|
608 | k8->push_back(make_pair(157,1.)); |
---|
609 | natEl.push_back(k8); |
---|
610 | vector<pair<G4int,G4double> >*k9=new vector<pair<G4int,G4double> >; |
---|
611 | k9->push_back(make_pair(157,1.)); |
---|
612 | natEl.push_back(k9); |
---|
613 | // Now fill natElements and natIsoCrossS |
---|
614 | G4int nona=natEl.size(); |
---|
615 | #ifdef cdebug |
---|
616 | G4cout<<"G4QIsotope::Constructor natEl filling is finished nE="<<nona<<G4endl; |
---|
617 | #endif |
---|
618 | for(G4int i=0; i<nona; i++) |
---|
619 | { |
---|
620 | vector<pair<G4int,G4double> >* is=natEl[i]; // Pointer to theElement |
---|
621 | G4int n=is->size(); |
---|
622 | #ifdef cdebug |
---|
623 | G4cout<<"G4QIsotope::Constructor: Element # "<<i<<", nOfIsotopes="<<n<<G4endl; |
---|
624 | #endif |
---|
625 | vector<pair<G4int,G4double>*>*a=new vector<pair<G4int,G4double>*>; |
---|
626 | vector<pair<G4int,G4double>*>*s=new vector<pair<G4int,G4double>*>; |
---|
627 | G4double last=0.; |
---|
628 | if(n) for(G4int j=0; j<n; j++) |
---|
629 | { |
---|
630 | G4int nn =(*is)[j].first; // #ofNeutrons in the isotope |
---|
631 | G4double cur=(*is)[j].second;// value of the summed abundancy |
---|
632 | pair<G4int,G4double>* aP = new pair<G4int,G4double>(nn,cur-last); |
---|
633 | last=cur; // Update the summed value |
---|
634 | pair<G4int,G4double>* sP = new pair<G4int,G4double>((*is)[j]); |
---|
635 | a->push_back(aP); |
---|
636 | s->push_back(sP); |
---|
637 | #ifdef cdebug |
---|
638 | G4cout<<"G4QIsotope::Constructor:Element# "<<i<<", Pair # "<<j<<" is filled"<<G4endl; |
---|
639 | #endif |
---|
640 | } |
---|
641 | natElements.push_back(a); // Fill abundancies for the particular isotope |
---|
642 | natSumAbund.push_back(s); // Fill summes abundancies up to this isotope |
---|
643 | #ifdef cdebug |
---|
644 | G4cout<<"G4QIsotope::Constructor: natElements is filled"<<G4endl; |
---|
645 | #endif |
---|
646 | vector<pair<G4int,G4double>*>*c=new vector<pair<G4int,G4double>*>; |
---|
647 | if(n) for(G4int j=0; j<n; j++) // Cross sections are 0. by default |
---|
648 | { |
---|
649 | pair<G4int,G4double>* cP = new pair<G4int,G4double>((*is)[j].first,0.); |
---|
650 | c->push_back(cP); |
---|
651 | #ifdef cdebug |
---|
652 | G4cout<<"G4QIsotope::Constructor:CrosSecPair i="<<i<<", j="<<j<<" is filled"<<G4endl; |
---|
653 | #endif |
---|
654 | } |
---|
655 | natIsoCrosS.push_back(c); // FillPrototypeCrossSec's (0) for theParticularIsotope |
---|
656 | #ifdef cdebug |
---|
657 | G4cout<<"G4QIsotope::Constructor: natIsoCrosS is filled"<<G4endl; |
---|
658 | #endif |
---|
659 | delete is; |
---|
660 | } |
---|
661 | #ifdef cdebug |
---|
662 | G4cout<<"G4QIsotope::Constructor: is finished"<<G4endl; |
---|
663 | #endif |
---|
664 | } |
---|
665 | |
---|
666 | G4QIsotope::~G4QIsotope() // The QIsotopes are destructed only in theEnd of theJob |
---|
667 | { |
---|
668 | #ifdef debug |
---|
669 | G4cout<<"G4QIsotope::Destructor is called"<<G4endl; |
---|
670 | #endif |
---|
671 | G4int uP=natElements.size(); |
---|
672 | if(uP) for(G4int i=0; i<uP; i++) |
---|
673 | { |
---|
674 | vector<pair<G4int,G4double>*>* curA=natElements[i]; |
---|
675 | G4int nn=curA->size(); // Can not be 0 by definition |
---|
676 | if(nn) for(G4int n=0; n<nn; n++) delete (*curA)[n]; // Delete pair(N,Ab) |
---|
677 | delete curA; // Delet abundancy vector |
---|
678 | vector<pair<G4int,G4double>*>* curS=natSumAbund[i]; |
---|
679 | G4int ns=curS->size(); // Can not be 0 by definition |
---|
680 | if(ns) for(G4int m=0; m<ns; m++) delete (*curS)[m]; // Delete pair(N,Ab) |
---|
681 | delete curS; // Delet abundancy vector |
---|
682 | vector<pair<G4int,G4double>*>* curC=natIsoCrosS[i]; |
---|
683 | G4int nc=curC->size(); // Can not be 0 by definition |
---|
684 | if(nc) for(G4int k=0; k<nc; k++) delete (*curC)[k]; // Delete pair(N,CS) |
---|
685 | delete curC; // Delete cross section vector |
---|
686 | } |
---|
687 | G4int nP=newElems.size(); |
---|
688 | if(nP) for(G4int j=0; j<nP; j++) // LOOP over new UserDefinedElements |
---|
689 | { |
---|
690 | pair<G4int, vector<pair<G4int,G4double>*>* >* nEl= newElems[j]; |
---|
691 | G4int nEn=nEl->second->size(); |
---|
692 | if(nEn) for(G4int k=0; k<nEn; k++) |
---|
693 | { |
---|
694 | pair<G4int,G4double>* curA=(*(nEl->second))[k]; |
---|
695 | delete curA; // Delete vect<pair(N,Abundancy)*> |
---|
696 | } |
---|
697 | delete nEl; // Delete vect<IndZ,vect<pair(N,Ab)*>*> newElementVector |
---|
698 | // |
---|
699 | pair<G4int, vector<pair<G4int,G4double>*>* >* nSA= newSumAb[j]; |
---|
700 | G4int nSn=nSA->second->size(); |
---|
701 | if(nSn) for(G4int n=0; n<nSn; n++) |
---|
702 | { |
---|
703 | pair<G4int,G4double>* curS=(*(nSA->second))[n]; |
---|
704 | delete curS; // Delete vect<pair(N,SumAbund)*> |
---|
705 | } |
---|
706 | delete nSA; // Delete vect<IndZ,vect<pair(N,SA)*>*> newSumAbunVector |
---|
707 | // |
---|
708 | pair<G4int, vector<pair<G4int,G4double>*>* >* nCS= newIsoCS[j]; |
---|
709 | G4int nCn=nCS->second->size(); |
---|
710 | if(nCn) for(G4int m=0; m<nCn; m++) |
---|
711 | { |
---|
712 | pair<G4int,G4double>* curC = (*(nCS->second))[m]; |
---|
713 | delete curC; // Delete vect<pair(N,CrossSect)*> |
---|
714 | } |
---|
715 | delete nCS; // Delete vect<IndZ,vect<pair(N,CS)*>*> newIsoCroSVector |
---|
716 | // |
---|
717 | if(nEn!=nCn) G4cerr<<"*G4QIsotope-WORNING-:#El="<<j<<":nE="<<nEn<<"!=nC="<<nCn<<G4endl; |
---|
718 | if(nEn!=nSn) G4cerr<<"*G4QIsotope-WORNING-:#El="<<j<<":nE="<<nEn<<"!=nS="<<nSn<<G4endl; |
---|
719 | } |
---|
720 | } |
---|
721 | |
---|
722 | // Returns Pointer to the G4QIsotope singletone |
---|
723 | G4QIsotope* G4QIsotope::Get() |
---|
724 | // ================= |
---|
725 | { |
---|
726 | #ifdef pdebug |
---|
727 | G4cout<<"G4QIsotope::Get is called"<<G4endl; |
---|
728 | #endif |
---|
729 | static G4QIsotope theIsotopes; // *** Static body of the G4QIsotope class *** |
---|
730 | return &theIsotopes; |
---|
731 | } |
---|
732 | |
---|
733 | // #ofProtons in stable isotopes with fixed A=Z+N. Returns length and fils VectOfIsotopes |
---|
734 | G4int G4QIsotope::GetProtons(G4int A, vector<G4int>& isoV) |
---|
735 | // ========================================================= |
---|
736 | { |
---|
737 | const G4int nAZ=270; // Dimension of the table |
---|
738 | // Best Z for the given A |
---|
739 | const G4int bestZ[nAZ] = { |
---|
740 | 0, 1, 1, 2, 2, 0, 3, 3, 4, 4, //0 |
---|
741 | 5, 5, 6, 6, 7, 7, 8, 8, 8, 9, //10 |
---|
742 | 10, 10, 10, 11, 12, 12, 12, 13, 14, 14, //20 |
---|
743 | 14, 15, 16, 16, 16, 17, 18, 17, 18, 19, //30 |
---|
744 | 18, 19, 20, 20, 20, 21, 22, 22, 23, 23, //40 |
---|
745 | 22, 23, 24, 24, 26, 25, 26, 26, 28, 27, //50 |
---|
746 | 28, 28, 28, 29, 30, 29, 30, 30, 30, 31, //60 |
---|
747 | 32, 31, 32, 32, 32, 33, 34, 34, 34, 35, //70 |
---|
748 | 34, 35, 36, 36, 36, 37, 39, 36, 38, 39, //80 |
---|
749 | 40, 40, 41, 40, 40, 42, 42, 42, 42, 44, //90 |
---|
750 | 44, 44, 44, 45, 44, 46, 46, 47, 46, 47, //100 |
---|
751 | 48, 48, 48, 48, 48, 49, 50, 50, 50, 50, //110 |
---|
752 | 50, 51, 50, 51, 50, 52, 52, 53, 52, 54, //120 |
---|
753 | 52, 54, 54, 55, 54, 56, 54, 56, 56, 57, //130 |
---|
754 | 58, 59, 60, 60, 60, 60, 60, 62, 62, 62, //140 |
---|
755 | 62, 63, 62, 63, 62, 64, 64, 64, 64, 65, //150 |
---|
756 | 64, 66, 66, 66, 66, 67, 68, 68, 68, 69, //160 |
---|
757 | 68, 70, 70, 70, 70, 71, 70, 72, 72, 72, //170 |
---|
758 | 72, 73, 74, 74, 74, 75, 74, 75, 76, 76, //180 |
---|
759 | 76, 77, 76, 77, 78, 78, 78, 79, 80, 80, //190 |
---|
760 | 80, 80, 80, 81, 80, 81, 82, 82, 82, 83, //200 |
---|
761 | 82, 0, 82, 0, 82, 0, 84, 0, 0, 0, //210 |
---|
762 | 86, 0, 86, 87, 88, 0, 88, 89, 88, 89, //220 |
---|
763 | 89, 91, 90, 0, 92, 92, 0, 93, 92, 94, //230 |
---|
764 | 0, 0, 0, 95, 94, 0, 0, 96, 0, 0, //240 |
---|
765 | 0, 98, 99, 0, 0, 0, 0,100,101,102, //250 |
---|
766 | 103,104,105,106, 0,108,109, 0, 0, 0}; //260 |
---|
767 | // 0 1 2 3 4 5 6 7 8 9 |
---|
768 | // Second candidate |
---|
769 | const G4int secoZ[nAZ] = { |
---|
770 | 0, 0, 0, 1, 0, 0, 0, 4, 0, 0, //0 |
---|
771 | 4, 6, 5, 7, 0, 8, 0, 0, 0, 8, //10 |
---|
772 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //20 |
---|
773 | 0, 0, 15, 0, 0, 16, 0, 0, 0, 0, //30 |
---|
774 | 20, 20, 0, 0, 0, 0, 20, 0, 0, 0, //40 |
---|
775 | 24, 0, 0, 0, 24, 0, 0, 0, 26, 0, //50 |
---|
776 | 27, 0, 0, 0, 28, 0, 0, 0, 0, 0, //60 |
---|
777 | 30, 0, 0, 0, 34, 0, 32, 0, 0, 0, //70 |
---|
778 | 36, 0, 34, 0, 38, 0, 38, 38, 0, 0, //80 |
---|
779 | 0, 0, 42, 0, 42, 0, 44, 0, 44, 0, //90 |
---|
780 | 42, 0, 46, 0, 46, 0, 48, 0, 48, 0, //100 |
---|
781 | 46, 0, 50, 49, 50, 50, 48, 0, 0, 0, //110 |
---|
782 | 52, 0, 52, 0, 52, 0, 54, 0, 54, 0, //120 |
---|
783 | 54, 0, 56, 0, 56, 0, 56, 0, 58, 0, //130 |
---|
784 | 54, 0, 58, 0, 62, 61, 0, 0, 60, 0, //140 |
---|
785 | 60, 0, 64, 0, 64, 0, 66, 0, 66, 0, //150 |
---|
786 | 66, 0, 68, 0, 68, 0, 0, 0, 70, 0, //160 |
---|
787 | 70, 0, 0, 0, 72, 0, 72, 0, 0, 0, //170 |
---|
788 | 74, 0, 0, 0, 76, 0, 76, 76, 0, 0, //180 |
---|
789 | 78, 0, 78, 0, 0, 0, 80, 0, 78, 0, //190 |
---|
790 | 0, 0, 0, 0, 82, 0, 0, 0, 0, 84, //200 |
---|
791 | 84, 0, 83, 0, 83, 0, 0, 0, 0, 0, //210 |
---|
792 | 0, 0, 0, 0, 0, 0, 0, 0, 89, 0, //220 |
---|
793 | 0, 0, 0, 0, 93, 0, 0, 0, 93, 0, //230 |
---|
794 | 0, 0, 0, 0, 0, 0, 0, 97, 0, 0, //240 |
---|
795 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //250 |
---|
796 | 0, 0,107, 0, 0, 0, 0, 0, 0, 0}; //260 |
---|
797 | // 0 1 2 3 4 5 6 7 8 9 |
---|
798 | // Third candidate |
---|
799 | const G4int thrdZ[nAZ] = { |
---|
800 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //0 |
---|
801 | 0, 0, 7, 0, 0, 0, 0, 0, 0, 0, //10 |
---|
802 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //20 |
---|
803 | 0, 0, 0, 0, 0,20, 0, 0, 0, 0, //30 |
---|
804 | 19, 0, 0, 0, 0, 0, 0, 0, 0, 0, //40 |
---|
805 | 23, 0, 0, 0, 0, 0, 0, 0, 0, 0, //50 |
---|
806 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //60 |
---|
807 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //70 |
---|
808 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //80 |
---|
809 | 0, 0,36, 0,38, 0,40, 0, 0, 0, //90 |
---|
810 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //100 |
---|
811 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //110 |
---|
812 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //120 |
---|
813 | 56, 0,50, 0, 0, 0,58, 0,57, 0, //130 |
---|
814 | 0, 0, 0, 0, 0, 0, 0, 0,65, 0, //140 |
---|
815 | 0, 0,66, 0, 0, 0, 0, 0, 0, 0, //150 |
---|
816 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //160 |
---|
817 | 0, 0, 0, 0, 0, 0,71, 0, 0, 0, //170 |
---|
818 | 73, 0, 0, 0, 0, 0, 0, 0, 0, 0, //180 |
---|
819 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //190 |
---|
820 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //200 |
---|
821 | 83, 0,84, 0,84, 0, 0, 0, 0, 0, //210 |
---|
822 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //220 |
---|
823 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //230 |
---|
824 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //240 |
---|
825 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //250 |
---|
826 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; //260 |
---|
827 | //0 1 2 3 4 5 6 7 8 9 |
---|
828 | // Fourth candidate (only two isotopes) |
---|
829 | const G4int quadZ[nAZ] = { |
---|
830 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //0 |
---|
831 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //10 |
---|
832 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //20 |
---|
833 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //30 |
---|
834 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //40 |
---|
835 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //50 |
---|
836 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //60 |
---|
837 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //70 |
---|
838 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //80 |
---|
839 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //90 |
---|
840 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //100 |
---|
841 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //110 |
---|
842 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //120 |
---|
843 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //130 |
---|
844 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //140 |
---|
845 | 0, 0,67, 0, 0, 0, 0, 0, 0, 0, //150 |
---|
846 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //160 |
---|
847 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //170 |
---|
848 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //180 |
---|
849 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //190 |
---|
850 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //200 |
---|
851 | 85, 0, 0, 0, 0, 0, 0, 0, 0, 0, //210 |
---|
852 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //220 |
---|
853 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //230 |
---|
854 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //240 |
---|
855 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, //250 |
---|
856 | 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; //260 |
---|
857 | //0 1 2 3 4 5 6 7 8 9 |
---|
858 | isoV.clear(); // Always cleans up before filling |
---|
859 | if(A>=nAZ) return 0; |
---|
860 | G4int fZ=bestZ[A]; |
---|
861 | if(fZ) |
---|
862 | { |
---|
863 | isoV.push_back(fZ); |
---|
864 | G4int sZ=secoZ[A]; |
---|
865 | if(sZ) |
---|
866 | { |
---|
867 | isoV.push_back(sZ); |
---|
868 | G4int tZ=thrdZ[A]; |
---|
869 | if(tZ) |
---|
870 | { |
---|
871 | isoV.push_back(tZ); |
---|
872 | G4int qZ=quadZ[A]; |
---|
873 | if(qZ) |
---|
874 | { |
---|
875 | isoV.push_back(qZ); |
---|
876 | return 4; |
---|
877 | } |
---|
878 | else return 3; |
---|
879 | } |
---|
880 | else return 2; |
---|
881 | } |
---|
882 | else return 1; |
---|
883 | } |
---|
884 | else return 0; |
---|
885 | } |
---|
886 | |
---|
887 | // Randomize a#ofNeutrons in the Element (independent function @@ can be used for initial) |
---|
888 | G4int G4QIsotope::RandomizeNeutrons(G4int i) |
---|
889 | // ====================================== |
---|
890 | { |
---|
891 | static const G4int nElements = 110; // Max=Meitnerium(Mt)Z=99(starts with Z=0 - neuteron) |
---|
892 | // If an error is found in this simple tree, please, correct the initialization above |
---|
893 | static G4int dN[nElements]= |
---|
894 | { |
---|
895 | // n Be F Al P |
---|
896 | 1, -1, -1, -1, 5, -1, -1, -1, -1, 10, -1, 12, -1, 14, -1, 16, -1, -1, -1, -1, |
---|
897 | // Sc Mn Co As Y |
---|
898 | -1, 24, -1, -1, -1, 30, -1, 32, -1, -1, -1, -1, -1, 42, -1, -1, -1, -1, -1, 50, |
---|
899 | // Nb Tc Rh I Cs Pr |
---|
900 | -1, 52, -1, 55, -1, 58, -1, -1, -1, -1, -1, -1, -1, 74, -1, 78, -1, -1, -1, 82, |
---|
901 | // Pm Tb Ho Tm Au |
---|
902 | -1,-85, -1, -1, -1, 94, -1, 98, -1,100, -1, -1, -1, -1, -1, -1, -1, -1, -1,118, |
---|
903 | // Po At Rn Fr Ra Ac |
---|
904 | -1, -1, -1, 126,-125,-136,-136,-138,-138,-142, |
---|
905 | // Th Pa U Np Pu Am Cm Bk Cf Es |
---|
906 | 142,-140, -1,-144,-150,-148,-151,-150,-153,-153 |
---|
907 | // Fm Md No Lr Rf Db Sg Bh Hs Mt |
---|
908 | -157,-157,-157,-157,-157,-157,-157,-155,-157,-157 |
---|
909 | }; |
---|
910 | #ifdef debug |
---|
911 | G4cout<<"G4QIsotope::RandomizeNeutrons is called Z="<<i<<G4endl; |
---|
912 | #endif |
---|
913 | G4int N=dN[i]; |
---|
914 | if(N==-1) |
---|
915 | { |
---|
916 | G4double rnd=G4UniformRand(); |
---|
917 | if (i<44) // ====== H - Mo |
---|
918 | { |
---|
919 | if (i<23) // ------ H - Ti |
---|
920 | { |
---|
921 | if (i<12) // ______ H - Ne |
---|
922 | { |
---|
923 | if (i<6) // ...... H - B |
---|
924 | { |
---|
925 | if (i<3) |
---|
926 | { |
---|
927 | if(i==1) // H |
---|
928 | { |
---|
929 | if(rnd>.00015) N=0; |
---|
930 | else N=1; |
---|
931 | } |
---|
932 | else // He (2) |
---|
933 | { |
---|
934 | if(rnd>1.37e-6) N=2; |
---|
935 | else N=1; |
---|
936 | } |
---|
937 | } |
---|
938 | else |
---|
939 | { |
---|
940 | if(i==3) // Li |
---|
941 | { |
---|
942 | if(rnd>.075) N=4; |
---|
943 | else N=3; |
---|
944 | } |
---|
945 | else // B (5) |
---|
946 | { |
---|
947 | if(rnd>.199) N=6; |
---|
948 | else N=5; |
---|
949 | } |
---|
950 | } |
---|
951 | } |
---|
952 | else // ...... C - Ne |
---|
953 | { |
---|
954 | if (i<8) |
---|
955 | { |
---|
956 | if(i==6) // C |
---|
957 | { |
---|
958 | if(rnd>.011) N=6; |
---|
959 | else N=7; |
---|
960 | } |
---|
961 | else // N (7) |
---|
962 | { |
---|
963 | if(rnd>.0037) N=7; |
---|
964 | else N=8; |
---|
965 | } |
---|
966 | } |
---|
967 | else |
---|
968 | { |
---|
969 | if(i==8) // O |
---|
970 | { |
---|
971 | if (rnd<.9976) N=8; |
---|
972 | else if(rnd<.9996) N=10; |
---|
973 | else N=9; |
---|
974 | } |
---|
975 | else // Ne (10) |
---|
976 | { |
---|
977 | if (rnd<.9948) N=10; |
---|
978 | else if(rnd<.9975) N=11; |
---|
979 | else N=12; |
---|
980 | } |
---|
981 | } |
---|
982 | } |
---|
983 | } |
---|
984 | else // ______ Mg - Ti |
---|
985 | { |
---|
986 | if (i<18) // ...... Mg - Cl |
---|
987 | { |
---|
988 | if (i<16) |
---|
989 | { |
---|
990 | if(i==12) // Mg |
---|
991 | { |
---|
992 | if (rnd<.7899) N=12; |
---|
993 | else if(rnd<.8899) N=13; |
---|
994 | else N=14; |
---|
995 | } |
---|
996 | else // Si (14) |
---|
997 | { |
---|
998 | if (rnd<.9223) N=14; |
---|
999 | else if(rnd<.969) N=15; |
---|
1000 | else N=16; |
---|
1001 | } |
---|
1002 | } |
---|
1003 | else |
---|
1004 | { |
---|
1005 | if(i==16) // S |
---|
1006 | { |
---|
1007 | if (rnd<.9502) N=16; |
---|
1008 | else if(rnd<.9923) N=18; |
---|
1009 | else if(rnd<.9998) N=17; |
---|
1010 | else N=20; |
---|
1011 | } |
---|
1012 | else // Cl (17) |
---|
1013 | { |
---|
1014 | if (rnd>.7577) N=18; |
---|
1015 | else N=20; |
---|
1016 | } |
---|
1017 | } |
---|
1018 | } |
---|
1019 | else // ...... Ar - Ti |
---|
1020 | { |
---|
1021 | if (i<20) |
---|
1022 | { |
---|
1023 | if(i==18) // Ar |
---|
1024 | { |
---|
1025 | if (rnd<.996) N=22; |
---|
1026 | else if(rnd<.99937) N=18; |
---|
1027 | else N=20; |
---|
1028 | } |
---|
1029 | else // K (19) |
---|
1030 | { |
---|
1031 | if (rnd<.932581) N=20; |
---|
1032 | else if(rnd<.999883) N=22; |
---|
1033 | else N=21; |
---|
1034 | } |
---|
1035 | } |
---|
1036 | else |
---|
1037 | { |
---|
1038 | if(i==20) // Ca |
---|
1039 | { |
---|
1040 | if (rnd<.96941) N=20; |
---|
1041 | else if(rnd<.99027) N=24; |
---|
1042 | else if(rnd<.99674) N=22; |
---|
1043 | else if(rnd<.99861) N=28; |
---|
1044 | else if(rnd<.99996) N=23; |
---|
1045 | else N=26; |
---|
1046 | } |
---|
1047 | else // Ti (22) |
---|
1048 | { |
---|
1049 | if (rnd<.738) N=26; |
---|
1050 | else if(rnd<.818) N=24; |
---|
1051 | else if(rnd<.891) N=25; |
---|
1052 | else if(rnd<.946) N=27; |
---|
1053 | else N=28; |
---|
1054 | } |
---|
1055 | } |
---|
1056 | } |
---|
1057 | } |
---|
1058 | } |
---|
1059 | else // ------ V - Mo |
---|
1060 | { |
---|
1061 | if (i<32) // ______ V - Ga |
---|
1062 | { |
---|
1063 | if (i<28) // ...... H - Fe |
---|
1064 | { |
---|
1065 | if (i<26) |
---|
1066 | { |
---|
1067 | if(i==23) // V |
---|
1068 | { |
---|
1069 | if (rnd<.9975) N=28; |
---|
1070 | else N=27; |
---|
1071 | } |
---|
1072 | else // Cr (24) |
---|
1073 | { |
---|
1074 | if (rnd<.8379) N=28; |
---|
1075 | else if(rnd<.9329) N=29; |
---|
1076 | else if(rnd<.97635) N=26; |
---|
1077 | else N=30; |
---|
1078 | } |
---|
1079 | } |
---|
1080 | else // Fe (26) |
---|
1081 | { |
---|
1082 | if (rnd<.9172) N=30; |
---|
1083 | else if(rnd<.9762) N=28; |
---|
1084 | else if(rnd<.9972) N=31; |
---|
1085 | else N=32; |
---|
1086 | } |
---|
1087 | } |
---|
1088 | else // ...... Ni - Ga |
---|
1089 | { |
---|
1090 | if (i<30) |
---|
1091 | { |
---|
1092 | if(i==28) // Ni |
---|
1093 | { |
---|
1094 | if (rnd<.68077) N=30; |
---|
1095 | else if(rnd<.943) N=32; |
---|
1096 | else if(rnd<.97934) N=34; |
---|
1097 | else if(rnd<.99074) N=33; |
---|
1098 | else N=36; |
---|
1099 | } |
---|
1100 | else // Cu (29) |
---|
1101 | { |
---|
1102 | if (rnd<.6917) N=34; |
---|
1103 | else N=36; |
---|
1104 | } |
---|
1105 | } |
---|
1106 | else |
---|
1107 | { |
---|
1108 | if(i==30) // Zn |
---|
1109 | { |
---|
1110 | if (rnd<.486) N=34; |
---|
1111 | else if(rnd<.765) N=36; |
---|
1112 | else if(rnd<.953) N=38; |
---|
1113 | else if(rnd<.994) N=37; |
---|
1114 | else N=40; |
---|
1115 | } |
---|
1116 | else // Ga (31) |
---|
1117 | { |
---|
1118 | if (rnd<.60108) N=38; |
---|
1119 | else N=40; |
---|
1120 | } |
---|
1121 | } |
---|
1122 | } |
---|
1123 | } |
---|
1124 | else // ______ Ge - Mo |
---|
1125 | { |
---|
1126 | if (i<37) // ...... H - B |
---|
1127 | { |
---|
1128 | if (i<35) |
---|
1129 | { |
---|
1130 | if(i==32) // Ge |
---|
1131 | { |
---|
1132 | if (rnd<.3594) N=42; |
---|
1133 | else if(rnd<.6360) N=40; |
---|
1134 | else if(rnd<.8484) N=38; |
---|
1135 | else if(rnd<.9256) N=41; |
---|
1136 | else N=44; |
---|
1137 | } |
---|
1138 | else // Se (34) |
---|
1139 | { |
---|
1140 | if (rnd>.4961) N=46; |
---|
1141 | else if(rnd<.7378) N=44; |
---|
1142 | else if(rnd<.8274) N=42; |
---|
1143 | else if(rnd<.9148) N=48; |
---|
1144 | else if(rnd<.9911) N=43; |
---|
1145 | else N=40; |
---|
1146 | } |
---|
1147 | } |
---|
1148 | else |
---|
1149 | { |
---|
1150 | if(i==35) // Br |
---|
1151 | { |
---|
1152 | if (rnd<.5069) N=44; |
---|
1153 | else N=46; |
---|
1154 | } |
---|
1155 | else // Kr (36) |
---|
1156 | { |
---|
1157 | if (rnd<.57) N=48; |
---|
1158 | else if(rnd<.743) N=50; |
---|
1159 | else if(rnd<.859) N=46; |
---|
1160 | else if(rnd<.974) N=47; |
---|
1161 | else if(rnd<.9965) N=44; |
---|
1162 | else N=42; |
---|
1163 | } |
---|
1164 | } |
---|
1165 | } |
---|
1166 | else // ...... Rb - Mo |
---|
1167 | { |
---|
1168 | if (i<40) |
---|
1169 | { |
---|
1170 | if(i==37) // Rb |
---|
1171 | { |
---|
1172 | if (rnd<.7217) N=48; |
---|
1173 | else N=50; |
---|
1174 | } |
---|
1175 | else // SR (38) |
---|
1176 | { |
---|
1177 | if (rnd<.8258) N=50; |
---|
1178 | else if(rnd<.9244) N=48; |
---|
1179 | else if(rnd<.9944) N=49; |
---|
1180 | else N=46; |
---|
1181 | } |
---|
1182 | } |
---|
1183 | else |
---|
1184 | { |
---|
1185 | if(i==40) // Zr |
---|
1186 | { |
---|
1187 | if (rnd<.5145) N=50; |
---|
1188 | else if(rnd<.6883) N=54; |
---|
1189 | else if(rnd<.8598) N=53; |
---|
1190 | else if(rnd<.972) N=51; |
---|
1191 | else N=56; |
---|
1192 | } |
---|
1193 | else // Mo (42) |
---|
1194 | { |
---|
1195 | if (rnd<.2413) N=56; |
---|
1196 | else if(rnd<.4081) N=54; |
---|
1197 | else if(rnd<.5673) N=53; |
---|
1198 | else if(rnd<.7157) N=50; |
---|
1199 | else if(rnd<.8120) N=58; |
---|
1200 | else if(rnd<.9075) N=55; |
---|
1201 | else N=52; |
---|
1202 | } |
---|
1203 | } |
---|
1204 | } |
---|
1205 | } |
---|
1206 | } |
---|
1207 | } |
---|
1208 | else // ====== Ru - U |
---|
1209 | { |
---|
1210 | if (i<66) // ------ Ru - Gd |
---|
1211 | { |
---|
1212 | if (i<54) // ______ Ru - Te |
---|
1213 | { |
---|
1214 | if (i<49) // ...... Ru - Cd |
---|
1215 | { |
---|
1216 | if (i<47) |
---|
1217 | { |
---|
1218 | if(i==44) // Ru |
---|
1219 | { |
---|
1220 | if (rnd<.316) N=58; |
---|
1221 | else if(rnd<.502) N=60; |
---|
1222 | else if(rnd<.673) N=57; |
---|
1223 | else if(rnd<.8) N=55; |
---|
1224 | else if(rnd<.926) N=56; |
---|
1225 | else if(rnd<.9814) N=52; |
---|
1226 | else N=54; |
---|
1227 | } |
---|
1228 | else // Pd (46) |
---|
1229 | { |
---|
1230 | if (rnd<.2733) N=60; |
---|
1231 | else if(rnd<.5379) N=62; |
---|
1232 | else if(rnd<.7612) N=59; |
---|
1233 | else if(rnd<.8784) N=55; |
---|
1234 | else if(rnd<.9898) N=58; |
---|
1235 | else N=56; |
---|
1236 | } |
---|
1237 | } |
---|
1238 | else |
---|
1239 | { |
---|
1240 | if(i==47) // Ag |
---|
1241 | { |
---|
1242 | if(rnd<.51839) N=60; |
---|
1243 | else N=62; |
---|
1244 | } |
---|
1245 | else // Cd (48) |
---|
1246 | { |
---|
1247 | if (rnd<.2873) N=66; |
---|
1248 | else if(rnd<.5286) N=64; |
---|
1249 | else if(rnd<.6566) N=59; |
---|
1250 | else if(rnd<.7815) N=62; |
---|
1251 | else if(rnd<.9037) N=65; |
---|
1252 | else if(rnd<.9786) N=68; |
---|
1253 | else if(rnd<.9911) N=58; |
---|
1254 | else N=60; |
---|
1255 | } |
---|
1256 | } |
---|
1257 | } |
---|
1258 | else // ...... In - Te |
---|
1259 | { |
---|
1260 | if (i<51) |
---|
1261 | { |
---|
1262 | if(i==49) // In |
---|
1263 | { |
---|
1264 | if (rnd<.9577) N=66; |
---|
1265 | else N=64; |
---|
1266 | } |
---|
1267 | else // Sn (50) |
---|
1268 | { |
---|
1269 | if (rnd<.3259) N=70; |
---|
1270 | else if(rnd<.5681) N=68; |
---|
1271 | else if(rnd<.7134) N=66; |
---|
1272 | else if(rnd<.7992) N=69; |
---|
1273 | else if(rnd<.8760) N=67; |
---|
1274 | else if(rnd<.9339) N=74; |
---|
1275 | else if(rnd<.9802) N=72; |
---|
1276 | else if(rnd<.9899) N=62; |
---|
1277 | else N=64; |
---|
1278 | //else if(rnd<.9964) N=64; |
---|
1279 | //else N=65; |
---|
1280 | } |
---|
1281 | } |
---|
1282 | else |
---|
1283 | { |
---|
1284 | if(i==51) // Sb |
---|
1285 | { |
---|
1286 | if (rnd<.5736) N=70; |
---|
1287 | else N=72; |
---|
1288 | } |
---|
1289 | else // Te (52) |
---|
1290 | { |
---|
1291 | if (rnd<.3387) N=78; |
---|
1292 | else if(rnd<.6557) N=76; |
---|
1293 | else if(rnd<.8450) N=74; |
---|
1294 | else if(rnd<.9162) N=73; |
---|
1295 | else if(rnd<.9641) N=72; |
---|
1296 | else if(rnd<.9900) N=70; |
---|
1297 | else if(rnd<.99905) N=71; |
---|
1298 | else N=68; |
---|
1299 | } |
---|
1300 | } |
---|
1301 | } |
---|
1302 | } |
---|
1303 | else // ______ Xe - Gd |
---|
1304 | { |
---|
1305 | if (i<60) // ...... Xe - B |
---|
1306 | { |
---|
1307 | if (i<57) |
---|
1308 | { |
---|
1309 | if(i==54) // Xe |
---|
1310 | { |
---|
1311 | if (rnd<.269) N=78; |
---|
1312 | else if(rnd<.533) N=75; |
---|
1313 | else if(rnd<.745) N=77; |
---|
1314 | else if(rnd<.849) N=80; |
---|
1315 | else if(rnd<.938) N=82; |
---|
1316 | else if(rnd<.979) N=76; |
---|
1317 | else if(rnd<.9981) N=74; |
---|
1318 | else if(rnd<.9991) N=70; |
---|
1319 | else N=72; |
---|
1320 | } |
---|
1321 | else // Ba (56) |
---|
1322 | { |
---|
1323 | if (rnd<.717) N=82; |
---|
1324 | else if(rnd<.8293) N=81; |
---|
1325 | else if(rnd<.9078) N=80; |
---|
1326 | else if(rnd<.97373) N=79; |
---|
1327 | else if(rnd<.99793) N=78; |
---|
1328 | else if(rnd<.99899) N=74; |
---|
1329 | else N=76; |
---|
1330 | } |
---|
1331 | } |
---|
1332 | else |
---|
1333 | { |
---|
1334 | if(i==57) // La |
---|
1335 | { |
---|
1336 | if (rnd<.999098)N=82; |
---|
1337 | else N=81; |
---|
1338 | } |
---|
1339 | else // Ce (58) |
---|
1340 | { |
---|
1341 | if (rnd<.8843) N=82; |
---|
1342 | else if(rnd<.9956) N=84; |
---|
1343 | else if(rnd<.9981) N=80; |
---|
1344 | else N=78; |
---|
1345 | } |
---|
1346 | } |
---|
1347 | } |
---|
1348 | else // ...... Nd - Gd |
---|
1349 | { |
---|
1350 | if (i<63) |
---|
1351 | { |
---|
1352 | if(i==60) // Nd |
---|
1353 | { |
---|
1354 | if (rnd<.2713) N=82; |
---|
1355 | else if(rnd<.5093) N=84; |
---|
1356 | else if(rnd<.6812) N=86; |
---|
1357 | else if(rnd<.8030) N=83; |
---|
1358 | else if(rnd<.8860) N=85; |
---|
1359 | else if(rnd<.9436) N=88; |
---|
1360 | else N=90; |
---|
1361 | } |
---|
1362 | else // Sm (62) |
---|
1363 | { |
---|
1364 | if (rnd<.267) N=90; |
---|
1365 | else if(rnd<.494) N=92; |
---|
1366 | else if(rnd<.644) N=85; |
---|
1367 | else if(rnd<.782) N=87; |
---|
1368 | else if(rnd<.895) N=86; |
---|
1369 | else if(rnd<.969) N=88; |
---|
1370 | else N=82; |
---|
1371 | } |
---|
1372 | } |
---|
1373 | else |
---|
1374 | { |
---|
1375 | if(i==63) // Eu |
---|
1376 | { |
---|
1377 | if (rnd<.522) N=90; |
---|
1378 | else N=89; |
---|
1379 | } |
---|
1380 | else // Gd (64) |
---|
1381 | { |
---|
1382 | if (rnd<.2484) N=94; |
---|
1383 | else if(rnd<.4670) N=96; |
---|
1384 | else if(rnd<.6717) N=92; |
---|
1385 | else if(rnd<.8282) N=93; |
---|
1386 | else if(rnd<.9762) N=91; |
---|
1387 | else if(rnd<.9980) N=90; |
---|
1388 | else N=88; |
---|
1389 | } |
---|
1390 | } |
---|
1391 | } |
---|
1392 | } |
---|
1393 | } |
---|
1394 | else // ------ Dy - U |
---|
1395 | { |
---|
1396 | if (i<76) // ______ Dy - Re |
---|
1397 | { |
---|
1398 | if (i<72) // ...... Dy - Lu |
---|
1399 | { |
---|
1400 | if (i<70) |
---|
1401 | { |
---|
1402 | if(i==66) // Dy |
---|
1403 | { |
---|
1404 | if (rnd<.282) N=98; |
---|
1405 | else if(rnd<.537) N=96; |
---|
1406 | else if(rnd<.786) N=97; |
---|
1407 | else if(rnd<.975) N=95; |
---|
1408 | else if(rnd<.9984) N=94; |
---|
1409 | else if(rnd<.9994) N=92; |
---|
1410 | else N=90; |
---|
1411 | } |
---|
1412 | else // Er (68) |
---|
1413 | { |
---|
1414 | if (rnd<.3360) N= 98; |
---|
1415 | else if(rnd<.6040) N=100; |
---|
1416 | else if(rnd<.8335) N= 99; |
---|
1417 | else if(rnd<.9825) N=102; |
---|
1418 | else if(rnd<.9986) N= 96; |
---|
1419 | else N= 94; |
---|
1420 | } |
---|
1421 | } |
---|
1422 | else |
---|
1423 | { |
---|
1424 | if(i==70) // Yb |
---|
1425 | { |
---|
1426 | if (rnd<.3180) N=104; |
---|
1427 | else if(rnd<.5370) N=102; |
---|
1428 | else if(rnd<.6982) N=103; |
---|
1429 | else if(rnd<.8412) N=101; |
---|
1430 | else if(rnd<.9682) N=106; |
---|
1431 | else if(rnd<.9987) N=100; |
---|
1432 | else N= 98; |
---|
1433 | } |
---|
1434 | else // Lu (71) |
---|
1435 | { |
---|
1436 | if (rnd<.9741) N=104; |
---|
1437 | else N=105; |
---|
1438 | } |
---|
1439 | } |
---|
1440 | } |
---|
1441 | else // ...... Hf - Re |
---|
1442 | { |
---|
1443 | if (i<74) |
---|
1444 | { |
---|
1445 | if(i==72) // Hf |
---|
1446 | { |
---|
1447 | if (rnd<.35100) N=108; |
---|
1448 | else if(rnd<.62397) N=106; |
---|
1449 | else if(rnd<.81003) N=105; |
---|
1450 | else if(rnd<.94632) N=107; |
---|
1451 | else if(rnd<.99838) N=104; |
---|
1452 | else N=102; |
---|
1453 | } |
---|
1454 | else // Ta (73) |
---|
1455 | { |
---|
1456 | if(rnd<.99988) N=108; |
---|
1457 | else N=107; |
---|
1458 | } |
---|
1459 | } |
---|
1460 | else |
---|
1461 | { |
---|
1462 | if(i==74) // W |
---|
1463 | { |
---|
1464 | if (rnd<.307) N=110; |
---|
1465 | else if(rnd<.593) N=112; |
---|
1466 | else if(rnd<.856) N=108; |
---|
1467 | else if(rnd<.9988) N=109; |
---|
1468 | else N=106; |
---|
1469 | } |
---|
1470 | else // Re (75) |
---|
1471 | { |
---|
1472 | if (rnd<.626) N=112; |
---|
1473 | else N=110; |
---|
1474 | } |
---|
1475 | } |
---|
1476 | } |
---|
1477 | } |
---|
1478 | else // ______ Os - U |
---|
1479 | { |
---|
1480 | if (i<81) // ...... Os - Hg |
---|
1481 | { |
---|
1482 | if (i<78) |
---|
1483 | { |
---|
1484 | if(i==76) // Os |
---|
1485 | { |
---|
1486 | if (rnd<.410) N=116; |
---|
1487 | else if(rnd<.674) N=114; |
---|
1488 | else if(rnd<.835) N=113; |
---|
1489 | else if(rnd<.968) N=112; |
---|
1490 | else if(rnd<.984) N=111; |
---|
1491 | else if(rnd<.9998) N=110; |
---|
1492 | else N=108; |
---|
1493 | } |
---|
1494 | else // Ir (77) |
---|
1495 | { |
---|
1496 | if (rnd<.627) N=116; |
---|
1497 | else N=114; |
---|
1498 | } |
---|
1499 | } |
---|
1500 | else |
---|
1501 | { |
---|
1502 | if(i==78) // Pt |
---|
1503 | { |
---|
1504 | if (rnd<.338) N=117; |
---|
1505 | else if(rnd<.667) N=116; |
---|
1506 | else if(rnd<.920) N=118; |
---|
1507 | else if(rnd<.992) N=120; |
---|
1508 | else if(rnd<.9999) N=114; |
---|
1509 | else N=112; |
---|
1510 | } |
---|
1511 | else // Hg (80) |
---|
1512 | { |
---|
1513 | if (rnd<.2986) N=122; |
---|
1514 | else if(rnd<.5296) N=120; |
---|
1515 | else if(rnd<.6983) N=119; |
---|
1516 | else if(rnd<.8301) N=121; |
---|
1517 | else if(rnd<.9298) N=118; |
---|
1518 | else if(rnd<.9985) N=124; |
---|
1519 | else N=116; |
---|
1520 | } |
---|
1521 | } |
---|
1522 | } |
---|
1523 | else // ...... Tl - U |
---|
1524 | { |
---|
1525 | if (i<92) |
---|
1526 | { |
---|
1527 | if (i==81) // Tl |
---|
1528 | { |
---|
1529 | if (rnd<.70476) N=124; |
---|
1530 | else N=122; |
---|
1531 | } |
---|
1532 | else // Pb (82) |
---|
1533 | { |
---|
1534 | if (rnd<.524) N=126; |
---|
1535 | else if(rnd<.765) N=124; |
---|
1536 | else if(rnd<.986) N=125; |
---|
1537 | else N=122; |
---|
1538 | } |
---|
1539 | } |
---|
1540 | else // U (92) |
---|
1541 | { |
---|
1542 | if (rnd<.992745)N=146; |
---|
1543 | else if(rnd<.999945)N=143; |
---|
1544 | else N=142; |
---|
1545 | } |
---|
1546 | } |
---|
1547 | } |
---|
1548 | } |
---|
1549 | } |
---|
1550 | } |
---|
1551 | else if(N<0) |
---|
1552 | { |
---|
1553 | N=-N; |
---|
1554 | G4cerr<<"---Worn---G4QIsotope::RandomizeNeutrons:UnstableElement's used Z="<<i<<G4endl; |
---|
1555 | } |
---|
1556 | #ifdef debug |
---|
1557 | G4cout<<"G4QIsotope::RandomizeNeutrons: Z="<<i<<", N="<<N<<G4endl; |
---|
1558 | #endif |
---|
1559 | return N; |
---|
1560 | } |
---|
1561 | |
---|
1562 | // Returns the input index (if it is >0 & unique) or theFirstFreeIndex (<=0 or nonunique) |
---|
1563 | G4int G4QIsotope::InitElement(G4int Z, G4int index, // Ret: -1 - Empty, -2 - Wrong (sum>1) |
---|
1564 | vector<pair<G4int,G4double>*>* abund) |
---|
1565 | // ======================================================================= |
---|
1566 | { |
---|
1567 | G4int I=abund->size(); |
---|
1568 | #ifdef debug |
---|
1569 | G4cout<<"G4QIsotope::InitElement: called with I="<<I<<" pairs,Z="<<Z<<",i="<<indexG4endl; |
---|
1570 | #endif |
---|
1571 | if(I<=0) |
---|
1572 | { |
---|
1573 | G4cerr<<"--Worning--G4QIsotope::InitEl:(-1)0VectorOfNewEl,Z="<<Z<<",i="<<index<<G4endl; |
---|
1574 | return -2; |
---|
1575 | } |
---|
1576 | if(IsDefined(Z,index)) // This index is already defined |
---|
1577 | { |
---|
1578 | G4cerr<<"-Worning-G4QIsotope::InitEl:VONewEl,Z="<<Z<<",ind="<<index<<" exists"<<G4endl; |
---|
1579 | return index; |
---|
1580 | } |
---|
1581 | G4int ZInd=1000*index+Z; // Fake Z increased by the UserDefinedIndex |
---|
1582 | vector<pair<G4int,G4double>*>*A = new vector<pair<G4int,G4double>*>; |
---|
1583 | vector<pair<G4int,G4double>*>*S = new vector<pair<G4int,G4double>*>; |
---|
1584 | vector<pair<G4int,G4double>*>*C = new vector<pair<G4int,G4double>*>; |
---|
1585 | #ifdef debug |
---|
1586 | G4cout<<"G4QIsotope::InitElement: A & S & C vectors are alocated"<<G4endl; |
---|
1587 | #endif |
---|
1588 | G4double sumAbu=0; // Summ of abbundancies |
---|
1589 | for(G4int j=0; j<I; j++) |
---|
1590 | { |
---|
1591 | G4int N=(*abund)[j]->first; |
---|
1592 | G4double abu=(*abund)[j]->second; |
---|
1593 | #ifdef debug |
---|
1594 | G4cout<<"G4QIsotope::InitElement: pair#"<<j<<", N="<<N<<", abund="<<abu<<G4endl; |
---|
1595 | #endif |
---|
1596 | sumAbu+=abu; |
---|
1597 | if(j==I-1.) |
---|
1598 | { |
---|
1599 | if(fabs(sumAbu-1.)>.00001) |
---|
1600 | { |
---|
1601 | G4cerr<<"--Worning--G4QIsotope::InitEl:maxSum="<<sumAbu<<" is fixed to 1."<<G4endl; |
---|
1602 | abu+=1.-sumAbu; |
---|
1603 | sumAbu=1.; |
---|
1604 | } |
---|
1605 | else if(sumAbu-abu>1.) |
---|
1606 | { |
---|
1607 | G4cerr<<"--Worning--G4QIsotope::InitEl:(-2)WrongAbund,Z="<<Z<<",i="<<index<<G4endl; |
---|
1608 | for(G4int k=0; k<I-1; k++) |
---|
1609 | { |
---|
1610 | delete (*A)[k]; |
---|
1611 | delete (*S)[k]; |
---|
1612 | delete (*C)[k]; |
---|
1613 | } |
---|
1614 | delete A; |
---|
1615 | delete S; |
---|
1616 | delete C; |
---|
1617 | return -2; |
---|
1618 | } |
---|
1619 | #ifdef debug |
---|
1620 | G4cout<<"G4QIsotope::InitElement:TheLastIsChecked it isOK or coredTo "<<sumAbu<<G4endl; |
---|
1621 | #endif |
---|
1622 | } |
---|
1623 | pair<G4int,G4double>* abP= new pair<G4int,G4double>(N,abu); |
---|
1624 | A->push_back(abP); |
---|
1625 | pair<G4int,G4double>* saP= new pair<G4int,G4double>(N,sumAbu); |
---|
1626 | S->push_back(saP); |
---|
1627 | pair<G4int,G4double>* csP= new pair<G4int,G4double>(N,0.); |
---|
1628 | C->push_back(csP); |
---|
1629 | #ifdef debug |
---|
1630 | G4cout<<"G4QIsotope::InitElement: A & S & C are filled nP="<<C->size()<<G4endl; |
---|
1631 | #endif |
---|
1632 | } |
---|
1633 | pair<G4int,vector<pair<G4int,G4double>*>*>* newAP= |
---|
1634 | new pair<G4int,vector<pair<G4int,G4double>*>*>(ZInd,A); |
---|
1635 | newElems.push_back(newAP); |
---|
1636 | pair<G4int,vector<pair<G4int,G4double>*>*>* newSA= |
---|
1637 | new pair<G4int,vector<pair<G4int,G4double>*>*>(ZInd,S); |
---|
1638 | newSumAb.push_back(newSA); |
---|
1639 | pair<G4int,vector<pair<G4int,G4double>*>*>* newCP= |
---|
1640 | new pair<G4int,vector<pair<G4int,G4double>*>*>(ZInd,C); |
---|
1641 | newIsoCS.push_back(newCP); |
---|
1642 | #ifdef debug |
---|
1643 | G4cout<<"G4QIsotope::InitElement: newElems & newSumAb & newIsoCS are filled "<<G4endl; |
---|
1644 | #endif |
---|
1645 | return index; |
---|
1646 | } |
---|
1647 | |
---|
1648 | // The highest index defined for Element with Z (Index>0 correspondToUserDefinedElements) |
---|
1649 | G4int G4QIsotope::GetLastIndex(G4int Z) // Returns theLastDefinedIndex (if onlyNatural: =0) |
---|
1650 | // ================================= |
---|
1651 | { |
---|
1652 | #ifdef debug |
---|
1653 | G4cout<<"G4QIsotope::GetLastIndex is called Z="<<Z<<G4endl; |
---|
1654 | #endif |
---|
1655 | G4int mind=0; // Prototype of the maximum existing index for this Z |
---|
1656 | G4int nE=newElems.size(); // A number of definitions in the newElements vector |
---|
1657 | if(nE) for(G4int i=0; i<nE; i++) // LOOP over new UserDefinedElements |
---|
1658 | { |
---|
1659 | G4int zin=newElems[i]->first; |
---|
1660 | G4int zi=zin%1000; // Existing Z |
---|
1661 | G4int in=zin/1000; // Existing index |
---|
1662 | if(Z==zi && in>mind) mind=in; // maximum index for this Z |
---|
1663 | } |
---|
1664 | return mind; |
---|
1665 | } |
---|
1666 | |
---|
1667 | // Indices can have differen numbers (not 1,2,3,...) & in different sequences (9,3,7,...) |
---|
1668 | G4bool G4QIsotope::IsDefined(G4int Z, G4int Ind) // Ind is an index to be found (true) |
---|
1669 | // ========================================== |
---|
1670 | { |
---|
1671 | #ifdef debug |
---|
1672 | G4cout<<"G4QIsotope::IsDefined is called Z="<<Z<<", I="<<Ind<<G4endl; |
---|
1673 | #endif |
---|
1674 | if(Ind<=0) |
---|
1675 | { |
---|
1676 | if(Ind<0) G4cerr<<"-W-G4QIsotope::IsDefined: Z="<<Z<<", Ind="<<Ind<<" < 0->=0"<<G4endl; |
---|
1677 | return true; // to avoid definition with the negative index |
---|
1678 | } |
---|
1679 | G4int nE=newElems.size(); // A number of definitions in the newElements vector |
---|
1680 | if(nE) for(G4int i=0; i<nE; i++) // LOOP over new UserDefinedElements |
---|
1681 | { |
---|
1682 | G4int zin=newElems[i]->first; |
---|
1683 | G4int zi=zin%1000; // Existing Z |
---|
1684 | G4int in=zin/1000; // Existing index |
---|
1685 | if(Z==zi && Ind==in) return true; // The index for the element Z is found |
---|
1686 | } |
---|
1687 | return false; // The index for the element Z is not found |
---|
1688 | } |
---|
1689 | |
---|
1690 | // A#ofNeutrons in theElement with Z & UserDefIndex. Universal for Nat(index=0) & UserDefEl |
---|
1691 | G4int G4QIsotope::GetNeutrons(G4int Z, G4int index) // If theElem doesn't exist, returns <0 |
---|
1692 | // ============================================= |
---|
1693 | { |
---|
1694 | #ifdef debug |
---|
1695 | G4cout<<"G4QIsotope::GetNeutrons is called Z="<<Z<<", index="<<index<<G4endl; |
---|
1696 | #endif |
---|
1697 | // To reduce the code, but make the member function a bit slower, one can use for natural |
---|
1698 | // isotopes the same algorithm as for the newElements, splitting the natElements Vector |
---|
1699 | if(!index) return RandomizeNeutrons(Z); // @@ Fast decision for the natural isotopes |
---|
1700 | else if(index<0) |
---|
1701 | { |
---|
1702 | G4cerr<<"---Worning---G4QIsotope::GetNeutrons:(-2) Negative Index i="<<index<<G4endl; |
---|
1703 | return -2; |
---|
1704 | } |
---|
1705 | // For the positive index tries to randomize the newUserDefinedElement |
---|
1706 | G4bool found=false; // Prototype of the"ZWithTheSameIndex is found" event |
---|
1707 | G4int nE=newElems.size(); // A number of definitions in the newElements Vector |
---|
1708 | G4int i=0; |
---|
1709 | if(nE) for(i=0; i<nE; i++) |
---|
1710 | { |
---|
1711 | G4int zin=newElems[i]->first; |
---|
1712 | G4int in=zin/1000; // Existing index |
---|
1713 | G4int zi=zin%1000; // Existing Z |
---|
1714 | if(Z==zi && in==index) |
---|
1715 | { |
---|
1716 | found=true; // The newElement with the same Z & index is found |
---|
1717 | break; // Finish the search and quit the loop |
---|
1718 | } |
---|
1719 | } |
---|
1720 | if(!found) |
---|
1721 | { |
---|
1722 | G4cerr<<"--Worning--G4QIsotope::GetNeutrons:(-1) NotFound Z="<<Z<<",i="<<index<<G4endl; |
---|
1723 | return -1; |
---|
1724 | } |
---|
1725 | vector<pair<G4int,G4double>*>* abu = newSumAb[i]->second; |
---|
1726 | G4int nn = abu->size(); // A#Of UserDefinedIsotopes for the newElement |
---|
1727 | if(nn>0) |
---|
1728 | { |
---|
1729 | if(nn==1) return (*abu)[0]->first; |
---|
1730 | else |
---|
1731 | { |
---|
1732 | G4double rnd=G4UniformRand(); |
---|
1733 | G4int j=0; |
---|
1734 | for(j=0; j<nn; j++) if ( rnd < (*abu)[j]->second ) break; |
---|
1735 | if(j>=nn) j=nn-1; |
---|
1736 | return (*abu)[j]->first; |
---|
1737 | } |
---|
1738 | } |
---|
1739 | else |
---|
1740 | { |
---|
1741 | G4cerr<<"--Worning--G4QIsotope::GetNeutrons:(-3) Empty Z="<<Z<<",i="<<index<<G4endl; |
---|
1742 | return -3; |
---|
1743 | } |
---|
1744 | } |
---|
1745 | |
---|
1746 | // Get a pointer to the vector of pairs(N,CrosSec), where N is used to calculate CrosSec |
---|
1747 | vector<pair<G4int,G4double>*>* G4QIsotope::GetCSVector(G4int Z, G4int index) |
---|
1748 | // ============================================= |
---|
1749 | { |
---|
1750 | #ifdef debug |
---|
1751 | G4cout<<"G4QIsotope::GetCSVector is called"<<G4endl; |
---|
1752 | #endif |
---|
1753 | if(index<0) |
---|
1754 | { |
---|
1755 | G4cerr<<"---Worning---G4QIsotope::GetSCVector:(-1) Negative Index i="<<index<<G4endl; |
---|
1756 | return 0; |
---|
1757 | } |
---|
1758 | else if(!index) return natIsoCrosS[Z]; |
---|
1759 | // For the positive index tries to find the newUserDefinedElement |
---|
1760 | G4bool found=false; // Prototype of the"ZWithTheSameIndex is found" event |
---|
1761 | G4int nE=newIsoCS.size(); // A number of definitions in the newElements Vector |
---|
1762 | G4int i=0; |
---|
1763 | if(nE) for(i=0; i<nE; i++) |
---|
1764 | { |
---|
1765 | G4int zin=newIsoCS[i]->first; |
---|
1766 | G4int in=zin/1000; // Existing index |
---|
1767 | G4int zi=zin%1000; // Existing Z |
---|
1768 | if(Z==zi && in==index) |
---|
1769 | { |
---|
1770 | found=true; // The newElement with the same Z & index is found |
---|
1771 | break; // Finish the search and quit the loop |
---|
1772 | } |
---|
1773 | } |
---|
1774 | if(!found) |
---|
1775 | { |
---|
1776 | G4cerr<<"--Worning--G4QIsotope::GetSCVector:(-2) NotFound Z="<<Z<<",i="<<index<<G4endl; |
---|
1777 | return 0; |
---|
1778 | } |
---|
1779 | return newIsoCS[i]->second; |
---|
1780 | } |
---|
1781 | |
---|
1782 | // Get a pointer to the vector of pairs(N,IntAbundancy) for the element with Z |
---|
1783 | vector<pair<G4int,G4double>*>* G4QIsotope::GetAbuVector(G4int Z, G4int index) |
---|
1784 | // ============================================== |
---|
1785 | { |
---|
1786 | #ifdef debug |
---|
1787 | G4cout<<"G4QIsotope::GetAbuVector is called"<<G4endl; |
---|
1788 | #endif |
---|
1789 | if(index<0) |
---|
1790 | { |
---|
1791 | G4cerr<<"---Worning---G4QIsotope::GetAbuVector:(-1) Negative Index i="<<index<<G4endl; |
---|
1792 | return 0; |
---|
1793 | } |
---|
1794 | else if(!index) return natElements[Z]; |
---|
1795 | // For the positive index tries to find the newUserDefinedElement |
---|
1796 | G4bool found=false; // Prototype of the"ZWithTheSameIndex is found" event |
---|
1797 | G4int nE=newElems.size(); // A number of definitions in the newElements Vector |
---|
1798 | G4int i=0; |
---|
1799 | if(nE) for(i=0; i<nE; i++) |
---|
1800 | { |
---|
1801 | G4int zin=newElems[i]->first; |
---|
1802 | G4int in=zin/1000; // Existing index |
---|
1803 | G4int zi=zin%1000; // Existing Z |
---|
1804 | if(Z==zi && in==index) |
---|
1805 | { |
---|
1806 | found=true; // The newElement with the same Z & index is found |
---|
1807 | break; // Finish the search and quit the loop |
---|
1808 | } |
---|
1809 | } |
---|
1810 | if(!found) |
---|
1811 | { |
---|
1812 | G4cerr<<"--Worning--G4QIsotope::GetAbuVector:(-2)NotFound Z="<<Z<<",i="<<index<<G4endl; |
---|
1813 | return 0; |
---|
1814 | } |
---|
1815 | return newElems[i]->second; |
---|
1816 | } |
---|
1817 | |
---|
1818 | // Get a pointer to the vector of pairs(N,SumAbundancy) for the element with Z |
---|
1819 | vector<pair<G4int,G4double>*>* G4QIsotope::GetSumAVector(G4int Z, G4int index) |
---|
1820 | // =============================================== |
---|
1821 | { |
---|
1822 | #ifdef debug |
---|
1823 | G4cout<<"G4QIsotope::GetSumAVector is called"<<G4endl; |
---|
1824 | #endif |
---|
1825 | if(index<0) |
---|
1826 | { |
---|
1827 | G4cerr<<"---Worning---G4QIsotope::GetSumAVector:(-1) Negative Index i="<<index<<G4endl; |
---|
1828 | return 0; |
---|
1829 | } |
---|
1830 | else if(!index) return natSumAbund[Z]; |
---|
1831 | // For the positive index tries to find the newUserDefinedElement |
---|
1832 | G4bool found=false; // Prototype of the"ZWithTheSameIndex is found" event |
---|
1833 | G4int nE=newSumAb.size(); // A number of definitions in the newElements Vector |
---|
1834 | G4int i=0; |
---|
1835 | if(nE) for(i=0; i<nE; i++) |
---|
1836 | { |
---|
1837 | G4int zin=newSumAb[i]->first; |
---|
1838 | G4int in=zin/1000; // Existing index |
---|
1839 | G4int zi=zin%1000; // Existing Z |
---|
1840 | if(Z==zi && in==index) |
---|
1841 | { |
---|
1842 | found=true; // The newElement with the same Z & index is found |
---|
1843 | break; // Finish the search and quit the loop |
---|
1844 | } |
---|
1845 | } |
---|
1846 | if(!found) |
---|
1847 | { |
---|
1848 | G4cerr<<"-Worning-G4QIsotope::GetSumAVector:(-2)Not Found Z="<<Z<<",i="<<index<<G4endl; |
---|
1849 | return 0; |
---|
1850 | } |
---|
1851 | return newSumAb[i]->second; |
---|
1852 | } |
---|
1853 | |
---|
1854 | // Calculates the mean Cross Section for the initialized Element(ind=0 Nat,ind>0 UserDef) |
---|
1855 | G4double G4QIsotope::GetMeanCrossSection(G4int Z, G4int index) |
---|
1856 | // ========================================= |
---|
1857 | { |
---|
1858 | vector<pair<G4int,G4double>*>* ab; |
---|
1859 | vector<pair<G4int,G4double>*>* cs; |
---|
1860 | #ifdef ppdebug |
---|
1861 | G4cout<<"G4QIsotope::GetMeanCrossSection is called"<<G4endl; |
---|
1862 | #endif |
---|
1863 | if(index<0) |
---|
1864 | { |
---|
1865 | G4cerr<<"---Worning---G4QIsotope::GetMeanCS:(-1) Negative Index i="<<index<<G4endl; |
---|
1866 | return -1.; |
---|
1867 | } |
---|
1868 | else if(!index) // =========> Natural Abundancies for Isotopes of the Element |
---|
1869 | { |
---|
1870 | #ifdef ppdebug |
---|
1871 | G4cout<<"G4QIsotope::GetMeanCrossSection: Nat Abundance, Z="<<Z<<G4endl; |
---|
1872 | #endif |
---|
1873 | ab=natElements[Z]; |
---|
1874 | cs=natIsoCrosS[Z]; |
---|
1875 | } |
---|
1876 | else // =========> UserDefinedAbundancies for Isotopes of theElement |
---|
1877 | { |
---|
1878 | #ifdef ppdebug |
---|
1879 | G4cout<<"G4QIsotope::GetMeanCrossSection: Art Abund, Z="<<Z<<",ind="<<index<<G4endl; |
---|
1880 | #endif |
---|
1881 | // For the positive index tries to find the newUserDefinedElement |
---|
1882 | G4bool found=false; // Prototype of the"ZWithTheSameIndex is found" event |
---|
1883 | G4int nE=newIsoCS.size(); // A number of definitions in the newElements Vector |
---|
1884 | G4int i=0; |
---|
1885 | if(nE) for(i=0; i<nE; i++) |
---|
1886 | { |
---|
1887 | G4int zin=newIsoCS[i]->first; |
---|
1888 | G4int in=zin/1000; // Existing index |
---|
1889 | G4int zi=zin%1000; // Existing Z |
---|
1890 | if(Z==zi && in==index) |
---|
1891 | { |
---|
1892 | found=true; // The newElement with the same Z & index is found |
---|
1893 | break; // Finish the search and quit the loop |
---|
1894 | } |
---|
1895 | } |
---|
1896 | if(!found) |
---|
1897 | { |
---|
1898 | G4cerr<<"--Worning--G4QIsotope::GetMeanCS:(-2) NotFound Z="<<Z<<",i="<<index<<G4endl; |
---|
1899 | return -2.; |
---|
1900 | } |
---|
1901 | ab=newElems[i]->second; |
---|
1902 | cs=newIsoCS[i]->second; |
---|
1903 | } |
---|
1904 | G4int nis=ab->size(); |
---|
1905 | //G4double last=0.; |
---|
1906 | if(!nis) |
---|
1907 | { |
---|
1908 | G4cerr<<"--Worning--G4QIsotope::GetMeanCS:(-3) Empty Z="<<Z<<",i="<<index<<G4endl; |
---|
1909 | return -3.; |
---|
1910 | } |
---|
1911 | else |
---|
1912 | { |
---|
1913 | G4double sum=0.; |
---|
1914 | for(G4int j=0; j<nis; j++) |
---|
1915 | { |
---|
1916 | G4double cur=(*ab)[j]->second; |
---|
1917 | //G4double abunda=cur-last; |
---|
1918 | //last=cur; |
---|
1919 | #ifdef ppdebug |
---|
1920 | G4cout<<"G4QIsot::GetMeanCS:j="<<j<<",ab="<<cur<<",CS="<<(*cs)[j]->second<<G4endl; |
---|
1921 | #endif |
---|
1922 | //sum+=abunda * (*cs)[j]->second; |
---|
1923 | sum+=cur * (*cs)[j]->second; |
---|
1924 | } |
---|
1925 | return sum; |
---|
1926 | } |
---|
1927 | } |
---|
1928 | |
---|
1929 | // Randomize A#OfNeutrons in the Isotope weighted by theAbubdancies and theCrossSections |
---|
1930 | G4int G4QIsotope::GetCSNeutrons(G4int Z, G4int index) |
---|
1931 | // =================================== |
---|
1932 | { |
---|
1933 | vector<pair<G4int,G4double>*>* ab; |
---|
1934 | vector<pair<G4int,G4double>*>* cs; |
---|
1935 | #ifdef debug |
---|
1936 | G4cout<<"G4QIsotope::GetCSNeutrons is called"<<G4endl; |
---|
1937 | #endif |
---|
1938 | if(index<0) |
---|
1939 | { |
---|
1940 | G4cerr<<"---Worning---G4QIsotope::GetCSNeutrons:(-1) Negative Index i="<<index<<G4endl; |
---|
1941 | return -1; |
---|
1942 | } |
---|
1943 | else if(!index) // =========> Natural Abundancies for Isotopes of the Element |
---|
1944 | { |
---|
1945 | ab=natElements[Z]; |
---|
1946 | cs=natIsoCrosS[Z]; |
---|
1947 | } |
---|
1948 | else // =========> UserDefinedAbundancies for Isotopes of theElement |
---|
1949 | { |
---|
1950 | // For the positive index tries to find the newUserDefinedElement |
---|
1951 | G4bool found=false; // Prototype of the"ZWithTheSameIndex is found" event |
---|
1952 | G4int nE=newIsoCS.size(); // A number of definitions in the newElements Vector |
---|
1953 | G4int i=0; |
---|
1954 | if(nE) for(i=0; i<nE; i++) |
---|
1955 | { |
---|
1956 | G4int zin=newIsoCS[i]->first; |
---|
1957 | G4int in=zin/1000; // Existing index |
---|
1958 | G4int zi=zin%1000; // Existing Z |
---|
1959 | if(Z==zi && in==index) |
---|
1960 | { |
---|
1961 | found=true; // The newElement with the same Z & index is found |
---|
1962 | break; // Finish the search and quit the loop |
---|
1963 | } |
---|
1964 | } |
---|
1965 | if(!found) |
---|
1966 | { |
---|
1967 | G4cerr<<"--Worning--G4QIsotope::GetCSNeut:(-2) NotFound Z="<<Z<<",i="<<index<<G4endl; |
---|
1968 | return -2; |
---|
1969 | } |
---|
1970 | ab=newElems[i]->second; |
---|
1971 | cs=newIsoCS[i]->second; |
---|
1972 | } |
---|
1973 | G4int nis=ab->size(); |
---|
1974 | G4double last=0.; |
---|
1975 | if(!nis) |
---|
1976 | { |
---|
1977 | G4cerr<<"--Worning--G4QIsotope::GetCSNeutrons:(-3) Empty Z="<<Z<<",i="<<index<<G4endl; |
---|
1978 | return -3; |
---|
1979 | } |
---|
1980 | else |
---|
1981 | { |
---|
1982 | G4double sum=0.; |
---|
1983 | vector<G4double> scs(nis); |
---|
1984 | for(G4int j=0; j<nis; j++) |
---|
1985 | { |
---|
1986 | G4double cur=(*ab)[j]->second; |
---|
1987 | G4double abunda=cur-last; |
---|
1988 | last=cur; |
---|
1989 | sum+=abunda * (*cs)[j]->second;; |
---|
1990 | scs.push_back(sum); |
---|
1991 | } |
---|
1992 | G4double rnd=sum*G4UniformRand(); |
---|
1993 | sum=0; |
---|
1994 | G4int k=0; |
---|
1995 | if(nis>1) for(k=0; k<nis; k++) if(rnd<scs[k]) break; |
---|
1996 | return (*ab)[k]->first; |
---|
1997 | } |
---|
1998 | } |
---|