Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions
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
 
- Public Member Functions inherited from G4VEmModel
 G4VEmModel (const G4String &nam)
 
virtual ~G4VEmModel ()
 
virtual void InitialiseLocal (const G4ParticleDefinition *, G4VEmModel *masterModel)
 
virtual void InitialiseForMaterial (const G4ParticleDefinition *, const G4Material *)
 
virtual void InitialiseForElement (const G4ParticleDefinition *, G4int Z)
 
virtual G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
 
virtual G4double ChargeSquareRatio (const G4Track &)
 
virtual G4double GetChargeSquareRatio (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual G4double GetParticleCharge (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual void StartTracking (G4Track *)
 
virtual void CorrectionsAlongStep (const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double &eloss, G4double &niel, G4double length)
 
virtual G4double Value (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
 
virtual G4double MinPrimaryEnergy (const G4Material *, const G4ParticleDefinition *, G4double cut=0.0)
 
virtual G4double MinEnergyCut (const G4ParticleDefinition *, const G4MaterialCutsCouple *)
 
virtual void SetupForMaterial (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual void DefineForRegion (const G4Region *)
 
void InitialiseElementSelectors (const G4ParticleDefinition *, const G4DataVector &)
 
std::vector
< G4EmElementSelector * > * 
GetElementSelectors ()
 
void SetElementSelectors (std::vector< G4EmElementSelector * > *)
 
G4double ComputeDEDX (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
 
G4double CrossSection (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeMeanFreePath (const G4ParticleDefinition *, G4double kineticEnergy, const G4Material *, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, const G4Element *, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4int SelectIsotopeNumber (const G4Element *)
 
const G4ElementSelectRandomAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementSelectRandomAtom (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4int SelectRandomAtomNumber (const G4Material *)
 
void SetParticleChange (G4VParticleChange *, G4VEmFluctuationModel *f=0)
 
void SetCrossSectionTable (G4PhysicsTable *, G4bool isLocal)
 
G4ElementDataGetElementData ()
 
G4PhysicsTableGetCrossSectionTable ()
 
G4VEmFluctuationModelGetModelOfFluctuations ()
 
G4VEmAngularDistributionGetAngularDistribution ()
 
void SetAngularDistribution (G4VEmAngularDistribution *)
 
G4double HighEnergyLimit () const
 
G4double LowEnergyLimit () const
 
G4double HighEnergyActivationLimit () const
 
G4double LowEnergyActivationLimit () const
 
G4double PolarAngleLimit () const
 
G4double SecondaryThreshold () const
 
G4bool LPMFlag () const
 
G4bool DeexcitationFlag () const
 
G4bool ForceBuildTableFlag () const
 
G4bool UseAngularGeneratorFlag () const
 
void SetAngularGeneratorFlag (G4bool)
 
void SetHighEnergyLimit (G4double)
 
void SetLowEnergyLimit (G4double)
 
void SetActivationHighEnergyLimit (G4double)
 
void SetActivationLowEnergyLimit (G4double)
 
G4bool IsActive (G4double kinEnergy)
 
void SetPolarAngleLimit (G4double)
 
void SetSecondaryThreshold (G4double)
 
void SetLPMFlag (G4bool val)
 
void SetDeexcitationFlag (G4bool val)
 
void SetForceBuildTable (G4bool val)
 
void SetMasterThread (G4bool val)
 
G4bool IsMaster () const
 
G4double MaxSecondaryKinEnergy (const G4DynamicParticle *dynParticle)
 
const G4StringGetName () const
 
void SetCurrentCouple (const G4MaterialCutsCouple *)
 
const G4ElementGetCurrentElement () const
 

Additional Inherited Members

- Protected Member Functions inherited from G4VEmModel
G4ParticleChangeForLossGetParticleChangeForLoss ()
 
G4ParticleChangeForGammaGetParticleChangeForGamma ()
 
virtual G4double MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kineticEnergy)
 
const G4MaterialCutsCoupleCurrentCouple () const
 
void SetCurrentElement (const G4Element *)
 
- Protected Attributes inherited from G4VEmModel
G4ElementDatafElementData
 
G4VParticleChangepParticleChange
 
G4PhysicsTablexSectionTable
 
const std::vector< G4double > * theDensityFactor
 
const std::vector< G4int > * theDensityIdx
 
size_t idxTable
 

Detailed Description

Definition at line 59 of file G4eeToHadronsModel.hh.

Constructor & Destructor Documentation

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

Definition at line 68 of file G4eeToHadronsModel.cc.

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

70  : G4VEmModel(nam),
71  model(mod),
72  crossPerElectron(0),
73  crossBornPerElectron(0),
74  isInitialised(false),
75  nbins(100),
76  verbose(ver)
77 {
78  theGamma = G4Gamma::Gamma();
79  highKinEnergy = HighEnergyLimit();
80  lowKinEnergy = LowEnergyLimit();
81  emin = lowKinEnergy;
82  emax = highKinEnergy;
83  peakKinEnergy = highKinEnergy;
84  epeak = emax;
85 }
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:599
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:592
G4VEmModel(const G4String &nam)
Definition: G4VEmModel.cc:65
const XML_Char XML_Content * model
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
G4eeToHadronsModel::~G4eeToHadronsModel ( )
virtual

Definition at line 89 of file G4eeToHadronsModel.cc.

90 {
91  delete model;
92  delete crossPerElectron;
93  delete crossBornPerElectron;
94 }
const XML_Char XML_Content * model

Member Function Documentation

G4double G4eeToHadronsModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition p,
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().

200 {
201  return Z*ComputeCrossSectionPerElectron(p, kineticEnergy);
202 }
virtual G4double ComputeCrossSectionPerElectron(const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
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 test::b, python.hepunit::electron_mass_c2, and G4PhysicsVector::GetValue().

Referenced by ComputeCrossSectionPerAtom(), and CrossSectionPerVolume().

210 {
211  G4double cross = 0.0;
212  if(crossPerElectron) {
213  G4bool b;
214  G4double e = 2.0*electron_mass_c2*
215  sqrt(1.0 + 0.5*kineticEnergy/electron_mass_c2);
216  cross = crossPerElectron->GetValue(e, b);
217  }
218  // G4cout << "e= " << kineticEnergy << " cross= " << cross << G4endl;
219  return cross;
220 }
G4double GetValue(G4double theEnergy, G4bool &isOutRange) const
bool G4bool
Definition: G4Types.hh:79
float electron_mass_c2
Definition: hepunit.py:274
double G4double
Definition: G4Types.hh:76
G4double G4eeToHadronsModel::CrossSectionPerVolume ( const G4Material mat,
const G4ParticleDefinition p,
G4double  kineticEnergy,
G4double  cutEnergy,
G4double  maxEnergy 
)
virtual

Reimplemented from G4VEmModel.

Definition at line 183 of file G4eeToHadronsModel.cc.

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

188 {
189  return mat->GetElectronDensity()*
190  ComputeCrossSectionPerElectron(p, kineticEnergy);
191 }
G4double GetElectronDensity() const
Definition: G4Material.hh:215
virtual G4double ComputeCrossSectionPerElectron(const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
G4DynamicParticle * G4eeToHadronsModel::GenerateCMPhoton ( G4double  e)

Definition at line 303 of file G4eeToHadronsModel.cc.

References test::b, python.hepunit::electron_mass_c2, python.hepunit::fine_structure_const, G4cout, G4endl, G4UniformRand, G4PhysicsVector::GetValue(), G4INCL::Math::max(), G4INCL::Math::min(), python.hepunit::pi, G4InuclParticleNames::s0, and test::x.

Referenced by SampleSecondaries().

304 {
305  G4bool b;
306  G4double x;
307  G4DynamicParticle* gamma = 0;
308  G4double L = 2.0*log(e/electron_mass_c2);
309  G4double bt = 2.0*fine_structure_const*(L - 1.)/pi;
310  G4double btm1= bt - 1.0;
311  G4double del = 1. + fine_structure_const*(1.5*L + pi*pi/3. -2.)/pi;
312 
313  G4double s0 = crossBornPerElectron->GetValue(e, b);
314  G4double de = (emax - emin)/(G4double)nbins;
315  G4double x0 = min(de,e - emin)/e;
316  G4double ds = crossBornPerElectron->GetValue(e, b)
317  *(del*pow(x0,bt) - bt*(x0 - 0.25*x0*x0));
318  G4double e1 = e*(1. - x0);
319 
320  if(e1 < emax && s0*G4UniformRand()<ds) {
321  x = x0*pow(G4UniformRand(),1./bt);
322  } else {
323 
324  x = 1. - e1/e;
325  G4double s1 = crossBornPerElectron->GetValue(e1, b);
326  G4double w1 = bt*(del*pow(x,btm1) - 1.0 + 0.5*x);
327  G4double grej = s1*w1;
328  G4double f;
329  // G4cout << "e= " << e/GeV << " epeak= " << epeak/GeV
330  // << " s1= " << s1 << " w1= " << w1
331  // << " grej= " << grej << G4endl;
332  // Above emax cross section is 0
333  if(e1 > emax) {
334  x = 1. - emax/e;
335  G4double s2 = crossBornPerElectron->GetValue(emax, b);
336  G4double w2 = bt*(del*pow(x,btm1) - 1.0 + 0.5*x);
337  grej = s2*w2;
338  // G4cout << "emax= " << emax << " s2= " << s2 << " w2= " << w2
339  // << " grej= " << grej << G4endl;
340  }
341 
342  if(e1 > epeak) {
343  x = 1. - epeak/e;
344  G4double s2 = crossBornPerElectron->GetValue(epeak, b);
345  G4double w2 = bt*(del*pow(x,btm1) - 1.0 + 0.5*x);
346  grej = max(grej,s2*w2);
347  //G4cout << "epeak= " << epeak << " s2= " << s2 << " w2= " << w2
348  // << " grej= " << grej << G4endl;
349  }
350  G4double xmin = 1. - e1/e;
351  if(e1 > emax) xmin = 1. - emax/e;
352  G4double xmax = 1. - emin/e;
353  do {
354  x = xmin + G4UniformRand()*(xmax - xmin);
355  G4double s2 = crossBornPerElectron->GetValue((1.0 - x)*e, b);
356  G4double w2 = bt*(del*pow(x,btm1) - 1.0 + 0.5*x);
357  //G4cout << "x= " << x << " xmin= " << xmin << " xmax= " << xmax
358  // << " s2= " << s2 << " w2= " << w2
359  // << G4endl;
360  f = s2*w2;
361  if(f > grej) {
362  G4cout << "G4DynamicParticle* G4eeToHadronsModel:WARNING "
363  << f << " > " << grej << " majorant is`small!"
364  << G4endl;
365  }
366  } while (f < grej*G4UniformRand());
367  }
368 
369  G4ThreeVector dir(0.0,0.0,1.0);
370  gamma = new G4DynamicParticle(theGamma,dir,x*e);
371  return gamma;
372 }
G4double GetValue(G4double theEnergy, G4bool &isOutRange) const
#define G4UniformRand()
Definition: Randomize.hh:87
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
float electron_mass_c2
Definition: hepunit.py:274
T max(const T t1, const T t2)
brief Return the largest of the two arguments
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
void G4eeToHadronsModel::Initialise ( const G4ParticleDefinition ,
const G4DataVector  
)
virtual

Implements G4VEmModel.

Definition at line 98 of file G4eeToHadronsModel.cc.

References test::b, python.hepunit::electron_mass_c2, G4cout, G4endl, G4PhysicsVector::GetLowEdgeEnergy(), G4PhysicsVector::GetValue(), G4PhysicsVector::GetVectorLength(), python.hepunit::GeV, G4VEmModel::HighEnergyLimit(), G4VEmModel::LowEnergyLimit(), python.hepunit::MeV, G4INCL::Math::min(), python.hepunit::nanobarn, G4PhysicsVector::PutValue(), G4VEmModel::SetHighEnergyLimit(), and G4VEmModel::SetLowEnergyLimit().

100 {
101  if(isInitialised) { return; }
102  isInitialised = true;
103 
104  // Lab system
105  highKinEnergy = HighEnergyLimit();
106  lowKinEnergy = LowEnergyLimit();
107 
108  // CM system
109  emin = model->LowEnergy();
110  emax = model->HighEnergy();
111 
112  G4double emin0 =
113  2.0*electron_mass_c2*sqrt(1.0 + 0.5*lowKinEnergy/electron_mass_c2);
114  G4double emax0 =
115  2.0*electron_mass_c2*sqrt(1.0 + 0.5*highKinEnergy/electron_mass_c2);
116 
117  // recompute low energy
118  if(emin0 > emax) {
119  emin0 = emax;
120  model->SetLowEnergy(emin0);
121  }
122  if(emin > emin0) {
123  emin0 = emin;
124  lowKinEnergy = 0.5*emin*emin/electron_mass_c2 - 2.0*electron_mass_c2;
125  SetLowEnergyLimit(lowKinEnergy);
126  }
127 
128  // recompute high energy
129  if(emax < emax0) {
130  emax0 = emax;
131  highKinEnergy = 0.5*emax*emax/electron_mass_c2 - 2.0*electron_mass_c2;
132  SetHighEnergyLimit(highKinEnergy);
133  }
134 
135  // peak energy
136  epeak = std::min(model->PeakEnergy(), emax);
137  peakKinEnergy = 0.5*epeak*epeak/electron_mass_c2 - 2.0*electron_mass_c2;
138 
139  if(verbose>0) {
140  G4cout << "G4eeToHadronsModel::Initialise: " << G4endl;
141  G4cout << "LabSystem: emin(GeV)= " << lowKinEnergy/GeV
142  << " epeak(GeV)= " << peakKinEnergy/GeV
143  << " emax(GeV)= " << highKinEnergy/GeV
144  << G4endl;
145  G4cout << "SM System: emin(MeV)= " << emin/MeV
146  << " epeak(MeV)= " << epeak/MeV
147  << " emax(MeV)= " << emax/MeV
148  << G4endl;
149  }
150 
151  if(lowKinEnergy < peakKinEnergy) {
152  crossBornPerElectron = model->PhysicsVector(emin, emax);
153  crossPerElectron = model->PhysicsVector(emin, emax);
154  nbins = crossPerElectron->GetVectorLength();
155  for(G4int i=0; i<nbins; i++) {
156  G4double e = crossPerElectron->GetLowEdgeEnergy(i);
157  G4double cs = model->ComputeCrossSection(e);
158  crossBornPerElectron->PutValue(i, cs);
159  }
160  ComputeCMCrossSectionPerElectron();
161  }
162  if(verbose>1) {
163  G4cout << "G4eeToHadronsModel: Cross secsions per electron"
164  << " nbins= " << nbins
165  << " emin(MeV)= " << emin/MeV
166  << " emax(MeV)= " << emax/MeV
167  << G4endl;
168  G4bool b;
169  for(G4int i=0; i<nbins; i++) {
170  G4double e = crossPerElectron->GetLowEdgeEnergy(i);
171  G4double s1 = crossPerElectron->GetValue(e, b);
172  G4double s2 = crossBornPerElectron->GetValue(e, b);
173  G4cout << "E(MeV)= " << e/MeV
174  << " cross(nb)= " << s1/nanobarn
175  << " crossBorn(nb)= " << s2/nanobarn
176  << G4endl;
177  }
178  }
179 }
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:599
G4double GetValue(G4double theEnergy, G4bool &isOutRange) const
int nanobarn
Definition: hepunit.py:42
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:592
size_t GetVectorLength() const
G4double GetLowEdgeEnergy(size_t binNumber) const
int G4int
Definition: G4Types.hh:78
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:683
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
void PutValue(size_t index, G4double theValue)
const XML_Char XML_Content * model
float electron_mass_c2
Definition: hepunit.py:274
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
void SetLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:690
G4double G4eeToHadronsModel::PeakEnergy ( ) const
inline

Definition at line 127 of file G4eeToHadronsModel.hh.

128 {
129  return peakKinEnergy;
130 }
void G4eeToHadronsModel::SampleSecondaries ( std::vector< G4DynamicParticle * > *  newp,
const G4MaterialCutsCouple ,
const G4DynamicParticle dParticle,
G4double  tmin = 0.0,
G4double  maxEnergy = DBL_MAX 
)
virtual

Implements G4VEmModel.

Definition at line 224 of file G4eeToHadronsModel.cc.

References CLHEP::HepLorentzVector::boost(), CLHEP::HepLorentzVector::boostVector(), python.hepunit::electron_mass_c2, GenerateCMPhoton(), G4DynamicParticle::Get4Momentum(), G4DynamicParticle::GetKineticEnergy(), G4DynamicParticle::GetMomentumDirection(), CLHEP::HepLorentzVector::m(), G4DynamicParticle::Set4Momentum(), and test::v.

229 {
230  if(crossPerElectron) {
231  G4double t = dParticle->GetKineticEnergy();
232  G4double e = 2.0*electron_mass_c2*sqrt(1.0 + 0.5*t/electron_mass_c2);
233  G4LorentzVector inlv = dParticle->Get4Momentum();
234  G4ThreeVector inBoost = inlv.boostVector();
235  if(e > emin) {
237  G4LorentzVector gLv = gamma->Get4Momentum();
238  G4LorentzVector lv(0.0,0.0,0.0,e);
239  lv -= gLv;
240  G4double mass = lv.m();
241  G4ThreeVector boost = lv.boostVector();
242  const G4ThreeVector dir = gamma->GetMomentumDirection();
243  model->SampleSecondaries(newp, mass, dir);
244  G4int np = newp->size();
245  for(G4int j=0; j<np; j++) {
246  G4DynamicParticle* dp = (*newp)[j];
248  v.boost(boost);
249  v.boost(inBoost);
250  dp->Set4Momentum(v);
251  }
252  gLv.boost(inBoost);
253  gamma->Set4Momentum(gLv);
254  newp->push_back(gamma);
255  }
256  }
257 }
Hep3Vector boostVector() const
G4double GetKineticEnergy() const
int G4int
Definition: G4Types.hh:78
const G4ThreeVector & GetMomentumDirection() const
HepLorentzVector & boost(double, double, double)
const XML_Char XML_Content * model
float electron_mass_c2
Definition: hepunit.py:274
G4LorentzVector Get4Momentum() const
void Set4Momentum(const G4LorentzVector &momentum)
G4DynamicParticle * GenerateCMPhoton(G4double)
double G4double
Definition: G4Types.hh:76

The documentation for this class was generated from the following files: