CMS 3D CMS Logo

Public Member Functions | Private Member Functions | Private Attributes

PairProductionSimulator Class Reference

#include <PairProductionSimulator.h>

Inheritance diagram for PairProductionSimulator:
MaterialEffectsSimulator

List of all members.

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.

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.

{}

Member Function Documentation

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

Member Data Documentation

The minimal photon energy for possible conversion.

Definition at line 38 of file PairProductionSimulator.h.

Referenced by compute(), and PairProductionSimulator().