#include <PairProductionSimulator.h>
Public Member Functions | |
PairProductionSimulator (double photonEnergyCut, const RandomEngine *engine) | |
Constructor. | |
~PairProductionSimulator () | |
Default Destructor. | |
Private Member Functions | |
void | compute (ParticlePropagator &Particle) |
Generate an e+e- pair according to the probability that it happens. | |
double | gbteth (double ener, double partm, double efrac) |
A universal angular distribution - still from GEANT. | |
Private Attributes | |
double | photonEnergy |
The minimal photon energy for possible conversion. |
Definition at line 24 of file PairProductionSimulator.h.
PairProductionSimulator::PairProductionSimulator | ( | double | photonEnergyCut, |
const RandomEngine * | engine | ||
) |
Constructor.
Definition at line 6 of file PairProductionSimulator.cc.
References max(), and photonEnergy.
: MaterialEffectsSimulator(engine) { // Set the minimal photon energy for possible conversion photonEnergy = std::max(0.100,photonEnergyCut); }
PairProductionSimulator::~PairProductionSimulator | ( | ) | [inline] |
void PairProductionSimulator::compute | ( | ParticlePropagator & | Particle | ) | [private, virtual] |
Generate an e+e- pair according to the probability that it happens.
Implements MaterialEffectsSimulator.
Definition at line 18 of file PairProductionSimulator.cc.
References MaterialEffectsSimulator::_theUpdatedState, funct::cos(), MaterialEffectsSimulator::eMass(), RandomEngine::flatShoot(), gbteth(), funct::log(), M_PI, max(), phi, photonEnergy, MaterialEffectsSimulator::radLengths, MaterialEffectsSimulator::random, funct::sin(), mathSSE::sqrt(), and CommonMethods::weight().
{ double eGamma = Particle.e(); // The photon has enough energy to create a pair if ( eGamma>=photonEnergy ) { // This is a simple version (a la PDG) of a photon conversion generator. // It replaces the buggy GEANT3 -> C++ former version. // Author : Patrick Janot - 7-Jan-2004 // Probability to convert is 7/9*(dx/X0) if ( -std::log(random->flatShoot()) <= (7./9.)*radLengths ) { double xe=0; double xm=eMass()/eGamma; double weight = 0.; // Generate electron energy between emass and eGamma-emass do { xe = random->flatShoot()*(1.-2.*xm) + xm; weight = 1. - 4./3.*xe*(1.-xe); } while ( weight < random->flatShoot() ); double eElectron = xe * eGamma; double tElectron = eElectron-eMass(); double pElectron = std::sqrt(std::max((eElectron+eMass())*tElectron,0.)); double ePositron = eGamma-eElectron; double tPositron = ePositron-eMass(); double pPositron = std::sqrt((ePositron+eMass())*tPositron); // Generate angles double phi = random->flatShoot()*2.*M_PI; double sphi = std::sin(phi); double cphi = std::cos(phi); double stheta1, stheta2, ctheta1, ctheta2; if ( eElectron > ePositron ) { double theta1 = gbteth(eElectron,eMass(),xe)*eMass()/eElectron; stheta1 = std::sin(theta1); ctheta1 = std::cos(theta1); stheta2 = stheta1*pElectron/pPositron; ctheta2 = std::sqrt(std::max(0.,1.0-(stheta2*stheta2))); } else { double theta2 = gbteth(ePositron,eMass(),xe)*eMass()/ePositron; stheta2 = std::sin(theta2); ctheta2 = std::cos(theta2); stheta1 = stheta2*pPositron/pElectron; ctheta1 = std::sqrt(std::max(0.,1.0-(stheta1*stheta1))); } double chi = Particle.theta(); double psi = Particle.phi(); RawParticle::RotationZ rotZ(psi); RawParticle::RotationY rotY(chi); _theUpdatedState.resize(2,RawParticle()); // The eletron _theUpdatedState[0].SetXYZT(pElectron*stheta1*cphi, pElectron*stheta1*sphi, pElectron*ctheta1, eElectron); _theUpdatedState[0].setID(+11); _theUpdatedState[0].rotate(rotY); _theUpdatedState[0].rotate(rotZ); // The positron _theUpdatedState[1].SetXYZT(-pPositron*stheta2*cphi, -pPositron*stheta2*sphi, pPositron*ctheta2, ePositron); _theUpdatedState[1].setID(-11); _theUpdatedState[1].rotate(rotY); _theUpdatedState[1].rotate(rotZ); } } }
double PairProductionSimulator::gbteth | ( | double | ener, |
double | partm, | ||
double | efrac | ||
) | [private] |
A universal angular distribution - still from GEANT.
Definition at line 103 of file PairProductionSimulator.cc.
References beta, RandomEngine::flatShoot(), funct::log(), M_PI, MaterialEffectsSimulator::random, and MaterialEffectsSimulator::theZ().
Referenced by compute().
{ const double alfa = 0.625; double d = 0.13*(0.8+1.3/theZ())*(100.0+(1.0/ener))*(1.0+efrac); double w1 = 9.0/(9.0+d); double umax = ener*M_PI/partm; double u; do { double beta; if (random->flatShoot()<=w1) beta = alfa; else beta = 3.0*alfa; u = -(std::log(random->flatShoot()*random->flatShoot()))/beta; } while (u>=umax); return u; }
double PairProductionSimulator::photonEnergy [private] |
The minimal photon energy for possible conversion.
Definition at line 38 of file PairProductionSimulator.h.
Referenced by compute(), and PairProductionSimulator().