Geant4  10.02.p03
G4MicroElecInelasticModel Class Reference

#include <G4MicroElecInelasticModel.hh>

Inheritance diagram for G4MicroElecInelasticModel:
Collaboration diagram for G4MicroElecInelasticModel:

Public Member Functions

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

Protected Attributes

G4ParticleChangeForGamma * fParticleChangeForGamma
 
- Protected Attributes inherited from G4VEmModel
G4ElementDatafElementData
 
G4VParticleChange * pParticleChange
 
G4PhysicsTablexSectionTable
 
const std::vector< G4double > * theDensityFactor
 
const std::vector< G4int > * theDensityIdx
 
size_t idxTable
 

Private Types

typedef std::map< G4String, G4String, std::less< G4String > > MapFile
 
typedef std::map< G4String, G4MicroElecCrossSectionDataSet *, std::less< G4String > > MapData
 
typedef std::map< double, std::map< double, double > > TriDimensionMap
 
typedef std::map< double, std::vector< double > > VecMap
 

Private Member Functions

G4double RandomizeEjectedElectronEnergy (G4ParticleDefinition *aParticleDefinition, G4double incomingParticleEnergy, G4int shell)
 
void RandomizeEjectedElectronDirection (G4ParticleDefinition *aParticleDefinition, G4double incomingParticleEnergy, G4double outgoingParticleEnergy, G4double &cosTheta, G4double &phi)
 
G4double LogLogInterpolate (G4double e1, G4double e2, G4double e, G4double xs1, G4double xs2)
 
G4double QuadInterpolator (G4double e11, G4double e12, G4double e21, G4double e22, G4double x11, G4double x12, G4double x21, G4double x22, G4double t1, G4double t2, G4double t, G4double e)
 
G4int RandomSelect (G4double energy, const G4String &particle)
 
G4MicroElecInelasticModeloperator= (const G4MicroElecInelasticModel &right)
 
 G4MicroElecInelasticModel (const G4MicroElecInelasticModel &)
 

Private Attributes

G4VAtomDeexcitationfAtomDeexcitation
 
G4MaterialnistSi
 
std::map< G4String, G4double, std::less< G4String > > lowEnergyLimit
 
std::map< G4String, G4double, std::less< G4String > > highEnergyLimit
 
G4bool isInitialised
 
G4int verboseLevel
 
MapFile tableFile
 
MapData tableData
 
G4MicroElecSiStructure SiStructure
 
TriDimensionMap eDiffCrossSectionData [7]
 
TriDimensionMap pDiffCrossSectionData [7]
 
std::vector< double > eTdummyVec
 
std::vector< double > pTdummyVec
 
VecMap eVecm
 
VecMap pVecm
 

Additional Inherited Members

- Protected Member Functions inherited from G4VEmModel
G4ParticleChangeForLoss * GetParticleChangeForLoss ()
 
G4ParticleChangeForGamma * GetParticleChangeForGamma ()
 
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 64 of file G4MicroElecInelasticModel.hh.

Member Typedef Documentation

◆ MapData

Definition at line 112 of file G4MicroElecInelasticModel.hh.

◆ MapFile

typedef std::map<G4String,G4String,std::less<G4String> > G4MicroElecInelasticModel::MapFile
private

Definition at line 109 of file G4MicroElecInelasticModel.hh.

◆ TriDimensionMap

typedef std::map<double, std::map<double, double> > G4MicroElecInelasticModel::TriDimensionMap
private

Definition at line 139 of file G4MicroElecInelasticModel.hh.

◆ VecMap

typedef std::map<double, std::vector<double> > G4MicroElecInelasticModel::VecMap
private

Definition at line 145 of file G4MicroElecInelasticModel.hh.

Constructor & Destructor Documentation

◆ G4MicroElecInelasticModel() [1/2]

G4MicroElecInelasticModel::G4MicroElecInelasticModel ( const G4ParticleDefinition p = 0,
const G4String nam = "MicroElecInelasticModel" 
)

Definition at line 62 of file G4MicroElecInelasticModel.cc.

65 {
67 
68  verboseLevel= 0;
69  // Verbosity scale:
70  // 0 = nothing
71  // 1 = warning for energy non-conservation
72  // 2 = details of energy budget
73  // 3 = calculation of cross sections, file openings, sampling of atoms
74  // 4 = entering in methods
75 
76  if( verboseLevel>0 )
77  {
78  G4cout << "MicroElec inelastic model is constructed " << G4endl;
79  }
80 
81  //Mark this model as "applicable" for atomic deexcitation
82  SetDeexcitationFlag(true);
84 
85  // default generator
87 }
G4Material * FindOrBuildMaterial(const G4String &name, G4bool isotopes=true, G4bool warning=false)
static G4NistManager * Instance()
G4VEmModel(const G4String &nam)
Definition: G4VEmModel.cc:69
G4GLOB_DLL std::ostream G4cout
G4ParticleChangeForGamma * fParticleChangeForGamma
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:624
#define G4endl
Definition: G4ios.hh:61
void SetDeexcitationFlag(G4bool val)
Definition: G4VEmModel.hh:781
G4VAtomDeexcitation * fAtomDeexcitation
Here is the call graph for this function:

◆ ~G4MicroElecInelasticModel()

G4MicroElecInelasticModel::~G4MicroElecInelasticModel ( )
virtual

Definition at line 91 of file G4MicroElecInelasticModel.cc.

92 {
93  // Cross section
94 
95  std::map< G4String,G4MicroElecCrossSectionDataSet*,std::less<G4String> >::iterator pos;
96  for (pos = tableData.begin(); pos != tableData.end(); ++pos)
97  {
98  G4MicroElecCrossSectionDataSet* table = pos->second;
99  delete table;
100  }
101 
102  // Final state
103 
104  eVecm.clear();
105  pVecm.clear();
106 
107 }
static const G4double pos

◆ G4MicroElecInelasticModel() [2/2]

G4MicroElecInelasticModel::G4MicroElecInelasticModel ( const G4MicroElecInelasticModel )
private

Member Function Documentation

◆ CrossSectionPerVolume()

G4double G4MicroElecInelasticModel::CrossSectionPerVolume ( const G4Material material,
const G4ParticleDefinition p,
G4double  ekin,
G4double  emin,
G4double  emax 
)
virtual

Reimplemented from G4VEmModel.

Definition at line 278 of file G4MicroElecInelasticModel.cc.

283 {
284  if (verboseLevel > 3)
285  G4cout << "Calling CrossSectionPerVolume() of G4MicroElecInelasticModel" << G4endl;
286 
288 
289  /* if (
290  particleDefinition != G4Proton::ProtonDefinition()
291  &&
292  particleDefinition != G4Electron::ElectronDefinition()
293  &&
294  particleDefinition != G4GenericIon::GenericIonDefinition()
295  )
296 
297  return 0;*/
298 
299  // Calculate total cross section for model
300 
301  G4double lowLim = 0;
302  G4double highLim = 0;
303  G4double sigma=0;
304 
305  const G4String& particleName = particleDefinition->GetParticleName();
306  G4String nameLocal = particleName ;
307 
308  G4double Zeff2 = 1.0;
309  G4double Mion_c2 = particleDefinition->GetPDGMass();
310 
311  if (Mion_c2 > proton_mass_c2)
312  {
313  G4ionEffectiveCharge EffCharge ;
314  G4double Zeff = EffCharge.EffectiveCharge(particleDefinition, material,ekin);
315  Zeff2 = Zeff*Zeff;
316 
317  if (verboseLevel > 3)
318  G4cout << "Before scaling : " << G4endl
319  << "Particle : " << nameLocal << ", mass : " << Mion_c2/proton_mass_c2 << "*mp, charge " << Zeff
320  << ", Ekin (eV) = " << ekin/eV << G4endl ;
321 
322  ekin *= proton_mass_c2/Mion_c2 ;
323  nameLocal = "proton" ;
324 
325  if (verboseLevel > 3)
326  G4cout << "After scaling : " << G4endl
327  << "Particle : " << nameLocal << ", Ekin (eV) = " << ekin/eV << G4endl ;
328  }
329 
330  if (material == nistSi || material->GetBaseMaterial() == nistSi)
331  {
332 
333  std::map< G4String,G4double,std::less<G4String> >::iterator pos1;
334  pos1 = lowEnergyLimit.find(nameLocal);
335  if (pos1 != lowEnergyLimit.end())
336  {
337  lowLim = pos1->second;
338  }
339 
340  std::map< G4String,G4double,std::less<G4String> >::iterator pos2;
341  pos2 = highEnergyLimit.find(nameLocal);
342  if (pos2 != highEnergyLimit.end())
343  {
344  highLim = pos2->second;
345  }
346 
347  if (ekin >= lowLim && ekin < highLim)
348  {
349  std::map< G4String,G4MicroElecCrossSectionDataSet*,std::less<G4String> >::iterator pos;
350  pos = tableData.find(nameLocal);
351 
352  if (pos != tableData.end())
353  {
354  G4MicroElecCrossSectionDataSet* table = pos->second;
355  if (table != 0)
356  {
357  sigma = table->FindValue(ekin);
358  }
359  }
360  else
361  {
362  G4Exception("G4MicroElecInelasticModel::CrossSectionPerVolume","em0002",FatalException,"Model not applicable to particle type.");
363  }
364  }
365  else
366  {
367  if (nameLocal!="e-")
368  {
369  // G4cout << "Particle : " << nameLocal << ", Ekin (eV) = " << ekin/eV << G4endl;
370  // G4cout << "### Warning: particle energy out of bounds! ###" << G4endl;
371  }
372  }
373 
374  if (verboseLevel > 3)
375  {
376  G4cout << "---> Kinetic energy (eV)=" << ekin/eV << G4endl;
377  G4cout << " - Cross section per Si atom (cm^2)=" << sigma*Zeff2/cm2 << G4endl;
378  G4cout << " - Cross section per Si atom (cm^-1)=" << sigma*density*Zeff2/(1./cm) << G4endl;
379  }
380 
381  } // if (SiMaterial)
382  return sigma*density*Zeff2;
383 
384 
385 }
static const double cm
Definition: G4SIunits.hh:118
static const double cm2
Definition: G4SIunits.hh:119
G4double GetTotNbOfAtomsPerVolume() const
Definition: G4Material.hh:209
G4double density
Definition: TRTMaterials.hh:39
G4GLOB_DLL std::ostream G4cout
float proton_mass_c2
Definition: hepunit.py:275
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4double EffectiveCharge(const G4ParticleDefinition *p, const G4Material *material, G4double kineticEnergy)
std::map< G4String, G4double, std::less< G4String > > highEnergyLimit
static const double eV
Definition: G4SIunits.hh:212
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
const G4Material * GetBaseMaterial() const
Definition: G4Material.hh:233
std::map< G4String, G4double, std::less< G4String > > lowEnergyLimit
static const G4double pos
virtual G4double FindValue(G4double e, G4int componentId=0) const
Here is the call graph for this function:

◆ DifferentialCrossSection()

double G4MicroElecInelasticModel::DifferentialCrossSection ( G4ParticleDefinition aParticleDefinition,
G4double  k,
G4double  energyTransfer,
G4int  shell 
)

Definition at line 631 of file G4MicroElecInelasticModel.cc.

635 {
636  G4double sigma = 0.;
637 
638  if (energyTransfer >= SiStructure.Energy(LevelIndex))
639  {
640  G4double valueT1 = 0;
641  G4double valueT2 = 0;
642  G4double valueE21 = 0;
643  G4double valueE22 = 0;
644  G4double valueE12 = 0;
645  G4double valueE11 = 0;
646 
647  G4double xs11 = 0;
648  G4double xs12 = 0;
649  G4double xs21 = 0;
650  G4double xs22 = 0;
651 
652  if (particleDefinition == G4Electron::ElectronDefinition())
653  {
654  // k should be in eV and energy transfer eV also
655 
656  std::vector<double>::iterator t2 = std::upper_bound(eTdummyVec.begin(),eTdummyVec.end(), k);
657  std::vector<double>::iterator t1 = t2-1;
658  // SI : the following condition avoids situations where energyTransfer >last vector element
659  if (energyTransfer <= eVecm[(*t1)].back() && energyTransfer <= eVecm[(*t2)].back() )
660  {
661  std::vector<double>::iterator e12 = std::upper_bound(eVecm[(*t1)].begin(),eVecm[(*t1)].end(), energyTransfer);
662  std::vector<double>::iterator e11 = e12-1;
663 
664  std::vector<double>::iterator e22 = std::upper_bound(eVecm[(*t2)].begin(),eVecm[(*t2)].end(), energyTransfer);
665  std::vector<double>::iterator e21 = e22-1;
666 
667  valueT1 =*t1;
668  valueT2 =*t2;
669  valueE21 =*e21;
670  valueE22 =*e22;
671  valueE12 =*e12;
672  valueE11 =*e11;
673 
674  xs11 = eDiffCrossSectionData[LevelIndex][valueT1][valueE11];
675  xs12 = eDiffCrossSectionData[LevelIndex][valueT1][valueE12];
676  xs21 = eDiffCrossSectionData[LevelIndex][valueT2][valueE21];
677  xs22 = eDiffCrossSectionData[LevelIndex][valueT2][valueE22];
678  }
679 
680  }
681 
682  if (particleDefinition == G4Proton::ProtonDefinition())
683  {
684  // k should be in eV and energy transfer eV also
685  std::vector<double>::iterator t2 = std::upper_bound(pTdummyVec.begin(),pTdummyVec.end(), k);
686  std::vector<double>::iterator t1 = t2-1;
687  if (energyTransfer <= pVecm[(*t1)].back() && energyTransfer <= pVecm[(*t2)].back() )
688  {
689  std::vector<double>::iterator e12 = std::upper_bound(pVecm[(*t1)].begin(),pVecm[(*t1)].end(), energyTransfer);
690  std::vector<double>::iterator e11 = e12-1;
691 
692  std::vector<double>::iterator e22 = std::upper_bound(pVecm[(*t2)].begin(),pVecm[(*t2)].end(), energyTransfer);
693  std::vector<double>::iterator e21 = e22-1;
694 
695  valueT1 =*t1;
696  valueT2 =*t2;
697  valueE21 =*e21;
698  valueE22 =*e22;
699  valueE12 =*e12;
700  valueE11 =*e11;
701 
702  xs11 = pDiffCrossSectionData[LevelIndex][valueT1][valueE11];
703  xs12 = pDiffCrossSectionData[LevelIndex][valueT1][valueE12];
704  xs21 = pDiffCrossSectionData[LevelIndex][valueT2][valueE21];
705  xs22 = pDiffCrossSectionData[LevelIndex][valueT2][valueE22];
706  }
707  }
708 
709  G4double xsProduct = xs11 * xs12 * xs21 * xs22;
710  if (xsProduct != 0.)
711  {
712  sigma = QuadInterpolator( valueE11, valueE12,
713  valueE21, valueE22,
714  xs11, xs12,
715  xs21, xs22,
716  valueT1, valueT2,
717  k, energyTransfer);
718  }
719 
720  }
721 
722  return sigma;
723 }
static G4Electron * ElectronDefinition()
Definition: G4Electron.cc:89
TTree * t1
Definition: plottest35.C:26
static G4Proton * ProtonDefinition()
Definition: G4Proton.cc:88
G4double Energy(G4int level)
G4double QuadInterpolator(G4double e11, G4double e12, G4double e21, G4double e22, G4double x11, G4double x12, G4double x21, G4double x22, G4double t1, G4double t2, G4double t, G4double e)
TTree * t2
Definition: plottest35.C:36
double G4double
Definition: G4Types.hh:76
Here is the call graph for this function:
Here is the caller graph for this function:

◆ Initialise()

void G4MicroElecInelasticModel::Initialise ( const G4ParticleDefinition particle,
const G4DataVector  
)
virtual

Implements G4VEmModel.

Definition at line 111 of file G4MicroElecInelasticModel.cc.

113 {
114 
115  if (verboseLevel > 3)
116  G4cout << "Calling G4MicroElecInelasticModel::Initialise()" << G4endl;
117 
118  // Energy limits
119 
120  G4String fileElectron("microelec/sigma_inelastic_e_Si");
121  G4String fileProton("microelec/sigma_inelastic_p_Si");
122 
125 
128 
129  G4double scaleFactor = 1e-18 * cm *cm;
130 
131  char *path = getenv("G4LEDATA");
132 
133  // *** ELECTRON
134  electron = electronDef->GetParticleName();
135 
136  tableFile[electron] = fileElectron;
137 
138  lowEnergyLimit[electron] = 16.7 * eV;
139  highEnergyLimit[electron] = 100.0 * MeV;
140 
141  // Cross section
142 
144  tableE->LoadData(fileElectron);
145 
146  tableData[electron] = tableE;
147 
148  // Final state
149 
150  std::ostringstream eFullFileName;
151  eFullFileName << path << "/microelec/sigmadiff_inelastic_e_Si.dat";
152  std::ifstream eDiffCrossSection(eFullFileName.str().c_str());
153 
154  if (!eDiffCrossSection)
155  {
156  G4Exception("G4MicroElecInelasticModel::Initialise","em0003",FatalException,"Missing data file:/microelec/sigmadiff_inelastic_e_Si.dat");
157  }
158 
159  //
160 
161  // Clear the arrays for re-initialization case (MT mode)
162  // Octobre 22nd, 2014 - Melanie Raine
163 
164  eTdummyVec.clear();
165  pTdummyVec.clear();
166 
167  eVecm.clear();
168  pVecm.clear();
169 
170  eDiffCrossSectionData->clear();
171  pDiffCrossSectionData->clear();
172 
173  //
174 
175 
176  eTdummyVec.push_back(0.);
177  while(!eDiffCrossSection.eof())
178  {
179  double tDummy;
180  double eDummy;
181  eDiffCrossSection>>tDummy>>eDummy;
182  if (tDummy != eTdummyVec.back()) eTdummyVec.push_back(tDummy);
183  for (int j=0; j<6; j++)
184  {
185  eDiffCrossSection>>eDiffCrossSectionData[j][tDummy][eDummy];
186 
187  // SI - only if eof is not reached !
188  if (!eDiffCrossSection.eof()) eDiffCrossSectionData[j][tDummy][eDummy]*=scaleFactor;
189 
190  eVecm[tDummy].push_back(eDummy);
191 
192  }
193  }
194  //
195 
196  // *** PROTON
197 
198  proton = protonDef->GetParticleName();
199 
200  tableFile[proton] = fileProton;
201 
202  lowEnergyLimit[proton] = 50. * keV;
203  highEnergyLimit[proton] = 10. * GeV;
204 
205  // Cross section
206 
208  tableP->LoadData(fileProton);
209 
210  tableData[proton] = tableP;
211 
212  // Final state
213 
214  std::ostringstream pFullFileName;
215  pFullFileName << path << "/microelec/sigmadiff_inelastic_p_Si.dat";
216  std::ifstream pDiffCrossSection(pFullFileName.str().c_str());
217 
218  if (!pDiffCrossSection)
219  {
220  G4Exception("G4MicroElecInelasticModel::Initialise","em0003",FatalException,"Missing data file:/microelec/sigmadiff_inelastic_p_Si.dat");
221  }
222 
223  pTdummyVec.push_back(0.);
224  while(!pDiffCrossSection.eof())
225  {
226  double tDummy;
227  double eDummy;
228  pDiffCrossSection>>tDummy>>eDummy;
229  if (tDummy != pTdummyVec.back()) pTdummyVec.push_back(tDummy);
230  for (int j=0; j<6; j++)
231  {
232  pDiffCrossSection>>pDiffCrossSectionData[j][tDummy][eDummy];
233 
234  // SI - only if eof is not reached !
235  if (!pDiffCrossSection.eof()) pDiffCrossSectionData[j][tDummy][eDummy]*=scaleFactor;
236 
237  pVecm[tDummy].push_back(eDummy);
238  }
239  }
240 
241 
242  if (particle==electronDef)
243  {
246  }
247 
248  if (particle==protonDef)
249  {
252  }
253 
254  if( verboseLevel>0 )
255  {
256  G4cout << "MicroElec Inelastic model is initialized " << G4endl
257  << "Energy range: "
258  << LowEnergyLimit() / keV << " keV - "
259  << HighEnergyLimit() / MeV << " MeV for "
260  << particle->GetParticleName()
261  << " with mass (amu) " << particle->GetPDGMass()/proton_mass_c2
262  << " and charge " << particle->GetPDGCharge()
263  << G4endl << G4endl ;
264  }
265 
266  //
267 
269 
270  if (isInitialised) { return; }
272  isInitialised = true;
273 
274 }
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:641
static const double cm
Definition: G4SIunits.hh:118
static G4Electron * ElectronDefinition()
Definition: G4Electron.cc:89
static const double MeV
Definition: G4SIunits.hh:211
static G4LossTableManager * Instance()
static G4Proton * ProtonDefinition()
Definition: G4Proton.cc:88
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:725
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:634
float proton_mass_c2
Definition: hepunit.py:275
static const double GeV
Definition: G4SIunits.hh:214
G4ParticleChangeForGamma * fParticleChangeForGamma
virtual G4bool LoadData(const G4String &argFileName)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
std::map< G4String, G4double, std::less< G4String > > highEnergyLimit
static const double eV
Definition: G4SIunits.hh:212
#define G4endl
Definition: G4ios.hh:61
G4VAtomDeexcitation * AtomDeexcitation()
static const double keV
Definition: G4SIunits.hh:213
double G4double
Definition: G4Types.hh:76
void SetLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:732
G4VAtomDeexcitation * fAtomDeexcitation
G4double GetPDGCharge() const
std::map< G4String, G4double, std::less< G4String > > lowEnergyLimit
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition: G4VEmModel.cc:134
Here is the call graph for this function:

◆ LogLogInterpolate()

G4double G4MicroElecInelasticModel::LogLogInterpolate ( G4double  e1,
G4double  e2,
G4double  e,
G4double  xs1,
G4double  xs2 
)
private

Definition at line 727 of file G4MicroElecInelasticModel.cc.

732 {
733  G4double a = (std::log10(xs2)-std::log10(xs1)) / (std::log10(e2)-std::log10(e1));
734  G4double b = std::log10(xs2) - a*std::log10(e2);
735  G4double sigma = a*std::log10(e) + b;
736  G4double value = (std::pow(10.,sigma));
737  return value;
738 }
static const G4double e2
static const G4double e1
double G4double
Definition: G4Types.hh:76
Here is the caller graph for this function:

◆ operator=()

G4MicroElecInelasticModel& G4MicroElecInelasticModel::operator= ( const G4MicroElecInelasticModel right)
private

◆ QuadInterpolator()

G4double G4MicroElecInelasticModel::QuadInterpolator ( G4double  e11,
G4double  e12,
G4double  e21,
G4double  e22,
G4double  x11,
G4double  x12,
G4double  x21,
G4double  x22,
G4double  t1,
G4double  t2,
G4double  t,
G4double  e 
)
private

Definition at line 742 of file G4MicroElecInelasticModel.cc.

748 {
749  G4double interpolatedvalue1 = LogLogInterpolate(e11, e12, e, xs11, xs12);
750  G4double interpolatedvalue2 = LogLogInterpolate(e21, e22, e, xs21, xs22);
751  G4double value = LogLogInterpolate(t1, t2, t, interpolatedvalue1, interpolatedvalue2);
752  return value;
753 }
TTree * t1
Definition: plottest35.C:26
TTree * t2
Definition: plottest35.C:36
double G4double
Definition: G4Types.hh:76
G4double LogLogInterpolate(G4double e1, G4double e2, G4double e, G4double xs1, G4double xs2)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ RandomizeEjectedElectronDirection()

void G4MicroElecInelasticModel::RandomizeEjectedElectronDirection ( G4ParticleDefinition aParticleDefinition,
G4double  incomingParticleEnergy,
G4double  outgoingParticleEnergy,
G4double cosTheta,
G4double phi 
)
private

◆ RandomizeEjectedElectronEnergy()

G4double G4MicroElecInelasticModel::RandomizeEjectedElectronEnergy ( G4ParticleDefinition aParticleDefinition,
G4double  incomingParticleEnergy,
G4int  shell 
)
private

Definition at line 525 of file G4MicroElecInelasticModel.cc.

527 {
528  if (particleDefinition == G4Electron::ElectronDefinition())
529  {
530  G4double maximumEnergyTransfer=0.;
531  if ((k+SiStructure.Energy(shell))/2. > k) maximumEnergyTransfer=k;
532  else maximumEnergyTransfer = (k+SiStructure.Energy(shell))/2.;
533 
534  G4double crossSectionMaximum = 0.;
535 
536  G4double minEnergy = SiStructure.Energy(shell);
537  G4double maxEnergy = maximumEnergyTransfer;
538  G4int nEnergySteps = 100;
539 
540  G4double value(minEnergy);
541  G4double stpEnergy(std::pow(maxEnergy/value, 1./static_cast<G4double>(nEnergySteps-1)));
542  G4int step(nEnergySteps);
543  while (step>0)
544  {
545  step--;
546  G4double differentialCrossSection = DifferentialCrossSection(particleDefinition, k/eV, value/eV, shell);
547  if(differentialCrossSection >= crossSectionMaximum) crossSectionMaximum = differentialCrossSection;
548  value*=stpEnergy;
549  }
550 
551 
552  G4double secondaryElectronKineticEnergy=0.;
553  do
554  {
555  secondaryElectronKineticEnergy = G4UniformRand() * (maximumEnergyTransfer-SiStructure.Energy(shell));
556  } while(G4UniformRand()*crossSectionMaximum >
557  DifferentialCrossSection(particleDefinition, k/eV,(secondaryElectronKineticEnergy+SiStructure.Energy(shell))/eV,shell));
558 
559  return secondaryElectronKineticEnergy;
560 
561  }
562 
563  if (particleDefinition == G4Proton::ProtonDefinition())
564  {
565  G4double maximumEnergyTransfer = 4.* (electron_mass_c2 / proton_mass_c2) * k;
566  G4double crossSectionMaximum = 0.;
567 
568  G4double minEnergy = SiStructure.Energy(shell);
569  G4double maxEnergy = maximumEnergyTransfer;
570  G4int nEnergySteps = 100;
571 
572  G4double value(minEnergy);
573  G4double stpEnergy(std::pow(maxEnergy/value, 1./static_cast<G4double>(nEnergySteps-1)));
574  G4int step(nEnergySteps);
575  while (step>0)
576  {
577  step--;
578  G4double differentialCrossSection = DifferentialCrossSection(particleDefinition, k/eV, value/eV, shell);
579  if(differentialCrossSection >= crossSectionMaximum) crossSectionMaximum = differentialCrossSection;
580  value*=stpEnergy;
581  }
582  G4double secondaryElectronKineticEnergy = 0.;
583  do
584  {
585  secondaryElectronKineticEnergy = G4UniformRand() * (maximumEnergyTransfer-SiStructure.Energy(shell));
586 
587  } while(G4UniformRand()*crossSectionMaximum >
588  DifferentialCrossSection(particleDefinition, k/eV,(secondaryElectronKineticEnergy+SiStructure.Energy(shell))/eV,shell));
589  return secondaryElectronKineticEnergy;
590  }
591 
592  return 0;
593 }
static G4Electron * ElectronDefinition()
Definition: G4Electron.cc:89
static G4Proton * ProtonDefinition()
Definition: G4Proton.cc:88
double DifferentialCrossSection(G4ParticleDefinition *aParticleDefinition, G4double k, G4double energyTransfer, G4int shell)
int G4int
Definition: G4Types.hh:78
G4double Energy(G4int level)
#define G4UniformRand()
Definition: Randomize.hh:97
float proton_mass_c2
Definition: hepunit.py:275
float electron_mass_c2
Definition: hepunit.py:274
static const double eV
Definition: G4SIunits.hh:212
double G4double
Definition: G4Types.hh:76
Here is the call graph for this function:
Here is the caller graph for this function:

◆ RandomSelect()

G4int G4MicroElecInelasticModel::RandomSelect ( G4double  energy,
const G4String particle 
)
private

Definition at line 757 of file G4MicroElecInelasticModel.cc.

758 {
759  G4int level = 0;
760 
761  std::map< G4String,G4MicroElecCrossSectionDataSet*,std::less<G4String> >::iterator pos;
762  pos = tableData.find(particle);
763 
764  if (pos != tableData.end())
765  {
766  G4MicroElecCrossSectionDataSet* table = pos->second;
767 
768  if (table != 0)
769  {
770  G4double* valuesBuffer = new G4double[table->NumberOfComponents()];
771  const size_t n(table->NumberOfComponents());
772  size_t i(n);
773  G4double value = 0.;
774 
775  while (i>0)
776  {
777  i--;
778  valuesBuffer[i] = table->GetComponent(i)->FindValue(k);
779  value += valuesBuffer[i];
780  }
781 
782  value *= G4UniformRand();
783 
784  i = n;
785 
786  while (i > 0)
787  {
788  i--;
789 
790  if (valuesBuffer[i] > value)
791  {
792  delete[] valuesBuffer;
793  return i;
794  }
795  value -= valuesBuffer[i];
796  }
797 
798  if (valuesBuffer) delete[] valuesBuffer;
799 
800  }
801  }
802  else
803  {
804  G4Exception("G4MicroElecInelasticModel::RandomSelect","em0002",FatalException,"Model not applicable to particle type.");
805  }
806 
807  return level;
808 }
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
int G4int
Definition: G4Types.hh:78
Char_t n[5]
#define G4UniformRand()
Definition: Randomize.hh:97
virtual const G4VEMDataSet * GetComponent(G4int componentId) const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
double G4double
Definition: G4Types.hh:76
static const G4double pos
Here is the call graph for this function:
Here is the caller graph for this function:

◆ SampleSecondaries()

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

Implements G4VEmModel.

Definition at line 389 of file G4MicroElecInelasticModel.cc.

394 {
395 
396  if (verboseLevel > 3)
397  G4cout << "Calling SampleSecondaries() of G4MicroElecInelasticModel" << G4endl;
398 
399  G4double lowLim = 0;
400  G4double highLim = 0;
401 
402  G4double ekin = particle->GetKineticEnergy();
403  G4double k = ekin ;
404 
405  G4ParticleDefinition* PartDef = particle->GetDefinition();
406  const G4String& particleName = PartDef->GetParticleName();
407  G4String nameLocal2 = particleName ;
408  G4double particleMass = particle->GetDefinition()->GetPDGMass();
409 
410  if (particleMass > proton_mass_c2)
411  {
412  k *= proton_mass_c2/particleMass ;
413  PartDef = G4Proton::ProtonDefinition();
414  nameLocal2 = "proton" ;
415  }
416 
417  std::map< G4String,G4double,std::less<G4String> >::iterator pos1;
418  pos1 = lowEnergyLimit.find(nameLocal2);
419 
420  if (pos1 != lowEnergyLimit.end())
421  {
422  lowLim = pos1->second;
423  }
424 
425  std::map< G4String,G4double,std::less<G4String> >::iterator pos2;
426  pos2 = highEnergyLimit.find(nameLocal2);
427 
428  if (pos2 != highEnergyLimit.end())
429  {
430  highLim = pos2->second;
431  }
432 
433  if (k >= lowLim && k < highLim)
434  {
435  G4ParticleMomentum primaryDirection = particle->GetMomentumDirection();
436  G4double totalEnergy = ekin + particleMass;
437  G4double pSquare = ekin * (totalEnergy + particleMass);
438  G4double totalMomentum = std::sqrt(pSquare);
439 
440  G4int Shell = RandomSelect(k,nameLocal2);
441  G4double bindingEnergy = SiStructure.Energy(Shell);
442 
443  if (verboseLevel > 3)
444  {
445  G4cout << "---> Kinetic energy (eV)=" << k/eV << G4endl ;
446  G4cout << "Shell: " << Shell << ", energy: " << bindingEnergy/eV << G4endl;
447  }
448 
449  // sample deexcitation
450 
451  G4int secNumberInit = 0; // need to know at a certain point the energy of secondaries
452  G4int secNumberFinal = 0; // So I'll make the difference and then sum the energies
453 
454  G4int Z = 14;
455 
456  if(fAtomDeexcitation && Shell > 2) {
457 
459 
460  if (Shell == 4)
461  {
462  as = G4AtomicShellEnumerator(1);
463  }
464  else if (Shell == 3)
465  {
466  as = G4AtomicShellEnumerator(3);
467  }
468 
469  const G4AtomicShell* shell = fAtomDeexcitation->GetAtomicShell(Z, as);
470  secNumberInit = fvect->size();
471  fAtomDeexcitation->GenerateParticles(fvect, shell, Z, 0, 0);
472  secNumberFinal = fvect->size();
473  }
474 
475  G4double secondaryKinetic = RandomizeEjectedElectronEnergy(PartDef,k,Shell);
476 
477  if (verboseLevel > 3)
478  {
479  G4cout << "Ionisation process" << G4endl;
480  G4cout << "Shell: " << Shell << " Kin. energy (eV)=" << k/eV
481  << " Sec. energy (eV)=" << secondaryKinetic/eV << G4endl;
482  }
483 
484  G4ThreeVector deltaDirection =
485  GetAngularDistribution()->SampleDirectionForShell(particle, secondaryKinetic,
486  Z, Shell,
487  couple->GetMaterial());
488 
489  //if (particle->GetDefinition() == G4Electron::ElectronDefinition())
490  //{
491  G4double deltaTotalMomentum = std::sqrt(secondaryKinetic*(secondaryKinetic + 2.*electron_mass_c2 ));
492 
493  G4double finalPx = totalMomentum*primaryDirection.x() - deltaTotalMomentum*deltaDirection.x();
494  G4double finalPy = totalMomentum*primaryDirection.y() - deltaTotalMomentum*deltaDirection.y();
495  G4double finalPz = totalMomentum*primaryDirection.z() - deltaTotalMomentum*deltaDirection.z();
496  G4double finalMomentum = std::sqrt(finalPx*finalPx + finalPy*finalPy + finalPz*finalPz);
497  finalPx /= finalMomentum;
498  finalPy /= finalMomentum;
499  finalPz /= finalMomentum;
500 
501  G4ThreeVector direction;
502  direction.set(finalPx,finalPy,finalPz);
503 
504  fParticleChangeForGamma->ProposeMomentumDirection(direction.unit()) ;
505  //}
506  //else fParticleChangeForGamma->ProposeMomentumDirection(primaryDirection) ;
507 
508  // note that secondaryKinetic is the energy of the delta ray, not of all secondaries.
509  G4double deexSecEnergy = 0;
510  for (G4int j=secNumberInit; j < secNumberFinal; j++) {
511  deexSecEnergy = deexSecEnergy + (*fvect)[j]->GetKineticEnergy();}
512 
513  fParticleChangeForGamma->SetProposedKineticEnergy(ekin-bindingEnergy-secondaryKinetic);
514  fParticleChangeForGamma->ProposeLocalEnergyDeposit(bindingEnergy-deexSecEnergy);
515 
516  G4DynamicParticle* dp = new G4DynamicParticle (G4Electron::Electron(),deltaDirection,secondaryKinetic) ;
517  fvect->push_back(dp);
518 
519  }
520 
521 }
void set(double x, double y, double z)
const G4Material * GetMaterial() const
static G4Proton * ProtonDefinition()
Definition: G4Proton.cc:88
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:617
int G4int
Definition: G4Types.hh:78
G4double Energy(G4int level)
virtual G4ThreeVector & SampleDirectionForShell(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, G4int shellID, const G4Material *)
virtual const G4AtomicShell * GetAtomicShell(G4int Z, G4AtomicShellEnumerator shell)=0
G4double GetKineticEnergy() const
G4double RandomizeEjectedElectronEnergy(G4ParticleDefinition *aParticleDefinition, G4double incomingParticleEnergy, G4int shell)
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
Float_t Z
Hep3Vector unit() const
G4int RandomSelect(G4double energy, const G4String &particle)
float proton_mass_c2
Definition: hepunit.py:275
float electron_mass_c2
Definition: hepunit.py:274
G4ParticleChangeForGamma * fParticleChangeForGamma
double x() const
std::map< G4String, G4double, std::less< G4String > > highEnergyLimit
static const double eV
Definition: G4SIunits.hh:212
double y() const
const G4ThreeVector & GetMomentumDirection() const
void GenerateParticles(std::vector< G4DynamicParticle *> *secVect, const G4AtomicShell *, G4int Z, G4int coupleIndex)
double z() const
static G4Electron * Electron()
Definition: G4Electron.cc:94
G4ParticleDefinition * GetDefinition() const
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
G4VAtomDeexcitation * fAtomDeexcitation
std::map< G4String, G4double, std::less< G4String > > lowEnergyLimit
G4AtomicShellEnumerator
Here is the call graph for this function:

Member Data Documentation

◆ eDiffCrossSectionData

TriDimensionMap G4MicroElecInelasticModel::eDiffCrossSectionData[7]
private

Definition at line 140 of file G4MicroElecInelasticModel.hh.

◆ eTdummyVec

std::vector<double> G4MicroElecInelasticModel::eTdummyVec
private

Definition at line 142 of file G4MicroElecInelasticModel.hh.

◆ eVecm

VecMap G4MicroElecInelasticModel::eVecm
private

Definition at line 146 of file G4MicroElecInelasticModel.hh.

◆ fAtomDeexcitation

G4VAtomDeexcitation* G4MicroElecInelasticModel::fAtomDeexcitation
private

Definition at line 97 of file G4MicroElecInelasticModel.hh.

◆ fParticleChangeForGamma

G4ParticleChangeForGamma* G4MicroElecInelasticModel::fParticleChangeForGamma
protected

Definition at line 92 of file G4MicroElecInelasticModel.hh.

◆ highEnergyLimit

std::map<G4String,G4double,std::less<G4String> > G4MicroElecInelasticModel::highEnergyLimit
private

Definition at line 102 of file G4MicroElecInelasticModel.hh.

◆ isInitialised

G4bool G4MicroElecInelasticModel::isInitialised
private

Definition at line 104 of file G4MicroElecInelasticModel.hh.

◆ lowEnergyLimit

std::map<G4String,G4double,std::less<G4String> > G4MicroElecInelasticModel::lowEnergyLimit
private

Definition at line 101 of file G4MicroElecInelasticModel.hh.

◆ nistSi

G4Material* G4MicroElecInelasticModel::nistSi
private

Definition at line 99 of file G4MicroElecInelasticModel.hh.

◆ pDiffCrossSectionData

TriDimensionMap G4MicroElecInelasticModel::pDiffCrossSectionData[7]
private

Definition at line 141 of file G4MicroElecInelasticModel.hh.

◆ pTdummyVec

std::vector<double> G4MicroElecInelasticModel::pTdummyVec
private

Definition at line 143 of file G4MicroElecInelasticModel.hh.

◆ pVecm

VecMap G4MicroElecInelasticModel::pVecm
private

Definition at line 147 of file G4MicroElecInelasticModel.hh.

◆ SiStructure

G4MicroElecSiStructure G4MicroElecInelasticModel::SiStructure
private

Definition at line 117 of file G4MicroElecInelasticModel.hh.

◆ tableData

MapData G4MicroElecInelasticModel::tableData
private

Definition at line 113 of file G4MicroElecInelasticModel.hh.

◆ tableFile

MapFile G4MicroElecInelasticModel::tableFile
private

Definition at line 110 of file G4MicroElecInelasticModel.hh.

◆ verboseLevel

G4int G4MicroElecInelasticModel::verboseLevel
private

Definition at line 105 of file G4MicroElecInelasticModel.hh.


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