#include <G4DipBustGenerator.hh>
Inheritance diagram for G4DipBustGenerator:
Public Member Functions | |
G4DipBustGenerator (const G4String &name="") | |
virtual | ~G4DipBustGenerator () |
virtual G4ThreeVector & | SampleDirection (const G4DynamicParticle *dp, G4double out_energy, G4int Z, const G4Material *mat=0) |
G4double | PolarAngle (const G4double initial_energy, const G4double final_energy, const G4int Z) |
void | PrintGeneratorInformation () const |
Definition at line 55 of file G4DipBustGenerator.hh.
G4DipBustGenerator::G4DipBustGenerator | ( | const G4String & | name = "" |
) |
Definition at line 57 of file G4DipBustGenerator.cc.
00058 : G4VEmAngularDistribution("DipBustGen") 00059 {}
G4DipBustGenerator::~G4DipBustGenerator | ( | ) | [virtual] |
G4double G4DipBustGenerator::PolarAngle | ( | const G4double | initial_energy, | |
const G4double | final_energy, | |||
const G4int | Z | |||
) |
Definition at line 103 of file G4DipBustGenerator.cc.
References G4UniformRand, and G4INCL::Math::pi.
00106 { 00107 G4double c, cosTheta, delta, cofA, signc = 1., a, power = 1./3.; 00108 G4double gamma, beta, theta; 00109 00110 c = 4. - 8.*G4UniformRand(); 00111 a = c; 00112 00113 if( c < 0. ) 00114 { 00115 signc = -1.; 00116 a = -c; 00117 } 00118 delta = std::sqrt(a*a+4.); 00119 delta += a; 00120 delta *= 0.5; 00121 00122 cofA = -signc*std::pow(delta, power); 00123 00124 cosTheta = cofA - 1./cofA; 00125 00126 gamma = 1. + eTkin/electron_mass_c2; 00127 beta = std::sqrt(1. - 1./gamma/gamma); 00128 00129 cosTheta = (cosTheta + beta)/(1 + cosTheta*beta); 00130 00131 theta = std::acos(cosTheta); 00132 00133 if( theta < 0. ) theta = 0.; 00134 if( theta > pi ) theta = pi; 00135 // G4cout <<"theta = "<<theta<<"; "; 00136 00137 return theta; 00138 }
void G4DipBustGenerator::PrintGeneratorInformation | ( | ) | const |
Definition at line 140 of file G4DipBustGenerator.cc.
References G4cout, and G4endl.
00141 { 00142 G4cout << "\n" << G4endl; 00143 G4cout << "Angular Generator based on classical formula from" << G4endl; 00144 G4cout << "J.D. Jackson, Classical Electrodynamics, Wiley, New York 1975" 00145 << G4endl; 00146 }
G4ThreeVector & G4DipBustGenerator::SampleDirection | ( | const G4DynamicParticle * | dp, | |
G4double | out_energy, | |||
G4int | Z, | |||
const G4Material * | mat = 0 | |||
) | [virtual] |
Implements G4VEmAngularDistribution.
Definition at line 65 of file G4DipBustGenerator.cc.
References G4VEmAngularDistribution::fLocalDirection, G4UniformRand, G4DynamicParticle::GetKineticEnergy(), and G4DynamicParticle::GetMomentumDirection().
00067 { 00068 G4double c, cosTheta, delta, cofA, signc = 1., a, power = 1./3.; 00069 00070 G4double eTkin = dp->GetKineticEnergy(); 00071 00072 c = 4. - 8.*G4UniformRand(); 00073 a = c; 00074 00075 if( c < 0. ) 00076 { 00077 signc = -1.; 00078 a = -c; 00079 } 00080 delta = std::sqrt(a*a+4.); 00081 delta += a; 00082 delta *= 0.5; 00083 00084 cofA = -signc*std::pow(delta, power); 00085 00086 cosTheta = cofA - 1./cofA; 00087 00088 G4double tau = eTkin/electron_mass_c2; 00089 G4double beta = std::sqrt(tau*(tau + 2.))/(tau + 1.); 00090 00091 cosTheta = (cosTheta + beta)/(1 + cosTheta*beta); 00092 00093 G4double sinTheta = std::sqrt((1 - cosTheta)*(1 + cosTheta)); 00094 G4double phi = twopi*G4UniformRand(); 00095 00096 fLocalDirection.set(sinTheta*std::cos(phi), sinTheta*std::sin(phi),cosTheta); 00097 fLocalDirection.rotateUz(dp->GetMomentumDirection()); 00098 00099 return fLocalDirection; 00100 00101 }