diff --git a/Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/CalculateChiSquared.h b/Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/CalculateChiSquared.h index 29f5386b2cea..43694114c164 100644 --- a/Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/CalculateChiSquared.h +++ b/Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/CalculateChiSquared.h @@ -49,4 +49,4 @@ class DLLExport CalculateChiSquared : public IFittingAlgorithm { } // namespace CurveFitting } // namespace Mantid -#endif /* MANTID_CURVEFITTING_CALCULATECHISQUARED_H_ */ \ No newline at end of file +#endif /* MANTID_CURVEFITTING_CALCULATECHISQUARED_H_ */ diff --git a/Code/Mantid/Framework/CurveFitting/src/CalculateChiSquared.cpp b/Code/Mantid/Framework/CurveFitting/src/CalculateChiSquared.cpp index 292b48bc8a4e..73a74b987b59 100644 --- a/Code/Mantid/Framework/CurveFitting/src/CalculateChiSquared.cpp +++ b/Code/Mantid/Framework/CurveFitting/src/CalculateChiSquared.cpp @@ -43,8 +43,31 @@ void CalculateChiSquared::initConcrete() { //---------------------------------------------------------------------------------------------- /// Execute the algorithm. void CalculateChiSquared::execConcrete() { + + // Function may need some preparation. + m_function->setUpForFit(); + + API::FunctionDomain_sptr domain; + API::FunctionValues_sptr values; + //m_domainCreator->ignoreInvalidData(getProperty("IgnoreInvalidData")); + m_domainCreator->createDomain(domain, values); + + // Do something with the function which may depend on workspace. + m_domainCreator->initFunction(m_function); + // Calculate function values. + m_function->function(*domain, *values); + + // Calculate the chi squared. + double chiSquared = 0.0; + for(size_t i = 0; i < values->size(); ++i) { + double tmp = values->getFitData(i) - values->getCalculated(i); + chiSquared += tmp * tmp; + } + + // Store the result. + setProperty("ChiSquared", chiSquared); } } // namespace CurveFitting -} // namespace Mantid \ No newline at end of file +} // namespace Mantid diff --git a/Code/Mantid/Framework/CurveFitting/test/CalculateChiSquaredTest.h b/Code/Mantid/Framework/CurveFitting/test/CalculateChiSquaredTest.h index 75765d79c0e3..39701361558a 100644 --- a/Code/Mantid/Framework/CurveFitting/test/CalculateChiSquaredTest.h +++ b/Code/Mantid/Framework/CurveFitting/test/CalculateChiSquaredTest.h @@ -4,41 +4,116 @@ #include #include "MantidCurveFitting/CalculateChiSquared.h" +#include "MantidAPI/FunctionDomain1D.h" #include "MantidAPI/FunctionFactory.h" +#include "MantidAPI/FunctionValues.h" +#include "MantidAPI/IFunction1D.h" #include "MantidKernel/EmptyValues.h" +#include + using Mantid::CurveFitting::CalculateChiSquared; using namespace Mantid; using namespace Mantid::API; -class CalculateChiSquaredTest : public CxxTest::TestSuite -{ +class CalculateChiSquaredTest : public CxxTest::TestSuite { public: // This pair of boilerplate methods prevent the suite being created statically // This means the constructor isn't called when running other tests - static CalculateChiSquaredTest *createSuite() { return new CalculateChiSquaredTest(); } - static void destroySuite( CalculateChiSquaredTest *suite ) { delete suite; } - + static CalculateChiSquaredTest *createSuite() { + return new CalculateChiSquaredTest(); + } + static void destroySuite(CalculateChiSquaredTest *suite) { delete suite; } - void test_Init() - { + void test_Init() { CalculateChiSquared alg; - TS_ASSERT_THROWS_NOTHING( alg.initialize() ) - TS_ASSERT( alg.isInitialized() ) + TS_ASSERT_THROWS_NOTHING(alg.initialize()) + TS_ASSERT(alg.isInitialized()) } - void test_exec() - { + void test_1D_empty_defaults() { Tester tester; tester.set1DFunction(); tester.set1DSpectrumEmpty(); tester.runAlgorithm(); + tester.check1DSpectrum(); + TS_ASSERT_DELTA(tester.chiSquared, 20338.0, 1.0); + } + + void test_1D_empty_all_x_range() { + Tester tester; + tester.set1DFunction(); + tester.set1DSpectrumEmpty(); + tester.setXRange_All(); + tester.runAlgorithm(); + tester.check1DSpectrum(); + TS_ASSERT_DELTA(tester.chiSquared, 20338.0, 1.0); + } + + void test_1D_empty_smaller() { + Tester tester; + tester.set1DFunction(); + tester.set1DSpectrumEmpty(); + tester.setXRange_Smaller(); + tester.runAlgorithm(); + tester.check1DSpectrum(); + TS_ASSERT_DELTA(tester.chiSquared, 1189.0, 1.0); + } + + void test_1D_empty_smaller1() { + Tester tester; + tester.set1DFunction(); + tester.set1DSpectrumEmpty(); + tester.setXRange_Smaller_bin_boundaries(); + tester.runAlgorithm(); + tester.check1DSpectrum(); + TS_ASSERT_DELTA(tester.chiSquared, 1189.0, 1.0); } - -private: - class Tester{ + void test_1D_values() { + Tester tester; + tester.set1DFunction(); + tester.set1DSpectrumValues(); + tester.runAlgorithm(); + tester.check1DSpectrum(); + TS_ASSERT_DELTA(tester.chiSquared, 1655.0, 1.0); + } + void test_1D_values_smaller() { + Tester tester; + tester.set1DFunction(); + tester.set1DSpectrumValues(); + tester.setXRange_Smaller(); + tester.runAlgorithm(); + tester.check1DSpectrum(); + TS_ASSERT_DELTA(tester.chiSquared, 153.0, 1.0); + std::cerr << tester.chiSquared << std::endl; + } + + void test_1D_values_point_data() { + Tester tester(3,10,false); + tester.set1DFunction(); + tester.set1DSpectrumValues(); + tester.runAlgorithm(); + tester.check1DSpectrum(); + TS_ASSERT_DELTA(tester.chiSquared, 307.0, 1.0); + } + + void test_1D_workspace_index() { + Tester tester; + tester.set1DFunction(); + tester.set1DSpectrumValues(5); + tester.setXRange_Smaller(); + tester.setWorkspaceIndex(); + tester.runAlgorithm(); + tester.check1DSpectrum(); + TS_ASSERT_DELTA(tester.chiSquared, 153.0, 1.0); + std::cerr << tester.chiSquared << std::endl; + } + +private: + class Tester { + public: // input parameters size_t nParams; size_t nData; @@ -47,6 +122,7 @@ class CalculateChiSquaredTest : public CxxTest::TestSuite double xMax; std::vector xBins; std::vector xValues; + std::vector yValues; // values for algorithm input properties IFunction_sptr function; @@ -55,60 +131,134 @@ class CalculateChiSquaredTest : public CxxTest::TestSuite double StartX; double EndX; - // algorithm output - double chiSquared; - - void makeXValues(){ + void makeXValues() { size_t dlt = isHisto ? 1 : 0; xBins.resize(nData + dlt); double dx = (xMax - xMin) / double(xBins.size() - 1); - for(size_t i = 0; i < xBins.size(); ++i) - { + for (size_t i = 0; i < xBins.size(); ++i) { xBins[i] = xMin + i * dx; } - if (isHisto){ + if (isHisto) { xValues.resize(nData); - std::transform(xBins.begin(),xBins.end()-1, xValues.begin(), std::bind2nd(std::plus(),dx/2)); + std::transform(xBins.begin(), xBins.end() - 1, xValues.begin(), + std::bind2nd(std::plus(), dx / 2)); } else { xValues = xBins; } } + void setDefaultXRange() { + if (StartX == EMPTY_DBL()) + StartX = xMin; + if (EndX == EMPTY_DBL()) { + EndX = xMax; + } else { + auto ix = std::upper_bound(xBins.begin(), xBins.end(), EndX); + if (ix != xBins.end()) EndX = *ix; + else + EndX = xMax; + } + } + public: + // algorithm output + double chiSquared; + Tester(size_t np = 3, size_t nd = 10, bool histo = true) : nParams(np), nData(nd), isHisto(histo), workspaceIndex(0), xMin(-10), - xMax(10), StartX(EMPTY_DBL()), EndX(EMPTY_DBL()) {} - - void runAlgorithm(){ + xMax(10), StartX(EMPTY_DBL()), EndX(EMPTY_DBL()) { + makeXValues(); + } + void runAlgorithm() { CalculateChiSquared alg; - TS_ASSERT_THROWS_NOTHING( alg.initialize() ) - TS_ASSERT( alg.isInitialized() ) - TS_ASSERT_THROWS_NOTHING( alg.setProperty("Function", function) ); - TS_ASSERT_THROWS_NOTHING( alg.setProperty("InputWorkspace", workspace) ); - TS_ASSERT_THROWS_NOTHING( alg.execute(); ); - TS_ASSERT( alg.isExecuted() ); - + TS_ASSERT_THROWS_NOTHING(alg.initialize()) + TS_ASSERT(alg.isInitialized()) + TS_ASSERT_THROWS_NOTHING(alg.setProperty("Function", function)); + TS_ASSERT_THROWS_NOTHING(alg.setProperty("InputWorkspace", workspace)); + if (dynamic_cast(function.get())) { + TS_ASSERT_THROWS_NOTHING(alg.setProperty("StartX", StartX)); + TS_ASSERT_THROWS_NOTHING(alg.setProperty("EndX", EndX)); + } + TS_ASSERT_THROWS_NOTHING(alg.execute()); + TS_ASSERT(alg.isExecuted()); chiSquared = alg.getProperty("ChiSquared"); } - void set1DFunction(){ - const std::string fun = "name=UserFunction,Formula=a+b*x+c*x^2,a=1,b=1,c=1"; + void setXRange_All() { + StartX = xMin; + EndX = xMax; + } + + void setXRange_Smaller_bin_boundaries() { + StartX = xBins[3]; + EndX = xBins[7]; + } + + void setXRange_Smaller() { + StartX = xBins[3] - 0.3; + EndX = xBins[7] + 0.7; + } + + void setWorkspaceIndex() { + workspaceIndex = 3; + auto ws = dynamic_cast(workspace.get()); + if (ws) + yValues = ws->readY(workspaceIndex); + else + TS_FAIL("Not a matrix workspace"); + } + + void set1DFunction() { + const std::string fun = + "name=UserFunction,Formula=a+b*x+c*x^2,a=1,b=1,c=1"; function = FunctionFactory::Instance().createInitialized(fun); TS_ASSERT_EQUALS(function->nParams(), nParams); } - void set1DSpectrumEmpty(){ - makeXValues(); + void set1DSpectrumEmpty() { + const size_t nSpec = 1; - auto space = WorkspaceFactory::Instance().create("Workspace2D", nSpec, nData+1, nData); - space->dataX(0).assign(xBins.begin(),xBins.end()); + size_t dn = isHisto ? 1 : 0; + auto space = WorkspaceFactory::Instance().create("Workspace2D", nSpec, + nData + dn, nData); + space->dataX(0).assign(xBins.begin(), xBins.end()); workspace = space; + yValues = space->readY(0); } - }; + void set1DSpectrumValues(const size_t nSpec = 1) { + size_t dn = isHisto ? 1 : 0; + auto space = WorkspaceFactory::Instance().create("Workspace2D", nSpec, + nData + dn, nData); + space->dataX(0).assign(xBins.begin(), xBins.end()); + for (size_t spec = 0; spec < nSpec; ++spec) { + for (size_t i = 0; i < nData; ++i) { + const double x = space->readX(0)[i]; + space->dataY(spec)[i] = (1.1 + 0.1 * double(spec)) * (1.0 + x + x * x); + } + } + workspace = space; + yValues = space->readY(0); + } + void check1DSpectrum() { + setDefaultXRange(); + double sum2 = 0.0; + for (size_t i = 0; i < xValues.size(); ++i) { + const double xValue = xValues[i]; + if (xValue >= StartX && xValue <= EndX) { + FunctionDomain1DVector x(xValue); + FunctionValues y(x); + function->function(x, y); + const double tmp = yValues[i] - y[0]; + sum2 += tmp * tmp; + } + } + TS_ASSERT_DIFFERS(sum2, 0); + TS_ASSERT_DELTA(sum2, chiSquared, 1e-10); + } + }; }; - -#endif /* MANTID_CURVEFITTING_CALCULATECHISQUAREDTEST_H_ */ \ No newline at end of file +#endif /* MANTID_CURVEFITTING_CALCULATECHISQUAREDTEST_H_ */