Geant4 11.1.1
Toolkit for the simulation of the passage of particles through matter
•All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4ICRU73QOModel.cc
Go to the documentation of this file.
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// -------------------------------------------------------------------
28//
29// GEANT4 Class file
30//
31//
32// File name: G4ICRU73QOModel
33//
34// Author: Alexander Bagulya
35//
36// Creation date: 21.05.2010
37//
38// Modifications:
39//
40//
41// -------------------------------------------------------------------
42//
43//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
44//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
45
46#include "G4ICRU73QOModel.hh"
47#include "G4PhysicalConstants.hh"
48#include "G4SystemOfUnits.hh"
49#include "Randomize.hh"
50#include "G4Electron.hh"
51#include "G4ParticleChangeForLoss.hh"
52#include "G4LossTableManager.hh"
53#include "G4AntiProton.hh"
54#include "G4DeltaAngle.hh"
55#include "G4Log.hh"
56#include "G4Exp.hh"
57
58//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
59
60using namespace std;
61
62//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
63
64G4ICRU73QOModel::G4ICRU73QOModel(const G4ParticleDefinition* p, const G4String& nam)
65 : G4VEmModel(nam),
66 particle(nullptr),
67 isInitialised(false)
68{
69 mass = charge = chargeSquare = massRate = ratio = 0.0;
70 if(p) { SetParticle(p); }
71 SetHighEnergyLimit(10.0*MeV);
72
73 lowestKinEnergy = 5.0*keV;
74
75 sizeL0 = 67;
76 sizeL1 = 22;
77 sizeL2 = 14;
78
79 theElectron = G4Electron::Electron();
80
81 for (G4int i = 0; i < 100; ++i)
82 {
83 indexZ[i] = -1;
84 }
85 for(G4int i = 0; i < NQOELEM; ++i)
86 {
87 if(ZElementAvailable[i] > 0) {
88 indexZ[ZElementAvailable[i]] = i;
89 }
90 }
91 fParticleChange = nullptr;
92 denEffData = nullptr;
93}
94
95//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
96
97void G4ICRU73QOModel::Initialise(const G4ParticleDefinition* p,
98 const G4DataVector&)
99{
100 if(p != particle) SetParticle(p);
101
102 // always false before the run
103 SetDeexcitationFlag(false);
104
105 if(!isInitialised) {
106 isInitialised = true;
107
108 if(UseAngularGeneratorFlag() && !GetAngularDistribution()) {
109 SetAngularDistribution(new G4DeltaAngle());
110 }
111
112 G4String pname = particle->GetParticleName();
113 fParticleChange = GetParticleChangeForLoss();
114 const G4MaterialTable* mtab = G4Material::GetMaterialTable();
115 denEffData = (*mtab)[0]->GetIonisation()->GetDensityEffectData();
116 }
117}
118
119//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
120
121G4double G4ICRU73QOModel::ComputeCrossSectionPerElectron(
122 const G4ParticleDefinition* p,
123 G4double kineticEnergy,
124 G4double cut,
125 G4double maxKinEnergy)
126{
127 G4double cross = 0.0;
128 const G4double tmax = MaxSecondaryEnergy(p, kineticEnergy);
129 const G4double maxEnergy = std::min(tmax, maxKinEnergy);
130 const G4double cutEnergy = std::max(cut, lowestKinEnergy*massRate);
131 if(cutEnergy < maxEnergy) {
132
133 const G4double energy = kineticEnergy + mass;
134 const G4double energy2 = energy*energy;
135 const G4double beta2 = kineticEnergy*(kineticEnergy + 2.0*mass)/energy2;
136 cross = (maxEnergy - cutEnergy)/(cutEnergy*maxEnergy)
137 - beta2*G4Log(maxEnergy/cutEnergy)/tmax;
138
139 cross *= CLHEP::twopi_mc2_rcl2*chargeSquare/beta2;
140 }
141
142 return cross;
143}
144
145//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
146
147G4double G4ICRU73QOModel::ComputeCrossSectionPerAtom(
148 const G4ParticleDefinition* p,
149 G4double kineticEnergy,
150 G4double Z, G4double,
151 G4double cutEnergy,
152 G4double maxEnergy)
153{
154 G4double cross = Z*ComputeCrossSectionPerElectron
155 (p,kineticEnergy,cutEnergy,maxEnergy);
156 return cross;
157}
158
159//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
160
161G4double G4ICRU73QOModel::CrossSectionPerVolume(
162 const G4Material* material,
163 const G4ParticleDefinition* p,
164 G4double kineticEnergy,
165 G4double cutEnergy,
166 G4double maxEnergy)
167{
168 G4double eDensity = material->GetElectronDensity();
169 G4double cross = eDensity*ComputeCrossSectionPerElectron
170 (p,kineticEnergy,cutEnergy,maxEnergy);
171 return cross;
172}
173
174//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
175
176G4double G4ICRU73QOModel::ComputeDEDXPerVolume(const G4Material* material,
177 const G4ParticleDefinition* p,
178 G4double kineticEnergy,
179 G4double cut)
180{
181 SetParticle(p);
182 const G4double tmax = MaxSecondaryEnergy(p, kineticEnergy);
183 const G4double tkin = kineticEnergy/massRate;
184 const G4double cutEnergy = std::max(cut, lowestKinEnergy*massRate);
185 G4double dedx = 0.0;
186 if(tkin > lowestKinEnergy) { dedx = DEDX(material, tkin); }
187 else { dedx = DEDX(material, lowestKinEnergy)*sqrt(tkin/lowestKinEnergy); }
188
189 if (cutEnergy < tmax) {
190
191 const G4double tau = kineticEnergy/mass;
192 const G4double x = cutEnergy/tmax;
193 dedx += (G4Log(x)*(tau + 1.)*(tau + 1.)/(tau * (tau + 2.0)) + 1.0 - x) *
194 CLHEP::twopi_mc2_rcl2 *chargeSquare * material->GetElectronDensity();
195 }
196 dedx = std::max(dedx, 0.0);
197 return dedx;
198}
199
200//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
201
202G4double G4ICRU73QOModel::DEDX(const G4Material* material,
203 G4double kineticEnergy)
204{
205 G4double eloss = 0.0;
206 const std::size_t numberOfElements = material->GetNumberOfElements();
207 const G4double* theAtomicNumDensityVector =
208 material->GetAtomicNumDensityVector();
209
210 // Bragg's rule calculation
211 const G4ElementVector* theElementVector =
212 material->GetElementVector() ;
213
214 // loop for the elements in the material
215 for (std::size_t i=0; i<numberOfElements; ++i)
216 {
217 const G4Element* element = (*theElementVector)[i] ;
218 eloss += DEDXPerElement(element->GetZasInt(), kineticEnergy)
219 * theAtomicNumDensityVector[i] * element->GetZ();
220 }
221 return eloss;
222}
223
224//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
225
226G4double G4ICRU73QOModel::DEDXPerElement(G4int AtomicNumber,
227 G4double kineticEnergy)
228{
229 G4int Z = std::min(AtomicNumber, 97);
230 G4int nbOfShells = std::max(GetNumberOfShells(Z), 1);
231
232 G4double v = CLHEP::c_light * std::sqrt( 2.0*kineticEnergy/proton_mass_c2 );
233
234 G4double fBetheVelocity = CLHEP::fine_structure_const*CLHEP::c_light/v;
235
236 G4double tau = kineticEnergy/proton_mass_c2;
237 G4double gam = tau + 1.0;
238 G4double bg2 = tau * (tau+2.0);
239 G4double beta2 = bg2/(gam*gam);
240
241 G4double l0Term = 0, l1Term = 0, l2Term = 0;
242
243 for (G4int nos = 0; nos < nbOfShells; ++nos){
244
245 G4double NormalizedEnergy = (2.0*CLHEP::electron_mass_c2*beta2) /
246 GetShellEnergy(Z,nos);
247
248 G4double shStrength = GetShellStrength(Z,nos);
249
250 G4double l0 = GetL0(NormalizedEnergy);
251 l0Term += shStrength * l0;
252
253 G4double l1 = GetL1(NormalizedEnergy);
254 l1Term += shStrength * l1;
255
256 G4double l2 = GetL2(NormalizedEnergy);
257 l2Term += shStrength * l2;
258
259 }
260 G4double dedx = 2*CLHEP::twopi_mc2_rcl2*chargeSquare*factorBethe[Z]*
261 (l0Term + charge*fBetheVelocity*l1Term
262 + chargeSquare*fBetheVelocity*fBetheVelocity*l2Term)/beta2;
263 return dedx;
264}
265
266
267//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
268
269G4double G4ICRU73QOModel::GetOscillatorEnergy(G4int Z,
270 G4int nbOfTheShell) const
271{
272 G4int idx = denEffData->GetElementIndex(Z, kStateUndefined);
273 if(idx == -1) { idx = denEffData->GetElementIndex(Z-1, kStateUndefined); }
274 G4double PlasmaEnergy = denEffData->GetPlasmaEnergy(idx);
275
276 G4double PlasmaEnergy2 = PlasmaEnergy * PlasmaEnergy;
277
278 G4double plasmonTerm = 0.66667
279 * G4AtomicShells::GetNumberOfElectrons(Z,nbOfTheShell)
280 * PlasmaEnergy2 / (Z*Z) ;
281
282 static const G4double exphalf = G4Exp(0.5);
283 G4double ionTerm = exphalf *
284 (G4AtomicShells::GetBindingEnergy(Z,nbOfTheShell)) ;
285 G4double ionTerm2 = ionTerm*ionTerm ;
286
287 G4double oscShellEnergy = std::sqrt( ionTerm2 + plasmonTerm );
288
289 return oscShellEnergy;
290}
291
292//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
293
294G4int G4ICRU73QOModel::GetNumberOfShells(G4int Z) const
295{
296 G4int nShell = 0;
297
298 if(indexZ[Z] >= 0) {
299 nShell = nbofShellsForElement[indexZ[Z]];
300 } else {
301 nShell = G4AtomicShells::GetNumberOfShells(Z);
302 }
303 return nShell;
304}
305
306//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
307
308G4double G4ICRU73QOModel::GetShellEnergy(G4int Z, G4int nbOfTheShell) const
309{
310 G4double shellEnergy = 0.;
311
312 G4int idx = indexZ[Z];
313
314 if(idx >= 0) {
315 shellEnergy = ShellEnergy[startElemIndex[idx] + nbOfTheShell]*CLHEP::eV;
316 } else {
317 shellEnergy = GetOscillatorEnergy(Z, nbOfTheShell);
318 }
319
320 return shellEnergy;
321}
322
323//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
324
325G4double G4ICRU73QOModel::GetShellStrength(G4int Z, G4int nbOfTheShell) const
326{
327 G4double shellStrength = 0.;
328
329 G4int idx = indexZ[Z];
330
331 if(idx >= 0) {
332 shellStrength = SubShellOccupation[startElemIndex[idx] + nbOfTheShell] / Z;
333 } else {
334 shellStrength = G4double(G4AtomicShells::GetNumberOfElectrons(Z,nbOfTheShell))/Z;
335 }
336
337 return shellStrength;
338}
339
340//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
341
342G4double G4ICRU73QOModel::GetL0(G4double normEnergy) const
343{
344 G4int n;
345
346 for(n = 0; n < sizeL0; n++) {
347 if( normEnergy < L0[n][0] ) break;
348 }
349 if(0 == n) { n = 1; }
350 if(n >= sizeL0) { n = sizeL0 - 1; }
351
352 G4double l0 = L0[n][1];
353 G4double l0p = L0[n-1][1];
354 G4double bethe = l0p + (l0 - l0p) * ( normEnergy - L0[n-1][0]) /
355 (L0[n][0] - L0[n-1][0]);
356
357 return bethe ;
358}
359
360//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
361
362G4double G4ICRU73QOModel::GetL1(G4double normEnergy) const
363{
364 G4int n;
365
366 for(n = 0; n < sizeL1; n++) {
367 if( normEnergy < L1[n][0] ) break;
368 }
369 if(0 == n) n = 1 ;
370 if(n >= sizeL1) n = sizeL1 - 1 ;
371
372 G4double l1 = L1[n][1];
373 G4double l1p = L1[n-1][1];
374 G4double barkas= l1p + (l1 - l1p) * ( normEnergy - L1[n-1][0]) /
375 (L1[n][0] - L1[n-1][0]);
376
377 return barkas;
378}
379
380//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
381
382G4double G4ICRU73QOModel::GetL2(G4double normEnergy) const
383{
384 G4int n;
385 for(n = 0; n < sizeL2; n++) {
386 if( normEnergy < L2[n][0] ) break;
387 }
388 if(0 == n) n = 1 ;
389 if(n >= sizeL2) n = sizeL2 - 1 ;
390
391 G4double l2 = L2[n][1];
392 G4double l2p = L2[n-1][1];
393 G4double bloch = l2p + (l2 - l2p) * ( normEnergy - L2[n-1][0]) /
394 (L2[n][0] - L2[n-1][0]);
395
396 return bloch;
397}
398
399//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
400
401void G4ICRU73QOModel::CorrectionsAlongStep(const G4MaterialCutsCouple*,
402 const G4DynamicParticle*,
403 const G4double&,
404 G4double&)
405{}
406
407//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
408
409void G4ICRU73QOModel::SampleSecondaries(std::vector<G4DynamicParticle*>* vdp,
410 const G4MaterialCutsCouple* couple,
411 const G4DynamicParticle* dp,
412 G4double minEnergy,
413 G4double maxEnergy)
414{
415 const G4double tmax = MaxSecondaryKinEnergy(dp);
416 const G4double xmax = std::min(tmax, maxEnergy);
417 const G4double xmin = std::max(minEnergy, lowestKinEnergy*massRate);
418 if(xmin >= xmax) { return; }
419
420 G4double kineticEnergy = dp->GetKineticEnergy();
421 const G4double energy = kineticEnergy + mass;
422 const G4double energy2 = energy*energy;
423 const G4double beta2 = kineticEnergy*(kineticEnergy + 2.0*mass)/energy2;
424 G4double grej = 1.0;
425 G4double deltaKinEnergy, f;
426
427 G4ThreeVector direction = dp->GetMomentumDirection();
428
429 // sampling follows ...
430 do {
431 G4double x = G4UniformRand();
432 deltaKinEnergy = xmin*xmax/(xmin*(1.0 - x) + xmax*x);
433
434 f = 1.0 - beta2*deltaKinEnergy/tmax;
435
436 if(f > grej) {
437 G4cout << "G4ICRU73QOModel::SampleSecondary Warning! "
438 << "Majorant " << grej << " < "
439 << f << " for e= " << deltaKinEnergy
440 << G4endl;
441 }
442
443 // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
444 } while( grej*G4UniformRand() >= f );
445
446 G4ThreeVector deltaDirection;
447
448 if(UseAngularGeneratorFlag()) {
449 const G4Material* mat = couple->GetMaterial();
450 G4int Z = SelectRandomAtomNumber(mat);
451
452 deltaDirection =
453 GetAngularDistribution()->SampleDirection(dp, deltaKinEnergy, Z, mat);
454
455 } else {
456
457 G4double deltaMomentum =
458 sqrt(deltaKinEnergy * (deltaKinEnergy + 2.0*electron_mass_c2));
459 G4double totMomentum = energy*sqrt(beta2);
460 G4double cost = deltaKinEnergy * (energy + electron_mass_c2) /
461 (deltaMomentum * totMomentum);
462 if(cost > 1.0) { cost = 1.0; }
463 G4double sint = sqrt((1.0 - cost)*(1.0 + cost));
464
465 G4double phi = twopi * G4UniformRand() ;
466
467 deltaDirection.set(sint*cos(phi),sint*sin(phi), cost) ;
468 deltaDirection.rotateUz(direction);
469 }
470 // create G4DynamicParticle object for delta ray
471 auto delta = new G4DynamicParticle(theElectron,deltaDirection,deltaKinEnergy);
472
473 // Change kinematics of primary particle
474 kineticEnergy -= deltaKinEnergy;
475 G4ThreeVector finalP = dp->GetMomentum() - delta->GetMomentum();
476 finalP = finalP.unit();
477
478 fParticleChange->SetProposedKineticEnergy(kineticEnergy);
479 fParticleChange->SetProposedMomentumDirection(finalP);
480
481 vdp->push_back(delta);
482}
483
484//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
485
486G4double G4ICRU73QOModel::MaxSecondaryEnergy(const G4ParticleDefinition* pd,
487 G4double kinEnergy)
488{
489 if(pd != particle) { SetParticle(pd); }
490 G4double tau = kinEnergy/mass;
491 G4double tmax = 2.0*electron_mass_c2*tau*(tau + 2.) /
492 (1. + 2.0*(tau + 1.)*ratio + ratio*ratio);
493 return tmax;
494}
495
496//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
497
498const G4int G4ICRU73QOModel::ZElementAvailable[NQOELEM] =
499 {1,2,4,6,7,8,10,13,14,-18,
500 22,26,28,29,32,36,42,47,
501 50,54,73,74,78,79,82,92};
502
503const G4int G4ICRU73QOModel::nbofShellsForElement[NQOELEM] =
504 {1,1,2,3,3,3,3,4,5,4,
505 5,5,5,5,6,4,6,6,
506 7,6,6,8,7,7,9,9};
507
508const G4int G4ICRU73QOModel::startElemIndex[NQOELEM] =
509 {0,1,2,4,7,10,13,16,20,25,
510 29,34,39,44,49,55,59,65,
511 71,78,84,90,98,105,112,121};
512
513//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
514
515// SubShellOccupation = Z * ShellStrength
516const G4double G4ICRU73QOModel::SubShellOccupation[NQODATA] =
517 {
518 1.000, // H 0
519 2.000, // He 1
520 1.930, 2.070, // Be 2-3
521 1.992, 1.841, 2.167, // C 4-6
522 1.741, 1.680, 3.579, // N 7-9
523 1.802, 1.849, 4.349, // O 10-12
524 1.788, 2.028, 6.184, // Ne 13-15
525 1.623, 2.147, 6.259, 2.971, // Al 16-19
526 1.631, 2.094, 6.588, 2.041, 1.646, // Si 20-24
527 1.535, 8.655, 1.706, 6.104, // Ar 25-28
528 1.581, 8.358, 8.183, 2.000, 1.878, // Ti 29-33
529 1.516, 8.325, 8.461, 6.579, 1.119, // Fe 34-38
530 1.422, 7.81, 8.385, 8.216, 2.167, // Ni 39-43
531 1.458, 8.049, 8.79, 9.695, 1.008, // Cu 44-48
532 1.442, 7.791, 7.837, 10.122, 2.463, 2.345, // Ge 49-54
533 1.645, 7.765, 19.192, 7.398, // Kr 55-58
534 1.313, 6.409, 19.229, 8.633, 5.036, 1.380, // Mo 59-64
535 1.295, 6.219, 18.751, 8.748, 10.184, 1.803, // Ag 65-70
536 1.277, 6.099, 20.386, 8.011, 10.007, 2.272, 1.948, // Sn 71-77
537 1.563, 6.312, 21.868, 5.762, 11.245, 7.250, // Xe 78-83
538 0.9198, 6.5408, 18.9727, 24.9149, 15.0161, 6.6284, // Ta 84-89
539 1.202, 5.582, 19.527, 18.741, 8.411, 14.387, 4.042, 2.108, // W 90-97
540 1.159, 5.467, 18.802, 33.905, 8.300, 9.342, 1.025, // Pt 98-104
541 1.124, 5.331, 18.078, 34.604, 8.127, 10.414, 1.322, // Au 105-111
542 2.000, 8.000, 18.000, 18.000, 14.000, 8.000, 10.000, 2.000, 2.000, // Pb 112-120
543 2.000, 8.000, 18.000, 32.000, 18.000, 8.000, 2.000, 1.000, 3.000 // U 121-129
544};
545
546//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
547
548// ShellEnergy in eV
549const G4double G4ICRU73QOModel::ShellEnergy[NQODATA] =
550 {
551 19.2, // H
552 41.8, // He
553 209.11, 21.68, // Be
554 486.2, 60.95, 23.43, // C
555 732.61, 100.646, 23.550, // N
556 965.1, 129.85, 31.60, // O
557 1525.9, 234.9, 56.18, // Ne
558 2701, 476.5, 150.42, 16.89, // Al
559 3206.1, 586.4, 186.8, 23.52, 14.91, // Si
560 5551.6, 472.43, 124.85, 22.332, // Ar
561 8554.6, 850.58, 93.47, 39.19, 19.46, // Ti
562 12254.7, 1279.29, 200.35, 49.19, 17.66, // Fe
563 14346.9, 1532.28, 262.71, 74.37, 23.03, // Ni
564 15438.5, 1667.96, 294.1, 70.69, 16.447, // Cu
565 19022.1, 2150.79, 455.79, 179.87, 57.89, 20.95, // Ge
566 24643, 2906.4, 366.85, 22.24, // Kr
567 34394, 4365.3, 589.36, 129.42, 35.59, 18.42, // Mo
568 43664.3, 5824.91, 909.79, 175.47, 54.89, 19.63, // Ag
569 49948, 6818.2, 1036.1, 172.65, 70.89, 33.87, 14.54, // Sn
570 58987, 8159, 1296.6, 356.75, 101.03, 16.52, // Xe
571 88926, 18012, 3210, 575, 108.7, 30.8, // Ta
572 115025.9, 17827.44, 3214.36, 750.41, 305.21, 105.50, 38.09, 21.25, // W
573 128342, 20254, 3601.8, 608.1, 115.0, 42.75, 17.04, // Pt
574 131872, 20903, 3757.4, 682.1, 105.2, 44.89, 17.575, // Au
575 154449, 25067, 5105.0, 987.44, 247.59, 188.1, 40.61, 19.2, 15.17, // Pb
576 167282, 27868, 6022.7, 1020.4, 244.81, 51.33, 13, 11.06, 14.43 // U
577};
578
579//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
580
581// Data for L0 from: Sigmund P., Haagerup U. Phys. Rev. A34 (1986) 892-910
582const G4double G4ICRU73QOModel::L0[67][2] =
583{
584 {0.00, 0.000001},
585 {0.10, 0.000001},
586 {0.12, 0.00001},
587 {0.14, 0.00005},
588 {0.16, 0.00014},
589 {0.18, 0.00030},
590 {0.20, 0.00057},
591 {0.25, 0.00189},
592 {0.30, 0.00429},
593 {0.35, 0.00784},
594 {0.40, 0.01248},
595 {0.45, 0.01811},
596 {0.50, 0.02462},
597 {0.60, 0.03980},
598 {0.70, 0.05731},
599 {0.80, 0.07662},
600 {0.90, 0.09733},
601 {1.00, 0.11916},
602 {1.20, 0.16532},
603 {1.40, 0.21376},
604 {1.60, 0.26362},
605 {1.80, 0.31428},
606 {2.00, 0.36532},
607 {2.50, 0.49272},
608 {3.00, 0.61765},
609 {3.50, 0.73863},
610 {4.00, 0.85496},
611 {4.50, 0.96634},
612 {5.00, 1.07272},
613 {6.00, 1.27086},
614 {7.00, 1.45075},
615 {8.00, 1.61412},
616 {9.00, 1.76277},
617 {10.00, 1.89836},
618 {12.00, 2.13625},
619 {14.00, 2.33787},
620 {16.00, 2.51093},
621 {18.00, 2.66134},
622 {20.00, 2.79358},
623 {25.00, 3.06539},
624 {30.00, 3.27902},
625 {35.00, 3.45430},
626 {40.00, 3.60281},
627 {45.00, 3.73167},
628 {50.00, 3.84555},
629 {60.00, 4.04011},
630 {70.00, 4.20264},
631 {80.00, 4.34229},
632 {90.00, 4.46474},
633 {100.00, 4.57378},
634 {120.00, 4.76155},
635 {140.00, 4.91953},
636 {160.00, 5.05590},
637 {180.00, 5.17588},
638 {200.00, 5.28299},
639 {250.00, 5.50925},
640 {300.00, 5.69364},
641 {350.00, 5.84926},
642 {400.00, 5.98388},
643 {450.00, 6.10252},
644 {500.00, 6.20856},
645 {600.00, 6.39189},
646 {700.00, 6.54677},
647 {800.00, 6.68084},
648 {900.00, 6.79905},
649 {1000.00, 6.90474}
650};
651
652//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
653
654// Data for L1 from: Mikkelsen H.H., Sigmund P. Phys. Rev. A40 (1989) 101-116
655const G4double G4ICRU73QOModel::L1[22][2] =
656{
657 {0.00, -0.000001},
658 {0.10, -0.00001},
659 {0.20, -0.00049},
660 {0.30, -0.00084},
661 {0.40, 0.00085},
662 {0.50, 0.00519},
663 {0.60, 0.01198},
664 {0.70, 0.02074},
665 {0.80, 0.03133},
666 {0.90, 0.04369},
667 {1.00, 0.06035},
668 {2.00, 0.24023},
669 {3.00, 0.44284},
670 {4.00, 0.62012},
671 {5.00, 0.77031},
672 {6.00, 0.90390},
673 {7.00, 1.02705},
674 {8.00, 1.10867},
675 {9.00, 1.17546},
676 {10.00, 1.21599},
677 {15.00, 1.24349},
678 {20.00, 1.16752}
679};
680
681//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
682
683// Data for L2 from: Mikkelsen H.H. Nucl. Instr. Meth. B58 (1991) 136-148
684const G4double G4ICRU73QOModel::L2[14][2] =
685{
686 {0.00, 0.000001},
687 {0.10, 0.00001},
688 {0.20, 0.00000},
689 {0.40, -0.00120},
690 {0.60, -0.00036},
691 {0.80, 0.00372},
692 {1.00, 0.01298},
693 {2.00, 0.08296},
694 {4.00, 0.21953},
695 {6.00, 0.23903},
696 {8.00, 0.20893},
697 {10.00, 0.10879},
698 {20.00, -0.88409},
699 {40.00, -1.13902}
700};
701
702//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
703
704// Correction obtained by V.Ivanchenko using G4BetheBlochModel
705const G4double G4ICRU73QOModel::factorBethe[99] = { 1.0,
7060.9637, 0.9872, 0.9469, 0.9875, 0.91, 0.989, 0.9507, 0.9773, 0.8621, 0.979, // 1 - 10
7070.8357, 0.868, 0.9417, 0.9466, 0.8911, 0.905, 0.944, 0.9607, 0.928, 0.96, // 11 - 20
7080.9098, 0.976, 0.8425, 0.8099, 0.7858, 0.947, 0.7248, 0.9106, 0.9246, 0.6821, // 21 - 30
7090.7223, 0.9784, 0.774, 0.7953, 0.829, 0.9405, 0.8318, 0.8583, 0.8563, 0.8481, // 31 - 40
7100.8207, 0.9033, 0.8063, 0.7837, 0.7818, 0.744, 0.875, 0.7693, 0.7871, 0.8459, // 41 - 50
7110.8231, 0.8462, 0.853, 0.8736, 0.856, 0.8762, 0.8629, 0.8323, 0.8064, 0.7828, // 51 - 60
7120.7533, 0.7273, 0.7093, 0.7157, 0.6823, 0.6612, 0.6418, 0.6395, 0.6323, 0.6221, // 61 - 70
7130.6497, 0.6746, 0.8568, 0.8541, 0.6958, 0.6962, 0.7051, 0.863, 0.8588, 0.7226, // 71 - 80
7140.7454, 0.78, 0.7783, 0.7996, 0.8216, 0.8632, 0.8558, 0.8792, 0.8745, 0.8676, // 81 - 90
7150.8321, 0.8272, 0.7999, 0.7934, 0.7787, 0.7851, 0.7692, 0.7598};
G4ElementVector
std::vector< const G4Element * > G4ElementVector
Definition: G4ElementVector.hh:44
G4Exp
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:180
G4Log
G4double G4Log(G4double x)
Definition: G4Log.hh:227
G4MaterialTable
std::vector< G4Material * > G4MaterialTable
Definition: G4MaterialTable.hh:44
kStateUndefined
@ kStateUndefined
Definition: G4Material.hh:110
G4double
double G4double
Definition: G4Types.hh:83
G4int
int G4int
Definition: G4Types.hh:85
Z
const G4int Z[17]
Definition: G4WaterStopping.cc:51
G4endl
#define G4endl
Definition: G4ios.hh:57
G4cout
G4GLOB_DLL std::ostream G4cout
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:52
CLHEP::Hep3Vector::unit
Hep3Vector unit() const
CLHEP::Hep3Vector::set
void set(double x, double y, double z)
CLHEP::Hep3Vector::rotateUz
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:33
G4AtomicShells::GetNumberOfElectrons
static G4int GetNumberOfElectrons(G4int Z, G4int SubshellNb)
Definition: G4AtomicShells.cc:793
G4AtomicShells::GetBindingEnergy
static G4double GetBindingEnergy(G4int Z, G4int SubshellNb)
Definition: G4AtomicShells.cc:779
G4AtomicShells::GetNumberOfShells
static G4int GetNumberOfShells(G4int Z)
Definition: G4AtomicShells.cc:768
G4DensityEffectData::GetPlasmaEnergy
G4double GetPlasmaEnergy(G4int idx) const
Definition: G4DensityEffectData.hh:107
G4DensityEffectData::GetElementIndex
G4int GetElementIndex(G4int Z, G4State st=kStateUndefined) const
Definition: G4DensityEffectData.cc:1293
G4DynamicParticle::GetMomentumDirection
const G4ThreeVector & GetMomentumDirection() const
G4DynamicParticle::GetKineticEnergy
G4double GetKineticEnergy() const
G4DynamicParticle::GetMomentum
G4ThreeVector GetMomentum() const
G4Electron::Electron
static G4Electron * Electron()
Definition: G4Electron.cc:93
G4Element::GetZ
G4double GetZ() const
Definition: G4Element.hh:131
G4Element::GetZasInt
G4int GetZasInt() const
Definition: G4Element.hh:132
G4ICRU73QOModel::MaxSecondaryEnergy
G4double MaxSecondaryEnergy(const G4ParticleDefinition *, G4double kinEnergy) final
Definition: G4ICRU73QOModel.cc:486
G4ICRU73QOModel::G4ICRU73QOModel
G4ICRU73QOModel(const G4ParticleDefinition *p=nullptr, const G4String &nam="ICRU73QO")
Definition: G4ICRU73QOModel.cc:64
G4ICRU73QOModel::Initialise
void Initialise(const G4ParticleDefinition *, const G4DataVector &) override
Definition: G4ICRU73QOModel.cc:97
G4ICRU73QOModel::CrossSectionPerVolume
G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy) override
Definition: G4ICRU73QOModel.cc:161
G4ICRU73QOModel::CorrectionsAlongStep
void CorrectionsAlongStep(const G4MaterialCutsCouple *, const G4DynamicParticle *, const G4double &length, G4double &eloss) override
Definition: G4ICRU73QOModel.cc:401
G4ICRU73QOModel::ComputeCrossSectionPerAtom
G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kineticEnergy, G4double Z, G4double A, G4double cutEnergy, G4double maxEnergy) override
Definition: G4ICRU73QOModel.cc:147
G4ICRU73QOModel::SampleSecondaries
void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy) override
Definition: G4ICRU73QOModel.cc:409
G4ICRU73QOModel::ComputeCrossSectionPerElectron
G4double ComputeCrossSectionPerElectron(const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
Definition: G4ICRU73QOModel.cc:121
G4ICRU73QOModel::ComputeDEDXPerVolume
G4double ComputeDEDXPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double) override
Definition: G4ICRU73QOModel.cc:176
G4MaterialCutsCouple::GetMaterial
const G4Material * GetMaterial() const
Definition: G4MaterialCutsCouple.hh:166
G4Material::GetElementVector
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:185
G4Material::GetNumberOfElements
size_t GetNumberOfElements() const
Definition: G4Material.hh:181
G4Material::GetAtomicNumDensityVector
const G4double * GetAtomicNumDensityVector() const
Definition: G4Material.hh:211
G4Material::GetElectronDensity
G4double GetElectronDensity() const
Definition: G4Material.hh:212
G4Material::GetMaterialTable
static G4MaterialTable * GetMaterialTable()
Definition: G4Material.cc:677
G4ParticleChangeForLoss::SetProposedKineticEnergy
void SetProposedKineticEnergy(G4double proposedKinEnergy)
Definition: G4ParticleChangeForLoss.hh:108
G4ParticleChangeForLoss::SetProposedMomentumDirection
void SetProposedMomentumDirection(const G4ThreeVector &dir)
Definition: G4ParticleChangeForLoss.hh:139
G4ParticleDefinition::GetParticleName
const G4String & GetParticleName() const
Definition: G4ParticleDefinition.hh:109
G4VEmAngularDistribution::SampleDirection
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
G4VEmModel::SetHighEnergyLimit
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:746
G4VEmModel::GetAngularDistribution
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:600
G4VEmModel::SelectRandomAtomNumber
G4int SelectRandomAtomNumber(const G4Material *) const
Definition: G4VEmModel.cc:253
G4VEmModel::SetDeexcitationFlag
void SetDeexcitationFlag(G4bool val)
Definition: G4VEmModel.hh:802
G4VEmModel::SetAngularDistribution
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:607
G4VEmModel::UseAngularGeneratorFlag
G4bool UseAngularGeneratorFlag() const
Definition: G4VEmModel.hh:697
G4VEmModel::MaxSecondaryKinEnergy
G4double MaxSecondaryKinEnergy(const G4DynamicParticle *dynParticle)
Definition: G4VEmModel.hh:501
G4VEmModel::GetParticleChangeForLoss
G4ParticleChangeForLoss * GetParticleChangeForLoss()
Definition: G4VEmModel.cc:109