Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4eBremsstrahlungModel Class Reference

#include <G4eBremsstrahlungModel.hh>

+ Inheritance diagram for G4eBremsstrahlungModel:

Public Member Functions

 G4eBremsstrahlungModel (const G4ParticleDefinition *p=0, const G4String &nam="eBrem")
 
virtual ~G4eBremsstrahlungModel ()
 
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)
 
virtual G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy)
 
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double tkin, G4double Z, G4double, G4double cut, G4double maxE=DBL_MAX)
 
virtual G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
 
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
 
- Public Member Functions inherited from G4VEmModel
 G4VEmModel (const G4String &nam)
 
virtual ~G4VEmModel ()
 
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)=0
 
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin=0.0, G4double tmax=DBL_MAX)=0
 
virtual G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
 
virtual G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=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 *)
 
virtual void SetupForMaterial (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual void DefineForRegion (const G4Region *)
 
void InitialiseElementSelectors (const G4ParticleDefinition *, const G4DataVector &)
 
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)
 
void SetParticleChange (G4VParticleChange *, G4VEmFluctuationModel *f=0)
 
void SetCrossSectionTable (G4PhysicsTable *)
 
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
 
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 ForceBuildTable (G4bool val)
 
G4double MaxSecondaryKinEnergy (const G4DynamicParticle *dynParticle)
 
const G4StringGetName () const
 
void SetCurrentCouple (const G4MaterialCutsCouple *)
 
const G4ElementGetCurrentElement () const
 

Protected Member Functions

const G4ElementSelectRandomAtom (const G4MaterialCutsCouple *couple)
 
- Protected Member Functions inherited from G4VEmModel
G4ParticleChangeForLossGetParticleChangeForLoss ()
 
G4ParticleChangeForGammaGetParticleChangeForGamma ()
 
virtual G4double MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kineticEnergy)
 
const G4MaterialCutsCoupleCurrentCouple () const
 
void SetCurrentElement (const G4Element *)
 

Protected Attributes

const G4ParticleDefinitionparticle
 
G4ParticleDefinitiontheGamma
 
G4ParticleChangeForLossfParticleChange
 
G4bool isElectron
 
- Protected Attributes inherited from G4VEmModel
G4VParticleChangepParticleChange
 
G4PhysicsTablexSectionTable
 
const std::vector< G4double > * theDensityFactor
 
const std::vector< G4int > * theDensityIdx
 

Detailed Description

Definition at line 64 of file G4eBremsstrahlungModel.hh.

Constructor & Destructor Documentation

◆ G4eBremsstrahlungModel()

G4eBremsstrahlungModel::G4eBremsstrahlungModel ( const G4ParticleDefinition p = 0,
const G4String nam = "eBrem" 
)

Definition at line 85 of file G4eBremsstrahlungModel.cc.

87 : G4VEmModel(nam),
88 particle(0),
89 isElectron(true),
90 probsup(1.0),
91 MigdalConstant(classic_electr_radius*electron_Compton_length*electron_Compton_length*4.0*pi),
92 LPMconstant(fine_structure_const*electron_mass_c2*electron_mass_c2/(4.*pi*hbarc)),
93 isInitialised(false)
94{
95 if(p) { SetParticle(p); }
98 highKinEnergy = HighEnergyLimit();
99 lowKinEnergy = LowEnergyLimit();
100 fParticleChange = 0;
101}
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:529
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:522
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:515
G4ParticleChangeForLoss * fParticleChange
G4ParticleDefinition * theGamma
const G4ParticleDefinition * particle

◆ ~G4eBremsstrahlungModel()

G4eBremsstrahlungModel::~G4eBremsstrahlungModel ( )
virtual

Definition at line 105 of file G4eBremsstrahlungModel.cc.

106{
107 size_t n = partialSumSigma.size();
108 if(n > 0) {
109 for(size_t i=0; i<n; i++) {
110 delete partialSumSigma[i];
111 }
112 }
113}

Member Function Documentation

◆ ComputeCrossSectionPerAtom()

G4double G4eBremsstrahlungModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition ,
G4double  tkin,
G4double  Z,
G4double  ,
G4double  cut,
G4double  maxE = DBL_MAX 
)
virtual

Reimplemented from G4VEmModel.

Definition at line 484 of file G4eBremsstrahlungModel.cc.

493{
494 G4double cross = 0.0 ;
495 if ( kineticEnergy < keV || kineticEnergy < cut) { return cross; }
496
497 static const G4double ksi=2.0, alfa=1.00;
498 static const G4double csigh = 0.127, csiglow = 0.25, asiglow = 0.020*MeV ;
499 static const G4double Tlim = 10.*MeV ;
500
501 static const G4double xlim = 1.2 ;
502 static const G4int NZ = 8 ;
503 static const G4int Nsig = 11 ;
504 static const G4double ZZ[NZ] =
505 {2.,4.,6.,14.,26.,50.,82.,92.} ;
506 static const G4double coefsig[NZ][Nsig] = {
507 // Z=2
508 { 0.4638, 0.37748, 0.32249, -0.060362, -0.065004,
509 -0.033457, -0.004583, 0.011954, 0.0030404, -0.0010077,
510 -0.00028131},
511
512 // Z=4
513 { 0.50008, 0.33483, 0.34364, -0.086262, -0.055361,
514 -0.028168, -0.0056172, 0.011129, 0.0027528, -0.00092265,
515 -0.00024348},
516
517 // Z=6
518 { 0.51587, 0.31095, 0.34996, -0.11623, -0.056167,
519 -0.0087154, 0.00053943, 0.0054092, 0.00077685, -0.00039635,
520 -6.7818e-05},
521
522 // Z=14
523 { 0.55058, 0.25629, 0.35854, -0.080656, -0.054308,
524 -0.049933, -0.00064246, 0.016597, 0.0021789, -0.001327,
525 -0.00025983},
526
527 // Z=26
528 { 0.5791, 0.26152, 0.38953, -0.17104, -0.099172,
529 0.024596, 0.023718, -0.0039205, -0.0036658, 0.00041749,
530 0.00023408},
531
532 // Z=50
533 { 0.62085, 0.27045, 0.39073, -0.37916, -0.18878,
534 0.23905, 0.095028, -0.068744, -0.023809, 0.0062408,
535 0.0020407},
536
537 // Z=82
538 { 0.66053, 0.24513, 0.35404, -0.47275, -0.22837,
539 0.35647, 0.13203, -0.1049, -0.034851, 0.0095046,
540 0.0030535},
541
542 // Z=92
543 { 0.67143, 0.23079, 0.32256, -0.46248, -0.20013,
544 0.3506, 0.11779, -0.1024, -0.032013, 0.0092279,
545 0.0028592}
546
547 } ;
548
549 G4int iz = 0 ;
550 G4double delz = 1.e6 ;
551 for (G4int ii=0; ii<NZ; ii++)
552 {
553 G4double absdelz = std::abs(Z-ZZ[ii]);
554 if(absdelz < delz)
555 {
556 iz = ii ;
557 delz = absdelz;
558 }
559 }
560
561 G4double xx = log10(kineticEnergy/MeV);
562 G4double fs = 1. ;
563
564 if (xx <= xlim) {
565
566 fs = coefsig[iz][Nsig-1] ;
567 for (G4int j=Nsig-2; j>=0; j--) {
568
569 fs = fs*xx+coefsig[iz][j] ;
570 }
571 if(fs < 0.) { fs = 0.; }
572 }
573
574 cross = Z*(Z+ksi)*(1.-csigh*exp(log(Z)/4.))*pow(log(kineticEnergy/cut),alfa);
575
576 if (kineticEnergy <= Tlim)
577 cross *= exp(csiglow*log(Tlim/kineticEnergy))
578 *(1.+asiglow/(sqrt(Z)*kineticEnergy));
579
580 if (!isElectron)
581 cross *= PositronCorrFactorSigma(Z, kineticEnergy, cut);
582
583 cross *= fs/Avogadro ;
584
585 if (cross < 0.) { cross = 0.; }
586 return cross;
587}
double G4double
Definition: G4Types.hh:64
int G4int
Definition: G4Types.hh:66

Referenced by CrossSectionPerVolume().

◆ ComputeDEDXPerVolume()

G4double G4eBremsstrahlungModel::ComputeDEDXPerVolume ( const G4Material material,
const G4ParticleDefinition p,
G4double  kineticEnergy,
G4double  cutEnergy 
)
virtual

Reimplemented from G4VEmModel.

Definition at line 166 of file G4eBremsstrahlungModel.cc.

171{
172 if(!particle) { SetParticle(p); }
173 if(kineticEnergy < lowKinEnergy) { return 0.0; }
174
175 const G4double thigh = 100.*GeV;
176
177 G4double cut = std::min(cutEnergy, kineticEnergy);
178
179 G4double rate, loss;
180 const G4double factorHigh = 36./(1450.*GeV);
181 const G4double coef1 = -0.5;
182 const G4double coef2 = 2./9.;
183
184 const G4ElementVector* theElementVector = material->GetElementVector();
185 const G4double* theAtomicNumDensityVector = material->GetAtomicNumDensityVector();
186
187 G4double totalEnergy = kineticEnergy + electron_mass_c2;
188 G4double dedx = 0.0;
189
190 // loop for elements in the material
191 for (size_t i=0; i<material->GetNumberOfElements(); i++) {
192
193 G4double Z = (*theElementVector)[i]->GetZ();
194 G4double natom = theAtomicNumDensityVector[i];
195
196 // loss for MinKinEnergy<KineticEnergy<=100 GeV
197 if (kineticEnergy <= thigh) {
198
199 // x = log(totalEnergy/electron_mass_c2);
200 loss = ComputeBremLoss(Z, kineticEnergy, cut) ;
201 if (!isElectron) loss *= PositronCorrFactorLoss(Z, kineticEnergy, cut);
202
203 // extrapolation for KineticEnergy>100 GeV
204 } else {
205
206 // G4double xhigh = log(thigh/electron_mass_c2);
207 G4double cuthigh = thigh*0.5;
208
209 if (cut < thigh) {
210
211 loss = ComputeBremLoss(Z, thigh, cut) ;
212 if (!isElectron) loss *= PositronCorrFactorLoss(Z, thigh, cut) ;
213 rate = cut/totalEnergy;
214 loss *= (1. + coef1*rate + coef2*rate*rate);
215 rate = cut/thigh;
216 loss /= (1.+coef1*rate+coef2*rate*rate);
217
218 } else {
219
220 loss = ComputeBremLoss(Z, thigh, cuthigh) ;
221 if (!isElectron) loss *= PositronCorrFactorLoss(Z, thigh, cuthigh) ;
222 rate = cut/totalEnergy;
223 loss *= (1. + coef1*rate + coef2*rate*rate);
224 loss *= cut*factorHigh;
225 }
226 }
227 loss *= natom;
228
229 G4double kp2 = MigdalConstant*totalEnergy*totalEnergy
230 * (material->GetElectronDensity()) ;
231
232 // now compute the correction due to the supression(s)
233 G4double kmin = 1.*eV;
234 G4double kmax = cut;
235
236 if (kmax > kmin) {
237
238 G4double floss = 0.;
239 G4int nmax = 100;
240
241 G4double vmin=log(kmin);
242 G4double vmax=log(kmax) ;
243 G4int nn = (G4int)(nmax*(vmax-vmin)/(log(highKinEnergy)-vmin)) ;
244 G4double u,fac,c,v,dv ;
245 if(nn > 0) {
246
247 dv = (vmax-vmin)/nn ;
248 v = vmin-dv ;
249
250 for(G4int n=0; n<=nn; n++) {
251
252 v += dv;
253 u = exp(v);
254 fac = u*SupressionFunction(material,kineticEnergy,u);
255 fac *= probsup*(u*u/(u*u+kp2))+1.-probsup;
256 if ((n==0)||(n==nn)) c=0.5;
257 else c=1. ;
258 fac *= c ;
259 floss += fac ;
260 }
261 floss *=dv/(kmax-kmin);
262
263 } else {
264 floss = 1.;
265 }
266 if(floss > 1.) floss = 1.;
267 // correct the loss
268 loss *= floss;
269 }
270 dedx += loss;
271 }
272 if(dedx < 0.) { dedx = 0.; }
273 return dedx;
274}
std::vector< G4Element * > G4ElementVector
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:189
size_t GetNumberOfElements() const
Definition: G4Material.hh:185
const G4double * GetAtomicNumDensityVector() const
Definition: G4Material.hh:215
G4double GetElectronDensity() const
Definition: G4Material.hh:216

◆ CrossSectionPerVolume()

G4double G4eBremsstrahlungModel::CrossSectionPerVolume ( const G4Material material,
const G4ParticleDefinition p,
G4double  kineticEnergy,
G4double  cutEnergy,
G4double  maxEnergy 
)
virtual

Reimplemented from G4VEmModel.

Definition at line 406 of file G4eBremsstrahlungModel.cc.

412{
413 if(!particle) { SetParticle(p); }
414 G4double cross = 0.0;
415 G4double tmax = min(maxEnergy, kineticEnergy);
416 G4double cut = min(cutEnergy, kineticEnergy);
417 if(cut >= tmax) { return cross; }
418
419 const G4ElementVector* theElementVector = material->GetElementVector();
420 const G4double* theAtomNumDensityVector = material->GetAtomicNumDensityVector();
421 G4double dum=0.;
422
423 for (size_t i=0; i<material->GetNumberOfElements(); i++) {
424
425 cross += theAtomNumDensityVector[i] * ComputeCrossSectionPerAtom(p,
426 kineticEnergy, (*theElementVector)[i]->GetZ(), dum, cut);
427 if (tmax < kineticEnergy) {
428 cross -= theAtomNumDensityVector[i] * ComputeCrossSectionPerAtom(p,
429 kineticEnergy, (*theElementVector)[i]->GetZ(), dum, tmax);
430 }
431 }
432
433 // now compute the correction due to the supression(s)
434
435 G4double kmax = tmax;
436 G4double kmin = cut;
437
438 G4double totalEnergy = kineticEnergy+electron_mass_c2 ;
439 G4double kp2 = MigdalConstant*totalEnergy*totalEnergy
440 *(material->GetElectronDensity());
441
442 G4double fsig = 0.;
443 G4int nmax = 100;
444 G4double vmin=log(kmin);
445 G4double vmax=log(kmax) ;
446 G4int nn = (G4int)(nmax*(vmax-vmin)/(log(highKinEnergy)-vmin));
447 G4double u,fac,c,v,dv,y ;
448 if(nn > 0) {
449
450 dv = (vmax-vmin)/nn ;
451 v = vmin-dv ;
452 for(G4int n=0; n<=nn; n++) {
453
454 v += dv;
455 u = exp(v);
456 fac = SupressionFunction(material, kineticEnergy, u);
457 y = u/kmax;
458 fac *= (4.-4.*y+3.*y*y)/3.;
459 fac *= probsup*(u*u/(u*u+kp2))+1.-probsup;
460
461 if ((n==0)||(n==nn)) c=0.5;
462 else c=1. ;
463
464 fac *= c;
465 fsig += fac;
466 }
467 y = kmin/kmax ;
468 fsig *=dv/(-4.*log(y)/3.-4.*(1.-y)/3.+0.5*(1.-y*y));
469
470 } else {
471
472 fsig = 1.;
473 }
474 if (fsig > 1.) fsig = 1.;
475
476 // correct the cross section
477 cross *= fsig;
478
479 return cross;
480}
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double tkin, G4double Z, G4double, G4double cut, G4double maxE=DBL_MAX)

◆ Initialise()

void G4eBremsstrahlungModel::Initialise ( const G4ParticleDefinition p,
const G4DataVector cuts 
)
virtual

Implements G4VEmModel.

Reimplemented in G4ePolarizedBremsstrahlungModel.

Definition at line 126 of file G4eBremsstrahlungModel.cc.

128{
129 if(p) { SetParticle(p); }
130 highKinEnergy = HighEnergyLimit();
131 lowKinEnergy = LowEnergyLimit();
132 const G4ProductionCutsTable* theCoupleTable=
134
135 if(theCoupleTable) {
136 G4int numOfCouples = theCoupleTable->GetTableSize();
137
138 G4int nn = partialSumSigma.size();
139 G4int nc = cuts.size();
140 if(nn > 0) {
141 for (G4int ii=0; ii<nn; ii++){
142 G4DataVector* a=partialSumSigma[ii];
143 if ( a ) delete a;
144 }
145 partialSumSigma.clear();
146 }
147 if(numOfCouples>0) {
148 for (G4int i=0; i<numOfCouples; i++) {
149 G4double cute = DBL_MAX;
150 if(i < nc) cute = cuts[i];
151 const G4MaterialCutsCouple* couple = theCoupleTable->GetMaterialCutsCouple(i);
152 const G4Material* material = couple->GetMaterial();
153 G4DataVector* dv = ComputePartialSumSigma(material, 0.5*highKinEnergy,
154 std::min(cute, 0.25*highKinEnergy));
155 partialSumSigma.push_back(dv);
156 }
157 }
158 }
159 if(isInitialised) return;
161 isInitialised = true;
162}
const G4Material * GetMaterial() const
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
static G4ProductionCutsTable * GetProductionCutsTable()
G4ParticleChangeForLoss * GetParticleChangeForLoss()
Definition: G4VEmModel.cc:95
#define DBL_MAX
Definition: templates.hh:83

Referenced by G4ePolarizedBremsstrahlungModel::Initialise().

◆ SampleSecondaries()

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

Implements G4VEmModel.

Reimplemented in G4ePolarizedBremsstrahlungModel.

Definition at line 643 of file G4eBremsstrahlungModel.cc.

661{
662 G4double kineticEnergy = dp->GetKineticEnergy();
663 G4double tmax = min(maxEnergy, kineticEnergy);
664 if(tmin >= tmax) { return; }
665
666//
667// GEANT4 internal units.
668//
669 static const G4double
670 ah10 = 4.67733E+00, ah11 =-6.19012E-01, ah12 = 2.02225E-02,
671 ah20 =-7.34101E+00, ah21 = 1.00462E+00, ah22 =-3.20985E-02,
672 ah30 = 2.93119E+00, ah31 =-4.03761E-01, ah32 = 1.25153E-02;
673
674 static const G4double
675 bh10 = 4.23071E+00, bh11 =-6.10995E-01, bh12 = 1.95531E-02,
676 bh20 =-7.12527E+00, bh21 = 9.69160E-01, bh22 =-2.74255E-02,
677 bh30 = 2.69925E+00, bh31 =-3.63283E-01, bh32 = 9.55316E-03;
678
679 static const G4double
680 al00 =-2.05398E+00, al01 = 2.38815E-02, al02 = 5.25483E-04,
681 al10 =-7.69748E-02, al11 =-6.91499E-02, al12 = 2.22453E-03,
682 al20 = 4.06463E-02, al21 =-1.01281E-02, al22 = 3.40919E-04;
683
684 static const G4double
685 bl00 = 1.04133E+00, bl01 =-9.43291E-03, bl02 =-4.54758E-04,
686 bl10 = 1.19253E-01, bl11 = 4.07467E-02, bl12 =-1.30718E-03,
687 bl20 =-1.59391E-02, bl21 = 7.27752E-03, bl22 =-1.94405E-04;
688
689 static const G4double tlow = 1.*MeV;
690
691 G4double gammaEnergy;
692 G4bool LPMOK = false;
693 const G4Material* material = couple->GetMaterial();
694
695 // select randomly one element constituing the material
696 const G4Element* anElement = SelectRandomAtom(couple);
697
698 // Extract Z factors for this Element
699 G4double lnZ = 3.*(anElement->GetIonisation()->GetlogZ3());
700 G4double FZ = lnZ* (4.- 0.55*lnZ);
701 G4double ZZ = anElement->GetIonisation()->GetZZ3();
702
703 // limits of the energy sampling
704 G4double totalEnergy = kineticEnergy + electron_mass_c2;
705
706 G4double xmin = tmin/kineticEnergy;
707 G4double xmax = tmax/kineticEnergy;
708 G4double kappa = 0.0;
709 if(xmax >= 1.) { xmax = 1.; }
710 else { kappa = log(xmax)/log(xmin); }
711 G4double epsilmin = tmin/totalEnergy;
712 G4double epsilmax = tmax/totalEnergy;
713
714 // Migdal factor
715 G4double MigdalFactor = (material->GetElectronDensity())*MigdalConstant
716 / (epsilmax*epsilmax);
717
718 G4double x, epsil, greject, migdal, grejmax, q;
719 G4double U = log(kineticEnergy/electron_mass_c2);
720 G4double U2 = U*U;
721
722 // precalculated parameters
723 G4double ah, bh;
724 G4double screenfac = 0.0;
725
726 if (kineticEnergy > tlow) {
727
728 G4double ah1 = ah10 + ZZ* (ah11 + ZZ* ah12);
729 G4double ah2 = ah20 + ZZ* (ah21 + ZZ* ah22);
730 G4double ah3 = ah30 + ZZ* (ah31 + ZZ* ah32);
731
732 G4double bh1 = bh10 + ZZ* (bh11 + ZZ* bh12);
733 G4double bh2 = bh20 + ZZ* (bh21 + ZZ* bh22);
734 G4double bh3 = bh30 + ZZ* (bh31 + ZZ* bh32);
735
736 ah = 1. + (ah1*U2 + ah2*U + ah3) / (U2*U);
737 bh = 0.75 + (bh1*U2 + bh2*U + bh3) / (U2*U);
738
739 // limit of the screening variable
740 screenfac =
741 136.*electron_mass_c2/((anElement->GetIonisation()->GetZ3())*totalEnergy);
742 G4double screenmin = screenfac*epsilmin/(1.-epsilmin);
743
744 // Compute the maximum of the rejection function
745 G4double F1 = max(ScreenFunction1(screenmin) - FZ ,0.);
746 G4double F2 = max(ScreenFunction2(screenmin) - FZ ,0.);
747 grejmax = (F1 - epsilmin* (F1*ah - bh*epsilmin*F2))/(42.392 - FZ);
748
749 } else {
750
751 G4double al0 = al00 + ZZ* (al01 + ZZ* al02);
752 G4double al1 = al10 + ZZ* (al11 + ZZ* al12);
753 G4double al2 = al20 + ZZ* (al21 + ZZ* al22);
754
755 G4double bl0 = bl00 + ZZ* (bl01 + ZZ* bl02);
756 G4double bl1 = bl10 + ZZ* (bl11 + ZZ* bl12);
757 G4double bl2 = bl20 + ZZ* (bl21 + ZZ* bl22);
758
759 ah = al0 + al1*U + al2*U2;
760 bh = bl0 + bl1*U + bl2*U2;
761
762 // Compute the maximum of the rejection function
763 grejmax = max(1. + xmin* (ah + bh*xmin), 1.+ah+bh);
764 G4double xm = -ah/(2.*bh);
765 if ( xmin < xm && xm < xmax) grejmax = max(grejmax, 1.+ xm* (ah + bh*xm));
766 }
767
768 //
769 // sample the energy rate of the emitted gamma for e- kin energy > 1 MeV
770 //
771
772 do {
773 if (kineticEnergy > tlow) {
774 do {
775 q = G4UniformRand();
776 x = pow(xmin, q + kappa*(1.0 - q));
777 epsil = x*kineticEnergy/totalEnergy;
778 G4double screenvar = screenfac*epsil/(1.0-epsil);
779 G4double F1 = max(ScreenFunction1(screenvar) - FZ ,0.);
780 G4double F2 = max(ScreenFunction2(screenvar) - FZ ,0.);
781 migdal = (1. + MigdalFactor)/(1. + MigdalFactor/(x*x));
782 greject = migdal*(F1 - epsil* (ah*F1 - bh*epsil*F2))/(42.392 - FZ);
783 /*
784 if ( greject > grejmax ) {
785 G4cout << "### G4eBremsstrahlungModel Warning: Majoranta exceeded! "
786 << greject << " > " << grejmax
787 << " x= " << x
788 << " e= " << kineticEnergy
789 << G4endl;
790 }
791 */
792 } while( greject < G4UniformRand()*grejmax );
793
794 } else {
795
796 do {
797 q = G4UniformRand();
798 x = pow(xmin, q + kappa*(1.0 - q));
799 migdal = (1. + MigdalFactor)/(1. + MigdalFactor/(x*x));
800 greject = migdal*(1. + x* (ah + bh*x));
801 /*
802 if ( greject > grejmax ) {
803 G4cout << "### G4eBremsstrahlungModel Warning: Majoranta exceeded! "
804 << greject << " > " << grejmax
805 << " x= " << x
806 << " e= " << kineticEnergy
807 << G4endl;
808 }
809 */
810 } while( greject < G4UniformRand()*grejmax );
811 }
812 gammaEnergy = x*kineticEnergy;
813
814 if (LPMFlag()) {
815 // take into account the supression due to the LPM effect
816 if (G4UniformRand() <= SupressionFunction(material,kineticEnergy,
817 gammaEnergy))
818 LPMOK = true;
819 }
820 else LPMOK = true;
821
822 } while (!LPMOK);
823
824 //
825 // angles of the emitted gamma. ( Z - axis along the parent particle)
826 // use general interface
827 //
828
829 G4ThreeVector gammaDirection =
830 GetAngularDistribution()->SampleDirection(dp, totalEnergy-gammaEnergy,
831 G4lrint(anElement->GetZ()),
832 couple->GetMaterial());
833
834 // create G4DynamicParticle object for the Gamma
835 G4DynamicParticle* gamma = new G4DynamicParticle(theGamma,gammaDirection,
836 gammaEnergy);
837 vdp->push_back(gamma);
838
839 G4double totMomentum = sqrt(kineticEnergy*(totalEnergy + electron_mass_c2));
840
841 G4ThreeVector direction = (totMomentum*dp->GetMomentumDirection()
842 - gammaEnergy*gammaDirection).unit();
843
844 // energy of primary
845 G4double finalE = kineticEnergy - gammaEnergy;
846
847 // stop tracking and create new secondary instead of primary
848 if(gammaEnergy > SecondaryThreshold()) {
851 G4DynamicParticle* el =
853 direction, finalE);
854 vdp->push_back(el);
855
856 // continue tracking
857 } else {
860 }
861}
@ fStopAndKill
bool G4bool
Definition: G4Types.hh:67
#define G4UniformRand()
Definition: Randomize.hh:53
const G4ThreeVector & GetMomentumDirection() const
G4double GetKineticEnergy() const
G4double GetZ() const
Definition: G4Element.hh:131
G4IonisParamElm * GetIonisation() const
Definition: G4Element.hh:209
G4double GetlogZ3() const
G4double GetZ3() const
G4double GetZZ3() const
void SetProposedKineticEnergy(G4double proposedKinEnergy)
void SetProposedMomentumDirection(const G4ThreeVector &dir)
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:508
G4bool LPMFlag() const
Definition: G4VEmModel.hh:564
G4double SecondaryThreshold() const
Definition: G4VEmModel.hh:557
void ProposeTrackStatus(G4TrackStatus status)
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *couple)
int G4lrint(double ad)
Definition: templates.hh:163

Referenced by G4ePolarizedBremsstrahlungModel::SampleSecondaries().

◆ SelectRandomAtom()

const G4Element * G4eBremsstrahlungModel::SelectRandomAtom ( const G4MaterialCutsCouple couple)
protected

Definition at line 865 of file G4eBremsstrahlungModel.cc.

867{
868 // select randomly 1 element within the material
869
870 const G4Material* material = couple->GetMaterial();
871 G4int nElements = material->GetNumberOfElements();
872 const G4ElementVector* theElementVector = material->GetElementVector();
873
874 const G4Element* elm = 0;
875
876 if(1 < nElements) {
877
878 --nElements;
879 G4DataVector* dv = partialSumSigma[couple->GetIndex()];
880 G4double rval = G4UniformRand()*((*dv)[nElements]);
881
882 elm = (*theElementVector)[nElements];
883 for (G4int i=0; i<nElements; ++i) {
884 if (rval <= (*dv)[i]) {
885 elm = (*theElementVector)[i];
886 break;
887 }
888 }
889 } else { elm = (*theElementVector)[0]; }
890
892 return elm;
893}
void SetCurrentElement(const G4Element *)
Definition: G4VEmModel.hh:384

Referenced by SampleSecondaries().

Member Data Documentation

◆ fParticleChange

G4ParticleChangeForLoss* G4eBremsstrahlungModel::fParticleChange
protected

◆ isElectron

G4bool G4eBremsstrahlungModel::isElectron
protected

Definition at line 133 of file G4eBremsstrahlungModel.hh.

Referenced by ComputeCrossSectionPerAtom(), and ComputeDEDXPerVolume().

◆ particle

const G4ParticleDefinition* G4eBremsstrahlungModel::particle
protected

◆ theGamma

G4ParticleDefinition* G4eBremsstrahlungModel::theGamma
protected

Definition at line 130 of file G4eBremsstrahlungModel.hh.

Referenced by G4eBremsstrahlungModel(), and SampleSecondaries().


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