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

#include <G4LivermoreIonisationModel.hh>

+ Inheritance diagram for G4LivermoreIonisationModel:

Public Member Functions

 G4LivermoreIonisationModel (const G4ParticleDefinition *p=0, const G4String &processName="LowEnergyIoni")
 
virtual ~G4LivermoreIonisationModel ()
 
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)
 
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0, G4double cut=0, G4double emax=DBL_MAX)
 
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
 
virtual G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy)
 
void SetVerboseLevel (G4int vl)
 
G4int GetVerboseLevel ()
 
- Public Member Functions inherited from G4VEmModel
 G4VEmModel (const G4String &nam)
 
virtual ~G4VEmModel ()
 
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)=0
 
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin=0.0, G4double tmax=DBL_MAX)=0
 
virtual G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
 
virtual G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
virtual G4double ChargeSquareRatio (const G4Track &)
 
virtual G4double GetChargeSquareRatio (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual G4double GetParticleCharge (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual void StartTracking (G4Track *)
 
virtual void CorrectionsAlongStep (const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double &eloss, G4double &niel, G4double length)
 
virtual G4double Value (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
 
virtual G4double MinPrimaryEnergy (const G4Material *, const G4ParticleDefinition *)
 
virtual void SetupForMaterial (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual void DefineForRegion (const G4Region *)
 
void InitialiseElementSelectors (const G4ParticleDefinition *, const G4DataVector &)
 
G4double ComputeDEDX (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
 
G4double CrossSection (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeMeanFreePath (const G4ParticleDefinition *, G4double kineticEnergy, const G4Material *, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, const G4Element *, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4int SelectIsotopeNumber (const G4Element *)
 
const G4ElementSelectRandomAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementSelectRandomAtom (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
void SetParticleChange (G4VParticleChange *, G4VEmFluctuationModel *f=0)
 
void SetCrossSectionTable (G4PhysicsTable *)
 
G4PhysicsTableGetCrossSectionTable ()
 
G4VEmFluctuationModelGetModelOfFluctuations ()
 
G4VEmAngularDistributionGetAngularDistribution ()
 
void SetAngularDistribution (G4VEmAngularDistribution *)
 
G4double HighEnergyLimit () const
 
G4double LowEnergyLimit () const
 
G4double HighEnergyActivationLimit () const
 
G4double LowEnergyActivationLimit () const
 
G4double PolarAngleLimit () const
 
G4double SecondaryThreshold () const
 
G4bool LPMFlag () const
 
G4bool DeexcitationFlag () const
 
G4bool ForceBuildTableFlag () const
 
void SetHighEnergyLimit (G4double)
 
void SetLowEnergyLimit (G4double)
 
void SetActivationHighEnergyLimit (G4double)
 
void SetActivationLowEnergyLimit (G4double)
 
G4bool IsActive (G4double kinEnergy)
 
void SetPolarAngleLimit (G4double)
 
void SetSecondaryThreshold (G4double)
 
void SetLPMFlag (G4bool val)
 
void SetDeexcitationFlag (G4bool val)
 
void ForceBuildTable (G4bool val)
 
G4double MaxSecondaryKinEnergy (const G4DynamicParticle *dynParticle)
 
const G4StringGetName () const
 
void SetCurrentCouple (const G4MaterialCutsCouple *)
 
const G4ElementGetCurrentElement () const
 

Protected Attributes

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

Additional Inherited Members

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

Detailed Description

Definition at line 58 of file G4LivermoreIonisationModel.hh.

Constructor & Destructor Documentation

◆ G4LivermoreIonisationModel()

G4LivermoreIonisationModel::G4LivermoreIonisationModel ( const G4ParticleDefinition p = 0,
const G4String processName = "LowEnergyIoni" 
)

Definition at line 76 of file G4LivermoreIonisationModel.cc.

77 :
79 isInitialised(false),
80 crossSectionHandler(0), energySpectrum(0)
81{
82 fIntrinsicLowEnergyLimit = 10.0*eV;
83 fIntrinsicHighEnergyLimit = 100.0*GeV;
84
85 verboseLevel = 0;
86
87 transitionManager = G4AtomicTransitionManager::Instance();
88}
static G4AtomicTransitionManager * Instance()
G4ParticleChangeForLoss * fParticleChange

◆ ~G4LivermoreIonisationModel()

G4LivermoreIonisationModel::~G4LivermoreIonisationModel ( )
virtual

Definition at line 92 of file G4LivermoreIonisationModel.cc.

93{
94 delete energySpectrum;
95 delete crossSectionHandler;
96}

Member Function Documentation

◆ ComputeCrossSectionPerAtom()

G4double G4LivermoreIonisationModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition ,
G4double  kinEnergy,
G4double  Z,
G4double  A = 0,
G4double  cut = 0,
G4double  emax = DBL_MAX 
)
virtual

Reimplemented from G4VEmModel.

Definition at line 170 of file G4LivermoreIonisationModel.cc.

175{
176 G4int iZ = (G4int) Z;
177 if (!crossSectionHandler)
178 {
179 G4Exception("G4LivermoreIonisationModel::ComputeCrossSectionPerAtom",
180 "em1007",FatalException,"The cross section handler is not correctly initialized");
181 return 0;
182 }
183
184 //The cut is already included in the crossSectionHandler
185 G4double cs =
186 crossSectionHandler->GetCrossSectionAboveThresholdForElement(energy,cutEnergy,iZ);
187
188 if (verboseLevel > 1)
189 {
190 G4cout << "G4LivermoreIonisationModel " << G4endl;
191 G4cout << "Cross section for delta emission > " << cutEnergy/keV << " keV at " <<
192 energy/keV << " keV and Z = " << iZ << " --> " << cs/barn << " barn" << G4endl;
193 }
194 return cs;
195}
@ FatalException
double G4double
Definition: G4Types.hh:64
int G4int
Definition: G4Types.hh:66
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
G4double GetCrossSectionAboveThresholdForElement(G4double energy, G4double cutEnergy, G4int Z)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41

◆ ComputeDEDXPerVolume()

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

Reimplemented from G4VEmModel.

Definition at line 200 of file G4LivermoreIonisationModel.cc.

204{
205 G4double sPower = 0.0;
206
207 const G4ElementVector* theElementVector = material->GetElementVector();
208 size_t NumberOfElements = material->GetNumberOfElements() ;
209 const G4double* theAtomicNumDensityVector =
210 material->GetAtomicNumDensityVector();
211
212 // loop for elements in the material
213 for (size_t iel=0; iel<NumberOfElements; iel++ )
214 {
215 G4int iZ = (G4int)((*theElementVector)[iel]->GetZ());
216 G4int nShells = transitionManager->NumberOfShells(iZ);
217 for (G4int n=0; n<nShells; n++)
218 {
219 G4double e = energySpectrum->AverageEnergy(iZ, 0.0,cutEnergy,
220 kineticEnergy, n);
221 G4double cs= crossSectionHandler->FindValue(iZ,kineticEnergy, n);
222 sPower += e * cs * theAtomicNumDensityVector[iel];
223 }
224 G4double esp = energySpectrum->Excitation(iZ,kineticEnergy);
225 sPower += esp * theAtomicNumDensityVector[iel];
226 }
227
228 if (verboseLevel > 2)
229 {
230 G4cout << "G4LivermoreIonisationModel " << G4endl;
231 G4cout << "Stopping power < " << cutEnergy/keV << " keV at " <<
232 kineticEnergy/keV << " keV = " << sPower/(keV/mm) << " keV/mm" << G4endl;
233 }
234
235 return sPower;
236}
std::vector< G4Element * > G4ElementVector
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:189
size_t GetNumberOfElements() const
Definition: G4Material.hh:185
const G4double * GetAtomicNumDensityVector() const
Definition: G4Material.hh:215
G4double FindValue(G4int Z, G4double e) const
virtual G4double Excitation(G4int Z, G4double kineticEnergy) const =0
virtual G4double AverageEnergy(G4int Z, G4double minKineticEnergy, G4double maxKineticEnergy, G4double kineticEnergy, G4int shell=0, const G4ParticleDefinition *pd=0) const =0

◆ GetVerboseLevel()

G4int G4LivermoreIonisationModel::GetVerboseLevel ( )
inline

Definition at line 90 of file G4LivermoreIonisationModel.hh.

90{return verboseLevel;};

◆ Initialise()

void G4LivermoreIonisationModel::Initialise ( const G4ParticleDefinition particle,
const G4DataVector cuts 
)
virtual

Implements G4VEmModel.

Definition at line 100 of file G4LivermoreIonisationModel.cc.

102{
103 //Check that the Livermore Ionisation is NOT attached to e+
104 if (particle != G4Electron::Electron())
105 {
106 G4Exception("G4LivermoreIonisationModel::Initialise",
107 "em0002",FatalException,"Livermore Ionisation Model is applicable only to electrons");
108 }
109
110 //Read energy spectrum
111 if (energySpectrum)
112 {
113 delete energySpectrum;
114 energySpectrum = 0;
115 }
116 energySpectrum = new G4eIonisationSpectrum();
117 if (verboseLevel > 3)
118 G4cout << "G4VEnergySpectrum is initialized" << G4endl;
119
120 //Initialize cross section handler
121 if (crossSectionHandler)
122 {
123 delete crossSectionHandler;
124 crossSectionHandler = 0;
125 }
126
127 const size_t nbins = 20;
128 G4double emin = LowEnergyLimit();
129 G4double emax = HighEnergyLimit();
130 G4int ndec = G4int(std::log10(emax/emin) + 0.5);
131 if(ndec <= 0) { ndec = 1; }
132
133 G4VDataSetAlgorithm* interpolation = new G4SemiLogInterpolation();
134 crossSectionHandler = new G4eIonisationCrossSectionHandler(energySpectrum,interpolation,
135 emin,emax,nbins*ndec);
136 crossSectionHandler->Clear();
137 crossSectionHandler->LoadShellData("ioni/ion-ss-cs-");
138 //This is used to retrieve cross section values later on
139 G4VEMDataSet* emdata =
140 crossSectionHandler->BuildMeanFreePathForMaterials(&cuts);
141 //The method BuildMeanFreePathForMaterials() is required here only to force
142 //the building of an internal table: the output pointer can be deleted
143 delete emdata;
144
145 if (verboseLevel > 0)
146 {
147 G4cout << "Livermore Ionisation model is initialized " << G4endl
148 << "Energy range: "
149 << LowEnergyLimit() / keV << " keV - "
150 << HighEnergyLimit() / GeV << " GeV"
151 << G4endl;
152 }
153
154 if (verboseLevel > 3)
155 {
156 G4cout << "Cross section data: " << G4endl;
157 crossSectionHandler->PrintData();
158 G4cout << "Parameters: " << G4endl;
159 energySpectrum->PrintData();
160 }
161
162 if(isInitialised) { return; }
164 isInitialised = true;
165}
static G4Electron * Electron()
Definition: G4Electron.cc:94
G4VEMDataSet * BuildMeanFreePathForMaterials(const G4DataVector *energyCuts=0)
void LoadShellData(const G4String &dataFile)
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:529
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:522
G4ParticleChangeForLoss * GetParticleChangeForLoss()
Definition: G4VEmModel.cc:95
virtual void PrintData() const =0

◆ SampleSecondaries()

void G4LivermoreIonisationModel::SampleSecondaries ( std::vector< G4DynamicParticle * > *  fvect,
const G4MaterialCutsCouple couple,
const G4DynamicParticle aDynamicParticle,
G4double  tmin,
G4double  maxEnergy 
)
virtual

Implements G4VEmModel.

Definition at line 240 of file G4LivermoreIonisationModel.cc.

245{
246
247 G4double kineticEnergy = aDynamicParticle->GetKineticEnergy();
248
249 if (kineticEnergy <= fIntrinsicLowEnergyLimit)
250 {
253 return;
254 }
255
256 // Select atom and shell
257 G4int Z = crossSectionHandler->SelectRandomAtom(couple, kineticEnergy);
258 G4int shellIndex = crossSectionHandler->SelectRandomShell(Z, kineticEnergy);
259 const G4AtomicShell* shell = transitionManager->Shell(Z,shellIndex);
261
262 // Sample delta energy using energy interval for delta-electrons
263 G4double energyMax =
264 std::min(maxE,energySpectrum->MaxEnergyOfSecondaries(kineticEnergy));
265 G4double energyDelta = energySpectrum->SampleEnergy(Z, cutE, energyMax,
266 kineticEnergy, shellIndex);
267
268 if (energyDelta == 0.) //nothing happens
269 { return; }
270
271 // Transform to shell potential
272 G4double deltaKinE = energyDelta + 2.0*bindingEnergy;
273 G4double primaryKinE = kineticEnergy + 2.0*bindingEnergy;
274
275 // sampling of scattering angle neglecting atomic motion
276 G4double deltaMom = std::sqrt(deltaKinE*(deltaKinE + 2.0*electron_mass_c2));
277 G4double primaryMom = std::sqrt(primaryKinE*(primaryKinE + 2.0*electron_mass_c2));
278
279 G4double cost = deltaKinE * (primaryKinE + 2.0*electron_mass_c2)
280 / (deltaMom * primaryMom);
281 if (cost > 1.) { cost = 1.; }
282 G4double sint = std::sqrt((1. - cost)*(1. + cost));
283 G4double phi = twopi * G4UniformRand();
284 G4double dirx = sint * std::cos(phi);
285 G4double diry = sint * std::sin(phi);
286 G4double dirz = cost;
287
288 // Rotate to incident electron direction
289 G4ThreeVector primaryDirection = aDynamicParticle->GetMomentumDirection();
290 G4ThreeVector deltaDir(dirx,diry,dirz);
291 deltaDir.rotateUz(primaryDirection);
292 //Updated components
293 dirx = deltaDir.x();
294 diry = deltaDir.y();
295 dirz = deltaDir.z();
296
297 // Take into account atomic motion del is relative momentum of the motion
298 // kinetic energy of the motion == bindingEnergy in V.Ivanchenko model
299 cost = 2.0*G4UniformRand() - 1.0;
300 sint = std::sqrt(1. - cost*cost);
301 phi = twopi * G4UniformRand();
302 G4double del = std::sqrt(bindingEnergy *(bindingEnergy + 2.0*electron_mass_c2))
303 / deltaMom;
304 dirx += del* sint * std::cos(phi);
305 diry += del* sint * std::sin(phi);
306 dirz += del* cost;
307
308 // Find out new primary electron direction
309 G4double finalPx = primaryMom*primaryDirection.x() - deltaMom*dirx;
310 G4double finalPy = primaryMom*primaryDirection.y() - deltaMom*diry;
311 G4double finalPz = primaryMom*primaryDirection.z() - deltaMom*dirz;
312
313 //Ok, ready to create the delta ray
314 G4DynamicParticle* theDeltaRay = new G4DynamicParticle();
315 theDeltaRay->SetKineticEnergy(energyDelta);
316 G4double norm = 1.0/std::sqrt(dirx*dirx + diry*diry + dirz*dirz);
317 dirx *= norm;
318 diry *= norm;
319 dirz *= norm;
320 theDeltaRay->SetMomentumDirection(dirx, diry, dirz);
321 theDeltaRay->SetDefinition(G4Electron::Electron());
322 fvect->push_back(theDeltaRay);
323
324 //This is the amount of energy available for fluorescence
325 G4double theEnergyDeposit = bindingEnergy;
326
327 // fill ParticleChange
328 // changed energy and momentum of the actual particle
329 G4double finalKinEnergy = kineticEnergy - energyDelta - theEnergyDeposit;
330 if(finalKinEnergy < 0.0)
331 {
332 theEnergyDeposit += finalKinEnergy;
333 finalKinEnergy = 0.0;
334 }
335 else
336 {
337 G4double normLocal = 1.0/std::sqrt(finalPx*finalPx+finalPy*finalPy+finalPz*finalPz);
338 finalPx *= normLocal;
339 finalPy *= normLocal;
340 finalPz *= normLocal;
341 fParticleChange->ProposeMomentumDirection(finalPx, finalPy, finalPz);
342 }
344
345 if (theEnergyDeposit < 0)
346 {
347 G4cout << "G4LivermoreIonisationModel: Negative energy deposit: "
348 << theEnergyDeposit/eV << " eV" << G4endl;
349 theEnergyDeposit = 0.0;
350 }
351
352 //Assign local energy deposit
354
355 if (verboseLevel > 1)
356 {
357 G4cout << "-----------------------------------------------------------" << G4endl;
358 G4cout << "Energy balance from G4LivermoreIonisation" << G4endl;
359 G4cout << "Incoming primary energy: " << kineticEnergy/keV << " keV" << G4endl;
360 G4cout << "-----------------------------------------------------------" << G4endl;
361 G4cout << "Outgoing primary energy: " << finalKinEnergy/keV << " keV" << G4endl;
362 G4cout << "Delta ray " << energyDelta/keV << " keV" << G4endl;
363 G4cout << "Fluorescence: " << (bindingEnergy-theEnergyDeposit)/keV << " keV" << G4endl;
364 G4cout << "Local energy deposit " << theEnergyDeposit/keV << " keV" << G4endl;
365 G4cout << "Total final state: " << (finalKinEnergy+energyDelta+bindingEnergy)
366 << " keV" << G4endl;
367 G4cout << "-----------------------------------------------------------" << G4endl;
368 }
369 return;
370}
#define G4UniformRand()
Definition: Randomize.hh:53
double z() const
double x() const
double y() const
G4double BindingEnergy() const
G4AtomicShell * Shell(G4int Z, size_t shellIndex) const
void SetMomentumDirection(const G4ThreeVector &aDirection)
const G4ThreeVector & GetMomentumDirection() const
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4double GetKineticEnergy() const
void SetKineticEnergy(G4double aEnergy)
void SetProposedKineticEnergy(G4double proposedKinEnergy)
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
G4int SelectRandomAtom(const G4MaterialCutsCouple *couple, G4double e) const
G4int SelectRandomShell(G4int Z, G4double e) const
virtual G4double SampleEnergy(G4int Z, G4double minKineticEnergy, G4double maxKineticEnergy, G4double kineticEnergy, G4int shell=0, const G4ParticleDefinition *pd=0) const =0
virtual G4double MaxEnergyOfSecondaries(G4double kineticEnergy, G4int Z=0, const G4ParticleDefinition *pd=0) const =0
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
G4double bindingEnergy(G4int A, G4int Z)

◆ SetVerboseLevel()

void G4LivermoreIonisationModel::SetVerboseLevel ( G4int  vl)
inline

Definition at line 89 of file G4LivermoreIonisationModel.hh.

89{verboseLevel = vl;};

Member Data Documentation

◆ fParticleChange

G4ParticleChangeForLoss* G4LivermoreIonisationModel::fParticleChange
protected

Definition at line 94 of file G4LivermoreIonisationModel.hh.

Referenced by Initialise(), and SampleSecondaries().


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