G4eBremParametrizedModel.cc File Reference

#include "G4eBremParametrizedModel.hh"
#include "G4PhysicalConstants.hh"
#include "G4SystemOfUnits.hh"
#include "G4Electron.hh"
#include "G4Positron.hh"
#include "G4Gamma.hh"
#include "Randomize.hh"
#include "G4Material.hh"
#include "G4Element.hh"
#include "G4ElementVector.hh"
#include "G4ProductionCutsTable.hh"
#include "G4ParticleChangeForLoss.hh"
#include "G4LossTableManager.hh"
#include "G4ModifiedTsai.hh"

Go to the source code of this file.

Functions

G4double ScreenFunction1 (G4double ScreenVariable)
G4double ScreenFunction2 (G4double ScreenVariable)
G4double ComputeParametrizedDXSectionPerAtom (G4double kineticEnergy, G4double gammaEnergy, G4double Z)


Function Documentation

G4double ComputeParametrizedDXSectionPerAtom ( G4double  kineticEnergy,
G4double  gammaEnergy,
G4double  Z 
)

Definition at line 357 of file G4eBremParametrizedModel.cc.

References ScreenFunction1(), and ScreenFunction2().

00357                                                                                                        {
00358   G4double lnZ = std::log(Z); // 3.*(anElement->GetIonisation()->GetlogZ3());
00359   G4double FZ = lnZ* (4.- 0.55*lnZ);
00360   G4double ZZ = std::pow (Z*(Z+1.),1./3.); // anElement->GetIonisation()->GetZZ3();
00361   G4double Z3 = std::pow (Z,1./3.); // (anElement->GetIonisation()->GetZ3())
00362 
00363   G4double totalEnergy = kineticEnergy + electron_mass_c2;
00364 
00365   //  G4double x, epsil, greject, migdal, grejmax, q;
00366   G4double epsil, greject;
00367   G4double U  = log(kineticEnergy/electron_mass_c2);
00368   G4double U2 = U*U;
00369 
00370   // precalculated parameters
00371   G4double ah, bh;
00372 
00373 if (kineticEnergy > tlow) {
00374        
00375     G4double ah1 = ah10 + ZZ* (ah11 + ZZ* ah12);
00376     G4double ah2 = ah20 + ZZ* (ah21 + ZZ* ah22);
00377     G4double ah3 = ah30 + ZZ* (ah31 + ZZ* ah32);
00378 
00379     G4double bh1 = bh10 + ZZ* (bh11 + ZZ* bh12);
00380     G4double bh2 = bh20 + ZZ* (bh21 + ZZ* bh22);
00381     G4double bh3 = bh30 + ZZ* (bh31 + ZZ* bh32);
00382 
00383     ah = 1.   + (ah1*U2 + ah2*U + ah3) / (U2*U);
00384     bh = 0.75 + (bh1*U2 + bh2*U + bh3) / (U2*U);
00385 
00386     // limit of the screening variable
00387     G4double screenfac =
00388       136.*electron_mass_c2/(Z3*totalEnergy);
00389 
00390     epsil = gammaEnergy/totalEnergy; //         epsil = x*kineticEnergy/totalEnergy;
00391         G4double screenvar = screenfac*epsil/(1.0-epsil);
00392         G4double F1 = max(ScreenFunction1(screenvar) - FZ ,0.);
00393         G4double F2 = max(ScreenFunction2(screenvar) - FZ ,0.);
00394 
00395 
00396         greject = (F1 - epsil* (ah*F1 - bh*epsil*F2))/8.; //  1./(42.392 - FZ);
00397 
00398     std::cout << " yy = "<<epsil<<std::endl;
00399     std::cout << " F1/(...) "<<F1/(42.392 - FZ)<<std::endl;
00400     std::cout << " F2/(...) "<<F2/(42.392 - FZ)<<std::endl;
00401     std::cout << " (42.392 - FZ) " << (42.392 - FZ) <<std::endl;
00402 
00403   } else {  
00404 
00405     G4double al0 = al00 + ZZ* (al01 + ZZ* al02);
00406     G4double al1 = al10 + ZZ* (al11 + ZZ* al12);
00407     G4double al2 = al20 + ZZ* (al21 + ZZ* al22);
00408  
00409     G4double bl0 = bl00 + ZZ* (bl01 + ZZ* bl02);
00410     G4double bl1 = bl10 + ZZ* (bl11 + ZZ* bl12);
00411     G4double bl2 = bl20 + ZZ* (bl21 + ZZ* bl22);
00412  
00413     ah = al0 + al1*U + al2*U2;
00414     bh = bl0 + bl1*U + bl2*U2;
00415 
00416     G4double x=gammaEnergy/kineticEnergy;
00417     greject=(1. + x* (ah + bh*x));
00418 
00419     /*
00420     // Compute the maximum of the rejection function
00421     grejmax = max(1. + xmin* (ah + bh*xmin), 1.+ah+bh);
00422     G4double xm = -ah/(2.*bh);
00423     if ( xmin < xm && xm < xmax) grejmax = max(grejmax, 1.+ xm* (ah + bh*xm));
00424     */
00425   }
00426 
00427  return greject;
00428 }

G4double ScreenFunction1 ( G4double  ScreenVariable  ) 

Definition at line 323 of file G4eBremParametrizedModel.cc.

Referenced by ComputeParametrizedDXSectionPerAtom().

00327 {
00328   G4double screenVal;
00329 
00330   if (ScreenVariable > 1.)
00331     screenVal = 42.24 - 8.368*std::log(ScreenVariable+0.952);
00332   else
00333     screenVal = 42.392 - ScreenVariable* (7.796 - 1.961*ScreenVariable);
00334 
00335   return screenVal;
00336 } 

G4double ScreenFunction2 ( G4double  ScreenVariable  ) 

Definition at line 340 of file G4eBremParametrizedModel.cc.

Referenced by ComputeParametrizedDXSectionPerAtom().

00344 {
00345   G4double screenVal;
00346 
00347   if (ScreenVariable > 1.)
00348     screenVal = 42.24 - 8.368*std::log(ScreenVariable+0.952);
00349   else
00350     screenVal = 41.734 - ScreenVariable* (6.484 - 1.250*ScreenVariable);
00351 
00352   return screenVal;
00353 } 


Generated on Mon May 27 17:50:41 2013 for Geant4 by  doxygen 1.4.7