#include <G4IonParametrisedLossModel.hh>

Inheritance diagram for G4IonParametrisedLossModel:

Collaboration diagram for G4IonParametrisedLossModel:

Public Member Functions
	G4IonParametrisedLossModel (const G4ParticleDefinition *particle=0, const G4String &name="ParamICRU73")

virtual	~G4IonParametrisedLossModel ()

virtual void	Initialise (const G4ParticleDefinition *, const G4DataVector &)

virtual G4double	MinEnergyCut (const G4ParticleDefinition , const G4MaterialCutsCouple )

virtual G4double	ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double, G4double, G4double, G4double, G4double)

virtual G4double	CrossSectionPerVolume (const G4Material , const G4ParticleDefinition , G4double, G4double, G4double)

virtual G4double	ComputeDEDXPerVolume (const G4Material , const G4ParticleDefinition , G4double, G4double)

G4double	ComputeLossForStep (const G4MaterialCutsCouple , const G4ParticleDefinition , G4double, G4double)

G4double	DeltaRayMeanEnergyTransferRate (const G4Material , const G4ParticleDefinition , G4double, G4double)

virtual void	SampleSecondaries (std::vector< G4DynamicParticle > , const G4MaterialCutsCouple , const G4DynamicParticle , G4double, G4double)

virtual G4double	GetChargeSquareRatio (const G4ParticleDefinition , const G4Material , G4double)

virtual G4double	GetParticleCharge (const G4ParticleDefinition , const G4Material , G4double)

virtual void	CorrectionsAlongStep (const G4MaterialCutsCouple , const G4DynamicParticle , G4double &, G4double &, G4double)

G4bool	AddDEDXTable (const G4String &name, G4VIonDEDXTable table, G4VIonDEDXScalingAlgorithm algorithm=0)

G4bool	RemoveDEDXTable (const G4String &name)

void	DeactivateICRU73Scaling ()

LossTableList::iterator	IsApplicable (const G4ParticleDefinition , const G4Material )

void	PrintDEDXTable (const G4ParticleDefinition , const G4Material , G4double, G4double, G4int, G4bool)

void	PrintDEDXTableHandlers (const G4ParticleDefinition , const G4Material , G4double, G4double, G4int, G4bool)

void	SetEnergyLossLimit (G4double ionEnergyLossLimit)

Public Member Functions inherited from G4VEmModel
	G4VEmModel (const G4String &nam)

virtual	~G4VEmModel ()

virtual void	InitialiseLocal (const G4ParticleDefinition , G4VEmModel masterModel)

virtual void	InitialiseForMaterial (const G4ParticleDefinition , const G4Material )

virtual void	InitialiseForElement (const G4ParticleDefinition *, G4int Z)

virtual G4double	GetPartialCrossSection (const G4Material , G4int, const G4ParticleDefinition , G4double)

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

virtual G4double	ChargeSquareRatio (const G4Track &)

virtual void	StartTracking (G4Track *)

virtual G4double	Value (const G4MaterialCutsCouple , const G4ParticleDefinition , G4double kineticEnergy)

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

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

virtual void	DefineForRegion (const G4Region *)

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

void	InitialiseElementSelectors (const G4ParticleDefinition *, const G4DataVector &)

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

void	SetElementSelectors (std::vector< G4EmElementSelector > )

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 G4Element *	SelectRandomAtom (const G4MaterialCutsCouple , const G4ParticleDefinition , G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)

const G4Element *	SelectRandomAtom (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)

G4ElementData *	GetElementData ()

G4PhysicsTable *	GetCrossSectionTable ()

G4VEmFluctuationModel *	GetModelOfFluctuations ()

G4VEmAngularDistribution *	GetAngularDistribution ()

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 G4String &	GetName () const

void	SetCurrentCouple (const G4MaterialCutsCouple *)

const G4Element *	GetCurrentElement () const

const G4Isotope *	GetCurrentIsotope () const

G4bool	IsLocked () const

void	SetLocked (G4bool)

Protected Member Functions
virtual G4double	MaxSecondaryEnergy (const G4ParticleDefinition *, G4double)

Protected Member Functions inherited from G4VEmModel
G4ParticleChangeForLoss *	GetParticleChangeForLoss ()

G4ParticleChangeForGamma *	GetParticleChangeForGamma ()

const G4MaterialCutsCouple *	CurrentCouple () const

void	SetCurrentElement (const G4Element *)

Private Types
typedef std::map< IonMatCouple, G4LPhysicsFreeVector * >	RangeEnergyTable

typedef std::map< IonMatCouple, G4LPhysicsFreeVector * >	EnergyRangeTable

Private Member Functions
void	UpdateDEDXCache (const G4ParticleDefinition , const G4Material , G4double cutEnergy)

void	UpdateRangeCache (const G4ParticleDefinition , const G4MaterialCutsCouple )

void	UpdateCache (const G4ParticleDefinition *)

void	BuildRangeVector (const G4ParticleDefinition , const G4MaterialCutsCouple )

G4IonParametrisedLossModel &	operator= (const G4IonParametrisedLossModel &right)

	G4IonParametrisedLossModel (const G4IonParametrisedLossModel &)

Private Attributes
G4BraggIonModel *	braggIonModel

G4BetheBlochModel *	betheBlochModel

LossTableList	lossTableList

RangeEnergyTable	r

EnergyRangeTable	E

G4double	lowerEnergyEdgeIntegr

G4double	upperEnergyEdgeIntegr

size_t	nmbBins

size_t	nmbSubBins

G4ParticleChangeForLoss *	particleChangeLoss

G4EmCorrections *	corrections

G4double	corrFactor

G4double	energyLossLimit

G4DataVector	cutEnergies

G4ParticleDefinition *	genericIon

G4double	genericIonPDGMass

const G4ParticleDefinition *	cacheParticle

G4double	cacheMass

G4double	cacheElecMassRatio

G4double	cacheChargeSquare

const G4ParticleDefinition *	rangeCacheParticle

const G4MaterialCutsCouple *	rangeCacheMatCutsCouple

G4PhysicsVector *	rangeCacheEnergyRange

G4PhysicsVector *	rangeCacheRangeEnergy

const G4ParticleDefinition *	dedxCacheParticle

const G4Material *	dedxCacheMaterial

G4double	dedxCacheEnergyCut

LossTableList::iterator	dedxCacheIter

G4double	dedxCacheTransitionEnergy

G4double	dedxCacheTransitionFactor

G4double	dedxCacheGenIonMassRatio

Additional Inherited Members
Protected Attributes inherited from G4VEmModel
G4ElementData *	fElementData

G4VParticleChange *	pParticleChange

G4PhysicsTable *	xSectionTable

const std::vector< G4double > *	theDensityFactor

const std::vector< G4int > *	theDensityIdx

size_t	idxTable

Static Protected Attributes inherited from G4VEmModel
static const G4double	inveplus = 1.0/CLHEP::eplus

Detailed Description

Definition at line 93 of file G4IonParametrisedLossModel.hh.

Member Typedef Documentation

◆ EnergyRangeTable

typedef std::map<IonMatCouple, G4LPhysicsFreeVector*> G4IonParametrisedLossModel::EnergyRangeTable

private

Definition at line 273 of file G4IonParametrisedLossModel.hh.

◆ RangeEnergyTable

typedef std::map<IonMatCouple, G4LPhysicsFreeVector*> G4IonParametrisedLossModel::RangeEnergyTable

private

Definition at line 270 of file G4IonParametrisedLossModel.hh.

Constructor & Destructor Documentation

◆ G4IonParametrisedLossModel() [1/2]

G4IonParametrisedLossModel::G4IonParametrisedLossModel	(	const G4ParticleDefinition *	particle = `0`,
		const G4String &	name = `"ParamICRU73"`
	)

Definition at line 105 of file G4IonParametrisedLossModel.cc.

   : G4VEmModel(nam),
     braggIonModel(0),
     betheBlochModel(0),
     nmbBins(90),
     nmbSubBins(100),
     particleChangeLoss(0),
     corrFactor(1.0),
     energyLossLimit(0.01),
     cutEnergies(0)
 {
   genericIon = G4GenericIon::Definition();
   genericIonPDGMass = genericIon -> GetPDGMass();
   corrections = G4LossTableManager::Instance() -> EmCorrections();
 
   // The upper limit of the current model is set to 100 TeV
   SetHighEnergyLimit(100.0 * TeV);
 
   // The Bragg ion and Bethe Bloch models are instantiated
   braggIonModel = new G4BraggIonModel();
   betheBlochModel = new G4BetheBlochModel();
 
   // By default ICRU 73 stopping power tables are loaded:
   AddDEDXTable("ICRU73",
            new G4IonStoppingData("ion_stopping_data/icru73"),
            new G4IonDEDXScalingICRU73());
 
   // The boundaries for the range tables are set
   lowerEnergyEdgeIntegr = 0.025 * MeV;
   upperEnergyEdgeIntegr = betheBlochModel -> HighEnergyLimit();
 
   // Cache parameters are set
   cacheParticle = 0;
   cacheMass = 0;
   cacheElecMassRatio = 0;
   cacheChargeSquare = 0;
 
   // Cache parameters are set
   rangeCacheParticle = 0;
   rangeCacheMatCutsCouple = 0;
   rangeCacheEnergyRange = 0;
   rangeCacheRangeEnergy = 0;
 
   // Cache parameters are set
   dedxCacheParticle = 0;
   dedxCacheMaterial = 0;
   dedxCacheEnergyCut = 0;
   dedxCacheIter = lossTableList.end();
   dedxCacheTransitionEnergy = 0.0;
   dedxCacheTransitionFactor = 0.0;
   dedxCacheGenIonMassRatio = 0.0;
 
   // default generator
   SetAngularDistribution(new G4DeltaAngle());
 }

Here is the call graph for this function:

◆ ~G4IonParametrisedLossModel()

G4IonParametrisedLossModel::~G4IonParametrisedLossModel ( )

virtual

Definition at line 165 of file G4IonParametrisedLossModel.cc.

                                                         {
 
   // dE/dx table objects are deleted and the container is cleared
   LossTableList::iterator iterTables = lossTableList.begin();
   LossTableList::iterator iterTables_end = lossTableList.end();
 
   for(;iterTables != iterTables_end; ++iterTables) { delete *iterTables; }
   lossTableList.clear();
 
   // range table
   RangeEnergyTable::iterator itr = r.begin();
   RangeEnergyTable::iterator itr_end = r.end();
   for(;itr != itr_end; ++itr) { delete itr->second; }
   r.clear();
 
   // inverse range
   EnergyRangeTable::iterator ite = E.begin();
   EnergyRangeTable::iterator ite_end = E.end();
   for(;ite != ite_end; ++ite) { delete ite->second; }
   E.clear();
 
 }

◆ G4IonParametrisedLossModel() [2/2]

G4IonParametrisedLossModel::G4IonParametrisedLossModel ( const G4IonParametrisedLossModel & )

private

Member Function Documentation

◆ AddDEDXTable()

G4bool G4IonParametrisedLossModel::AddDEDXTable	(	const G4String &	name,
		G4VIonDEDXTable *	table,
		G4VIonDEDXScalingAlgorithm *	algorithm = `0`
	)

Definition at line 1243 of file G4IonParametrisedLossModel.cc.

                                                                        {
 
   if(table == 0) {
      G4cerr << "G4IonParametrisedLossModel::AddDEDXTable() Cannot "
             << " add table: Invalid pointer."
             << G4endl;
 
      return false;
   }
 
   // Checking uniqueness of name
   LossTableList::iterator iter = lossTableList.begin();
   LossTableList::iterator iter_end = lossTableList.end();
 
   for(;iter != iter_end; iter++) {
      G4String tableName = (*iter) -> GetName();
 
      if(tableName == nam) {
         G4cerr << "G4IonParametrisedLossModel::AddDEDXTable() Cannot "
                << " add table: Name already exists."
                << G4endl;
 
         return false;
      }
   }
 
   G4VIonDEDXScalingAlgorithm* scalingAlgorithm = algorithm;
   if(scalingAlgorithm == 0)
      scalingAlgorithm = new G4VIonDEDXScalingAlgorithm;
 
   G4IonDEDXHandler* handler =
                       new G4IonDEDXHandler(table, scalingAlgorithm, nam);
 
   lossTableList.push_front(handler);
 
   return true;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ BuildRangeVector()

void G4IonParametrisedLossModel::BuildRangeVector	(	const G4ParticleDefinition *	particle,
		const G4MaterialCutsCouple *	matCutsCouple
	)

private

Definition at line 1065 of file G4IonParametrisedLossModel.cc.

                                                                 {
 
   G4double cutEnergy = DBL_MAX;
   size_t cutIndex = matCutsCouple -> GetIndex();
   cutEnergy = cutEnergies[cutIndex];
 
   const G4Material* material = matCutsCouple -> GetMaterial();
 
   G4double massRatio = genericIonPDGMass / particle -> GetPDGMass();
 
   G4double lowerEnergy = lowerEnergyEdgeIntegr / massRatio;
   G4double upperEnergy = upperEnergyEdgeIntegr / massRatio;
 
   G4double logLowerEnergyEdge = std::log(lowerEnergy);
   G4double logUpperEnergyEdge = std::log(upperEnergy);
 
   G4double logDeltaEnergy = (logUpperEnergyEdge - logLowerEnergyEdge) /
                                                         G4double(nmbBins);
 
   G4double logDeltaIntegr = logDeltaEnergy / G4double(nmbSubBins);
 
   G4LPhysicsFreeVector* energyRangeVector =
                               new G4LPhysicsFreeVector(nmbBins+1,
                                                        lowerEnergy,
                                                        upperEnergy);
 
   G4double dedxLow = ComputeDEDXPerVolume(material,
                                           particle,
                                           lowerEnergy,
                                           cutEnergy);
 
   G4double range = 2.0 * lowerEnergy / dedxLow;
 
   energyRangeVector -> PutValues(0, lowerEnergy, range);
 
   G4double logEnergy = std::log(lowerEnergy);
   for(size_t i = 1; i < nmbBins+1; i++) {
 
       G4double logEnergyIntegr = logEnergy;
 
       for(size_t j = 0; j < nmbSubBins; j++) {
 
           G4double binLowerBoundary = G4Exp(logEnergyIntegr);
           logEnergyIntegr += logDeltaIntegr;
 
           G4double binUpperBoundary = G4Exp(logEnergyIntegr);
           G4double deltaIntegr = binUpperBoundary - binLowerBoundary;
 
           G4double energyIntegr = binLowerBoundary + 0.5 * deltaIntegr;
 
           G4double dedxValue = ComputeDEDXPerVolume(material,
                                                     particle,
                                                     energyIntegr,
                             cutEnergy);
 
           if(dedxValue > 0.0) range += deltaIntegr / dedxValue;
 
 #ifdef PRINT_DEBUG_DETAILS
               G4cout << "   E = "<< energyIntegr/MeV
                      << " MeV -> dE = " << deltaIntegr/MeV
                      << " MeV -> dE/dx = " << dedxValue/MeV*mm
                      << " MeV/mm -> dE/(dE/dx) = " << deltaIntegr /
                                                                dedxValue / mm
                      << " mm -> range = " << range / mm
                      << " mm " << G4endl;
 #endif
       }
 
       logEnergy += logDeltaEnergy;
 
       G4double energy = G4Exp(logEnergy);
 
       energyRangeVector -> PutValues(i, energy, range);
 
 #ifdef PRINT_DEBUG_DETAILS
       G4cout << "G4IonParametrisedLossModel::BuildRangeVector() bin = "
              << i <<", E = "
              << energy / MeV << " MeV, R = "
              << range / mm << " mm"
              << G4endl;
 #endif
 
   }
   energyRangeVector -> SetSpline(true);
 
   G4double lowerRangeEdge =
                     energyRangeVector -> Value(lowerEnergy);
   G4double upperRangeEdge =
                     energyRangeVector -> Value(upperEnergy);
 
   G4LPhysicsFreeVector* rangeEnergyVector
                       = new G4LPhysicsFreeVector(nmbBins+1,
                                                  lowerRangeEdge,
                                                  upperRangeEdge);
 
   for(size_t i = 0; i < nmbBins+1; i++) {
       G4double energy = energyRangeVector -> Energy(i);
       rangeEnergyVector ->
              PutValues(i, energyRangeVector -> Value(energy), energy);
   }
 
   rangeEnergyVector -> SetSpline(true);
 
 #ifdef PRINT_DEBUG_TABLES
   G4cout << *energyLossVector
          << *energyRangeVector
          << *rangeEnergyVector << G4endl;
 #endif
 
   IonMatCouple ionMatCouple = std::make_pair(particle, matCutsCouple);
 
   E[ionMatCouple] = energyRangeVector;
   r[ionMatCouple] = rangeEnergyVector;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ ComputeCrossSectionPerAtom()

G4double G4IonParametrisedLossModel::ComputeCrossSectionPerAtom	(	const G4ParticleDefinition *	particle,
		G4double	kineticEnergy,
		G4double	atomicNumber,
		G4double	,
		G4double	cutEnergy,
		G4double	maxKinEnergy
	)

virtual

Reimplemented from G4VEmModel.

Definition at line 374 of file G4IonParametrisedLossModel.cc.

                                                           {
 
   // ############## Cross section per atom  ################################
   // Function computes ionization cross section per atom
   //
   // See Geant4 physics reference manual (version 9.1), section 9.1.3
   //
   // Ref.: W.M. Yao et al, Jour. of Phys. G 33 (2006) 1.
   //       B. Rossi, High energy particles, New York, NY: Prentice-Hall (1952).
   //
   // (Implementation adapted from G4BraggIonModel)
 
   G4double crosssection     = 0.0;
   G4double tmax      = MaxSecondaryEnergy(particle, kineticEnergy);
   G4double maxEnergy = std::min(tmax, maxKinEnergy);
 
   if(cutEnergy < tmax) {
 
      G4double energy  = kineticEnergy + cacheMass;
      G4double betaSquared  = kineticEnergy *
                                      (energy + cacheMass) / (energy * energy);
 
      crosssection = 1.0 / cutEnergy - 1.0 / maxEnergy -
                          betaSquared * std::log(maxEnergy / cutEnergy) / tmax;
 
      crosssection *= twopi_mc2_rcl2 * cacheChargeSquare / betaSquared;
   }
 
 #ifdef PRINT_DEBUG_CS
   G4cout << "########################################################"
          << G4endl
          << "# G4IonParametrisedLossModel::ComputeCrossSectionPerAtom"
          << G4endl
          << "# particle =" << particle -> GetParticleName()
          <<  G4endl
          << "# cut(MeV) = " << cutEnergy/MeV
          << G4endl;
 
   G4cout << "#"
          << std::setw(13) << std::right << "E(MeV)"
          << std::setw(14) << "CS(um)"
          << std::setw(14) << "E_max_sec(MeV)"
          << G4endl
          << "# ------------------------------------------------------"
          << G4endl;
 
   G4cout << std::setw(14) << std::right << kineticEnergy / MeV
          << std::setw(14) << crosssection / (um * um)
          << std::setw(14) << tmax / MeV
          << G4endl;
 #endif
 
   crosssection *= atomicNumber;
 
   return crosssection;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ ComputeDEDXPerVolume()

G4double G4IonParametrisedLossModel::ComputeDEDXPerVolume	(	const G4Material *	material,
		const G4ParticleDefinition *	particle,
		G4double	kineticEnergy,
		G4double	cutEnergy
	)

virtual

Reimplemented from G4VEmModel.

Definition at line 458 of file G4IonParametrisedLossModel.cc.

                                           {
 
   // ############## dE/dx ##################################################
   // Function computes dE/dx values, where following rules are adopted:
   //   A. If the ion-material pair is covered by any native ion data
   //      parameterisation, then:
   //      * This parameterization is used for energies below a given energy
   //        limit,
   //      * whereas above the limit the Bethe-Bloch model is applied, in
   //        combination with an effective charge estimate and high order
   //        correction terms.
   //      A smoothing procedure is applied to dE/dx values computed with
   //      the second approach. The smoothing factor is based on the dE/dx
   //      values of both approaches at the transition energy (high order
   //      correction terms are included in the calculation of the transition
   //      factor).
   //   B. If the particle is a generic ion, the BraggIon and Bethe-Bloch
   //      models are used and a smoothing procedure is applied to values
   //      obtained with the second approach.
   //   C. If the ion-material is not covered by any ion data parameterization
   //      then:
   //      * The BraggIon model is used for energies below a given energy
   //        limit,
   //      * whereas above the limit the Bethe-Bloch model is applied, in
   //        combination with an effective charge estimate and high order
   //        correction terms.
   //      Also in this case, a smoothing procedure is applied to dE/dx values
   //      computed with the second model.
 
   G4double dEdx = 0.0;
   UpdateDEDXCache(particle, material, cutEnergy);
 
   LossTableList::iterator iter = dedxCacheIter;
 
   if(iter != lossTableList.end()) {
 
      G4double transitionEnergy = dedxCacheTransitionEnergy;
 
      if(transitionEnergy > kineticEnergy) {
 
         dEdx = (*iter) -> GetDEDX(particle, material, kineticEnergy);
 
         G4double dEdxDeltaRays = DeltaRayMeanEnergyTransferRate(material,
                                            particle,
                                            kineticEnergy,
                                            cutEnergy);
         dEdx -= dEdxDeltaRays;
      }
      else {
         G4double massRatio = dedxCacheGenIonMassRatio;
 
         G4double chargeSquare =
                        GetChargeSquareRatio(particle, material, kineticEnergy);
 
         G4double scaledKineticEnergy = kineticEnergy * massRatio;
         G4double scaledTransitionEnergy = transitionEnergy * massRatio;
 
         G4double lowEnergyLimit = betheBlochModel -> LowEnergyLimit();
 
         if(scaledTransitionEnergy >= lowEnergyLimit) {
 
            dEdx = betheBlochModel -> ComputeDEDXPerVolume(
                                       material, genericIon,
                           scaledKineticEnergy, cutEnergy);
 
            dEdx *= chargeSquare;
 
            dEdx += corrections -> ComputeIonCorrections(particle,
                                                  material, kineticEnergy);
 
            G4double factor = 1.0 + dedxCacheTransitionFactor /
                                                            kineticEnergy;
 
        dEdx *= factor;
     }
      }
   }
   else {
      G4double massRatio = 1.0;
      G4double chargeSquare = 1.0;
 
      if(particle != genericIon) {
 
         chargeSquare = GetChargeSquareRatio(particle, material, kineticEnergy);
         massRatio = genericIonPDGMass / particle -> GetPDGMass();
      }
 
      G4double scaledKineticEnergy = kineticEnergy * massRatio;
 
      G4double lowEnergyLimit = betheBlochModel -> LowEnergyLimit();
      if(scaledKineticEnergy < lowEnergyLimit) {
         dEdx = braggIonModel -> ComputeDEDXPerVolume(
                                       material, genericIon,
                           scaledKineticEnergy, cutEnergy);
 
         dEdx *= chargeSquare;
      }
      else {
         G4double dEdxLimitParam = braggIonModel -> ComputeDEDXPerVolume(
                                       material, genericIon,
                           lowEnergyLimit, cutEnergy);
 
         G4double dEdxLimitBetheBloch = betheBlochModel -> ComputeDEDXPerVolume(
                                       material, genericIon,
                           lowEnergyLimit, cutEnergy);
 
         if(particle != genericIon) {
            G4double chargeSquareLowEnergyLimit =
                GetChargeSquareRatio(particle, material,
                                     lowEnergyLimit / massRatio);
 
            dEdxLimitParam *= chargeSquareLowEnergyLimit;
            dEdxLimitBetheBloch *= chargeSquareLowEnergyLimit;
 
            dEdxLimitBetheBloch +=
                     corrections -> ComputeIonCorrections(particle,
                                       material, lowEnergyLimit / massRatio);
         }
 
         G4double factor = (1.0 + (dEdxLimitParam/dEdxLimitBetheBloch - 1.0)
                                * lowEnergyLimit / scaledKineticEnergy);
 
         dEdx = betheBlochModel -> ComputeDEDXPerVolume(
                                       material, genericIon,
                           scaledKineticEnergy, cutEnergy);
 
         dEdx *= chargeSquare;
 
         if(particle != genericIon) {
            dEdx += corrections -> ComputeIonCorrections(particle,
                                       material, kineticEnergy);
         }
 
         dEdx *= factor;
      }
 
   }
 
   if (dEdx < 0.0) dEdx = 0.0;
 
   return dEdx;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ ComputeLossForStep()

G4double G4IonParametrisedLossModel::ComputeLossForStep	(	const G4MaterialCutsCouple *	matCutsCouple,
		const G4ParticleDefinition *	particle,
		G4double	kineticEnergy,
		G4double	stepLength
	)

Definition at line 1184 of file G4IonParametrisedLossModel.cc.

                                           {
 
   G4double loss = 0.0;
 
   UpdateRangeCache(particle, matCutsCouple);
 
   G4PhysicsVector* energyRange = rangeCacheEnergyRange;
   G4PhysicsVector* rangeEnergy = rangeCacheRangeEnergy;
 
   if(energyRange != 0 && rangeEnergy != 0) {
 
      G4double lowerEnEdge = energyRange -> Energy( 0 );
      G4double lowerRangeEdge = rangeEnergy -> Energy( 0 );
 
      // Computing range for pre-step kinetic energy:
      G4double range = energyRange -> Value(kineticEnergy);
 
      // Energy below vector boundary:
      if(kineticEnergy < lowerEnEdge) {
 
         range =  energyRange -> Value(lowerEnEdge);
         range *= std::sqrt(kineticEnergy / lowerEnEdge);
      }
 
 #ifdef PRINT_DEBUG
      G4cout << "G4IonParametrisedLossModel::ComputeLossForStep() range = "
             << range / mm << " mm, step = " << stepLength / mm << " mm"
             << G4endl;
 #endif
 
      // Remaining range:
      G4double remRange = range - stepLength;
 
      // If range is smaller than step length, the loss is set to kinetic
      // energy
      if(remRange < 0.0) loss = kineticEnergy;
      else if(remRange < lowerRangeEdge) {
 
         G4double ratio = remRange / lowerRangeEdge;
         loss = kineticEnergy - ratio * ratio * lowerEnEdge;
      }
      else {
 
         G4double energy = rangeEnergy -> Value(range - stepLength);
         loss = kineticEnergy - energy;
      }
   }
 
   if(loss < 0.0) loss = 0.0;
 
   return loss;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ CorrectionsAlongStep()

void G4IonParametrisedLossModel::CorrectionsAlongStep	(	const G4MaterialCutsCouple *	couple,
		const G4DynamicParticle *	dynamicParticle,
		G4double &	eloss,
		G4double &	,
		G4double	length
	)

virtual

Reimplemented from G4VEmModel.

Definition at line 918 of file G4IonParametrisedLossModel.cc.

                                               {
 
   // ############## Corrections for along step energy loss calculation ######
   // The computed energy loss (due to electronic stopping) is overwritten
   // by this function if an ion data parameterization is available for the
   // current ion-material pair.
   // No action on the energy loss (due to electronic stopping) is performed
   // if no parameterization is available. In this case the original
   // generic ion tables (in combination with the effective charge) are used
   // in the along step DoIt function.
   //
   //
   // (Implementation partly adapted from G4BraggIonModel/G4BetheBlochModel)
 
   const G4ParticleDefinition* particle = dynamicParticle -> GetDefinition();
   const G4Material* material = couple -> GetMaterial();
 
   G4double kineticEnergy = dynamicParticle -> GetKineticEnergy();
 
   if(kineticEnergy == eloss) { return; }
 
   G4double cutEnergy = DBL_MAX;
   size_t cutIndex = couple -> GetIndex();
   cutEnergy = cutEnergies[cutIndex];
 
   UpdateDEDXCache(particle, material, cutEnergy);
 
   LossTableList::iterator iter = dedxCacheIter;
 
   // If parameterization for ions is available the electronic energy loss
   // is overwritten
   if(iter != lossTableList.end()) {
 
      // The energy loss is calculated using the ComputeDEDXPerVolume function
      // and the step length (it is assumed that dE/dx does not change
      // considerably along the step)
      eloss =
         length * ComputeDEDXPerVolume(material, particle,
                                       kineticEnergy, cutEnergy);
 
 #ifdef PRINT_DEBUG
   G4cout.precision(6);
   G4cout << "########################################################"
          << G4endl
          << "# G4IonParametrisedLossModel::CorrectionsAlongStep"
          << G4endl
          << "# cut(MeV) = " << cutEnergy/MeV
          << G4endl;
 
   G4cout << "#"
          << std::setw(13) << std::right << "E(MeV)"
          << std::setw(14) << "l(um)"
          << std::setw(14) << "l*dE/dx(MeV)"
          << std::setw(14) << "(l*dE/dx)/E"
          << G4endl
          << "# ------------------------------------------------------"
          << G4endl;
 
   G4cout << std::setw(14) << std::right << kineticEnergy / MeV
          << std::setw(14) << length / um
          << std::setw(14) << eloss / MeV
          << std::setw(14) << eloss / kineticEnergy * 100.0
          << G4endl;
 #endif
 
      // If the energy loss exceeds a certain fraction of the kinetic energy
      // (the fraction is indicated by the parameter "energyLossLimit") then
      // the range tables are used to derive a more accurate value of the
      // energy loss
      if(eloss > energyLossLimit * kineticEnergy) {
 
         eloss = ComputeLossForStep(couple, particle,
                                    kineticEnergy,length);
 
 #ifdef PRINT_DEBUG
   G4cout << "# Correction applied:"
          << G4endl;
 
   G4cout << std::setw(14) << std::right << kineticEnergy / MeV
          << std::setw(14) << length / um
          << std::setw(14) << eloss / MeV
          << std::setw(14) << eloss / kineticEnergy * 100.0
          << G4endl;
 #endif
 
      }
   }
 
   // For all corrections below a kinetic energy between the Pre- and
   // Post-step energy values is used
   G4double energy = kineticEnergy - eloss * 0.5;
   if(energy < 0.0) energy = kineticEnergy * 0.5;
 
   G4double chargeSquareRatio = corrections ->
                                      EffectiveChargeSquareRatio(particle,
                                                                 material,
                                                                 energy);
   GetModelOfFluctuations() -> SetParticleAndCharge(particle,
                                                    chargeSquareRatio);
 
   // A correction is applied considering the change of the effective charge
   // along the step (the parameter "corrFactor" refers to the effective
   // charge at the beginning of the step). Note: the correction is not
   // applied for energy loss values deriving directly from parameterized
   // ion stopping power tables
   G4double transitionEnergy = dedxCacheTransitionEnergy;
 
   if(iter != lossTableList.end() && transitionEnergy < kineticEnergy) {
      chargeSquareRatio *= corrections -> EffectiveChargeCorrection(particle,
                                                                 material,
                                                                 energy);
 
      G4double chargeSquareRatioCorr = chargeSquareRatio/corrFactor;
      eloss *= chargeSquareRatioCorr;
   }
   else if (iter == lossTableList.end()) {
 
      chargeSquareRatio *= corrections -> EffectiveChargeCorrection(particle,
                                                                 material,
                                                                 energy);
 
      G4double chargeSquareRatioCorr = chargeSquareRatio/corrFactor;
      eloss *= chargeSquareRatioCorr;
   }
 
   // Ion high order corrections are applied if the current model does not
   // overwrite the energy loss (i.e. when the effective charge approach is
   // used)
   if(iter == lossTableList.end()) {
 
      G4double scaledKineticEnergy = kineticEnergy * dedxCacheGenIonMassRatio;
      G4double lowEnergyLimit = betheBlochModel -> LowEnergyLimit();
 
      // Corrections are only applied in the Bethe-Bloch energy region
      if(scaledKineticEnergy > lowEnergyLimit)
        eloss += length *
             corrections -> IonHighOrderCorrections(particle, couple, energy);
   }
 }

Here is the call graph for this function:

◆ CrossSectionPerVolume()

G4double G4IonParametrisedLossModel::CrossSectionPerVolume	(	const G4Material *	material,
		const G4ParticleDefinition *	particle,
		G4double	kineticEnergy,
		G4double	cutEnergy,
		G4double	maxEnergy
	)

virtual

Reimplemented from G4VEmModel.

Definition at line 439 of file G4IonParametrisedLossModel.cc.

                                                        {
 
   G4double nbElecPerVolume = material -> GetTotNbOfElectPerVolume();
   G4double cross = ComputeCrossSectionPerAtom(particle,
                                               kineticEnergy,
                                               nbElecPerVolume, 0,
                                               cutEnergy,
                                               maxEnergy);
 
   return cross;
 }

Here is the call graph for this function:

◆ DeactivateICRU73Scaling()

void G4IonParametrisedLossModel::DeactivateICRU73Scaling ( )

Definition at line 1325 of file G4IonParametrisedLossModel.cc.

                                                          {
 
   RemoveDEDXTable("ICRU73");
   AddDEDXTable("ICRU73", new G4IonStoppingData("ion_stopping_data/icru73"));
 }

Here is the call graph for this function:

◆ DeltaRayMeanEnergyTransferRate()

G4double G4IonParametrisedLossModel::DeltaRayMeanEnergyTransferRate	(	const G4Material *	,
		const G4ParticleDefinition *	,
		G4double	,
		G4double
	)

inline

Here is the caller graph for this function:

◆ GetChargeSquareRatio()

G4double G4IonParametrisedLossModel::GetChargeSquareRatio	(	const G4ParticleDefinition *	particle,
		const G4Material *	material,
		G4double	kineticEnergy
	)

virtual

Reimplemented from G4VEmModel.

Definition at line 228 of file G4IonParametrisedLossModel.cc.

                                                      {    // Kinetic energy
 
   G4double chargeSquareRatio = corrections ->
                                      EffectiveChargeSquareRatio(particle,
                                                                 material,
                                                                 kineticEnergy);
   corrFactor = chargeSquareRatio *
                        corrections -> EffectiveChargeCorrection(particle,
                                                                 material,
                                                                 kineticEnergy);
   return corrFactor;
 }

Here is the caller graph for this function:

◆ GetParticleCharge()

G4double G4IonParametrisedLossModel::GetParticleCharge	(	const G4ParticleDefinition *	particle,
		const G4Material *	material,
		G4double	kineticEnergy
	)

virtual

Reimplemented from G4VEmModel.

Definition at line 246 of file G4IonParametrisedLossModel.cc.

                                                      {   // Kinetic energy
 
   return corrections -> GetParticleCharge(particle, material, kineticEnergy);
 }

◆ Initialise()

void G4IonParametrisedLossModel::Initialise	(	const G4ParticleDefinition *	particle,
		const G4DataVector &	cuts
	)

virtual

Implements G4VEmModel.

Definition at line 256 of file G4IonParametrisedLossModel.cc.

                                                                  {
 
   // Cached parameters are reset
   cacheParticle = 0;
   cacheMass = 0;
   cacheElecMassRatio = 0;
   cacheChargeSquare = 0;
 
   // Cached parameters are reset
   rangeCacheParticle = 0;
   rangeCacheMatCutsCouple = 0;
   rangeCacheEnergyRange = 0;
   rangeCacheRangeEnergy = 0;
 
   // Cached parameters are reset
   dedxCacheParticle = 0;
   dedxCacheMaterial = 0;
   dedxCacheEnergyCut = 0;
   dedxCacheIter = lossTableList.end();
   dedxCacheTransitionEnergy = 0.0;
   dedxCacheTransitionFactor = 0.0;
   dedxCacheGenIonMassRatio = 0.0;
 
   // The cache of loss tables is cleared
   LossTableList::iterator iterTables = lossTableList.begin();
   LossTableList::iterator iterTables_end = lossTableList.end();
 
   for(;iterTables != iterTables_end; iterTables++)
                                        (*iterTables) -> ClearCache();
 
   // Range vs energy and energy vs range vectors from previous runs are
   // cleared
   RangeEnergyTable::iterator iterRange = r.begin();
   RangeEnergyTable::iterator iterRange_end = r.end();
 
   for(;iterRange != iterRange_end; iterRange++) {
     delete iterRange->second;
   }
   r.clear();
 
   EnergyRangeTable::iterator iterEnergy = E.begin();
   EnergyRangeTable::iterator iterEnergy_end = E.end();
 
   for(;iterEnergy != iterEnergy_end; iterEnergy++) {
     delete iterEnergy->second;
   }
   E.clear();
 
   // The cut energies are (re)loaded
   size_t size = cuts.size();
   cutEnergies.clear();
   for(size_t i = 0; i < size; i++) cutEnergies.push_back(cuts[i]);
 
   // All dE/dx vectors are built
   const G4ProductionCutsTable* coupleTable=
                      G4ProductionCutsTable::GetProductionCutsTable();
   size_t nmbCouples = coupleTable -> GetTableSize();
 
 #ifdef PRINT_TABLE_BUILT
     G4cout << "G4IonParametrisedLossModel::Initialise():"
            << " Building dE/dx vectors:"
            << G4endl;
 #endif
 
   for (size_t i = 0; i < nmbCouples; i++) {
 
     const G4MaterialCutsCouple* couple =
                                      coupleTable -> GetMaterialCutsCouple(i);
 
     const G4Material* material = couple -> GetMaterial();
     //    G4ProductionCuts* productionCuts = couple -> GetProductionCuts();
 
     for(G4int atomicNumberIon = 3; atomicNumberIon < 102; atomicNumberIon++) {
 
        LossTableList::iterator iter = lossTableList.begin();
        LossTableList::iterator iter_end = lossTableList.end();
 
        for(;iter != iter_end; iter++) {
 
           if(*iter == 0) {
               G4cout << "G4IonParametrisedLossModel::Initialise():"
                      << " Skipping illegal table."
                      << G4endl;
           }
 
           G4bool isApplicable =
                     (*iter) -> BuildDEDXTable(atomicNumberIon, material);
           if(isApplicable) {
 
 #ifdef PRINT_TABLE_BUILT
              G4cout << "  Atomic Number Ion = " << atomicNumberIon
                     << ", Material = " << material -> GetName()
                     << ", Table = " << (*iter) -> GetName()
                     << G4endl;
 #endif
              break;
       }
        }
     }
   }
 
   // The particle change object
   if(! particleChangeLoss) {
     particleChangeLoss = GetParticleChangeForLoss();
     braggIonModel -> SetParticleChange(particleChangeLoss, 0);
     betheBlochModel -> SetParticleChange(particleChangeLoss, 0);
   }
 
   // The G4BraggIonModel and G4BetheBlochModel instances are initialised with
   // the same settings as the current model:
   braggIonModel -> Initialise(particle, cuts);
   betheBlochModel -> Initialise(particle, cuts);
 }

Here is the call graph for this function:

◆ IsApplicable()

LossTableList::iterator G4IonParametrisedLossModel::IsApplicable	(	const G4ParticleDefinition *	,
		const G4Material *
	)

inline

Here is the caller graph for this function:

◆ MaxSecondaryEnergy()

G4double G4IonParametrisedLossModel::MaxSecondaryEnergy	(	const G4ParticleDefinition *	particle,
		G4double	kineticEnergy
	)

protectedvirtual

Reimplemented from G4VEmModel.

Definition at line 200 of file G4IonParametrisedLossModel.cc.

                                                      {
 
   // ############## Maximum energy of secondaries ##########################
   // Function computes maximum energy of secondary electrons which are
   // released by an ion
   //
   // See Geant4 physics reference manual (version 9.1), section 9.1.1
   //
   // Ref.: W.M. Yao et al, Jour. of Phys. G 33 (2006) 1.
   //       C.Caso et al. (Part. Data Group), Europ. Phys. Jour. C 3 1 (1998).
   //       B. Rossi, High energy particles, New York, NY: Prentice-Hall (1952).
   //
   // (Implementation adapted from G4BraggIonModel)
 
   if(particle != cacheParticle) UpdateCache(particle);
 
   G4double tau  = kineticEnergy/cacheMass;
   G4double tmax = 2.0 * electron_mass_c2 * tau * (tau + 2.) /
                   (1. + 2.0 * (tau + 1.) * cacheElecMassRatio +
                   cacheElecMassRatio * cacheElecMassRatio);
 
   return tmax;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ MinEnergyCut()

G4double G4IonParametrisedLossModel::MinEnergyCut	(	const G4ParticleDefinition *	,
		const G4MaterialCutsCouple *	couple
	)

virtual

Reimplemented from G4VEmModel.

Definition at line 190 of file G4IonParametrisedLossModel.cc.

                                                                            {
 
   return couple -> GetMaterial() -> GetIonisation() ->
                                                   GetMeanExcitationEnergy();
 }

◆ operator=()

G4IonParametrisedLossModel& G4IonParametrisedLossModel::operator= ( const G4IonParametrisedLossModel & right )

private

◆ PrintDEDXTable()

void G4IonParametrisedLossModel::PrintDEDXTable	(	const G4ParticleDefinition *	particle,
		const G4Material *	material,
		G4double	lowerBoundary,
		G4double	upperBoundary,
		G4int	numBins,
		G4bool	logScaleEnergy
	)

Definition at line 607 of file G4IonParametrisedLossModel.cc.

                                           {     // Logarithmic scaling of energy
 
   G4double atomicMassNumber = particle -> GetAtomicMass();
   G4double materialDensity = material -> GetDensity();
 
   G4cout << "# dE/dx table for " << particle -> GetParticleName()
          << " in material " << material -> GetName()
          << " of density " << materialDensity / g * cm3
          << " g/cm3"
          << G4endl
          << "# Projectile mass number A1 = " << atomicMassNumber
          << G4endl
          << "# ------------------------------------------------------"
          << G4endl;
   G4cout << "#"
          << std::setw(13) << std::right << "E"
          << std::setw(14) << "E/A1"
          << std::setw(14) << "dE/dx"
          << std::setw(14) << "1/rho*dE/dx"
          << G4endl;
   G4cout << "#"
          << std::setw(13) << std::right << "(MeV)"
          << std::setw(14) << "(MeV)"
          << std::setw(14) << "(MeV/cm)"
          << std::setw(14) << "(MeV*cm2/mg)"
          << G4endl
          << "# ------------------------------------------------------"
          << G4endl;
 
   G4double energyLowerBoundary = lowerBoundary * atomicMassNumber;
   G4double energyUpperBoundary = upperBoundary * atomicMassNumber;
 
   if(logScaleEnergy) {
 
      energyLowerBoundary = std::log(energyLowerBoundary);
      energyUpperBoundary = std::log(energyUpperBoundary);
   }
 
   G4double deltaEnergy = (energyUpperBoundary - energyLowerBoundary) /
                                                            G4double(nmbBins);
 
   for(int i = 0; i < numBins + 1; i++) {
 
       G4double energy = energyLowerBoundary + i * deltaEnergy;
       if(logScaleEnergy) energy = G4Exp(energy);
 
       G4double dedx = ComputeDEDXPerVolume(material, particle, energy, DBL_MAX);
       G4cout.precision(6);
       G4cout << std::setw(14) << std::right << energy / MeV
              << std::setw(14) << energy / atomicMassNumber / MeV
              << std::setw(14) << dedx / MeV * cm
              << std::setw(14) << dedx / materialDensity / (MeV*cm2/(0.001*g))
              << G4endl;
   }
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ PrintDEDXTableHandlers()

void G4IonParametrisedLossModel::PrintDEDXTableHandlers	(	const G4ParticleDefinition *	particle,
		const G4Material *	material,
		G4double	lowerBoundary,
		G4double	upperBoundary,
		G4int	numBins,
		G4bool	logScaleEnergy
	)

Definition at line 671 of file G4IonParametrisedLossModel.cc.

                                           {     // Logarithmic scaling of energy
 
   LossTableList::iterator iter = lossTableList.begin();
   LossTableList::iterator iter_end = lossTableList.end();
 
   for(;iter != iter_end; iter++) {
       G4bool isApplicable =  (*iter) -> IsApplicable(particle, material);
       if(isApplicable) {
     (*iter) -> PrintDEDXTable(particle, material,
                                   lowerBoundary, upperBoundary,
                                   numBins,logScaleEnergy);
         break;
       }
   }
 }

Here is the call graph for this function:

◆ RemoveDEDXTable()

G4bool G4IonParametrisedLossModel::RemoveDEDXTable ( const G4String & name )

Definition at line 1286 of file G4IonParametrisedLossModel.cc.

                                       {
 
   LossTableList::iterator iter = lossTableList.begin();
   LossTableList::iterator iter_end = lossTableList.end();
 
   for(;iter != iter_end; iter++) {
      G4String tableName = (*iter) -> GetName();
 
      if(tableName == nam) {
         delete (*iter);
 
         // Remove from table list
         lossTableList.erase(iter);
 
         // Range vs energy and energy vs range vectors are cleared
         RangeEnergyTable::iterator iterRange = r.begin();
         RangeEnergyTable::iterator iterRange_end = r.end();
 
         for(;iterRange != iterRange_end; iterRange++)
                                             delete iterRange -> second;
         r.clear();
 
         EnergyRangeTable::iterator iterEnergy = E.begin();
         EnergyRangeTable::iterator iterEnergy_end = E.end();
 
         for(;iterEnergy != iterEnergy_end; iterEnergy++)
                                             delete iterEnergy -> second;
         E.clear();
 
         return true;
      }
   }
 
   return false;
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ SampleSecondaries()

void G4IonParametrisedLossModel::SampleSecondaries	(	std::vector< G4DynamicParticle >	secondaries,
		const G4MaterialCutsCouple *	couple,
		const G4DynamicParticle *	particle,
		G4double	cutKinEnergySec,
		G4double	userMaxKinEnergySec
	)

virtual

Implements G4VEmModel.

Definition at line 695 of file G4IonParametrisedLossModel.cc.

                                                {
 
 
   // ############## Sampling of secondaries #################################
   // The probability density function (pdf) of the kinetic energy T of a
   // secondary electron may be written as:
   //    pdf(T) = f(T) * g(T)
   // where
   //    f(T) = (Tmax - Tcut) / (Tmax * Tcut) * (1 / T^2)
   //    g(T) = 1 - beta^2 * T / Tmax
   // where Tmax is the maximum kinetic energy of the secondary, Tcut
   // is the lower energy cut and beta is the kinetic energy of the
   // projectile.
   //
   // Sampling of the kinetic energy of a secondary electron:
   //  1) T0 is sampled from f(T) using the cumulated distribution function
   //     F(T) = int_Tcut^T f(T')dT'
   //  2) T is accepted or rejected by evaluating the rejection function g(T)
   //     at the sampled energy T0 against a randomly sampled value
   //
   //
   // See Geant4 physics reference manual (version 9.1), section 9.1.4
   //
   //
   // Reference pdf: W.M. Yao et al, Jour. of Phys. G 33 (2006) 1.
   //
   // (Implementation adapted from G4BraggIonModel)
 
   G4double rossiMaxKinEnergySec = MaxSecondaryKinEnergy(particle);
   G4double maxKinEnergySec =
                         std::min(rossiMaxKinEnergySec, userMaxKinEnergySec);
 
   if(cutKinEnergySec >= maxKinEnergySec) return;
 
   G4double kineticEnergy = particle -> GetKineticEnergy();
 
   G4double energy  = kineticEnergy + cacheMass;
   G4double betaSquared = kineticEnergy * (energy + cacheMass)
     / (energy * energy);
 
   G4double kinEnergySec;
   G4double grej;
 
   do {
 
     // Sampling kinetic energy from f(T) (using F(T)):
     G4double xi = G4UniformRand();
     kinEnergySec = cutKinEnergySec * maxKinEnergySec /
                         (maxKinEnergySec * (1.0 - xi) + cutKinEnergySec * xi);
 
     // Deriving the value of the rejection function at the obtained kinetic
     // energy:
     grej = 1.0 - betaSquared * kinEnergySec / rossiMaxKinEnergySec;
 
     if(grej > 1.0) {
         G4cout << "G4IonParametrisedLossModel::SampleSecondary Warning: "
                << "Majorant 1.0 < "
                << grej << " for e= " << kinEnergySec
                << G4endl;
     }
 
   } while( G4UniformRand() >= grej );
 
   const G4Material* mat =  couple->GetMaterial();
   G4int Z = SelectRandomAtomNumber(mat);
 
   const G4ParticleDefinition* electron = G4Electron::Electron();
 
   G4DynamicParticle* delta = new G4DynamicParticle(electron,
     GetAngularDistribution()->SampleDirection(particle, kinEnergySec,
                                               Z, mat),
                                                    kinEnergySec);
 
 
   secondaries->push_back(delta);
 
   // Change kinematics of primary particle
   G4ThreeVector direction = particle ->GetMomentumDirection();
   G4double totalMomentum = std::sqrt(kineticEnergy*(energy + cacheMass));
 
   G4ThreeVector finalP = totalMomentum*direction - delta->GetMomentum();
   finalP               = finalP.unit();
 
   kineticEnergy       -= kinEnergySec;
 
   particleChangeLoss->SetProposedKineticEnergy(kineticEnergy);
   particleChangeLoss->SetProposedMomentumDirection(finalP);
 }

Here is the call graph for this function:

◆ SetEnergyLossLimit()

void G4IonParametrisedLossModel::SetEnergyLossLimit ( G4double ionEnergyLossLimit )

inline

◆ UpdateCache()

void G4IonParametrisedLossModel::UpdateCache ( const G4ParticleDefinition * )

inlineprivate

Here is the caller graph for this function:

◆ UpdateDEDXCache()

void G4IonParametrisedLossModel::UpdateDEDXCache	(	const G4ParticleDefinition *	particle,
		const G4Material *	material,
		G4double	cutEnergy
	)

private

Definition at line 830 of file G4IonParametrisedLossModel.cc.

                                          {
 
   // ############## Caching ##################################################
   // If the ion-material combination is covered by any native ion data
   // parameterisation (for low energies), a transition factor is computed
   // which is applied to Bethe-Bloch results at higher energies to
   // guarantee a smooth transition.
   // This factor only needs to be calculated for the first step an ion
   // performs inside a certain material.
 
   if(particle == dedxCacheParticle &&
      material == dedxCacheMaterial &&
      cutEnergy == dedxCacheEnergyCut) {
   }
   else {
 
      dedxCacheParticle = particle;
      dedxCacheMaterial = material;
      dedxCacheEnergyCut = cutEnergy;
 
      G4double massRatio = genericIonPDGMass / particle -> GetPDGMass();
      dedxCacheGenIonMassRatio = massRatio;
 
      LossTableList::iterator iter = IsApplicable(particle, material);
      dedxCacheIter = iter;
 
      // If any table is applicable, the transition factor is computed:
      if(iter != lossTableList.end()) {
 
         // Retrieving the transition energy from the parameterisation table
         G4double transitionEnergy =
                  (*iter) -> GetUpperEnergyEdge(particle, material);
         dedxCacheTransitionEnergy = transitionEnergy;
 
         // Computing dE/dx from low-energy parameterisation at
         // transition energy
         G4double dEdxParam = (*iter) -> GetDEDX(particle, material,
                                            transitionEnergy);
 
         G4double dEdxDeltaRays = DeltaRayMeanEnergyTransferRate(material,
                                            particle,
                                            transitionEnergy,
                                            cutEnergy);
         dEdxParam -= dEdxDeltaRays;
 
         // Computing dE/dx from Bethe-Bloch formula at transition
         // energy
         G4double transitionChargeSquare =
               GetChargeSquareRatio(particle, material, transitionEnergy);
 
         G4double scaledTransitionEnergy = transitionEnergy * massRatio;
 
         G4double dEdxBetheBloch =
                            betheBlochModel -> ComputeDEDXPerVolume(
                                         material, genericIon,
                             scaledTransitionEnergy, cutEnergy);
         dEdxBetheBloch *= transitionChargeSquare;
 
         // Additionally, high order corrections are added
         dEdxBetheBloch +=
             corrections -> ComputeIonCorrections(particle,
                                                  material, transitionEnergy);
 
         // Computing transition factor from both dE/dx values
         dedxCacheTransitionFactor =
                      (dEdxParam - dEdxBetheBloch)/dEdxBetheBloch
                              * transitionEnergy;
      }
      else {
 
         dedxCacheParticle = particle;
         dedxCacheMaterial = material;
         dedxCacheEnergyCut = cutEnergy;
 
         dedxCacheGenIonMassRatio =
                              genericIonPDGMass / particle -> GetPDGMass();
 
         dedxCacheTransitionEnergy = 0.0;
         dedxCacheTransitionFactor = 0.0;
      }
   }
 }

Here is the call graph for this function:

Here is the caller graph for this function:

◆ UpdateRangeCache()

void G4IonParametrisedLossModel::UpdateRangeCache	(	const G4ParticleDefinition *	particle,
		const G4MaterialCutsCouple *	matCutsCouple
	)

private

Definition at line 791 of file G4IonParametrisedLossModel.cc.

                                                                 {
 
   // ############## Caching ##################################################
   // If the ion-material-cut combination is covered by any native ion data
   // parameterisation (for low energies), range vectors are computed
 
   if(particle == rangeCacheParticle &&
      matCutsCouple == rangeCacheMatCutsCouple) {
   }
   else{
      rangeCacheParticle = particle;
      rangeCacheMatCutsCouple = matCutsCouple;
 
      const G4Material* material = matCutsCouple -> GetMaterial();
      LossTableList::iterator iter = IsApplicable(particle, material);
 
      // If any table is applicable, the transition factor is computed:
      if(iter != lossTableList.end()) {
 
         // Build range-energy and energy-range vectors if they don't exist
         IonMatCouple ionMatCouple = std::make_pair(particle, matCutsCouple);
         RangeEnergyTable::iterator iterRange = r.find(ionMatCouple);
 
         if(iterRange == r.end()) BuildRangeVector(particle, matCutsCouple);
 
         rangeCacheEnergyRange = E[ionMatCouple];
         rangeCacheRangeEnergy = r[ionMatCouple];
      }
      else {
         rangeCacheEnergyRange = 0;
         rangeCacheRangeEnergy = 0;
      }
   }
 }