84 for (
G4int i = 0; i < 100; ++i)
123 denEffData = (*mtab)[0]->GetIonisation()->GetDensityEffectData();
138 if(cutEnergy < maxEnergy) {
142 G4double beta2 = kineticEnergy*(kineticEnergy + 2.0*
mass)/energy2;
143 cross = (maxEnergy - cutEnergy)/(cutEnergy*maxEnergy)
144 - beta2*
G4Log(maxEnergy/cutEnergy)/tmax;
162 (p,kineticEnergy,cutEnergy,maxEnergy);
177 (p,kineticEnergy,cutEnergy,maxEnergy);
195 if (cutEnergy < tmax) {
206 if(dedx < 0.0) { dedx = 0.0; }
217 const G4double* theAtomicNumDensityVector =
225 for (
G4int i=0; i<numberOfElements; ++i)
227 const G4Element* element = (*theElementVector)[i] ;
229 * theAtomicNumDensityVector[i] *
G4int(element->
GetZ());
239 G4int Z = AtomicNumber;
240 if(Z > 97) { Z = 97; }
242 if(nbOfShells < 1) { nbOfShells = 1; }
244 G4double v = CLHEP::c_light * std::sqrt( 2.0*kineticEnergy/proton_mass_c2 );
246 G4double fBetheVelocity = CLHEP::fine_structure_const*CLHEP::c_light/v;
248 G4double tau = kineticEnergy/proton_mass_c2;
253 G4double l0Term = 0, l1Term = 0, l2Term = 0;
255 for (
G4int nos = 0; nos < nbOfShells; ++nos){
257 G4double NormalizedEnergy = (2.0*CLHEP::electron_mass_c2*beta2) /
263 l0Term += shStrength * l0;
266 l1Term += shStrength * l1;
269 l2Term += shStrength * l2;
273 (l0Term +
charge*fBetheVelocity*l1Term
274 +
chargeSquare*fBetheVelocity*fBetheVelocity*l2Term)/beta2;
282 G4int nbOfTheShell)
const
288 G4double PlasmaEnergy2 = PlasmaEnergy * PlasmaEnergy;
292 * PlasmaEnergy2 / (Z*Z) ;
296 G4double ionTerm2 = ionTerm*ionTerm ;
298 G4double oscShellEnergy = std::sqrt( ionTerm2 + plasmonTerm );
300 return oscShellEnergy;
309 for(n = 0; n <
sizeL0; n++) {
310 if( normEnergy <
L0[n][0] )
break;
312 if(0 == n) { n = 1; }
313 if(n >= sizeL0) { n = sizeL0 - 1; }
317 G4double bethe = l0p + (l0 - l0p) * ( normEnergy -
L0[n-1][0]) /
318 (
L0[
n][0] -
L0[n-1][0]);
329 for(n = 0; n <
sizeL1; n++) {
330 if( normEnergy <
L1[n][0] )
break;
333 if(n >= sizeL1) n = sizeL1 - 1 ;
337 G4double barkas= l1p + (l1 - l1p) * ( normEnergy -
L1[n-1][0]) /
338 (
L1[
n][0] -
L1[n-1][0]);
348 for(n = 0; n <
sizeL2; n++) {
349 if( normEnergy <
L2[n][0] )
break;
352 if(n >= sizeL2) n = sizeL2 - 1 ;
356 G4double bloch = l2p + (l2 - l2p) * ( normEnergy -
L2[n-1][0]) /
357 (
L2[
n][0] -
L2[n-1][0]);
381 if(xmin >= xmax) {
return; }
386 G4double beta2 = kineticEnergy*(kineticEnergy + 2.0*
mass)/energy2;
395 deltaKinEnergy = xmin*xmax/(xmin*(1.0 -
x) + xmax*x);
397 f = 1.0 - beta2*deltaKinEnergy/tmax;
400 G4cout <<
"G4ICRU73QOModel::SampleSecondary Warning! "
401 <<
"Majorant " << grej <<
" < "
402 << f <<
" for e= " << deltaKinEnergy
421 sqrt(deltaKinEnergy * (deltaKinEnergy + 2.0*electron_mass_c2));
422 G4double totMomentum = energy*sqrt(beta2);
423 G4double cost = deltaKinEnergy * (energy + electron_mass_c2) /
424 (deltaMomentum * totMomentum);
425 if(cost > 1.0) { cost = 1.0; }
426 G4double sint = sqrt((1.0 - cost)*(1.0 + cost));
430 deltaDirection.set(sint*cos(phi),sint*sin(phi), cost) ;
431 deltaDirection.rotateUz(direction);
438 kineticEnergy -= deltaKinEnergy;
440 finalP = finalP.unit();
445 vdp->push_back(delta);
455 G4double tmax = 2.0*electron_mass_c2*tau*(tau + 2.) /
462 const G4int G4ICRU73QOModel::ZElementAvailable[NQOELEM] = {1,2,4,6,7,8,10,13,14,-18,
463 22,26,28,29,32,36,42,47,
464 50,54,73,74,78,79,82,92};
466 const G4int G4ICRU73QOModel::nbofShellsForElement[NQOELEM] = {1,1,2,3,3,3,3,4,5,4,
470 const G4int G4ICRU73QOModel::startElemIndex[NQOELEM] = {0,1,2,4,7,10,13,16,20,25,
471 29,34,39,44,49,55,59,65,
472 71,78,84,90,98,105,112,121};
486 1.623, 2.147, 6.259, 2.971,
487 1.631, 2.094, 6.588, 2.041, 1.646,
488 1.535, 8.655, 1.706, 6.104,
489 1.581, 8.358, 8.183, 2.000, 1.878,
490 1.516, 8.325, 8.461, 6.579, 1.119,
491 1.422, 7.81, 8.385, 8.216, 2.167,
492 1.458, 8.049, 8.79, 9.695, 1.008,
493 1.442, 7.791, 7.837, 10.122, 2.463, 2.345,
494 1.645, 7.765, 19.192, 7.398,
495 1.313, 6.409, 19.229, 8.633, 5.036, 1.380,
496 1.295, 6.219, 18.751, 8.748, 10.184, 1.803,
497 1.277, 6.099, 20.386, 8.011, 10.007, 2.272, 1.948,
498 1.563, 6.312, 21.868, 5.762, 11.245, 7.250,
499 0.9198, 6.5408, 18.9727, 24.9149, 15.0161, 6.6284,
500 1.202, 5.582, 19.527, 18.741, 8.411, 14.387, 4.042, 2.108,
501 1.159, 5.467, 18.802, 33.905, 8.300, 9.342, 1.025,
502 1.124, 5.331, 18.078, 34.604, 8.127, 10.414, 1.322,
503 2.000, 8.000, 18.000, 18.000, 14.000, 8.000, 10.000, 2.000, 2.000,
504 2.000, 8.000, 18.000, 32.000, 18.000, 8.000, 2.000, 1.000, 3.000
516 732.61, 100.646, 23.550,
517 965.1, 129.85, 31.60,
518 1525.9, 234.9, 56.18,
519 2701, 476.5, 150.42, 16.89,
520 3206.1, 586.4, 186.8, 23.52, 14.91,
521 5551.6, 472.43, 124.85, 22.332,
522 8554.6, 850.58, 93.47, 39.19, 19.46,
523 12254.7, 1279.29, 200.35, 49.19, 17.66,
524 14346.9, 1532.28, 262.71, 74.37, 23.03,
525 15438.5, 1667.96, 294.1, 70.69, 16.447,
526 19022.1, 2150.79, 455.79, 179.87, 57.89, 20.95,
527 24643, 2906.4, 366.85, 22.24,
528 34394, 4365.3, 589.36, 129.42, 35.59, 18.42,
529 43664.3, 5824.91, 909.79, 175.47, 54.89, 19.63,
530 49948, 6818.2, 1036.1, 172.65, 70.89, 33.87, 14.54,
531 58987, 8159, 1296.6, 356.75, 101.03, 16.52,
532 88926, 18012, 3210, 575, 108.7, 30.8,
533 115025.9, 17827.44, 3214.36, 750.41, 305.21, 105.50, 38.09, 21.25,
534 128342, 20254, 3601.8, 608.1, 115.0, 42.75, 17.04,
535 131872, 20903, 3757.4, 682.1, 105.2, 44.89, 17.575,
536 154449, 25067, 5105.0, 987.44, 247.59, 188.1, 40.61, 19.2, 15.17,
537 167282, 27868, 6022.7, 1020.4, 244.81, 51.33, 13, 11.06, 14.43
667 0.9637, 0.9872, 0.9469, 0.9875, 0.91, 0.989, 0.9507, 0.9773, 0.8621, 0.979,
668 0.8357, 0.868, 0.9417, 0.9466, 0.8911, 0.905, 0.944, 0.9607, 0.928, 0.96,
669 0.9098, 0.976, 0.8425, 0.8099, 0.7858, 0.947, 0.7248, 0.9106, 0.9246, 0.6821,
670 0.7223, 0.9784, 0.774, 0.7953, 0.829, 0.9405, 0.8318, 0.8583, 0.8563, 0.8481,
671 0.8207, 0.9033, 0.8063, 0.7837, 0.7818, 0.744, 0.875, 0.7693, 0.7871, 0.8459,
672 0.8231, 0.8462, 0.853, 0.8736, 0.856, 0.8762, 0.8629, 0.8323, 0.8064, 0.7828,
673 0.7533, 0.7273, 0.7093, 0.7157, 0.6823, 0.6612, 0.6418, 0.6395, 0.6323, 0.6221,
674 0.6497, 0.6746, 0.8568, 0.8541, 0.6958, 0.6962, 0.7051, 0.863, 0.8588, 0.7226,
675 0.7454, 0.78, 0.7783, 0.7996, 0.8216, 0.8632, 0.8558, 0.8792, 0.8745, 0.8676,
676 0.8321, 0.8272, 0.7999, 0.7934, 0.7787, 0.7851, 0.7692, 0.7598};
virtual G4double ComputeDEDXPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double)
G4double MaxSecondaryKinEnergy(const G4DynamicParticle *dynParticle)
static const G4int nbofShellsForElement[NQOELEM]
static const G4double factorBethe[99]
G4double DEDXPerElement(G4int Z, G4double kineticEnergy)
G4double GetL1(G4double normEnergy) const
G4ParticleChangeForLoss * GetParticleChangeForLoss()
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
std::vector< G4Element * > G4ElementVector
static const G4double SubShellOccupation[NQODATA]
G4double DEDX(const G4Material *material, G4double kineticEnergy)
G4double GetKineticEnergy() const
CLHEP::Hep3Vector G4ThreeVector
static const G4double L0[67][2]
G4int GetElementIndex(G4int Z, G4State mState)
static const G4double L2[14][2]
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
static G4MaterialTable * GetMaterialTable()
std::vector< G4Material * > G4MaterialTable
G4VEmAngularDistribution * GetAngularDistribution()
G4double GetL2(G4double normEnergy) const
static const G4double L1[22][2]
G4ParticleChangeForLoss * fParticleChange
const G4ElementVector * GetElementVector() const
G4ICRU73QOModel(const G4ParticleDefinition *p=0, const G4String &nam="ICRU73QO")
const G4String & GetParticleName() const
G4double GetShellStrength(G4int Z, G4int nbOfTheShell) const
void SetHighEnergyLimit(G4double)
G4ParticleDefinition * theElectron
G4double GetShellEnergy(G4int Z, G4int nbOfTheShell) const
G4GLOB_DLL std::ostream G4cout
G4DensityEffectData * denEffData
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
virtual G4double ComputeCrossSectionPerElectron(const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
G4double GetElectronDensity() const
G4bool UseAngularGeneratorFlag() const
static G4double GetBindingEnergy(G4int Z, G4int SubshellNb)
const G4ThreeVector & GetMomentumDirection() const
static const double twopi
void SetProposedKineticEnergy(G4double proposedKinEnergy)
G4double GetPlasmaEnergy(G4int idx)
void SetParticle(const G4ParticleDefinition *p)
virtual ~G4ICRU73QOModel()
void SetProposedMomentumDirection(const G4ThreeVector &dir)
G4double GetOscillatorEnergy(G4int Z, G4int nbOfTheShell) const
const G4double * GetAtomicNumDensityVector() const
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)
G4double G4Log(G4double x)
static const G4int NQOELEM
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
G4double energy(const ThreeVector &p, const G4double m)
void SetAngularDistribution(G4VEmAngularDistribution *)
const G4double x[NPOINTSGL]
static G4int GetNumberOfElectrons(G4int Z, G4int SubshellNb)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
static G4Electron * Electron()
size_t GetNumberOfElements() const
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kineticEnergy, G4double Z, G4double A, G4double cutEnergy, G4double maxEnergy)
virtual G4double MaxSecondaryEnergy(const G4ParticleDefinition *, G4double kinEnergy)
void SetDeexcitationFlag(G4bool val)
static const G4double ShellEnergy[NQODATA]
G4int GetNumberOfShells(G4int Z) const
static const G4int ZElementAvailable[NQOELEM]
const G4ParticleDefinition * particle
static const G4int startElemIndex[NQOELEM]
G4double GetL0(G4double normEnergy) const
virtual void CorrectionsAlongStep(const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double &eloss, G4double &niel, G4double length)
G4ThreeVector GetMomentum() const
const G4Material * GetMaterial() const
G4int SelectRandomAtomNumber(const G4Material *)