Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4DNARuddIonisationExtendedModel Class Reference

#include <G4DNARuddIonisationExtendedModel.hh>

Inheritance diagram for G4DNARuddIonisationExtendedModel:
Collaboration diagram for G4DNARuddIonisationExtendedModel:

Public Member Functions

 G4DNARuddIonisationExtendedModel (const G4ParticleDefinition *p=0, const G4String &nam="DNARuddIonisationExtendedModel")
 
virtual ~G4DNARuddIonisationExtendedModel ()
 
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)
 
void SelectStationary (G4bool input)
 
- 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 level, const G4ParticleDefinition *, G4double kineticEnergy)
 
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 * > *)
 
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 SetFluctuationFlag (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

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

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 46 of file G4DNARuddIonisationExtendedModel.hh.

Constructor & Destructor Documentation

G4DNARuddIonisationExtendedModel::G4DNARuddIonisationExtendedModel ( const G4ParticleDefinition p = 0,
const G4String nam = "DNARuddIonisationExtendedModel" 
)

Definition at line 53 of file G4DNARuddIonisationExtendedModel.cc.

55  :G4VEmModel(nam),isInitialised(false)
56 {
57  // nistwater = G4NistManager::Instance()->FindOrBuildMaterial("G4_WATER");
58  fpWaterDensity = 0;
59 
60  slaterEffectiveCharge[0]=0.;
61  slaterEffectiveCharge[1]=0.;
62  slaterEffectiveCharge[2]=0.;
63  sCoefficient[0]=0.;
64  sCoefficient[1]=0.;
65  sCoefficient[2]=0.;
66 
67  lowEnergyLimitForA[1] = 0 * eV;
68  lowEnergyLimitForA[2] = 0 * eV;
69  lowEnergyLimitForA[3] = 0 * eV;
70  lowEnergyLimitOfModelForA[1] = 100 * eV;
71  lowEnergyLimitOfModelForA[4] = 1 * keV;
72  lowEnergyLimitOfModelForA[5] = 0.5 * MeV; // For A = 3 or above, limit is MeV/uma
73  killBelowEnergyForA[1] = lowEnergyLimitOfModelForA[1];
74  killBelowEnergyForA[4] = lowEnergyLimitOfModelForA[4];
75  killBelowEnergyForA[5] = lowEnergyLimitOfModelForA[5];
76 
77  verboseLevel= 0;
78  // Verbosity scale:
79  // 0 = nothing
80  // 1 = warning for energy non-conservation
81  // 2 = details of energy budget
82  // 3 = calculation of cross sections, file openings, sampling of atoms
83  // 4 = entering in methods
84 
85  if( verboseLevel>0 )
86  {
87  G4cout << "Rudd ionisation model is constructed " << G4endl;
88  }
89 
90  // Define default angular generator
92 
93  // Mark this model as "applicable" for atomic deexcitation
94  SetDeexcitationFlag(true);
95  fAtomDeexcitation = 0;
97 
98  // Selection of stationary mode
99 
100  statCode = false;
101 }
G4VEmModel(const G4String &nam)
Definition: G4VEmModel.cc:68
G4GLOB_DLL std::ostream G4cout
static constexpr double eV
Definition: G4SIunits.hh:215
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:626
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
void SetDeexcitationFlag(G4bool val)
Definition: G4VEmModel.hh:788
static constexpr double keV
Definition: G4SIunits.hh:216

Here is the call graph for this function:

G4DNARuddIonisationExtendedModel::~G4DNARuddIonisationExtendedModel ( )
virtual

Definition at line 105 of file G4DNARuddIonisationExtendedModel.cc.

106 {
107  // Cross section
108 
109  std::map< G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator pos;
110  for (pos = tableData.begin(); pos != tableData.end(); ++pos)
111  {
112  G4DNACrossSectionDataSet* table = pos->second;
113  delete table;
114  }
115 
116  // The following removal is forbidden G4VEnergyLossModel takes care of deletion
117  // however coverity will signal this as an error
118  // if (fAtomDeexcitation) {delete fAtomDeexcitation;}
119 
120 }
static const G4double pos

Member Function Documentation

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

Reimplemented from G4VEmModel.

Definition at line 540 of file G4DNARuddIonisationExtendedModel.cc.

545 {
546  //SEB: particleDefinition->GetParticleName() is for eg. Fe56
547  // particleDefinition->GetPDGMass() is correct
548  // particleDefinition->GetAtomicNumber() is correct
549 
550  if (verboseLevel > 3)
551  G4cout << "Calling CrossSectionPerVolume() of G4DNARuddIonisationExtendedModel" << G4endl;
552 
553  // Calculate total cross section for model
554 
557 
558  if (
559  particleDefinition != G4Proton::ProtonDefinition()
560  &&
561  particleDefinition != instance->GetIon("hydrogen")
562  &&
563  particleDefinition != instance->GetIon("alpha++")
564  &&
565  particleDefinition != instance->GetIon("alpha+")
566  &&
567  particleDefinition != instance->GetIon("helium")
568  &&
569  //SEB
570  //particleDefinition != instance->GetIon("carbon")
571  //&&
572  //particleDefinition != instance->GetIon("nitrogen")
573  //&&
574  //particleDefinition != instance->GetIon("oxygen")
575  //&&
576  //particleDefinition != instance->GetIon("iron")
577  particleDefinition != G4IonTable::GetIonTable()->GetIon(3,7)
578  &&
579  particleDefinition != G4IonTable::GetIonTable()->GetIon(4,9)
580  &&
581  particleDefinition != G4IonTable::GetIonTable()->GetIon(5,11)
582  &&
583  particleDefinition != G4IonTable::GetIonTable()->GetIon(6,12)
584  &&
585  particleDefinition != G4IonTable::GetIonTable()->GetIon(7,14)
586  &&
587  particleDefinition != G4IonTable::GetIonTable()->GetIon(8,16)
588  &&
589  particleDefinition != G4IonTable::GetIonTable()->GetIon(14,28)
590  &&
591  particleDefinition != G4IonTable::GetIonTable()->GetIon(26,56)
592  //
593  )
594 
595  return 0;
596 
597  G4double lowLim = 0;
598 
599  if ( particleDefinition == G4Proton::ProtonDefinition()
600  || particleDefinition == instance->GetIon("hydrogen")
601  )
602 
603  lowLim = lowEnergyLimitOfModelForA[1];
604 
605  else if ( particleDefinition == instance->GetIon("alpha++")
606  || particleDefinition == instance->GetIon("alpha+")
607  || particleDefinition == instance->GetIon("helium")
608  )
609 
610  lowLim = lowEnergyLimitOfModelForA[4];
611 
612  else lowLim = lowEnergyLimitOfModelForA[5];
613 
614  G4double highLim = 0;
615  G4double sigma=0;
616 
617 
618  G4double waterDensity = (*fpWaterDensity)[material->GetIndex()];
619 
620  if(waterDensity!= 0.0)
621 // if (material == nistwater || material->GetBaseMaterial() == nistwater)
622  {
623  const G4String& particleName = particleDefinition->GetParticleName();
624 
625  std::map< G4String,G4double,std::less<G4String> >::iterator pos2;
626  pos2 = highEnergyLimit.find(particleName);
627 
628  if (pos2 != highEnergyLimit.end())
629  {
630  highLim = pos2->second;
631  }
632 
633  if (k <= highLim)
634  {
635 
636  //SI : XS must not be zero otherwise sampling of secondaries method ignored
637 
638  if (k < lowLim) k = lowLim;
639 
640  //
641 
642  std::map< G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator pos;
643  pos = tableData.find(particleName);
644 
645  if (pos != tableData.end())
646  {
647  G4DNACrossSectionDataSet* table = pos->second;
648  if (table != 0)
649  {
650  sigma = table->FindValue(k);
651  }
652  }
653  else
654  {
655  G4Exception("G4DNARuddIonisationExtendedModel::CrossSectionPerVolume","em0002",
656  FatalException,"Model not applicable to particle type.");
657  }
658 
659  } // if (k >= lowLim && k < highLim)
660 
661  if (verboseLevel > 2)
662  {
663  G4cout << "__________________________________" << G4endl;
664  G4cout << "G4DNARuddIonisationExtendedModel - XS INFO START" << G4endl;
665  G4cout << "Kinetic energy(eV)=" << k/eV << " particle : " << particleDefinition->GetParticleName() << G4endl;
666  G4cout << "Cross section per water molecule (cm^2)=" << sigma/cm/cm << G4endl;
667  G4cout << "Cross section per water molecule (cm^-1)=" << sigma*waterDensity/(1./cm) << G4endl;
668  //G4cout << " - Cross section per water molecule (cm^-1)=" << sigma*material->GetAtomicNumDensityVector()[1]/(1./cm) << G4endl;
669  G4cout << "G4DNARuddIonisationExtendedModel - XS INFO END" << G4endl;
670 
671  }
672 
673  } // if (waterMaterial)
674 
675  return sigma*waterDensity;
676 // return sigma*material->GetAtomicNumDensityVector()[1];
677 
678 }
size_t GetIndex() const
Definition: G4Material.hh:262
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:503
static G4Proton * ProtonDefinition()
Definition: G4Proton.cc:88
G4GLOB_DLL std::ostream G4cout
static constexpr double cm
Definition: G4SIunits.hh:119
static constexpr double eV
Definition: G4SIunits.hh:215
static G4DNAGenericIonsManager * Instance(void)
virtual G4double FindValue(G4double e, G4int componentId=0) const
static G4IonTable * GetIonTable()
Definition: G4IonTable.hh:78
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define G4endl
Definition: G4ios.hh:61
static MCTruthManager * instance
double G4double
Definition: G4Types.hh:76
static const G4double pos
G4ParticleDefinition * GetIon(const G4String &name)

Here is the call graph for this function:

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

Implements G4VEmModel.

Definition at line 124 of file G4DNARuddIonisationExtendedModel.cc.

126 {
127  if (verboseLevel > 3)
128  G4cout << "Calling G4DNARuddIonisationExtendedModel::Initialise()" << G4endl;
129 
130  // Energy limits
131 
132  G4String fileProton("dna/sigma_ionisation_p_rudd");
133  G4String fileHydrogen("dna/sigma_ionisation_h_rudd");
134  G4String fileAlphaPlusPlus("dna/sigma_ionisation_alphaplusplus_rudd");
135  G4String fileAlphaPlus("dna/sigma_ionisation_alphaplus_rudd");
136  G4String fileHelium("dna/sigma_ionisation_he_rudd");
137  G4String fileLithium("dna/sigma_ionisation_li_rudd");
138  G4String fileBeryllium("dna/sigma_ionisation_be_rudd");
139  G4String fileBoron("dna/sigma_ionisation_b_rudd");
140  G4String fileCarbon("dna/sigma_ionisation_c_rudd");
141  G4String fileNitrogen("dna/sigma_ionisation_n_rudd");
142  G4String fileOxygen("dna/sigma_ionisation_o_rudd");
143  G4String fileSilicon("dna/sigma_ionisation_si_rudd");
144  G4String fileIron("dna/sigma_ionisation_fe_rudd");
145 
149  G4ParticleDefinition* hydrogenDef = instance->GetIon("hydrogen");
150  G4ParticleDefinition* alphaPlusPlusDef = instance->GetIon("alpha++");
151  G4ParticleDefinition* alphaPlusDef = instance->GetIon("alpha+");
152  G4ParticleDefinition* heliumDef = instance->GetIon("helium");
153 
154  //SEB
155  //G4ParticleDefinition* carbonDef = instance->GetIon("carbon");
156  //G4ParticleDefinition* nitrogenDef = instance->GetIon("nitrogen");
157  //G4ParticleDefinition* oxygenDef = instance->GetIon("oxygen");
158  //G4ParticleDefinition* siliconDef = instance->GetIon("silicon");
159  //G4ParticleDefinition* ironDef = instance->GetIon("iron");
161  G4ParticleDefinition* berylliumDef = G4IonTable::GetIonTable()->GetIon(4,9);
164  G4ParticleDefinition* nitrogenDef = G4IonTable::GetIonTable()->GetIon(7,14);
166  G4ParticleDefinition* siliconDef = G4IonTable::GetIonTable()->GetIon(14,28);
168  //
169 
171  G4String hydrogen;
172  G4String alphaPlusPlus;
173  G4String alphaPlus;
174  G4String helium;
175  G4String lithium;
176  G4String beryllium;
177  G4String boron;
178  G4String carbon;
179  G4String nitrogen;
180  G4String oxygen;
181  G4String silicon;
182  G4String iron;
183 
184  G4double scaleFactor = 1 * m*m;
185 
186  // LIMITS AND DATA
187 
188  // **********************************************************************************************
189 
190  proton = protonDef->GetParticleName();
191  tableFile[proton] = fileProton;
192  lowEnergyLimit[proton] = lowEnergyLimitForA[1];
193  highEnergyLimit[proton] = 500. * keV;
194 
195  // Cross section
196 
198  eV,
199  scaleFactor );
200  tableProton->LoadData(fileProton);
201  tableData[proton] = tableProton;
202 
203  // **********************************************************************************************
204 
205  hydrogen = hydrogenDef->GetParticleName();
206  tableFile[hydrogen] = fileHydrogen;
207 
208  lowEnergyLimit[hydrogen] = lowEnergyLimitForA[1];
209  highEnergyLimit[hydrogen] = 100. * MeV;
210 
211  // Cross section
212 
214  eV,
215  scaleFactor );
216  tableHydrogen->LoadData(fileHydrogen);
217 
218  tableData[hydrogen] = tableHydrogen;
219 
220  // **********************************************************************************************
221 
222  alphaPlusPlus = alphaPlusPlusDef->GetParticleName();
223  tableFile[alphaPlusPlus] = fileAlphaPlusPlus;
224 
225  lowEnergyLimit[alphaPlusPlus] = lowEnergyLimitForA[4];
226  highEnergyLimit[alphaPlusPlus] = 400. * MeV;
227 
228  // Cross section
229 
231  eV,
232  scaleFactor );
233  tableAlphaPlusPlus->LoadData(fileAlphaPlusPlus);
234 
235  tableData[alphaPlusPlus] = tableAlphaPlusPlus;
236 
237  // **********************************************************************************************
238 
239  alphaPlus = alphaPlusDef->GetParticleName();
240  tableFile[alphaPlus] = fileAlphaPlus;
241 
242  lowEnergyLimit[alphaPlus] = lowEnergyLimitForA[4];
243  highEnergyLimit[alphaPlus] = 400. * MeV;
244 
245  // Cross section
246 
248  eV,
249  scaleFactor );
250  tableAlphaPlus->LoadData(fileAlphaPlus);
251  tableData[alphaPlus] = tableAlphaPlus;
252 
253  // **********************************************************************************************
254 
255  helium = heliumDef->GetParticleName();
256  tableFile[helium] = fileHelium;
257 
258  lowEnergyLimit[helium] = lowEnergyLimitForA[4];
259  highEnergyLimit[helium] = 400. * MeV;
260 
261  // Cross section
262 
264  eV,
265  scaleFactor );
266  tableHelium->LoadData(fileHelium);
267  tableData[helium] = tableHelium;
268 
269  // **********************************************************************************************
270 
271  lithium = lithiumDef->GetParticleName();
272  tableFile[lithium] = fileLithium;
273 
274  //SEB
275  //lowEnergyLimit[carbon] = lowEnergyLimitForA[5] * particle->GetAtomicMass();
276  //highEnergyLimit[carbon] = 1e6* particle->GetAtomicMass() * MeV;
277  lowEnergyLimit[lithium] = 0.5*7*MeV;
278  highEnergyLimit[lithium] = 1e6*7*MeV;
279  //
280 
281  // Cross section
282 
284  eV,
285  scaleFactor );
286  tableLithium->LoadData(fileLithium);
287  tableData[lithium] = tableLithium;
288 
289  // **********************************************************************************************
290 
291  beryllium = berylliumDef->GetParticleName();
292  tableFile[beryllium] = fileBeryllium;
293 
294  //SEB
295  //lowEnergyLimit[carbon] = lowEnergyLimitForA[5] * particle->GetAtomicMass();
296  //highEnergyLimit[carbon] = 1e6* particle->GetAtomicMass() * MeV;
297  lowEnergyLimit[beryllium] = 0.5*9*MeV;
298  highEnergyLimit[beryllium] = 1e6*9*MeV;
299  //
300 
301  // Cross section
302 
304  eV,
305  scaleFactor );
306  tableBeryllium->LoadData(fileBeryllium);
307  tableData[beryllium] = tableBeryllium;
308 
309  // **********************************************************************************************
310 
311  boron = boronDef->GetParticleName();
312  tableFile[boron] = fileBoron;
313 
314  //SEB
315  //lowEnergyLimit[carbon] = lowEnergyLimitForA[5] * particle->GetAtomicMass();
316  //highEnergyLimit[carbon] = 1e6* particle->GetAtomicMass() * MeV;
317  lowEnergyLimit[boron] = 0.5*11*MeV;
318  highEnergyLimit[boron] = 1e6*11*MeV;
319  //
320 
321  // Cross section
322 
324  eV,
325  scaleFactor );
326  tableBoron->LoadData(fileBoron);
327  tableData[boron] = tableBoron;
328 
329  // **********************************************************************************************
330 
331  carbon = carbonDef->GetParticleName();
332  tableFile[carbon] = fileCarbon;
333 
334  //SEB
335  //lowEnergyLimit[carbon] = lowEnergyLimitForA[5] * particle->GetAtomicMass();
336  //highEnergyLimit[carbon] = 1e6* particle->GetAtomicMass() * MeV;
337  lowEnergyLimit[carbon] = 0.5*12*MeV;
338  highEnergyLimit[carbon] = 1e6*12*MeV;
339  //
340 
341  // Cross section
342 
344  eV,
345  scaleFactor );
346  tableCarbon->LoadData(fileCarbon);
347  tableData[carbon] = tableCarbon;
348 
349  // **********************************************************************************************
350 
351  oxygen = oxygenDef->GetParticleName();
352  tableFile[oxygen] = fileOxygen;
353 
354  //SEB
355  //lowEnergyLimit[oxygen] = lowEnergyLimitForA[5]* particle->GetAtomicMass();
356  //highEnergyLimit[oxygen] = 1e6* particle->GetAtomicMass()* MeV;
357  lowEnergyLimit[oxygen] = 0.5*16*MeV;
358  highEnergyLimit[oxygen] = 1e6*16*MeV;
359  //
360 
361  // Cross section
362 
364  eV,
365  scaleFactor );
366  tableOxygen->LoadData(fileOxygen);
367  tableData[oxygen] = tableOxygen;
368 
369  // **********************************************************************************************
370 
371  nitrogen = nitrogenDef->GetParticleName();
372  tableFile[nitrogen] = fileNitrogen;
373 
374  //SEB
375  //lowEnergyLimit[nitrogen] = lowEnergyLimitForA[5]* particle->GetAtomicMass();
376  //highEnergyLimit[nitrogen] = 1e6* particle->GetAtomicMass()* MeV;
377  lowEnergyLimit[nitrogen] = 0.5*14*MeV;
378  highEnergyLimit[nitrogen] = 1e6*14*MeV;
379  //
380 
381  // Cross section
382 
384  eV,
385  scaleFactor );
386  tableNitrogen->LoadData(fileNitrogen);
387  tableData[nitrogen] = tableNitrogen;
388 
389  // **********************************************************************************************
390 
391  silicon = siliconDef->GetParticleName();
392  tableFile[silicon] = fileSilicon;
393 
394  //lowEnergyLimit[silicon] = lowEnergyLimitForA[5]* particle->GetAtomicMass();
395  //highEnergyLimit[silicon] = 1e6* particle->GetAtomicMass()* MeV;
396  lowEnergyLimit[silicon] = 0.5*28*MeV;
397  highEnergyLimit[silicon] = 1e6*28*MeV;
398  //
399 
400  // Cross section
401 
403  eV,
404  scaleFactor );
405  tableSilicon->LoadData(fileSilicon);
406  tableData[silicon] = tableSilicon;
407 
408  // **********************************************************************************************
409 
410  iron = ironDef->GetParticleName();
411  tableFile[iron] = fileIron;
412 
413  //SEB
414  //lowEnergyLimit[iron] = lowEnergyLimitForA[5]* particle->GetAtomicMass();
415  //highEnergyLimit[iron] = 1e6* particle->GetAtomicMass()* MeV;
416  lowEnergyLimit[iron] = 0.5*56*MeV;
417  highEnergyLimit[iron] = 1e6*56*MeV;
418  //
419 
420  // Cross section
421 
423  eV,
424  scaleFactor );
425  tableIron->LoadData(fileIron);
426  tableData[iron] = tableIron;
427 
428  // **********************************************************************************************
429 
430  //SEB: not anymore
431  // ZF Following lines can be replaced by:
432  //SetLowEnergyLimit(lowEnergyLimit[particle->GetParticleName()]);
433  //SetHighEnergyLimit(highEnergyLimit[particle->GetParticleName()]);
434  // at least for HZE
435 
436  if (particle==protonDef)
437  {
438  SetLowEnergyLimit(lowEnergyLimit[proton]);
439  SetHighEnergyLimit(highEnergyLimit[proton]);
440  }
441 
442  if (particle==hydrogenDef)
443  {
444  SetLowEnergyLimit(lowEnergyLimit[hydrogen]);
445  SetHighEnergyLimit(highEnergyLimit[hydrogen]);
446  }
447 
448  if (particle==heliumDef)
449  {
450  SetLowEnergyLimit(lowEnergyLimit[helium]);
451  SetHighEnergyLimit(highEnergyLimit[helium]);
452  }
453 
454  if (particle==alphaPlusDef)
455  {
456  SetLowEnergyLimit(lowEnergyLimit[alphaPlus]);
457  SetHighEnergyLimit(highEnergyLimit[alphaPlus]);
458  }
459 
460  if (particle==alphaPlusPlusDef)
461  {
462  SetLowEnergyLimit(lowEnergyLimit[alphaPlusPlus]);
463  SetHighEnergyLimit(highEnergyLimit[alphaPlusPlus]);
464  }
465 
466  if (particle==lithiumDef)
467  {
468  SetLowEnergyLimit(lowEnergyLimit[lithium]);
469  SetHighEnergyLimit(highEnergyLimit[lithium]);
470  }
471 
472  if (particle==berylliumDef)
473  {
474  SetLowEnergyLimit(lowEnergyLimit[beryllium]);
475  SetHighEnergyLimit(highEnergyLimit[beryllium]);
476  }
477 
478  if (particle==boronDef)
479  {
480  SetLowEnergyLimit(lowEnergyLimit[boron]);
481  SetHighEnergyLimit(highEnergyLimit[boron]);
482  }
483 
484  if (particle==carbonDef)
485  {
486  SetLowEnergyLimit(lowEnergyLimit[carbon]);
487  SetHighEnergyLimit(highEnergyLimit[carbon]);
488  }
489 
490  if (particle==nitrogenDef)
491  {
492  SetLowEnergyLimit(lowEnergyLimit[nitrogen]);
493  SetHighEnergyLimit(highEnergyLimit[nitrogen]);
494  }
495 
496  if (particle==oxygenDef)
497  {
498  SetLowEnergyLimit(lowEnergyLimit[oxygen]);
499  SetHighEnergyLimit(highEnergyLimit[oxygen]);
500  }
501 
502  if (particle==siliconDef)
503  {
504  SetLowEnergyLimit(lowEnergyLimit[silicon]);
505  SetHighEnergyLimit(highEnergyLimit[silicon]);
506  }
507 
508  if (particle==ironDef)
509  {
510  SetLowEnergyLimit(lowEnergyLimit[iron]);
511  SetHighEnergyLimit(highEnergyLimit[iron]);
512  }
513 
514  //----------------------------------------------------------------------
515 
516  if( verboseLevel>0 )
517  {
518  G4cout << "Rudd ionisation model is initialized " << G4endl
519  << "Energy range: "
520  << LowEnergyLimit() / eV << " eV - "
521  << HighEnergyLimit() / keV << " keV for "
522  << particle->GetParticleName()
523  << G4endl;
524  }
525 
526  // Initialize water density pointer
528 
529  //
530 
531  fAtomDeexcitation = G4LossTableManager::Instance()->AtomDeexcitation();
532 
533  if (isInitialised) { return; }
535  isInitialised = true;
536 }
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:643
static G4LossTableManager * Instance()
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:636
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
Definition: G4Material.cc:602
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:503
static G4Proton * ProtonDefinition()
Definition: G4Proton.cc:88
virtual G4bool LoadData(const G4String &argFileName)
const G4String & GetParticleName() const
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:732
G4GLOB_DLL std::ostream G4cout
const std::vector< double > * GetNumMolPerVolTableFor(const G4Material *) const
static constexpr double m
Definition: G4SIunits.hh:129
static constexpr double eV
Definition: G4SIunits.hh:215
static G4DNAGenericIonsManager * Instance(void)
static G4IonTable * GetIonTable()
Definition: G4IonTable.hh:78
static G4DNAMolecularMaterial * Instance()
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
static MCTruthManager * instance
G4VAtomDeexcitation * AtomDeexcitation()
double G4double
Definition: G4Types.hh:76
void SetLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:739
static constexpr double keV
Definition: G4SIunits.hh:216
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition: G4VEmModel.cc:132
G4ParticleDefinition * GetIon(const G4String &name)

Here is the call graph for this function:

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

Implements G4VEmModel.

Definition at line 682 of file G4DNARuddIonisationExtendedModel.cc.

687 {
688  //SEB: particle->GetDefinition()->GetParticleName() is for eg. Fe56
689  // particle->GetDefinition()->GetPDGMass() is correct
690  // particle->GetDefinition()->GetAtomicNumber() is correct
691  // particle->GetDefinition()->GetAtomicMass() is correct
692 
693  if (verboseLevel > 3)
694  G4cout << "Calling SampleSecondaries() of G4DNARuddIonisationExtendedModel" << G4endl;
695 
696  G4double lowLim = 0;
697  G4double highLim = 0;
698 
699  // ZF: the following line summarizes the commented part
700 
701  if(particle->GetDefinition()->GetAtomicMass() <= 4) lowLim = killBelowEnergyForA[particle->GetDefinition()->GetAtomicMass()];
702 
703  else lowLim = killBelowEnergyForA[5]*particle->GetDefinition()->GetAtomicMass();
704 
705  /*
706 
707  if(particle->GetDefinition()->GetAtomicMass() >= 5) lowLim = killBelowEnergyForA[5]*particle->GetDefinition()->GetAtomicMass();
708 
709  if ( particle->GetDefinition() == G4Proton::ProtonDefinition()
710  || particle->GetDefinition() == instance->GetIon("hydrogen")
711  )
712 
713  lowLim = killBelowEnergyForA[1];
714 
715  if ( particle->GetDefinition() == instance->GetIon("alpha++")
716  || particle->GetDefinition() == instance->GetIon("alpha+")
717  || particle->GetDefinition() == instance->GetIon("helium")
718  )
719 
720  lowLim = killBelowEnergyForA[4];
721 
722  */
723  //
724 
725  G4double k = particle->GetKineticEnergy();
726 
727  const G4String& particleName = particle->GetDefinition()->GetParticleName();
728 
729  // SI - the following is useless since lowLim is already defined
730  /*
731  std::map< G4String,G4double,std::less<G4String> >::iterator pos1;
732  pos1 = lowEnergyLimit.find(particleName);
733 
734  if (pos1 != lowEnergyLimit.end())
735  {
736  lowLim = pos1->second;
737  }
738  */
739 
740  std::map< G4String,G4double,std::less<G4String> >::iterator pos2;
741  pos2 = highEnergyLimit.find(particleName);
742 
743  if (pos2 != highEnergyLimit.end())highLim = pos2->second;
744 
745  if (k >= lowLim && k < highLim)
746  {
747  G4ParticleDefinition* definition = particle->GetDefinition();
748  G4ParticleMomentum primaryDirection = particle->GetMomentumDirection();
749  /*
750  G4double particleMass = definition->GetPDGMass();
751  G4double totalEnergy = k + particleMass;
752  G4double pSquare = k*(totalEnergy+particleMass);
753  G4double totalMomentum = std::sqrt(pSquare);
754  */
755 
756  G4int ionizationShell = RandomSelect(k,particleName);
757 
758  // sample deexcitation
759  // here we assume that H_{2}O electronic levels are the same as Oxygen.
760  // this can be considered true with a rough 10% error in energy on K-shell,
761 
762  G4int secNumberInit = 0; // need to know at a certain point the energy of secondaries
763  G4int secNumberFinal = 0; // So I'll make the diference and then sum the energies
765  bindingEnergy = waterStructure.IonisationEnergy(ionizationShell);
766 
767  //SI: additional protection if tcs interpolation method is modified
768  if (k<bindingEnergy) return;
769  //
770 
771  G4int Z = 8;
772 
773  if(fAtomDeexcitation) {
775 
776  if (ionizationShell <5 && ionizationShell >1)
777  {
778  as = G4AtomicShellEnumerator(4-ionizationShell);
779  }
780  else if (ionizationShell <2)
781  {
782  as = G4AtomicShellEnumerator(3);
783  }
784 
785  // DEBUG
786  // if (ionizationShell == 4) {
787  //
788  // G4cout << "Z: " << Z << " as: " << as
789  // << " ionizationShell: " << ionizationShell << " bindingEnergy: "<< bindingEnergy/eV << G4endl;
790  // G4cout << "Press <Enter> key to continue..." << G4endl;
791  // G4cin.ignore();
792  // }
793 
794 
795  const G4AtomicShell* shell = fAtomDeexcitation->GetAtomicShell(Z, as);
796  secNumberInit = fvect->size();
797  fAtomDeexcitation->GenerateParticles(fvect, shell, Z, 0, 0);
798  secNumberFinal = fvect->size();
799  }
800 
801  G4double secondaryKinetic = RandomizeEjectedElectronEnergy(definition,k,ionizationShell);
802 
803  G4ThreeVector deltaDirection =
804  GetAngularDistribution()->SampleDirectionForShell(particle, secondaryKinetic,
805  Z, ionizationShell,
806  couple->GetMaterial());
807 
808  // SI: the following lines are not needed anymore
809  /*
810  G4double cosTheta = 0.;
811  G4double phi = 0.;
812  RandomizeEjectedElectronDirection(definition, k,secondaryKinetic, cosTheta, phi, ionizationShell);
813 
814  G4double sinTheta = std::sqrt(1.-cosTheta*cosTheta);
815  G4double dirX = sinTheta*std::cos(phi);
816  G4double dirY = sinTheta*std::sin(phi);
817  G4double dirZ = cosTheta;
818  G4ThreeVector deltaDirection(dirX,dirY,dirZ);
819  deltaDirection.rotateUz(primaryDirection);
820  */
821 
822  // Ignored for ions on electrons
823  /*
824  G4double deltaTotalMomentum = std::sqrt(secondaryKinetic*(secondaryKinetic + 2.*electron_mass_c2 ));
825 
826  G4double finalPx = totalMomentum*primaryDirection.x() - deltaTotalMomentum*deltaDirection.x();
827  G4double finalPy = totalMomentum*primaryDirection.y() - deltaTotalMomentum*deltaDirection.y();
828  G4double finalPz = totalMomentum*primaryDirection.z() - deltaTotalMomentum*deltaDirection.z();
829  G4double finalMomentum = std::sqrt(finalPx*finalPx+finalPy*finalPy+finalPz*finalPz);
830  finalPx /= finalMomentum;
831  finalPy /= finalMomentum;
832  finalPz /= finalMomentum;
833 
834  G4ThreeVector direction;
835  direction.set(finalPx,finalPy,finalPz);
836 
837  fParticleChangeForGamma->ProposeMomentumDirection(direction.unit()) ;
838  */
839 
841  G4double scatteredEnergy = k-bindingEnergy-secondaryKinetic;
842  G4double deexSecEnergy = 0;
843  for (G4int j=secNumberInit; j < secNumberFinal; j++) {
844 
845  deexSecEnergy = deexSecEnergy + (*fvect)[j]->GetKineticEnergy();
846 
847  }
848 
849  if (!statCode)
850  {
852  fParticleChangeForGamma->ProposeLocalEnergyDeposit(k-scatteredEnergy-secondaryKinetic-deexSecEnergy);
853  }
854  else
855  {
858  }
859 
860  G4DynamicParticle* dp = new G4DynamicParticle (G4Electron::Electron(),deltaDirection,secondaryKinetic) ;
861  fvect->push_back(dp);
862 
863  const G4Track * theIncomingTrack = fParticleChangeForGamma->GetCurrentTrack();
865  ionizationShell,
866  theIncomingTrack);
867  }
868 
869  // SI - not useful since low energy of model is 0 eV
870 
871  if (k < lowLim)
872  {
876  }
877 
878 }
G4double GetKineticEnergy() const
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:619
if(nIso!=0)
G4ParticleDefinition * GetDefinition() const
int G4int
Definition: G4Types.hh:78
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
const G4String & GetParticleName() const
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
virtual G4ThreeVector & SampleDirectionForShell(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, G4int shellID, const G4Material *)
virtual const G4AtomicShell * GetAtomicShell(G4int Z, G4AtomicShellEnumerator shell)=0
G4GLOB_DLL std::ostream G4cout
const G4ThreeVector & GetMomentumDirection() const
G4int GetAtomicMass() const
static G4DNAChemistryManager * Instance()
void CreateWaterMolecule(ElectronicModification, G4int, const G4Track *)
const G4Track * GetCurrentTrack() const
static G4Electron * Electron()
Definition: G4Electron.cc:94
void SetProposedKineticEnergy(G4double proposedKinEnergy)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
void GenerateParticles(std::vector< G4DynamicParticle * > *secVect, const G4AtomicShell *, G4int Z, G4int coupleIndex)
G4double bindingEnergy(G4int A, G4int Z)
G4AtomicShellEnumerator
const G4Material * GetMaterial() const

Here is the call graph for this function:

void G4DNARuddIonisationExtendedModel::SelectStationary ( G4bool  input)
inline

Definition at line 172 of file G4DNARuddIonisationExtendedModel.hh.

173 {
174  statCode = input;
175 }

Member Data Documentation

G4ParticleChangeForGamma* G4DNARuddIonisationExtendedModel::fParticleChangeForGamma
protected

Definition at line 74 of file G4DNARuddIonisationExtendedModel.hh.


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