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

#include <G4BetheBlochModel.hh>

+ Inheritance diagram for G4BetheBlochModel:

Public Member Functions

 G4BetheBlochModel (const G4ParticleDefinition *p=nullptr, const G4String &nam="BetheBloch")
 
 ~G4BetheBlochModel () override
 
void Initialise (const G4ParticleDefinition *, const G4DataVector &) override
 
G4double MinEnergyCut (const G4ParticleDefinition *, const G4MaterialCutsCouple *couple) override
 
virtual G4double ComputeCrossSectionPerElectron (const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
 
G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kineticEnergy, G4double Z, G4double A, G4double cutEnergy, G4double maxEnergy) override
 
G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy) override
 
G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy) override
 
G4double GetChargeSquareRatio (const G4ParticleDefinition *p, const G4Material *mat, G4double kineticEnergy) override
 
G4double GetParticleCharge (const G4ParticleDefinition *p, const G4Material *mat, G4double kineticEnergy) override
 
void CorrectionsAlongStep (const G4MaterialCutsCouple *couple, const G4DynamicParticle *dp, const G4double &length, G4double &eloss) override
 
void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy) override
 
G4BetheBlochModeloperator= (const G4BetheBlochModel &right)=delete
 
 G4BetheBlochModel (const G4BetheBlochModel &)=delete
 
- 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 FillNumberOfSecondaries (G4int &numberOfTriplets, G4int &numberOfRecoil)
 
virtual void ModelDescription (std::ostream &outFile) const
 
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)
 
const G4ElementSelectRandomAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementSelectTargetAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double logKineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementSelectRandomAtom (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementGetCurrentElement (const G4Material *mat=nullptr) const
 
G4int SelectRandomAtomNumber (const G4Material *) const
 
const G4IsotopeGetCurrentIsotope (const G4Element *elm=nullptr) const
 
G4int SelectIsotopeNumber (const G4Element *) const
 
void SetParticleChange (G4VParticleChange *, G4VEmFluctuationModel *f=nullptr)
 
void SetCrossSectionTable (G4PhysicsTable *, G4bool isLocal)
 
G4ElementDataGetElementData ()
 
G4PhysicsTableGetCrossSectionTable ()
 
G4VEmFluctuationModelGetModelOfFluctuations ()
 
G4VEmAngularDistributionGetAngularDistribution ()
 
G4VEmModelGetTripletModel ()
 
void SetTripletModel (G4VEmModel *)
 
void SetAngularDistribution (G4VEmAngularDistribution *)
 
G4double HighEnergyLimit () const
 
G4double LowEnergyLimit () const
 
G4double HighEnergyActivationLimit () const
 
G4double LowEnergyActivationLimit () const
 
G4double PolarAngleLimit () const
 
G4double SecondaryThreshold () 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) const
 
void SetPolarAngleLimit (G4double)
 
void SetSecondaryThreshold (G4double)
 
void SetDeexcitationFlag (G4bool val)
 
void SetForceBuildTable (G4bool val)
 
void SetFluctuationFlag (G4bool val)
 
void SetMasterThread (G4bool val)
 
G4bool IsMaster () const
 
void SetUseBaseMaterials (G4bool val)
 
G4bool UseBaseMaterials () const
 
G4double MaxSecondaryKinEnergy (const G4DynamicParticle *dynParticle)
 
const G4StringGetName () const
 
void SetCurrentCouple (const G4MaterialCutsCouple *)
 
G4bool IsLocked () const
 
void SetLocked (G4bool)
 
void SetLPMFlag (G4bool)
 
G4VEmModeloperator= (const G4VEmModel &right)=delete
 
 G4VEmModel (const G4VEmModel &)=delete
 

Protected Member Functions

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 = nullptr
 
G4VParticleChangepParticleChange = nullptr
 
G4PhysicsTablexSectionTable = nullptr
 
const G4MaterialpBaseMaterial = nullptr
 
const std::vector< G4double > * theDensityFactor = nullptr
 
const std::vector< G4int > * theDensityIdx = nullptr
 
G4double inveplus
 
G4double pFactor = 1.0
 
std::size_t currentCoupleIndex = 0
 
std::size_t basedCoupleIndex = 0
 
G4bool lossFlucFlag = true
 

Detailed Description

Definition at line 60 of file G4BetheBlochModel.hh.

Constructor & Destructor Documentation

◆ G4BetheBlochModel() [1/2]

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

Definition at line 74 of file G4BetheBlochModel.cc.

76 : G4VEmModel(nam),
77 twoln10(2.0*G4Log(10.0)),
78 fAlphaTlimit(1*CLHEP::GeV),
79 fProtonTlimit(10*CLHEP::GeV)
80{
81 theElectron = G4Electron::Electron();
84 SetLowEnergyLimit(2.0*CLHEP::MeV);
85}
G4double G4Log(G4double x)
Definition G4Log.hh:227
static G4Electron * Electron()
Definition G4Electron.cc:91
static G4LossTableManager * Instance()
G4EmCorrections * EmCorrections()
static G4NistManager * Instance()
void SetLowEnergyLimit(G4double)
G4VEmModel(const G4String &nam)
Definition G4VEmModel.cc:67

◆ ~G4BetheBlochModel()

G4BetheBlochModel::~G4BetheBlochModel ( )
overridedefault

◆ G4BetheBlochModel() [2/2]

G4BetheBlochModel::G4BetheBlochModel ( const G4BetheBlochModel & )
delete

Member Function Documentation

◆ ComputeCrossSectionPerAtom()

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{
223 return Z*ComputeCrossSectionPerElectron(p,kinEnergy,cutEnergy,maxEnergy);
224}
virtual G4double ComputeCrossSectionPerElectron(const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)

◆ ComputeCrossSectionPerElectron()

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

Definition at line 184 of file G4BetheBlochModel.cc.

188{
189 G4double cross = 0.0;
190 const G4double tmax = MaxSecondaryEnergy(p, kineticEnergy);
191 const G4double cutEnergy = std::min(std::min(cut,tmax), tlimit);
192 const G4double maxEnergy = std::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 *= CLHEP::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}
double G4double
Definition G4Types.hh:83
G4double MaxSecondaryEnergy(const G4ParticleDefinition *, G4double kinEnergy) override

Referenced by ComputeCrossSectionPerAtom(), G4LindhardSorensenIonModel::ComputeCrossSectionPerElectron(), and CrossSectionPerVolume().

◆ ComputeDEDXPerVolume()

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 const G4double tmax = MaxSecondaryEnergy(p, kineticEnergy);
251 // projectile formfactor limit energy loss
252 const G4double cutEnergy = std::min(std::min(cut,tmax), tlimit);
253
254 G4double tau = kineticEnergy/mass;
255 G4double gam = tau + 1.0;
256 G4double bg2 = tau * (tau+2.0);
257 G4double beta2 = bg2/(gam*gam);
258 G4double xc = cutEnergy/tmax;
259
260 G4double eexc = material->GetIonisation()->GetMeanExcitationEnergy();
261 G4double eexc2 = eexc*eexc;
262
263 G4double eDensity = material->GetElectronDensity();
264
265 // added ICRU90 stopping data for limited list of materials
266 /*
267 G4cout << "### DEDX ICRI90:" << (nullptr != fICRU90)
268 << " Ekin=" << kineticEnergy
269 << " " << p->GetParticleName()
270 << " q2=" << chargeSquare
271 << " inside " << material->GetName() << G4endl;
272 */
273 if(nullptr != fICRU90 && kineticEnergy < fProtonTlimit) {
274 if(material != currentMaterial) {
275 currentMaterial = material;
276 baseMaterial = material->GetBaseMaterial()
277 ? material->GetBaseMaterial() : material;
278 iICRU90 = fICRU90->GetIndex(baseMaterial);
279 }
280 if(iICRU90 >= 0) {
281 G4double dedx = 0.0;
282 // only for alpha
283 if(isAlpha) {
284 if(kineticEnergy <= fAlphaTlimit) {
285 dedx = fICRU90->GetElectronicDEDXforAlpha(iICRU90, kineticEnergy);
286 } else {
287 const G4double e = kineticEnergy*CLHEP::proton_mass_c2/mass;
288 dedx = fICRU90->GetElectronicDEDXforProton(iICRU90, e)*chargeSquare;
289 }
290 } else {
291 dedx = fICRU90->GetElectronicDEDXforProton(iICRU90, kineticEnergy)
292 *chargeSquare;
293 }
294 dedx *= material->GetDensity();
295 if(cutEnergy < tmax) {
296 dedx += (G4Log(xc) + (1.0 - xc)*beta2)*CLHEP::twopi_mc2_rcl2
297 *(eDensity*chargeSquare/beta2);
298 }
299 //G4cout << " iICRU90=" << iICRU90 << " dedx=" << dedx << G4endl;
300 if(dedx > 0.0) { return dedx; }
301 }
302 }
303 // general Bethe-Bloch formula
304 G4double dedx = G4Log(2.0*CLHEP::electron_mass_c2*bg2*cutEnergy/eexc2)
305 - (1.0 + xc)*beta2;
306
307 if(0.0 < spin) {
308 G4double del = 0.5*cutEnergy/(kineticEnergy + mass);
309 dedx += del*del;
310 }
311
312 // density correction
313 G4double x = G4Log(bg2)/twoln10;
314 dedx -= material->GetIonisation()->DensityCorrection(x);
315
316 // shell correction
317 dedx -= 2.0*corr->ShellCorrection(p,material,kineticEnergy);
318
319 // now compute the total ionization loss
320 dedx *= CLHEP::twopi_mc2_rcl2*chargeSquare*eDensity/beta2;
321
322 //High order correction different for hadrons and ions
323 if(isIon) {
324 dedx += corr->IonBarkasCorrection(p,material,kineticEnergy);
325 } else {
326 dedx += corr->HighOrderCorrections(p,material,kineticEnergy,cutEnergy);
327 }
328
329 dedx = std::max(dedx, 0.0);
330 /*
331 G4cout << "E(MeV)= " << kineticEnergy/CLHEP::MeV << " dedx= " << dedx
332 << " " << material->GetName() << G4endl;
333 */
334 return dedx;
335}
G4double HighOrderCorrections(const G4ParticleDefinition *, const G4Material *, const G4double kineticEnergy, const G4double cutEnergy)
G4double ShellCorrection(const G4ParticleDefinition *, const G4Material *, const G4double kineticEnergy)
G4double IonBarkasCorrection(const G4ParticleDefinition *, const G4Material *, const G4double kineticEnergy)
G4double GetElectronicDEDXforProton(const G4Material *, G4double kinEnergy) const
G4int GetIndex(const G4Material *) const
G4double GetElectronicDEDXforAlpha(const G4Material *, G4double scaledKinEnergy) const
G4double DensityCorrection(G4double x) const
G4double GetMeanExcitationEnergy() const
G4double GetDensity() const
const G4Material * GetBaseMaterial() const
G4IonisParamMat * GetIonisation() const
G4double GetElectronDensity() const

Referenced by G4BetheBlochNoDeltaModel::ComputeDEDXPerVolume(), and G4LindhardSorensenIonModel::ComputeDEDXPerVolume().

◆ CorrectionsAlongStep()

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

Reimplemented from G4VEmModel.

Definition at line 339 of file G4BetheBlochModel.cc.

343{
344 // no correction for alpha
345 if(isAlpha) { return; }
346
347 // no correction at the last step or at small step
348 const G4double preKinEnergy = dp->GetKineticEnergy();
349 if(eloss >= preKinEnergy || eloss < preKinEnergy*0.05) { return; }
350
351 // corrections for all charged particles with Q > 1
352 const G4ParticleDefinition* p = dp->GetDefinition();
353 if(p != particle) { SetupParameters(p); }
354 if(!isIon) { return; }
355
356 // effective energy and charge at a step
357 const G4double e = std::max(preKinEnergy - eloss*0.5, preKinEnergy*0.5);
358 const G4Material* mat = couple->GetMaterial();
359 const G4double q20 = corr->EffectiveChargeSquareRatio(p, mat, preKinEnergy);
360 const G4double q2 = corr->EffectiveChargeSquareRatio(p, mat, e);
361 const G4double qfactor = q2/q20;
362
363 /*
364 G4cout << "G4BetheBlochModel::CorrectionsAlongStep: Epre(MeV)="
365 << preKinEnergy << " Eeff(MeV)=" << e
366 << " eloss=" << eloss << " elossnew=" << eloss*qfactor
367 << " qfactor=" << qfactor << " Qpre=" << q20
368 << p->GetParticleName() <<G4endl;
369 */
370 eloss *= qfactor;
371}
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
G4double EffectiveChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, const G4double kineticEnergy)
const G4Material * GetMaterial() const

◆ CrossSectionPerVolume()

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

Reimplemented from G4VEmModel.

Reimplemented in G4BetheBlochNoDeltaModel.

Definition at line 228 of file G4BetheBlochModel.cc.

234{
235 G4double sigma = mat->GetElectronDensity()
236 *ComputeCrossSectionPerElectron(p,kinEnergy,cutEnergy,maxEnergy);
237 if(isAlpha) {
238 sigma *= corr->EffectiveChargeSquareRatio(p,mat,kinEnergy)/chargeSquare;
239 }
240 return sigma;
241}

◆ GetChargeSquareRatio() [1/2]

G4double G4BetheBlochModel::GetChargeSquareRatio ( ) const
inlineprotected

Definition at line 162 of file G4BetheBlochModel.hh.

163{
164 return chargeSquare;
165}

◆ GetChargeSquareRatio() [2/2]

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

Reimplemented from G4VEmModel.

Definition at line 125 of file G4BetheBlochModel.cc.

128{
129 // this method is called only for ions, so no check if it is an ion
130 if(isAlpha) { return 1.0; }
131 chargeSquare = corr->EffectiveChargeSquareRatio(p, mat, kinEnergy);
132 return chargeSquare;
133}

◆ GetParticleCharge()

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

Reimplemented from G4VEmModel.

Definition at line 137 of file G4BetheBlochModel.cc.

140{
141 // this method is called only for ions, so no check if it is an ion
142 return corr->GetParticleCharge(p, mat, kineticEnergy);
143}
G4double GetParticleCharge(const G4ParticleDefinition *, const G4Material *, const G4double kineticEnergy)

◆ Initialise()

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

Implements G4VEmModel.

Definition at line 93 of file G4BetheBlochModel.cc.

95{
96 if(p != particle) { SetupParameters(p); }
97
98 // always false before the run
100
101 // initialisation once
102 if(nullptr == fParticleChange) {
103 const G4String& pname = particle->GetParticleName();
104 if(G4EmParameters::Instance()->UseICRU90Data() &&
105 (pname == "proton" || pname == "GenericIon" || pname == "alpha")) {
106 fICRU90 = nist->GetICRU90StoppingData();
107 }
108 if(particle->GetPDGCharge() > CLHEP::eplus ||
109 pname == "GenericIon") { isIon = true; }
110 if(pname == "alpha") { isAlpha = true; }
111
112 fParticleChange = GetParticleChangeForLoss();
113 if(UseAngularGeneratorFlag() && nullptr == GetAngularDistribution()) {
115 }
116 }
117 // initialisation for each new run
118 if(IsMaster() && nullptr != fICRU90) {
119 fICRU90->Initialise();
120 }
121}
static G4EmParameters * Instance()
G4ICRU90StoppingData * GetICRU90StoppingData()
const G4String & GetParticleName() const
G4VEmAngularDistribution * GetAngularDistribution()
G4bool IsMaster() const
void SetDeexcitationFlag(G4bool val)
void SetAngularDistribution(G4VEmAngularDistribution *)
G4bool UseAngularGeneratorFlag() const
G4ParticleChangeForLoss * GetParticleChangeForLoss()

Referenced by G4LindhardSorensenIonModel::Initialise().

◆ MaxSecondaryEnergy()

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

Reimplemented from G4VEmModel.

Definition at line 487 of file G4BetheBlochModel.cc.

489{
490 // here particle type is checked for the case,
491 // when this model is shared between particles
492 if(pd != particle) { SetupParameters(pd); }
493 G4double tau = kinEnergy/mass;
494 return 2.0*CLHEP::electron_mass_c2*tau*(tau + 2.) /
495 (1. + 2.0*(tau + 1.)*ratio + ratio*ratio);
496}

Referenced by ComputeCrossSectionPerElectron(), ComputeDEDXPerVolume(), and SampleSecondaries().

◆ MinEnergyCut()

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

Reimplemented from G4VEmModel.

Definition at line 175 of file G4BetheBlochModel.cc.

177{
179}

◆ operator=()

G4BetheBlochModel & G4BetheBlochModel::operator= ( const G4BetheBlochModel & right)
delete

◆ SampleSecondaries()

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

Implements G4VEmModel.

Definition at line 375 of file G4BetheBlochModel.cc.

380{
381 G4double kinEnergy = dp->GetKineticEnergy();
382 const G4double tmax = MaxSecondaryEnergy(dp->GetDefinition(), kinEnergy);
383 const G4double minKinEnergy = std::min(cut, tmax);
384 const G4double maxKinEnergy = std::min(maxEnergy, tmax);
385 if(minKinEnergy >= maxKinEnergy) { return; }
386
387 //G4cout << "G4BetheBlochModel::SampleSecondaries Emin= " << minKinEnergy
388 // << " Emax= " << maxKinEnergy << G4endl;
389
390 const G4double totEnergy = kinEnergy + mass;
391 const G4double etot2 = totEnergy*totEnergy;
392 const G4double beta2 = kinEnergy*(kinEnergy + 2.0*mass)/etot2;
393
394 G4double deltaKinEnergy, f;
395 G4double f1 = 0.0;
396 G4double fmax = 1.0;
397 if( 0.0 < spin ) { fmax += 0.5*maxKinEnergy*maxKinEnergy/etot2; }
398
399 CLHEP::HepRandomEngine* rndmEngineMod = G4Random::getTheEngine();
400 G4double rndm[2];
401
402 // sampling without nuclear size effect
403 do {
404 rndmEngineMod->flatArray(2, rndm);
405 deltaKinEnergy = minKinEnergy*maxKinEnergy
406 /(minKinEnergy*(1.0 - rndm[0]) + maxKinEnergy*rndm[0]);
407
408 f = 1.0 - beta2*deltaKinEnergy/tmax;
409 if( 0.0 < spin ) {
410 f1 = 0.5*deltaKinEnergy*deltaKinEnergy/etot2;
411 f += f1;
412 }
413
414 // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
415 } while( fmax*rndm[1] > f);
416
417 // projectile formfactor - suppresion of high energy
418 // delta-electron production at high energy
419
420 G4double x = formfact*deltaKinEnergy;
421 if(x > 1.e-6) {
422
423 G4double x1 = 1.0 + x;
424 G4double grej = 1.0/(x1*x1);
425 if( 0.0 < spin ) {
426 G4double x2 = 0.5*electron_mass_c2*deltaKinEnergy/(mass*mass);
427 grej *= (1.0 + magMoment2*(x2 - f1/f)/(1.0 + x2));
428 }
429 if(grej > 1.1) {
430 G4cout << "### G4BetheBlochModel WARNING: grej= " << grej
431 << " " << dp->GetDefinition()->GetParticleName()
432 << " Ekin(MeV)= " << kinEnergy
433 << " delEkin(MeV)= " << deltaKinEnergy
434 << G4endl;
435 }
436 if(rndmEngineMod->flat() > grej) { return; }
437 }
438
439 G4ThreeVector deltaDirection;
440
442 const G4Material* mat = couple->GetMaterial();
443 deltaDirection =
444 GetAngularDistribution()->SampleDirection(dp, deltaKinEnergy,
446 mat);
447 } else {
448
449 G4double deltaMomentum =
450 std::sqrt(deltaKinEnergy * (deltaKinEnergy + 2.0*electron_mass_c2));
451 G4double cost = deltaKinEnergy * (totEnergy + electron_mass_c2) /
452 (deltaMomentum * dp->GetTotalMomentum());
453 cost = std::min(cost, 1.0);
454 const G4double sint = std::sqrt((1.0 - cost)*(1.0 + cost));
455 const G4double phi = twopi*rndmEngineMod->flat();
456
457 deltaDirection.set(sint*std::cos(phi),sint*std::sin(phi), cost) ;
458 deltaDirection.rotateUz(dp->GetMomentumDirection());
459 }
460 /*
461 G4cout << "### G4BetheBlochModel "
462 << dp->GetDefinition()->GetParticleName()
463 << " Ekin(MeV)= " << kinEnergy
464 << " delEkin(MeV)= " << deltaKinEnergy
465 << " tmin(MeV)= " << minKinEnergy
466 << " tmax(MeV)= " << maxKinEnergy
467 << " dir= " << dp->GetMomentumDirection()
468 << " dirDelta= " << deltaDirection
469 << G4endl;
470 */
471 // create G4DynamicParticle object for delta ray
472 auto delta = new G4DynamicParticle(theElectron,deltaDirection,deltaKinEnergy);
473
474 vdp->push_back(delta);
475
476 // Change kinematics of primary particle
477 kinEnergy -= deltaKinEnergy;
478 G4ThreeVector finalP = dp->GetMomentum() - delta->GetMomentum();
479 finalP = finalP.unit();
480
481 fParticleChange->SetProposedKineticEnergy(kinEnergy);
482 fParticleChange->SetProposedMomentumDirection(finalP);
483}
#define G4endl
Definition G4ios.hh:67
G4GLOB_DLL std::ostream G4cout
Hep3Vector unit() const
void set(double x, double y, double z)
Hep3Vector & rotateUz(const Hep3Vector &)
virtual double flat()=0
virtual void flatArray(const int size, double *vect)=0
const G4ThreeVector & GetMomentumDirection() const
G4ThreeVector GetMomentum() const
G4double GetTotalMomentum() const
void SetProposedKineticEnergy(G4double proposedKinEnergy)
void SetProposedMomentumDirection(const G4ThreeVector &dir)
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
G4int SelectRandomAtomNumber(const G4Material *) const

◆ SetChargeSquareRatio()

void G4BetheBlochModel::SetChargeSquareRatio ( G4double val)
inlineprotected

Definition at line 169 of file G4BetheBlochModel.hh.

170{
171 chargeSquare = val;
172}

Referenced by G4BetheBlochIonGasModel::ChargeSquareRatio().


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