Geant4  10.02.p03
G4VEmModel Class Referenceabstract

#include <G4VEmModel.hh>

Inherited by G4BetheBlochModel, G4BetheHeitlerModel, G4BoldyshevTripletModel, G4BraggIonModel, G4BraggModel, G4DiscreteScatteringModel, G4DNABornExcitationModel1, G4DNABornExcitationModel2, G4DNABornIonisationModel1, G4DNABornIonisationModel2, G4DNAChampionElasticModel, G4DNADingfelderChargeDecreaseModel, G4DNADingfelderChargeIncreaseModel, G4DNAEmfietzoglouExcitationModel, G4DNAEmfietzoglouIonisationModel, G4DNAIonElasticModel, G4DNAMeltonAttachmentModel, G4DNAMillerGreenExcitationModel, G4DNAOneStepThermalizationModel, G4DNARuddIonisationExtendedModel, G4DNARuddIonisationModel, G4DNASancheExcitationModel, G4DNAScreenedRutherfordElasticModel, G4DNATransformElectronModel, G4DNAUeharaScreenedRutherfordElasticModel, G4eBremParametrizedModel, G4eBremsstrahlungRelModel, G4eCoulombScatteringModel, G4eeToHadronsModel, G4eeToHadronsMultiModel, G4eeToTwoGammaModel, G4EmMultiModel, G4eSingleCoulombScatteringModel, G4hCoulombScatteringModel, G4ICRU49NuclearStoppingModel, G4ICRU73QOModel, G4IonCoulombScatteringModel, G4IonParametrisedLossModel, G4KleinNishinaCompton, G4KleinNishinaModel, G4LivermoreComptonModel, G4LivermoreComptonModifiedModel, G4LivermoreGammaConversionModel, G4LivermoreGammaConversionModelRC, G4LivermoreIonisationModel, G4LivermoreNuclearGammaConversionModel, G4LivermorePolarizedComptonModel, G4LivermorePolarizedGammaConversionModel, G4LivermorePolarizedPhotoElectricGDModel, G4LivermorePolarizedPhotoElectricModel, G4LivermorePolarizedRayleighModel, G4LivermoreRayleighModel, G4LowEPComptonModel, G4LowEPPolarizedComptonModel, G4MicroElecElasticModel, G4MicroElecInelasticModel, G4MollerBhabhaModel, G4mplIonisationModel, G4mplIonisationWithDeltaModel, G4MuBetheBlochModel, G4MuBremsstrahlungModel, G4MuElecElasticModel, G4MuElecInelasticModel, G4MuPairProductionModel, G4PAIModel, G4PAIPhotModel, G4PairProductionRelModel, G4PEEffectFluoModel, G4PenelopeAnnihilationModel, G4PenelopeComptonModel, G4PenelopeGammaConversionModel, G4PenelopeIonisationModel, G4PenelopePhotoElectricModel, G4PenelopeRayleighModel, G4VDNAPTBModel, G4VLEPTSModel, G4VMscModel, and G4XrayRayleighModel.

Collaboration diagram for G4VEmModel:

Public Member Functions

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

Protected Member Functions

G4ParticleChangeForLoss * GetParticleChangeForLoss ()
 
G4ParticleChangeForGamma * GetParticleChangeForGamma ()
 
virtual G4double MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kineticEnergy)
 
const G4MaterialCutsCoupleCurrentCouple () const
 
void SetCurrentElement (const G4Element *)
 

Protected Attributes

G4ElementDatafElementData
 
G4VParticleChange * pParticleChange
 
G4PhysicsTablexSectionTable
 
const std::vector< G4double > * theDensityFactor
 
const std::vector< G4int > * theDensityIdx
 
size_t idxTable
 

Static Protected Attributes

static const G4double inveplus = 1.0/CLHEP::eplus
 

Private Member Functions

G4VEmModeloperator= (const G4VEmModel &right)
 
 G4VEmModel (const G4VEmModel &)
 

Private Attributes

G4VEmFluctuationModelflucModel
 
G4VEmAngularDistributionanglModel
 
const G4String name
 
G4double lowLimit
 
G4double highLimit
 
G4double eMinActive
 
G4double eMaxActive
 
G4double polarAngleLimit
 
G4double secondaryThreshold
 
G4bool theLPMflag
 
G4bool flagDeexcitation
 
G4bool flagForceBuildTable
 
G4bool isMaster
 
G4bool localTable
 
G4bool localElmSelectors
 
G4bool useAngularGenerator
 
G4bool isLocked
 
G4int nSelectors
 
std::vector< G4EmElementSelector * > * elmSelectors
 
G4LossTableManagerfManager
 
const G4MaterialCutsCouplefCurrentCouple
 
const G4ElementfCurrentElement
 
const G4IsotopefCurrentIsotope
 
G4int nsec
 
std::vector< G4doublexsec
 

Detailed Description

Definition at line 104 of file G4VEmModel.hh.

Constructor & Destructor Documentation

◆ G4VEmModel() [1/2]

G4VEmModel::G4VEmModel ( const G4String nam)

Definition at line 69 of file G4VEmModel.cc.

69  :
70  flucModel(nullptr),anglModel(nullptr), name(nam), lowLimit(0.1*CLHEP::keV),
73  theLPMflag(false),flagDeexcitation(false),flagForceBuildTable(false),
74  isMaster(true),fElementData(nullptr),pParticleChange(nullptr),xSectionTable(nullptr),
75  theDensityFactor(nullptr),theDensityIdx(nullptr),fCurrentCouple(nullptr),
76  fCurrentElement(nullptr),fCurrentIsotope(nullptr),nsec(5)
77 {
78  xsec.resize(nsec);
79  nSelectors = 0;
80  elmSelectors = nullptr;
81  localElmSelectors = true;
82  localTable = true;
83  useAngularGenerator = false;
84  isLocked = false;
85  idxTable = 0;
86 
88  fManager->Register(this);
89 }
G4bool localTable
Definition: G4VEmModel.hh:412
size_t idxTable
Definition: G4VEmModel.hh:426
G4bool flagForceBuildTable
Definition: G4VEmModel.hh:409
static G4LossTableManager * Instance()
G4bool flagDeexcitation
Definition: G4VEmModel.hh:408
G4bool theLPMflag
Definition: G4VEmModel.hh:407
G4bool isMaster
Definition: G4VEmModel.hh:410
const G4String name
Definition: G4VEmModel.hh:397
const G4Element * fCurrentElement
Definition: G4VEmModel.hh:435
G4VEmFluctuationModel * flucModel
Definition: G4VEmModel.hh:395
G4int nSelectors
Definition: G4VEmModel.hh:416
G4bool isLocked
Definition: G4VEmModel.hh:415
G4LossTableManager * fManager
Definition: G4VEmModel.hh:433
const std::vector< G4int > * theDensityIdx
Definition: G4VEmModel.hh:425
G4double eMaxActive
Definition: G4VEmModel.hh:404
G4double secondaryThreshold
Definition: G4VEmModel.hh:406
const G4MaterialCutsCouple * fCurrentCouple
Definition: G4VEmModel.hh:434
G4double polarAngleLimit
Definition: G4VEmModel.hh:405
std::vector< G4EmElementSelector * > * elmSelectors
Definition: G4VEmModel.hh:417
static const double pi
Definition: SystemOfUnits.h:53
G4VEmAngularDistribution * anglModel
Definition: G4VEmModel.hh:396
void Register(G4VEnergyLossProcess *p)
G4int nsec
Definition: G4VEmModel.hh:438
const std::vector< G4double > * theDensityFactor
Definition: G4VEmModel.hh:424
G4double eMinActive
Definition: G4VEmModel.hh:403
std::vector< G4double > xsec
Definition: G4VEmModel.hh:439
G4bool localElmSelectors
Definition: G4VEmModel.hh:413
G4double highLimit
Definition: G4VEmModel.hh:402
G4bool useAngularGenerator
Definition: G4VEmModel.hh:414
static const double TeV
G4double lowLimit
Definition: G4VEmModel.hh:401
G4VParticleChange * pParticleChange
Definition: G4VEmModel.hh:422
G4PhysicsTable * xSectionTable
Definition: G4VEmModel.hh:423
const G4Isotope * fCurrentIsotope
Definition: G4VEmModel.hh:436
#define DBL_MAX
Definition: templates.hh:83
static const double keV
G4ElementData * fElementData
Definition: G4VEmModel.hh:421
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ~G4VEmModel()

G4VEmModel::~G4VEmModel ( )
virtual

Definition at line 93 of file G4VEmModel.cc.

94 {
95  if(localElmSelectors) {
96  if(nSelectors > 0) {
97  for(G4int i=0; i<nSelectors; ++i) {
98  delete (*elmSelectors)[i];
99  }
100  }
101  delete elmSelectors;
102  }
103  delete anglModel;
104 
105  if(localTable && xSectionTable) {
107  delete xSectionTable;
108  xSectionTable = nullptr;
109  }
110  if(isMaster && fElementData) {
111  delete fElementData;
112  fElementData = nullptr;
113  }
114 
115  fManager->DeRegister(this);
116 }
G4bool localTable
Definition: G4VEmModel.hh:412
void DeRegister(G4VEnergyLossProcess *p)
G4bool isMaster
Definition: G4VEmModel.hh:410
int G4int
Definition: G4Types.hh:78
G4int nSelectors
Definition: G4VEmModel.hh:416
G4LossTableManager * fManager
Definition: G4VEmModel.hh:433
std::vector< G4EmElementSelector * > * elmSelectors
Definition: G4VEmModel.hh:417
G4VEmAngularDistribution * anglModel
Definition: G4VEmModel.hh:396
G4bool localElmSelectors
Definition: G4VEmModel.hh:413
G4PhysicsTable * xSectionTable
Definition: G4VEmModel.hh:423
G4ElementData * fElementData
Definition: G4VEmModel.hh:421
void clearAndDestroy()
Here is the call graph for this function:

◆ G4VEmModel() [2/2]

G4VEmModel::G4VEmModel ( const G4VEmModel )
private

Member Function Documentation

◆ ChargeSquareRatio()

G4double G4VEmModel::ChargeSquareRatio ( const G4Track &  track)
virtual

Reimplemented in G4BraggIonGasModel, and G4BetheBlochIonGasModel.

Definition at line 337 of file G4VEmModel.cc.

338 {
339  return GetChargeSquareRatio(track.GetParticleDefinition(),
340  track.GetMaterial(), track.GetKineticEnergy());
341 }
virtual G4double GetChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:345
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ComputeCrossSectionPerAtom() [1/2]

G4double G4VEmModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition ,
G4double  kinEnergy,
G4double  Z,
G4double  A = 0.,
G4double  cutEnergy = 0.0,
G4double  maxEnergy = DBL_MAX 
)
virtual

Reimplemented in G4IonParametrisedLossModel, G4LowEPPolarizedComptonModel, G4LowEPComptonModel, G4BetheBlochModel, G4eCoulombScatteringModel, G4MuBetheBlochModel, G4MuBremsstrahlungModel, G4BraggIonModel, G4BraggModel, G4ICRU73QOModel, G4MuPairProductionModel, G4UrbanMscModel, G4PenelopeComptonModel, G4eBremParametrizedModel, G4PenelopeIonisationModel, G4eBremsstrahlungRelModel, G4WentzelVIModel, G4WentzelVIRelModel, G4eeToHadronsMultiModel, G4mplIonisationWithDeltaModel, G4hCoulombScatteringModel, G4IonCoulombScatteringModel, G4MollerBhabhaModel, G4EmMultiModel, G4eSingleCoulombScatteringModel, G4eeToTwoGammaModel, G4eeToHadronsModel, G4PairProductionRelModel, G4BetheHeitlerModel, G4PenelopePhotoElectricModel, G4KleinNishinaCompton, G4KleinNishinaModel, G4LivermoreIonisationModel, G4PenelopeRayleighModel, G4PEEffectFluoModel, G4LivermorePolarizedComptonModel, G4PenelopeAnnihilationModel, G4PenelopeGammaConversionModel, G4PolarizedComptonModel, G4LivermorePolarizedPhotoElectricModel, G4LivermorePolarizedRayleighModel, G4LivermoreComptonModel, G4LivermoreGammaConversionModel, G4LivermoreComptonModifiedModel, G4LivermorePolarizedGammaConversionModel, G4LivermoreNuclearGammaConversionModel, G4BoldyshevTripletModel, G4LivermoreGammaConversionModelRC, G4LivermorePolarizedPhotoElectricGDModel, G4LivermoreRayleighModel, G4XrayRayleighModel, and G4DiscreteScatteringModel.

Definition at line 313 of file G4VEmModel.cc.

316 {
317  return 0.0;
318 }
Here is the caller graph for this function:

◆ ComputeCrossSectionPerAtom() [2/2]

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

Definition at line 530 of file G4VEmModel.hh.

535 {
536  SetCurrentElement(elm);
537  return ComputeCrossSectionPerAtom(part,kinEnergy,elm->GetZ(),elm->GetN(),
538  cutEnergy,maxEnergy);
539 }
G4double GetN() const
Definition: G4Element.hh:134
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:313
void SetCurrentElement(const G4Element *)
Definition: G4VEmModel.hh:459
G4double GetZ() const
Definition: G4Element.hh:131
Here is the call graph for this function:

◆ ComputeCrossSectionPerShell()

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

Definition at line 323 of file G4VEmModel.cc.

326 {
327  return 0.0;
328 }
Here is the caller graph for this function:

◆ ComputeDEDX()

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

Definition at line 490 of file G4VEmModel.hh.

494 {
495  SetCurrentCouple(couple);
496  return ComputeDEDXPerVolume(couple->GetMaterial(),part,kinEnergy,cutEnergy);
497 }
const G4Material * GetMaterial() const
TString part[npart]
virtual G4double ComputeDEDXPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
Definition: G4VEmModel.cc:249
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:445
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ComputeDEDXPerVolume()

G4double G4VEmModel::ComputeDEDXPerVolume ( const G4Material ,
const G4ParticleDefinition ,
G4double  kineticEnergy,
G4double  cutEnergy = DBL_MAX 
)
virtual

◆ ComputeMeanFreePath()

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

Definition at line 515 of file G4VEmModel.hh.

520 {
521  G4double mfp = DBL_MAX;
522  G4double cross = CrossSectionPerVolume(material,part,ekin,emin,emax);
523  if (cross > DBL_MIN) { mfp = 1./cross; }
524  return mfp;
525 }
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:258
static const G4double emax
#define DBL_MIN
Definition: templates.hh:75
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83
Here is the call graph for this function:
Here is the caller graph for this function:

◆ CorrectionsAlongStep()

void G4VEmModel::CorrectionsAlongStep ( const G4MaterialCutsCouple ,
const G4DynamicParticle ,
G4double eloss,
G4double niel,
G4double  length 
)
virtual

Reimplemented in G4IonParametrisedLossModel, G4BraggIonModel, G4BetheBlochModel, and G4ICRU73QOModel.

Definition at line 362 of file G4VEmModel.cc.

365 {}
Here is the caller graph for this function:

◆ CrossSection()

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

Definition at line 501 of file G4VEmModel.hh.

506 {
507  SetCurrentCouple(couple);
508  return CrossSectionPerVolume(couple->GetMaterial(),part,kinEnergy,
509  cutEnergy,maxEnergy);
510 }
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:258
const G4Material * GetMaterial() const
TString part[npart]
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:445
Here is the call graph for this function:
Here is the caller graph for this function:

◆ CrossSectionPerVolume()

G4double G4VEmModel::CrossSectionPerVolume ( const G4Material material,
const G4ParticleDefinition p,
G4double  kineticEnergy,
G4double  cutEnergy = 0.0,
G4double  maxEnergy = DBL_MAX 
)
virtual

Reimplemented in G4IonParametrisedLossModel, G4BetheBlochModel, G4MuBetheBlochModel, G4BraggIonModel, G4BraggModel, G4ICRU73QOModel, G4PenelopeIonisationModel, G4MollerBhabhaModel, G4PAIModel, G4eeToTwoGammaModel, G4PAIPhotModel, G4PenelopeComptonModel, G4PEEffectFluoModel, G4MicroElecInelasticModel, G4MuElecInelasticModel, G4DNAOneStepThermalizationModel, G4eeToHadronsMultiModel, G4eeToHadronsModel, G4DNATransformElectronModel, G4BetheBlochNoDeltaModel, G4BraggNoDeltaModel, G4ICRU73NoDeltaModel, G4DNAEmfietzoglouIonisationModel, G4VDNAPTBModel, G4DNAEmfietzoglouExcitationModel, G4MicroElecElasticModel, G4MuElecElasticModel, G4DNAIonElasticModel, G4LivermorePolarizedPhotoElectricModel, G4DNABornIonisationModel1, G4DNABornIonisationModel2, G4DNARuddIonisationExtendedModel, G4DNARuddIonisationModel, G4DNASancheExcitationModel, G4DNABornExcitationModel2, G4DNAMeltonAttachmentModel, G4DNAMillerGreenExcitationModel, MyKleinNishinaCompton, G4DNABornExcitationModel1, G4DNAChampionElasticModel, G4DNADingfelderChargeIncreaseModel, G4DNADingfelderChargeDecreaseModel, G4DNAScreenedRutherfordElasticModel, G4LEPTSAttachmentModel, G4LEPTSExcitationModel, G4DNAUeharaScreenedRutherfordElasticModel, G4LEPTSDissociationModel, G4LEPTSElasticModel, G4LEPTSIonisationModel, G4LEPTSPositroniumModel, G4LEPTSRotExcitationModel, and G4LEPTSVibExcitationModel.

Definition at line 258 of file G4VEmModel.cc.

263 {
264  SetupForMaterial(p, material, ekin);
265  G4double cross = 0.0;
266  const G4ElementVector* theElementVector = material->GetElementVector();
267  const G4double* theAtomNumDensityVector =
268  material->GetVecNbOfAtomsPerVolume();
269  G4int nelm = material->GetNumberOfElements();
270  if(nelm > nsec) {
271  xsec.resize(nelm);
272  nsec = nelm;
273  }
274  for (G4int i=0; i<nelm; ++i) {
275  cross += theAtomNumDensityVector[i]*
276  ComputeCrossSectionPerAtom(p,(*theElementVector)[i],ekin,emin,emax);
277  xsec[i] = cross;
278  }
279  return cross;
280 }
std::vector< G4Element * > G4ElementVector
virtual void SetupForMaterial(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:403
int G4int
Definition: G4Types.hh:78
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:206
G4int nsec
Definition: G4VEmModel.hh:438
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:313
static const G4double emax
std::vector< G4double > xsec
Definition: G4VEmModel.hh:439
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
double G4double
Definition: G4Types.hh:76
Here is the call graph for this function:
Here is the caller graph for this function:

◆ CurrentCouple()

const G4MaterialCutsCouple * G4VEmModel::CurrentCouple ( ) const
inlineprotected

Definition at line 452 of file G4VEmModel.hh.

453 {
454  return fCurrentCouple;
455 }
const G4MaterialCutsCouple * fCurrentCouple
Definition: G4VEmModel.hh:434
Here is the caller graph for this function:

◆ DeexcitationFlag()

G4bool G4VEmModel::DeexcitationFlag ( ) const
inline

Definition at line 683 of file G4VEmModel.hh.

684 {
685  return flagDeexcitation;
686 }
G4bool flagDeexcitation
Definition: G4VEmModel.hh:408
Here is the caller graph for this function:

◆ DefineForRegion()

void G4VEmModel::DefineForRegion ( const G4Region )
virtual

Reimplemented in G4PAIModel, and G4PAIPhotModel.

Definition at line 332 of file G4VEmModel.cc.

333 {}
Here is the caller graph for this function:

◆ ForceBuildTableFlag()

G4bool G4VEmModel::ForceBuildTableFlag ( ) const
inline

Definition at line 690 of file G4VEmModel.hh.

691 {
692  return flagForceBuildTable;
693 }
G4bool flagForceBuildTable
Definition: G4VEmModel.hh:409
Here is the caller graph for this function:

◆ GetAngularDistribution()

G4VEmAngularDistribution * G4VEmModel::GetAngularDistribution ( )
inline

Definition at line 617 of file G4VEmModel.hh.

618 {
619  return anglModel;
620 }
G4VEmAngularDistribution * anglModel
Definition: G4VEmModel.hh:396
Here is the caller graph for this function:

◆ GetChargeSquareRatio()

G4double G4VEmModel::GetChargeSquareRatio ( const G4ParticleDefinition p,
const G4Material ,
G4double  kineticEnergy 
)
virtual

Reimplemented in G4IonParametrisedLossModel, G4BraggIonModel, G4BraggModel, and G4BetheBlochModel.

Definition at line 345 of file G4VEmModel.cc.

347 {
348  G4double q = p->GetPDGCharge()*inveplus;
349  return q*q;
350 }
static const G4double inveplus
Definition: G4VEmModel.hh:427
double G4double
Definition: G4Types.hh:76
G4double GetPDGCharge() const
Here is the call graph for this function:
Here is the caller graph for this function:

◆ GetCrossSectionTable()

G4PhysicsTable * G4VEmModel::GetCrossSectionTable ( )
inline

Definition at line 826 of file G4VEmModel.hh.

827 {
828  return xSectionTable;
829 }
G4PhysicsTable * xSectionTable
Definition: G4VEmModel.hh:423
Here is the caller graph for this function:

◆ GetCurrentElement()

const G4Element * G4VEmModel::GetCurrentElement ( ) const
inline

Definition at line 467 of file G4VEmModel.hh.

468 {
469  return fCurrentElement;
470 }
const G4Element * fCurrentElement
Definition: G4VEmModel.hh:435
Here is the caller graph for this function:

◆ GetCurrentIsotope()

const G4Isotope * G4VEmModel::GetCurrentIsotope ( ) const
inline

Definition at line 474 of file G4VEmModel.hh.

475 {
476  return fCurrentIsotope;
477 }
const G4Isotope * fCurrentIsotope
Definition: G4VEmModel.hh:436
Here is the caller graph for this function:

◆ GetElementData()

G4ElementData * G4VEmModel::GetElementData ( )
inline

Definition at line 819 of file G4VEmModel.hh.

820 {
821  return fElementData;
822 }
G4ElementData * fElementData
Definition: G4VEmModel.hh:421
Here is the caller graph for this function:

◆ GetElementSelectors()

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

Definition at line 802 of file G4VEmModel.hh.

803 {
804  return elmSelectors;
805 }
std::vector< G4EmElementSelector * > * elmSelectors
Definition: G4VEmModel.hh:417
Here is the caller graph for this function:

◆ GetModelOfFluctuations()

G4VEmFluctuationModel * G4VEmModel::GetModelOfFluctuations ( )
inline

Definition at line 610 of file G4VEmModel.hh.

611 {
612  return flucModel;
613 }
G4VEmFluctuationModel * flucModel
Definition: G4VEmModel.hh:395
Here is the caller graph for this function:

◆ GetName()

const G4String & G4VEmModel::GetName ( void  ) const
inline

Definition at line 795 of file G4VEmModel.hh.

796 {
797  return name;
798 }
const G4String name
Definition: G4VEmModel.hh:397
Here is the caller graph for this function:

◆ GetPartialCrossSection()

virtual G4double G4VEmModel::GetPartialCrossSection ( const G4Material ,
G4int  ,
const G4ParticleDefinition ,
G4double   
)
inlinevirtual

◆ GetParticleChangeForGamma()

G4ParticleChangeForGamma * G4VEmModel::GetParticleChangeForGamma ( )
protected

Definition at line 134 of file G4VEmModel.cc.

135 {
136  G4ParticleChangeForGamma* p = nullptr;
137  if (pParticleChange) {
138  p = static_cast<G4ParticleChangeForGamma*>(pParticleChange);
139  } else {
140  p = new G4ParticleChangeForGamma();
141  pParticleChange = p;
142  }
143  return p;
144 }
G4VParticleChange * pParticleChange
Definition: G4VEmModel.hh:422

◆ GetParticleChangeForLoss()

G4ParticleChangeForLoss * G4VEmModel::GetParticleChangeForLoss ( )
protected

Definition at line 120 of file G4VEmModel.cc.

121 {
122  G4ParticleChangeForLoss* p = nullptr;
123  if (pParticleChange) {
124  p = static_cast<G4ParticleChangeForLoss*>(pParticleChange);
125  } else {
126  p = new G4ParticleChangeForLoss();
127  pParticleChange = p;
128  }
129  return p;
130 }
G4VParticleChange * pParticleChange
Definition: G4VEmModel.hh:422
Here is the caller graph for this function:

◆ GetParticleCharge()

G4double G4VEmModel::GetParticleCharge ( const G4ParticleDefinition p,
const G4Material ,
G4double  kineticEnergy 
)
virtual

Reimplemented in G4IonParametrisedLossModel, G4BraggIonModel, G4BraggModel, G4BetheBlochModel, G4BraggIonGasModel, and G4BetheBlochIonGasModel.

Definition at line 354 of file G4VEmModel.cc.

356 {
357  return p->GetPDGCharge();
358 }
G4double GetPDGCharge() const
Here is the call graph for this function:
Here is the caller graph for this function:

◆ HighEnergyActivationLimit()

G4double G4VEmModel::HighEnergyActivationLimit ( ) const
inline

Definition at line 648 of file G4VEmModel.hh.

649 {
650  return eMaxActive;
651 }
G4double eMaxActive
Definition: G4VEmModel.hh:404
Here is the caller graph for this function:

◆ HighEnergyLimit()

G4double G4VEmModel::HighEnergyLimit ( ) const
inline

Definition at line 634 of file G4VEmModel.hh.

635 {
636  return highLimit;
637 }
G4double highLimit
Definition: G4VEmModel.hh:402

◆ Initialise()

virtual void G4VEmModel::Initialise ( const G4ParticleDefinition ,
const G4DataVector  
)
pure virtual

Implemented in G4GoudsmitSaundersonMscModel, G4IonParametrisedLossModel, G4LowEPPolarizedComptonModel, G4LowEPComptonModel, G4eCoulombScatteringModel, G4UrbanMscModel, G4MuBremsstrahlungModel, G4MuPairProductionModel, G4BetheBlochModel, G4BraggModel, G4ICRU73QOModel, G4MuBetheBlochModel, G4IonCoulombScatteringModel, G4BraggIonModel, G4PAIModel, G4WentzelVIRelModel, G4PenelopeIonisationModel, G4hCoulombScatteringModel, G4eSingleCoulombScatteringModel, G4WentzelVIModel, G4MicroElecInelasticModel, G4MuElecInelasticModel, G4DNAOneStepThermalizationModel, G4eeToHadronsMultiModel, G4PenelopeComptonModel, G4PolarizedAnnihilationModel, G4MollerBhabhaModel, G4PAIPhotModel, G4PolarizedPEEffectModel, G4eBremsstrahlungRelModel, G4LivermoreBremsstrahlungModel, G4eBremParametrizedModel, G4ICRU49NuclearStoppingModel, G4PairProductionRelModel, G4SeltzerBergerModel, G4BetheHeitlerModel, G4eeToTwoGammaModel, G4EmMultiModel, G4eeToHadronsModel, G4PenelopePhotoElectricModel, G4mplIonisationModel, G4mplIonisationWithDeltaModel, G4LivermoreIonisationModel, G4ePolarizedBremsstrahlungModel, G4PolarizedGammaConversionModel, G4KleinNishinaCompton, G4KleinNishinaModel, G4PEEffectFluoModel, G4DNATransformElectronModel, G4PenelopeGammaConversionModel, G4PenelopeRayleighModel, G4DNAEmfietzoglouIonisationModel, G4PenelopeAnnihilationModel, G4LivermorePolarizedComptonModel, G4DummyModel, G4DNAEmfietzoglouExcitationModel, G4MicroElecElasticModel, G4MuElecElasticModel, G4LivermorePolarizedRayleighModel, G4LivermorePolarizedPhotoElectricModel, G4DNABornIonisationModel1, G4DNABornIonisationModel2, G4DNAIonElasticModel, G4LivermoreComptonModifiedModel, G4LivermoreComptonModel, G4DNARuddIonisationExtendedModel, G4DNARuddIonisationModel, G4DNASancheExcitationModel, G4DNAMeltonAttachmentModel, G4DNAMillerGreenExcitationModel, G4DNABornExcitationModel2, G4LivermorePolarizedGammaConversionModel, G4LivermorePolarizedPhotoElectricGDModel, G4XrayRayleighModel, G4DNABornExcitationModel1, G4DNAChampionElasticModel, G4DNADingfelderChargeIncreaseModel, G4LivermoreGammaConversionModel, G4DNADingfelderChargeDecreaseModel, G4LivermoreGammaConversionModelRC, G4DNAScreenedRutherfordElasticModel, G4LivermoreNuclearGammaConversionModel, G4LivermoreRayleighModel, G4VDNAPTBModel, G4BoldyshevTripletModel, G4DNAUeharaScreenedRutherfordElasticModel, G4LEPTSAttachmentModel, G4LEPTSDissociationModel, G4LEPTSElasticModel, G4LEPTSExcitationModel, G4LEPTSIonisationModel, G4LEPTSPositroniumModel, G4LEPTSRotExcitationModel, G4LEPTSVibExcitationModel, and G4DiscreteScatteringModel.

Here is the caller graph for this function:

◆ InitialiseElementSelectors()

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

Definition at line 148 of file G4VEmModel.cc.

150 {
151  // using spline for element selectors should be investigated in details
152  // because small number of points may provide biased results
153  // large number of points requires significant increase of memory
154  //G4bool spline = fManager->SplineFlag();
155  G4bool spline = false;
156 
157  //G4cout << "IES: for " << GetName() << " Emin(MeV)= " << lowLimit/MeV
158  // << " Emax(MeV)= " << highLimit/MeV << G4endl;
159 
160  // two times less bins because probability functon is normalized
161  // so correspondingly is more smooth
162  if(highLimit <= lowLimit) { return; }
163 
164  G4ProductionCutsTable* theCoupleTable=
166  G4int numOfCouples = theCoupleTable->GetTableSize();
167 
168  // prepare vector
169  if(!elmSelectors) {
170  elmSelectors = new std::vector<G4EmElementSelector*>;
171  }
172  if(numOfCouples > nSelectors) {
173  for(G4int i=nSelectors; i<numOfCouples; ++i) {
174  elmSelectors->push_back(nullptr);
175  }
176  nSelectors = numOfCouples;
177  }
178 
179  // initialise vector
180  for(G4int i=0; i<numOfCouples; ++i) {
181 
182  // no need in element selectors for infionite cuts
183  if(cuts[i] == DBL_MAX) { continue; }
184 
185  fCurrentCouple = theCoupleTable->GetMaterialCutsCouple(i);
187 
188  // selector already exist check if should be deleted
189  G4bool create = true;
190  if((*elmSelectors)[i]) {
191  if(material == ((*elmSelectors)[i])->GetMaterial()) { create = false; }
192  else { delete (*elmSelectors)[i]; }
193  }
194  if(create) {
195  G4double emin = std::max(lowLimit,
196  MinPrimaryEnergy(material, part, cuts[i]));
197  G4double emax = std::max(highLimit, 10*emin);
199  *G4Log(emax/emin)/log106);
200  nbins = std::max(nbins, 3);
201 
202  (*elmSelectors)[i] = new G4EmElementSelector(this,material,nbins,
203  emin,emax,spline);
204  }
205  ((*elmSelectors)[i])->Initialise(part, cuts[i]);
206  /*
207  G4cout << "G4VEmModel::InitialiseElmSelectors i= " << i
208  << " idx= " << fCurrentCouple->GetIndex()
209  << " " << part->GetParticleName()
210  << " for " << GetName() << " cut= " << cuts[i]
211  << " " << (*elmSelectors)[i] << G4endl;
212  ((*elmSelectors)[i])->Dump(part);
213  */
214  }
215 }
const G4Material * GetMaterial() const
virtual G4double MinPrimaryEnergy(const G4Material *, const G4ParticleDefinition *, G4double cut=0.0)
Definition: G4VEmModel.cc:378
int G4int
Definition: G4Types.hh:78
G4int nSelectors
Definition: G4VEmModel.hh:416
G4LossTableManager * fManager
Definition: G4VEmModel.hh:433
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)=0
string material
Definition: eplot.py:19
const G4MaterialCutsCouple * fCurrentCouple
Definition: G4VEmModel.hh:434
std::vector< G4EmElementSelector * > * elmSelectors
Definition: G4VEmModel.hh:417
bool G4bool
Definition: G4Types.hh:79
G4int GetNumberOfBinsPerDecade() const
static const G4double emax
G4double G4Log(G4double x)
Definition: G4Log.hh:230
static G4ProductionCutsTable * GetProductionCutsTable()
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
G4double highLimit
Definition: G4VEmModel.hh:402
G4double lowLimit
Definition: G4VEmModel.hh:401
const G4double log106
Definition: G4VEmModel.cc:67
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83
Here is the call graph for this function:
Here is the caller graph for this function:

◆ InitialiseForElement()

◆ InitialiseForMaterial()

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

Definition at line 225 of file G4VEmModel.cc.

227 {
228  if(material) {
229  const G4ElementVector* theElementVector = material->GetElementVector();
230  G4int n = material->GetNumberOfElements();
231  for(G4int i=0; i<n; ++i) {
232  G4int Z = G4lrint(((*theElementVector)[i])->GetZ());
233  InitialiseForElement(part, Z);
234  }
235  } else {
236  //G4cout << "G4VEmModel::InitialiseForMaterial for " << GetName();
237  //if(part) { G4cout << " and " << part->GetParticleName(); }
238  //G4cout << " with no material" << G4endl;
239  }
240 }
virtual void InitialiseForElement(const G4ParticleDefinition *, G4int Z)
Definition: G4VEmModel.cc:244
std::vector< G4Element * > G4ElementVector
int G4int
Definition: G4Types.hh:78
Char_t n[5]
Float_t Z
int G4lrint(double ad)
Definition: templates.hh:163
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
Here is the call graph for this function:
Here is the caller graph for this function:

◆ InitialiseLocal()

◆ IsActive()

G4bool G4VEmModel::IsActive ( G4double  kinEnergy)
inline

Definition at line 753 of file G4VEmModel.hh.

754 {
755  return (kinEnergy >= eMinActive && kinEnergy <= eMaxActive);
756 }
G4double eMaxActive
Definition: G4VEmModel.hh:404
G4double eMinActive
Definition: G4VEmModel.hh:403
Here is the caller graph for this function:

◆ IsLocked()

G4bool G4VEmModel::IsLocked ( ) const
inline

Definition at line 833 of file G4VEmModel.hh.

834 {
835  return isLocked;
836 }
G4bool isLocked
Definition: G4VEmModel.hh:415
Here is the caller graph for this function:

◆ IsMaster()

G4bool G4VEmModel::IsMaster ( ) const
inline

Definition at line 718 of file G4VEmModel.hh.

719 {
720  return isMaster;
721 }
G4bool isMaster
Definition: G4VEmModel.hh:410

◆ LowEnergyActivationLimit()

G4double G4VEmModel::LowEnergyActivationLimit ( ) const
inline

Definition at line 655 of file G4VEmModel.hh.

656 {
657  return eMinActive;
658 }
G4double eMinActive
Definition: G4VEmModel.hh:403
Here is the caller graph for this function:

◆ LowEnergyLimit()

G4double G4VEmModel::LowEnergyLimit ( ) const
inline

Definition at line 641 of file G4VEmModel.hh.

642 {
643  return lowLimit;
644 }
G4double lowLimit
Definition: G4VEmModel.hh:401

◆ LPMFlag()

G4bool G4VEmModel::LPMFlag ( ) const
inline

Definition at line 676 of file G4VEmModel.hh.

677 {
678  return theLPMflag;
679 }
G4bool theLPMflag
Definition: G4VEmModel.hh:407
Here is the caller graph for this function:

◆ MaxSecondaryEnergy()

G4double G4VEmModel::MaxSecondaryEnergy ( const G4ParticleDefinition ,
G4double  kineticEnergy 
)
protectedvirtual

Reimplemented in G4IonParametrisedLossModel, G4BetheBlochModel, G4BraggIonModel, G4BraggModel, G4PAIModel, G4ICRU73QOModel, G4PAIPhotModel, G4MuBetheBlochModel, G4mplIonisationWithDeltaModel, and G4MollerBhabhaModel.

Definition at line 395 of file G4VEmModel.cc.

397 {
398  return kineticEnergy;
399 }
Here is the caller graph for this function:

◆ MaxSecondaryKinEnergy()

G4double G4VEmModel::MaxSecondaryKinEnergy ( const G4DynamicParticle dynParticle)
inline

Definition at line 482 of file G4VEmModel.hh.

483 {
484  return MaxSecondaryEnergy(dynPart->GetParticleDefinition(),
485  dynPart->GetKineticEnergy());
486 }
virtual G4double MaxSecondaryEnergy(const G4ParticleDefinition *, G4double kineticEnergy)
Definition: G4VEmModel.cc:395
Here is the call graph for this function:
Here is the caller graph for this function:

◆ MinEnergyCut()

G4double G4VEmModel::MinEnergyCut ( const G4ParticleDefinition ,
const G4MaterialCutsCouple  
)
virtual

Reimplemented in G4PenelopeIonisationModel, G4IonParametrisedLossModel, G4MuBremsstrahlungModel, G4BetheBlochModel, G4MuBetheBlochModel, G4BraggIonModel, and G4eBremParametrizedModel.

Definition at line 387 of file G4VEmModel.cc.

389 {
390  return 0.0;
391 }
Here is the caller graph for this function:

◆ MinPrimaryEnergy()

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

Reimplemented in G4MuBremsstrahlungModel, G4eCoulombScatteringModel, G4MuPairProductionModel, G4eBremsstrahlungRelModel, G4hCoulombScatteringModel, G4LivermoreGammaConversionModel, G4LivermorePolarizedGammaConversionModel, G4LivermoreNuclearGammaConversionModel, G4BoldyshevTripletModel, and G4LivermoreGammaConversionModelRC.

Definition at line 378 of file G4VEmModel.cc.

381 {
382  return 0.0;
383 }
Here is the caller graph for this function:

◆ ModelDescription()

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

Definition at line 432 of file G4VEmModel.cc.

433 {
434  outFile << "The description for this model has not been written yet.\n";
435 }
Here is the caller graph for this function:

◆ operator=()

G4VEmModel& G4VEmModel::operator= ( const G4VEmModel right)
private
Here is the caller graph for this function:

◆ PolarAngleLimit()

G4double G4VEmModel::PolarAngleLimit ( ) const
inline

Definition at line 662 of file G4VEmModel.hh.

663 {
664  return polarAngleLimit;
665 }
G4double polarAngleLimit
Definition: G4VEmModel.hh:405
Here is the caller graph for this function:

◆ SampleSecondaries()

virtual void G4VEmModel::SampleSecondaries ( std::vector< G4DynamicParticle *> *  ,
const G4MaterialCutsCouple ,
const G4DynamicParticle ,
G4double  tmin = 0.0,
G4double  tmax = DBL_MAX 
)
pure virtual

Implemented in G4IonParametrisedLossModel, G4BetheBlochModel, G4LowEPPolarizedComptonModel, G4MuBetheBlochModel, G4BraggIonModel, G4BraggModel, G4ICRU73QOModel, G4LowEPComptonModel, G4MuBremsstrahlungModel, G4eCoulombScatteringModel, G4MuPairProductionModel, G4MollerBhabhaModel, G4PenelopeIonisationModel, G4PAIModel, G4PenelopeComptonModel, G4eeToTwoGammaModel, G4eeToHadronsModel, G4eBremParametrizedModel, G4hCoulombScatteringModel, G4IonCoulombScatteringModel, G4PAIPhotModel, G4eeToHadronsMultiModel, G4mplIonisationWithDeltaModel, G4eBremsstrahlungRelModel, G4VMscModel, G4eSingleCoulombScatteringModel, G4EmMultiModel, G4PolarizedAnnihilationModel, G4PairProductionRelModel, G4BetheHeitlerModel, G4PEEffectFluoModel, G4MicroElecInelasticModel, G4MuElecInelasticModel, G4DNAOneStepThermalizationModel, G4PenelopePhotoElectricModel, G4KleinNishinaCompton, G4KleinNishinaModel, G4ICRU49NuclearStoppingModel, G4PolarizedComptonModel, G4LivermoreIonisationModel, G4PenelopeRayleighModel, G4LivermorePolarizedComptonModel, G4PenelopeGammaConversionModel, G4PenelopeAnnihilationModel, G4DNATransformElectronModel, G4mplIonisationModel, G4LivermoreBremsstrahlungModel, G4LivermorePolarizedPhotoElectricModel, G4LivermorePolarizedRayleighModel, G4PolarizedMollerBhabhaModel, G4PolarizedPEEffectModel, G4SeltzerBergerModel, G4DNAEmfietzoglouIonisationModel, G4VDNAPTBModel, G4LivermoreComptonModel, G4DNAEmfietzoglouExcitationModel, G4ePolarizedBremsstrahlungModel, G4PolarizedGammaConversionModel, G4LivermoreGammaConversionModel, G4MicroElecElasticModel, G4MuElecElasticModel, G4LivermorePolarizedGammaConversionModel, G4DNABornExcitationModel2, G4DNAIonElasticModel, G4DNAMillerGreenExcitationModel, G4LivermoreComptonModifiedModel, G4HeatedKleinNishinaCompton, G4DNABornIonisationModel1, G4DNABornIonisationModel2, G4LivermoreNuclearGammaConversionModel, G4DNABornExcitationModel1, G4BoldyshevTripletModel, G4LivermoreGammaConversionModelRC, G4DummyModel, G4DNARuddIonisationExtendedModel, G4DNARuddIonisationModel, G4DNASancheExcitationModel, G4LivermorePolarizedPhotoElectricGDModel, G4LivermoreRayleighModel, G4DNAMeltonAttachmentModel, G4XrayRayleighModel, MyKleinNishinaCompton, G4DNAChampionElasticModel, G4DNADingfelderChargeIncreaseModel, G4DNADingfelderChargeDecreaseModel, G4DNAScreenedRutherfordElasticModel, G4DNAUeharaScreenedRutherfordElasticModel, G4LEPTSAttachmentModel, G4LEPTSExcitationModel, G4LEPTSDissociationModel, G4LEPTSElasticModel, G4LEPTSIonisationModel, G4LEPTSPositroniumModel, G4LEPTSRotExcitationModel, G4LEPTSVibExcitationModel, and G4DiscreteScatteringModel.

Here is the caller graph for this function:

◆ SecondaryThreshold()

G4double G4VEmModel::SecondaryThreshold ( ) const
inline

Definition at line 669 of file G4VEmModel.hh.

670 {
671  return secondaryThreshold;
672 }
G4double secondaryThreshold
Definition: G4VEmModel.hh:406
Here is the caller graph for this function:

◆ SelectIsotopeNumber()

G4int G4VEmModel::SelectIsotopeNumber ( const G4Element elm)
inline

Definition at line 585 of file G4VEmModel.hh.

586 {
587  SetCurrentElement(elm);
588  G4int N = G4lrint(elm->GetN());
589  G4int ni = elm->GetNumberOfIsotopes();
590  fCurrentIsotope = 0;
591  if(ni > 0) {
592  G4int idx = 0;
593  if(ni > 1) {
596  for(; idx<ni; ++idx) {
597  x -= ab[idx];
598  if (x <= 0.0) { break; }
599  }
600  if(idx >= ni) { idx = ni - 1; }
601  }
602  fCurrentIsotope = elm->GetIsotope(idx);
603  N = fCurrentIsotope->GetN();
604  }
605  return N;
606 }
const G4Isotope * GetIsotope(G4int iso) const
Definition: G4Element.hh:169
int G4int
Definition: G4Types.hh:78
size_t GetNumberOfIsotopes() const
Definition: G4Element.hh:158
G4double * GetRelativeAbundanceVector() const
Definition: G4Element.hh:166
G4double GetN() const
Definition: G4Element.hh:134
#define G4UniformRand()
Definition: Randomize.hh:97
static const G4double ab
int G4lrint(double ad)
Definition: templates.hh:163
G4int GetN() const
Definition: G4Isotope.hh:94
**D E S C R I P T I O N
double G4double
Definition: G4Types.hh:76
const G4Isotope * fCurrentIsotope
Definition: G4VEmModel.hh:436
void SetCurrentElement(const G4Element *)
Definition: G4VEmModel.hh:459
Here is the call graph for this function:
Here is the caller graph for this function:

◆ SelectRandomAtom() [1/2]

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

Definition at line 544 of file G4VEmModel.hh.

549 {
550  fCurrentCouple = couple;
551  if(nSelectors > 0) {
552  fCurrentElement =
553  ((*elmSelectors)[couple->GetIndex()])->SelectRandomAtom(kinEnergy);
554  } else {
555  fCurrentElement = SelectRandomAtom(couple->GetMaterial(),part,kinEnergy,
556  cutEnergy,maxEnergy);
557  }
558  fCurrentIsotope = 0;
559  return fCurrentElement;
560 }
const G4Material * GetMaterial() const
const G4Element * fCurrentElement
Definition: G4VEmModel.hh:435
G4int nSelectors
Definition: G4VEmModel.hh:416
const G4MaterialCutsCouple * fCurrentCouple
Definition: G4VEmModel.hh:434
TString part[npart]
const G4Isotope * fCurrentIsotope
Definition: G4VEmModel.hh:436
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:544
Here is the call graph for this function:
Here is the caller graph for this function:

◆ SelectRandomAtom() [2/2]

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

Definition at line 289 of file G4VEmModel.cc.

294 {
295  const G4ElementVector* theElementVector = material->GetElementVector();
296  G4int n = material->GetNumberOfElements() - 1;
297  fCurrentElement = (*theElementVector)[n];
298  if (n > 0) {
300  G4VEmModel::CrossSectionPerVolume(material,pd,kinEnergy,tcut,tmax);
301  for(G4int i=0; i<n; ++i) {
302  if (x <= xsec[i]) {
303  fCurrentElement = (*theElementVector)[i];
304  break;
305  }
306  }
307  }
308  return fCurrentElement;
309 }
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:258
std::vector< G4Element * > G4ElementVector
int G4int
Definition: G4Types.hh:78
const G4Element * fCurrentElement
Definition: G4VEmModel.hh:435
Char_t n[5]
#define G4UniformRand()
Definition: Randomize.hh:97
std::vector< G4double > xsec
Definition: G4VEmModel.hh:439
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
double G4double
Definition: G4Types.hh:76
Here is the call graph for this function:

◆ SelectRandomAtomNumber()

G4int G4VEmModel::SelectRandomAtomNumber ( const G4Material mat)
inline

Definition at line 564 of file G4VEmModel.hh.

565 {
566  // this algorith assumes that cross section is proportional to
567  // number electrons multiplied by number of atoms
568  size_t nn = mat->GetNumberOfElements();
569  const G4ElementVector* elmv = mat->GetElementVector();
570  G4int Z = G4lrint((*elmv)[0]->GetZ());
571  if(1 < nn) {
572  const G4double* at = mat->GetVecNbOfAtomsPerVolume();
574  for( size_t i=0; i<nn; ++i) {
575  Z = G4lrint((*elmv)[0]->GetZ());
576  tot -= Z*at[i];
577  if(tot <= 0.0) { break; }
578  }
579  }
580  return Z;
581 }
std::vector< G4Element * > G4ElementVector
G4double GetTotNbOfAtomsPerVolume() const
Definition: G4Material.hh:209
int G4int
Definition: G4Types.hh:78
#define G4UniformRand()
Definition: Randomize.hh:97
Float_t Z
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:206
int G4lrint(double ad)
Definition: templates.hh:163
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
double G4double
Definition: G4Types.hh:76
Here is the call graph for this function:
Here is the caller graph for this function:

◆ SetActivationHighEnergyLimit()

void G4VEmModel::SetActivationHighEnergyLimit ( G4double  val)
inline

Definition at line 739 of file G4VEmModel.hh.

740 {
741  eMaxActive = val;
742 }
G4double eMaxActive
Definition: G4VEmModel.hh:404
Here is the caller graph for this function:

◆ SetActivationLowEnergyLimit()

void G4VEmModel::SetActivationLowEnergyLimit ( G4double  val)
inline

Definition at line 746 of file G4VEmModel.hh.

747 {
748  eMinActive = val;
749 }
G4double eMinActive
Definition: G4VEmModel.hh:403
Here is the caller graph for this function:

◆ SetAngularDistribution()

void G4VEmModel::SetAngularDistribution ( G4VEmAngularDistribution p)
inline

Definition at line 624 of file G4VEmModel.hh.

625 {
626  if(p != anglModel) {
627  delete anglModel;
628  anglModel = p;
629  }
630 }
G4VEmAngularDistribution * anglModel
Definition: G4VEmModel.hh:396
Here is the caller graph for this function:

◆ SetAngularGeneratorFlag()

void G4VEmModel::SetAngularGeneratorFlag ( G4bool  val)
inline

Definition at line 704 of file G4VEmModel.hh.

705 {
706  useAngularGenerator = val;
707 }
G4bool useAngularGenerator
Definition: G4VEmModel.hh:414
Here is the caller graph for this function:

◆ SetCrossSectionTable()

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

Definition at line 418 of file G4VEmModel.cc.

419 {
420  if(p != xSectionTable) {
421  if(xSectionTable && localTable) {
423  delete xSectionTable;
424  }
425  xSectionTable = p;
426  }
427  localTable = isLocal;
428 }
G4bool localTable
Definition: G4VEmModel.hh:412
G4PhysicsTable * xSectionTable
Definition: G4VEmModel.hh:423
void clearAndDestroy()
Here is the call graph for this function:
Here is the caller graph for this function:

◆ SetCurrentCouple()

void G4VEmModel::SetCurrentCouple ( const G4MaterialCutsCouple p)
inline

Definition at line 445 of file G4VEmModel.hh.

446 {
447  fCurrentCouple = p;
448 }
const G4MaterialCutsCouple * fCurrentCouple
Definition: G4VEmModel.hh:434
Here is the caller graph for this function:

◆ SetCurrentElement()

void G4VEmModel::SetCurrentElement ( const G4Element elm)
inlineprotected

Definition at line 459 of file G4VEmModel.hh.

460 {
461  fCurrentElement = elm;
462  fCurrentIsotope = 0;
463 }
const G4Element * fCurrentElement
Definition: G4VEmModel.hh:435
const G4Isotope * fCurrentIsotope
Definition: G4VEmModel.hh:436
Here is the caller graph for this function:

◆ SetDeexcitationFlag()

void G4VEmModel::SetDeexcitationFlag ( G4bool  val)
inline

Definition at line 781 of file G4VEmModel.hh.

782 {
783  flagDeexcitation = val;
784 }
G4bool flagDeexcitation
Definition: G4VEmModel.hh:408
Here is the caller graph for this function:

◆ SetElementSelectors()

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

Definition at line 810 of file G4VEmModel.hh.

811 {
812  elmSelectors = p;
813  if(elmSelectors) { nSelectors = elmSelectors->size(); }
814  localElmSelectors = false;
815 }
G4int nSelectors
Definition: G4VEmModel.hh:416
std::vector< G4EmElementSelector * > * elmSelectors
Definition: G4VEmModel.hh:417
G4bool localElmSelectors
Definition: G4VEmModel.hh:413
Here is the caller graph for this function:

◆ SetForceBuildTable()

void G4VEmModel::SetForceBuildTable ( G4bool  val)
inline

Definition at line 788 of file G4VEmModel.hh.

789 {
790  flagForceBuildTable = val;
791 }
G4bool flagForceBuildTable
Definition: G4VEmModel.hh:409
Here is the caller graph for this function:

◆ SetHighEnergyLimit()

void G4VEmModel::SetHighEnergyLimit ( G4double  val)
inline

Definition at line 725 of file G4VEmModel.hh.

726 {
727  highLimit = val;
728 }
G4double highLimit
Definition: G4VEmModel.hh:402

◆ SetLocked()

void G4VEmModel::SetLocked ( G4bool  val)
inline

Definition at line 840 of file G4VEmModel.hh.

841 {
842  isLocked = val;
843 }
G4bool isLocked
Definition: G4VEmModel.hh:415
Here is the caller graph for this function:

◆ SetLowEnergyLimit()

void G4VEmModel::SetLowEnergyLimit ( G4double  val)
inline

Definition at line 732 of file G4VEmModel.hh.

733 {
734  lowLimit = val;
735 }
G4double lowLimit
Definition: G4VEmModel.hh:401

◆ SetLPMFlag()

void G4VEmModel::SetLPMFlag ( G4bool  val)
inline

Definition at line 774 of file G4VEmModel.hh.

775 {
776  theLPMflag = val;
777 }
G4bool theLPMflag
Definition: G4VEmModel.hh:407
Here is the caller graph for this function:

◆ SetMasterThread()

void G4VEmModel::SetMasterThread ( G4bool  val)
inline

Definition at line 711 of file G4VEmModel.hh.

712 {
713  isMaster = val;
714 }
G4bool isMaster
Definition: G4VEmModel.hh:410
Here is the caller graph for this function:

◆ SetParticleChange()

void G4VEmModel::SetParticleChange ( G4VParticleChange *  p,
G4VEmFluctuationModel f = 0 
)

Definition at line 410 of file G4VEmModel.cc.

411 {
412  if(p && pParticleChange != p) { pParticleChange = p; }
413  flucModel = f;
414 }
G4VEmFluctuationModel * flucModel
Definition: G4VEmModel.hh:395
G4VParticleChange * pParticleChange
Definition: G4VEmModel.hh:422
Here is the caller graph for this function:

◆ SetPolarAngleLimit()

void G4VEmModel::SetPolarAngleLimit ( G4double  val)
inline

Definition at line 760 of file G4VEmModel.hh.

761 {
762  if(!isLocked) { polarAngleLimit = val; }
763 }
G4bool isLocked
Definition: G4VEmModel.hh:415
G4double polarAngleLimit
Definition: G4VEmModel.hh:405
Here is the caller graph for this function:

◆ SetSecondaryThreshold()

void G4VEmModel::SetSecondaryThreshold ( G4double  val)
inline

Definition at line 767 of file G4VEmModel.hh.

768 {
769  secondaryThreshold = val;
770 }
G4double secondaryThreshold
Definition: G4VEmModel.hh:406
Here is the caller graph for this function:

◆ SetupForMaterial()

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

Reimplemented in G4eBremParametrizedModel, G4eBremsstrahlungRelModel, and G4PairProductionRelModel.

Definition at line 403 of file G4VEmModel.cc.

405 {}
Here is the caller graph for this function:

◆ StartTracking()

void G4VEmModel::StartTracking ( G4Track *  )
virtual

Reimplemented in G4GoudsmitSaundersonMscModel, G4UrbanMscModel, G4WentzelVIModel, and G4WentzelVIRelModel.

Definition at line 284 of file G4VEmModel.cc.

285 {}
Here is the caller graph for this function:

◆ UseAngularGeneratorFlag()

G4bool G4VEmModel::UseAngularGeneratorFlag ( ) const
inline

Definition at line 697 of file G4VEmModel.hh.

698 {
699  return useAngularGenerator;
700 }
G4bool useAngularGenerator
Definition: G4VEmModel.hh:414
Here is the caller graph for this function:

◆ Value()

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

Definition at line 369 of file G4VEmModel.cc.

371 {
372  SetCurrentCouple(couple);
373  return e*e*CrossSectionPerVolume(couple->GetMaterial(),p,e,0.0,DBL_MAX);
374 }
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:258
const G4Material * GetMaterial() const
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:445
#define DBL_MAX
Definition: templates.hh:83
Here is the call graph for this function:
Here is the caller graph for this function:

Member Data Documentation

◆ anglModel

G4VEmAngularDistribution* G4VEmModel::anglModel
private

Definition at line 396 of file G4VEmModel.hh.

◆ elmSelectors

std::vector<G4EmElementSelector*>* G4VEmModel::elmSelectors
private

Definition at line 417 of file G4VEmModel.hh.

◆ eMaxActive

G4double G4VEmModel::eMaxActive
private

Definition at line 404 of file G4VEmModel.hh.

◆ eMinActive

G4double G4VEmModel::eMinActive
private

Definition at line 403 of file G4VEmModel.hh.

◆ fCurrentCouple

const G4MaterialCutsCouple* G4VEmModel::fCurrentCouple
private

Definition at line 434 of file G4VEmModel.hh.

◆ fCurrentElement

const G4Element* G4VEmModel::fCurrentElement
private

Definition at line 435 of file G4VEmModel.hh.

◆ fCurrentIsotope

const G4Isotope* G4VEmModel::fCurrentIsotope
private

Definition at line 436 of file G4VEmModel.hh.

◆ fElementData

G4ElementData* G4VEmModel::fElementData
protected

Definition at line 421 of file G4VEmModel.hh.

◆ flagDeexcitation

G4bool G4VEmModel::flagDeexcitation
private

Definition at line 408 of file G4VEmModel.hh.

◆ flagForceBuildTable

G4bool G4VEmModel::flagForceBuildTable
private

Definition at line 409 of file G4VEmModel.hh.

◆ flucModel

G4VEmFluctuationModel* G4VEmModel::flucModel
private

Definition at line 395 of file G4VEmModel.hh.

◆ fManager

G4LossTableManager* G4VEmModel::fManager
private

Definition at line 433 of file G4VEmModel.hh.

◆ highLimit

G4double G4VEmModel::highLimit
private

Definition at line 402 of file G4VEmModel.hh.

◆ idxTable

size_t G4VEmModel::idxTable
protected

Definition at line 426 of file G4VEmModel.hh.

◆ inveplus

const G4double G4VEmModel::inveplus = 1.0/CLHEP::eplus
staticprotected

Definition at line 427 of file G4VEmModel.hh.

◆ isLocked

G4bool G4VEmModel::isLocked
private

Definition at line 415 of file G4VEmModel.hh.

◆ isMaster

G4bool G4VEmModel::isMaster
private

Definition at line 410 of file G4VEmModel.hh.

◆ localElmSelectors

G4bool G4VEmModel::localElmSelectors
private

Definition at line 413 of file G4VEmModel.hh.

◆ localTable

G4bool G4VEmModel::localTable
private

Definition at line 412 of file G4VEmModel.hh.

◆ lowLimit

G4double G4VEmModel::lowLimit
private

Definition at line 401 of file G4VEmModel.hh.

◆ name

const G4String G4VEmModel::name
private

Definition at line 397 of file G4VEmModel.hh.

◆ nsec

G4int G4VEmModel::nsec
private

Definition at line 438 of file G4VEmModel.hh.

◆ nSelectors

G4int G4VEmModel::nSelectors
private

Definition at line 416 of file G4VEmModel.hh.

◆ polarAngleLimit

G4double G4VEmModel::polarAngleLimit
private

Definition at line 405 of file G4VEmModel.hh.

◆ pParticleChange

G4VParticleChange* G4VEmModel::pParticleChange
protected

Definition at line 422 of file G4VEmModel.hh.

◆ secondaryThreshold

G4double G4VEmModel::secondaryThreshold
private

Definition at line 406 of file G4VEmModel.hh.

◆ theDensityFactor

const std::vector<G4double>* G4VEmModel::theDensityFactor
protected

Definition at line 424 of file G4VEmModel.hh.

◆ theDensityIdx

const std::vector<G4int>* G4VEmModel::theDensityIdx
protected

Definition at line 425 of file G4VEmModel.hh.

◆ theLPMflag

G4bool G4VEmModel::theLPMflag
private

Definition at line 407 of file G4VEmModel.hh.

◆ useAngularGenerator

G4bool G4VEmModel::useAngularGenerator
private

Definition at line 414 of file G4VEmModel.hh.

◆ xsec

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

Definition at line 439 of file G4VEmModel.hh.

◆ xSectionTable

G4PhysicsTable* G4VEmModel::xSectionTable
protected

Definition at line 423 of file G4VEmModel.hh.


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