Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4LivermoreIonisationModel Class Reference

#include <G4LivermoreIonisationModel.hh>

Inheritance diagram for G4LivermoreIonisationModel:
Collaboration 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 InitialiseLocal (const G4ParticleDefinition *, G4VEmModel *masterModel)
 
virtual void InitialiseForMaterial (const G4ParticleDefinition *, const G4Material *)
 
virtual void InitialiseForElement (const G4ParticleDefinition *, G4int Z)
 
virtual G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
virtual G4double GetPartialCrossSection (const G4Material *, G4int level, const G4ParticleDefinition *, G4double kineticEnergy)
 
virtual G4double ComputeCrossSectionPerShell (const G4ParticleDefinition *, G4int Z, G4int shellIdx, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
virtual G4double ChargeSquareRatio (const G4Track &)
 
virtual 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 * > *)
 
virtual 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=nullptr)
 
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 Attributes

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

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 *)
 
- Static Protected Attributes inherited from G4VEmModel
static const G4double inveplus = 1.0/CLHEP::eplus
 

Detailed Description

Definition at line 58 of file G4LivermoreIonisationModel.hh.

Constructor & Destructor Documentation

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

Definition at line 77 of file G4LivermoreIonisationModel.cc.

78  :
79  G4VEmModel(nam), fParticleChange(0),
80  isInitialised(false),
81  crossSectionHandler(0), energySpectrum(0)
82 {
83  fIntrinsicLowEnergyLimit = 12.*eV;
84  fIntrinsicHighEnergyLimit = 100.0*GeV;
85 
86  verboseLevel = 0;
88 
89  transitionManager = G4AtomicTransitionManager::Instance();
90 }
G4VEmModel(const G4String &nam)
Definition: G4VEmModel.cc:68
G4ParticleChangeForLoss * fParticleChange
static constexpr double eV
Definition: G4SIunits.hh:215
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:623
static constexpr double GeV
Definition: G4SIunits.hh:217
static G4AtomicTransitionManager * Instance()

Here is the call graph for this function:

G4LivermoreIonisationModel::~G4LivermoreIonisationModel ( )
virtual

Definition at line 94 of file G4LivermoreIonisationModel.cc.

95 {
96  delete energySpectrum;
97  delete crossSectionHandler;
98 }

Member Function Documentation

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 176 of file G4LivermoreIonisationModel.cc.

182 {
183  G4int iZ = (G4int) Z;
184  if (!crossSectionHandler)
185  {
186  G4Exception("G4LivermoreIonisationModel::ComputeCrossSectionPerAtom",
187  "em1007",FatalException,
188  "The cross section handler is not correctly initialized");
189  return 0;
190  }
191 
192  //The cut is already included in the crossSectionHandler
193  G4double cs =
194  crossSectionHandler->GetCrossSectionAboveThresholdForElement(energy,
195  cutEnergy,
196  iZ);
197 
198  if (verboseLevel > 1)
199  {
200  G4cout << "G4LivermoreIonisationModel " << G4endl;
201  G4cout << "Cross section for delta emission > "
202  << cutEnergy/keV << " keV at "
203  << energy/keV << " keV and Z = " << iZ << " --> "
204  << cs/barn << " barn" << G4endl;
205  }
206  return cs;
207 }
G4double GetCrossSectionAboveThresholdForElement(G4double energy, G4double cutEnergy, G4int Z)
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4double energy(const ThreeVector &p, const G4double m)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
static constexpr double barn
Definition: G4SIunits.hh:105
static constexpr double keV
Definition: G4SIunits.hh:216

Here is the call graph for this function:

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

Reimplemented from G4VEmModel.

Definition at line 213 of file G4LivermoreIonisationModel.cc.

217 {
218  G4double sPower = 0.0;
219 
220  const G4ElementVector* theElementVector = material->GetElementVector();
221  size_t NumberOfElements = material->GetNumberOfElements() ;
222  const G4double* theAtomicNumDensityVector =
223  material->GetAtomicNumDensityVector();
224 
225  // loop for elements in the material
226  for (size_t iel=0; iel<NumberOfElements; iel++ )
227  {
228  G4int iZ = (G4int)((*theElementVector)[iel]->GetZ());
229  G4int nShells = transitionManager->NumberOfShells(iZ);
230  for (G4int n=0; n<nShells; n++)
231  {
232  G4double e = energySpectrum->AverageEnergy(iZ, 0.0,cutEnergy,
233  kineticEnergy, n);
234  G4double cs= crossSectionHandler->FindValue(iZ,kineticEnergy, n);
235  sPower += e * cs * theAtomicNumDensityVector[iel];
236  }
237  G4double esp = energySpectrum->Excitation(iZ,kineticEnergy);
238  sPower += esp * theAtomicNumDensityVector[iel];
239  }
240 
241  if (verboseLevel > 2)
242  {
243  G4cout << "G4LivermoreIonisationModel " << G4endl;
244  G4cout << "Stopping power < " << cutEnergy/keV
245  << " keV at " << kineticEnergy/keV << " keV = "
246  << sPower/(keV/mm) << " keV/mm" << G4endl;
247  }
248 
249  return sPower;
250 }
static constexpr double mm
Definition: G4SIunits.hh:115
std::vector< G4Element * > G4ElementVector
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
int G4int
Definition: G4Types.hh:78
G4double FindValue(G4int Z, G4double e) const
virtual G4double AverageEnergy(G4int Z, G4double minKineticEnergy, G4double maxKineticEnergy, G4double kineticEnergy, G4int shell=0, const G4ParticleDefinition *pd=0) const =0
G4GLOB_DLL std::ostream G4cout
const G4double * GetAtomicNumDensityVector() const
Definition: G4Material.hh:216
virtual G4double Excitation(G4int Z, G4double kineticEnergy) const =0
#define G4endl
Definition: G4ios.hh:61
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
double G4double
Definition: G4Types.hh:76
static constexpr double keV
Definition: G4SIunits.hh:216

Here is the call graph for this function:

G4int G4LivermoreIonisationModel::GetVerboseLevel ( )
inline

Definition at line 90 of file G4LivermoreIonisationModel.hh.

90 {return verboseLevel;};
void G4LivermoreIonisationModel::Initialise ( const G4ParticleDefinition particle,
const G4DataVector cuts 
)
virtual

Implements G4VEmModel.

Definition at line 102 of file G4LivermoreIonisationModel.cc.

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

Here is the call graph for this function:

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

Implements G4VEmModel.

Definition at line 254 of file G4LivermoreIonisationModel.cc.

260 {
261 
262  G4double kineticEnergy = aDynamicParticle->GetKineticEnergy();
263 
264  if (kineticEnergy <= fIntrinsicLowEnergyLimit)
265  {
268  return;
269  }
270 
271  // Select atom and shell
272  G4int Z = crossSectionHandler->SelectRandomAtom(couple, kineticEnergy);
273  G4int shellIndex = crossSectionHandler->SelectRandomShell(Z, kineticEnergy);
274  const G4AtomicShell* shell = transitionManager->Shell(Z,shellIndex);
276 
277  // Sample delta energy using energy interval for delta-electrons
278  G4double energyMax =
279  std::min(maxE,energySpectrum->MaxEnergyOfSecondaries(kineticEnergy));
280  G4double energyDelta = energySpectrum->SampleEnergy(Z, cutE, energyMax,
281  kineticEnergy, shellIndex);
282 
283  if (energyDelta == 0.) //nothing happens
284  { return; }
285 
287  G4DynamicParticle* delta = new G4DynamicParticle(electron,
288  GetAngularDistribution()->SampleDirectionForShell(aDynamicParticle, energyDelta,
289  Z, shellIndex,
290  couple->GetMaterial()),
291  energyDelta);
292 
293  fvect->push_back(delta);
294 
295  // Change kinematics of primary particle
296  G4ThreeVector direction = aDynamicParticle->GetMomentumDirection();
297  G4double totalMomentum = std::sqrt(kineticEnergy*(kineticEnergy + 2*electron_mass_c2));
298 
299  G4ThreeVector finalP = totalMomentum*direction - delta->GetMomentum();
300  finalP = finalP.unit();
301 
302  //This is the amount of energy available for fluorescence
303  G4double theEnergyDeposit = bindingEnergy;
304 
305  // fill ParticleChange
306  // changed energy and momentum of the actual particle
307  G4double finalKinEnergy = kineticEnergy - energyDelta - theEnergyDeposit;
308  if(finalKinEnergy < 0.0)
309  {
310  theEnergyDeposit += finalKinEnergy;
311  finalKinEnergy = 0.0;
312  }
313  else
314  {
316  }
317  fParticleChange->SetProposedKineticEnergy(finalKinEnergy);
318 
319  if (theEnergyDeposit < 0)
320  {
321  G4cout << "G4LivermoreIonisationModel: Negative energy deposit: "
322  << theEnergyDeposit/eV << " eV" << G4endl;
323  theEnergyDeposit = 0.0;
324  }
325 
326  //Assign local energy deposit
327  fParticleChange->ProposeLocalEnergyDeposit(theEnergyDeposit);
328 
329  if (verboseLevel > 1)
330  {
331  G4cout << "-----------------------------------------------------------" << G4endl;
332  G4cout << "Energy balance from G4LivermoreIonisation" << G4endl;
333  G4cout << "Incoming primary energy: " << kineticEnergy/keV << " keV" << G4endl;
334  G4cout << "-----------------------------------------------------------" << G4endl;
335  G4cout << "Outgoing primary energy: " << finalKinEnergy/keV << " keV" << G4endl;
336  G4cout << "Delta ray " << energyDelta/keV << " keV" << G4endl;
337  G4cout << "Fluorescence: " << (bindingEnergy-theEnergyDeposit)/keV << " keV" << G4endl;
338  G4cout << "Local energy deposit " << theEnergyDeposit/keV << " keV" << G4endl;
339  G4cout << "Total final state: " << (finalKinEnergy+energyDelta+bindingEnergy)
340  << " keV" << G4endl;
341  G4cout << "-----------------------------------------------------------" << G4endl;
342  }
343  return;
344 }
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
G4double GetKineticEnergy() const
G4int SelectRandomShell(G4int Z, G4double e) const
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:616
G4double BindingEnergy() const
int G4int
Definition: G4Types.hh:78
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
static constexpr double electron_mass_c2
G4ParticleChangeForLoss * fParticleChange
virtual G4double SampleEnergy(G4int Z, G4double minKineticEnergy, G4double maxKineticEnergy, G4double kineticEnergy, G4int shell=0, const G4ParticleDefinition *pd=0) const =0
G4GLOB_DLL std::ostream G4cout
const G4ThreeVector & GetMomentumDirection() const
void SetProposedKineticEnergy(G4double proposedKinEnergy)
static constexpr double eV
Definition: G4SIunits.hh:215
virtual G4double MaxEnergyOfSecondaries(G4double kineticEnergy, G4int Z=0, const G4ParticleDefinition *pd=0) const =0
Hep3Vector unit() const
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
static G4Electron * Electron()
Definition: G4Electron.cc:94
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
G4double bindingEnergy(G4int A, G4int Z)
static constexpr double keV
Definition: G4SIunits.hh:216
G4AtomicShell * Shell(G4int Z, size_t shellIndex) const
G4int SelectRandomAtom(const G4MaterialCutsCouple *couple, G4double e) const
G4ThreeVector GetMomentum() const
const G4Material * GetMaterial() const

Here is the call graph for this function:

void G4LivermoreIonisationModel::SetVerboseLevel ( G4int  vl)
inline

Definition at line 89 of file G4LivermoreIonisationModel.hh.

89 {verboseLevel = vl;};

Member Data Documentation

G4ParticleChangeForLoss* G4LivermoreIonisationModel::fParticleChange
protected

Definition at line 90 of file G4LivermoreIonisationModel.hh.


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