46 fpMolWaterDensity = 0;
47 numberOfPartialCrossSections[0] = 0;
48 numberOfPartialCrossSections[1] = 0;
49 numberOfPartialCrossSections[2] = 0;
61 G4cout <<
"Dingfelder charge decrease model is constructed " <<
G4endl;
84 G4cout <<
"Calling G4DNADingfelderChargeDecreaseModel::Initialise()"
107 lowEnergyLimit[alphaPlusPlus] = 1. *
keV;
108 highEnergyLimit[alphaPlusPlus] = 400. *
MeV;
111 lowEnergyLimit[alphaPlus] = 1. *
keV;
112 highEnergyLimit[alphaPlus] = 400. *
MeV;
116 if (particle==protonDef)
122 if (particle==alphaPlusPlusDef)
128 if (particle==alphaPlusDef)
147 numberOfPartialCrossSections[0] = 1;
150 f0[0][1]=1.; a0[0][1]=0.95;
171 numberOfPartialCrossSections[1] = 2;
185 numberOfPartialCrossSections[2] = 1;
191 G4cout <<
"Dingfelder charge decrease model is initialized " <<
G4endl
205 isInitialised =
true;
217 if (verboseLevel > 3)
220 <<
"Calling CrossSectionPerVolume() of G4DNADingfelderChargeDecreaseModel"
232 particleDefinition != instance->
GetIon(
"alpha++")
234 particleDefinition != instance->
GetIon(
"alpha+")
245 if(waterDensity!= 0.0)
250 std::map< G4String,G4double,std::less<G4String> >::iterator pos1;
251 pos1 = lowEnergyLimit.find(particleName);
253 if (pos1 != lowEnergyLimit.end())
255 lowLim = pos1->second;
258 std::map< G4String,G4double,std::less<G4String> >::iterator pos2;
259 pos2 = highEnergyLimit.find(particleName);
261 if (pos2 != highEnergyLimit.end())
263 highLim = pos2->second;
266 if (k >= lowLim && k < highLim)
268 crossSection = Sum(k,particleDefinition);
271 if (verboseLevel > 2)
273 G4cout <<
"_______________________________________" <<
G4endl;
274 G4cout <<
"G4DNADingfelderChargeDecreaeModel" <<
G4endl;
275 G4cout <<
"Kinetic energy(eV)=" << k/
eV <<
"particle :" << particleName <<
G4endl;
276 G4cout <<
"Cross section per water molecule (cm^2)=" << crossSection/
cm/
cm <<
G4endl;
277 G4cout <<
"Cross section per water molecule (cm^-1)=" << crossSection*
278 waterDensity/(1./
cm) << G4endl;
283 return crossSection*waterDensity;
297 if (verboseLevel > 3)
300 <<
"Calling SampleSecondaries() of G4DNADingfelderChargeDecreaseModel"
310 G4int finalStateIndex = RandomSelect(inK,definition);
312 G4int n = NumberOfFinalStates(definition, finalStateIndex);
313 G4double waterBindingEnergy = WaterBindingEnergyConstant(definition, finalStateIndex);
314 G4double outgoingParticleBindingEnergy = OutgoingParticleBindingEnergyConstant(definition, finalStateIndex);
321 - waterBindingEnergy + outgoingParticleBindingEnergy;
328 - waterBindingEnergy + outgoingParticleBindingEnergy;
334 G4Exception(
"G4DNADingfelderChargeDecreaseModel::SampleSecondaries",
"em0004",
347 + waterBindingEnergy - outgoingParticleBindingEnergy);
350 + waterBindingEnergy - outgoingParticleBindingEnergy);
356 fvect->push_back(dp);
368 G4int finalStateIndex)
377 if (particleDefinition == instance->
GetIon(
"alpha++"))
379 if (finalStateIndex == 0)
384 if (particleDefinition == instance->
GetIon(
"alpha+"))
393 G4int finalStateIndex)
398 return instance->
GetIon(
"hydrogen");
400 if (particleDefinition == instance->
GetIon(
"alpha++"))
402 if (finalStateIndex == 0)
403 return instance->
GetIon(
"alpha+");
404 return instance->
GetIon(
"helium");
407 if (particleDefinition == instance->
GetIon(
"alpha+"))
408 return instance->
GetIon(
"helium");
416 G4int finalStateIndex)
427 if (particleDefinition == instance->
GetIon(
"alpha++"))
435 if (finalStateIndex == 0)
438 return 10.79 * 2 *
eV;
441 if (particleDefinition == instance->
GetIon(
"alpha+"))
458 G4int finalStateIndex)
465 if (particleDefinition == instance->
GetIon(
"alpha++"))
470 if (finalStateIndex == 0)
473 return (54.509 + 24.587) *
eV;
476 if (particleDefinition == instance->
GetIon(
"alpha+"))
489 G4double G4DNADingfelderChargeDecreaseModel::PartialCrossSection(
G4double k,
493 G4int particleTypeIndex = 0;
498 particleTypeIndex = 0;
499 if (particleDefinition == instance->
GetIon(
"alpha++"))
500 particleTypeIndex = 1;
501 if (particleDefinition == instance->
GetIon(
"alpha+"))
502 particleTypeIndex = 2;
522 if (x1[index][particleTypeIndex] < x0[index][particleTypeIndex])
532 x1[index][particleTypeIndex] = x0[index][particleTypeIndex]
533 + std::pow((a0[index][particleTypeIndex] - a1[index][particleTypeIndex])
534 / (c0[index][particleTypeIndex]
535 * d0[index][particleTypeIndex]),
536 1. / (d0[index][particleTypeIndex] - 1.));
537 b1[index][particleTypeIndex] = (a0[index][particleTypeIndex]
538 - a1[index][particleTypeIndex]) * x1[index][particleTypeIndex]
539 + b0[index][particleTypeIndex]
540 - c0[index][particleTypeIndex]
541 * std::pow(x1[index][particleTypeIndex]
542 - x0[index][particleTypeIndex],
543 d0[index][particleTypeIndex]);
549 if (x < x0[index][particleTypeIndex])
550 y = a0[index][particleTypeIndex] * x + b0[index][particleTypeIndex];
551 else if (x < x1[index][particleTypeIndex])
552 y = a0[index][particleTypeIndex] * x + b0[index][particleTypeIndex]
553 - c0[index][particleTypeIndex]
554 * std::pow(x - x0[index][particleTypeIndex],
555 d0[index][particleTypeIndex]);
557 y = a1[index][particleTypeIndex] * x + b1[index][particleTypeIndex];
559 return f0[index][particleTypeIndex] * std::pow(10., y) *
m *
m;
563 G4int G4DNADingfelderChargeDecreaseModel::RandomSelect(
G4double k,
566 G4int particleTypeIndex = 0;
571 particleTypeIndex = 0;
572 if (particleDefinition == instance->
GetIon(
"alpha++"))
573 particleTypeIndex = 1;
574 if (particleDefinition == instance->
GetIon(
"alpha+"))
575 particleTypeIndex = 2;
577 const G4int n = numberOfPartialCrossSections[particleTypeIndex];
585 values[i] = PartialCrossSection(k, i, particleDefinition);
596 if (values[i] >
value)
612 G4int particleTypeIndex = 0;
617 particleTypeIndex = 0;
618 if (particleDefinition == instance->
GetIon(
"alpha++"))
619 particleTypeIndex = 1;
620 if (particleDefinition == instance->
GetIon(
"alpha+"))
621 particleTypeIndex = 2;
625 for (
G4int i = 0; i < numberOfPartialCrossSections[particleTypeIndex]; i++)
627 totalCrossSection += PartialCrossSection(k, i, particleDefinition);
629 return totalCrossSection;
G4DNADingfelderChargeDecreaseModel(const G4ParticleDefinition *p=0, const G4String &nam="DNADingfelderChargeDecreaseModel")
G4double LowEnergyLimit() const
static constexpr double proton_mass_c2
G4double GetKineticEnergy() const
G4double HighEnergyLimit() const
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
virtual G4double CrossSectionPerVolume(const G4Material *material, const G4ParticleDefinition *p, G4double ekin, G4double emin, G4double emax)
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
static G4Proton * ProtonDefinition()
G4ParticleDefinition * GetDefinition() const
virtual ~G4DNADingfelderChargeDecreaseModel()
G4ParticleChangeForGamma * fParticleChangeForGamma
const G4String & GetParticleName() const
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
static constexpr double electron_mass_c2
void SetHighEnergyLimit(G4double)
G4GLOB_DLL std::ostream G4cout
const std::vector< double > * GetNumMolPerVolTableFor(const G4Material *) const
static constexpr double m
const XML_Char int const XML_Char * value
const G4ThreeVector & GetMomentumDirection() const
static constexpr double cm
static G4Proton * Proton()
static constexpr double eV
static G4DNAGenericIonsManager * Instance(void)
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
static G4DNAChemistryManager * Instance()
static G4DNAMolecularMaterial * Instance()
G4double GetPDGMass() const
void CreateWaterMolecule(ElectronicModification, G4int, const G4Track *)
const G4Track * GetCurrentTrack() const
static constexpr double MeV
void ProposeTrackStatus(G4TrackStatus status)
void SetLowEnergyLimit(G4double)
static constexpr double keV
G4ParticleChangeForGamma * GetParticleChangeForGamma()
G4ParticleDefinition * GetIon(const G4String &name)