#include <G4PhotoElectricAngularGeneratorSauterGavrila.hh>
Inheritance diagram for G4PhotoElectricAngularGeneratorSauterGavrila:
Public Member Functions | |
G4PhotoElectricAngularGeneratorSauterGavrila () | |
~G4PhotoElectricAngularGeneratorSauterGavrila () | |
virtual G4ThreeVector & | SampleDirection (const G4DynamicParticle *dp, G4double e=0.0, G4int shellId=0, const G4Material *mat=0) |
void | PrintGeneratorInformation () const |
Definition at line 57 of file G4PhotoElectricAngularGeneratorSauterGavrila.hh.
G4PhotoElectricAngularGeneratorSauterGavrila::G4PhotoElectricAngularGeneratorSauterGavrila | ( | ) |
Definition at line 54 of file G4PhotoElectricAngularGeneratorSauterGavrila.cc.
00054 : 00055 G4VEmAngularDistribution("AngularGenSauterGavrilaLowE") 00056 {}
G4PhotoElectricAngularGeneratorSauterGavrila::~G4PhotoElectricAngularGeneratorSauterGavrila | ( | ) |
void G4PhotoElectricAngularGeneratorSauterGavrila::PrintGeneratorInformation | ( | ) | const |
Definition at line 104 of file G4PhotoElectricAngularGeneratorSauterGavrila.cc.
References G4cout, and G4endl.
00105 { 00106 G4cout << "\n" << G4endl; 00107 G4cout << "" << G4endl; 00108 G4cout << "Re-implementation of the photolectric angular distribution" << G4endl; 00109 G4cout << "developed my M. Maire for the Standard EM Physics G4PhotoElectricEffect" << G4endl; 00110 G4cout << "It computes the theta distribution of the emitted electron, with respect to the" << G4endl; 00111 G4cout << "incident Gamma, using the Sauter-Gavrila distribution for the K-shell\n" << G4endl; 00112 }
G4ThreeVector & G4PhotoElectricAngularGeneratorSauterGavrila::SampleDirection | ( | const G4DynamicParticle * | dp, | |
G4double | e = 0.0 , |
|||
G4int | shellId = 0 , |
|||
const G4Material * | mat = 0 | |||
) | [virtual] |
Implements G4VEmAngularDistribution.
Definition at line 62 of file G4PhotoElectricAngularGeneratorSauterGavrila.cc.
References G4VEmAngularDistribution::fLocalDirection, G4UniformRand, G4DynamicParticle::GetKineticEnergy(), and G4DynamicParticle::GetMomentumDirection().
00065 { 00066 00067 // Compute Theta distribution of the emitted electron, with respect to the 00068 // incident Gamma. 00069 // The Sauter-Gavrila distribution for the K-shell is used. 00070 00071 G4double costeta = 1.; 00072 G4double Phi = twopi * G4UniformRand(); 00073 G4double cosphi = std::cos(Phi); 00074 G4double sinphi = std::sin(Phi); 00075 G4double sinteta = 0; 00076 G4double gamma = 1. + dp->GetKineticEnergy()/electron_mass_c2; 00077 00078 if (gamma > 5.) { 00079 fLocalDirection = dp->GetMomentumDirection(); 00080 return fLocalDirection; 00081 // Bugzilla 1120 00082 // SI on 05/09/2010 as suggested by JG 04/09/10 00083 } 00084 00085 G4double beta = std::sqrt((gamma - 1)*(gamma + 1))/gamma; 00086 G4double b = 0.5*gamma*(gamma - 1)*(gamma - 2); 00087 00088 G4double rndm,term,greject,grejsup; 00089 if (gamma < 2.) grejsup = gamma*gamma*(1.+b-beta*b); 00090 else grejsup = gamma*gamma*(1.+b+beta*b); 00091 00092 do { rndm = 1.-2*G4UniformRand(); 00093 costeta = (rndm+beta)/(rndm*beta+1.); 00094 term = 1.-beta*costeta; 00095 greject = (1.-costeta*costeta)*(1.+b*term)/(term*term); 00096 } while(greject < G4UniformRand()*grejsup); 00097 00098 sinteta = std::sqrt((1 - costeta)*(1 + costeta)); 00099 fLocalDirection.set(sinteta*cosphi, sinteta*sinphi, costeta); 00100 fLocalDirection.rotateUz(dp->GetMomentumDirection()); 00101 return fLocalDirection; 00102 }