Commit 91481fe2 authored by Hahn, Steven's avatar Hahn, Steven
Browse files

When possible use math constants.

parent 58d30f52
......@@ -30,10 +30,10 @@ public:
void setCentre(const double c) override { setParameter("Center", c); }
double fwhm() const override {
return getParameter("Sigma") * (2.0 * sqrt(2.0 * log(2.0)));
return getParameter("Sigma") * (2.0 * sqrt(2.0 * M_LN2));
}
void setFwhm(const double w) override {
setParameter("Sigma", w / (2.0 * sqrt(2.0 * log(2.0))));
setParameter("Sigma", w / (2.0 * sqrt(2.0 * M_LN2)));
}
double height() const override { return getParameter("Height"); }
......@@ -80,8 +80,8 @@ private:
}
double getGaussianAnalyticalInfiniteIntegral(IPeakFunction_sptr gaussian) {
return gaussian->height() * gaussian->fwhm() /
(2.0 * sqrt(2.0 * log(2.0))) * sqrt(2.0 * M_PI);
return gaussian->height() * gaussian->fwhm() / (2.0 * sqrt(M_LN2)) *
sqrt(M_PI);
}
public:
......@@ -169,7 +169,7 @@ public:
* -integral from -3 to 3 should give approx. 0.997
*/
IPeakFunction_sptr gaussian =
getGaussian(0.0, 2.0 * sqrt(2.0 * log(2.0)), 1.0 / sqrt(2.0 * M_PI));
getGaussian(0.0, 2.0 * sqrt(2.0 * M_LN2), 1.0 / sqrt(2.0 * M_PI));
PeakFunctionIntegrator integrator(1e-10);
IntegrationResult rOneSigma = integrator.integrate(*gaussian, -1.0, 1.0);
......
......@@ -300,7 +300,7 @@ void GetAllEi::exec() {
// preprocess first monitors peaks;
g_log.debug() << "*Looking for real energy peaks on first monitor\n";
findBinRanges(monitorWS->readX(0), monitorWS->readY(0), guess_ei,
this->m_min_Eresolution / (2 * std::sqrt(2 * std::log(2.))),
this->m_min_Eresolution / (2. * std::sqrt(2. * M_LN2)),
irange_min, irange_max, guessValid);
// remove invalid guess values
......@@ -311,7 +311,7 @@ void GetAllEi::exec() {
if (!this->getProperty("IgnoreSecondMonitor")) {
g_log.debug() << "*Looking for real energy peaks on second monitor\n";
findBinRanges(monitorWS->readX(1), monitorWS->readY(1), guess_ei,
this->m_min_Eresolution / (2 * std::sqrt(2 * std::log(2.))),
this->m_min_Eresolution / (2. * std::sqrt(2. * M_LN2)),
irange1_min, irange1_max, guessValid);
removeInvalidValues<double>(guessValid, guess_ei);
removeInvalidValues<size_t>(guessValid, irange_min);
......@@ -471,7 +471,7 @@ bool GetAllEi::peakGuess(const API::MatrixWorkspace_sptr &inputWS, size_t index,
double &peakTwoSigma) {
// calculate sigma from half-width parameters
double maxSigma = Ei * m_min_Eresolution / (2 * std::sqrt(2 * std::log(2.)));
double maxSigma = Ei * m_min_Eresolution / (2. * std::sqrt(2. * M_LN2));
double sMin(std::numeric_limits<double>::max());
double sMax(-sMin);
......@@ -621,8 +621,8 @@ bool GetAllEi::findMonitorPeak(const API::MatrixWorkspace_sptr &inputWS,
double &position, double &height,
double &twoSigma) {
// calculate sigma from half-width parameters
double maxSigma = Ei * m_min_Eresolution / (2 * std::sqrt(2 * std::log(2.)));
double minSigma = Ei * m_max_Eresolution / (2 * std::sqrt(2 * std::log(2.)));
double maxSigma = Ei * m_min_Eresolution / (2. * std::sqrt(2. * M_LN2));
double minSigma = Ei * m_max_Eresolution / (2. * std::sqrt(2. * M_LN2));
//--------------------------------------------------------------------
double peak1Pos, peak1TwoSigma, peak1Height;
if (!peakGuess(inputWS, 0, Ei, monsRangeMin, monsRangeMax, peak1Pos,
......
......@@ -26,7 +26,7 @@ namespace Mantid {
namespace Algorithms {
namespace {
/// Factor to convert full width half max to sigma for calculations of I/sigma.
const double FWHM_TO_SIGMA = 2.0 * sqrt(2.0 * std::log(2.0));
const double FWHM_TO_SIGMA = 2.0 * sqrt(2.0 * M_LN2);
const double BAD_OFFSET(1000.); // mark things that didn't work with this
//--------------------------------------------------------------------------------------------
......@@ -1087,7 +1087,7 @@ void GetDetOffsetsMultiPeaks::generatePeaksList(
vec_widthDivPos.push_back(widthdevpos);
// g_log.debug() << " h:" << height << " c:" << centre << " w:" <<
// (width/(2.*std::sqrt(2.*std::log(2.))))
// (width/(2.*std::sqrt(2.*M_LN2)))
// << " b:" << background << " chisq:" << chi2 << "\n";
// Add peak to vectors
......
......@@ -249,7 +249,7 @@ IFunction_sptr GetDetectorOffsets::createFunction(const double peakHeight,
peak->setHeight(peakHeight);
peak->setCentre(peakLoc);
const double sigma(10.0);
peak->setFwhm(2.0 * std::sqrt(2.0 * std::log(2.0)) * sigma);
peak->setFwhm(2.0 * std::sqrt(2.0 * M_LN2) * sigma);
auto fitFunc = new CompositeFunction(); // Takes ownership of the functions
fitFunc->addFunction(background);
......
......@@ -53,7 +53,7 @@ public:
for (size_t i = 0; i < Y.size(); ++i) {
TS_ASSERT_EQUALS(X[i], static_cast<double>(i) / 2.0)
TS_ASSERT_EQUALS(Y[i], 1)
TS_ASSERT_EQUALS(E[i], sqrt(2.0) / 2.0)
TS_ASSERT_EQUALS(E[i], M_SQRT2 / 2.0)
}
TS_ASSERT(!output->isDistribution())
......
......@@ -55,7 +55,7 @@ public:
for (size_t i = 0; i < Y.size(); ++i) {
TS_ASSERT_EQUALS(X[i], static_cast<double>(i) / 2.0)
TS_ASSERT_EQUALS(Y[i], 4)
TS_ASSERT_EQUALS(E[i], sqrt(2.0) / 0.5)
TS_ASSERT_EQUALS(E[i], M_SQRT2 / 0.5)
}
TS_ASSERT(output->isDistribution())
}
......
......@@ -47,7 +47,7 @@ public:
inWS->dataY(0)[2] = 0.0;
inWS->dataE(0)[2] = 0.0;
inWS->dataY(0)[3] = 2.0;
inWS->dataE(0)[3] = sqrt(2.0);
inWS->dataE(0)[3] = M_SQRT2;
inWS->dataY(0)[4] = 10000.0;
inWS->dataE(0)[4] = 100.0;
......
......@@ -88,7 +88,7 @@ public:
// {
// // Now the data. Y and E unchanged
// TS_ASSERT_EQUALS(yValues[j], 2.0);
// TS_ASSERT_EQUALS(eValues[j], sqrt(2.0));
// TS_ASSERT_EQUALS(eValues[j], M_SQRT2);
// // X data originally was 0->10 in steps of 1. Now it should be the
// centre of each bin which is
......
......@@ -78,7 +78,7 @@ public:
for (int j = 0; j < numBins; ++j) {
// Now the data. Y and E unchanged
TS_ASSERT_EQUALS(yValues[j], 2.0);
TS_ASSERT_EQUALS(eValues[j], sqrt(2.0));
TS_ASSERT_EQUALS(eValues[j], M_SQRT2);
// X data originally was 0->10 in steps of 1. Now it should be the
// centre of each bin which is
// 1.0 away from the last centre
......@@ -128,7 +128,7 @@ public:
for (size_t j = 0; j < numBins; ++j) {
// Now the data. Y and E unchanged
TS_ASSERT_EQUALS(yValues[j], 2.0);
TS_ASSERT_EQUALS(eValues[j], sqrt(2.0));
TS_ASSERT_EQUALS(eValues[j], M_SQRT2);
// X data originally was 0->10 in steps of 1. Now it should be the
// centre of each bin which is
// 1.0 away from the last centre
......
......@@ -633,7 +633,7 @@ private:
TS_ASSERT(ws.hasDx(0));
TS_ASSERT_EQUALS(ws.readDx(0)[0], 0.0);
TS_ASSERT_EQUALS(ws.readDx(0)[1], 1.0);
TS_ASSERT_EQUALS(ws.readDx(0)[2], sqrt(2.0));
TS_ASSERT_EQUALS(ws.readDx(0)[2], M_SQRT2);
TS_ASSERT_EQUALS(ws.readDx(0)[3], sqrt(3.0));
// Check that the length of x and dx is the same
auto x = ws.readX(0);
......@@ -644,7 +644,7 @@ private:
TS_ASSERT(ws.hasDx(0));
TS_ASSERT_EQUALS(ws.readDx(0)[0], 0.0 + 1.0);
TS_ASSERT_EQUALS(ws.readDx(0)[1], 1.0 + 1.0);
TS_ASSERT_EQUALS(ws.readDx(0)[2], sqrt(2.0) + 1.0);
TS_ASSERT_EQUALS(ws.readDx(0)[2], M_SQRT2 + 1.0);
TS_ASSERT_EQUALS(ws.readDx(0)[3], sqrt(3.0) + 1.0);
} else {
TSM_ASSERT("Should never reach here", false);
......
......@@ -72,7 +72,7 @@ public:
TS_ASSERT_EQUALS(outputWS->readE(i)[j], 0.0);
} else {
TS_ASSERT_EQUALS(outputWS->readY(i)[j], 2.0);
TS_ASSERT_DELTA(outputWS->readE(i)[j], sqrt(2.0), 0.0001);
TS_ASSERT_DELTA(outputWS->readE(i)[j], M_SQRT2, 0.0001);
}
TS_ASSERT_EQUALS(outputWS->readX(i)[j], j);
}
......
......@@ -130,7 +130,7 @@ public:
void testExec() {
AnalysisDataService::Instance().add(
"normIn", WorkspaceCreationHelper::Create2DWorkspaceBinned(10, 3, 1));
doTest("normIn", "normOut", 1.0, sqrt(2.0) / 2.0);
doTest("normIn", "normOut", 1.0, 0.5 * M_SQRT2);
AnalysisDataService::Instance().remove("normIn");
AnalysisDataService::Instance().remove("normOut");
}
......@@ -268,7 +268,7 @@ public:
void testExec_InPlace() {
AnalysisDataService::Instance().add(
"normIn", WorkspaceCreationHelper::Create2DWorkspaceBinned(10, 3, 1));
doTest("normIn", "normIn", 1.0, sqrt(2.0) / 2.0);
doTest("normIn", "normIn", 1.0, 0.5 * M_SQRT2);
AnalysisDataService::Instance().remove("normIn");
}
......@@ -279,7 +279,7 @@ public:
EventWorkspace_const_sptr outputEvent;
outputEvent = boost::dynamic_pointer_cast<const EventWorkspace>(
doTest("normInEvent", "normOutEvent", 1.0, sqrt(2.0) / 2.0));
doTest("normInEvent", "normOutEvent", 1.0, 0.5 * M_SQRT2));
// Output is an event workspace
TS_ASSERT(outputEvent);
......@@ -294,7 +294,7 @@ public:
EventWorkspace_const_sptr outputEvent;
outputEvent = boost::dynamic_pointer_cast<const EventWorkspace>(
doTest("normInEvent", "normInEvent", 1.0, sqrt(2.0) / 2.0));
doTest("normInEvent", "normInEvent", 1.0, 0.5 * M_SQRT2));
// Output is an event workspace
TS_ASSERT(outputEvent);
......
......@@ -47,7 +47,7 @@ public:
return;
// Check the results
TS_ASSERT_DELTA(ws->readY(0)[0], sqrt(2.0), 1e-5);
TS_ASSERT_DELTA(ws->readY(0)[0], M_SQRT2, 1e-5);
TS_ASSERT_DELTA(ws->readE(0)[0], 0.0, 1e-5);
// Remove workspace from the data service.
......
......@@ -62,7 +62,7 @@ public:
TS_ASSERT_DELTA(result->readE(0)[1176], 104.841321, 0.000001)
// Now one where first is zero
TS_ASSERT_EQUALS(result->readY(0)[2], 2.0)
TS_ASSERT_EQUALS(result->readE(0)[2], std::sqrt(2.0))
TS_ASSERT_EQUALS(result->readE(0)[2], M_SQRT2)
// And one where second is zero
TS_ASSERT_EQUALS(result->readY(0)[113], 97.0)
TS_ASSERT_EQUALS(result->readE(0)[113], std::sqrt(97.0))
......
......@@ -123,8 +123,8 @@ void TOFExtinction::exec() {
sigfsq_ys = getSigFsqr(EgLaueI, cell, wl, twoth, tbar, fsq, sigfsq);
} else if (cType.compare("Type I Gaussian") == 0) {
// Apply correction to fsq with Type-I BCG for testing
double EgLaueI = std::sqrt(2.0) *
getEgLaue(Eg, twoth, wl, divBeam, betaBeam) * 2.0 / 3.0;
double EgLaueI =
M_SQRT2 * getEgLaue(Eg, twoth, wl, divBeam, betaBeam) * 2.0 / 3.0;
double Xqt = getXqt(EgLaueI, cell, wl, twoth, tbar, fsq);
y_corr = getGaussian(Xqt, twoth);
sigfsq_ys = getSigFsqr(EgLaueI, cell, wl, twoth, tbar, fsq, sigfsq);
......
......@@ -226,7 +226,7 @@ void Bk2BkExpConvPV::geneatePeak(double *out, const double *xValues,
void Bk2BkExpConvPV::calHandEta(double sigma2, double gamma, double &H,
double &eta) const {
// 1. Calculate H
double H_G = sqrt(8.0 * sigma2 * log(2.0));
double H_G = sqrt(8.0 * sigma2 * M_LN2);
double H_L = gamma;
double temp1 = std::pow(H_L, 5) + 0.07842 * H_G * std::pow(H_L, 4) +
......
......@@ -662,8 +662,8 @@ void calculateEigesystem(DoubleFortranVector &eigenvalues,
// -------------------------------------------------------------------
ComplexFortranMatrix bex(1, 1, -1, 1);
bex(1, 1) = -(bext(1) - i * bext(2)) / sqrt(2.0);
bex(1, -1) = (bext(1) + i * bext(2)) / sqrt(2.0);
bex(1, 1) = -(bext(1) - i * bext(2)) * M_SQRT1_2;
bex(1, -1) = (bext(1) + i * bext(2)) * M_SQRT1_2;
bex(1, 0) = bext(3);
// calculates Bex(1,q) for a canted moment
......@@ -683,8 +683,8 @@ void calculateEigesystem(DoubleFortranVector &eigenvalues,
}
ComplexType rbextp, rbextm, rbextz;
rbextp = rbex(1, -1) * sqrt(2.0);
rbextm = rbex(1, 1) * (-sqrt(2.0));
rbextp = rbex(1, -1) * M_SQRT2;
rbextm = rbex(1, 1) * (-M_SQRT2);
rbextz = rbex(1, 0);
auto facmol = 2 * (gj - 1) * myb;
......
......@@ -77,7 +77,7 @@ double Gaussian::activeParameter(size_t i) const {
double Gaussian::centre() const { return getParameter("PeakCentre"); }
double Gaussian::height() const { return getParameter("Height"); }
double Gaussian::fwhm() const {
return 2.0 * sqrt(2.0 * std::log(2.0)) * getParameter("Sigma");
return 2.0 * sqrt(2.0 * M_LN2) * getParameter("Sigma");
}
double Gaussian::intensity() const {
auto sigma = getParameter("Sigma");
......@@ -96,7 +96,7 @@ double Gaussian::intensity() const {
void Gaussian::setCentre(const double c) { setParameter("PeakCentre", c); }
void Gaussian::setHeight(const double h) { setParameter("Height", h); }
void Gaussian::setFwhm(const double w) {
setParameter("Sigma", w / (2.0 * sqrt(2.0 * std::log(2.0))));
setParameter("Sigma", w / (2.0 * sqrt(2.0 * M_LN2)));
}
void Gaussian::setIntensity(const double i) {
m_intensityCache = i;
......
......@@ -293,7 +293,7 @@ void GramCharlierComptonProfile::addMassProfile(
const double denom = ((std::pow(2.0, static_cast<int>(npoly))) * factorial);
for (int j = 0; j < NFINE_Y; ++j) {
const double y = m_yfine[j] / std::sqrt(2.) / wg;
const double y = m_yfine[j] * M_SQRT1_2 / wg;
const double hermiteI = Math::hermitePoly(npoly, y);
result[j] += ampNorm * std::exp(-y * y) * hermiteI * hermiteCoeff / denom;
}
......@@ -316,7 +316,7 @@ void GramCharlierComptonProfile::addFSETerm(std::vector<double> &lhs) const {
kfse *= getParameter("C_0");
for (int j = 0; j < NFINE_Y; ++j) {
const double y = m_yfine[j] / std::sqrt(2.) / wg;
const double y = m_yfine[j] * M_SQRT1_2 / wg;
const double he3 = Math::hermitePoly(3, y);
lhs[j] += ampNorm * std::exp(-y * y) * he3 * (kfse / m_qfine[j]);
}
......
......@@ -398,7 +398,7 @@ void NeutronBk2BkExpConvPVoigt::function1D(double *out, const double *xValues,
void NeutronBk2BkExpConvPVoigt::calHandEta(double sigma2, double gamma,
double &H, double &eta) const {
// 1. Calculate H
double H_G = sqrt(8.0 * sigma2 * log(2.0));
double H_G = sqrt(8.0 * sigma2 * M_LN2);
double H_L = gamma;
double temp1 = std::pow(H_L, 5) + 0.07842 * H_G * std::pow(H_L, 4) +
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment