69 0.5917, 0.7628, 0.8983, 0.9801 };
71 0.1813, 0.1569, 0.1112, 0.0506 };
73 static const G4double tlow = 1.*CLHEP::MeV;
79 ah10 = 4.67733E+00, ah11 =-6.19012E-01, ah12 = 2.02225E-02,
80 ah20 =-7.34101E+00, ah21 = 1.00462E+00, ah22 =-3.20985E-02,
81 ah30 = 2.93119E+00, ah31 =-4.03761E-01, ah32 = 1.25153E-02;
84 bh10 = 4.23071E+00, bh11 =-6.10995E-01, bh12 = 1.95531E-02,
85 bh20 =-7.12527E+00, bh21 = 9.69160E-01, bh22 =-2.74255E-02,
86 bh30 = 2.69925E+00, bh31 =-3.63283E-01, bh32 = 9.55316E-03;
89 al00 =-2.05398E+00, al01 = 2.38815E-02, al02 = 5.25483E-04,
90 al10 =-7.69748E-02, al11 =-6.91499E-02, al12 = 2.22453E-03,
91 al20 = 4.06463E-02, al21 =-1.01281E-02, al22 = 3.40919E-04;
94 bl00 = 1.04133E+00, bl01 =-9.43291E-03, bl02 =-4.54758E-04,
95 bl10 = 1.19253E-01, bl11 = 4.07467E-02, bl12 =-1.30718E-03,
96 bl20 =-1.59391E-02, bl21 = 7.27752E-03, bl22 =-1.94405E-04;
112 lowKinEnergy = 10.*
MeV;
122 InitialiseConstants();
123 if(p) { SetParticle(p); }
128 void G4eBremParametrizedModel::InitialiseConstants()
178 if(p) { SetParticle(p); }
186 if(isInitialised) {
return; }
188 isInitialised =
true;
208 if(kineticEnergy < lowKinEnergy) {
return 0.0; }
210 if(cut == 0.0) {
return 0.0; }
223 SetCurrentElement((*theElementVector)[i]->GetZ());
225 dedx += theAtomicNumDensityVector[i]*
currentZ*
currentZ*ComputeBremLoss(cut);
247 for(
G4int l=0; l<
n; l++) {
249 for(
G4int i=0; i<8; i++) {
253 xs = ComputeDXSectionPerAtom(eg);
275 if(kineticEnergy < lowKinEnergy) {
return 0.0; }
279 if(cut >= tmax) {
return 0.0; }
281 SetCurrentElement(Z);
283 G4double cross = ComputeXSectionPerAtom(cut);
286 if(tmax <
kinEnergy) { cross -= ComputeXSectionPerAtom(tmax); }
288 cross *= Z*Z*bremFactor;
310 for(
G4int l=0; l<
n; l++) {
312 for(
G4int i=0; i<8; i++) {
316 xs = ComputeDXSectionPerAtom(eg);
336 if (ScreenVariable > 1.)
337 screenVal = 42.24 - 8.368*
G4Log(ScreenVariable+0.952);
339 screenVal = 42.392 - ScreenVariable* (7.796 - 1.961*ScreenVariable);
352 if (ScreenVariable > 1.)
353 screenVal = 42.24 - 8.368*
G4Log(ScreenVariable+0.952);
355 screenVal = 41.734 - ScreenVariable* (6.484 - 1.250*ScreenVariable);
363 G4double G4eBremParametrizedModel::ComputeParametrizedDXSectionPerAtom(
367 SetCurrentElement(Z);
382 if (kineticEnergy > tlow) {
384 G4double ah1 = ah10 + ZZ* (ah11 + ZZ* ah12);
385 G4double ah2 = ah20 + ZZ* (ah21 + ZZ* ah22);
386 G4double ah3 = ah30 + ZZ* (ah31 + ZZ* ah32);
388 G4double bh1 = bh10 + ZZ* (bh11 + ZZ* bh12);
389 G4double bh2 = bh20 + ZZ* (bh21 + ZZ* bh22);
390 G4double bh3 = bh30 + ZZ* (bh31 + ZZ* bh32);
392 ah = 1. + (ah1*U2 + ah2*U + ah3) / (U2*U);
393 bh = 0.75 + (bh1*U2 + bh2*U + bh3) / (U2*U);
400 G4double screenvar = screenfac*epsil/(1.0-epsil);
401 G4double F1 =
max(ScreenFunction1(screenvar) - FZ ,0.);
402 G4double F2 =
max(ScreenFunction2(screenvar) - FZ ,0.);
405 greject = (F1 - epsil* (ah*F1 - bh*epsil*F2))/8.;
407 std::cout <<
" yy = "<<epsil<<std::endl;
408 std::cout <<
" F1/(...) "<<F1/(42.392 - FZ)<<std::endl;
409 std::cout <<
" F2/(...) "<<F2/(42.392 - FZ)<<std::endl;
410 std::cout <<
" (42.392 - FZ) " << (42.392 - FZ) <<std::endl;
414 G4double al0 = al00 + ZZ* (al01 + ZZ* al02);
415 G4double al1 = al10 + ZZ* (al11 + ZZ* al12);
416 G4double al2 = al20 + ZZ* (al21 + ZZ* al22);
418 G4double bl0 = bl00 + ZZ* (bl01 + ZZ* bl02);
419 G4double bl1 = bl10 + ZZ* (bl11 + ZZ* bl12);
420 G4double bl2 = bl20 + ZZ* (bl21 + ZZ* bl22);
422 ah = al0 + al1*U + al2*U2;
423 bh = bl0 + bl1*U + bl2*U2;
426 greject=(1. + x* (ah + bh*
x));
441 G4double G4eBremParametrizedModel::ComputeDXSectionPerAtom(
G4double gammaEnergy)
444 if(gammaEnergy < 0.0) {
return 0.0; }
456 std::cout<<
" F1(0) "<<ScreenFunction1(0.) <<std::endl;
457 std::cout<<
" F1(0) "<<ScreenFunction2(0.) <<std::endl;
458 std::cout<<
"Ekin = "<<
kinEnergy<<std::endl;
459 std::cout<<
"Z = "<<
currentZ<<std::endl;
460 std::cout<<
"main = "<<main<<std::endl;
461 std::cout<<
" y = "<<y<<std::endl;
462 std::cout<<
" Fel-fCoulomb "<< (
Fel-
fCoulomb) <<std::endl;
465 std::cout<<
"main2 = "<<main2<<std::endl;
475 std::vector<G4DynamicParticle*>* vdp,
482 if(kineticEnergy < lowKinEnergy) {
return; }
485 if(cut >= emax) {
return; }
491 SetCurrentElement(elm->GetZ());
504 gammaEnergy = sqrt(x);
505 f = ComputeDXSectionPerAtom(gammaEnergy);
508 G4cout <<
"### G4eBremParametrizedModel Warning: Majoranta exceeded! "
509 << f <<
" > " <<
fMax
510 <<
" Egamma(MeV)= " << gammaEnergy
511 <<
" E(mEV)= " << kineticEnergy
529 vdp->push_back(gamma);
531 G4double totMomentum = sqrt(kineticEnergy*(totalEnergy + electron_mass_c2));
533 - gammaEnergy*gammaDirection).unit();
536 G4double finalE = kineticEnergy - gammaEnergy;
virtual void InitialiseLocal(const G4ParticleDefinition *, G4VEmModel *masterModel)
G4double LowEnergyLimit() const
G4ParticleChangeForLoss * GetParticleChangeForLoss()
G4double SecondaryThreshold() const
G4bool isElectron(G4int ityp)
std::vector< G4Element * > G4ElementVector
G4double GetKineticEnergy() const
void InitialiseElementSelectors(const G4ParticleDefinition *, const G4DataVector &)
G4ParticleChangeForLoss * fParticleChange
G4double GetZ13(G4double Z)
virtual G4double ComputeDEDXPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy)
G4VEmAngularDistribution * GetAngularDistribution()
G4eBremParametrizedModel(const G4ParticleDefinition *p=0, const G4String &nam="eBremParam")
int main(int argc, char **argv)
static const G4double xgi[8]
const G4ElementVector * GetElementVector() const
static G4NistManager * Instance()
G4GLOB_DLL std::ostream G4cout
G4ParticleDefinition * theGamma
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
G4double GetElectronDensity() const
const G4ThreeVector & GetMomentumDirection() const
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double cutEnergy, G4double maxEnergy)
void SetProposedKineticEnergy(G4double proposedKinEnergy)
std::vector< G4EmElementSelector * > * GetElementSelectors()
void SetProposedMomentumDirection(const G4ThreeVector &dir)
const G4double * GetAtomicNumDensityVector() const
virtual ~G4eBremParametrizedModel()
G4double G4Log(G4double x)
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
void SetElementSelectors(std::vector< G4EmElementSelector * > *)
G4double GetPDGMass() const
T max(const T t1, const T t2)
brief Return the largest of the two arguments
void SetAngularDistribution(G4VEmAngularDistribution *)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
const G4ParticleDefinition * particle
static G4Electron * Electron()
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double tkin, G4double Z, G4double, G4double cutEnergy, G4double maxEnergy=DBL_MAX)
size_t GetNumberOfElements() const
virtual void SetupForMaterial(const G4ParticleDefinition *, const G4Material *, G4double)
void ProposeTrackStatus(G4TrackStatus status)
void SetLowEnergyLimit(G4double)
virtual G4double MinEnergyCut(const G4ParticleDefinition *, const G4MaterialCutsCouple *)
static const G4double wgi[8]
void SetCurrentElement(const G4Element *)
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
const G4Material * GetMaterial() const
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)