From d12cf0441f0d17dbfabe23bb9530643f836ec6ea Mon Sep 17 00:00:00 2001 From: Roman Tolchenov Date: Thu, 5 Mar 2015 13:05:00 +0000 Subject: [PATCH] Re #11131. Unit tests for numeric parameter calculations. --- .../BackToBackExponential.h | 1 - .../src/PeakParametersNumeric.cpp | 6 + .../CurveFitting/test/ChebfunBaseTest.h | 2 - .../test/IPeakFunctionIntensityTest.h | 2 +- .../test/PeakParametersNumericTest.h | 283 ++++++++++++++++++ 5 files changed, 290 insertions(+), 4 deletions(-) create mode 100644 Code/Mantid/Framework/CurveFitting/test/PeakParametersNumericTest.h diff --git a/Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/BackToBackExponential.h b/Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/BackToBackExponential.h index bdf4e79ca415..0a29fef9099d 100644 --- a/Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/BackToBackExponential.h +++ b/Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/BackToBackExponential.h @@ -65,7 +65,6 @@ class DLLExport BackToBackExponential : public PeakParametersNumeric { virtual void functionDeriv1D(API::Jacobian *out, const double *xValues, const size_t nData); - //double operator()(double) const; std::pair getExtent() const; protected: diff --git a/Code/Mantid/Framework/CurveFitting/src/PeakParametersNumeric.cpp b/Code/Mantid/Framework/CurveFitting/src/PeakParametersNumeric.cpp index fbdaae0d728c..e596643ed063 100644 --- a/Code/Mantid/Framework/CurveFitting/src/PeakParametersNumeric.cpp +++ b/Code/Mantid/Framework/CurveFitting/src/PeakParametersNumeric.cpp @@ -95,6 +95,8 @@ double PeakParametersNumeric::fwhm() const { /// to defineWidthParameter(...). /// @param w :: New value for the width. void PeakParametersNumeric::setFwhm(const double w) { + const double c = centre(); + const double h = height(); double wOld = fwhm(); double factor = w / wOld; for (size_t i = 0; i < m_widthIndices.size(); ++i) { @@ -113,6 +115,10 @@ void PeakParametersNumeric::setFwhm(const double w) { } setParameter(index, value); } + // The width parameters can shift the centre and height, + // restore them. + setCentre(c); + setHeight(h); } /// Calculate function value for a single argument. diff --git a/Code/Mantid/Framework/CurveFitting/test/ChebfunBaseTest.h b/Code/Mantid/Framework/CurveFitting/test/ChebfunBaseTest.h index eafa66f878c2..6e09f3c68f77 100644 --- a/Code/Mantid/Framework/CurveFitting/test/ChebfunBaseTest.h +++ b/Code/Mantid/Framework/CurveFitting/test/ChebfunBaseTest.h @@ -6,8 +6,6 @@ #include "MantidCurveFitting/ChebfunBase.h" #include -#include "C:/Users/hqs74821/Work/Mantid_stuff/Testing/class/MyTestDef.h" - using namespace Mantid; using namespace Mantid::API; using namespace Mantid::CurveFitting; diff --git a/Code/Mantid/Framework/CurveFitting/test/IPeakFunctionIntensityTest.h b/Code/Mantid/Framework/CurveFitting/test/IPeakFunctionIntensityTest.h index 513f22227ae1..c75e7594d56b 100644 --- a/Code/Mantid/Framework/CurveFitting/test/IPeakFunctionIntensityTest.h +++ b/Code/Mantid/Framework/CurveFitting/test/IPeakFunctionIntensityTest.h @@ -75,7 +75,7 @@ class IPeakFunctionIntensityTest : public CxxTest::TestSuite { "), but intensity changed from " + DBL2STR(oldIntensity) + " to " + DBL2STR(newIntensity) + " (ratio " + DBL2STR(intensityRatio) + ").", - intensityRatio, heightRatio, 1e-10); + intensityRatio, heightRatio, 1e-8); } initialIntensities = newIntensities; diff --git a/Code/Mantid/Framework/CurveFitting/test/PeakParametersNumericTest.h b/Code/Mantid/Framework/CurveFitting/test/PeakParametersNumericTest.h new file mode 100644 index 000000000000..ed9d296f087e --- /dev/null +++ b/Code/Mantid/Framework/CurveFitting/test/PeakParametersNumericTest.h @@ -0,0 +1,283 @@ +#ifndef PEAKPARAMETERSNUMERICTEST_H_ +#define PEAKPARAMETERSNUMERICTEST_H_ + +#include + +#include "MantidCurveFitting/PeakParametersNumeric.h" +#include "MantidCurveFitting/BackToBackExponential.h" +#include + +using namespace Mantid; +using namespace Mantid::API; +using namespace Mantid::CurveFitting; + +class GaussFun : public PeakParametersNumeric { +public: + std::pair getExtent() const + { + double c = getParameter("c"); + double w = getTrueFwhm(); + return std::make_pair(c - 2*w, c + 2*w); + } + + virtual double getTrueFwhm() const = 0; +}; + + +class GaussLinearW : public GaussFun { +public: + /// Default constructor. + GaussLinearW() : GaussFun() { + declareParameter("h",1.0); + declareParameter("s",1.0); + declareParameter("c",0.0); + + defineCentreParameter("c"); + defineHeightParameter("h"); + defineWidthParameter("s",Linear); + } + + std::string name() const { return "GaussLinearW"; } + virtual const std::string category() const { return "Peak"; } + virtual void function1D(double *out, const double *xValues, + const size_t nData) const + { + double h = getParameter("h"); + double s = getParameter("s"); + double c = getParameter("c"); + + for(size_t i = 0; i < nData; ++i) + { + double tmp = (xValues[i] - c ) / s; + out[i] = h * exp(- tmp * tmp / 2 ); + } + } + + double getTrueFwhm() const + { + double s = getParameter("s"); + return 2 * sqrt(2.0 * log(2.0)) * s; + } + +protected: + virtual void functionLocal(double *, const double *, const size_t) const {} + virtual void functionDerivLocal(API::Jacobian *, const double *, + const size_t) {} + double expWidth() const; +}; + +class GaussInverseW : public GaussFun { +public: + /// Default constructor. + GaussInverseW() : GaussFun() { + declareParameter("h",1.0); + declareParameter("s",1.0); + declareParameter("c",0.0); + + defineCentreParameter("c"); + defineHeightParameter("h"); + defineWidthParameter("s",Inverse); + } + + std::string name() const { return "GaussInverseW"; } + virtual const std::string category() const { return "Peak"; } + virtual void function1D(double *out, const double *xValues, + const size_t nData) const + { + double h = getParameter("h"); + double s = getParameter("s"); + double c = getParameter("c"); + + for(size_t i = 0; i < nData; ++i) + { + double tmp = (xValues[i] - c ) * s; + out[i] = h * exp(- tmp * tmp / 2 ); + } + } + + double getTrueFwhm() const + { + double s = getParameter("s"); + return 2 * sqrt(2.0 * log(2.0)) / s; + } + +protected: + virtual void functionLocal(double *, const double *, const size_t) const {} + virtual void functionDerivLocal(API::Jacobian *, const double *, + const size_t) {} + double expWidth() const; +}; + +class GaussSquaredW : public GaussFun { +public: + /// Default constructor. + GaussSquaredW() : GaussFun() { + declareParameter("h",1.0); + declareParameter("s",1.0); + declareParameter("c",0.0); + + defineCentreParameter("c"); + defineHeightParameter("h"); + defineWidthParameter("s",Square); + } + + std::string name() const { return "GaussSquaredW"; } + virtual const std::string category() const { return "Peak"; } + virtual void function1D(double *out, const double *xValues, + const size_t nData) const + { + double h = getParameter("h"); + double s = getParameter("s"); + double c = getParameter("c"); + + for(size_t i = 0; i < nData; ++i) + { + double tmp = (xValues[i] - c ); + out[i] = h * exp(- tmp * tmp / 2 / s ); + } + } + + double getTrueFwhm() const + { + double s = getParameter("s"); + return 2 * sqrt(2.0 * log(2.0) * s); + } + +protected: + virtual void functionLocal(double *, const double *, const size_t) const {} + virtual void functionDerivLocal(API::Jacobian *, const double *, + const size_t) {} + double expWidth() const; +}; + + +class PeakParametersNumericTest : public CxxTest::TestSuite { +public: + + void test_GaussLinearW() + { + do_test_Gauss(GaussLinearW(),1e-7); + } + + void test_GaussInverseW() + { + do_test_Gauss(GaussInverseW(),1e-4); + } + + void test_GaussSquaredW() + { + do_test_Gauss(GaussSquaredW(),1e-7); + } + + void test_Back2Back() + { + BackToBackExponential fun; + fun.initialize(); + double tol = 1e-4; + TS_ASSERT_DELTA(fun.centre(), 0.0335, tol); + TS_ASSERT_DELTA(fun.height(), 2.0953, tol); + TS_ASSERT_DELTA(fun.fwhm(), 0.4027, tol); + + double dh; + double x = fun.centre() - tol; + fun.function1D(&dh,&x,1); + TS_ASSERT( dh < fun.height() ); + + x = fun.centre() + tol; + fun.function1D(&dh,&x,1); + TS_ASSERT( dh < fun.height() ); + + auto range = fun.getExtent(); + double left = 0.0; + double right = 0.0; + dh = fun.height() / 2; + for(double x = range.first; x < range.second; x += tol) + { + double y; + fun.function1D(&y,&x,1); + if (left == 0.0 && y >= dh ) + { + left = x; + } + if (left != 0.0 && right == 0.0 && y <= dh ) + { + right = x; + break; + } + } + TS_ASSERT_DELTA(right - left, fun.fwhm(), tol); + + fun.setCentre(0.0); + TS_ASSERT_DELTA(fun.centre(), 0.0, tol); + TS_ASSERT_DELTA(fun.height(), 2.0953, tol); + TS_ASSERT_DELTA(fun.fwhm(), 0.4027, tol); + + fun.setCentre(-1.0); + TS_ASSERT_DELTA(fun.centre(), -1.0, tol); + TS_ASSERT_DELTA(fun.height(), 2.0953, tol); + TS_ASSERT_DELTA(fun.fwhm(), 0.4027, tol); + + fun.setHeight(1.0); + TS_ASSERT_DELTA(fun.centre(), -1.0, tol); + TS_ASSERT_DELTA(fun.height(), 1.0, tol); + TS_ASSERT_DELTA(fun.fwhm(), 0.4027, tol); + std::cerr << fun.intensity() << std::endl; + + fun.setHeight(10.0); + TS_ASSERT_DELTA(fun.centre(), -1.0, tol); + TS_ASSERT_DELTA(fun.height(), 10.0, tol); + TS_ASSERT_DELTA(fun.fwhm(), 0.4027, tol); + std::cerr << fun.intensity() << std::endl; + + fun.setFwhm(1.0); + TS_ASSERT_DELTA(fun.centre(), -1.0, tol); + TS_ASSERT_DELTA(fun.height(), 10.0, tol); + TS_ASSERT_DELTA(fun.fwhm(), 1.0, tol); + + } + +private: + + void do_test_Gauss(GaussFun& fun, double tol) + { + //std::cerr << fun.centre() << ' ' << fun.height() << ' ' << fun.fwhm() - fun.getTrueFwhm() << std::endl; + TS_ASSERT_DELTA(fun.centre(), 0.0, tol); + TS_ASSERT_DELTA(fun.height(), 1.0, tol); + TS_ASSERT_DELTA(fun.fwhm(), fun.getTrueFwhm(), tol); + + fun.setHeight(2.1); + TS_ASSERT_DELTA(fun.centre(), 0.0, tol); + TS_ASSERT_DELTA(fun.height(), 2.1, tol); + TS_ASSERT_DELTA(fun.fwhm(), fun.getTrueFwhm(), tol); + + fun.setHeight(0.3); + TS_ASSERT_DELTA(fun.centre(), 0.0, tol); + TS_ASSERT_DELTA(fun.height(), 0.3, tol); + TS_ASSERT_DELTA(fun.fwhm(), fun.getTrueFwhm(), tol); + + fun.setCentre(1.3); + TS_ASSERT_DELTA(fun.centre(), 1.3, tol); + TS_ASSERT_DELTA(fun.height(), 0.3, tol); + TS_ASSERT_DELTA(fun.fwhm(), fun.getTrueFwhm(), tol); + + fun.setCentre(-1.3); + TS_ASSERT_DELTA(fun.centre(), -1.3, tol); + TS_ASSERT_DELTA(fun.height(), 0.3, tol); + TS_ASSERT_DELTA(fun.fwhm(), fun.getTrueFwhm(), tol); + + fun.setFwhm(2.0); + TS_ASSERT_DELTA(fun.centre(), -1.3, tol); + TS_ASSERT_DELTA(fun.height(), 0.3, tol); + TS_ASSERT_DELTA(fun.fwhm(), fun.getTrueFwhm(), tol); + TS_ASSERT_DELTA(fun.fwhm(), 2.0, tol); + + fun.setFwhm(0.001); + TS_ASSERT_DELTA(fun.centre(), -1.3, tol); + TS_ASSERT_DELTA(fun.height(), 0.3, tol); + TS_ASSERT_DELTA(fun.fwhm(), fun.getTrueFwhm(), tol); + TS_ASSERT_DELTA(fun.fwhm(), 0.001, tol); + } + +}; + +#endif /*PEAKPARAMETERSNUMERICTEST_H_*/