Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4BetheBlochModel Class Reference

#include <G4BetheBlochModel.hh>

Inheritance diagram for G4BetheBlochModel:
Collaboration diagram for G4BetheBlochModel:

Public Member Functions

 G4BetheBlochModel (const G4ParticleDefinition *p=nullptr, const G4String &nam="BetheBloch")
 
virtual ~G4BetheBlochModel ()
 
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &) override
 
virtual G4double MinEnergyCut (const G4ParticleDefinition *, const G4MaterialCutsCouple *couple) override
 
virtual G4double ComputeCrossSectionPerElectron (const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
 
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kineticEnergy, G4double Z, G4double A, G4double cutEnergy, G4double maxEnergy) override
 
virtual G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy) override
 
virtual G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy) override
 
virtual G4double GetChargeSquareRatio (const G4ParticleDefinition *p, const G4Material *mat, G4double kineticEnergy) override
 
virtual G4double GetParticleCharge (const G4ParticleDefinition *p, const G4Material *mat, G4double kineticEnergy) override
 
virtual void CorrectionsAlongStep (const G4MaterialCutsCouple *couple, const G4DynamicParticle *dp, G4double &eloss, G4double &, G4double length) override
 
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy) override
 
- 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 GetPartialCrossSection (const G4Material *, G4int level, const G4ParticleDefinition *, G4double kineticEnergy)
 
virtual G4double ComputeCrossSectionPerShell (const G4ParticleDefinition *, G4int Z, G4int shellIdx, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
virtual G4double ChargeSquareRatio (const G4Track &)
 
virtual void StartTracking (G4Track *)
 
virtual G4double Value (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
 
virtual G4double MinPrimaryEnergy (const G4Material *, const G4ParticleDefinition *, G4double cut=0.0)
 
virtual void SetupForMaterial (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual void DefineForRegion (const G4Region *)
 
virtual void ModelDescription (std::ostream &outFile) const
 
void InitialiseElementSelectors (const G4ParticleDefinition *, const G4DataVector &)
 
std::vector
< G4EmElementSelector * > * 
GetElementSelectors ()
 
void SetElementSelectors (std::vector< G4EmElementSelector * > *)
 
virtual 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=nullptr)
 
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 SetFluctuationFlag (G4bool val)
 
void SetMasterThread (G4bool val)
 
G4bool IsMaster () const
 
G4double MaxSecondaryKinEnergy (const G4DynamicParticle *dynParticle)
 
const G4StringGetName () const
 
void SetCurrentCouple (const G4MaterialCutsCouple *)
 
const G4ElementGetCurrentElement () const
 
const G4IsotopeGetCurrentIsotope () const
 
G4bool IsLocked () const
 
void SetLocked (G4bool)
 

Protected Member Functions

virtual G4double MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kinEnergy) override
 
G4double GetChargeSquareRatio () const
 
void SetChargeSquareRatio (G4double val)
 
- Protected Member Functions inherited from G4VEmModel
G4ParticleChangeForLossGetParticleChangeForLoss ()
 
G4ParticleChangeForGammaGetParticleChangeForGamma ()
 
const G4MaterialCutsCoupleCurrentCouple () const
 
void SetCurrentElement (const G4Element *)
 

Additional Inherited Members

- Protected Attributes inherited from G4VEmModel
G4ElementDatafElementData
 
G4VParticleChangepParticleChange
 
G4PhysicsTablexSectionTable
 
const std::vector< G4double > * theDensityFactor
 
const std::vector< G4int > * theDensityIdx
 
size_t idxTable
 
G4bool lossFlucFlag
 
- Static Protected Attributes inherited from G4VEmModel
static const G4double inveplus = 1.0/CLHEP::eplus
 

Detailed Description

Definition at line 72 of file G4BetheBlochModel.hh.

Constructor & Destructor Documentation

G4BetheBlochModel::G4BetheBlochModel ( const G4ParticleDefinition p = nullptr,
const G4String nam = "BetheBloch" 
)
explicit

Definition at line 75 of file G4BetheBlochModel.cc.

77  : G4VEmModel(nam),
78  particle(nullptr),
79  tlimit(DBL_MAX),
80  twoln10(2.0*G4Log(10.0)),
81  isIon(false)
82 {
83  fParticleChange = nullptr;
84  theElectron = G4Electron::Electron();
85  if(p) {
86  SetGenericIon(p);
87  SetParticle(p);
88  } else {
89  SetParticle(theElectron);
90  }
92  nist = G4NistManager::Instance();
94 }
static G4LossTableManager * Instance()
static G4NistManager * Instance()
G4VEmModel(const G4String &nam)
Definition: G4VEmModel.cc:68
G4EmCorrections * EmCorrections()
G4double G4Log(G4double x)
Definition: G4Log.hh:230
static G4Electron * Electron()
Definition: G4Electron.cc:94
static constexpr double MeV
Definition: G4SIunits.hh:214
void SetLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:739
#define DBL_MAX
Definition: templates.hh:83

Here is the call graph for this function:

G4BetheBlochModel::~G4BetheBlochModel ( )
virtual

Definition at line 98 of file G4BetheBlochModel.cc.

99 {}

Member Function Documentation

G4double G4BetheBlochModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition p,
G4double  kineticEnergy,
G4double  Z,
G4double  A,
G4double  cutEnergy,
G4double  maxEnergy 
)
overridevirtual

Reimplemented from G4VEmModel.

Definition at line 216 of file G4BetheBlochModel.cc.

222 {
224  (p,kineticEnergy,cutEnergy,maxEnergy);
225  return cross;
226 }
virtual G4double ComputeCrossSectionPerElectron(const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4BetheBlochModel::ComputeCrossSectionPerElectron ( const G4ParticleDefinition p,
G4double  kineticEnergy,
G4double  cutEnergy,
G4double  maxEnergy 
)
virtual

Definition at line 185 of file G4BetheBlochModel.cc.

189 {
190  G4double cross = 0.0;
191  G4double tmax = MaxSecondaryEnergy(p, kineticEnergy);
192  G4double maxEnergy = min(tmax,maxKinEnergy);
193  if(cutEnergy < maxEnergy) {
194 
195  G4double totEnergy = kineticEnergy + mass;
196  G4double energy2 = totEnergy*totEnergy;
197  G4double beta2 = kineticEnergy*(kineticEnergy + 2.0*mass)/energy2;
198 
199  cross = (maxEnergy - cutEnergy)/(cutEnergy*maxEnergy)
200  - beta2*G4Log(maxEnergy/cutEnergy)/tmax;
201 
202  // +term for spin=1/2 particle
203  if( 0.0 < spin ) { cross += 0.5*(maxEnergy - cutEnergy)/energy2; }
204 
205  cross *= twopi_mc2_rcl2*chargeSquare/beta2;
206  }
207 
208  // G4cout << "BB: e= " << kineticEnergy << " tmin= " << cutEnergy
209  // << " tmax= " << tmax << " cross= " << cross << G4endl;
210 
211  return cross;
212 }
G4double G4Log(G4double x)
Definition: G4Log.hh:230
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
double G4double
Definition: G4Types.hh:76
virtual G4double MaxSecondaryEnergy(const G4ParticleDefinition *, G4double kinEnergy) override

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4BetheBlochModel::ComputeDEDXPerVolume ( const G4Material material,
const G4ParticleDefinition p,
G4double  kineticEnergy,
G4double  cutEnergy 
)
overridevirtual

Reimplemented from G4VEmModel.

Reimplemented in G4BetheBlochNoDeltaModel.

Definition at line 245 of file G4BetheBlochModel.cc.

249 {
250  G4double tmax = MaxSecondaryEnergy(p, kineticEnergy);
251  G4double cutEnergy = std::min(cut,tmax);
252 
253  G4double tau = kineticEnergy/mass;
254  G4double gam = tau + 1.0;
255  G4double bg2 = tau * (tau+2.0);
256  G4double beta2 = bg2/(gam*gam);
257 
258  G4double eexc = material->GetIonisation()->GetMeanExcitationEnergy();
259  G4double eexc2 = eexc*eexc;
260 
261  G4double eDensity = material->GetElectronDensity();
262 
263  G4double dedx = G4Log(2.0*electron_mass_c2*bg2*cutEnergy/eexc2)
264  - (1.0 + cutEnergy/tmax)*beta2;
265 
266  if(0.0 < spin) {
267  G4double del = 0.5*cutEnergy/(kineticEnergy + mass);
268  dedx += del*del;
269  }
270 
271  // density correction
272  G4double x = G4Log(bg2)/twoln10;
273  dedx -= material->GetIonisation()->DensityCorrection(x);
274 
275  // shell correction
276  dedx -= 2.0*corr->ShellCorrection(p,material,kineticEnergy);
277 
278  // now compute the total ionization loss
279  dedx *= twopi_mc2_rcl2*chargeSquare*eDensity/beta2;
280 
281  //High order correction different for hadrons and ions
282  if(isIon) {
283  dedx += corr->IonBarkasCorrection(p,material,kineticEnergy);
284  } else {
285  dedx += corr->HighOrderCorrections(p,material,kineticEnergy,cutEnergy);
286  }
287 
288  dedx = std::max(dedx, 0.0);
289 
290  //G4cout << "E(MeV)= " << kineticEnergy/MeV << " dedx= " << dedx
291  // << " " << material->GetName() << G4endl;
292 
293  return dedx;
294 }
G4IonisParamMat * GetIonisation() const
Definition: G4Material.hh:226
G4double HighOrderCorrections(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy, G4double cutEnergy)
tuple x
Definition: test.py:50
G4double GetElectronDensity() const
Definition: G4Material.hh:217
G4double IonBarkasCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
float electron_mass_c2
Definition: hepunit.py:274
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4double ShellCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4double DensityCorrection(G4double x)
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
G4double GetMeanExcitationEnergy() const
double G4double
Definition: G4Types.hh:76
virtual G4double MaxSecondaryEnergy(const G4ParticleDefinition *, G4double kinEnergy) override

Here is the call graph for this function:

Here is the caller graph for this function:

void G4BetheBlochModel::CorrectionsAlongStep ( const G4MaterialCutsCouple couple,
const G4DynamicParticle dp,
G4double eloss,
G4double ,
G4double  length 
)
overridevirtual

Reimplemented from G4VEmModel.

Definition at line 298 of file G4BetheBlochModel.cc.

303 {
304  if(isIon) {
305  const G4ParticleDefinition* p = dp->GetDefinition();
306  const G4Material* mat = couple->GetMaterial();
307  G4double preKinEnergy = dp->GetKineticEnergy();
308  G4double e = preKinEnergy - eloss*0.5;
309  if(e < preKinEnergy*0.75) { e = preKinEnergy*0.75; }
310 
311  G4double q2 = corr->EffectiveChargeSquareRatio(p,mat,e);
313  G4double qfactor = q2*corr->EffectiveChargeCorrection(p,mat,e)/corrFactor;
314  G4double highOrder = length*corr->IonHighOrderCorrections(p,couple,e);
315  G4double elossnew = eloss*qfactor + highOrder;
316  if(elossnew > preKinEnergy) { elossnew = preKinEnergy; }
317  else if(elossnew < eloss*0.5) { elossnew = eloss*0.5; }
318  eloss = elossnew;
319  //G4cout << "G4BetheBlochModel::CorrectionsAlongStep: e= " << preKinEnergy
320  // << " qfactor= " << qfactor
321  // << " highOrder= " << highOrder << " ("
322  // << highOrder/eloss << ")" << G4endl;
323  }
324 }
G4double GetKineticEnergy() const
G4double EffectiveChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
const char * p
Definition: xmltok.h:285
G4double EffectiveChargeCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4ParticleDefinition * GetDefinition() const
G4VEmFluctuationModel * GetModelOfFluctuations()
Definition: G4VEmModel.hh:612
G4double IonHighOrderCorrections(const G4ParticleDefinition *, const G4MaterialCutsCouple *, G4double kineticEnergy)
double G4double
Definition: G4Types.hh:76
virtual void SetParticleAndCharge(const G4ParticleDefinition *, G4double q2)
const G4Material * GetMaterial() const

Here is the call graph for this function:

G4double G4BetheBlochModel::CrossSectionPerVolume ( const G4Material material,
const G4ParticleDefinition p,
G4double  kineticEnergy,
G4double  cutEnergy,
G4double  maxEnergy 
)
overridevirtual

Reimplemented from G4VEmModel.

Reimplemented in G4BetheBlochNoDeltaModel.

Definition at line 230 of file G4BetheBlochModel.cc.

236 {
237  G4double eDensity = material->GetElectronDensity();
239  (p,kineticEnergy,cutEnergy,maxEnergy);
240  return cross;
241 }
G4double GetElectronDensity() const
Definition: G4Material.hh:217
virtual G4double ComputeCrossSectionPerElectron(const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4BetheBlochModel::GetChargeSquareRatio ( const G4ParticleDefinition p,
const G4Material mat,
G4double  kineticEnergy 
)
overridevirtual

Reimplemented from G4VEmModel.

Definition at line 126 of file G4BetheBlochModel.cc.

129 {
130  // this method is called only for ions
131  G4double q2 = corr->EffectiveChargeSquareRatio(p,mat,kineticEnergy);
132  corrFactor = q2*corr->EffectiveChargeCorrection(p,mat,kineticEnergy);
133  return corrFactor;
134 }
G4double EffectiveChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4double EffectiveChargeCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4BetheBlochModel::GetChargeSquareRatio ( ) const
inlineprotected

Definition at line 194 of file G4BetheBlochModel.hh.

195 {
196  return chargeSquare;
197 }
G4double G4BetheBlochModel::GetParticleCharge ( const G4ParticleDefinition p,
const G4Material mat,
G4double  kineticEnergy 
)
overridevirtual

Reimplemented from G4VEmModel.

Reimplemented in G4BetheBlochIonGasModel.

Definition at line 138 of file G4BetheBlochModel.cc.

141 {
142  //G4cout<<"G4BetheBlochModel::GetParticleCharge e= "<<kineticEnergy <<
143  // " q= " << corr->GetParticleCharge(p,mat,kineticEnergy) <<G4endl;
144  // this method is called only for ions, so no check if it is an ion
145  return corr->GetParticleCharge(p,mat,kineticEnergy);
146 }
G4double GetParticleCharge(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)

Here is the call graph for this function:

void G4BetheBlochModel::Initialise ( const G4ParticleDefinition p,
const G4DataVector  
)
overridevirtual

Implements G4VEmModel.

Definition at line 103 of file G4BetheBlochModel.cc.

105 {
106  SetGenericIon(p);
107  SetParticle(p);
108 
109  //G4cout << "G4BetheBlochModel::Initialise for " << p->GetParticleName()
110  // << " isIon= " << isIon
111  // << G4endl;
112 
113  // always false before the run
114  SetDeexcitationFlag(false);
115 
116  if(nullptr == fParticleChange) {
117  fParticleChange = GetParticleChangeForLoss();
120  }
121  }
122 }
G4ParticleChangeForLoss * GetParticleChangeForLoss()
Definition: G4VEmModel.cc:118
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:619
G4bool UseAngularGeneratorFlag() const
Definition: G4VEmModel.hh:699
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:626
void SetDeexcitationFlag(G4bool val)
Definition: G4VEmModel.hh:788

Here is the call graph for this function:

G4double G4BetheBlochModel::MaxSecondaryEnergy ( const G4ParticleDefinition pd,
G4double  kinEnergy 
)
overrideprotectedvirtual

Reimplemented from G4VEmModel.

Definition at line 444 of file G4BetheBlochModel.cc.

446 {
447  // here particle type is checked for any method
448  SetParticle(pd);
449  G4double tau = kinEnergy/mass;
450  G4double tmax = 2.0*electron_mass_c2*tau*(tau + 2.) /
451  (1. + 2.0*(tau + 1.)*ratio + ratio*ratio);
452  return std::min(tmax,tlimit);
453 }
float electron_mass_c2
Definition: hepunit.py:274
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4BetheBlochModel::MinEnergyCut ( const G4ParticleDefinition ,
const G4MaterialCutsCouple couple 
)
overridevirtual

Reimplemented from G4VEmModel.

Definition at line 176 of file G4BetheBlochModel.cc.

178 {
179  return couple->GetMaterial()->GetIonisation()->GetMeanExcitationEnergy();
180 }
G4IonisParamMat * GetIonisation() const
Definition: G4Material.hh:226
G4double GetMeanExcitationEnergy() const
const G4Material * GetMaterial() const

Here is the call graph for this function:

void G4BetheBlochModel::SampleSecondaries ( std::vector< G4DynamicParticle * > *  vdp,
const G4MaterialCutsCouple couple,
const G4DynamicParticle dp,
G4double  tmin,
G4double  maxEnergy 
)
overridevirtual

Implements G4VEmModel.

Definition at line 328 of file G4BetheBlochModel.cc.

333 {
334  G4double kineticEnergy = dp->GetKineticEnergy();
335  G4double tmax = MaxSecondaryEnergy(dp->GetDefinition(),kineticEnergy);
336 
337  G4double maxKinEnergy = std::min(maxEnergy,tmax);
338  if(minKinEnergy >= maxKinEnergy) { return; }
339 
340  //G4cout << "G4BetheBlochModel::SampleSecondaries Emin= " << minKinEnergy
341  // << " Emax= " << maxKinEnergy << G4endl;
342 
343  G4double totEnergy = kineticEnergy + mass;
344  G4double etot2 = totEnergy*totEnergy;
345  G4double beta2 = kineticEnergy*(kineticEnergy + 2.0*mass)/etot2;
346 
347  G4double deltaKinEnergy, f;
348  G4double f1 = 0.0;
349  G4double fmax = 1.0;
350  if( 0.0 < spin ) { fmax += 0.5*maxKinEnergy*maxKinEnergy/etot2; }
351 
352  CLHEP::HepRandomEngine* rndmEngineMod = G4Random::getTheEngine();
353  G4double rndm[2];
354 
355  // sampling without nuclear size effect
356  do {
357  rndmEngineMod->flatArray(2, rndm);
358  deltaKinEnergy = minKinEnergy*maxKinEnergy
359  /(minKinEnergy*(1.0 - rndm[0]) + maxKinEnergy*rndm[0]);
360 
361  f = 1.0 - beta2*deltaKinEnergy/tmax;
362  if( 0.0 < spin ) {
363  f1 = 0.5*deltaKinEnergy*deltaKinEnergy/etot2;
364  f += f1;
365  }
366 
367  // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
368  } while( fmax*rndm[1] > f);
369 
370  // projectile formfactor - suppresion of high energy
371  // delta-electron production at high energy
372 
373  G4double x = formfact*deltaKinEnergy;
374  if(x > 1.e-6) {
375 
376  G4double x1 = 1.0 + x;
377  G4double grej = 1.0/(x1*x1);
378  if( 0.0 < spin ) {
379  G4double x2 = 0.5*electron_mass_c2*deltaKinEnergy/(mass*mass);
380  grej *= (1.0 + magMoment2*(x2 - f1/f)/(1.0 + x2));
381  }
382  if(grej > 1.1) {
383  G4cout << "### G4BetheBlochModel WARNING: grej= " << grej
384  << " " << dp->GetDefinition()->GetParticleName()
385  << " Ekin(MeV)= " << kineticEnergy
386  << " delEkin(MeV)= " << deltaKinEnergy
387  << G4endl;
388  }
389  if(rndmEngineMod->flat() > grej) { return; }
390  }
391 
392  G4ThreeVector deltaDirection;
393 
395 
396  const G4Material* mat = couple->GetMaterial();
398 
399  deltaDirection =
400  GetAngularDistribution()->SampleDirection(dp, deltaKinEnergy, Z, mat);
401 
402  } else {
403 
404  G4double deltaMomentum =
405  sqrt(deltaKinEnergy * (deltaKinEnergy + 2.0*electron_mass_c2));
406  G4double cost = deltaKinEnergy * (totEnergy + electron_mass_c2) /
407  (deltaMomentum * dp->GetTotalMomentum());
408  if(cost > 1.0) { cost = 1.0; }
409  G4double sint = sqrt((1.0 - cost)*(1.0 + cost));
410 
411  G4double phi = twopi*rndmEngineMod->flat();
412 
413  deltaDirection.set(sint*cos(phi),sint*sin(phi), cost) ;
414  deltaDirection.rotateUz(dp->GetMomentumDirection());
415  }
416  /*
417  G4cout << "### G4BetheBlochModel "
418  << dp->GetDefinition()->GetParticleName()
419  << " Ekin(MeV)= " << kineticEnergy
420  << " delEkin(MeV)= " << deltaKinEnergy
421  << " tmin(MeV)= " << minKinEnergy
422  << " tmax(MeV)= " << maxKinEnergy
423  << " dir= " << dp->GetMomentumDirection()
424  << " dirDelta= " << deltaDirection
425  << G4endl;
426  */
427  // create G4DynamicParticle object for delta ray
428  G4DynamicParticle* delta =
429  new G4DynamicParticle(theElectron,deltaDirection,deltaKinEnergy);
430 
431  vdp->push_back(delta);
432 
433  // Change kinematics of primary particle
434  kineticEnergy -= deltaKinEnergy;
435  G4ThreeVector finalP = dp->GetMomentum() - delta->GetMomentum();
436  finalP = finalP.unit();
437 
438  fParticleChange->SetProposedKineticEnergy(kineticEnergy);
439  fParticleChange->SetProposedMomentumDirection(finalP);
440 }
void set(double x, double y, double z)
G4double GetKineticEnergy() const
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:619
G4ParticleDefinition * GetDefinition() const
tuple x
Definition: test.py:50
virtual double flat()=0
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
static constexpr double twopi
Definition: G4SIunits.hh:76
G4double GetTotalMomentum() const
G4GLOB_DLL std::ostream G4cout
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
G4bool UseAngularGeneratorFlag() const
Definition: G4VEmModel.hh:699
const G4ThreeVector & GetMomentumDirection() const
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:38
void SetProposedKineticEnergy(G4double proposedKinEnergy)
float electron_mass_c2
Definition: hepunit.py:274
void SetProposedMomentumDirection(const G4ThreeVector &dir)
Hep3Vector unit() const
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
virtual void flatArray(const int size, double *vect)=0
virtual G4double MaxSecondaryEnergy(const G4ParticleDefinition *, G4double kinEnergy) override
G4ThreeVector GetMomentum() const
const G4Material * GetMaterial() const
G4int SelectRandomAtomNumber(const G4Material *)
Definition: G4VEmModel.hh:564

Here is the call graph for this function:

void G4BetheBlochModel::SetChargeSquareRatio ( G4double  val)
inlineprotected

Definition at line 201 of file G4BetheBlochModel.hh.

202 {
203  chargeSquare = val;
204 }

Here is the caller graph for this function:


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