#include <PairProductionSimulator.h>

Inheritance diagram for PairProductionSimulator:

Public Member Functions
	PairProductionSimulator (double photonEnergyCut, const RandomEngine *engine)
	Constructor. More...

	~PairProductionSimulator ()
	Default Destructor. More...

Public Member Functions inherited from MaterialEffectsSimulator
RHEP_const_iter	beginDaughters () const
	Returns const iterator to the beginning of the daughters list. More...

int	closestDaughterId ()
	The id of the closest charged daughter (filled for nuclear interactions only) More...

double	eMass () const
	Electron mass in GeV/c2. More...

RHEP_const_iter	endDaughters () const
	Returns const iterator to the end of the daughters list. More...

double	excitE () const
	Mean excitation energy (in GeV) More...

	MaterialEffectsSimulator (const RandomEngine *engine, double A=28.0855, double Z=14.0000, double density=2.329, double radLen=9.360)

unsigned	nDaughters () const
	Returns the number of daughters. More...

XYZVector	orthogonal (const XYZVector &) const
	A vector orthogonal to another one (because it's not in XYZTLorentzVector) More...

double	radLenIncm () const
	One radiation length in cm. More...

double	rho () const
	Density in g/cm3. More...

void	setNormalVector (const GlobalVector &normal)
	Sets the vector normal to the surface traversed. More...

double	theA () const
	A. More...

double	theZ () const
	Z. More...

void	updateState (ParticlePropagator &myTrack, double radlen)
	Compute the material effect (calls the sub class) More...

virtual	~MaterialEffectsSimulator ()

Private Member Functions
void	compute (ParticlePropagator &Particle)
	Generate an e+e- pair according to the probability that it happens. More...

double	gbteth (double ener, double partm, double efrac)
	A universal angular distribution - still from GEANT. More...

Private Attributes
double	photonEnergy
	The minimal photon energy for possible conversion. More...

Additional Inherited Members
Public Types inherited from MaterialEffectsSimulator
typedef std::vector < RawParticle > ::const_iterator	RHEP_const_iter

Protected Attributes inherited from MaterialEffectsSimulator
std::vector< RawParticle >	_theUpdatedState

double	A

double	density

double	radLen

double	radLengths

const RandomEngine *	random

int	theClosestChargedDaughterId

GlobalVector	theNormalVector

double	Z

Detailed Description

Definition at line 24 of file PairProductionSimulator.h.

Constructor & Destructor Documentation

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

Default Destructor.

Definition at line 33 of file PairProductionSimulator.h.

33 {}

Member Function Documentation

void PairProductionSimulator::compute ( ParticlePropagator & Particle )

privatevirtual

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;
 }

Member Data Documentation

double PairProductionSimulator::photonEnergy

private

The minimal photon energy for possible conversion.

Definition at line 38 of file PairProductionSimulator.h.

Referenced by compute(), and PairProductionSimulator().

Public Member Functions

Private Member Functions

Private Attributes

Additional Inherited Members

Detailed Description

Constructor & Destructor Documentation

Member Function Documentation

Member Data Documentation