87 std::istringstream& dataStream )
88 : Isotope_(WhichIsotope),
89 MetaState_(WhichMetaState),
91 YieldType_(WhichYieldType),
100 }
catch (std::exception e)
115 std::istringstream& dataStream)
116 : Isotope_(WhichIsotope),
117 MetaState_(WhichMetaState),
119 YieldType_(WhichYieldType),
120 Verbosity_(Verbosity)
128 }
catch (std::exception e)
179 }
catch (std::exception& e)
204 std::vector< G4ReactionProduct* >* Alphas =
new std::vector< G4ReactionProduct* >;
205 std::vector< G4ReactionProduct* >* Neutrons =
new std::vector< G4ReactionProduct* >;
206 std::vector< G4ReactionProduct* >* Gammas =
new std::vector< G4ReactionProduct* >;
252 G4cout <<
" -- First daughter product sampled" <<
G4endl;
269 if(
Verbosity_ & G4FFGEnumerations::DAUGHTER_INFO)
274 G4cout <<
" -- Second daughter product sampled" <<
G4endl;
317 G4int icounter_max=1024;
321 if ( icounter > icounter_max ) {
322 G4cout <<
"Loop-counter exceeded the threshold value at " << __LINE__ <<
"th line of " << __FILE__ <<
"." <<
G4endl;
345 G4double NumeratorSqrt, NumeratorOther, Denominator;
350 if(Alphas->size() > 0)
375 MassRatio = 1 / MassRatio;
383 const G4double MeanAlphaAngle = 0.3644 * MassRatio * MassRatio * MassRatio
384 - 1.9766 * MassRatio * MassRatio
391 const G4double MeanAlphaAngleStdDev = 0.0523598776;
394 for(
unsigned int i = 0; i < Alphas->size(); ++i)
397 Theta = MeanAlphaAngle + PlusMinus;
401 }
else if(Theta >
pi)
403 Theta = (2 *
pi - Theta);
407 Direction.setRThetaPhi(1.0, Theta, Phi);
408 Magnitude = std::sqrt(2 * Alphas->at(i)->GetKineticEnergy() * Alphas->at(i)->GetDefinition()->GetPDGMass());
409 Alphas->at(i)->SetMomentum(Direction * Magnitude);
410 ResultantVector += Alphas->at(i)->GetMomentum();
415 if(Neutrons->size() != 0)
417 for(
unsigned int i = 0; i < Neutrons->size(); ++i)
422 Direction.setRThetaPhi(1.0, Theta, Phi);
423 Magnitude = std::sqrt(2 * Neutrons->at(i)->GetKineticEnergy() * Neutrons->at(i)->GetDefinition()->GetPDGMass());
424 Neutrons->at(i)->SetMomentum(Direction * Magnitude);
425 ResultantVector += Neutrons->at(i)->GetMomentum();
430 if(Gammas->size() != 0)
432 for(
unsigned int i = 0; i < Gammas->size(); ++i)
437 Direction.setRThetaPhi(1.0, Theta, Phi);
438 Magnitude = Gammas->at(i)->GetKineticEnergy() / CLHEP::c_light;
439 Gammas->at(i)->SetMomentum(Direction * Magnitude);
440 ResultantVector += Gammas->at(i)->GetMomentum();
449 G4double ResultantX, ResultantY, ResultantZ;
450 ResultantX = ResultantVector.getX();
451 ResultantY = ResultantVector.getY();
452 ResultantZ = ResultantVector.getZ();
456 LightFragment = FirstDaughter;
457 HeavyFragment = SecondDaughter;
460 LightFragment = SecondDaughter;
461 HeavyFragment = FirstDaughter;
466 LightFragmentDirection.setRThetaPhi(1.0, 0, 0);
501 NumeratorSqrt = std::sqrt(LightFragmentMass *
502 (LightFragmentMass * (2 * HeavyFragmentMass * FragmentsKE
503 - ResultantX * ResultantX
504 - ResultantY * ResultantY)
505 + HeavyFragmentMass * (2 * HeavyFragmentMass * FragmentsKE
506 - ResultantX * ResultantX
507 - ResultantY * ResultantY
508 - ResultantZ - ResultantZ)));
511 NumeratorOther = LightFragmentMass * ResultantZ;
514 Denominator = LightFragmentMass + HeavyFragmentMass;
517 const G4double LightFragmentMomentum = (NumeratorSqrt + NumeratorOther) / Denominator;
518 const G4double HeavyFragmentMomentum = std::sqrt(ResultantX * ResultantX
519 + ResultantY * ResultantY
523 LightFragment->
SetMomentum(LightFragmentDirection * LightFragmentMomentum);
524 HeavyFragmentDirection.setX(-ResultantX);
525 HeavyFragmentDirection.setY(-ResultantY);
526 HeavyFragmentDirection.setZ((-LightFragmentMomentum - ResultantZ));
528 HeavyFragmentDirection.setR(1.0);
529 HeavyFragment->
SetMomentum(HeavyFragmentDirection * HeavyFragmentMomentum);
536 G4cout <<
" -- Daugher product momenta finalized" <<
G4endl;
548 for(
unsigned int i = 0; i < Neutrons->size(); i++)
553 for(
unsigned int i = 0; i < Gammas->size(); i++)
558 for(
unsigned int i = 0; i < Alphas->size(); i++)
570 RotationAxis.setRThetaPhi(1.0, Theta, Phi);
573 ResultantVector.set(0.0, 0.0, 0.0);
574 for(
unsigned int i = 0; i < FissionProducts->size(); i++)
576 Direction = FissionProducts->at(i)->GetMomentumDirection();
577 Direction.rotateUz(RotationAxis);
578 FissionProducts->at(i)->SetMomentumDirection(Direction);
579 ResultantVector += FissionProducts->at(i)->GetMomentum();
584 G4double PossibleImbalance = ResultantVector.mag();
585 if(PossibleImbalance > 0.01)
587 std::ostringstream Temp;
588 Temp <<
"Momenta imbalance of ";
589 Temp << PossibleImbalance / (
MeV / CLHEP::c_light);
590 Temp <<
" MeV/c in the system";
591 G4Exception(
"G4FissionProductYieldDist::G4GetFission()",
594 "Results may not be valid");
598 delete FirstDaughter;
599 delete SecondDaughter;
605 return FissionProducts;
692 G4bool lowerExists =
false;
693 G4bool higherExists =
false;
708 }
else if(energyGroup == YieldEnergyGroups_ - 1)
727 G4Ions* FoundParticle = NULL;
728 if(isExact || YieldEnergyGroups_ == 1)
735 if(RandomParticle <=
Trees_[tree].ProbabilityRangeEnd[energyGroup])
746 while((RangeIsSmaller = (RandomParticle < Branch->ProbabilityRangeBottom[energyGroup]))
752 Branch = Branch->
Left;
754 Branch = Branch->
Right;
759 }
else if(lowerExists && higherExists)
771 return FoundParticle;
776 G4bool LowerEnergyGroupExists )
780 G4Ions* FoundParticle = NULL;
782 G4int NextNearestEnergy;
785 if(LowerEnergyGroupExists ==
true)
788 NextNearestEnergy = NearestEnergy - 1;
792 NextNearestEnergy = 1;
795 for(
G4int Tree = 0; Tree <
TreeCount_ && FoundParticle == NULL; Tree++)
804 return FoundParticle;
809 G4int LowerEnergyGroup )
813 G4Ions* FoundParticle = NULL;
814 G4int HigherEnergyGroup = LowerEnergyGroup + 1;
816 for(
G4int Tree = 0; Tree <
TreeCount_ && FoundParticle == NULL; Tree++)
825 return FoundParticle;
844 || EnergyGroup1 == EnergyGroup2
852 G4Ions* FoundParticle = NULL;
864 if(RandomParticle < RangeAtIncidentEnergy)
878 if(RandomParticle > RangeAtIncidentEnergy)
891 Particle = FoundParticle;
907 G4int NumberOfAlphasToProduce;
922 for(
int i = 0; i < NumberOfAlphasToProduce; i++)
941 G4int NeutronProduction;
946 for(
int i = 0; i < NeutronProduction; i++)
970 G4int Z = (Product -
A) / 1000;
1042 std::ostringstream DirectoryName;
1047 return DirectoryName.str();
1058 std::ostringstream FileName;
1061 if(Isotope < 100000)
1070 return FileName.str();
1081 return DynamicParticle;
1091 G4int A = Isotope % 1000;
1092 G4int Z = (Isotope -
A) / 1000;
1095 std::ostringstream IsotopeName;
1097 IsotopeName << Z <<
"_" <<
A;
1114 return IsotopeName.str();
1163 for(
G4int i = 0; i < ProductCount; i++)
1226 TotalAlphaEnergy = 0;
1230 for(
unsigned int i = 0; i < Alphas->size(); i++)
1236 Alphas->at(i)->SetKineticEnergy(AlphaEnergy);
1239 TotalAlphaEnergy += AlphaEnergy;
1243 MeanAlphaEnergy -= 0.1;
1268 G4int icounter_max=1024;
1272 if ( icounter > icounter_max ) {
1273 G4cout <<
"Loop-counter exceeded the threshold value at " << __LINE__ <<
"th line of " << __FILE__ <<
"." <<
G4endl;
1294 Gammas->back()->SetTotalEnergy(SampleEnergy);
1311 Gammas->back()->SetTotalEnergy(SampleEnergy);
1334 G4int icounter_max=1024;
1338 if ( icounter > icounter_max ) {
1339 G4cout <<
"Loop-counter exceeded the threshold value at " << __LINE__ <<
"th line of " << __FILE__ <<
"." <<
G4endl;
1342 TotalNeutronEnergy = 0;
1347 for(
unsigned int i = 0; i < Neutrons->size(); i++)
1351 Neutrons->at(i)->SetKineticEnergy(NeutronEnergy);
1354 TotalNeutronEnergy +=NeutronEnergy;
1386 + *(WhichNubar + 2) * BFactor;
1387 while(*WhichNubar != -1)
1392 + *(WhichNubar + 2) * BFactor;
1401 while(*WhichNubar != -1)
1423 NewBranch->Left = NULL;
1424 NewBranch->Right = NULL;
1477 while(BranchPosition > 1)
1479 if(BranchPosition & 1)
1482 WhichBranch = &((*WhichBranch)->Right);
1486 WhichBranch = &((*WhichBranch)->Left);
1489 BranchPosition >>= 1;
1492 *WhichBranch = NewBranch;
1504 G4int WhichTree = 0;
1535 delete Branch->
Left;
1537 delete Branch->
Right;
#define G4FFG_RECURSIVE_FUNCTIONLEAVE__
virtual ~G4FissionProductYieldDist(void)
Default deconstructor.
static G4Pow * GetInstance()
const G4FFGEnumerations::YieldType YieldType_
The type of yield to be used: INDEPENDET or CUMULATIVE.
G4ENDFYieldDataContainer * G4GetYield(G4int WhichYield)
Returns the data for the requested fission product.
G4FFGEnumerations::MetaState GetMetaState(void)
Get the meta state.
G4double powA(G4double A, G4double y) const
MetaState
ENDF format provides for 3 isomers - 1 ground state and 2 meta states.
G4double G4SampleGaussian(G4double Mean, G4double StdDev)
Returns a double value taken from a Gaussian distribution about Mean and with a standard deviation of...
void Copy(G4int Elements, T *To, T *From)
Copy values from one array to another.
G4double G4SampleWatt(G4int WhatIsotope, G4FFGEnumerations::FissionCause WhatCause, G4double WhatEnergy)
Samples the Watt fission spectrum for the selected isotope, using an algorithm adopted from Ref...
G4int IncidentEnergiesCount
Number of different incident energies are available.
G4ParticleHPNames * ElementNames_
Pointer to G4NeutronHPNames Provides access to the list of element names included in Geant4...
void Divide(G4int Elements, T *To, T *Numerator, T *Denominator=NULL)
Divide an array by another.
G4double * ProbabilityRangeEnd
The value of the highest probability segment stored in this tree.
G4double powN(G4double x, G4int n) const
G4double G4SampleUniform(void)
Returns a double value evenly distributed in the range (0, 1].
CLHEP::Hep3Vector G4ThreeVector
void CheckAlphaSanity(void)
Checks to make sure that alpha overpopulation will not occur, which could result in an unsolvable zer...
G4Ions * LargestA_
Defines the largest Z particle in the field of trees.
G4double * GetYieldProbability(void)
Get the yield probability.
G4FissionProductYieldDist(G4int WhichIsotope, G4FFGEnumerations::MetaState WhichMetaState, G4FFGEnumerations::FissionCause WhichCause, G4FFGEnumerations::YieldType WhichYieldType, std::istringstream &dataStream)
Default constructor.
G4Ions * FindParticle(G4double RandomParticle)
Returns the G4Ions definitions pointer for the particle whose probability segment contains the (0...
G4Ions * Particle
Pointer to the G4Ions definition that contains the data for the isotope/isomer for this ProbabilityBr...
G4FPYSamplingOps * RandomEngine_
Pointer to the CLHEP library random engine.
G4double * YieldEnergies_
Energy values of each energy.
void G4SetTernaryProbability(G4double TernaryProbability)
Sets the probability of ternary fission.
void SetMomentum(const G4double x, const G4double y, const G4double z)
#define G4FFG_LOCATION__
G4FFG_LOCATION__ outputs the current location in the code.
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
static G4Alpha * Definition()
ProbabilityBranch * Right
Pointer to the branch that has a higher probability.
G4int Verbosity_
Verbosity level.
static const G4int SpontaneousNubar_[][3]
Evaluated nubar values for neutron induced fission.
G4DynamicParticleVector * G4GetFission(void)
Generates a fission event using default sampling and returns the pointer to that fission event...
G4String MakeDirectoryName(void)
Generates the directory location for the data file referenced by G4FissionProductYieldDist.
G4int GetPDGEncoding() const
G4ENDFTapeRead * ENDFData_
Name of the fission yield product data file that G4FissionProductYieldDist references.
void G4SetAlphaProduction(G4double WhatAlphaProduction)
Set the alpha production behavior for fission event generation.
void Initialize(std::istringstream &dataStream)
Initialize is a common function called by all constructors.
G4String MakeFileName(G4int Isotope, G4FFGEnumerations::MetaState MetaState)
Generates the appropriate file name for the isotope requested.
G4double RemainingEnergy_
Container for the energy remaining to be assigned in the fission generation.
YieldType
The two types of fission data available.
static const G4String theString[100]
G4double * ProbabilityRangeBottom
The upper limit of the probability segment, indexed by incident energy.
#define G4FFG_DATA_FUNCTIONENTER__
virtual void SortProbability(G4ENDFYieldDataContainer *YieldData)
Sorts information for a potential new particle into the correct tree.
void G4SetVerbosity(G4int WhatVerbosity)
Sets the verbosity levels.
G4double NubarWidth_
Width of the gaussian distribution that samples nubar for the isotope and incident neutron energy tha...
const G4String & GetParticleName() const
G4int RemainingA_
Counter for the number of nucleons available to the fission event.
G4Gamma * GammaDefinition_
Contains the g4ParticleDefinition pointer to a gamma particle.
G4double * ProbabilityRangeTop
The lower limit of the probability segment, indexed by incident energy.
FissionCause
Causes of fission.
G4int GetAtomicNumber() const
G4Ions * NeutronDefinition_
Contains the G4ParticleDefinition pointer to a neutron, cast as a G4Ion for compatibility.
Verbosity
These are the verbosity levels.
G4double TernaryProbability_
Sets the ternary fission probability.
static const G4int SpontaneousNubarWidth_[][2]
Recommended Gaussian widths for neutron induced fission.
virtual void GenerateNeutrons(std::vector< G4ReactionProduct * > *Neutrons)
Generate a linked chain of neutrons and return the pointer to the last neutron in the chain...
virtual void GenerateAlphas(std::vector< G4ReactionProduct * > *Alphas)
Generates a G4DynamicParticleVector with the fission alphas.
G4DynamicParticle * MakeG4DynamicParticle(G4ReactionProduct *)
Creates a G4DynamicParticle from an existing G4ReactionProduct.
const G4ParticleDefinition * GetDefinition() const
void SampleAlphaEnergies(std::vector< G4ReactionProduct * > *Alphas)
Sample the energy of the alpha particles.
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
virtual void MakeTrees(void)
Dynamically allocates and initializes the 'field' of 'trees' with the 'trunks'.
ProbabilityTree is the 'root' for each tree hierarchy, and contains pointer to the first ProbabilityB...
G4double Nubar_
Nubar for the isotope and incident neutron energy that G4FissionProductYieldDist references.
#define G4FFG_DATA_FUNCTIONLEAVE__
G4int YieldEnergyGroups_
Number of specific energy groups.
G4double AlphaProduction_
Controls whether alpha particles are emitted, and how many.
G4double * DataTotal_
A running total of all the probabilities.
G4FFGDefaultValues is a one-stop shop for storing the default values to variables that configure how ...
G4int BranchCount_
A run-time counter for the total number of branches stored.
static const G4int NeutronInducedNubar_[][3]
Evaluated nubar values for neutron induced fission.
void Multiply(G4int Elements, T *To, T *M1, T *M2=NULL)
Multiply two arrays together.
ProbabilityTree * Trees_
An array, or 'field', of the probability trees.
G4Ions * SmallestA_
Defines the smallest A particle in the field of trees.
std::vector< G4DynamicParticle * > G4DynamicParticleVector
const G4int Isotope_
Number in ZZZAAA format of the isotope that G4FissionProductYieldDist references. ...
G4Ions * GetParticleDefinition(G4int Product, G4FFGEnumerations::MetaState MetaState)
Returns the G4Ions definition pointer to the isotope defined by Product and MetaState.
G4Ions * FindParticleBranchSearch(ProbabilityBranch *Branch, G4double RandomParticle, G4int EnergyGroup1, G4int EnergyGroup2)
Returns the G4Ions definitions pointer for the particle whose probability segment contains the (0...
void G4SetEnergy(G4double WhatIncidentEnergy)
Sets the energy of the incident particle.
static const G4double MeanGammaEnergy
Default mean gamma energy for gamma sampling.
static const char ENDFFissionDataLocation[]
ENDF data tape location, reference against G4HPNEUTRONDATA.
void DeleteVectorOfPointers(std::vector< T > &Vector)
static const G4int Isotope
Default Isotope.
G4ENDFYieldDataContainer is a simple data storage class that handles the memory management internally...
G4int GetAtomicMass() const
static G4IonTable * GetIonTable()
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
void Renormalize(ProbabilityBranch *Branch)
Renormalizes the data in a ProbabilityTree.
G4int GetProduct(void)
Get the product.
G4int G4GetNumberOfEnergyGroups(void)
Returns the number of energy yield groups that were extracted from the ENDF tape file.
ProbabilityBranch * Trunk
First branch, or 'trunk', in the tree.
G4Ions * SmallestZ_
Defines the smallest Z particle in the field of trees.
G4double GetPDGMass() const
G4ENDFTapeRead is a class designed to read in data from unformatted ENDF data tapes for MT = 454 or M...
void SampleGammaEnergies(std::vector< G4ReactionProduct * > *Gammas)
Samples the energy of the gamma rays.
G4Ions * FindParticleExtrapolation(G4double RandomParticle, G4bool LowerEnergyGroupExists)
Returns the G4Ions definitions pointer for the particle whose probability segment contains the (0...
G4bool IsEnd
Variable to identify if this is the last tree in or not.
G4int BranchCount
Counter for the number of branches that this tree has.
void SetNubar(void)
Sets the nubar values for the isotope referenced by G4FissionProductYieldDistdefined from the data se...
G4double IncidentEnergy_
Kinetic energy, if any, of the incident particle in GeV.
G4Ions * FindParticleInterpolation(G4double RandomParticle, G4int LowerEnergyGroup)
Returns the G4Ions definitions pointer for the particle whose probability segment contains the (0...
virtual void ReadProbabilities(void)
Reads in the probability data from the data file.
const G4FFGEnumerations::FissionCause Cause_
The cause of fission: SPONTANEOUS or N_INDUCED.
static const G4int NeutronInducedNubarWidth_[][2]
Recommended Gaussian widths for neutron induced fission.
G4double * IncidentEnergies
The different energies available.
static G4Neutron * Definition()
ProbabilityBranch is a tree hierarchy storage method.
void Set(G4int Elements, T *To, T Value)
Set's all the values in an array to a constant.
void G4SetVerbosity(G4int WhatVerbosity)
Sets the verbosity levels.
G4ThreeVector GetMomentum() const
#define G4FFG_FUNCTIONLEAVE__
G4FPYSamplingOps performs all the uniform and Gaussian distribution sampling operations.
G4double * MaintainNormalizedData_
Variable for ensuring that the input data is normalized.
#define G4FFG_SPACING__
G4FFG_SPACING__ indents the debug messages according to the debugging depth.
void SampleNeutronEnergies(std::vector< G4ReactionProduct * > *Neutrons)
Sample the energy of the neutrons using the Watt fission spectrum.
#define G4FFG_RECURSIVE_FUNCTIONENTER__
G4int G4GetNumberOfFissionProducts(void)
Returns the number of fission products that were extracted from the ENDF tape file.
ProbabilityBranch * Left
Pointer to the branch that has a lower probability.
void BurnTree(ProbabilityBranch *Branch)
Recursively burns each branch in a probability tree.
G4IonTable * IonTable_
Pointer to G4IonTable All G4Ions are created using G4IonTable.
void G4SetVerbosity(G4int WhatVerbosity)
Sets the verbosity levels.
G4Ions * AlphaDefinition_
Contains the G4Ions pointer to an alpha particle.
G4ThreeVector G4ParticleMomentum
G4int TreeCount_
The number of trees in the field.
virtual G4Ions * GetFissionProduct(void)=0
Selects a fission product from the probability tree, limited by the number of nucleons available to t...
G4Ions * G4GetFissionProduct(void)
Selects a fission fragment at random from the probability tree and returns the G4Ions pointer...
#define G4FFG_FUNCTIONENTER__
G4String MakeIsotopeName(G4int Isotope, G4FFGEnumerations::MetaState MetaState)
Generates the unique name for an isotope/isomer defined by Isotope\ and MetaState in the following fo...
G4Ions * LargestZ_
Defines the largest Z particle in the field of trees.
void Add(G4int Elements, T *To, T *A1, T *A2=NULL)
Add two arrays together.
G4double * G4GetEnergyGroupValues(void)
Returns and array containing the values of each of the energy groups.
static G4Gamma * Definition()
G4int RemainingZ_
Counter for the number of protons available to the fission event.
G4int G4SampleIntegerGaussian(G4double Mean, G4double StdDev)
Returns an integer value taken from a Gaussian distribution.