#include <G4SauterGavrilaAngularDistribution.hh>
Inheritance diagram for G4SauterGavrilaAngularDistribution:
Public Member Functions | |
G4SauterGavrilaAngularDistribution () | |
virtual | ~G4SauterGavrilaAngularDistribution () |
virtual G4ThreeVector & | SampleDirection (const G4DynamicParticle *dp, G4double e=0.0, G4int shellId=0, const G4Material *mat=0) |
void | PrintGeneratorInformation () const |
Definition at line 55 of file G4SauterGavrilaAngularDistribution.hh.
G4SauterGavrilaAngularDistribution::G4SauterGavrilaAngularDistribution | ( | ) |
Definition at line 48 of file G4SauterGavrilaAngularDistribution.cc.
00049 : G4VEmAngularDistribution("AngularGenSauterGavrila") 00050 {}
G4SauterGavrilaAngularDistribution::~G4SauterGavrilaAngularDistribution | ( | ) | [virtual] |
void G4SauterGavrilaAngularDistribution::PrintGeneratorInformation | ( | ) | const |
Definition at line 95 of file G4SauterGavrilaAngularDistribution.cc.
References G4cout, and G4endl.
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 }
G4ThreeVector & G4SauterGavrilaAngularDistribution::SampleDirection | ( | const G4DynamicParticle * | dp, | |
G4double | e = 0.0 , |
|||
G4int | shellId = 0 , |
|||
const G4Material * | mat = 0 | |||
) | [virtual] |
Implements G4VEmAngularDistribution.
Definition at line 56 of file G4SauterGavrilaAngularDistribution.cc.
References G4VEmAngularDistribution::fLocalDirection, G4UniformRand, G4DynamicParticle::GetKineticEnergy(), and G4DynamicParticle::GetMomentumDirection().
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 // Bugzilla 1120 00065 // SI on 05/09/2010 as suggested by JG 04/09/10 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 }