Geant4-11
Public Member Functions | Protected Member Functions | Protected Attributes | Private Member Functions | Private Attributes | Static Private Attributes
G4LindhardSorensenIonModel Class Reference

#include <G4LindhardSorensenIonModel.hh>

Inheritance diagram for G4LindhardSorensenIonModel:
G4VEmModel

Public Member Functions

virtual G4double ChargeSquareRatio (const G4Track &)
 
G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, const G4Element *, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kineticEnergy, G4double Z, G4double A, G4double cutEnergy, G4double maxEnergy) override
 
G4double ComputeCrossSectionPerElectron (const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
 
virtual G4double ComputeCrossSectionPerShell (const G4ParticleDefinition *, G4int Z, G4int shellIdx, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeDEDX (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
 
G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy) override
 
G4double ComputeMeanFreePath (const G4ParticleDefinition *, G4double kineticEnergy, const G4Material *, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
void CorrectionsAlongStep (const G4MaterialCutsCouple *couple, const G4DynamicParticle *dp, const G4double &length, G4double &eloss) override
 
G4double CrossSection (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy) override
 
G4bool DeexcitationFlag () const
 
virtual void DefineForRegion (const G4Region *)
 
virtual void FillNumberOfSecondaries (G4int &numberOfTriplets, G4int &numberOfRecoil)
 
G4bool ForceBuildTableFlag () const
 
 G4LindhardSorensenIonModel (const G4LindhardSorensenIonModel &)=delete
 
 G4LindhardSorensenIonModel (const G4ParticleDefinition *p=nullptr, const G4String &nam="LindhardSorensen")
 
G4VEmAngularDistributionGetAngularDistribution ()
 
G4double GetChargeSquareRatio (const G4ParticleDefinition *p, const G4Material *mat, G4double kineticEnergy) override
 
G4PhysicsTableGetCrossSectionTable ()
 
const G4ElementGetCurrentElement () const
 
const G4IsotopeGetCurrentIsotope () const
 
G4ElementDataGetElementData ()
 
std::vector< G4EmElementSelector * > * GetElementSelectors ()
 
G4VEmFluctuationModelGetModelOfFluctuations ()
 
const G4StringGetName () const
 
virtual G4double GetPartialCrossSection (const G4Material *, G4int level, const G4ParticleDefinition *, G4double kineticEnergy)
 
G4double GetParticleCharge (const G4ParticleDefinition *p, const G4Material *mat, G4double kineticEnergy) override
 
G4VEmModelGetTripletModel ()
 
G4double HighEnergyActivationLimit () const
 
G4double HighEnergyLimit () const
 
void Initialise (const G4ParticleDefinition *, const G4DataVector &) override
 
void InitialiseElementSelectors (const G4ParticleDefinition *, const G4DataVector &)
 
virtual void InitialiseForElement (const G4ParticleDefinition *, G4int Z)
 
virtual void InitialiseForMaterial (const G4ParticleDefinition *, const G4Material *)
 
virtual void InitialiseLocal (const G4ParticleDefinition *, G4VEmModel *masterModel)
 
G4bool IsActive (G4double kinEnergy) const
 
G4bool IsLocked () const
 
G4bool IsMaster () const
 
G4double LowEnergyActivationLimit () const
 
G4double LowEnergyLimit () const
 
G4bool LPMFlag () const
 
G4double MaxSecondaryKinEnergy (const G4DynamicParticle *dynParticle)
 
G4double MinEnergyCut (const G4ParticleDefinition *, const G4MaterialCutsCouple *couple) override
 
virtual G4double MinPrimaryEnergy (const G4Material *, const G4ParticleDefinition *, G4double cut=0.0)
 
virtual void ModelDescription (std::ostream &outFile) const
 
G4LindhardSorensenIonModeloperator= (const G4LindhardSorensenIonModel &right)=delete
 
G4double PolarAngleLimit () const
 
void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy) override
 
G4double SecondaryThreshold () const
 
G4int SelectIsotopeNumber (const G4Element *)
 
const G4ElementSelectRandomAtom (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementSelectRandomAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4int SelectRandomAtomNumber (const G4Material *)
 
const G4ElementSelectTargetAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double logKineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
void SetActivationHighEnergyLimit (G4double)
 
void SetActivationLowEnergyLimit (G4double)
 
void SetAngularDistribution (G4VEmAngularDistribution *)
 
void SetAngularGeneratorFlag (G4bool)
 
void SetCrossSectionTable (G4PhysicsTable *, G4bool isLocal)
 
void SetCurrentCouple (const G4MaterialCutsCouple *)
 
void SetDeexcitationFlag (G4bool val)
 
void SetElementSelectors (std::vector< G4EmElementSelector * > *)
 
void SetFluctuationFlag (G4bool val)
 
void SetForceBuildTable (G4bool val)
 
void SetHighEnergyLimit (G4double)
 
void SetLocked (G4bool)
 
void SetLowEnergyLimit (G4double)
 
void SetLPMFlag (G4bool val)
 
void SetMasterThread (G4bool val)
 
void SetParticleChange (G4VParticleChange *, G4VEmFluctuationModel *f=nullptr)
 
void SetPolarAngleLimit (G4double)
 
void SetSecondaryThreshold (G4double)
 
void SetTripletModel (G4VEmModel *)
 
virtual void SetupForMaterial (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
void SetUseBaseMaterials (G4bool val)
 
virtual void StartTracking (G4Track *)
 
G4bool UseAngularGeneratorFlag () const
 
G4bool UseBaseMaterials () const
 
virtual G4double Value (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
 
 ~G4LindhardSorensenIonModel () override
 

Protected Member Functions

const G4MaterialCutsCoupleCurrentCouple () const
 
G4double GetChargeSquareRatio () const
 
G4ParticleChangeForGammaGetParticleChangeForGamma ()
 
G4ParticleChangeForLossGetParticleChangeForLoss ()
 
G4double MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kinEnergy) override
 
void SetChargeSquareRatio (G4double val)
 
void SetCurrentElement (const G4Element *)
 

Protected Attributes

size_t basedCoupleIndex = 0
 
size_t currentCoupleIndex = 0
 
G4ElementDatafElementData = nullptr
 
G4double inveplus
 
G4bool lossFlucFlag = true
 
const G4MaterialpBaseMaterial = nullptr
 
G4double pFactor = 1.0
 
G4VParticleChangepParticleChange = nullptr
 
const std::vector< G4double > * theDensityFactor = nullptr
 
const std::vector< G4int > * theDensityIdx = nullptr
 
G4PhysicsTablexSectionTable = nullptr
 

Private Member Functions

G4double ComputeDEDXPerVolumeLS (const G4Material *, const G4ParticleDefinition *, G4double kinEnergy, G4double cutEnergy)
 
void InitialiseLS ()
 
void SetParticle (const G4ParticleDefinition *p)
 
void SetupParameters ()
 

Private Attributes

G4VEmAngularDistributionanglModel = nullptr
 
G4double charge = 1.0
 
G4double chargeSquare = 1.0
 
G4EmCorrectionscorr
 
std::vector< G4EmElementSelector * > * elmSelectors = nullptr
 
G4double eMaxActive = DBL_MAX
 
G4double eMinActive = 0.0
 
G4double eRatio = 0.0
 
G4BetheBlochModelfBBModel
 
G4BraggIonModelfBraggModel
 
const G4MaterialCutsCouplefCurrentCouple = nullptr
 
const G4ElementfCurrentElement = nullptr
 
const G4IsotopefCurrentIsotope = nullptr
 
G4double fElimit
 
G4LossTableManagerfEmManager
 
G4bool flagDeexcitation = false
 
G4bool flagForceBuildTable = false
 
G4VEmFluctuationModelflucModel = nullptr
 
G4double formfact = 0.0
 
G4ParticleChangeForLossfParticleChange
 
G4VEmModelfTripletModel = nullptr
 
G4double highLimit
 
G4bool isLocked = false
 
G4bool isMaster = true
 
G4bool localElmSelectors = true
 
G4bool localTable = true
 
G4double lowLimit
 
G4double magMoment2 = 0.0
 
G4double mass = 0.0
 
const G4String name
 
G4NistManagernist
 
G4int nsec = 5
 
G4int nSelectors = 0
 
const G4ParticleDefinitionparticle
 
G4double polarAngleLimit
 
G4double pRatio = 1.0
 
G4double secondaryThreshold = DBL_MAX
 
G4double spin = 0.0
 
G4ParticleDefinitiontheElectron
 
G4bool theLPMflag = false
 
G4double tlimit = DBL_MAX
 
G4double twoln10
 
G4bool useAngularGenerator = false
 
G4bool useBaseMaterials = false
 
std::vector< G4doublexsec
 
G4int Zin = 1
 

Static Private Attributes

static std::vector< G4float > * fact [MAXZION] = {nullptr}
 
static G4IonICRU73DatafIonData = nullptr
 
static G4LindhardSorensenDatalsdata = nullptr
 
static const G4int MAXZION = 93
 

Detailed Description

Definition at line 63 of file G4LindhardSorensenIonModel.hh.

Constructor & Destructor Documentation

◆ G4LindhardSorensenIonModel() [1/2]

G4LindhardSorensenIonModel::G4LindhardSorensenIonModel ( const G4ParticleDefinition p = nullptr,
const G4String nam = "LindhardSorensen" 
)
explicit

Definition at line 67 of file G4LindhardSorensenIonModel.cc.

69 : G4VEmModel(nam),
70 particle(nullptr),
71 twoln10(2.0*G4Log(10.0))
72{
73 fParticleChange = nullptr;
79 fElimit = 2.0*CLHEP::MeV;
80}
G4double G4Log(G4double x)
Definition: G4Log.hh:226
static G4Electron * Electron()
Definition: G4Electron.cc:93
G4ParticleChangeForLoss * fParticleChange
const G4ParticleDefinition * particle
static G4LossTableManager * Instance()
G4EmCorrections * EmCorrections()
static G4NistManager * Instance()
G4VEmModel(const G4String &nam)
Definition: G4VEmModel.cc:66
static constexpr double MeV

References corr, G4Electron::Electron(), G4LossTableManager::EmCorrections(), fBBModel, fBraggModel, fElimit, fParticleChange, G4NistManager::Instance(), G4LossTableManager::Instance(), CLHEP::MeV, nist, and theElectron.

◆ ~G4LindhardSorensenIonModel()

G4LindhardSorensenIonModel::~G4LindhardSorensenIonModel ( )
override

Definition at line 84 of file G4LindhardSorensenIonModel.cc.

85{}

◆ G4LindhardSorensenIonModel() [2/2]

G4LindhardSorensenIonModel::G4LindhardSorensenIonModel ( const G4LindhardSorensenIonModel )
delete

Member Function Documentation

◆ ChargeSquareRatio()

G4double G4VEmModel::ChargeSquareRatio ( const G4Track track)
virtualinherited

Reimplemented in G4BraggIonGasModel, and G4BetheBlochIonGasModel.

Definition at line 374 of file G4VEmModel.cc.

375{
377 track.GetMaterial(), track.GetKineticEnergy());
378}
const G4ParticleDefinition * GetParticleDefinition() const
G4Material * GetMaterial() const
G4double GetKineticEnergy() const
virtual G4double GetChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:382

References G4VEmModel::GetChargeSquareRatio(), G4Track::GetKineticEnergy(), G4Track::GetMaterial(), and G4Track::GetParticleDefinition().

Referenced by G4VEnergyLossProcess::PostStepGetPhysicalInteractionLength().

◆ ComputeCrossSectionPerAtom() [1/2]

G4double G4VEmModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition part,
const G4Element elm,
G4double  kinEnergy,
G4double  cutEnergy = 0.0,
G4double  maxEnergy = DBL_MAX 
)
inlineinherited

Definition at line 566 of file G4VEmModel.hh.

571{
573 return ComputeCrossSectionPerAtom(part,kinEnergy,elm->GetZ(),elm->GetN(),
574 cutEnergy,maxEnergy);
575}
G4double GetZ() const
Definition: G4Element.hh:131
G4double GetN() const
Definition: G4Element.hh:135
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:341
void SetCurrentElement(const G4Element *)
Definition: G4VEmModel.hh:497

References G4VEmModel::ComputeCrossSectionPerAtom(), G4Element::GetN(), G4Element::GetZ(), and G4VEmModel::SetCurrentElement().

◆ ComputeCrossSectionPerAtom() [2/2]

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

Reimplemented from G4VEmModel.

Definition at line 192 of file G4LindhardSorensenIonModel.cc.

198{
199 return Z*ComputeCrossSectionPerElectron(p,kineticEnergy,cutEnergy,maxEnergy);
200}
const G4int Z[17]
G4double ComputeCrossSectionPerElectron(const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)

References ComputeCrossSectionPerElectron(), and Z.

◆ ComputeCrossSectionPerElectron()

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

Definition at line 172 of file G4LindhardSorensenIonModel.cc.

177{
178 // take into account formfactor
179 G4double tmax = MaxSecondaryEnergy(p, kinEnergy);
180 G4double emax = std::min(tmax, maxKinEnergy);
181 G4double escaled = kinEnergy*pRatio;
182 G4double cross = (escaled <= fElimit)
183 ? fBraggModel->ComputeCrossSectionPerElectron(p,kinEnergy,cutEnergy,emax)
184 : fBBModel->ComputeCrossSectionPerElectron(p,kinEnergy,cutEnergy,emax);
185 // G4cout << "LS: e= " << kinEnergy << " tmin= " << cutEnergy
186 // << " tmax= " << maxEnergy << " cross= " << cross << G4endl;
187 return cross;
188}
static const G4double emax
double G4double
Definition: G4Types.hh:83
virtual G4double ComputeCrossSectionPerElectron(const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
G4double ComputeCrossSectionPerElectron(const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
G4double MaxSecondaryEnergy(const G4ParticleDefinition *, G4double kinEnergy) override
T min(const T t1, const T t2)
brief Return the smallest of the two arguments

References G4BetheBlochModel::ComputeCrossSectionPerElectron(), G4BraggIonModel::ComputeCrossSectionPerElectron(), emax, fBBModel, fBraggModel, fElimit, MaxSecondaryEnergy(), G4INCL::Math::min(), and pRatio.

Referenced by ComputeCrossSectionPerAtom(), and CrossSectionPerVolume().

◆ ComputeCrossSectionPerShell()

G4double G4VEmModel::ComputeCrossSectionPerShell ( const G4ParticleDefinition ,
G4int  Z,
G4int  shellIdx,
G4double  kinEnergy,
G4double  cutEnergy = 0.0,
G4double  maxEnergy = DBL_MAX 
)
virtualinherited

Definition at line 351 of file G4VEmModel.cc.

354{
355 return 0.0;
356}

Referenced by G4EmCalculator::ComputeCrossSectionPerShell().

◆ ComputeDEDX()

G4double G4VEmModel::ComputeDEDX ( const G4MaterialCutsCouple couple,
const G4ParticleDefinition part,
G4double  kineticEnergy,
G4double  cutEnergy = DBL_MAX 
)
inlineinherited

Definition at line 528 of file G4VEmModel.hh.

532{
533 SetCurrentCouple(couple);
534 return pFactor*ComputeDEDXPerVolume(pBaseMaterial,part,kinEnergy,cutEnergy);
535}
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:472
const G4Material * pBaseMaterial
Definition: G4VEmModel.hh:427
G4double pFactor
Definition: G4VEmModel.hh:432
virtual G4double ComputeDEDXPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
Definition: G4VEmModel.cc:228

References G4VEmModel::ComputeDEDXPerVolume(), G4VEmModel::pBaseMaterial, G4VEmModel::pFactor, and G4VEmModel::SetCurrentCouple().

◆ ComputeDEDXPerVolume()

G4double G4LindhardSorensenIonModel::ComputeDEDXPerVolume ( const G4Material mat,
const G4ParticleDefinition p,
G4double  kineticEnergy,
G4double  cutEnergy 
)
overridevirtual

Reimplemented from G4VEmModel.

Definition at line 218 of file G4LindhardSorensenIonModel.cc.

222{
223 // formfactor is taken into account in CorrectionsAlongStep(..)
224 G4double tmax = MaxSecondaryEnergy(p, kinEnergy);
225 G4double cutEnergy = std::min(std::min(cut,tmax), tlimit);
226
227 G4double escaled = kinEnergy*pRatio;
228 G4double dedx = (escaled <= fElimit)
229 ? fBraggModel->ComputeDEDXPerVolume(mat, p, kinEnergy, cutEnergy)
230 : fBBModel->ComputeDEDXPerVolume(mat, p, kinEnergy, cutEnergy);
231
232 //G4cout << "E(MeV)=" << kinEnergy/MeV << " dedx=" << dedx
233 // << " " << material->GetName() << " Ecut(MeV)=" << cutEnergy << G4endl;
234 return dedx;
235}
G4double ComputeDEDXPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy) override
G4double ComputeDEDXPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy) override

References G4BetheBlochModel::ComputeDEDXPerVolume(), G4BraggIonModel::ComputeDEDXPerVolume(), fBBModel, fBraggModel, fElimit, MaxSecondaryEnergy(), G4INCL::Math::min(), pRatio, and tlimit.

◆ ComputeDEDXPerVolumeLS()

G4double G4LindhardSorensenIonModel::ComputeDEDXPerVolumeLS ( const G4Material ,
const G4ParticleDefinition ,
G4double  kinEnergy,
G4double  cutEnergy 
)
private

◆ ComputeMeanFreePath()

G4double G4VEmModel::ComputeMeanFreePath ( const G4ParticleDefinition part,
G4double  kineticEnergy,
const G4Material material,
G4double  cutEnergy = 0.0,
G4double  maxEnergy = DBL_MAX 
)
inlineinherited

Definition at line 553 of file G4VEmModel.hh.

558{
559 G4double cross = CrossSectionPerVolume(material,part,ekin,emin,emax);
560 return (cross > 0.0) ? 1./cross : DBL_MAX;
561}
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:237
string material
Definition: eplot.py:19
#define DBL_MAX
Definition: templates.hh:62

References G4VEmModel::CrossSectionPerVolume(), DBL_MAX, emax, and eplot::material.

◆ CorrectionsAlongStep()

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

Reimplemented from G4VEmModel.

Definition at line 239 of file G4LindhardSorensenIonModel.cc.

244{
245 // no correction at the last step
246 const G4double preKinEnergy = dp->GetKineticEnergy();
247 if(eloss >= preKinEnergy) { return; }
248
249 const G4ParticleDefinition* p = dp->GetDefinition();
250 const G4Material* mat = couple->GetMaterial();
251 const G4double eDensity = mat->GetElectronDensity();
252 const G4double e = std::max(preKinEnergy - eloss*0.5, preKinEnergy*0.75);
253 const G4double tmax = MaxSecondaryEnergy(p, e);
254 const G4double escaled = e*pRatio;
255 const G4double tau = e/mass;
256
259 const G4int Z = std::max(80, p->GetAtomicNumber());
260
261 G4double res;
262 if(escaled <= fElimit) {
263 res = fIonData->GetDEDX(mat, Z, escaled, G4Log(escaled));
264 if(res > 0.0) {
265 G4double cut = couple->GetProductionCuts()->GetProductionCut(1);
266 if(cut < tmax) {
267 const G4double x = cut/tmax;
268 res += (G4Log(x)*(tau + 1.)*(tau + 1.)/(tau * (tau + 2.0)) + 1.0 - x)
269 *q2*CLHEP::twopi_mc2_rcl2*eDensity;
270 }
271 } else {
272 res = eloss*q2*corr->EffectiveChargeCorrection(p,mat,e)/chargeSquare;
273 }
274 } else {
275 const G4double gam = tau + 1.0;
276 const G4double beta2 = tau * (tau+2.0)/(gam*gam);
277 G4double deltaL0 = 2.0*corr->BarkasCorrection(p, mat, e)*(charge-1.)/charge;
278 G4double deltaL = lsdata->GetDeltaL(Zin, gam);
279
280 res = eloss +
281 CLHEP::twopi_mc2_rcl2*q2*eDensity*(deltaL+deltaL0)*length/beta2;
282 /*
283 G4cout << "G4LindhardSorensenIonModel::CorrectionsAlongStep: E(GeV)= "
284 << preKinEnergy/GeV << " eloss(MeV)= " << eloss
285 << " L= " << eloss*beta2/(twopi_mc2_rcl2*chargeSquare*eDensity*length)
286 << " dL0= " << deltaL0
287 << " dL= " << deltaL << G4endl;
288 */
289 }
290 if(res > preKinEnergy) { res = preKinEnergy; }
291 else if(res < 0.0) { res = eloss; }
292 //G4cout << "G4LindhardSorensenIonModel::CorrectionsAlongStep: E(GeV)="
293 // << preKinEnergy/GeV << " eloss(MeV)=" << eloss
294 //<< " res(MeV)=" << res << G4endl;
295 eloss = res;
296}
int G4int
Definition: G4Types.hh:85
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
G4double EffectiveChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4double BarkasCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4double EffectiveChargeCorrection(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
G4double GetDEDX(const G4Material *, const G4int Z, const G4double e, const G4double loge) const
G4double GetDeltaL(G4int Z, G4double gamma) const
static G4LindhardSorensenData * lsdata
const G4Material * GetMaterial() const
G4ProductionCuts * GetProductionCuts() const
G4double GetElectronDensity() const
Definition: G4Material.hh:213
G4int GetAtomicNumber() const
G4double GetProductionCut(G4int index) const
virtual void SetParticleAndCharge(const G4ParticleDefinition *, G4double q2)
G4VEmFluctuationModel * GetModelOfFluctuations()
Definition: G4VEmModel.hh:614
static constexpr double twopi_mc2_rcl2
T max(const T t1, const T t2)
brief Return the largest of the two arguments

References G4EmCorrections::BarkasCorrection(), charge, chargeSquare, corr, G4EmCorrections::EffectiveChargeCorrection(), G4EmCorrections::EffectiveChargeSquareRatio(), fElimit, fIonData, G4Log(), G4InuclParticleNames::gam, G4ParticleDefinition::GetAtomicNumber(), G4IonICRU73Data::GetDEDX(), G4DynamicParticle::GetDefinition(), G4LindhardSorensenData::GetDeltaL(), G4Material::GetElectronDensity(), G4DynamicParticle::GetKineticEnergy(), G4MaterialCutsCouple::GetMaterial(), G4VEmModel::GetModelOfFluctuations(), G4ProductionCuts::GetProductionCut(), G4MaterialCutsCouple::GetProductionCuts(), lsdata, mass, G4INCL::Math::max(), MaxSecondaryEnergy(), pRatio, G4VEmFluctuationModel::SetParticleAndCharge(), CLHEP::twopi_mc2_rcl2, Z, and Zin.

◆ CrossSection()

G4double G4VEmModel::CrossSection ( const G4MaterialCutsCouple couple,
const G4ParticleDefinition part,
G4double  kineticEnergy,
G4double  cutEnergy = 0.0,
G4double  maxEnergy = DBL_MAX 
)
inlineinherited

◆ CrossSectionPerVolume()

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

Reimplemented from G4VEmModel.

Definition at line 204 of file G4LindhardSorensenIonModel.cc.

210{
211 return material->GetElectronDensity()
212 *ComputeCrossSectionPerElectron(p,kineticEnergy,cutEnergy,maxEnergy);
213}

References ComputeCrossSectionPerElectron(), and eplot::material.

◆ CurrentCouple()

const G4MaterialCutsCouple * G4VEmModel::CurrentCouple ( ) const
inlineprotectedinherited

◆ DeexcitationFlag()

G4bool G4VEmModel::DeexcitationFlag ( ) const
inlineinherited

Definition at line 704 of file G4VEmModel.hh.

705{
706 return flagDeexcitation;
707}
G4bool flagDeexcitation
Definition: G4VEmModel.hh:455

References G4VEmModel::flagDeexcitation.

Referenced by G4EmModelManager::DumpModelList().

◆ DefineForRegion()

void G4VEmModel::DefineForRegion ( const G4Region )
virtualinherited

Reimplemented in G4PAIModel, and G4PAIPhotModel.

Definition at line 360 of file G4VEmModel.cc.

361{}

Referenced by G4EmModelManager::AddEmModel().

◆ FillNumberOfSecondaries()

void G4VEmModel::FillNumberOfSecondaries ( G4int numberOfTriplets,
G4int numberOfRecoil 
)
virtualinherited

Definition at line 365 of file G4VEmModel.cc.

367{
368 numberOfTriplets = 0;
369 numberOfRecoil = 0;
370}

Referenced by G4VEmProcess::PostStepDoIt(), and G4VEnergyLossProcess::PostStepDoIt().

◆ ForceBuildTableFlag()

G4bool G4VEmModel::ForceBuildTableFlag ( ) const
inlineinherited

Definition at line 711 of file G4VEmModel.hh.

712{
713 return flagForceBuildTable;
714}
G4bool flagForceBuildTable
Definition: G4VEmModel.hh:456

References G4VEmModel::flagForceBuildTable.

Referenced by G4VMscModel::GetParticleChangeForMSC().

◆ GetAngularDistribution()

G4VEmAngularDistribution * G4VEmModel::GetAngularDistribution ( )
inlineinherited

Definition at line 621 of file G4VEmModel.hh.

622{
623 return anglModel;
624}
G4VEmAngularDistribution * anglModel
Definition: G4VEmModel.hh:414

References G4VEmModel::anglModel.

Referenced by G4EmModelManager::DumpModelList(), G4AtimaEnergyLossModel::Initialise(), G4BetheBlochModel::Initialise(), G4BraggIonModel::Initialise(), G4BraggModel::Initialise(), G4ICRU73QOModel::Initialise(), Initialise(), G4MollerBhabhaModel::Initialise(), G4AdjointBremsstrahlungModel::RapidSampleSecondaries(), G4LivermoreBremsstrahlungModel::SampleSecondaries(), G4eBremParametrizedModel::SampleSecondaries(), G4eBremsstrahlungRelModel::SampleSecondaries(), G4SeltzerBergerModel::SampleSecondaries(), G4DNABornIonisationModel1::SampleSecondaries(), G4DNABornIonisationModel2::SampleSecondaries(), G4DNAEmfietzoglouIonisationModel::SampleSecondaries(), G4DNARuddIonisationExtendedModel::SampleSecondaries(), G4DNARuddIonisationModel::SampleSecondaries(), G4PAIModel::SampleSecondaries(), G4PAIPhotModel::SampleSecondaries(), G4LivermoreIonisationModel::SampleSecondaries(), G4LivermorePhotoElectricModel::SampleSecondaries(), G4LivermoreRayleighModel::SampleSecondaries(), G4MicroElecInelasticModel::SampleSecondaries(), G4MicroElecInelasticModel_new::SampleSecondaries(), G4MuBremsstrahlungModel::SampleSecondaries(), G4MuPairProductionModel::SampleSecondaries(), G4AtimaEnergyLossModel::SampleSecondaries(), G4BetheBlochModel::SampleSecondaries(), G4BetheHeitlerModel::SampleSecondaries(), G4BraggIonModel::SampleSecondaries(), G4BraggModel::SampleSecondaries(), G4ICRU73QOModel::SampleSecondaries(), SampleSecondaries(), G4MollerBhabhaModel::SampleSecondaries(), G4PairProductionRelModel::SampleSecondaries(), G4PEEffectFluoModel::SampleSecondaries(), and G4IonParametrisedLossModel::SampleSecondaries().

◆ GetChargeSquareRatio() [1/2]

G4double G4LindhardSorensenIonModel::GetChargeSquareRatio ( ) const
inlineprotected

Definition at line 192 of file G4LindhardSorensenIonModel.hh.

193{
194 return chargeSquare;
195}

References chargeSquare.

◆ GetChargeSquareRatio() [2/2]

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

Reimplemented from G4VEmModel.

Definition at line 121 of file G4LindhardSorensenIonModel.cc.

124{
125 return corr->EffectiveChargeSquareRatio(p,mat,kinEnergy);
126}

References corr, and G4EmCorrections::EffectiveChargeSquareRatio().

◆ GetCrossSectionTable()

G4PhysicsTable * G4VEmModel::GetCrossSectionTable ( )
inlineinherited

◆ GetCurrentElement()

const G4Element * G4VEmModel::GetCurrentElement ( ) const
inlineinherited

◆ GetCurrentIsotope()

const G4Isotope * G4VEmModel::GetCurrentIsotope ( ) const
inlineinherited

Definition at line 512 of file G4VEmModel.hh.

513{
514 return fCurrentIsotope;
515}
const G4Isotope * fCurrentIsotope
Definition: G4VEmModel.hh:418

References G4VEmModel::fCurrentIsotope.

Referenced by G4VEmProcess::GetTargetIsotope().

◆ GetElementData()

G4ElementData * G4VEmModel::GetElementData ( )
inlineinherited

◆ GetElementSelectors()

std::vector< G4EmElementSelector * > * G4VEmModel::GetElementSelectors ( )
inlineinherited

◆ GetModelOfFluctuations()

G4VEmFluctuationModel * G4VEmModel::GetModelOfFluctuations ( )
inlineinherited

◆ GetName()

const G4String & G4VEmModel::GetName ( ) const
inlineinherited

◆ GetPartialCrossSection()

G4double G4VEmModel::GetPartialCrossSection ( const G4Material ,
G4int  level,
const G4ParticleDefinition ,
G4double  kineticEnergy 
)
virtualinherited

◆ GetParticleChangeForGamma()

G4ParticleChangeForGamma * G4VEmModel::GetParticleChangeForGamma ( )
protectedinherited

Definition at line 123 of file G4VEmModel.cc.

124{
125 G4ParticleChangeForGamma* p = nullptr;
126 if (pParticleChange != nullptr) {
127 p = static_cast<G4ParticleChangeForGamma*>(pParticleChange);
128 } else {
129 p = new G4ParticleChangeForGamma();
130 pParticleChange = p;
131 }
132 if(fTripletModel != nullptr) { fTripletModel->SetParticleChange(p); }
133 return p;
134}
void SetParticleChange(G4VParticleChange *, G4VEmFluctuationModel *f=nullptr)
Definition: G4VEmModel.cc:447
G4VParticleChange * pParticleChange
Definition: G4VEmModel.hh:425
G4VEmModel * fTripletModel
Definition: G4VEmModel.hh:415

References G4VEmModel::fTripletModel, G4VEmModel::pParticleChange, and G4VEmModel::SetParticleChange().

Referenced by G4MicroElecLOPhononModel::G4MicroElecLOPhononModel(), G4DNAChampionElasticModel::Initialise(), G4DNACPA100ElasticModel::Initialise(), G4DNADingfelderChargeDecreaseModel::Initialise(), G4DNADingfelderChargeIncreaseModel::Initialise(), G4DNAMeltonAttachmentModel::Initialise(), G4DNAMillerGreenExcitationModel::Initialise(), G4DNARuddIonisationExtendedModel::Initialise(), G4DNARuddIonisationModel::Initialise(), G4DNASancheExcitationModel::Initialise(), G4DNAScreenedRutherfordElasticModel::Initialise(), G4DNATransformElectronModel::Initialise(), G4DNAUeharaScreenedRutherfordElasticModel::Initialise(), G4LEPTSAttachmentModel::Initialise(), G4LEPTSDissociationModel::Initialise(), G4LEPTSElasticModel::Initialise(), G4LEPTSExcitationModel::Initialise(), G4LEPTSIonisationModel::Initialise(), G4LEPTSPositroniumModel::Initialise(), G4LEPTSRotExcitationModel::Initialise(), G4LEPTSVibExcitationModel::Initialise(), G4BoldyshevTripletModel::Initialise(), G4eplusTo3GammaOKVIModel::Initialise(), G4eSingleCoulombScatteringModel::Initialise(), G4IonCoulombScatteringModel::Initialise(), G4eeToHadronsMultiModel::Initialise(), G4JAEAElasticScatteringModel::Initialise(), G4JAEAPolarizedElasticScatteringModel::Initialise(), G4LivermoreComptonModel::Initialise(), G4LivermoreGammaConversion5DModel::Initialise(), G4LivermoreGammaConversionModel::Initialise(), G4LivermoreNuclearGammaConversionModel::Initialise(), G4LivermorePhotoElectricModel::Initialise(), G4LivermorePolarizedComptonModel::Initialise(), G4LivermorePolarizedGammaConversionModel::Initialise(), G4LivermorePolarizedRayleighModel::Initialise(), G4LivermoreRayleighModel::Initialise(), G4LowEPComptonModel::Initialise(), G4LowEPPolarizedComptonModel::Initialise(), G4MicroElecElasticModel::Initialise(), G4MicroElecElasticModel_new::Initialise(), G4MicroElecInelasticModel::Initialise(), G4MicroElecInelasticModel_new::Initialise(), G4MicroElecLOPhononModel::Initialise(), G4PenelopeAnnihilationModel::Initialise(), G4PenelopeComptonModel::Initialise(), G4PenelopeGammaConversionModel::Initialise(), G4PenelopePhotoElectricModel::Initialise(), G4PenelopeRayleighModel::Initialise(), G4PenelopeRayleighModelMI::Initialise(), G4PolarizedAnnihilationModel::Initialise(), G4BetheHeitlerModel::Initialise(), G4eCoulombScatteringModel::Initialise(), G4eDPWACoulombScatteringModel::Initialise(), G4eeToTwoGammaModel::Initialise(), G4eplusTo2GammaOKVIModel::Initialise(), G4hCoulombScatteringModel::Initialise(), G4KleinNishinaCompton::Initialise(), G4KleinNishinaModel::Initialise(), G4PairProductionRelModel::Initialise(), G4PEEffectFluoModel::Initialise(), G4XrayRayleighModel::Initialise(), G4DNABornExcitationModel1::Initialise(), G4DNABornExcitationModel2::Initialise(), G4DNABornIonisationModel1::Initialise(), G4DNABornIonisationModel2::Initialise(), G4DNACPA100ExcitationModel::Initialise(), G4DNACPA100IonisationModel::Initialise(), G4DNADiracRMatrixExcitationModel::Initialise(), G4DNAEmfietzoglouExcitationModel::Initialise(), G4DNAEmfietzoglouIonisationModel::Initialise(), G4DNAQuinnPlasmonExcitationModel::Initialise(), G4DNARelativisticIonisationModel::Initialise(), G4DNAELSEPAElasticModel::Initialise(), G4DNAModelInterface::Initialise(), and G4DNAIonElasticModel::Initialise().

◆ GetParticleChangeForLoss()

G4ParticleChangeForLoss * G4VEmModel::GetParticleChangeForLoss ( )
protectedinherited

◆ GetParticleCharge()

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

Reimplemented from G4VEmModel.

Definition at line 131 of file G4LindhardSorensenIonModel.cc.

134{
135 return corr->GetParticleCharge(p,mat,kinEnergy);
136}
G4double GetParticleCharge(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)

References corr, and G4EmCorrections::GetParticleCharge().

◆ GetTripletModel()

G4VEmModel * G4VEmModel::GetTripletModel ( )
inlineinherited

◆ HighEnergyActivationLimit()

G4double G4VEmModel::HighEnergyActivationLimit ( ) const
inlineinherited

◆ HighEnergyLimit()

G4double G4VEmModel::HighEnergyLimit ( ) const
inlineinherited

Definition at line 655 of file G4VEmModel.hh.

656{
657 return highLimit;
658}
G4double highLimit
Definition: G4VEmModel.hh:437

References G4VEmModel::highLimit.

Referenced by G4DNAChampionElasticModel::CrossSectionPerVolume(), G4DNACPA100ElasticModel::CrossSectionPerVolume(), G4DNACPA100ExcitationModel::CrossSectionPerVolume(), G4DNACPA100IonisationModel::CrossSectionPerVolume(), G4DNAEmfietzoglouExcitationModel::CrossSectionPerVolume(), G4DNAEmfietzoglouIonisationModel::CrossSectionPerVolume(), G4DNAMeltonAttachmentModel::CrossSectionPerVolume(), G4DNASancheExcitationModel::CrossSectionPerVolume(), G4DNAScreenedRutherfordElasticModel::CrossSectionPerVolume(), G4DNATransformElectronModel::CrossSectionPerVolume(), G4DNAELSEPAElasticModel::CrossSectionPerVolume(), G4DNAChampionElasticModel::G4DNAChampionElasticModel(), G4DNACPA100ElasticModel::G4DNACPA100ElasticModel(), G4DNAMeltonAttachmentModel::G4DNAMeltonAttachmentModel(), G4DNASancheExcitationModel::G4DNASancheExcitationModel(), G4DNAUeharaScreenedRutherfordElasticModel::G4DNAUeharaScreenedRutherfordElasticModel(), G4eeToHadronsModel::G4eeToHadronsModel(), G4IonParametrisedLossModel::G4IonParametrisedLossModel(), G4LivermorePolarizedRayleighModel::G4LivermorePolarizedRayleighModel(), G4ContinuousGainOfEnergy::GetContinuousStepLimit(), G4PenelopeBremsstrahlungModel::GetCrossSectionTableForCouple(), G4VMscModel::GetParticleChangeForMSC(), G4DNAChampionElasticModel::Initialise(), G4DNACPA100ElasticModel::Initialise(), G4DNADingfelderChargeDecreaseModel::Initialise(), G4DNADingfelderChargeIncreaseModel::Initialise(), G4DNAMeltonAttachmentModel::Initialise(), G4DNAMillerGreenExcitationModel::Initialise(), G4DNARuddIonisationExtendedModel::Initialise(), G4DNARuddIonisationModel::Initialise(), G4DNASancheExcitationModel::Initialise(), G4DNAScreenedRutherfordElasticModel::Initialise(), G4DNAUeharaScreenedRutherfordElasticModel::Initialise(), G4BoldyshevTripletModel::Initialise(), G4PAIModel::Initialise(), G4PAIPhotModel::Initialise(), G4JAEAElasticScatteringModel::Initialise(), G4JAEAPolarizedElasticScatteringModel::Initialise(), G4LivermoreComptonModel::Initialise(), G4LivermoreGammaConversion5DModel::Initialise(), G4LivermoreGammaConversionModel::Initialise(), G4LivermoreIonisationModel::Initialise(), G4LivermoreNuclearGammaConversionModel::Initialise(), G4LivermorePolarizedComptonModel::Initialise(), G4LivermorePolarizedGammaConversionModel::Initialise(), G4LivermoreRayleighModel::Initialise(), G4LowEPComptonModel::Initialise(), G4LowEPPolarizedComptonModel::Initialise(), G4MicroElecElasticModel::Initialise(), G4MicroElecElasticModel_new::Initialise(), G4MicroElecInelasticModel::Initialise(), G4MicroElecInelasticModel_new::Initialise(), G4PenelopeAnnihilationModel::Initialise(), G4PenelopeBremsstrahlungModel::Initialise(), G4PenelopeComptonModel::Initialise(), G4PenelopeGammaConversionModel::Initialise(), G4PenelopeIonisationModel::Initialise(), G4PenelopePhotoElectricModel::Initialise(), G4PenelopeRayleighModel::Initialise(), G4PenelopeRayleighModelMI::Initialise(), G4MuBremsstrahlungModel::Initialise(), G4MuPairProductionModel::Initialise(), G4eBremsstrahlungRelModel::Initialise(), G4eDPWACoulombScatteringModel::Initialise(), G4GoudsmitSaundersonMscModel::Initialise(), G4PairProductionRelModel::Initialise(), G4SeltzerBergerModel::Initialise(), G4WentzelVIModel::Initialise(), G4DNABornExcitationModel1::Initialise(), G4DNABornExcitationModel2::Initialise(), G4DNABornIonisationModel1::Initialise(), G4DNABornIonisationModel2::Initialise(), G4DNACPA100ExcitationModel::Initialise(), G4DNACPA100IonisationModel::Initialise(), G4DNADiracRMatrixExcitationModel::Initialise(), G4DNAEmfietzoglouExcitationModel::Initialise(), G4DNAEmfietzoglouIonisationModel::Initialise(), G4DNAQuinnPlasmonExcitationModel::Initialise(), G4DNARelativisticIonisationModel::Initialise(), G4EmModelManager::Initialise(), G4DNAELSEPAElasticModel::Initialise(), G4DNADummyModel::Initialise(), G4DNAIonElasticModel::Initialise(), G4mplIonisation::InitialiseEnergyLossProcess(), G4ionIonisation::InitialiseEnergyLossProcess(), G4MuBremsstrahlungModel::InitialiseLocal(), G4MuPairProductionModel::InitialiseLocal(), G4eBremsstrahlungRelModel::InitialiseLocal(), G4PairProductionRelModel::InitialiseLocal(), G4CoulombScattering::InitialiseProcess(), G4VEmProcess::PostStepDoIt(), G4VEmProcess::PreparePhysicsTable(), G4VEnergyLossProcess::PreparePhysicsTable(), G4VMultipleScattering::PreparePhysicsTable(), G4DNACPA100IonisationModel::SampleSecondaries(), G4DNAEmfietzoglouIonisationModel::SampleSecondaries(), G4DNASancheExcitationModel::SampleSecondaries(), G4EmConfigurator::SetExtraEmModel(), G4mplIonisationModel::SetParticle(), G4mplIonisationWithDeltaModel::SetParticle(), G4eBremsstrahlung::StreamProcessInfo(), and G4EmConfigurator::UpdateModelEnergyRange().

◆ Initialise()

void G4LindhardSorensenIonModel::Initialise ( const G4ParticleDefinition p,
const G4DataVector ptr 
)
overridevirtual

Implements G4VEmModel.

Definition at line 89 of file G4LindhardSorensenIonModel.cc.

91{
92 fBraggModel->Initialise(p, ptr);
93 fBBModel->Initialise(p, ptr);
94 SetParticle(p);
95 //G4cout << "G4LindhardSorensenIonModel::Initialise for "
96 // << p->GetParticleName() << G4endl;
97
98 // always false before the run
100
101 if(nullptr == fParticleChange) {
103 if(UseAngularGeneratorFlag() && nullptr == GetAngularDistribution()) {
105 }
106 }
107 if(IsMaster()) {
108 if(nullptr == lsdata) {
110 }
111 if(nullptr == fIonData) {
113 }
115 }
116}
void Initialise(const G4ParticleDefinition *, const G4DataVector &) override
void Initialise(const G4ParticleDefinition *, const G4DataVector &) override
void SetParticle(const G4ParticleDefinition *p)
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:621
G4bool IsMaster() const
Definition: G4VEmModel.hh:746
void SetDeexcitationFlag(G4bool val)
Definition: G4VEmModel.hh:823
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:628
G4bool UseAngularGeneratorFlag() const
Definition: G4VEmModel.hh:718
G4ParticleChangeForLoss * GetParticleChangeForLoss()
Definition: G4VEmModel.cc:108

References fBBModel, fBraggModel, fIonData, fParticleChange, G4VEmModel::GetAngularDistribution(), G4VEmModel::GetParticleChangeForLoss(), G4IonICRU73Data::Initialise(), G4BetheBlochModel::Initialise(), G4BraggIonModel::Initialise(), G4VEmModel::IsMaster(), lsdata, G4VEmModel::SetAngularDistribution(), G4VEmModel::SetDeexcitationFlag(), SetParticle(), and G4VEmModel::UseAngularGeneratorFlag().

◆ InitialiseElementSelectors()

void G4VEmModel::InitialiseElementSelectors ( const G4ParticleDefinition part,
const G4DataVector cuts 
)
inherited

Definition at line 138 of file G4VEmModel.cc.

140{
141 // using spline for element selectors should be investigated in details
142 // because small number of points may provide biased results
143 // large number of points requires significant increase of memory
144 G4bool spline = false;
145
146 //G4cout << "IES: for " << GetName() << " Emin(MeV)= " << lowLimit/MeV
147 // << " Emax(MeV)= " << highLimit/MeV << G4endl;
148
149 // two times less bins because probability functon is normalized
150 // so correspondingly is more smooth
151 if(highLimit <= lowLimit) { return; }
152
154
155 G4ProductionCutsTable* theCoupleTable=
157 G4int numOfCouples = theCoupleTable->GetTableSize();
158
159 // prepare vector
160 if(!elmSelectors) {
161 elmSelectors = new std::vector<G4EmElementSelector*>;
162 }
163 if(numOfCouples > nSelectors) {
164 for(G4int i=nSelectors; i<numOfCouples; ++i) {
165 elmSelectors->push_back(nullptr);
166 }
167 nSelectors = numOfCouples;
168 }
169
170 // initialise vector
171 for(G4int i=0; i<numOfCouples; ++i) {
172
173 // no need in element selectors for infinite cuts
174 if(cuts[i] == DBL_MAX) { continue; }
175
176 auto couple = theCoupleTable->GetMaterialCutsCouple(i);
177 auto material = couple->GetMaterial();
178 SetCurrentCouple(couple);
179
180 // selector already exist then delete
181 delete (*elmSelectors)[i];
182
183 G4double emin = std::max(lowLimit, MinPrimaryEnergy(material, part, cuts[i]));
184 G4double emax = std::max(highLimit, 10*emin);
185 static const G4double invlog106 = 1.0/(6*G4Log(10.));
186 G4int nbins = (G4int)(nbinsPerDec*G4Log(emax/emin)*invlog106);
187 nbins = std::max(nbins, 3);
188
189 (*elmSelectors)[i] = new G4EmElementSelector(this,material,nbins,
190 emin,emax,spline);
191 ((*elmSelectors)[i])->Initialise(part, cuts[i]);
192 /*
193 G4cout << "G4VEmModel::InitialiseElmSelectors i= " << i
194 << " " << part->GetParticleName()
195 << " for " << GetName() << " cut= " << cuts[i]
196 << " " << (*elmSelectors)[i] << G4endl;
197 ((*elmSelectors)[i])->Dump(part);
198 */
199 }
200}
bool G4bool
Definition: G4Types.hh:86
static G4EmParameters * Instance()
G4int NumberOfBinsPerDecade() const
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
std::size_t GetTableSize() const
static G4ProductionCutsTable * GetProductionCutsTable()
virtual G4double MinPrimaryEnergy(const G4Material *, const G4ParticleDefinition *, G4double cut=0.0)
Definition: G4VEmModel.cc:415
G4double lowLimit
Definition: G4VEmModel.hh:436
G4int nSelectors
Definition: G4VEmModel.hh:443
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)=0

References DBL_MAX, G4VEmModel::elmSelectors, emax, G4Log(), G4ProductionCutsTable::GetMaterialCutsCouple(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetTableSize(), G4VEmModel::highLimit, G4VEmModel::Initialise(), G4EmParameters::Instance(), G4VEmModel::lowLimit, eplot::material, G4INCL::Math::max(), G4VEmModel::MinPrimaryEnergy(), G4VEmModel::nSelectors, G4EmParameters::NumberOfBinsPerDecade(), and G4VEmModel::SetCurrentCouple().

Referenced by G4eSingleCoulombScatteringModel::Initialise(), G4PAIModel::Initialise(), G4PAIPhotModel::Initialise(), G4JAEAElasticScatteringModel::Initialise(), G4JAEAPolarizedElasticScatteringModel::Initialise(), G4LivermoreComptonModel::Initialise(), G4LivermoreGammaConversion5DModel::Initialise(), G4LivermoreGammaConversionModel::Initialise(), G4LivermoreNuclearGammaConversionModel::Initialise(), G4LivermorePolarizedComptonModel::Initialise(), G4LivermorePolarizedGammaConversionModel::Initialise(), G4LivermorePolarizedRayleighModel::Initialise(), G4LivermoreRayleighModel::Initialise(), G4LowEPComptonModel::Initialise(), G4LowEPPolarizedComptonModel::Initialise(), G4PenelopePhotoElectricModel::Initialise(), G4MuBremsstrahlungModel::Initialise(), G4MuPairProductionModel::Initialise(), G4BetheHeitlerModel::Initialise(), G4eBremParametrizedModel::Initialise(), G4eBremsstrahlungRelModel::Initialise(), G4eCoulombScatteringModel::Initialise(), G4eDPWACoulombScatteringModel::Initialise(), G4hCoulombScatteringModel::Initialise(), G4KleinNishinaCompton::Initialise(), G4KleinNishinaModel::Initialise(), G4PairProductionRelModel::Initialise(), G4SeltzerBergerModel::Initialise(), and G4XrayRayleighModel::Initialise().

◆ InitialiseForElement()

void G4VEmModel::InitialiseForElement ( const G4ParticleDefinition ,
G4int  Z 
)
virtualinherited

◆ InitialiseForMaterial()

void G4VEmModel::InitialiseForMaterial ( const G4ParticleDefinition part,
const G4Material material 
)
virtualinherited

Definition at line 209 of file G4VEmModel.cc.

211{
212 if(material != nullptr) {
213 size_t n = material->GetNumberOfElements();
214 for(size_t i=0; i<n; ++i) {
215 G4int Z = material->GetElement(i)->GetZasInt();
216 InitialiseForElement(part, Z);
217 }
218 }
219}
virtual void InitialiseForElement(const G4ParticleDefinition *, G4int Z)
Definition: G4VEmModel.cc:223

References G4VEmModel::InitialiseForElement(), eplot::material, CLHEP::detail::n, and Z.

Referenced by G4EmCalculator::FindEmModel().

◆ InitialiseLocal()

void G4VEmModel::InitialiseLocal ( const G4ParticleDefinition ,
G4VEmModel masterModel 
)
virtualinherited

◆ InitialiseLS()

void G4LindhardSorensenIonModel::InitialiseLS ( )
private

◆ IsActive()

G4bool G4VEmModel::IsActive ( G4double  kinEnergy) const
inlineinherited

◆ IsLocked()

G4bool G4VEmModel::IsLocked ( ) const
inlineinherited

◆ IsMaster()

G4bool G4VEmModel::IsMaster ( ) const
inlineinherited

Definition at line 746 of file G4VEmModel.hh.

747{
748 return isMaster;
749}
G4bool isMaster
Definition: G4VEmModel.hh:457

References G4VEmModel::isMaster.

Referenced by G4PenelopeBremsstrahlungModel::BuildXSTable(), G4PenelopeBremsstrahlungModel::ClearTables(), G4MuPairProductionModel::DataCorrupted(), G4PenelopePhotoElectricModel::GetNumberOfShellXS(), G4VMscModel::GetParticleChangeForMSC(), G4BoldyshevTripletModel::Initialise(), G4eSingleCoulombScatteringModel::Initialise(), G4PAIModel::Initialise(), G4PAIPhotModel::Initialise(), G4EmMultiModel::Initialise(), G4mplIonisationModel::Initialise(), G4mplIonisationWithDeltaModel::Initialise(), G4JAEAElasticScatteringModel::Initialise(), G4JAEAPolarizedElasticScatteringModel::Initialise(), G4LivermoreBremsstrahlungModel::Initialise(), G4LivermoreComptonModel::Initialise(), G4LivermoreGammaConversion5DModel::Initialise(), G4LivermoreGammaConversionModel::Initialise(), G4LivermoreNuclearGammaConversionModel::Initialise(), G4LivermorePhotoElectricModel::Initialise(), G4LivermorePolarizedComptonModel::Initialise(), G4LivermorePolarizedGammaConversionModel::Initialise(), G4LivermorePolarizedRayleighModel::Initialise(), G4LivermoreRayleighModel::Initialise(), G4LowEPComptonModel::Initialise(), G4LowEPPolarizedComptonModel::Initialise(), G4PenelopeAnnihilationModel::Initialise(), G4PenelopeBremsstrahlungModel::Initialise(), G4PenelopeComptonModel::Initialise(), G4PenelopeGammaConversionModel::Initialise(), G4PenelopeIonisationModel::Initialise(), G4PenelopePhotoElectricModel::Initialise(), G4PenelopeRayleighModel::Initialise(), G4PenelopeRayleighModelMI::Initialise(), G4MuBremsstrahlungModel::Initialise(), G4MuPairProductionModel::Initialise(), G4BetheBlochModel::Initialise(), G4BetheHeitlerModel::Initialise(), G4BraggIonModel::Initialise(), G4BraggModel::Initialise(), G4eBremParametrizedModel::Initialise(), G4eBremsstrahlungRelModel::Initialise(), G4eCoulombScatteringModel::Initialise(), G4eDPWACoulombScatteringModel::Initialise(), G4eeToTwoGammaModel::Initialise(), G4eplusTo2GammaOKVIModel::Initialise(), G4GoudsmitSaundersonMscModel::Initialise(), G4hCoulombScatteringModel::Initialise(), G4KleinNishinaCompton::Initialise(), G4KleinNishinaModel::Initialise(), Initialise(), G4PairProductionRelModel::Initialise(), G4SeltzerBergerModel::Initialise(), G4WentzelVIModel::Initialise(), G4PenelopeBremsstrahlungModel::InitialiseLocal(), G4PenelopeGammaConversionModel::ReadDataFile(), G4PenelopePhotoElectricModel::ReadDataFile(), G4BetheHeitlerModel::~G4BetheHeitlerModel(), G4BoldyshevTripletModel::~G4BoldyshevTripletModel(), G4BraggIonModel::~G4BraggIonModel(), G4BraggModel::~G4BraggModel(), G4eBremsstrahlungRelModel::~G4eBremsstrahlungRelModel(), G4eDPWACoulombScatteringModel::~G4eDPWACoulombScatteringModel(), G4GoudsmitSaundersonMscModel::~G4GoudsmitSaundersonMscModel(), G4JAEAElasticScatteringModel::~G4JAEAElasticScatteringModel(), G4JAEAPolarizedElasticScatteringModel::~G4JAEAPolarizedElasticScatteringModel(), G4LivermoreBremsstrahlungModel::~G4LivermoreBremsstrahlungModel(), G4LivermoreComptonModel::~G4LivermoreComptonModel(), G4LivermoreGammaConversion5DModel::~G4LivermoreGammaConversion5DModel(), G4LivermoreGammaConversionModel::~G4LivermoreGammaConversionModel(), G4LivermoreNuclearGammaConversionModel::~G4LivermoreNuclearGammaConversionModel(), G4LivermorePhotoElectricModel::~G4LivermorePhotoElectricModel(), G4LivermorePolarizedComptonModel::~G4LivermorePolarizedComptonModel(), G4LivermorePolarizedGammaConversionModel::~G4LivermorePolarizedGammaConversionModel(), G4LivermorePolarizedRayleighModel::~G4LivermorePolarizedRayleighModel(), G4LivermoreRayleighModel::~G4LivermoreRayleighModel(), G4LowEPComptonModel::~G4LowEPComptonModel(), G4LowEPPolarizedComptonModel::~G4LowEPPolarizedComptonModel(), G4mplIonisationModel::~G4mplIonisationModel(), G4mplIonisationWithDeltaModel::~G4mplIonisationWithDeltaModel(), G4PAIModel::~G4PAIModel(), G4PAIPhotModel::~G4PAIPhotModel(), G4PairProductionRelModel::~G4PairProductionRelModel(), G4PenelopeBremsstrahlungModel::~G4PenelopeBremsstrahlungModel(), G4PenelopeGammaConversionModel::~G4PenelopeGammaConversionModel(), G4PenelopeIonisationModel::~G4PenelopeIonisationModel(), G4PenelopePhotoElectricModel::~G4PenelopePhotoElectricModel(), G4PenelopeRayleighModel::~G4PenelopeRayleighModel(), G4PenelopeRayleighModelMI::~G4PenelopeRayleighModelMI(), G4SeltzerBergerModel::~G4SeltzerBergerModel(), and G4WentzelVIModel::~G4WentzelVIModel().

◆ LowEnergyActivationLimit()

G4double G4VEmModel::LowEnergyActivationLimit ( ) const
inlineinherited

◆ LowEnergyLimit()

G4double G4VEmModel::LowEnergyLimit ( ) const
inlineinherited

Definition at line 662 of file G4VEmModel.hh.

663{
664 return lowLimit;
665}

References G4VEmModel::lowLimit.

Referenced by G4eBremsstrahlungRelModel::ComputeCrossSectionPerAtom(), G4KleinNishinaCompton::ComputeCrossSectionPerAtom(), G4KleinNishinaModel::ComputeCrossSectionPerAtom(), G4LivermoreComptonModel::ComputeCrossSectionPerAtom(), G4LivermorePolarizedComptonModel::ComputeCrossSectionPerAtom(), G4LowEPComptonModel::ComputeCrossSectionPerAtom(), G4LowEPPolarizedComptonModel::ComputeCrossSectionPerAtom(), G4mplIonisationWithDeltaModel::ComputeCrossSectionPerElectron(), G4EmCalculator::ComputeDEDX(), G4eBremsstrahlungRelModel::ComputeDEDXPerVolume(), G4mplIonisationWithDeltaModel::ComputeDEDXPerVolume(), G4IonParametrisedLossModel::ComputeDEDXPerVolume(), G4IonParametrisedLossModel::CorrectionsAlongStep(), G4PenelopeRayleighModelMI::CrossSectionPerVolume(), G4PenelopeComptonModel::CrossSectionPerVolume(), G4DNAChampionElasticModel::CrossSectionPerVolume(), G4DNACPA100ElasticModel::CrossSectionPerVolume(), G4DNACPA100ExcitationModel::CrossSectionPerVolume(), G4DNACPA100IonisationModel::CrossSectionPerVolume(), G4DNAEmfietzoglouExcitationModel::CrossSectionPerVolume(), G4DNAEmfietzoglouIonisationModel::CrossSectionPerVolume(), G4DNAMeltonAttachmentModel::CrossSectionPerVolume(), G4DNASancheExcitationModel::CrossSectionPerVolume(), G4DNAScreenedRutherfordElasticModel::CrossSectionPerVolume(), G4DNAELSEPAElasticModel::CrossSectionPerVolume(), G4EmCalculator::FindEmModel(), G4DNAChampionElasticModel::G4DNAChampionElasticModel(), G4DNACPA100ElasticModel::G4DNACPA100ElasticModel(), G4DNAMeltonAttachmentModel::G4DNAMeltonAttachmentModel(), G4DNASancheExcitationModel::G4DNASancheExcitationModel(), G4DNAUeharaScreenedRutherfordElasticModel::G4DNAUeharaScreenedRutherfordElasticModel(), G4eeToHadronsModel::G4eeToHadronsModel(), G4LivermorePolarizedRayleighModel::G4LivermorePolarizedRayleighModel(), G4PenelopeBremsstrahlungModel::GetCrossSectionTableForCouple(), G4VMscModel::GetParticleChangeForMSC(), G4DNAChampionElasticModel::Initialise(), G4DNACPA100ElasticModel::Initialise(), G4DNADingfelderChargeDecreaseModel::Initialise(), G4DNADingfelderChargeIncreaseModel::Initialise(), G4DNAMeltonAttachmentModel::Initialise(), G4DNAMillerGreenExcitationModel::Initialise(), G4DNARuddIonisationExtendedModel::Initialise(), G4DNARuddIonisationModel::Initialise(), G4DNASancheExcitationModel::Initialise(), G4DNAScreenedRutherfordElasticModel::Initialise(), G4DNAUeharaScreenedRutherfordElasticModel::Initialise(), G4BoldyshevTripletModel::Initialise(), G4PAIModel::Initialise(), G4PAIPhotModel::Initialise(), G4JAEAElasticScatteringModel::Initialise(), G4JAEAPolarizedElasticScatteringModel::Initialise(), G4LivermoreComptonModel::Initialise(), G4LivermoreGammaConversion5DModel::Initialise(), G4LivermoreGammaConversionModel::Initialise(), G4LivermoreIonisationModel::Initialise(), G4LivermoreNuclearGammaConversionModel::Initialise(), G4LivermorePolarizedComptonModel::Initialise(), G4LivermorePolarizedGammaConversionModel::Initialise(), G4LivermoreRayleighModel::Initialise(), G4LowEPComptonModel::Initialise(), G4LowEPPolarizedComptonModel::Initialise(), G4MicroElecElasticModel::Initialise(), G4MicroElecElasticModel_new::Initialise(), G4MicroElecInelasticModel::Initialise(), G4MicroElecInelasticModel_new::Initialise(), G4PenelopeAnnihilationModel::Initialise(), G4PenelopeBremsstrahlungModel::Initialise(), G4PenelopeComptonModel::Initialise(), G4PenelopeGammaConversionModel::Initialise(), G4PenelopeIonisationModel::Initialise(), G4PenelopePhotoElectricModel::Initialise(), G4PenelopeRayleighModel::Initialise(), G4PenelopeRayleighModelMI::Initialise(), G4MuPairProductionModel::Initialise(), G4eBremParametrizedModel::Initialise(), G4eBremsstrahlungRelModel::Initialise(), G4eDPWACoulombScatteringModel::Initialise(), G4GoudsmitSaundersonMscModel::Initialise(), G4PairProductionRelModel::Initialise(), G4SeltzerBergerModel::Initialise(), G4WentzelVIModel::Initialise(), G4DNABornExcitationModel1::Initialise(), G4DNABornExcitationModel2::Initialise(), G4DNABornIonisationModel1::Initialise(), G4DNABornIonisationModel2::Initialise(), G4DNACPA100ExcitationModel::Initialise(), G4DNACPA100IonisationModel::Initialise(), G4DNADiracRMatrixExcitationModel::Initialise(), G4DNAEmfietzoglouExcitationModel::Initialise(), G4DNAEmfietzoglouIonisationModel::Initialise(), G4DNAQuinnPlasmonExcitationModel::Initialise(), G4DNARelativisticIonisationModel::Initialise(), G4EmModelManager::Initialise(), G4DNAELSEPAElasticModel::Initialise(), G4DNADummyModel::Initialise(), G4DNAIonElasticModel::Initialise(), G4mplIonisation::InitialiseEnergyLossProcess(), G4eBremsstrahlungRelModel::InitialiseLocal(), G4PairProductionRelModel::InitialiseLocal(), G4CoulombScattering::InitialiseProcess(), G4VEmProcess::PostStepDoIt(), G4eBremsstrahlungRelModel::SampleSecondaries(), G4DNACPA100IonisationModel::SampleSecondaries(), G4DNAEmfietzoglouIonisationModel::SampleSecondaries(), G4LivermoreComptonModel::SampleSecondaries(), G4LivermorePolarizedComptonModel::SampleSecondaries(), G4LowEPComptonModel::SampleSecondaries(), G4LowEPPolarizedComptonModel::SampleSecondaries(), G4PenelopeComptonModel::SampleSecondaries(), G4PolarizedComptonModel::SampleSecondaries(), G4KleinNishinaCompton::SampleSecondaries(), G4KleinNishinaModel::SampleSecondaries(), G4EmConfigurator::SetExtraEmModel(), G4mplIonisationModel::SetParticle(), G4mplIonisationWithDeltaModel::SetParticle(), and G4EmConfigurator::UpdateModelEnergyRange().

◆ LPMFlag()

G4bool G4VEmModel::LPMFlag ( ) const
inlineinherited

◆ MaxSecondaryEnergy()

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

Reimplemented from G4VEmModel.

Definition at line 419 of file G4LindhardSorensenIonModel.cc.

421{
422 // here particle type is checked for any method
423 SetParticle(pd);
424 G4double tau = kinEnergy/mass;
425 return 2.0*CLHEP::electron_mass_c2*tau*(tau + 2.) /
426 (1. + 2.0*(tau + 1.)*eRatio + eRatio*eRatio);
427}
static constexpr double electron_mass_c2

References CLHEP::electron_mass_c2, eRatio, mass, and SetParticle().

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

◆ MaxSecondaryKinEnergy()

G4double G4VEmModel::MaxSecondaryKinEnergy ( const G4DynamicParticle dynParticle)
inlineinherited

◆ MinEnergyCut()

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

Reimplemented from G4VEmModel.

Definition at line 163 of file G4LindhardSorensenIonModel.cc.

165{
167}
G4double GetMeanExcitationEnergy() const
G4IonisParamMat * GetIonisation() const
Definition: G4Material.hh:222

References G4Material::GetIonisation(), G4MaterialCutsCouple::GetMaterial(), and G4IonisParamMat::GetMeanExcitationEnergy().

◆ MinPrimaryEnergy()

G4double G4VEmModel::MinPrimaryEnergy ( const G4Material ,
const G4ParticleDefinition ,
G4double  cut = 0.0 
)
virtualinherited

◆ ModelDescription()

void G4VEmModel::ModelDescription ( std::ostream &  outFile) const
virtualinherited

Reimplemented in G4eeToHadronsMultiModel.

Definition at line 469 of file G4VEmModel.cc.

470{
471 outFile << "The description for this model has not been written yet.\n";
472}

◆ operator=()

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

◆ PolarAngleLimit()

G4double G4VEmModel::PolarAngleLimit ( ) const
inlineinherited

◆ SampleSecondaries()

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

Implements G4VEmModel.

Definition at line 300 of file G4LindhardSorensenIonModel.cc.

306{
307 G4double kineticEnergy = dp->GetKineticEnergy();
308 // take into account formfactor
309 G4double tmax = MaxSecondaryEnergy(dp->GetDefinition(),kineticEnergy);
310
311 G4double maxKinEnergy = std::min(maxEnergy,tmax);
312 if(minKinEnergy >= maxKinEnergy) { return; }
313
314 //G4cout << "G4LindhardSorensenIonModel::SampleSecondaries Emin= "
315 // << minKinEnergy << " Emax= " << maxKinEnergy << G4endl;
316
317 G4double totEnergy = kineticEnergy + mass;
318 G4double etot2 = totEnergy*totEnergy;
319 G4double beta2 = kineticEnergy*(kineticEnergy + 2.0*mass)/etot2;
320
321 G4double deltaKinEnergy, f;
322 G4double f1 = 0.0;
323 G4double fmax = 1.0;
324 if( 0.0 < spin ) { fmax += 0.5*maxKinEnergy*maxKinEnergy/etot2; }
325
326 CLHEP::HepRandomEngine* rndmEngineMod = G4Random::getTheEngine();
327 G4double rndm[2];
328
329 // sampling without nuclear size effect
330 do {
331 rndmEngineMod->flatArray(2, rndm);
332 deltaKinEnergy = minKinEnergy*maxKinEnergy
333 /(minKinEnergy*(1.0 - rndm[0]) + maxKinEnergy*rndm[0]);
334
335 f = 1.0 - beta2*deltaKinEnergy/tmax;
336 if( 0.0 < spin ) {
337 f1 = 0.5*deltaKinEnergy*deltaKinEnergy/etot2;
338 f += f1;
339 }
340
341 // Loop checking, 03-Aug-2015, Vladimir Ivanchenko
342 } while( fmax*rndm[1] > f);
343
344 // projectile formfactor - suppresion of high energy
345 // delta-electron production at high energy
346
347 G4double x = formfact*deltaKinEnergy;
348 if(x > 1.e-6) {
349
350 G4double x1 = 1.0 + x;
351 G4double grej = 1.0/(x1*x1);
352 if( 0.0 < spin ) {
353 G4double x2 = 0.5*electron_mass_c2*deltaKinEnergy/(mass*mass);
354 grej *= (1.0 + magMoment2*(x2 - f1/f)/(1.0 + x2));
355 }
356 if(grej > 1.1) {
357 G4cout << "### G4LindhardSorensenIonModel WARNING: grej= " << grej
358 << " " << dp->GetDefinition()->GetParticleName()
359 << " Ekin(MeV)= " << kineticEnergy
360 << " delEkin(MeV)= " << deltaKinEnergy
361 << G4endl;
362 }
363 if(rndmEngineMod->flat() > grej) { return; }
364 }
365
366 G4ThreeVector deltaDirection;
367
369
370 const G4Material* mat = couple->GetMaterial();
372
373 deltaDirection =
374 GetAngularDistribution()->SampleDirection(dp, deltaKinEnergy, Z, mat);
375
376 } else {
377
378 G4double deltaMomentum =
379 std::sqrt(deltaKinEnergy * (deltaKinEnergy + 2.0*electron_mass_c2));
380 G4double cost = deltaKinEnergy * (totEnergy + electron_mass_c2) /
381 (deltaMomentum * dp->GetTotalMomentum());
382 cost = std::min(cost, 1.0);
383 G4double sint = std::sqrt((1.0 - cost)*(1.0 + cost));
384
385 G4double phi = CLHEP::twopi*rndmEngineMod->flat();
386
387 deltaDirection.set(sint*std::cos(phi),sint*std::sin(phi), cost) ;
388 deltaDirection.rotateUz(dp->GetMomentumDirection());
389 }
390 /*
391 G4cout << "### G4LindhardSorensenIonModel "
392 << dp->GetDefinition()->GetParticleName()
393 << " Ekin(MeV)= " << kineticEnergy
394 << " delEkin(MeV)= " << deltaKinEnergy
395 << " tmin(MeV)= " << minKinEnergy
396 << " tmax(MeV)= " << maxKinEnergy
397 << " dir= " << dp->GetMomentumDirection()
398 << " dirDelta= " << deltaDirection
399 << G4endl;
400 */
401 // create G4DynamicParticle object for delta ray
402 G4DynamicParticle* delta =
403 new G4DynamicParticle(theElectron,deltaDirection,deltaKinEnergy);
404
405 vdp->push_back(delta);
406
407 // Change kinematics of primary particle
408 kineticEnergy -= deltaKinEnergy;
409 G4ThreeVector finalP = dp->GetMomentum() - delta->GetMomentum();
410 finalP = finalP.unit();
411
414}
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
Hep3Vector unit() const
void set(double x, double y, double z)
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:33
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)
const G4String & GetParticleName() const
virtual G4ThreeVector & SampleDirection(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, const G4Material *)=0
G4int SelectRandomAtomNumber(const G4Material *)
Definition: G4VEmModel.cc:297
static constexpr double twopi
Definition: SystemOfUnits.h:56
float electron_mass_c2
Definition: hepunit.py:273

References source.hepunit::electron_mass_c2, CLHEP::HepRandomEngine::flat(), CLHEP::HepRandomEngine::flatArray(), formfact, fParticleChange, G4cout, G4endl, G4VEmModel::GetAngularDistribution(), G4DynamicParticle::GetDefinition(), G4DynamicParticle::GetKineticEnergy(), G4MaterialCutsCouple::GetMaterial(), G4DynamicParticle::GetMomentum(), G4DynamicParticle::GetMomentumDirection(), G4ParticleDefinition::GetParticleName(), G4DynamicParticle::GetTotalMomentum(), magMoment2, mass, MaxSecondaryEnergy(), G4INCL::Math::min(), CLHEP::Hep3Vector::rotateUz(), G4VEmAngularDistribution::SampleDirection(), G4VEmModel::SelectRandomAtomNumber(), CLHEP::Hep3Vector::set(), G4ParticleChangeForLoss::SetProposedKineticEnergy(), G4ParticleChangeForLoss::SetProposedMomentumDirection(), spin, theElectron, CLHEP::twopi, CLHEP::Hep3Vector::unit(), G4VEmModel::UseAngularGeneratorFlag(), and Z.

◆ SecondaryThreshold()

G4double G4VEmModel::SecondaryThreshold ( ) const
inlineinherited

◆ SelectIsotopeNumber()

G4int G4VEmModel::SelectIsotopeNumber ( const G4Element elm)
inherited

Definition at line 319 of file G4VEmModel.cc.

320{
322 const size_t ni = elm->GetNumberOfIsotopes();
323 fCurrentIsotope = elm->GetIsotope(0);
324 size_t idx = 0;
325 if(ni > 1) {
326 const G4double* ab = elm->GetRelativeAbundanceVector();
328 for(; idx<ni; ++idx) {
329 x -= ab[idx];
330 if (x <= 0.0) {
331 fCurrentIsotope = elm->GetIsotope(idx);
332 break;
333 }
334 }
335 }
336 return fCurrentIsotope->GetN();
337}
static const G4double ab
#define G4UniformRand()
Definition: Randomize.hh:52
G4double * GetRelativeAbundanceVector() const
Definition: G4Element.hh:167
const G4Isotope * GetIsotope(G4int iso) const
Definition: G4Element.hh:170
size_t GetNumberOfIsotopes() const
Definition: G4Element.hh:159
G4int GetN() const
Definition: G4Isotope.hh:93

References ab, G4VEmModel::fCurrentIsotope, G4UniformRand, G4Element::GetIsotope(), G4Isotope::GetN(), G4Element::GetNumberOfIsotopes(), G4Element::GetRelativeAbundanceVector(), and G4VEmModel::SetCurrentElement().

Referenced by G4eSingleCoulombScatteringModel::SampleSecondaries(), G4IonCoulombScatteringModel::SampleSecondaries(), G4eCoulombScatteringModel::SampleSecondaries(), G4hCoulombScatteringModel::SampleSecondaries(), and G4BetheHeitler5DModel::SampleSecondaries().

◆ SelectRandomAtom() [1/2]

const G4Element * G4VEmModel::SelectRandomAtom ( const G4Material mat,
const G4ParticleDefinition pd,
G4double  kineticEnergy,
G4double  cutEnergy = 0.0,
G4double  maxEnergy = DBL_MAX 
)
inherited

Definition at line 275 of file G4VEmModel.cc.

280{
281 size_t n = mat->GetNumberOfElements();
282 fCurrentElement = mat->GetElement(0);
283 if (n > 1) {
284 const G4double x = G4UniformRand()*
285 G4VEmModel::CrossSectionPerVolume(mat,pd,kinEnergy,tcut,tmax);
286 for(size_t i=0; i<n; ++i) {
287 if (x <= xsec[i]) {
288 fCurrentElement = mat->GetElement(i);
289 break;
290 }
291 }
292 }
293 return fCurrentElement;
294}
const G4Element * GetElement(G4int iel) const
Definition: G4Material.hh:198
size_t GetNumberOfElements() const
Definition: G4Material.hh:182
std::vector< G4double > xsec
Definition: G4VEmModel.hh:466

References G4VEmModel::CrossSectionPerVolume(), G4VEmModel::fCurrentElement, G4UniformRand, G4Material::GetElement(), G4Material::GetNumberOfElements(), CLHEP::detail::n, and G4VEmModel::xsec.

◆ SelectRandomAtom() [2/2]

const G4Element * G4VEmModel::SelectRandomAtom ( const G4MaterialCutsCouple couple,
const G4ParticleDefinition part,
G4double  kineticEnergy,
G4double  cutEnergy = 0.0,
G4double  maxEnergy = DBL_MAX 
)
inlineinherited

Definition at line 580 of file G4VEmModel.hh.

585{
586 SetCurrentCouple(couple);
588 ((*elmSelectors)[couple->GetIndex()])->SelectRandomAtom(kinEnergy) :
589 SelectRandomAtom(pBaseMaterial,part,kinEnergy,cutEnergy,maxEnergy);
590 fCurrentIsotope = nullptr;
591 return fCurrentElement;
592}
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:580

References G4VEmModel::elmSelectors, G4VEmModel::fCurrentElement, G4VEmModel::fCurrentIsotope, G4MaterialCutsCouple::GetIndex(), G4VEmModel::nSelectors, G4VEmModel::pBaseMaterial, G4VEmModel::SelectRandomAtom(), and G4VEmModel::SetCurrentCouple().

Referenced by G4AdjointBremsstrahlungModel::RapidSampleSecondaries(), G4LivermoreBremsstrahlungModel::SampleSecondaries(), G4JAEAElasticScatteringModel::SampleSecondaries(), G4JAEAPolarizedElasticScatteringModel::SampleSecondaries(), G4LivermoreComptonModel::SampleSecondaries(), G4LivermoreNuclearGammaConversionModel::SampleSecondaries(), G4LivermorePhotoElectricModel::SampleSecondaries(), G4LivermorePolarizedComptonModel::SampleSecondaries(), G4LivermorePolarizedGammaConversionModel::SampleSecondaries(), G4LivermorePolarizedRayleighModel::SampleSecondaries(), G4LivermoreRayleighModel::SampleSecondaries(), G4LowEPComptonModel::SampleSecondaries(), G4LowEPPolarizedComptonModel::SampleSecondaries(), G4PenelopePhotoElectricModel::SampleSecondaries(), G4MuBremsstrahlungModel::SampleSecondaries(), G4MuPairProductionModel::SampleSecondaries(), G4hCoulombScatteringModel::SampleSecondaries(), G4KleinNishinaModel::SampleSecondaries(), G4PEEffectFluoModel::SampleSecondaries(), G4VEmModel::SelectRandomAtom(), and G4VEmModel::SelectTargetAtom().

◆ SelectRandomAtomNumber()

G4int G4VEmModel::SelectRandomAtomNumber ( const G4Material mat)
inherited

Definition at line 297 of file G4VEmModel.cc.

298{
299 // this algorith assumes that cross section is proportional to
300 // number electrons multiplied by number of atoms
301 const size_t nn = mat->GetNumberOfElements();
302 fCurrentElement = mat->GetElement(0);
303 if(1 < nn) {
304 const G4double* at = mat->GetVecNbOfAtomsPerVolume();
306 for(size_t i=0; i<nn; ++i) {
307 tot -= at[i];
308 if(tot <= 0.0) {
309 fCurrentElement = mat->GetElement(i);
310 break;
311 }
312 }
313 }
314 return fCurrentElement->GetZasInt();
315}
G4int GetZasInt() const
Definition: G4Element.hh:132
G4double GetTotNbOfAtomsPerVolume() const
Definition: G4Material.hh:205
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:202

References G4VEmModel::fCurrentElement, G4UniformRand, G4Material::GetElement(), G4Material::GetNumberOfElements(), G4Material::GetTotNbOfAtomsPerVolume(), G4Material::GetVecNbOfAtomsPerVolume(), G4Element::GetZasInt(), and G4InuclParticleNames::nn.

Referenced by G4AtimaEnergyLossModel::SampleSecondaries(), G4BetheBlochModel::SampleSecondaries(), G4BraggIonModel::SampleSecondaries(), G4BraggModel::SampleSecondaries(), G4ICRU73QOModel::SampleSecondaries(), SampleSecondaries(), G4MollerBhabhaModel::SampleSecondaries(), and G4IonParametrisedLossModel::SampleSecondaries().

◆ SelectTargetAtom()

const G4Element * G4VEmModel::SelectTargetAtom ( const G4MaterialCutsCouple couple,
const G4ParticleDefinition part,
G4double  kineticEnergy,
G4double  logKineticEnergy,
G4double  cutEnergy = 0.0,
G4double  maxEnergy = DBL_MAX 
)
inlineinherited

◆ SetActivationHighEnergyLimit()

void G4VEmModel::SetActivationHighEnergyLimit ( G4double  val)
inlineinherited

◆ SetActivationLowEnergyLimit()

void G4VEmModel::SetActivationLowEnergyLimit ( G4double  val)
inlineinherited

◆ SetAngularDistribution()

void G4VEmModel::SetAngularDistribution ( G4VEmAngularDistribution p)
inlineinherited

Definition at line 628 of file G4VEmModel.hh.

629{
630 if(p != anglModel) {
631 delete anglModel;
632 anglModel = p;
633 }
634}

References G4VEmModel::anglModel.

Referenced by G4EmLivermorePhysics::ConstructProcess(), G4EmLowEPPhysics::ConstructProcess(), G4EmStandardPhysics::ConstructProcess(), G4EmStandardPhysics_option3::ConstructProcess(), G4EmStandardPhysics_option4::ConstructProcess(), G4EmStandardPhysicsSS::ConstructProcess(), G4BetheHeitlerModel::G4BetheHeitlerModel(), G4DNABornIonisationModel1::G4DNABornIonisationModel1(), G4DNABornIonisationModel2::G4DNABornIonisationModel2(), G4DNAEmfietzoglouIonisationModel::G4DNAEmfietzoglouIonisationModel(), G4DNARuddIonisationExtendedModel::G4DNARuddIonisationExtendedModel(), G4DNARuddIonisationModel::G4DNARuddIonisationModel(), G4eBremParametrizedModel::G4eBremParametrizedModel(), G4eBremsstrahlungRelModel::G4eBremsstrahlungRelModel(), G4IonParametrisedLossModel::G4IonParametrisedLossModel(), G4LivermoreBremsstrahlungModel::G4LivermoreBremsstrahlungModel(), G4LivermoreIonisationModel::G4LivermoreIonisationModel(), G4LivermorePhotoElectricModel::G4LivermorePhotoElectricModel(), G4LivermoreRayleighModel::G4LivermoreRayleighModel(), G4MicroElecInelasticModel::G4MicroElecInelasticModel(), G4MicroElecInelasticModel_new::G4MicroElecInelasticModel_new(), G4MuBremsstrahlungModel::G4MuBremsstrahlungModel(), G4MuPairProductionModel::G4MuPairProductionModel(), G4PAIModel::G4PAIModel(), G4PAIPhotModel::G4PAIPhotModel(), G4PairProductionRelModel::G4PairProductionRelModel(), G4PEEffectFluoModel::G4PEEffectFluoModel(), G4SeltzerBergerModel::G4SeltzerBergerModel(), G4AtimaEnergyLossModel::Initialise(), G4BetheBlochModel::Initialise(), G4BraggIonModel::Initialise(), G4BraggModel::Initialise(), G4ICRU73QOModel::Initialise(), Initialise(), and G4MollerBhabhaModel::Initialise().

◆ SetAngularGeneratorFlag()

void G4VEmModel::SetAngularGeneratorFlag ( G4bool  val)
inlineinherited

Definition at line 725 of file G4VEmModel.hh.

726{
728}
G4bool useAngularGenerator
Definition: G4VEmModel.hh:461

References G4VEmModel::useAngularGenerator.

Referenced by G4VEnergyLossProcess::PreparePhysicsTable().

◆ SetChargeSquareRatio()

void G4LindhardSorensenIonModel::SetChargeSquareRatio ( G4double  val)
inlineprotected

Definition at line 199 of file G4LindhardSorensenIonModel.hh.

200{
201 chargeSquare = val;
202}

References chargeSquare.

◆ SetCrossSectionTable()

void G4VEmModel::SetCrossSectionTable ( G4PhysicsTable p,
G4bool  isLocal 
)
inherited

Definition at line 455 of file G4VEmModel.cc.

456{
457 if(p != xSectionTable) {
458 if(xSectionTable != nullptr && localTable) {
460 delete xSectionTable;
461 }
462 xSectionTable = p;
463 }
464 localTable = isLocal;
465}
void clearAndDestroy()
G4bool localTable
Definition: G4VEmModel.hh:459

References G4PhysicsTable::clearAndDestroy(), G4VEmModel::localTable, and G4VEmModel::xSectionTable.

Referenced by G4VMultipleScattering::BuildPhysicsTable().

◆ SetCurrentCouple()

void G4VEmModel::SetCurrentCouple ( const G4MaterialCutsCouple ptr)
inlineinherited

Definition at line 472 of file G4VEmModel.hh.

473{
474 if(fCurrentCouple != ptr) {
475 fCurrentCouple = ptr;
477 pBaseMaterial = ptr->GetMaterial();
478 pFactor = 1.0;
479 if(useBaseMaterials) {
480 basedCoupleIndex = (*theDensityIdx)[currentCoupleIndex];
481 if(nullptr != pBaseMaterial->GetBaseMaterial())
483 pFactor = (*theDensityFactor)[currentCoupleIndex];
484 }
485 }
486}
const G4Material * GetBaseMaterial() const
Definition: G4Material.hh:229
G4bool useBaseMaterials
Definition: G4VEmModel.hh:462
size_t currentCoupleIndex
Definition: G4VEmModel.hh:448
size_t basedCoupleIndex
Definition: G4VEmModel.hh:449

References G4VEmModel::basedCoupleIndex, G4VEmModel::currentCoupleIndex, G4VEmModel::fCurrentCouple, G4Material::GetBaseMaterial(), G4MaterialCutsCouple::GetIndex(), G4MaterialCutsCouple::GetMaterial(), G4VEmModel::pBaseMaterial, G4VEmModel::pFactor, and G4VEmModel::useBaseMaterials.

Referenced by G4VMultipleScattering::AlongStepGetPhysicalInteractionLength(), G4EmMultiModel::ComputeCrossSectionPerAtom(), G4VEmModel::ComputeDEDX(), G4TablesForExtrapolator::ComputeTrasportXS(), G4UrbanAdjointMscModel::ComputeTruePathLengthLimit(), G4GoudsmitSaundersonMscModel::ComputeTruePathLengthLimit(), G4UrbanMscModel::ComputeTruePathLengthLimit(), G4VEmModel::CrossSection(), G4AdjointPhotoElectricModel::DefineCurrentMaterialAndElectronEnergy(), G4WentzelVIModel::DefineMaterial(), G4WentzelVIRelModel::DefineMaterial(), G4EmCalculator::GetCrossSectionPerVolume(), G4LivermoreGammaConversion5DModel::Initialise(), G4VEmModel::InitialiseElementSelectors(), G4PEEffectFluoModel::SampleSecondaries(), G4EmMultiModel::SampleSecondaries(), G4VEnergyLossProcess::SelectModel(), G4VEmProcess::SelectModel(), G4VEmModel::SelectRandomAtom(), G4VEmModel::SelectTargetAtom(), and G4VEmModel::Value().

◆ SetCurrentElement()

void G4VEmModel::SetCurrentElement ( const G4Element elm)
inlineprotectedinherited

◆ SetDeexcitationFlag()

void G4VEmModel::SetDeexcitationFlag ( G4bool  val)
inlineinherited

Definition at line 823 of file G4VEmModel.hh.

824{
825 flagDeexcitation = val;
826}

References G4VEmModel::flagDeexcitation.

Referenced by G4DNABornIonisationModel1::G4DNABornIonisationModel1(), G4DNABornIonisationModel2::G4DNABornIonisationModel2(), G4DNACPA100IonisationModel::G4DNACPA100IonisationModel(), G4DNAEmfietzoglouIonisationModel::G4DNAEmfietzoglouIonisationModel(), G4DNARelativisticIonisationModel::G4DNARelativisticIonisationModel(), G4DNARuddIonisationExtendedModel::G4DNARuddIonisationExtendedModel(), G4DNARuddIonisationModel::G4DNARuddIonisationModel(), G4KleinNishinaModel::G4KleinNishinaModel(), G4LEPTSIonisationModel::G4LEPTSIonisationModel(), G4LivermoreComptonModel::G4LivermoreComptonModel(), G4LivermorePhotoElectricModel::G4LivermorePhotoElectricModel(), G4LivermorePolarizedComptonModel::G4LivermorePolarizedComptonModel(), G4LowEPComptonModel::G4LowEPComptonModel(), G4LowEPPolarizedComptonModel::G4LowEPPolarizedComptonModel(), G4MicroElecInelasticModel::G4MicroElecInelasticModel(), G4MicroElecInelasticModel_new::G4MicroElecInelasticModel_new(), G4PEEffectFluoModel::G4PEEffectFluoModel(), G4PenelopeBremsstrahlungModel::G4PenelopeBremsstrahlungModel(), G4PenelopeComptonModel::G4PenelopeComptonModel(), G4PenelopeIonisationModel::G4PenelopeIonisationModel(), G4PenelopePhotoElectricModel::G4PenelopePhotoElectricModel(), G4AtimaEnergyLossModel::Initialise(), G4BetheBlochModel::Initialise(), G4BraggIonModel::Initialise(), G4BraggModel::Initialise(), G4ICRU73QOModel::Initialise(), and Initialise().

◆ SetElementSelectors()

void G4VEmModel::SetElementSelectors ( std::vector< G4EmElementSelector * > *  p)
inlineinherited

Definition at line 852 of file G4VEmModel.hh.

853{
854 if(p != elmSelectors) {
855 elmSelectors = p;
856 nSelectors = (nullptr != elmSelectors) ? G4int(elmSelectors->size()) : 0;
857 localElmSelectors = false;
858 }
859}
G4bool localElmSelectors
Definition: G4VEmModel.hh:460

References G4VEmModel::elmSelectors, G4VEmModel::localElmSelectors, and G4VEmModel::nSelectors.

Referenced by G4eDPWACoulombScatteringModel::InitialiseLocal(), G4eSingleCoulombScatteringModel::InitialiseLocal(), G4PAIModel::InitialiseLocal(), G4PAIPhotModel::InitialiseLocal(), G4JAEAElasticScatteringModel::InitialiseLocal(), G4JAEAPolarizedElasticScatteringModel::InitialiseLocal(), G4LivermoreComptonModel::InitialiseLocal(), G4LivermoreNuclearGammaConversionModel::InitialiseLocal(), G4LivermorePolarizedComptonModel::InitialiseLocal(), G4LivermorePolarizedGammaConversionModel::InitialiseLocal(), G4LivermorePolarizedRayleighModel::InitialiseLocal(), G4LivermoreRayleighModel::InitialiseLocal(), G4LowEPComptonModel::InitialiseLocal(), G4LowEPPolarizedComptonModel::InitialiseLocal(), G4PenelopePhotoElectricModel::InitialiseLocal(), G4MuBremsstrahlungModel::InitialiseLocal(), G4MuPairProductionModel::InitialiseLocal(), G4BetheHeitlerModel::InitialiseLocal(), G4eBremParametrizedModel::InitialiseLocal(), G4eBremsstrahlungRelModel::InitialiseLocal(), G4eCoulombScatteringModel::InitialiseLocal(), G4hCoulombScatteringModel::InitialiseLocal(), G4KleinNishinaCompton::InitialiseLocal(), G4KleinNishinaModel::InitialiseLocal(), and G4PairProductionRelModel::InitialiseLocal().

◆ SetFluctuationFlag()

void G4VEmModel::SetFluctuationFlag ( G4bool  val)
inlineinherited

Definition at line 732 of file G4VEmModel.hh.

733{
734 lossFlucFlag = val;
735}
G4bool lossFlucFlag
Definition: G4VEmModel.hh:450

References G4VEmModel::lossFlucFlag.

Referenced by G4EmCalculator::ComputeNuclearDEDX().

◆ SetForceBuildTable()

void G4VEmModel::SetForceBuildTable ( G4bool  val)
inlineinherited

Definition at line 830 of file G4VEmModel.hh.

831{
833}

References G4VEmModel::flagForceBuildTable.

◆ SetHighEnergyLimit()

void G4VEmModel::SetHighEnergyLimit ( G4double  val)
inlineinherited

Definition at line 767 of file G4VEmModel.hh.

768{
769 highLimit = val;
770}

References G4VEmModel::highLimit.

Referenced by G4EmModelActivator::ActivateEmOptions(), G4EmModelActivator::ActivatePAI(), LBE::ConstructEM(), G4EmDNAPhysics_option2::ConstructProcess(), G4EmDNAPhysics_option5::ConstructProcess(), G4EmDNAPhysics_option7::ConstructProcess(), G4EmDNAPhysics_stationary::ConstructProcess(), G4EmLivermorePhysics::ConstructProcess(), G4EmLowEPPhysics::ConstructProcess(), G4EmPenelopePhysics::ConstructProcess(), G4EmStandardPhysics::ConstructProcess(), G4EmStandardPhysics_option1::ConstructProcess(), G4EmStandardPhysics_option2::ConstructProcess(), G4EmStandardPhysics_option4::ConstructProcess(), G4EmStandardPhysicsGS::ConstructProcess(), G4EmStandardPhysicsWVI::ConstructProcess(), G4BraggIonModel::G4BraggIonModel(), G4BraggModel::G4BraggModel(), G4DNAChampionElasticModel::G4DNAChampionElasticModel(), G4DNACPA100ElasticModel::G4DNACPA100ElasticModel(), G4DNACPA100ExcitationModel::G4DNACPA100ExcitationModel(), G4DNACPA100IonisationModel::G4DNACPA100IonisationModel(), G4DNAELSEPAElasticModel::G4DNAELSEPAElasticModel(), G4DNAEmfietzoglouExcitationModel::G4DNAEmfietzoglouExcitationModel(), G4DNAEmfietzoglouIonisationModel::G4DNAEmfietzoglouIonisationModel(), G4DNAIonElasticModel::G4DNAIonElasticModel(), G4DNAMeltonAttachmentModel::G4DNAMeltonAttachmentModel(), G4DNASancheExcitationModel::G4DNASancheExcitationModel(), G4DNAScreenedRutherfordElasticModel::G4DNAScreenedRutherfordElasticModel(), G4DNATransformElectronModel::G4DNATransformElectronModel(), G4DNAUeharaScreenedRutherfordElasticModel::G4DNAUeharaScreenedRutherfordElasticModel(), G4eDPWACoulombScatteringModel::G4eDPWACoulombScatteringModel(), G4ICRU73QOModel::G4ICRU73QOModel(), G4MicroElecElasticModel::G4MicroElecElasticModel(), G4MicroElecElasticModel_new::G4MicroElecElasticModel_new(), G4PenelopeAnnihilationModel::G4PenelopeAnnihilationModel(), G4PenelopeBremsstrahlungModel::G4PenelopeBremsstrahlungModel(), G4PenelopeComptonModel::G4PenelopeComptonModel(), G4PenelopeGammaConversionModel::G4PenelopeGammaConversionModel(), G4PenelopeIonisationModel::G4PenelopeIonisationModel(), G4PenelopePhotoElectricModel::G4PenelopePhotoElectricModel(), G4PenelopeRayleighModel::G4PenelopeRayleighModel(), G4PenelopeRayleighModelMI::G4PenelopeRayleighModelMI(), G4TDNAOneStepThermalizationModel< MODEL >::G4TDNAOneStepThermalizationModel(), G4XrayRayleighModel::G4XrayRayleighModel(), G4VLEPTSModel::Init(), G4DNAChampionElasticModel::Initialise(), G4DNACPA100ElasticModel::Initialise(), G4DNADingfelderChargeDecreaseModel::Initialise(), G4DNADingfelderChargeIncreaseModel::Initialise(), G4DNAMeltonAttachmentModel::Initialise(), G4DNAMillerGreenExcitationModel::Initialise(), G4DNARuddIonisationExtendedModel::Initialise(), G4DNARuddIonisationModel::Initialise(), G4DNASancheExcitationModel::Initialise(), G4MicroElecElasticModel::Initialise(), G4MicroElecInelasticModel::Initialise(), G4MicroElecInelasticModel_new::Initialise(), G4DNABornExcitationModel1::Initialise(), G4DNABornExcitationModel2::Initialise(), G4DNABornIonisationModel1::Initialise(), G4DNABornIonisationModel2::Initialise(), G4DNAELSEPAElasticModel::Initialise(), G4DNAModelInterface::Initialise(), G4DNAIonElasticModel::Initialise(), G4hhIonisation::InitialiseEnergyLossProcess(), G4mplIonisation::InitialiseEnergyLossProcess(), G4ePairProduction::InitialiseEnergyLossProcess(), G4MuBremsstrahlung::InitialiseEnergyLossProcess(), G4MuIonisation::InitialiseEnergyLossProcess(), G4MuPairProduction::InitialiseEnergyLossProcess(), G4PolarizedBremsstrahlung::InitialiseEnergyLossProcess(), G4PolarizedIonisation::InitialiseEnergyLossProcess(), G4eBremsstrahlung::InitialiseEnergyLossProcess(), G4eIonisation::InitialiseEnergyLossProcess(), G4hIonisation::InitialiseEnergyLossProcess(), G4ionIonisation::InitialiseEnergyLossProcess(), G4DNAAttachment::InitialiseProcess(), G4DNAChargeDecrease::InitialiseProcess(), G4DNAChargeIncrease::InitialiseProcess(), G4DNADissociation::InitialiseProcess(), G4DNAElastic::InitialiseProcess(), G4DNAExcitation::InitialiseProcess(), G4DNAIonisation::InitialiseProcess(), G4DNAPlasmonExcitation::InitialiseProcess(), G4DNAPositronium::InitialiseProcess(), G4DNARotExcitation::InitialiseProcess(), G4DNAVibExcitation::InitialiseProcess(), G4PolarizedCompton::InitialiseProcess(), G4PolarizedGammaConversion::InitialiseProcess(), G4PolarizedPhotoElectric::InitialiseProcess(), G4ComptonScattering::InitialiseProcess(), G4CoulombScattering::InitialiseProcess(), G4eplusAnnihilation::InitialiseProcess(), G4GammaConversion::InitialiseProcess(), G4PhotoElectricEffect::InitialiseProcess(), G4VEmProcess::PreparePhysicsTable(), G4VEnergyLossProcess::PreparePhysicsTable(), G4VMultipleScattering::PreparePhysicsTable(), G4DNAUeharaScreenedRutherfordElasticModel::SelectHighEnergyLimit(), G4VEmAdjointModel::SetHighEnergyLimit(), G4DNAELSEPAElasticModel::SetMaximumEnergy(), G4mplIonisationModel::SetParticle(), G4mplIonisationWithDeltaModel::SetParticle(), and G4EmConfigurator::UpdateModelEnergyRange().

◆ SetLocked()

void G4VEmModel::SetLocked ( G4bool  val)
inlineinherited

◆ SetLowEnergyLimit()

void G4VEmModel::SetLowEnergyLimit ( G4double  val)
inlineinherited

Definition at line 774 of file G4VEmModel.hh.

775{
776 lowLimit = val;
777}

References G4VEmModel::lowLimit.

Referenced by G4EmModelActivator::ActivatePAI(), G4EmDNAPhysics_option2::ConstructProcess(), G4EmDNAPhysics_option5::ConstructProcess(), G4EmDNAPhysics_stationary::ConstructProcess(), G4EmLivermorePhysics::ConstructProcess(), G4EmLowEPPhysics::ConstructProcess(), G4EmPenelopePhysics::ConstructProcess(), G4EmStandardPhysics::ConstructProcess(), G4EmStandardPhysics_option1::ConstructProcess(), G4EmStandardPhysics_option2::ConstructProcess(), G4EmStandardPhysics_option3::ConstructProcess(), G4EmStandardPhysics_option4::ConstructProcess(), G4EmStandardPhysicsGS::ConstructProcess(), G4EmStandardPhysicsWVI::ConstructProcess(), G4DNASancheExcitationModel::ExtendLowEnergyLimit(), G4AtimaEnergyLossModel::G4AtimaEnergyLossModel(), G4BetheBlochModel::G4BetheBlochModel(), G4BetheHeitler5DModel::G4BetheHeitler5DModel(), G4DNAChampionElasticModel::G4DNAChampionElasticModel(), G4DNACPA100ElasticModel::G4DNACPA100ElasticModel(), G4DNACPA100ExcitationModel::G4DNACPA100ExcitationModel(), G4DNACPA100IonisationModel::G4DNACPA100IonisationModel(), G4DNAELSEPAElasticModel::G4DNAELSEPAElasticModel(), G4DNAEmfietzoglouExcitationModel::G4DNAEmfietzoglouExcitationModel(), G4DNAEmfietzoglouIonisationModel::G4DNAEmfietzoglouIonisationModel(), G4DNAIonElasticModel::G4DNAIonElasticModel(), G4DNAMeltonAttachmentModel::G4DNAMeltonAttachmentModel(), G4DNASancheExcitationModel::G4DNASancheExcitationModel(), G4DNAScreenedRutherfordElasticModel::G4DNAScreenedRutherfordElasticModel(), G4DNATransformElectronModel::G4DNATransformElectronModel(), G4DNAUeharaScreenedRutherfordElasticModel::G4DNAUeharaScreenedRutherfordElasticModel(), G4DummyModel::G4DummyModel(), G4eBremParametrizedModel::G4eBremParametrizedModel(), G4eBremsstrahlungRelModel::G4eBremsstrahlungRelModel(), G4eDPWACoulombScatteringModel::G4eDPWACoulombScatteringModel(), G4LivermoreBremsstrahlungModel::G4LivermoreBremsstrahlungModel(), G4MicroElecElasticModel::G4MicroElecElasticModel(), G4MicroElecElasticModel_new::G4MicroElecElasticModel_new(), G4SeltzerBergerModel::G4SeltzerBergerModel(), G4TDNAOneStepThermalizationModel< MODEL >::G4TDNAOneStepThermalizationModel(), G4VLEPTSModel::Init(), G4DNAChampionElasticModel::Initialise(), G4DNACPA100ElasticModel::Initialise(), G4DNADingfelderChargeDecreaseModel::Initialise(), G4DNADingfelderChargeIncreaseModel::Initialise(), G4DNAMeltonAttachmentModel::Initialise(), G4DNAMillerGreenExcitationModel::Initialise(), G4DNARuddIonisationExtendedModel::Initialise(), G4DNARuddIonisationModel::Initialise(), G4MicroElecElasticModel::Initialise(), G4MicroElecInelasticModel::Initialise(), G4MicroElecInelasticModel_new::Initialise(), G4BetheHeitler5DModel::Initialise(), G4DNABornExcitationModel1::Initialise(), G4DNABornExcitationModel2::Initialise(), G4DNABornIonisationModel1::Initialise(), G4DNABornIonisationModel2::Initialise(), G4DNAELSEPAElasticModel::Initialise(), G4DNAModelInterface::Initialise(), G4DNAIonElasticModel::Initialise(), G4hhIonisation::InitialiseEnergyLossProcess(), G4mplIonisation::InitialiseEnergyLossProcess(), G4ePairProduction::InitialiseEnergyLossProcess(), G4MuBremsstrahlung::InitialiseEnergyLossProcess(), G4MuIonisation::InitialiseEnergyLossProcess(), G4MuPairProduction::InitialiseEnergyLossProcess(), G4PolarizedBremsstrahlung::InitialiseEnergyLossProcess(), G4PolarizedIonisation::InitialiseEnergyLossProcess(), G4eBremsstrahlung::InitialiseEnergyLossProcess(), G4eIonisation::InitialiseEnergyLossProcess(), G4hIonisation::InitialiseEnergyLossProcess(), G4ionIonisation::InitialiseEnergyLossProcess(), G4DNAAttachment::InitialiseProcess(), G4DNAChargeDecrease::InitialiseProcess(), G4DNAChargeIncrease::InitialiseProcess(), G4DNADissociation::InitialiseProcess(), G4DNAElastic::InitialiseProcess(), G4DNAExcitation::InitialiseProcess(), G4DNAIonisation::InitialiseProcess(), G4DNAPlasmonExcitation::InitialiseProcess(), G4DNAPositronium::InitialiseProcess(), G4DNARotExcitation::InitialiseProcess(), G4DNAVibExcitation::InitialiseProcess(), G4PolarizedCompton::InitialiseProcess(), G4PolarizedGammaConversion::InitialiseProcess(), G4PolarizedPhotoElectric::InitialiseProcess(), G4ComptonScattering::InitialiseProcess(), G4CoulombScattering::InitialiseProcess(), G4eplusAnnihilation::InitialiseProcess(), G4GammaConversion::InitialiseProcess(), G4PhotoElectricEffect::InitialiseProcess(), G4VEmAdjointModel::SetLowEnergyLimit(), G4mplIonisationModel::SetParticle(), G4mplIonisationWithDeltaModel::SetParticle(), and G4EmConfigurator::UpdateModelEnergyRange().

◆ SetLPMFlag()

void G4VEmModel::SetLPMFlag ( G4bool  val)
inlineinherited

◆ SetMasterThread()

void G4VEmModel::SetMasterThread ( G4bool  val)
inlineinherited

◆ SetParticle()

void G4LindhardSorensenIonModel::SetParticle ( const G4ParticleDefinition p)
inlineprivate

Definition at line 182 of file G4LindhardSorensenIonModel.hh.

183{
184 if(particle != p) {
185 particle = p;
187 }
188}

References particle, and SetupParameters().

Referenced by Initialise(), and MaxSecondaryEnergy().

◆ SetParticleChange()

void G4VEmModel::SetParticleChange ( G4VParticleChange p,
G4VEmFluctuationModel f = nullptr 
)
inherited

◆ SetPolarAngleLimit()

void G4VEmModel::SetPolarAngleLimit ( G4double  val)
inlineinherited

◆ SetSecondaryThreshold()

void G4VEmModel::SetSecondaryThreshold ( G4double  val)
inlineinherited

◆ SetTripletModel()

void G4VEmModel::SetTripletModel ( G4VEmModel p)
inlineinherited

Definition at line 645 of file G4VEmModel.hh.

646{
647 if(p != fTripletModel) {
648 delete fTripletModel;
649 fTripletModel = p;
650 }
651}

References G4VEmModel::fTripletModel.

Referenced by G4eplusTo2GammaOKVIModel::G4eplusTo2GammaOKVIModel().

◆ SetupForMaterial()

void G4VEmModel::SetupForMaterial ( const G4ParticleDefinition ,
const G4Material ,
G4double  kineticEnergy 
)
virtualinherited

◆ SetupParameters()

void G4LindhardSorensenIonModel::SetupParameters ( )
private

Definition at line 140 of file G4LindhardSorensenIonModel.cc.

141{
145 Zin = G4lrint(std::abs(charge));
149 const G4double aMag =
152 magMoment2 = magmom*magmom - 1.0;
153 G4double x = 0.8426*CLHEP::GeV;
154 if(spin == 0.0 && mass < GeV) { x = 0.736*CLHEP::GeV; }
155 else if (Zin > 1) { x /= nist->GetA27(Zin); }
156
158 tlimit = 2.0/formfact;
159}
static constexpr double GeV
Definition: G4SIunits.hh:203
G4double GetA27(G4int Z) const
G4double GetPDGMagneticMoment() const
G4double GetPDGCharge() const
G4double inveplus
Definition: G4VEmModel.hh:431
static constexpr double eplus
static constexpr double proton_mass_c2
static constexpr double GeV
static constexpr double hbar_Planck
static constexpr double c_squared
int G4lrint(double ad)
Definition: templates.hh:134

References CLHEP::c_squared, charge, chargeSquare, CLHEP::electron_mass_c2, CLHEP::eplus, eRatio, formfact, G4lrint(), G4NistManager::GetA27(), G4ParticleDefinition::GetPDGCharge(), G4ParticleDefinition::GetPDGMagneticMoment(), G4ParticleDefinition::GetPDGMass(), G4ParticleDefinition::GetPDGSpin(), CLHEP::GeV, GeV, CLHEP::hbar_Planck, G4VEmModel::inveplus, magMoment2, mass, nist, particle, pRatio, CLHEP::proton_mass_c2, spin, tlimit, and Zin.

Referenced by SetParticle().

◆ SetUseBaseMaterials()

void G4VEmModel::SetUseBaseMaterials ( G4bool  val)
inlineinherited

◆ StartTracking()

void G4VEmModel::StartTracking ( G4Track )
virtualinherited

◆ UseAngularGeneratorFlag()

G4bool G4VEmModel::UseAngularGeneratorFlag ( ) const
inlineinherited

◆ UseBaseMaterials()

G4bool G4VEmModel::UseBaseMaterials ( ) const
inlineinherited

Definition at line 760 of file G4VEmModel.hh.

761{
762 return useBaseMaterials;
763}

References G4VEmModel::useBaseMaterials.

◆ Value()

G4double G4VEmModel::Value ( const G4MaterialCutsCouple couple,
const G4ParticleDefinition p,
G4double  kineticEnergy 
)
virtualinherited

Field Documentation

◆ anglModel

G4VEmAngularDistribution* G4VEmModel::anglModel = nullptr
privateinherited

◆ basedCoupleIndex

size_t G4VEmModel::basedCoupleIndex = 0
protectedinherited

◆ charge

G4double G4LindhardSorensenIonModel::charge = 1.0
private

Definition at line 166 of file G4LindhardSorensenIonModel.hh.

Referenced by CorrectionsAlongStep(), and SetupParameters().

◆ chargeSquare

G4double G4LindhardSorensenIonModel::chargeSquare = 1.0
private

◆ corr

G4EmCorrections* G4LindhardSorensenIonModel::corr
private

◆ currentCoupleIndex

size_t G4VEmModel::currentCoupleIndex = 0
protectedinherited

Definition at line 448 of file G4VEmModel.hh.

Referenced by G4VEmModel::SetCurrentCouple().

◆ elmSelectors

std::vector<G4EmElementSelector*>* G4VEmModel::elmSelectors = nullptr
privateinherited

◆ eMaxActive

G4double G4VEmModel::eMaxActive = DBL_MAX
privateinherited

◆ eMinActive

G4double G4VEmModel::eMinActive = 0.0
privateinherited

◆ eRatio

G4double G4LindhardSorensenIonModel::eRatio = 0.0
private

Definition at line 167 of file G4LindhardSorensenIonModel.hh.

Referenced by MaxSecondaryEnergy(), and SetupParameters().

◆ fact

std::vector< G4float > * G4LindhardSorensenIonModel::fact = {nullptr}
staticprivate

Definition at line 150 of file G4LindhardSorensenIonModel.hh.

◆ fBBModel

G4BetheBlochModel* G4LindhardSorensenIonModel::fBBModel
private

◆ fBraggModel

G4BraggIonModel* G4LindhardSorensenIonModel::fBraggModel
private

◆ fCurrentCouple

const G4MaterialCutsCouple* G4VEmModel::fCurrentCouple = nullptr
privateinherited

Definition at line 416 of file G4VEmModel.hh.

Referenced by G4VEmModel::CurrentCouple(), and G4VEmModel::SetCurrentCouple().

◆ fCurrentElement

const G4Element* G4VEmModel::fCurrentElement = nullptr
privateinherited

◆ fCurrentIsotope

const G4Isotope* G4VEmModel::fCurrentIsotope = nullptr
privateinherited

◆ fElementData

G4ElementData* G4VEmModel::fElementData = nullptr
protectedinherited

◆ fElimit

G4double G4LindhardSorensenIonModel::fElimit
private

◆ fEmManager

G4LossTableManager* G4VEmModel::fEmManager
privateinherited

Definition at line 420 of file G4VEmModel.hh.

Referenced by G4VEmModel::G4VEmModel(), and G4VEmModel::~G4VEmModel().

◆ fIonData

G4IonICRU73Data * G4LindhardSorensenIonModel::fIonData = nullptr
staticprivate

Definition at line 148 of file G4LindhardSorensenIonModel.hh.

Referenced by CorrectionsAlongStep(), and Initialise().

◆ flagDeexcitation

G4bool G4VEmModel::flagDeexcitation = false
privateinherited

Definition at line 455 of file G4VEmModel.hh.

Referenced by G4VEmModel::DeexcitationFlag(), and G4VEmModel::SetDeexcitationFlag().

◆ flagForceBuildTable

G4bool G4VEmModel::flagForceBuildTable = false
privateinherited

◆ flucModel

G4VEmFluctuationModel* G4VEmModel::flucModel = nullptr
privateinherited

◆ formfact

G4double G4LindhardSorensenIonModel::formfact = 0.0
private

Definition at line 169 of file G4LindhardSorensenIonModel.hh.

Referenced by SampleSecondaries(), and SetupParameters().

◆ fParticleChange

G4ParticleChangeForLoss* G4LindhardSorensenIonModel::fParticleChange
private

◆ fTripletModel

G4VEmModel* G4VEmModel::fTripletModel = nullptr
privateinherited

◆ highLimit

G4double G4VEmModel::highLimit
privateinherited

◆ inveplus

G4double G4VEmModel::inveplus
protectedinherited

◆ isLocked

G4bool G4VEmModel::isLocked = false
privateinherited

◆ isMaster

G4bool G4VEmModel::isMaster = true
privateinherited

◆ localElmSelectors

G4bool G4VEmModel::localElmSelectors = true
privateinherited

Definition at line 460 of file G4VEmModel.hh.

Referenced by G4VEmModel::SetElementSelectors(), and G4VEmModel::~G4VEmModel().

◆ localTable

G4bool G4VEmModel::localTable = true
privateinherited

Definition at line 459 of file G4VEmModel.hh.

Referenced by G4VEmModel::SetCrossSectionTable(), and G4VEmModel::~G4VEmModel().

◆ lossFlucFlag

G4bool G4VEmModel::lossFlucFlag = true
protectedinherited

◆ lowLimit

G4double G4VEmModel::lowLimit
privateinherited

◆ lsdata

G4LindhardSorensenData * G4LindhardSorensenIonModel::lsdata = nullptr
staticprivate

Definition at line 149 of file G4LindhardSorensenIonModel.hh.

Referenced by CorrectionsAlongStep(), and Initialise().

◆ magMoment2

G4double G4LindhardSorensenIonModel::magMoment2 = 0.0
private

Definition at line 164 of file G4LindhardSorensenIonModel.hh.

Referenced by SampleSecondaries(), and SetupParameters().

◆ mass

G4double G4LindhardSorensenIonModel::mass = 0.0
private

◆ MAXZION

const G4int G4LindhardSorensenIonModel::MAXZION = 93
staticprivate

Definition at line 146 of file G4LindhardSorensenIonModel.hh.

◆ name

const G4String G4VEmModel::name
privateinherited

◆ nist

G4NistManager* G4LindhardSorensenIonModel::nist
private

Definition at line 156 of file G4LindhardSorensenIonModel.hh.

Referenced by G4LindhardSorensenIonModel(), and SetupParameters().

◆ nsec

G4int G4VEmModel::nsec = 5
privateinherited

Definition at line 444 of file G4VEmModel.hh.

Referenced by G4VEmModel::CrossSectionPerVolume(), and G4VEmModel::G4VEmModel().

◆ nSelectors

G4int G4VEmModel::nSelectors = 0
privateinherited

◆ particle

const G4ParticleDefinition* G4LindhardSorensenIonModel::particle
private

◆ pBaseMaterial

const G4Material* G4VEmModel::pBaseMaterial = nullptr
protectedinherited

◆ pFactor

G4double G4VEmModel::pFactor = 1.0
protectedinherited

◆ polarAngleLimit

G4double G4VEmModel::polarAngleLimit
privateinherited

Definition at line 441 of file G4VEmModel.hh.

Referenced by G4VEmModel::PolarAngleLimit(), and G4VEmModel::SetPolarAngleLimit().

◆ pParticleChange

G4VParticleChange* G4VEmModel::pParticleChange = nullptr
protectedinherited

◆ pRatio

G4double G4LindhardSorensenIonModel::pRatio = 1.0
private

◆ secondaryThreshold

G4double G4VEmModel::secondaryThreshold = DBL_MAX
privateinherited

◆ spin

G4double G4LindhardSorensenIonModel::spin = 0.0
private

Definition at line 163 of file G4LindhardSorensenIonModel.hh.

Referenced by SampleSecondaries(), and SetupParameters().

◆ theDensityFactor

const std::vector<G4double>* G4VEmModel::theDensityFactor = nullptr
protectedinherited

Definition at line 428 of file G4VEmModel.hh.

Referenced by G4VEmModel::G4VEmModel().

◆ theDensityIdx

const std::vector<G4int>* G4VEmModel::theDensityIdx = nullptr
protectedinherited

Definition at line 429 of file G4VEmModel.hh.

Referenced by G4VEmModel::G4VEmModel().

◆ theElectron

G4ParticleDefinition* G4LindhardSorensenIonModel::theElectron
private

Definition at line 153 of file G4LindhardSorensenIonModel.hh.

Referenced by G4LindhardSorensenIonModel(), and SampleSecondaries().

◆ theLPMflag

G4bool G4VEmModel::theLPMflag = false
privateinherited

Definition at line 454 of file G4VEmModel.hh.

Referenced by G4VEmModel::LPMFlag(), and G4VEmModel::SetLPMFlag().

◆ tlimit

G4double G4LindhardSorensenIonModel::tlimit = DBL_MAX
private

Definition at line 162 of file G4LindhardSorensenIonModel.hh.

Referenced by ComputeDEDXPerVolume(), and SetupParameters().

◆ twoln10

G4double G4LindhardSorensenIonModel::twoln10
private

Definition at line 170 of file G4LindhardSorensenIonModel.hh.

◆ useAngularGenerator

G4bool G4VEmModel::useAngularGenerator = false
privateinherited

◆ useBaseMaterials

G4bool G4VEmModel::useBaseMaterials = false
privateinherited

◆ xsec

std::vector<G4double> G4VEmModel::xsec
privateinherited

◆ xSectionTable

G4PhysicsTable* G4VEmModel::xSectionTable = nullptr
protectedinherited

◆ Zin

G4int G4LindhardSorensenIonModel::Zin = 1
private

Definition at line 160 of file G4LindhardSorensenIonModel.hh.

Referenced by CorrectionsAlongStep(), and SetupParameters().


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