00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044 #include "G4SauterGavrilaAngularDistribution.hh"
00045 #include "G4PhysicalConstants.hh"
00046 #include "Randomize.hh"
00047
00048 G4SauterGavrilaAngularDistribution::G4SauterGavrilaAngularDistribution()
00049 : G4VEmAngularDistribution("AngularGenSauterGavrila")
00050 {}
00051
00052 G4SauterGavrilaAngularDistribution::~G4SauterGavrilaAngularDistribution()
00053 {}
00054
00055 G4ThreeVector&
00056 G4SauterGavrilaAngularDistribution::SampleDirection(
00057 const G4DynamicParticle* dp, G4double, G4int, const G4Material*)
00058 {
00059 G4double tau = dp->GetKineticEnergy()/electron_mass_c2;
00060 const G4double taulimit = 30.0;
00061
00062 if (tau > taulimit) {
00063 fLocalDirection = dp->GetMomentumDirection();
00064
00065
00066 } else {
00067
00068 G4double invgamma = 1.0/(tau + 1.0);
00069 G4double beta = std::sqrt(tau*(tau + 2.0))*invgamma;
00070 G4double b = 0.5*tau*(tau*tau - 1.0);
00071 G4double invgamma2 = invgamma*invgamma;
00072
00073 G4double rndm,term,greject,grejsup,costeta,sint2;
00074 if (tau < 1.) { grejsup = (1.+b-beta*b)/invgamma2; }
00075 else { grejsup = (1.+b+beta*b)/invgamma2; }
00076
00077 do {
00078 rndm = 1 - 2*G4UniformRand();
00079 costeta = (rndm + beta)/(rndm*beta + 1);
00080 term = invgamma2/(1 + beta*rndm);
00081 sint2 = (1 - costeta)*(1 + costeta);
00082 greject = sint2*(1 + b*term)/(term*term);
00083
00084 } while(greject < G4UniformRand()*grejsup);
00085
00086 G4double sinteta = std::sqrt(sint2);
00087 G4double phi = CLHEP::twopi*G4UniformRand();
00088
00089 fLocalDirection.set(sinteta*std::cos(phi), sinteta*std::sin(phi), costeta);
00090 fLocalDirection.rotateUz(dp->GetMomentumDirection());
00091 }
00092 return fLocalDirection;
00093 }
00094
00095 void G4SauterGavrilaAngularDistribution::PrintGeneratorInformation() const
00096 {
00097 G4cout << "\n" << G4endl;
00098 G4cout << "Non-polarized photoelectric effect angular generator." << G4endl;
00099 G4cout << "The Sauter-Gavrila distribution for the K-shell is used."<<G4endl;
00100 G4cout << "Originally developed by M.Maire for Geant3"
00101 << G4endl;
00102 }