90 pi_rcl2(
pi*classic_electr_radius*classic_electr_radius),
124 G4double tau = ekin/electron_mass_c2;
174 CLHEP::HepRandomEngine* rndmEngine = G4Random::getTheEngine();
177 if(PositKinEnergy == 0.0) {
178 G4double cost = 2.*rndmEngine->flat()-1.;
179 G4double sint = sqrt((1. - cost)*(1. + cost));
182 phi =
twopi * rndmEngine->flat();
194 G4double tau = PositKinEnergy/electron_mass_c2;
197 G4double sqgrate = sqrt(tau/tau2)*0.5;
203 G4double epsilqot = epsilmax/epsilmin;
211 epsil = epsilmin*
G4Exp(
G4Log(epsilqot)*rndmEngine->flat());
212 greject = 1. - epsil + (2.*gam*epsil-1.)/(epsil*tau2*tau2);
214 }
while( greject < rndmEngine->
flat());
220 G4double cost = (epsil*tau2-1.)/(epsil*sqg2m1);
221 if(std::abs(cost) > 1.0) {
222 G4cout <<
"### G4eeToTwoGammaModel WARNING cost= " << cost
223 <<
" positron Ekin(MeV)= " << PositKinEnergy
224 <<
" gamma epsil= " << epsil
226 if(cost > 1.0) cost = 1.0;
229 G4double sint = sqrt((1.+cost)*(1.-cost));
236 G4double TotalAvailableEnergy = PositKinEnergy + 2.0*electron_mass_c2;
237 G4double Phot1Energy = epsil*TotalAvailableEnergy;
239 G4ThreeVector Phot1Direction(sint*cos(phi), sint*sin(phi), cost);
240 Phot1Direction.rotateUz(PositDirection);
244 pol1.rotateUz(Phot1Direction);
247 G4double Phot2Energy =(1.-epsil)*TotalAvailableEnergy;
248 G4double PositP= sqrt(PositKinEnergy*(PositKinEnergy+2.*electron_mass_c2));
249 G4ThreeVector dir = PositDirection*PositP - Phot1Direction*Phot1Energy;
267 vdp->push_back(aGamma1);
268 vdp->push_back(aGamma2);
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
G4double GetKineticEnergy() const
CLHEP::Hep3Vector G4ThreeVector
G4ParticleChangeForGamma * fParticleChange
G4GLOB_DLL std::ostream G4cout
virtual G4double ComputeCrossSectionPerElectron(const G4ParticleDefinition *, G4double kinEnergy, G4double cutEnergy=0., G4double maxEnergy=DBL_MAX)
G4double GetElectronDensity() const
const G4ThreeVector & GetMomentumDirection() const
static const double twopi
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0., G4double maxEnergy=DBL_MAX)
void SetPolarization(G4double polX, G4double polY, G4double polZ)
G4double G4Log(G4double x)
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4eeToTwoGammaModel(const G4ParticleDefinition *p=0, const G4String &nam="eplus2gg")
G4ParticleDefinition * theGamma
void SetProposedKineticEnergy(G4double proposedKinEnergy)
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)
void ProposeTrackStatus(G4TrackStatus status)
G4ParticleChangeForGamma * GetParticleChangeForGamma()
virtual ~G4eeToTwoGammaModel()