G4eeToHadronsModel Class Reference

#include <G4eeToHadronsModel.hh>

Inheritance diagram for G4eeToHadronsModel:

G4VEmModel

Public Member Functions

 G4eeToHadronsModel (G4Vee2hadrons *, G4int ver=0, const G4String &nam="eeToHadrons")
virtual ~G4eeToHadronsModel ()
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)
virtual G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kineticEnergy, G4double Z, G4double A, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
virtual G4double ComputeCrossSectionPerElectron (const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin=0.0, G4double maxEnergy=DBL_MAX)
G4DynamicParticleGenerateCMPhoton (G4double)
G4double PeakEnergy () const

Detailed Description

Definition at line 59 of file G4eeToHadronsModel.hh.


Constructor & Destructor Documentation

G4eeToHadronsModel::G4eeToHadronsModel ( G4Vee2hadrons ,
G4int  ver = 0,
const G4String nam = "eeToHadrons" 
)

Definition at line 68 of file G4eeToHadronsModel.cc.

References G4Gamma::Gamma(), G4VEmModel::HighEnergyLimit(), and G4VEmModel::LowEnergyLimit().

00070   : G4VEmModel(nam),
00071     model(mod),
00072     crossPerElectron(0),
00073     crossBornPerElectron(0),
00074     isInitialised(false),
00075     nbins(100),
00076     verbose(ver)
00077 {
00078   theGamma = G4Gamma::Gamma();
00079   highKinEnergy = HighEnergyLimit();
00080   lowKinEnergy  = LowEnergyLimit();
00081   emin = lowKinEnergy;
00082   emax = highKinEnergy;
00083   peakKinEnergy = highKinEnergy;
00084   epeak = emax;
00085 }

G4eeToHadronsModel::~G4eeToHadronsModel (  )  [virtual]

Definition at line 89 of file G4eeToHadronsModel.cc.

00090 {
00091   delete model;
00092   delete crossPerElectron;
00093   delete crossBornPerElectron;
00094 }


Member Function Documentation

G4double G4eeToHadronsModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition ,
G4double  kineticEnergy,
G4double  Z,
G4double  A,
G4double  cutEnergy = 0.0,
G4double  maxEnergy = DBL_MAX 
) [virtual]

Reimplemented from G4VEmModel.

Definition at line 195 of file G4eeToHadronsModel.cc.

References ComputeCrossSectionPerElectron().

00200 {
00201   return Z*ComputeCrossSectionPerElectron(p, kineticEnergy);
00202 }

G4double G4eeToHadronsModel::ComputeCrossSectionPerElectron ( const G4ParticleDefinition ,
G4double  kineticEnergy,
G4double  cutEnergy = 0.0,
G4double  maxEnergy = DBL_MAX 
) [virtual]

Definition at line 206 of file G4eeToHadronsModel.cc.

References G4PhysicsVector::GetValue().

Referenced by ComputeCrossSectionPerAtom(), and CrossSectionPerVolume().

00210 {
00211   G4double cross = 0.0;
00212   if(crossPerElectron) {
00213     G4bool b;
00214     G4double e = 2.0*electron_mass_c2*
00215                  sqrt(1.0 + 0.5*kineticEnergy/electron_mass_c2);
00216     cross = crossPerElectron->GetValue(e, b);
00217   }
00218   //  G4cout << "e= " << kineticEnergy << " cross= " << cross << G4endl;
00219   return cross;
00220 }

G4double G4eeToHadronsModel::CrossSectionPerVolume ( const G4Material ,
const G4ParticleDefinition ,
G4double  kineticEnergy,
G4double  cutEnergy,
G4double  maxEnergy 
) [virtual]

Reimplemented from G4VEmModel.

Definition at line 183 of file G4eeToHadronsModel.cc.

References ComputeCrossSectionPerElectron(), and G4Material::GetElectronDensity().

00188 {
00189   return mat->GetElectronDensity()*
00190     ComputeCrossSectionPerElectron(p, kineticEnergy);
00191 }

G4DynamicParticle * G4eeToHadronsModel::GenerateCMPhoton ( G4double   ) 

Definition at line 303 of file G4eeToHadronsModel.cc.

References G4cout, G4endl, G4UniformRand, G4PhysicsVector::GetValue(), G4INCL::Math::pi, and G4InuclParticleNames::s0.

Referenced by SampleSecondaries().

00304 {
00305   G4bool b;   
00306   G4double x;
00307   G4DynamicParticle* gamma = 0;
00308   G4double L   = 2.0*log(e/electron_mass_c2);
00309   G4double bt  = 2.0*fine_structure_const*(L - 1.)/pi;
00310   G4double btm1= bt - 1.0;
00311   G4double del = 1. + fine_structure_const*(1.5*L + pi*pi/3. -2.)/pi;
00312 
00313   G4double s0 = crossBornPerElectron->GetValue(e, b);
00314   G4double de = (emax - emin)/(G4double)nbins;
00315   G4double x0 = min(de,e - emin)/e;
00316   G4double ds = crossBornPerElectron->GetValue(e, b)
00317               *(del*pow(x0,bt) - bt*(x0 - 0.25*x0*x0));
00318   G4double e1 = e*(1. - x0);
00319 
00320   if(e1 < emax && s0*G4UniformRand()<ds) { 
00321     x = x0*pow(G4UniformRand(),1./bt);
00322   } else {    
00323 
00324     x  = 1. - e1/e;
00325     G4double s1 = crossBornPerElectron->GetValue(e1, b);
00326     G4double w1 = bt*(del*pow(x,btm1) - 1.0 + 0.5*x);
00327     G4double grej = s1*w1;
00328     G4double f;
00329     //    G4cout << "e= " << e/GeV << " epeak= " << epeak/GeV 
00330     //       << " s1= " << s1 << " w1= " << w1 
00331     //       << " grej= " << grej << G4endl;
00332     // Above emax cross section is 0
00333     if(e1 > emax) {
00334       x  = 1. - emax/e;
00335       G4double s2 = crossBornPerElectron->GetValue(emax, b);
00336       G4double w2 = bt*(del*pow(x,btm1) - 1.0 + 0.5*x);
00337       grej = s2*w2;
00338       //  G4cout << "emax= " << emax << " s2= " << s2 << " w2= " << w2 
00339       //   << " grej= " << grej << G4endl;
00340     }
00341 
00342     if(e1 > epeak) {
00343       x  = 1. - epeak/e;
00344       G4double s2 = crossBornPerElectron->GetValue(epeak, b);
00345       G4double w2 = bt*(del*pow(x,btm1) - 1.0 + 0.5*x);
00346       grej = max(grej,s2*w2);
00347       //G4cout << "epeak= " << epeak << " s2= " << s2 << " w2= " << w2 
00348       //     << " grej= " << grej << G4endl;
00349     }
00350     G4double xmin = 1. - e1/e;
00351     if(e1 > emax) xmin = 1. - emax/e;
00352     G4double xmax = 1. - emin/e;
00353     do {
00354       x = xmin + G4UniformRand()*(xmax - xmin);
00355       G4double s2 = crossBornPerElectron->GetValue((1.0 - x)*e, b);
00356       G4double w2 = bt*(del*pow(x,btm1) - 1.0 + 0.5*x);
00357       //G4cout << "x= " << x << " xmin= " << xmin << " xmax= " << xmax
00358       //     << " s2= " << s2 << " w2= " << w2 
00359       //           << G4endl;
00360       f = s2*w2;
00361       if(f > grej) {
00362         G4cout << "G4DynamicParticle* G4eeToHadronsModel:WARNING "
00363                << f << " > " << grej << " majorant is`small!" 
00364                << G4endl; 
00365       }
00366     } while (f < grej*G4UniformRand());
00367   }
00368 
00369   G4ThreeVector dir(0.0,0.0,1.0);
00370   gamma = new G4DynamicParticle(theGamma,dir,x*e);
00371   return gamma;
00372 }

void G4eeToHadronsModel::Initialise ( const G4ParticleDefinition ,
const G4DataVector  
) [virtual]

Implements G4VEmModel.

Definition at line 98 of file G4eeToHadronsModel.cc.

References G4Vee2hadrons::ComputeCrossSection(), G4cout, G4endl, G4PhysicsVector::GetLowEdgeEnergy(), G4PhysicsVector::GetValue(), G4PhysicsVector::GetVectorLength(), G4Vee2hadrons::HighEnergy(), G4VEmModel::HighEnergyLimit(), G4Vee2hadrons::LowEnergy(), G4VEmModel::LowEnergyLimit(), G4Vee2hadrons::PeakEnergy(), G4Vee2hadrons::PhysicsVector(), G4PhysicsVector::PutValue(), G4VEmModel::SetHighEnergyLimit(), G4Vee2hadrons::SetLowEnergy(), and G4VEmModel::SetLowEnergyLimit().

00100 {
00101   if(isInitialised) { return; }
00102   isInitialised  = true;
00103 
00104   // Lab system
00105   highKinEnergy = HighEnergyLimit();
00106   lowKinEnergy  = LowEnergyLimit();
00107 
00108   // CM system
00109   emin  = model->LowEnergy();
00110   emax  = model->HighEnergy();
00111 
00112   G4double emin0 = 
00113     2.0*electron_mass_c2*sqrt(1.0 + 0.5*lowKinEnergy/electron_mass_c2);
00114   G4double emax0 = 
00115     2.0*electron_mass_c2*sqrt(1.0 + 0.5*highKinEnergy/electron_mass_c2);
00116 
00117   // recompute low energy
00118   if(emin0 > emax) {
00119     emin0 = emax;
00120     model->SetLowEnergy(emin0);
00121   }
00122   if(emin > emin0) {
00123     emin0 = emin;
00124     lowKinEnergy  = 0.5*emin*emin/electron_mass_c2 - 2.0*electron_mass_c2;
00125     SetLowEnergyLimit(lowKinEnergy);
00126   }
00127 
00128   // recompute high energy
00129   if(emax < emax0) {
00130     emax0 = emax;
00131     highKinEnergy = 0.5*emax*emax/electron_mass_c2 - 2.0*electron_mass_c2;
00132     SetHighEnergyLimit(highKinEnergy);
00133   }
00134 
00135   // peak energy
00136   epeak = std::min(model->PeakEnergy(), emax);
00137   peakKinEnergy  = 0.5*epeak*epeak/electron_mass_c2 - 2.0*electron_mass_c2;
00138 
00139   if(verbose>0) {
00140     G4cout << "G4eeToHadronsModel::Initialise: " << G4endl;
00141     G4cout << "LabSystem: emin(GeV)= " << lowKinEnergy/GeV
00142            << " epeak(GeV)= " << peakKinEnergy/GeV
00143            << " emax(GeV)= " << highKinEnergy/GeV
00144            << G4endl;
00145     G4cout << "SM System: emin(MeV)= " << emin/MeV
00146            << " epeak(MeV)= " << epeak/MeV
00147            << " emax(MeV)= " << emax/MeV
00148            << G4endl;
00149   }
00150 
00151   if(lowKinEnergy < peakKinEnergy) {
00152     crossBornPerElectron = model->PhysicsVector(emin, emax);
00153     crossPerElectron     = model->PhysicsVector(emin, emax);
00154     nbins = crossPerElectron->GetVectorLength();
00155     for(G4int i=0; i<nbins; i++) {
00156       G4double e  = crossPerElectron->GetLowEdgeEnergy(i);
00157       G4double cs = model->ComputeCrossSection(e);
00158       crossBornPerElectron->PutValue(i, cs);
00159     }
00160     ComputeCMCrossSectionPerElectron();
00161   }
00162   if(verbose>1) {
00163     G4cout << "G4eeToHadronsModel: Cross secsions per electron"
00164            << " nbins= " << nbins
00165            << " emin(MeV)= " << emin/MeV
00166            << " emax(MeV)= " << emax/MeV
00167            << G4endl;
00168     G4bool b;
00169     for(G4int i=0; i<nbins; i++) {
00170       G4double e  = crossPerElectron->GetLowEdgeEnergy(i);
00171       G4double s1 = crossPerElectron->GetValue(e, b);
00172       G4double s2 = crossBornPerElectron->GetValue(e, b);
00173       G4cout << "E(MeV)= " << e/MeV
00174              << "  cross(nb)= " << s1/nanobarn
00175              << "  crossBorn(nb)= " << s2/nanobarn
00176              << G4endl;
00177     }
00178   }
00179 }

G4double G4eeToHadronsModel::PeakEnergy (  )  const [inline]

Definition at line 127 of file G4eeToHadronsModel.hh.

00128 {
00129   return peakKinEnergy;
00130 }

void G4eeToHadronsModel::SampleSecondaries ( std::vector< G4DynamicParticle * > *  ,
const G4MaterialCutsCouple ,
const G4DynamicParticle ,
G4double  tmin = 0.0,
G4double  maxEnergy = DBL_MAX 
) [virtual]

Implements G4VEmModel.

Definition at line 224 of file G4eeToHadronsModel.cc.

References GenerateCMPhoton(), G4DynamicParticle::Get4Momentum(), G4DynamicParticle::GetKineticEnergy(), G4DynamicParticle::GetMomentumDirection(), G4Vee2hadrons::SampleSecondaries(), and G4DynamicParticle::Set4Momentum().

00229 {
00230   if(crossPerElectron) {
00231     G4double t = dParticle->GetKineticEnergy();
00232     G4double e = 2.0*electron_mass_c2*sqrt(1.0 + 0.5*t/electron_mass_c2);
00233     G4LorentzVector inlv = dParticle->Get4Momentum();
00234     G4ThreeVector inBoost = inlv.boostVector();
00235     if(e > emin) {
00236       G4DynamicParticle* gamma = GenerateCMPhoton(e);
00237       G4LorentzVector gLv = gamma->Get4Momentum();
00238       G4LorentzVector lv(0.0,0.0,0.0,e);
00239       lv -= gLv;
00240       G4double mass = lv.m();
00241       G4ThreeVector boost = lv.boostVector();
00242       const G4ThreeVector dir = gamma->GetMomentumDirection();
00243       model->SampleSecondaries(newp, mass, dir);
00244       G4int np = newp->size();
00245       for(G4int j=0; j<np; j++) {
00246         G4DynamicParticle* dp = (*newp)[j];
00247         G4LorentzVector v = dp->Get4Momentum();
00248         v.boost(boost);
00249         v.boost(inBoost);
00250         dp->Set4Momentum(v);
00251       }
00252       gLv.boost(inBoost);
00253       gamma->Set4Momentum(gLv);
00254       newp->push_back(gamma);
00255     }
00256   }
00257 }


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:51:50 2013 for Geant4 by  doxygen 1.4.7