G4FPYSamplingOps Class Reference

#include <G4FPYSamplingOps.hh>

Public Member Functions
	G4FPYSamplingOps ()
	G4FPYSamplingOps (G4int Verbosity)
G4int	G4SampleIntegerGaussian (G4double Mean, G4double StdDev)
G4int	G4SampleIntegerGaussian (G4double Mean, G4double StdDev, G4FFGEnumerations::GaussianRange Range)
G4double	G4SampleGaussian (G4double Mean, G4double StdDev)
G4double	G4SampleGaussian (G4double Mean, G4double StdDev, G4FFGEnumerations::GaussianRange Range)
G4double	G4SampleUniform ()
G4double	G4SampleUniform (G4double Lower, G4double Upper)
G4double	G4SampleWatt (G4int WhatIsotope, G4FFGEnumerations::FissionCause WhatCause, G4double WhatEnergy)
void	G4SetVerbosity (G4int WhatVerbosity)
	~G4FPYSamplingOps ()

Protected Member Functions
void	Initialize ()
G4bool	CheckAndSetParameters ()
void	EvaluateWattConstants ()
G4double	SampleGaussian ()
void	ShiftParameters (G4FFGEnumerations::GaussianReturnType Type)

Protected Attributes
G4double	Mean_
G4double	StdDev_
G4ShiftedGaussian *	ShiftedGaussianValues_
G4int	Verbosity_
WattSpectrumConstants *	WattConstants_
CLHEP::HepRandomEngine *	RandomEngine_
G4bool	NextGaussianIsStoredInMemory_
G4double	GaussianOne_
G4double	GaussianTwo_
G4double	Tolerance_

Detailed Description

G4FPYSamplingOps performs all the uniform and Gaussian distribution sampling operations

Definition at line 51 of file G4FPYSamplingOps.hh.

Constructor & Destructor Documentation

G4FPYSamplingOps() [1/2]

G4FPYSamplingOps::G4FPYSamplingOps

(

)

Default constructor

• Usage: No arguments required
• Notes:

Definition at line 49 of file G4FPYSamplingOps.cc.

{
  // Set the default verbosity
  Verbosity_ = G4FFGDefaultValues::Verbosity;
 
  // Initialize the class
  Initialize();
}

G4FPYSamplingOps() [2/2]

G4FPYSamplingOps::G4FPYSamplingOps

(

G4int

Verbosity

)

Overloaded constructor

• Usage:
  • Verbosity: Verbosity level
• Notes:

Definition at line 58 of file G4FPYSamplingOps.cc.

{
  // Set the default verbosity
  Verbosity_ = Verbosity;
 
  // Initialize the class
  Initialize();
}

~G4FPYSamplingOps()

G4FPYSamplingOps::~G4FPYSamplingOps

(

)

Default deconstructor.

Definition at line 609 of file G4FPYSamplingOps.cc.

{
  G4FFG_FUNCTIONENTER__
 
  delete ShiftedGaussianValues_;
  delete WattConstants_;
 
  G4FFG_FUNCTIONLEAVE__
}

Member Function Documentation

CheckAndSetParameters()

G4bool G4FPYSamplingOps::CheckAndSetParameters ( )

protected

Check to see if the user requested parameters have already been calculated. If they have, it recalls the stored parameters and sets them as the current values.

Definition at line 308 of file G4FPYSamplingOps.cc.

{
  G4FFG_SAMPLING_FUNCTIONENTER__
 
  G4double ShiftedMean = ShiftedGaussianValues_->G4FindShiftedMean(Mean_, StdDev_);
  if (ShiftedMean == 0) {
    G4FFG_SAMPLING_FUNCTIONLEAVE__
    return FALSE;
  }
 
  Mean_ = ShiftedMean;
 
  G4FFG_SAMPLING_FUNCTIONLEAVE__
  return TRUE;
}

Referenced by ShiftParameters().

EvaluateWattConstants()

void G4FPYSamplingOps::EvaluateWattConstants ( )

protected

Evaluates the constants that are required for the Watt fission spectrum sampling.

Definition at line 324 of file G4FPYSamplingOps.cc.

{
  G4FFG_SAMPLING_FUNCTIONENTER__
 
  G4double A, K;
  A = K = 0;
  // Use the default values if IsotopeIndex is not reset
  G4int IsotopeIndex = 0;
 
  if (WattConstants_->Cause == G4FFGEnumerations::SPONTANEOUS) {
    // Determine if the isotope requested exists in the lookup tables
    for (G4int i = 0; SpontaneousWattIsotopesIndex[i] != -1; i++) {
      if (SpontaneousWattIsotopesIndex[i] == WattConstants_->Product) {
        IsotopeIndex = i;
 
        break;
      }
    }
 
    // Get A and B
    A = SpontaneousWattConstants[IsotopeIndex][0];
    WattConstants_->B = SpontaneousWattConstants[IsotopeIndex][1];
  }
  else if (WattConstants_->Cause == G4FFGEnumerations::NEUTRON_INDUCED) {
    // Determine if the isotope requested exists in the lookup tables
    for (G4int i = 0; NeutronInducedWattIsotopesIndex[i] != -1; i++) {
      if (NeutronInducedWattIsotopesIndex[i] == WattConstants_->Product) {
        IsotopeIndex = i;
        break;
      }
    }
 
    // Determine the Watt fission spectrum constants based on the energy of
    // the incident neutron
    if (WattConstants_->Energy == G4FFGDefaultValues::ThermalNeutronEnergy) {
      A = NeutronInducedWattConstants[IsotopeIndex][0][0];
      WattConstants_->B = NeutronInducedWattConstants[IsotopeIndex][0][1];
    }
    else if (WattConstants_->Energy > 14.0 * CLHEP::MeV) {
      G4Exception("G4FPYSamplingOps::G4SampleWatt()",
                  "Incident neutron energy above 14 MeV requested.", JustWarning,
                  "Using Watt fission constants for 14 Mev.");
 
      A = NeutronInducedWattConstants[IsotopeIndex][2][0];
      WattConstants_->B = NeutronInducedWattConstants[IsotopeIndex][2][1];
    }
    else {
      G4int EnergyIndex = 0;
      G4double EnergyDifference = 0;
      G4double RangeDifference, ConstantDifference;
 
      for (G4int i = 1; IncidentEnergyBins[i] != -1; i++) {
        if (WattConstants_->Energy <= IncidentEnergyBins[i]) {
          EnergyIndex = i;
          EnergyDifference = IncidentEnergyBins[EnergyIndex] - WattConstants_->Energy;
          if (EnergyDifference != 0) {
            std::ostringstream Temp;
            Temp << "Incident neutron energy of ";
            Temp << WattConstants_->Energy << " MeV is not ";
            Temp << "explicitly listed in the data tables";
            //                        G4Exception("G4FPYSamplingOps::G4SampleWatt()",
            //                                    Temp.str().c_str(),
            //                                    JustWarning,
            //                                    "Using linear interpolation.");
          }
          break;
        }
      }
 
      RangeDifference = IncidentEnergyBins[EnergyIndex] - IncidentEnergyBins[EnergyIndex - 1];
 
      // Interpolate the value for 'a'
      ConstantDifference = NeutronInducedWattConstants[IsotopeIndex][EnergyIndex][0]
                           - NeutronInducedWattConstants[IsotopeIndex][EnergyIndex - 1][0];
      A = (EnergyDifference / RangeDifference) * ConstantDifference
          + NeutronInducedWattConstants[IsotopeIndex][EnergyIndex - 1][0];
 
      // Interpolate the value for 'b'
      ConstantDifference = NeutronInducedWattConstants[IsotopeIndex][EnergyIndex][1]
                           - NeutronInducedWattConstants[IsotopeIndex][EnergyIndex - 1][1];
      WattConstants_->B = (EnergyDifference / RangeDifference) * ConstantDifference
                          + NeutronInducedWattConstants[IsotopeIndex][EnergyIndex - 1][1];
    }
  }
  else {
    // Produce an error since an unsupported fission type was requested and
    // abort the current run
    G4String Temp = "Watt fission spectra data not available for ";
    if (WattConstants_->Cause == G4FFGEnumerations::PROTON_INDUCED) {
      Temp += "proton induced fission.";
    }
    else if (WattConstants_->Cause == G4FFGEnumerations::GAMMA_INDUCED) {
      Temp += "gamma induced fission.";
    }
    else {
      Temp += "!Warning! unknown cause.";
    }
    G4Exception("G4FPYSamplingOps::G4SampleWatt()", Temp, RunMustBeAborted,
                "Fission events will not be sampled in this run.");
  }
 
  // Calculate the constants
  K = 1 + (WattConstants_->B / (8.0 * A));
  WattConstants_->L = (K + G4Pow::GetInstance()->powA((K * K - 1), 0.5)) / A;
  WattConstants_->M = A * WattConstants_->L - 1;
 
  G4FFG_SAMPLING_FUNCTIONLEAVE__
}

Referenced by G4SampleWatt().

G4SampleGaussian() [1/2]

G4double G4FPYSamplingOps::G4SampleGaussian	(	G4double	Mean,
		G4double	StdDev )

Returns a double value taken from a Gaussian distribution about Mean and with a standard deviation of StdDev.

• Usage:
  • Mean: Mean about which the Gaussian distribution will be sampled
  • StdDev: Standard deviation of the Gaussian distribution. 68.3% of the values will lie within the first standard deviation, 95.4% within the second standard deviation, etc...
• Notes:

Definition at line 89 of file G4FPYSamplingOps.cc.

{
  G4FFG_SAMPLING_FUNCTIONENTER__
 
  // Determine if the parameters have changed
  G4bool ParametersChanged = (Mean_ != Mean || StdDev_ != StdDev);
  if (static_cast<int>(ParametersChanged) == TRUE) {
    // Set the new values if the parameters have changed
    NextGaussianIsStoredInMemory_ = FALSE;
 
    Mean_ = Mean;
    StdDev_ = StdDev;
  }
 
  G4double Sample = SampleGaussian();
 
  G4FFG_SAMPLING_FUNCTIONLEAVE__
  return Sample;
}

Referenced by G4FissionProductYieldDist::G4GetFission(), G4SampleGaussian(), G4FissionProductYieldDist::SampleAlphaEnergies(), and G4FissionProductYieldDist::SampleGammaEnergies().

G4SampleGaussian() [2/2]

G4double G4FPYSamplingOps::G4SampleGaussian	(	G4double	Mean,
		G4double	StdDev,
		G4FFGEnumerations::GaussianRange	Range )

Returns a double value taken from a Gaussian distribution about Mean and with a standard deviation of StdDev.

• Usage:
  • Mean: Mean about which the Gaussian distribution will be sampled
  • StdDev: Standard deviation of the Gaussian distribution. 68.3% of the values will lie within the first standard deviation, 95.4% within the second standard deviation, etc...
  • Range: POSITIVE or ALL
• Notes:

Definition at line 109 of file G4FPYSamplingOps.cc.

{
  G4FFG_SAMPLING_FUNCTIONENTER__
 
  if (Range == G4FFGEnumerations::ALL) {
    // Call the overloaded function
    G4double Sample = G4SampleGaussian(Mean, StdDev);
 
    G4FFG_SAMPLING_FUNCTIONLEAVE__
    return Sample;
  }
 
  // Determine if the parameters have changed
  G4bool ParametersChanged = (Mean_ != Mean || StdDev_ != StdDev);
  if (static_cast<int>(ParametersChanged) == TRUE) {
    if (Mean <= 0) {
      std::ostringstream Temp;
      Temp << "Mean value of " << Mean << " out of range";
      G4Exception("G4FPYGaussianOps::G4SampleIntegerGaussian()", Temp.str().c_str(), JustWarning,
                  "A value of '0' will be used instead.");
 
      G4FFG_SAMPLING_FUNCTIONLEAVE__
      return 0;
    }
 
    // Set the new values if the parameters have changed and then perform
    // the shift
    Mean_ = Mean;
    StdDev_ = StdDev;
 
    ShiftParameters(G4FFGEnumerations::DOUBLE);
  }
 
  // Sample the Gaussian distribution
  G4double Rand;
  do {
    Rand = SampleGaussian();
  } while (Rand < 0);  // Loop checking, 11.05.2015, T. Koi
 
  G4FFG_SAMPLING_FUNCTIONLEAVE__
  return Rand;
}

G4SampleIntegerGaussian() [1/2]

G4int G4FPYSamplingOps::G4SampleIntegerGaussian	(	G4double	Mean,
		G4double	StdDev )

Returns an integer value taken from a Gaussian distribution. This overloaded version assumes that the range is not restricted to positive values only.

• Usage:
  • Mean: Mean about which the Gaussian distribution will be sampled
  • StdDev: Standard deviation of the Gaussian distribution. 68.3% of the values will lie within the first standard deviation, 95.4% within the second standard deviation, etc...
• Notes:

Definition at line 153 of file G4FPYSamplingOps.cc.

{
  G4FFG_SAMPLING_FUNCTIONENTER__
 
  // Determine if the parameters have changed
  G4bool ParametersChanged = (Mean_ != Mean || StdDev_ != StdDev);
  if (static_cast<int>(ParametersChanged) == TRUE) {
    // Set the new values if the parameters have changed
    NextGaussianIsStoredInMemory_ = FALSE;
 
    Mean_ = Mean;
    StdDev_ = StdDev;
  }
 
  // Return the sample integer value
  auto Sample = (G4int)(std::floor(SampleGaussian()));
 
  G4FFG_SAMPLING_FUNCTIONLEAVE__
  return Sample;
}

Referenced by G4SampleIntegerGaussian(), G4FissionProductYieldDist::GenerateAlphas(), and G4FissionProductYieldDist::GenerateNeutrons().

G4SampleIntegerGaussian() [2/2]

G4int G4FPYSamplingOps::G4SampleIntegerGaussian	(	G4double	Mean,
		G4double	StdDev,
		G4FFGEnumerations::GaussianRange	Range )

Returns an integer value taken from a Gaussian distribution about Mean and with a standard deviation of StdDev.

• Usage:
  • Mean: Mean about which the Gaussian distribution will be sampled
  • StdDev: Standard deviation of the Gaussian distribution. 68.3% of the values will lie within the first standard deviation, 95.4% within the second standard deviation, etc...
  • Range: POSITIVE or ALL
• Notes:

Definition at line 174 of file G4FPYSamplingOps.cc.

{
  G4FFG_SAMPLING_FUNCTIONENTER__
 
  if (Range == G4FFGEnumerations::ALL) {
    // Call the overloaded function
    G4int Sample = G4SampleIntegerGaussian(Mean, StdDev);
 
    G4FFG_SAMPLING_FUNCTIONLEAVE__
    return Sample;
  }
  // ParameterShift() locks up if the mean is less than 1.
  std::ostringstream Temp;
  if (Mean < 1) {
    //    Temp << "Mean value of " << Mean << " out of range";
    //    G4Exception("G4FPYGaussianOps::G4SampleIntegerGaussian()",
    //                Temp.str().c_str(),
    //                JustWarning,
    //                "A value of '0' will be used instead.");
 
    //    return 0;
  }
 
  if (Mean / StdDev < 2) {
    // Temp << "Non-ideal conditions:\tMean:" << Mean << "\tStdDev: "
    //         << StdDev;
    // G4Exception("G4FPYGaussianOps::G4SampleIntegerGaussian()",
    //             Temp.str().c_str(),
    //             JustWarning,
    //             "Results not guaranteed: Consider tightening the standard deviation");
  }
 
  // Determine if the parameters have changed
  G4bool ParametersChanged = (Mean_ != Mean || StdDev_ != StdDev);
  if (static_cast<int>(ParametersChanged) == TRUE) {
    // Set the new values if the parameters have changed and then perform
    // the shift
    Mean_ = Mean;
    StdDev_ = StdDev;
 
    ShiftParameters(G4FFGEnumerations::INT);
  }
 
  G4int RandInt;
  // Sample the Gaussian distribution - only non-negative values are
  // accepted
  do {
    RandInt = (G4int)floor(SampleGaussian());
  } while (RandInt < 0);  // Loop checking, 11.05.2015, T. Koi
 
  G4FFG_SAMPLING_FUNCTIONLEAVE__
  return RandInt;
}

G4SampleUniform() [1/2]

G4double G4FPYSamplingOps::G4SampleUniform

(

)

Returns a double value evenly distributed in the range (0, 1].

• Usage: No arguments required
• Notes:

Definition at line 229 of file G4FPYSamplingOps.cc.

{
  G4FFG_SAMPLING_FUNCTIONENTER__
 
  G4double Sample = RandomEngine_->flat();
 
  G4FFG_SAMPLING_FUNCTIONLEAVE__
  return Sample;
}

Referenced by G4FissionProductYieldDist::G4GetFission(), G4FissionProductYieldDist::G4GetFissionProduct(), G4SampleUniform(), G4SampleWatt(), G4FissionProductYieldDist::GenerateAlphas(), G4FPYBiasedLightFragmentDist::GetFissionProduct(), G4FPYNormalFragmentDist::GetFissionProduct(), and SampleGaussian().

G4SampleUniform() [2/2]

G4double G4FPYSamplingOps::G4SampleUniform	(	G4double	Lower,
		G4double	Upper )

Returns a double value evenly distributed in the range (Lower, Upper].

• Usage:
  • Lower: Lower bounds of the distribution
  • Upper: Upper bounds of the distribution
• Notes:

Definition at line 239 of file G4FPYSamplingOps.cc.

{
  G4FFG_SAMPLING_FUNCTIONENTER__
 
  // Calculate the difference
  G4double Difference = Upper - Lower;
 
  // Scale appropriately and return the value
  G4double Sample = G4SampleUniform() * Difference + Lower;
 
  G4FFG_SAMPLING_FUNCTIONLEAVE__
  return Sample;
}

G4SampleWatt()

G4double G4FPYSamplingOps::G4SampleWatt	(	G4int	WhatIsotope,
		G4FFGEnumerations::FissionCause	WhatCause,
		G4double	WhatEnergy )

Samples the Watt fission spectrum for the selected isotope, using an algorithm adopted from Ref. 1

• Usage:
  • WhatIsotope: The isotope that is to be sampled
  • WhatCause: The cause of the isotope to be sampled
  • WhatEnergy: The energy, in MeV of the incident particle
• Notes:
  • All variables needed for this function are grouped together in WattConstants_.

Definition at line 253 of file G4FPYSamplingOps.cc.

{
  G4FFG_SAMPLING_FUNCTIONENTER__
 
  // Determine if the parameters are different
  // TK modified 131108
  // if(WattConstants_->Product != WhatIsotope
  if (WattConstants_->Product != WhatIsotope / 10 || WattConstants_->Cause != WhatCause
      || WattConstants_->Energy != WhatEnergy)
  {
    // If the parameters are different or have not yet been defined then we
    // need to re-evaluate the constants
    // TK modified 131108
    // WattConstants_->Product = WhatIsotope;
    WattConstants_->Product = WhatIsotope / 10;
    WattConstants_->Cause = WhatCause;
    WattConstants_->Energy = WhatEnergy;
 
    EvaluateWattConstants();
  }
 
  G4double X = -G4Log(G4SampleUniform());
  G4double Y = -G4Log(G4SampleUniform());
  G4int icounter = 0;
  G4int icounter_max = 1024;
  while (G4Pow::GetInstance()->powN(Y - WattConstants_->M * (X + 1), 2)
         > WattConstants_->B * WattConstants_->L * X)  // Loop checking, 11.05.2015, T. Koi
  {
    icounter++;
    if (icounter > icounter_max) {
      G4cout << "Loop-counter exceeded the threshold value at " << __LINE__ << "th line of "
             << __FILE__ << "." << G4endl;
      break;
    }
    X = -G4Log(G4SampleUniform());
    Y = -G4Log(G4SampleUniform());
  }
 
  G4FFG_SAMPLING_FUNCTIONLEAVE__
  return WattConstants_->L * X;
}

Referenced by G4FissionProductYieldDist::SampleNeutronEnergies().

G4SetVerbosity()

void G4FPYSamplingOps::G4SetVerbosity

(

G4int

WhatVerbosity

)

Sets the verbosity levels

• Usage:
  • WhichVerbosity: Combination of levels
• Notes:
  • SILENT: All verbose output is repressed
  • UPDATES: Only high-level internal changes are reported
  • DAUGHTER_INFO: Displays information about daughter product sampling
  • NEUTRON_INFO: Displays information about neutron sampling
  • GAMMA_INFO: Displays information about gamma sampling
  • ALPHA_INFO: Displays information about alpha sampling
  • MOMENTUM_INFO: Displays information about momentum balancing
  • EXTRAPOLATION_INTERPOLATION_INFO: Displays information about any data extrapolation or interpolation that occurs
  • DEBUG: Reports program flow as it steps through functions
  • PRINT_ALL: Displays any and all output

Definition at line 297 of file G4FPYSamplingOps.cc.

{
  G4FFG_SAMPLING_FUNCTIONENTER__
 
  Verbosity_ = Verbosity;
 
  ShiftedGaussianValues_->G4SetVerbosity(Verbosity_);
 
  G4FFG_SAMPLING_FUNCTIONLEAVE__
}

Referenced by G4FissionProductYieldDist::G4SetVerbosity().

Initialize()

void G4FPYSamplingOps::Initialize ( )

protected

Initialize is a common function called by all constructors.

Definition at line 67 of file G4FPYSamplingOps.cc.

{
  G4FFG_FUNCTIONENTER__
 
  // Get the pointer to the random number generator
  // RandomEngine_ = CLHEP::HepRandom::getTheEngine();
  RandomEngine_ = G4Random::getTheEngine();
 
  // Initialize the data members
  ShiftedGaussianValues_ = new G4ShiftedGaussian;
  Mean_ = 0;
  StdDev_ = 0;
  NextGaussianIsStoredInMemory_ = FALSE;
  GaussianOne_ = 0;
  GaussianTwo_ = 0;
  Tolerance_ = 0.000001;
  WattConstants_ = new WattSpectrumConstants;
  WattConstants_->Product = 0;
 
  G4FFG_FUNCTIONLEAVE__
}

Referenced by G4FPYSamplingOps(), and G4FPYSamplingOps().

SampleGaussian()

G4double G4FPYSamplingOps::SampleGaussian ( )

protected

Samples a Gaussian distribution defined by the internal class variables NewMean_ and NewStdDev_.

Definition at line 433 of file G4FPYSamplingOps.cc.

{
  G4FFG_SAMPLING_FUNCTIONENTER__
 
  if (static_cast<int>(NextGaussianIsStoredInMemory_) == TRUE) {
    NextGaussianIsStoredInMemory_ = FALSE;
 
    G4FFG_SAMPLING_FUNCTIONLEAVE__
    return GaussianTwo_;
  }
 
  // Define the variables to be used
  G4double Radius;
  G4double MappingFactor;
 
  // Sample from the unit circle (21.4% rejection probability)
  do {
    GaussianOne_ = 2.0 * G4SampleUniform() - 1.0;
    GaussianTwo_ = 2.0 * G4SampleUniform() - 1.0;
    Radius = GaussianOne_ * GaussianOne_ + GaussianTwo_ * GaussianTwo_;
  } while (Radius > 1.0);  // Loop checking, 11.05.2015, T. Koi
 
  // Translate the values to Gaussian space
  MappingFactor = std::sqrt(-2.0 * G4Log(Radius) / Radius) * StdDev_;
  GaussianOne_ = Mean_ + GaussianOne_ * MappingFactor;
  GaussianTwo_ = Mean_ + GaussianTwo_ * MappingFactor;
 
  // Set the flag that a value is now stored in memory
  NextGaussianIsStoredInMemory_ = TRUE;
 
  G4FFG_SAMPLING_FUNCTIONLEAVE__
  return GaussianOne_;
}

Referenced by G4SampleGaussian(), G4SampleGaussian(), G4SampleIntegerGaussian(), and G4SampleIntegerGaussian().

ShiftParameters()

void G4FPYSamplingOps::ShiftParameters ( G4FFGEnumerations::GaussianReturnType Type )

protected

Sets the mean and standard deviation of the Gaussian distribution sampled by this class when POSITIVE values are requested. ShiftMean() performs two different operations based on the requested data type.

• INTEGER: Iteratively searches for an adjusted mean that produces the same result as the mean requested by the implementor. In this instance the standard deviation is not adjusted.
• DOUBLE: Adjusts the standard deviation of the Gaussian distribution so that the first seven standard deviations occur are all positive. The chance that a negative value will result using this method is 2.56E^-12

Definition at line 467 of file G4FPYSamplingOps.cc.

{
  G4FFG_SAMPLING_FUNCTIONENTER__
 
  // Set the flag that any second value stored is no longer valid
  NextGaussianIsStoredInMemory_ = FALSE;
 
  // Check if the requested parameters have already been calculated
  if (static_cast<int>(CheckAndSetParameters()) == TRUE) {
    G4FFG_SAMPLING_FUNCTIONLEAVE__
    return;
  }
 
  // If the requested type is INT, then perform an iterative refinement of the
  // mean that is going to be sampled
  if (Type == G4FFGEnumerations::INT) {
    // Return if the mean is greater than 7 standard deviations away from 0
    // since there is essentially 0 probability that a sampled number will
    // be less than 0
    if (Mean_ > 7 * StdDev_) {
      G4FFG_SAMPLING_FUNCTIONLEAVE__
      return;
    }
    // Variables that contain the area and weighted area information for
    // calculating the statistical mean of the Gaussian distribution when
    // converted to a step function
    G4double ErfContainer, AdjustedErfContainer, Container;
 
    // Variables that store lower and upper bounds information
    G4double LowErf, HighErf;
 
    // Values for determining the convergence of the solution
    G4double AdjMean = Mean_;
    G4double Delta = 1.0;
    G4bool HalfDelta = false;
    G4bool ToleranceCheck = false;
 
    // Normalization constant for use with erf()
    const G4double Normalization = StdDev_ * std::sqrt(2.0);
 
    // Determine the upper limit of the estimates
    const auto UpperLimit = (G4int)std::ceil(Mean_ + 9 * StdDev_);
 
    // Calculate the integral of the Gaussian distribution
 
    G4int icounter = 0;
    G4int icounter_max = 1024;
    do {
      icounter++;
      if (icounter > icounter_max) {
        G4cout << "Loop-counter exceeded the threshold value at " << __LINE__ << "th line of "
               << __FILE__ << "." << G4endl;
        break;
      }
      // Set the variables
      ErfContainer = 0;
      AdjustedErfContainer = 0;
 
      // Calculate the area and weighted area
      for (G4int i = 0; i <= UpperLimit; i++) {
        // Determine the lower and upper bounds
        LowErf = ((AdjMean - i) / Normalization);
        HighErf = ((AdjMean - (i + 1.0)) / Normalization);
 
        // Correct the bounds for how they lie on the x-axis with
        // respect to the mean
        if (LowErf <= 0) {
          LowErf *= -1;
          HighErf *= -1;
// Container = (erf(HighErf) - erf(LowErf))/2.0;
#if defined WIN32 - VC
          Container = (CLHEP::HepStat::erf(HighErf) - CLHEP::HepStat::erf(LowErf)) / 2.0;
#else
          Container = (erf(HighErf) - erf(LowErf)) / 2.0;
#endif
        }
        else if (HighErf < 0) {
          HighErf *= -1;
 
// Container = (erf(HighErf) + erf(LowErf))/2.0;
#if defined WIN32 - VC
          Container = (CLHEP::HepStat::erf(HighErf) + CLHEP::HepStat::erf(LowErf)) / 2.0;
#else
          Container = (erf(HighErf) + erf(LowErf)) / 2.0;
#endif
        }
        else {
// Container = (erf(LowErf) - erf(HighErf))/2.0;
#if defined WIN32 - VC
          Container = (CLHEP::HepStat::erf(LowErf) - CLHEP::HepStat::erf(HighErf)) / 2.0;
#else
          Container = (erf(LowErf) - erf(HighErf)) / 2.0;
#endif
        }
 
#if defined WIN32 - VC
        // TK modified to avoid problem caused by QNAN
        if (Container != Container) Container = 0;
#endif
 
        // Add up the weighted area
        ErfContainer += Container;
        AdjustedErfContainer += Container * i;
      }
 
      // Calculate the statistical mean
      Container = AdjustedErfContainer / ErfContainer;
 
      // Is it close enough to what we want?
      ToleranceCheck = (std::fabs(Mean_ - Container) < Tolerance_);
      if (static_cast<int>(ToleranceCheck) == TRUE) {
        break;
      }
 
      // Determine the step size
      if (static_cast<int>(HalfDelta) == TRUE) {
        Delta /= 2;
      }
 
      // Step in the appropriate direction
      if (Container > Mean_) {
        AdjMean -= Delta;
      }
      else {
        HalfDelta = TRUE;
        AdjMean += Delta;
      }
    } while (TRUE);  // Loop checking, 11.05.2015, T. Koi
 
    ShiftedGaussianValues_->G4InsertShiftedMean(AdjMean, Mean_, StdDev_);
    Mean_ = AdjMean;
  }
  else if (Mean_ / 7 < StdDev_) {
    // If the requested type is double, then just re-define the standard
    // deviation appropriately - chances are approximately 2.56E-12 that
    // the value will be negative using this sampling scheme
    StdDev_ = Mean_ / 7;
  }
 
  G4FFG_SAMPLING_FUNCTIONLEAVE__
}

Referenced by G4SampleGaussian(), and G4SampleIntegerGaussian().

Member Data Documentation

GaussianOne_

G4double G4FPYSamplingOps::GaussianOne_

protected

Contains the first of the two paired random numbers from the Gaussian distribution sampling.

Definition at line 198 of file G4FPYSamplingOps.hh.

Referenced by Initialize(), and SampleGaussian().

GaussianTwo_

G4double G4FPYSamplingOps::GaussianTwo_

protected

Contains the second of the two paired random numbers from the Gaussian distribution sampling.

Definition at line 202 of file G4FPYSamplingOps.hh.

Referenced by Initialize(), and SampleGaussian().

Mean_

G4double G4FPYSamplingOps::Mean_

protected

Mean for sampling a Gaussian distribution

Definition at line 172 of file G4FPYSamplingOps.hh.

Referenced by CheckAndSetParameters(), G4SampleGaussian(), G4SampleGaussian(), G4SampleIntegerGaussian(), G4SampleIntegerGaussian(), Initialize(), SampleGaussian(), and ShiftParameters().

NextGaussianIsStoredInMemory_

G4bool G4FPYSamplingOps::NextGaussianIsStoredInMemory_

protected

Declares whether the second paired random number has been already returned.

Definition at line 194 of file G4FPYSamplingOps.hh.

Referenced by G4SampleGaussian(), G4SampleIntegerGaussian(), Initialize(), SampleGaussian(), and ShiftParameters().

RandomEngine_

CLHEP::HepRandomEngine* G4FPYSamplingOps::RandomEngine_

protected

Pointer to the CLHEP random number generator.

Definition at line 188 of file G4FPYSamplingOps.hh.

Referenced by G4SampleUniform(), and Initialize().

ShiftedGaussianValues_

G4ShiftedGaussian* G4FPYSamplingOps::ShiftedGaussianValues_

protected

Structure chain that contains the all the previous values used for sampling a Gaussian distribution

Definition at line 178 of file G4FPYSamplingOps.hh.

Referenced by CheckAndSetParameters(), G4SetVerbosity(), Initialize(), ShiftParameters(), and ~G4FPYSamplingOps().

StdDev_

G4double G4FPYSamplingOps::StdDev_

protected

Standard deviation for sampling a GaussianDistribution

Definition at line 174 of file G4FPYSamplingOps.hh.

Referenced by CheckAndSetParameters(), G4SampleGaussian(), G4SampleGaussian(), G4SampleIntegerGaussian(), G4SampleIntegerGaussian(), Initialize(), SampleGaussian(), and ShiftParameters().

Tolerance_

G4double G4FPYSamplingOps::Tolerance_

protected

Defines the tolerance that ShiftParameters() must match.

Definition at line 204 of file G4FPYSamplingOps.hh.

Referenced by Initialize(), and ShiftParameters().

Verbosity_

G4int G4FPYSamplingOps::Verbosity_

protected

Verbosity level

Definition at line 180 of file G4FPYSamplingOps.hh.

Referenced by G4FPYSamplingOps(), G4FPYSamplingOps(), and G4SetVerbosity().

WattConstants_

WattSpectrumConstants* G4FPYSamplingOps::WattConstants_

protected

Structure that contains the values for sampling the Watt fission spectrum

Definition at line 184 of file G4FPYSamplingOps.hh.

Referenced by EvaluateWattConstants(), G4SampleWatt(), Initialize(), and ~G4FPYSamplingOps().

---

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

• G4FPYSamplingOps.hh
• G4FPYSamplingOps.cc