77 G4cout <<
"Rudd ionisation model is constructed " <<
G4endl;
99 std::map<G4String, G4DNACrossSectionDataSet*, std::less<G4String> >::iterator
pos;
120 G4cout <<
"Calling G4DNARuddIonisationModel::Initialise()" <<
G4endl;
125 G4String fileProton(
"dna/sigma_ionisation_p_rudd");
126 G4String fileHydrogen(
"dna/sigma_ionisation_h_rudd");
127 G4String fileAlphaPlusPlus(
"dna/sigma_ionisation_alphaplusplus_rudd");
128 G4String fileAlphaPlus(
"dna/sigma_ionisation_alphaplus_rudd");
129 G4String fileHelium(
"dna/sigma_ionisation_he_rudd");
178 tableHydrogen->
LoadData(fileHydrogen);
185 tableFile[alphaPlusPlus] = fileAlphaPlusPlus;
195 tableAlphaPlusPlus->
LoadData(fileAlphaPlusPlus);
197 tableData[alphaPlusPlus] = tableAlphaPlusPlus;
212 tableAlphaPlus->
LoadData(fileAlphaPlus);
233 if (particle==protonDef)
239 if (particle==hydrogenDef)
245 if (particle==heliumDef)
251 if (particle==alphaPlusDef)
257 if (particle==alphaPlusPlusDef)
265 G4cout <<
"Rudd ionisation model is initialized " <<
G4endl
297 G4cout <<
"Calling CrossSectionPerVolume() of G4DNARuddIonisationModel"
309 particleDefinition != instance->
GetIon(
"hydrogen")
311 particleDefinition != instance->
GetIon(
"alpha++")
313 particleDefinition != instance->
GetIon(
"alpha+")
315 particleDefinition != instance->
GetIon(
"helium")
323 || particleDefinition == instance->
GetIon(
"hydrogen")
328 if ( particleDefinition == instance->
GetIon(
"alpha++")
329 || particleDefinition == instance->
GetIon(
"alpha+")
330 || particleDefinition == instance->
GetIon(
"helium")
340 if(waterDensity!= 0.0)
356 std::map< G4String,G4double,std::less<G4String> >::iterator pos2;
361 highLim = pos2->second;
368 if (k < lowLim) k = lowLim;
372 std::map< G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator
pos;
385 G4Exception(
"G4DNARuddIonisationModel::CrossSectionPerVolume",
"em0002",
393 G4cout <<
"__________________________________" <<
G4endl;
394 G4cout <<
"G4DNARuddIonisationModel - XS INFO START" <<
G4endl;
396 G4cout <<
"Cross section per water molecule (cm^2)=" << sigma/
cm/
cm <<
G4endl;
397 G4cout <<
"Cross section per water molecule (cm^-1)=" << sigma*waterDensity/(1./
cm) << G4endl;
399 G4cout <<
"G4DNARuddIonisationModel - XS INFO END" <<
G4endl;
407 G4cout <<
"Warning : RuddIonisationModel: WATER DENSITY IS NULL" <<
G4endl;
411 return sigma*waterDensity;
426 G4cout <<
"Calling SampleSecondaries() of G4DNARuddIonisationModel"
464 std::map< G4String,G4double,std::less<G4String> >::iterator pos2;
469 highLim = pos2->second;
472 if (k >= lowLim && k <= highLim)
489 if (k<bindingEnergy)
return;
504 fvect->push_back(dp);
530 G4int secNumberInit = 0;
531 G4int secNumberFinal = 0;
537 if (ionizationShell <5 && ionizationShell >1)
541 else if (ionizationShell <2)
556 secNumberInit = fvect->size();
558 secNumberFinal = fvect->size();
561 G4double scatteredEnergy = k-bindingEnergy-secondaryKinetic;
563 for (
G4int j=secNumberInit; j < secNumberFinal; j++)
565 deexSecEnergy = deexSecEnergy + (*fvect)[j]->GetKineticEnergy();
615 G4double maximumKineticEnergyTransfer = 0.;
621 || particleDefinition == instance->
GetIon(
"hydrogen"))
626 else if (particleDefinition == instance->
GetIon(
"helium")
627 || particleDefinition == instance->
GetIon(
"alpha+")
628 || particleDefinition == instance->
GetIon(
"alpha++"))
630 maximumKineticEnergyTransfer = 4. * (0.511 / 3728) * k;
641 if (differentialCrossSection >= crossSectionMaximum)
642 crossSectionMaximum = differentialCrossSection;
649 secElecKinetic =
G4UniformRand()* maximumKineticEnergyTransfer;
655 return (secElecKinetic);
703 G4int ionizationLevelIndex)
721 const G4int j = ionizationLevelIndex;
736 const G4double Bj[5] = { 12.60 *
eV, 14.70 *
eV, 18.40 *
eV, 32.20 *
eV, 540
770 const G4double Gj[5] = { 0.99, 1.11, 1.11, 0.52, 1. };
780 G4double w = wBig / Bj[ionizationLevelIndex];
789 G4bool isProtonOrHydrogen =
false;
793 || particleDefinition == instance->
GetIon(
"hydrogen"))
795 isProtonOrHydrogen =
true;
799 else if (particleDefinition == instance->
GetIon(
"helium")
800 || particleDefinition == instance->
GetIon(
"alpha+")
801 || particleDefinition == instance->
GetIon(
"alpha++"))
804 tau = (0.511 / 3728.) * k;
808 * std::pow((Ry / Bj[ionizationLevelIndex]), 2);
814 G4double v2 = tau / Bj[ionizationLevelIndex];
819 G4double wc = 4. * v2 - 2. * v - (Ry / (4. * Bj[ionizationLevelIndex]));
821 wc = 4. * v2 - 2. * v
824 G4double L1 = (C1 * std::pow(v, (D1))) / (1. + E1 * std::pow(v, (D1 + 4.)));
825 G4double L2 = C2 * std::pow(v, (D2));
826 G4double H1 = (A1 * std::log(1. + v2)) / (v2 + (B1 / v2));
827 G4double H2 = (A2 / v2) + (B2 / (v2 * v2));
830 G4double F2 = (L2 * H2) / (L2 + H2);
834 * (S / Bj[ionizationLevelIndex])
836 / (std::pow((1. + w), 3)
837 * (1. +
G4Exp(alphaConst * (w - wc) / v))));
843 / (std::pow((1. + w), 3)
844 * (1. +
G4Exp(alphaConst * (w - wc) / v))));
846 if ((particleDefinition == instance->
GetIon(
"hydrogen"))
847 && (ionizationLevelIndex == 4))
852 / (std::pow((1. + w), 3)
853 * (1. +
G4Exp(alphaConst * (w - wc) / v))));
858 if (isProtonOrHydrogen)
863 if (particleDefinition == instance->
GetIon(
"alpha++"))
873 else if (particleDefinition == instance->
GetIon(
"alpha+"))
885 else if (particleDefinition == instance->
GetIon(
"helium"))
901 sigma = Gj[j] * (S / Bj[ionizationLevelIndex])
903 / (std::pow((1. + w), 3)
904 * (1. +
G4Exp(alphaConst * (w - wc) / v))));
910 / (std::pow((1. + w), 3)
911 * (1. +
G4Exp(alphaConst * (w - wc) / v))));
923 return zEff * zEff * sigma;
939 G4double r =
R(t, energyTransferred, slaterEffectiveChg, shellNumber);
955 G4double r =
R(t, energyTransferred, slaterEffectiveChg, shellNumber);
957 -
G4Exp(-2 * r) * (((2. * r * r + 2.) * r + 2.) * r + 1.);
973 G4double r =
R(t, energyTransferred, slaterEffectiveChg, shellNumber);
976 * ((((2. / 3. * r + 4. / 3.) * r + 2.) * r + 2.) * r + 1.);
991 G4double tElectron = 0.511 / 3728. * t;
994 G4double value = std::sqrt(2. * tElectron / H) / (energyTransferred / H)
995 * (slaterEffectiveChg / shellNumber);
1011 }
else if (particleDefinition == instance->
GetIon(
"hydrogen"))
1015 return ((0.6 / (1 +
G4Exp(value))) + 0.9);
1036 std::map<G4String, G4DNACrossSectionDataSet*, std::less<G4String> >::iterator
pos;
1055 value += valuesBuffer[i];
1066 if (valuesBuffer[i] > value)
1068 delete[] valuesBuffer;
1071 value -= valuesBuffer[i];
1075 delete[] valuesBuffer;
1080 G4Exception(
"G4DNARuddIonisationModel::RandomSelect",
1083 "Model not applicable to particle type.");
1103 std::map<G4String, G4double, std::less<G4String> >::iterator pos1;
1108 lowLim = pos1->second;
1111 std::map<G4String, G4double, std::less<G4String> >::iterator pos2;
1116 highLim = pos2->second;
1119 if (k >= lowLim && k <= highLim)
1121 std::map<G4String, G4DNACrossSectionDataSet*, std::less<G4String> >::iterator
pos;
1133 G4Exception(
"G4DNARuddIonisationModel::PartialCrossSection",
1136 "Model not applicable to particle type.");
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
G4double S_2p(G4double t, G4double energyTransferred, G4double slaterEffectiveChg, G4double shellNumber)
virtual G4bool LoadData(const G4String &argFileName)
static const G4double pos
virtual const G4AtomicShell * GetAtomicShell(G4int Z, G4AtomicShellEnumerator shell)=0
virtual G4double CrossSectionPerVolume(const G4Material *material, const G4ParticleDefinition *p, G4double ekin, G4double emin, G4double emax)
G4VAtomDeexcitation * AtomDeexcitation()
static constexpr double MeV
void SetHighEnergyLimit(G4double)
virtual size_t NumberOfComponents(void) const
static constexpr double keV
std::map< G4String, G4double, std::less< G4String > > lowEnergyLimit
const G4ThreeVector & GetMomentumDirection() const
const G4String & GetParticleName() const
G4double GetPDGCharge() const
void SetProposedKineticEnergy(G4double proposedKinEnergy)
G4double RandomizeEjectedElectronEnergy(G4ParticleDefinition *particleDefinition, G4double incomingParticleEnergy, G4int shell)
G4VAtomDeexcitation * fAtomDeexcitation
static G4Proton * Proton()
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
virtual ~G4DNARuddIonisationModel()
G4double IonisationEnergy(G4int level)
G4double S_2s(G4double t, G4double energyTransferred, G4double slaterEffectiveChg, G4double shellNumber)
G4double LowEnergyLimit() const
const std::vector< double > * GetNumMolPerVolTableFor(const G4Material *) const
Retrieve a table of molecular densities (number of molecules per unit volume) in the G4 unit system f...
static constexpr double proton_mass_c2
G4double PartialCrossSection(const G4Track &track)
static constexpr double m
G4double lowEnergyLimitOfModelForZ1
void SetDeexcitationFlag(G4bool val)
G4DNAWaterIonisationStructure waterStructure
G4ParticleDefinition * GetDefinition() const
void GenerateParticles(std::vector< G4DynamicParticle * > *secVect, const G4AtomicShell *, G4int Z, G4int coupleIndex)
G4double slaterEffectiveCharge[3]
const XML_Char int const XML_Char * value
static G4Proton * ProtonDefinition()
G4double killBelowEnergyForZ2
static constexpr double electron_mass_c2
static G4DNAMolecularMaterial * Instance()
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
G4double lowEnergyLimitForZ2
G4double S_1s(G4double t, G4double energyTransferred, G4double slaterEffectiveChg, G4double shellNumber)
void CreateWaterMolecule(ElectronicModification, G4int, const G4Track *)
virtual G4ThreeVector & SampleDirectionForShell(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, G4int shellID, const G4Material *)
static constexpr double eV
void SetAngularDistribution(G4VEmAngularDistribution *)
const G4Track * GetCurrentTrack() const
std::map< G4String, G4double, std::less< G4String > > highEnergyLimit
static G4Electron * Electron()
G4double DifferentialCrossSection(G4ParticleDefinition *particleDefinition, G4double k, G4double energyTransfer, G4int shell)
G4int RandomSelect(G4double energy, const G4String &particle)
const std::vector< G4double > * fpWaterDensity
G4ParticleChangeForGamma * fParticleChangeForGamma
static constexpr double eplus
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
G4ParticleDefinition * GetIon(const G4String &name)
void SetLowEnergyLimit(G4double)
static MCTruthManager * instance
G4int GetLeptonNumber() const
static G4DNAGenericIonsManager * Instance(void)
G4double GetKineticEnergy() const
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
virtual G4double FindValue(G4double e, G4int componentId=0) const
static constexpr double cm
G4GLOB_DLL std::ostream G4cout
virtual const G4VEMDataSet * GetComponent(G4int componentId) const
static G4LossTableManager * Instance()
G4VEmAngularDistribution * GetAngularDistribution()
const G4Material * GetMaterial() const
static constexpr double pi
G4ParticleChangeForGamma * GetParticleChangeForGamma()
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
G4DNARuddIonisationModel(const G4ParticleDefinition *p=0, const G4String &nam="DNARuddIonisationModel")
G4double R(G4double t, G4double energyTransferred, G4double slaterEffectiveChg, G4double shellNumber)
G4double bindingEnergy(G4int A, G4int Z)
G4double HighEnergyLimit() const
G4double lowEnergyLimitForZ1
G4double lowEnergyLimitOfModelForZ2
G4double CorrectionFactor(G4ParticleDefinition *particleDefinition, G4double k)
void ProposeTrackStatus(G4TrackStatus status)
static constexpr double Bohr_radius
G4double killBelowEnergyForZ1
const G4DynamicParticle * GetDynamicParticle() const
static G4DNAChemistryManager * Instance()
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
G4double Sum(G4double energy, const G4String &particle)