1 | // |
---|
2 | // ******************************************************************** |
---|
3 | // * License and Disclaimer * |
---|
4 | // * * |
---|
5 | // * The Geant4 software is copyright of the Copyright Holders of * |
---|
6 | // * the Geant4 Collaboration. It is provided under the terms and * |
---|
7 | // * conditions of the Geant4 Software License, included in the file * |
---|
8 | // * LICENSE and available at http://cern.ch/geant4/license . These * |
---|
9 | // * include a list of copyright holders. * |
---|
10 | // * * |
---|
11 | // * Neither the authors of this software system, nor their employing * |
---|
12 | // * institutes,nor the agencies providing financial support for this * |
---|
13 | // * work make any representation or warranty, express or implied, * |
---|
14 | // * regarding this software system or assume any liability for its * |
---|
15 | // * use. Please see the license in the file LICENSE and URL above * |
---|
16 | // * for the full disclaimer and the limitation of liability. * |
---|
17 | // * * |
---|
18 | // * This code implementation is the result of the scientific and * |
---|
19 | // * technical work of the GEANT4 collaboration. * |
---|
20 | // * By using, copying, modifying or distributing the software (or * |
---|
21 | // * any work based on the software) you agree to acknowledge its * |
---|
22 | // * use in resulting scientific publications, and indicate your * |
---|
23 | // * acceptance of all terms of the Geant4 Software license. * |
---|
24 | // ******************************************************************** |
---|
25 | // |
---|
26 | // $Id: G4ICRU73QOModel.cc,v 1.3 2010/06/04 09:08:09 vnivanch Exp $ |
---|
27 | // GEANT4 tag $Name: geant4-09-04-beta-01 $ |
---|
28 | // |
---|
29 | // ------------------------------------------------------------------- |
---|
30 | // |
---|
31 | // GEANT4 Class file |
---|
32 | // |
---|
33 | // |
---|
34 | // File name: G4ICRU73QOModel |
---|
35 | // |
---|
36 | // Author: Alexander Bagulya |
---|
37 | // |
---|
38 | // Creation date: 21.05.2010 |
---|
39 | // |
---|
40 | // Modifications: |
---|
41 | // |
---|
42 | // |
---|
43 | // ------------------------------------------------------------------- |
---|
44 | // |
---|
45 | |
---|
46 | |
---|
47 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
48 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
49 | |
---|
50 | #include "G4ICRU73QOModel.hh" |
---|
51 | #include "Randomize.hh" |
---|
52 | #include "G4Electron.hh" |
---|
53 | #include "G4ParticleChangeForLoss.hh" |
---|
54 | #include "G4LossTableManager.hh" |
---|
55 | #include "G4AntiProton.hh" |
---|
56 | |
---|
57 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
58 | |
---|
59 | using namespace std; |
---|
60 | |
---|
61 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
62 | |
---|
63 | G4ICRU73QOModel::G4ICRU73QOModel(const G4ParticleDefinition* p, const G4String& nam) |
---|
64 | : G4VEmModel(nam), |
---|
65 | particle(0), |
---|
66 | isInitialised(false) |
---|
67 | { |
---|
68 | if(p) SetParticle(p); |
---|
69 | SetHighEnergyLimit(10.0*MeV); |
---|
70 | |
---|
71 | lowestKinEnergy = 10.0*keV; |
---|
72 | |
---|
73 | sizeL0 = 67; |
---|
74 | sizeL1 = 22; |
---|
75 | sizeL2 = 14; |
---|
76 | |
---|
77 | theElectron = G4Electron::Electron(); |
---|
78 | |
---|
79 | for (G4int i = 0; i < 100; ++i) |
---|
80 | { |
---|
81 | indexZ[i] = -1; |
---|
82 | } |
---|
83 | for(G4int i = 0; i < NQOELEM; ++i) |
---|
84 | { |
---|
85 | if(ZElementAvailable[i] > 0) { |
---|
86 | indexZ[ZElementAvailable[i]] = i; |
---|
87 | } |
---|
88 | } |
---|
89 | } |
---|
90 | |
---|
91 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
92 | |
---|
93 | G4ICRU73QOModel::~G4ICRU73QOModel() |
---|
94 | {} |
---|
95 | |
---|
96 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
97 | |
---|
98 | G4double G4ICRU73QOModel::MinEnergyCut(const G4ParticleDefinition*, |
---|
99 | const G4MaterialCutsCouple* ) |
---|
100 | { |
---|
101 | // return couple->GetMaterial()->GetIonisation()->GetMeanExcitationEnergy(); |
---|
102 | return 100*keV; |
---|
103 | } |
---|
104 | |
---|
105 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
106 | |
---|
107 | void G4ICRU73QOModel::Initialise(const G4ParticleDefinition* p, |
---|
108 | const G4DataVector&) |
---|
109 | { |
---|
110 | if(p != particle) SetParticle(p); |
---|
111 | |
---|
112 | // always false before the run |
---|
113 | SetDeexcitationFlag(false); |
---|
114 | |
---|
115 | if(!isInitialised) { |
---|
116 | isInitialised = true; |
---|
117 | |
---|
118 | G4String pname = particle->GetParticleName(); |
---|
119 | fParticleChange = GetParticleChangeForLoss(); |
---|
120 | const G4MaterialTable* mtab = G4Material::GetMaterialTable(); |
---|
121 | denEffData = (*mtab)[0]->GetIonisation()->GetDensityEffectData(); |
---|
122 | } |
---|
123 | } |
---|
124 | |
---|
125 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
126 | |
---|
127 | G4double G4ICRU73QOModel::ComputeCrossSectionPerElectron( |
---|
128 | const G4ParticleDefinition* p, |
---|
129 | G4double kineticEnergy, |
---|
130 | G4double cutEnergy, |
---|
131 | G4double maxKinEnergy) |
---|
132 | { |
---|
133 | G4double cross = 0.0; |
---|
134 | G4double tmax = MaxSecondaryEnergy(p, kineticEnergy); |
---|
135 | G4double maxEnergy = std::min(tmax,maxKinEnergy); |
---|
136 | if(cutEnergy < maxEnergy) { |
---|
137 | |
---|
138 | G4double energy = kineticEnergy + mass; |
---|
139 | G4double energy2 = energy*energy; |
---|
140 | G4double beta2 = kineticEnergy*(kineticEnergy + 2.0*mass)/energy2; |
---|
141 | cross = 1.0/cutEnergy - 1.0/maxEnergy - beta2*log(maxEnergy/cutEnergy)/tmax; |
---|
142 | |
---|
143 | cross *= CLHEP::twopi_mc2_rcl2*chargeSquare/beta2; |
---|
144 | } |
---|
145 | |
---|
146 | return cross; |
---|
147 | } |
---|
148 | |
---|
149 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
150 | |
---|
151 | G4double G4ICRU73QOModel::ComputeCrossSectionPerAtom( |
---|
152 | const G4ParticleDefinition* p, |
---|
153 | G4double kineticEnergy, |
---|
154 | G4double Z, G4double, |
---|
155 | G4double cutEnergy, |
---|
156 | G4double maxEnergy) |
---|
157 | { |
---|
158 | G4double cross = Z*ComputeCrossSectionPerElectron |
---|
159 | (p,kineticEnergy,cutEnergy,maxEnergy); |
---|
160 | return cross; |
---|
161 | } |
---|
162 | |
---|
163 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
164 | |
---|
165 | G4double G4ICRU73QOModel::CrossSectionPerVolume( |
---|
166 | const G4Material* material, |
---|
167 | const G4ParticleDefinition* p, |
---|
168 | G4double kineticEnergy, |
---|
169 | G4double cutEnergy, |
---|
170 | G4double maxEnergy) |
---|
171 | { |
---|
172 | G4double eDensity = material->GetElectronDensity(); |
---|
173 | G4double cross = eDensity*ComputeCrossSectionPerElectron |
---|
174 | (p,kineticEnergy,cutEnergy,maxEnergy); |
---|
175 | return cross; |
---|
176 | } |
---|
177 | |
---|
178 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
179 | |
---|
180 | G4double G4ICRU73QOModel::ComputeDEDXPerVolume(const G4Material* material, |
---|
181 | const G4ParticleDefinition* p, |
---|
182 | G4double kineticEnergy, |
---|
183 | G4double cutEnergy) |
---|
184 | { |
---|
185 | SetParticle(p); |
---|
186 | G4double tmax = MaxSecondaryEnergy(p, kineticEnergy); |
---|
187 | G4double tkin = kineticEnergy/massRate; |
---|
188 | G4double dedx = 0.0; |
---|
189 | if(tkin > lowestKinEnergy) { dedx = DEDX(material, tkin); } |
---|
190 | else { dedx = DEDX(material, lowestKinEnergy)*sqrt(tkin/lowestKinEnergy); } |
---|
191 | |
---|
192 | if (cutEnergy < tmax) { |
---|
193 | |
---|
194 | G4double tau = kineticEnergy/mass; |
---|
195 | G4double gam = tau + 1.0; |
---|
196 | G4double bg2 = tau * (tau+2.0); |
---|
197 | G4double beta2 = bg2/(gam*gam); |
---|
198 | G4double x = cutEnergy/tmax; |
---|
199 | |
---|
200 | dedx += chargeSquare*( log(x) + (1.0 - x)*beta2 ) * twopi_mc2_rcl2 |
---|
201 | * material->GetElectronDensity()/beta2; |
---|
202 | } |
---|
203 | |
---|
204 | return dedx; |
---|
205 | } |
---|
206 | |
---|
207 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
208 | |
---|
209 | G4double G4ICRU73QOModel::DEDX(const G4Material* material, |
---|
210 | G4double kineticEnergy) |
---|
211 | { |
---|
212 | G4double eloss = 0.0; |
---|
213 | const G4int numberOfElements = material->GetNumberOfElements(); |
---|
214 | const G4double* theAtomicNumDensityVector = |
---|
215 | material->GetAtomicNumDensityVector(); |
---|
216 | |
---|
217 | // Bragg's rule calculation |
---|
218 | const G4ElementVector* theElementVector = |
---|
219 | material->GetElementVector() ; |
---|
220 | |
---|
221 | // loop for the elements in the material |
---|
222 | for (G4int i=0; i<numberOfElements; ++i) |
---|
223 | { |
---|
224 | const G4Element* element = (*theElementVector)[i] ; |
---|
225 | eloss += DEDXPerElement(G4int(element->GetZ()), kineticEnergy) |
---|
226 | * theAtomicNumDensityVector[i] * G4int(element->GetZ()); |
---|
227 | } |
---|
228 | return eloss; |
---|
229 | } |
---|
230 | |
---|
231 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
232 | |
---|
233 | G4double G4ICRU73QOModel::DEDXPerElement(G4int AtomicNumber, |
---|
234 | G4double kineticEnergy) |
---|
235 | { |
---|
236 | G4int Z = AtomicNumber; |
---|
237 | if(Z > 97) { Z = 97; } |
---|
238 | G4int nbOfShells = GetNumberOfShells(Z); |
---|
239 | if(nbOfShells < 1) nbOfShells = 1; |
---|
240 | |
---|
241 | G4double v = CLHEP::c_light * std::sqrt( 2.0*kineticEnergy/proton_mass_c2 ); |
---|
242 | |
---|
243 | G4double fBetheVelocity = CLHEP::fine_structure_const*CLHEP::c_light/v; |
---|
244 | |
---|
245 | G4double tau = kineticEnergy/proton_mass_c2; |
---|
246 | G4double gam = tau + 1.0; |
---|
247 | G4double bg2 = tau * (tau+2.0); |
---|
248 | G4double beta2 = bg2/(gam*gam); |
---|
249 | |
---|
250 | G4double l0Term = 0, l1Term = 0, l2Term = 0; |
---|
251 | |
---|
252 | for (G4int nos = 0; nos < nbOfShells; ++nos){ |
---|
253 | |
---|
254 | G4double NormalizedEnergy = (2.0*CLHEP::electron_mass_c2*beta2) / |
---|
255 | GetShellEnergy(Z,nos); |
---|
256 | |
---|
257 | G4double shStrength = GetShellStrength(Z,nos); |
---|
258 | |
---|
259 | G4double l0 = GetL0(NormalizedEnergy); |
---|
260 | l0Term += shStrength * l0; |
---|
261 | |
---|
262 | G4double l1 = GetL1(NormalizedEnergy); |
---|
263 | l1Term += shStrength * l1; |
---|
264 | |
---|
265 | G4double l2 = GetL2(NormalizedEnergy); |
---|
266 | l2Term += shStrength * l2; |
---|
267 | |
---|
268 | } |
---|
269 | G4double dedx = 2*CLHEP::twopi_mc2_rcl2*chargeSquare*factorBethe[Z]* |
---|
270 | (l0Term + q*fBetheVelocity*l1Term |
---|
271 | + chargeSquare*fBetheVelocity*fBetheVelocity*l2Term)/beta2; |
---|
272 | return dedx; |
---|
273 | } |
---|
274 | |
---|
275 | |
---|
276 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
277 | |
---|
278 | G4double G4ICRU73QOModel::GetOscillatorEnergy(G4int Z, |
---|
279 | G4int nbOfTheShell) const |
---|
280 | { |
---|
281 | G4int idx = denEffData->GetElementIndex(Z, kStateUndefined); |
---|
282 | if(idx == -1) { idx = denEffData->GetElementIndex(Z-1, kStateUndefined); } |
---|
283 | G4double PlasmaEnergy = denEffData->GetPlasmaEnergy(idx); |
---|
284 | |
---|
285 | G4double PlasmaEnergy2 = PlasmaEnergy * PlasmaEnergy; |
---|
286 | |
---|
287 | G4double plasmonTerm = 0.66667 * G4AtomicShells::GetNumberOfElectrons(Z,nbOfTheShell) |
---|
288 | * PlasmaEnergy2 / (Z*Z) ; |
---|
289 | |
---|
290 | G4double ionTerm = std::exp(0.5) * (G4AtomicShells::GetBindingEnergy(Z,nbOfTheShell)) ; |
---|
291 | G4double ionTerm2 = ionTerm*ionTerm ; |
---|
292 | |
---|
293 | G4double oscShellEnergy = std::sqrt( ionTerm2 + plasmonTerm ); |
---|
294 | |
---|
295 | return oscShellEnergy; |
---|
296 | } |
---|
297 | |
---|
298 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
299 | |
---|
300 | G4double G4ICRU73QOModel::GetL0(G4double normEnergy) const |
---|
301 | { |
---|
302 | G4int n; |
---|
303 | |
---|
304 | for(n = 0; n < sizeL0; n++) { |
---|
305 | if( normEnergy < L0[n][0] ) break; |
---|
306 | } |
---|
307 | if(0 == n) n = 1 ; |
---|
308 | if(n >= sizeL0) n = sizeL0 - 1 ; |
---|
309 | |
---|
310 | G4double l0 = L0[n][1]; |
---|
311 | G4double l0p = L0[n-1][1]; |
---|
312 | G4double bethe = l0p + (l0 - l0p) * ( normEnergy - L0[n-1][0]) / |
---|
313 | (L0[n][0] - L0[n-1][0]); |
---|
314 | |
---|
315 | return bethe ; |
---|
316 | } |
---|
317 | |
---|
318 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
319 | |
---|
320 | G4double G4ICRU73QOModel::GetL1(G4double normEnergy) const |
---|
321 | { |
---|
322 | G4int n; |
---|
323 | |
---|
324 | for(n = 0; n < sizeL1; n++) { |
---|
325 | if( normEnergy < L1[n][0] ) break; |
---|
326 | } |
---|
327 | if(0 == n) n = 1 ; |
---|
328 | if(n >= sizeL1) n = sizeL1 - 1 ; |
---|
329 | |
---|
330 | G4double l1 = L1[n][1]; |
---|
331 | G4double l1p = L1[n-1][1]; |
---|
332 | G4double barkas= l1p + (l1 - l1p) * ( normEnergy - L1[n-1][0]) / |
---|
333 | (L1[n][0] - L1[n-1][0]); |
---|
334 | |
---|
335 | return barkas; |
---|
336 | } |
---|
337 | |
---|
338 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
339 | |
---|
340 | G4double G4ICRU73QOModel::GetL2(G4double normEnergy) const |
---|
341 | { |
---|
342 | G4int n; |
---|
343 | for(n = 0; n < sizeL2; n++) { |
---|
344 | if( normEnergy < L2[n][0] ) break; |
---|
345 | } |
---|
346 | if(0 == n) n = 1 ; |
---|
347 | if(n >= sizeL2) n = sizeL2 - 1 ; |
---|
348 | |
---|
349 | G4double l2 = L2[n][1]; |
---|
350 | G4double l2p = L2[n-1][1]; |
---|
351 | G4double bloch = l2p + (l2 - l2p) * ( normEnergy - L2[n-1][0]) / |
---|
352 | (L2[n][0] - L2[n-1][0]); |
---|
353 | |
---|
354 | return bloch; |
---|
355 | } |
---|
356 | |
---|
357 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
358 | |
---|
359 | void G4ICRU73QOModel::CorrectionsAlongStep(const G4MaterialCutsCouple*, |
---|
360 | const G4DynamicParticle*, |
---|
361 | G4double&, |
---|
362 | G4double&, |
---|
363 | G4double) |
---|
364 | {} |
---|
365 | |
---|
366 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
367 | |
---|
368 | void G4ICRU73QOModel::SampleSecondaries(std::vector<G4DynamicParticle*>* vdp, |
---|
369 | const G4MaterialCutsCouple*, |
---|
370 | const G4DynamicParticle* dp, |
---|
371 | G4double xmin, |
---|
372 | G4double maxEnergy) |
---|
373 | { |
---|
374 | G4double tmax = MaxSecondaryKinEnergy(dp); |
---|
375 | G4double xmax = std::min(tmax, maxEnergy); |
---|
376 | if(xmin >= xmax) return; |
---|
377 | |
---|
378 | G4double kineticEnergy = dp->GetKineticEnergy(); |
---|
379 | G4double energy = kineticEnergy + mass; |
---|
380 | G4double energy2 = energy*energy; |
---|
381 | G4double beta2 = kineticEnergy*(kineticEnergy + 2.0*mass)/energy2; |
---|
382 | G4double grej = 1.0; |
---|
383 | G4double deltaKinEnergy, f; |
---|
384 | |
---|
385 | G4ThreeVector direction = dp->GetMomentumDirection(); |
---|
386 | |
---|
387 | // sampling follows ... |
---|
388 | do { |
---|
389 | G4double q = G4UniformRand(); |
---|
390 | deltaKinEnergy = xmin*xmax/(xmin*(1.0 - q) + xmax*q); |
---|
391 | |
---|
392 | f = 1.0 - beta2*deltaKinEnergy/tmax; |
---|
393 | |
---|
394 | if(f > grej) { |
---|
395 | G4cout << "G4ICRU73QOModel::SampleSecondary Warning! " |
---|
396 | << "Majorant " << grej << " < " |
---|
397 | << f << " for e= " << deltaKinEnergy |
---|
398 | << G4endl; |
---|
399 | } |
---|
400 | |
---|
401 | } while( grej*G4UniformRand() >= f ); |
---|
402 | |
---|
403 | G4double deltaMomentum = |
---|
404 | sqrt(deltaKinEnergy * (deltaKinEnergy + 2.0*electron_mass_c2)); |
---|
405 | G4double totMomentum = energy*sqrt(beta2); |
---|
406 | G4double cost = deltaKinEnergy * (energy + electron_mass_c2) / |
---|
407 | (deltaMomentum * totMomentum); |
---|
408 | if(cost > 1.0) cost = 1.0; |
---|
409 | G4double sint = sqrt((1.0 - cost)*(1.0 + cost)); |
---|
410 | |
---|
411 | G4double phi = twopi * G4UniformRand() ; |
---|
412 | |
---|
413 | G4ThreeVector deltaDirection(sint*cos(phi),sint*sin(phi), cost) ; |
---|
414 | deltaDirection.rotateUz(direction); |
---|
415 | |
---|
416 | // Change kinematics of primary particle |
---|
417 | kineticEnergy -= deltaKinEnergy; |
---|
418 | G4ThreeVector finalP = direction*totMomentum - deltaDirection*deltaMomentum; |
---|
419 | finalP = finalP.unit(); |
---|
420 | |
---|
421 | fParticleChange->SetProposedKineticEnergy(kineticEnergy); |
---|
422 | fParticleChange->SetProposedMomentumDirection(finalP); |
---|
423 | |
---|
424 | // create G4DynamicParticle object for delta ray |
---|
425 | G4DynamicParticle* delta = new G4DynamicParticle(theElectron,deltaDirection, |
---|
426 | deltaKinEnergy); |
---|
427 | |
---|
428 | vdp->push_back(delta); |
---|
429 | } |
---|
430 | |
---|
431 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
432 | |
---|
433 | G4double G4ICRU73QOModel::MaxSecondaryEnergy(const G4ParticleDefinition* pd, |
---|
434 | G4double kinEnergy) |
---|
435 | { |
---|
436 | if(pd != particle) SetParticle(pd); |
---|
437 | G4double tau = kinEnergy/mass; |
---|
438 | G4double tmax = 2.0*electron_mass_c2*tau*(tau + 2.) / |
---|
439 | (1. + 2.0*(tau + 1.)*ratio + ratio*ratio); |
---|
440 | return tmax; |
---|
441 | } |
---|
442 | |
---|
443 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
444 | |
---|
445 | const G4int G4ICRU73QOModel::ZElementAvailable[NQOELEM] = {1,2,4,6,7,8,10,13,14,-18, |
---|
446 | 22,26,28,29,32,36,42,47, |
---|
447 | 50,54,73,74,78,79,82,92}; |
---|
448 | |
---|
449 | const G4int G4ICRU73QOModel::nbofShellsForElement[NQOELEM] = {1,1,2,3,3,3,3,4,5,4, |
---|
450 | 5,5,5,5,6,4,6,6, |
---|
451 | 7,6,6,8,7,7,9,9}; |
---|
452 | |
---|
453 | const G4int G4ICRU73QOModel::startElemIndex[NQOELEM] = {0,1,2,4,7,10,13,16,20,25, |
---|
454 | 29,34,39,44,49,55,59,65, |
---|
455 | 71,78,84,90,98,105,112,121}; |
---|
456 | |
---|
457 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
458 | |
---|
459 | // SubShellOccupation = Z * ShellStrength |
---|
460 | const G4double G4ICRU73QOModel::SubShellOccupation[NQODATA] = |
---|
461 | { |
---|
462 | 1.000, // H 0 |
---|
463 | 2.000, // He 1 |
---|
464 | 1.930, 2.070, // Be 2-3 |
---|
465 | 1.992, 1.841, 2.167, // C 4-6 |
---|
466 | 1.741, 1.680, 3.579, // N 7-9 |
---|
467 | 1.802, 1.849, 4.349, // O 10-12 |
---|
468 | 1.788, 2.028, 6.184, // Ne 13-15 |
---|
469 | 1.623, 2.147, 6.259, 2.971, // Al 16-19 |
---|
470 | 1.631, 2.094, 6.588, 2.041, 1.646, // Si 20-24 |
---|
471 | 1.535, 8.655, 1.706, 6.104, // Ar 25-28 |
---|
472 | 1.581, 8.358, 8.183, 2.000, 1.878, // Ti 29-33 |
---|
473 | 1.516, 8.325, 8.461, 6.579, 1.119, // Fe 34-38 |
---|
474 | 1.422, 7.81, 8.385, 8.216, 2.167, // Ni 39-43 |
---|
475 | 1.458, 8.049, 8.79, 9.695, 1.008, // Cu 44-48 |
---|
476 | 1.442, 7.791, 7.837, 10.122, 2.463, 2.345, // Ge 49-54 |
---|
477 | 1.645, 7.765, 19.192, 7.398, // Kr 55-58 |
---|
478 | 1.313, 6.409, 19.229, 8.633, 5.036, 1.380, // Mo 59-64 |
---|
479 | 1.295, 6.219, 18.751, 8.748, 10.184, 1.803, // Ag 65-70 |
---|
480 | 1.277, 6.099, 20.386, 8.011, 10.007, 2.272, 1.948, // Sn 71-77 |
---|
481 | 1.563, 6.312, 21.868, 5.762, 11.245, 7.250, // Xe 78-83 |
---|
482 | 0.9198, 6.5408, 18.9727, 24.9149, 15.0161, 6.6284, // Ta 84-89 |
---|
483 | 1.202, 5.582, 19.527, 18.741, 8.411, 14.387, 4.042, 2.108, // W 90-97 |
---|
484 | 1.159, 5.467, 18.802, 33.905, 8.300, 9.342, 1.025, // Pt 98-104 |
---|
485 | 1.124, 5.331, 18.078, 34.604, 8.127, 10.414, 1.322, // Au 105-111 |
---|
486 | 2.000, 8.000, 18.000, 18.000, 14.000, 8.000, 10.000, 2.000, 2.000, // Pb 112-120 |
---|
487 | 2.000, 8.000, 18.000, 32.000, 18.000, 8.000, 2.000, 1.000, 3.000 // U 121-129 |
---|
488 | }; |
---|
489 | |
---|
490 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
491 | |
---|
492 | // ShellEnergy in eV |
---|
493 | const G4double G4ICRU73QOModel::ShellEnergy[NQODATA] = |
---|
494 | { |
---|
495 | 19.2, // H |
---|
496 | 41.8, // He |
---|
497 | 209.11, 21.68, // Be |
---|
498 | 486.2, 60.95, 23.43, // C |
---|
499 | 732.61, 100.646, 23.550, // N |
---|
500 | 965.1, 129.85, 31.60, // O |
---|
501 | 1525.9, 234.9, 56.18, // Ne |
---|
502 | 2701, 476.5, 150.42, 16.89, // Al |
---|
503 | 3206.1, 586.4, 186.8, 23.52, 14.91, // Si |
---|
504 | 5551.6, 472.43, 124.85, 22.332, // Ar |
---|
505 | 8554.6, 850.58, 93.47, 39.19, 19.46, // Ti |
---|
506 | 12254.7, 1279.29, 200.35, 49.19, 17.66, // Fe |
---|
507 | 14346.9, 1532.28, 262.71, 74.37, 23.03, // Ni |
---|
508 | 15438.5, 1667.96, 294.1, 70.69, 16.447, // Cu |
---|
509 | 19022.1, 2150.79, 455.79, 179.87, 57.89, 20.95, // Ge |
---|
510 | 24643, 2906.4, 366.85, 22.24, // Kr |
---|
511 | 34394, 4365.3, 589.36, 129.42, 35.59, 18.42, // Mo |
---|
512 | 43664.3, 5824.91, 909.79, 175.47, 54.89, 19.63, // Ag |
---|
513 | 49948, 6818.2, 1036.1, 172.65, 70.89, 33.87, 14.54, // Sn |
---|
514 | 58987, 8159, 1296.6, 356.75, 101.03, 16.52, // Xe |
---|
515 | 88926, 18012, 3210, 575, 108.7, 30.8, // Ta |
---|
516 | 115025.9, 17827.44, 3214.36, 750.41, 305.21, 105.50, 38.09, 21.25, // W |
---|
517 | 128342, 20254, 3601.8, 608.1, 115.0, 42.75, 17.04, // Pt |
---|
518 | 131872, 20903, 3757.4, 682.1, 105.2, 44.89, 17.575, // Au |
---|
519 | 154449, 25067, 5105.0, 987.44, 247.59, 188.1, 40.61, 19.2, 15.17, // Pb |
---|
520 | 167282, 27868, 6022.7, 1020.4, 244.81, 51.33, 13, 11.06, 14.43 // U |
---|
521 | }; |
---|
522 | |
---|
523 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
524 | |
---|
525 | // Data for L0 from: Sigmund P., Haagerup U. Phys. Rev. A34 (1986) 892-910 |
---|
526 | const G4double G4ICRU73QOModel::L0[67][2] = |
---|
527 | { |
---|
528 | {0.00, 0.000001}, |
---|
529 | {0.10, 0.000001}, |
---|
530 | {0.12, 0.00001}, |
---|
531 | {0.14, 0.00005}, |
---|
532 | {0.16, 0.00014}, |
---|
533 | {0.18, 0.00030}, |
---|
534 | {0.20, 0.00057}, |
---|
535 | {0.25, 0.00189}, |
---|
536 | {0.30, 0.00429}, |
---|
537 | {0.35, 0.00784}, |
---|
538 | {0.40, 0.01248}, |
---|
539 | {0.45, 0.01811}, |
---|
540 | {0.50, 0.02462}, |
---|
541 | {0.60, 0.03980}, |
---|
542 | {0.70, 0.05731}, |
---|
543 | {0.80, 0.07662}, |
---|
544 | {0.90, 0.09733}, |
---|
545 | {1.00, 0.11916}, |
---|
546 | {1.20, 0.16532}, |
---|
547 | {1.40, 0.21376}, |
---|
548 | {1.60, 0.26362}, |
---|
549 | {1.80, 0.31428}, |
---|
550 | {2.00, 0.36532}, |
---|
551 | {2.50, 0.49272}, |
---|
552 | {3.00, 0.61765}, |
---|
553 | {3.50, 0.73863}, |
---|
554 | {4.00, 0.85496}, |
---|
555 | {4.50, 0.96634}, |
---|
556 | {5.00, 1.07272}, |
---|
557 | {6.00, 1.27086}, |
---|
558 | {7.00, 1.45075}, |
---|
559 | {8.00, 1.61412}, |
---|
560 | {9.00, 1.76277}, |
---|
561 | {10.00, 1.89836}, |
---|
562 | {12.00, 2.13625}, |
---|
563 | {14.00, 2.33787}, |
---|
564 | {16.00, 2.51093}, |
---|
565 | {18.00, 2.66134}, |
---|
566 | {20.00, 2.79358}, |
---|
567 | {25.00, 3.06539}, |
---|
568 | {30.00, 3.27902}, |
---|
569 | {35.00, 3.45430}, |
---|
570 | {40.00, 3.60281}, |
---|
571 | {45.00, 3.73167}, |
---|
572 | {50.00, 3.84555}, |
---|
573 | {60.00, 4.04011}, |
---|
574 | {70.00, 4.20264}, |
---|
575 | {80.00, 4.34229}, |
---|
576 | {90.00, 4.46474}, |
---|
577 | {100.00, 4.57378}, |
---|
578 | {120.00, 4.76155}, |
---|
579 | {140.00, 4.91953}, |
---|
580 | {160.00, 5.05590}, |
---|
581 | {180.00, 5.17588}, |
---|
582 | {200.00, 5.28299}, |
---|
583 | {250.00, 5.50925}, |
---|
584 | {300.00, 5.69364}, |
---|
585 | {350.00, 5.84926}, |
---|
586 | {400.00, 5.98388}, |
---|
587 | {450.00, 6.10252}, |
---|
588 | {500.00, 6.20856}, |
---|
589 | {600.00, 6.39189}, |
---|
590 | {700.00, 6.54677}, |
---|
591 | {800.00, 6.68084}, |
---|
592 | {900.00, 6.79905}, |
---|
593 | {1000.00, 6.90474} |
---|
594 | }; |
---|
595 | |
---|
596 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
597 | |
---|
598 | // Data for L1 from: Mikkelsen H.H., Sigmund P. Phys. Rev. A40 (1989) 101-116 |
---|
599 | const G4double G4ICRU73QOModel::L1[22][2] = |
---|
600 | { |
---|
601 | {0.00, -0.000001}, |
---|
602 | {0.10, -0.00001}, |
---|
603 | {0.20, -0.00049}, |
---|
604 | {0.30, -0.00084}, |
---|
605 | {0.40, 0.00085}, |
---|
606 | {0.50, 0.00519}, |
---|
607 | {0.60, 0.01198}, |
---|
608 | {0.70, 0.02074}, |
---|
609 | {0.80, 0.03133}, |
---|
610 | {0.90, 0.04369}, |
---|
611 | {1.00, 0.06035}, |
---|
612 | {2.00, 0.24023}, |
---|
613 | {3.00, 0.44284}, |
---|
614 | {4.00, 0.62012}, |
---|
615 | {5.00, 0.77031}, |
---|
616 | {6.00, 0.90390}, |
---|
617 | {7.00, 1.02705}, |
---|
618 | {8.00, 1.10867}, |
---|
619 | {9.00, 1.17546}, |
---|
620 | {10.00, 1.21599}, |
---|
621 | {15.00, 1.24349}, |
---|
622 | {20.00, 1.16752} |
---|
623 | }; |
---|
624 | |
---|
625 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
626 | |
---|
627 | // Data for L2 from: Mikkelsen H.H. Nucl. Instr. Meth. B58 (1991) 136-148 |
---|
628 | const G4double G4ICRU73QOModel::L2[14][2] = |
---|
629 | { |
---|
630 | {0.00, 0.000001}, |
---|
631 | {0.10, 0.00001}, |
---|
632 | {0.20, 0.00000}, |
---|
633 | {0.40, -0.00120}, |
---|
634 | {0.60, -0.00036}, |
---|
635 | {0.80, 0.00372}, |
---|
636 | {1.00, 0.01298}, |
---|
637 | {2.00, 0.08296}, |
---|
638 | {4.00, 0.21953}, |
---|
639 | {6.00, 0.23903}, |
---|
640 | {8.00, 0.20893}, |
---|
641 | {10.00, 0.10879}, |
---|
642 | {20.00, -0.88409}, |
---|
643 | {40.00, -1.13902} |
---|
644 | }; |
---|
645 | |
---|
646 | //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... |
---|
647 | |
---|
648 | // Correction obtained by V.Ivanchenko using G4BetheBlochModel |
---|
649 | const G4double G4ICRU73QOModel::factorBethe[99] = { 1.0, |
---|
650 | 0.9637, 0.9877, 0.9467, 0.9856, 0.9066, 0.9886, 0.9256, 0.9822, 0.8729, 0.9965, // 1 - 10 |
---|
651 | 0.8564, 0.8879, 0.9577, 0.9606, 0.9104, 0.8962, 0.9316, 0.9378, 0.9189, 0.9543, // 11 - 20 |
---|
652 | 0.9114, 0.9846, 0.8497, 0.8311, 0.8081, 0.985, 0.76, 0.9696, 0.9964, 0.7361, // 21 - 30 |
---|
653 | 0.7787, 1.066, 0.8547, 0.8759, 0.9135, 1.042, 0.9256, 0.9555, 0.965, 0.9619, // 31 - 40 |
---|
654 | 0.9361, 0.9078, 0.9291, 0.9128, 0.9115, 0.876, 0.8701, 0.7955, 0.8197, 0.8551, // 41 - 50 |
---|
655 | 0.8318, 0.8646, 0.8802, 0.8871, 0.8881, 0.9103, 0.9113, 0.8879, 0.8588, 0.8279, // 51 - 60 |
---|
656 | 0.7953, 0.7685, 0.7433, 0.7219, 0.7003, 0.6732, 0.6485, 0.63, 0.6232, 0.611, // 61 - 70 |
---|
657 | 0.6385, 0.6642, 0.8452, 0.8127, 0.6982, 0.6971, 0.6889, 0.8416, 0.8341, 0.6956, // 71 - 80 |
---|
658 | 0.7308, 0.7777, 0.7692, 0.7924, 0.8054, 0.8203, 0.8386, 0.8618, 0.8654, 0.8558, // 81 - 90 |
---|
659 | 0.8323, 0.8362, 0.8044, 0.7811, 0.7685, 0.7714, 0.759, 0.7506}; |
---|