Skip to content

Commit

Permalink
Re #11617. Division by DOF
Browse files Browse the repository at this point in the history
  • Loading branch information
@mantid-roman
mantid-roman committed May 14, 2015
1 parent 2db38c1 commit 8f220b0
Show file tree
Hide file tree
Showing 2 changed files with 83 additions and 3 deletions.
22 changes: 22 additions & 0 deletions Code/Mantid/Framework/CurveFitting/src/CalculateChiSquared.cpp
View file
Expand Up @@ -37,6 +37,9 @@ const std::string CalculateChiSquared::summary() const {
//----------------------------------------------------------------------------------------------
/// Initialize the algorithm's properties.
void CalculateChiSquared::initConcrete() {
declareProperty("DivideByDOF", false, "Flag to divide the chi^2 by the "
"number of degrees of freedom (NofData "
"- nOfParams).");
declareProperty("ChiSquared", 0.0, "Output value of chi squared.");
}

Expand Down Expand Up @@ -65,6 +68,25 @@ void CalculateChiSquared::execConcrete() {
chiSquared += tmp * tmp;
}
}
g_log.debug() << "Chi squared " << chiSquared << std::endl;

// Divide by the DOF
bool divideByDOF = getProperty("DivideByDOF");
if (divideByDOF) {
// Get the number of free fitting parameters
size_t nParams = 0;
for(size_t i = 0; i < m_function->nParams(); ++i) {
if ( !m_function->isFixed(i) ) nParams += 1;
}
double dof = static_cast<double>(domain->size()) - static_cast<double>(nParams);
g_log.debug() << "DOF " << dof << std::endl;
if (dof <= 0.0) {
dof = 1.0;
g_log.warning() << "DOF has a non-positive value, changing to 1,0." << std::endl;
}
chiSquared /= dof;
g_log.debug() << "Chi squared / DOF " << chiSquared << std::endl;
}
// Store the result.
setProperty("ChiSquared", chiSquared);
}
Expand Down
64 changes: 61 additions & 3 deletions Code/Mantid/Framework/CurveFitting/test/CalculateChiSquaredTest.h
View file
Expand Up @@ -129,6 +129,46 @@ class CalculateChiSquaredTest : public CxxTest::TestSuite {
TS_ASSERT_DELTA(tester.chiSquared, 1450.39, 1.0);
}

void test_1D_values_divide_by_dof() {
Tester tester;
tester.set1DFunction();
tester.set1DSpectrumValues();
tester.setDivideByDOF();
tester.runAlgorithm();
tester.check1DSpectrum();
TS_ASSERT_DELTA(tester.chiSquared, 236.5, 1.0);
}

void test_1D_values_divide_by_dof_fixed_params() {
Tester tester(1);
tester.set1DFunction();
tester.set1DSpectrumValues();
tester.setDivideByDOF();
tester.runAlgorithm();
tester.check1DSpectrum();
TS_ASSERT_DELTA(tester.chiSquared, 184.0, 1.0);
}

void test_1D_values_divide_by_dof_zero() {
Tester tester(3,3);
tester.set1DFunction();
tester.set1DSpectrumValues();
tester.setDivideByDOF();
tester.runAlgorithm();
tester.check1DSpectrum();
TS_ASSERT_DELTA(tester.chiSquared, 5069.0, 1.0);
}

void test_1D_values_divide_by_dof_negative() {
Tester tester(3,2);
tester.set1DFunction();
tester.set1DSpectrumValues();
tester.setDivideByDOF();
tester.runAlgorithm();
tester.check1DSpectrum();
TS_ASSERT_DELTA(tester.chiSquared, 7151.0, 1.0);
}

private:
class Tester {
// input parameters
Expand All @@ -147,6 +187,7 @@ class CalculateChiSquaredTest : public CxxTest::TestSuite {
double StartX;
double EndX;
bool ignoreInvalidData;
bool divideByDOF;

void makeXValues() {
size_t dlt = isHisto ? 1 : 0;
Expand Down Expand Up @@ -191,7 +232,7 @@ class CalculateChiSquaredTest : public CxxTest::TestSuite {
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()),
ignoreInvalidData(false), chiSquared(-1), isExecuted(false) {
ignoreInvalidData(false), divideByDOF(false), chiSquared(-1), isExecuted(false) {
makeXValues();
}

Expand All @@ -202,6 +243,7 @@ class CalculateChiSquaredTest : public CxxTest::TestSuite {
TS_ASSERT_THROWS_NOTHING(alg.setProperty("Function", function));
TS_ASSERT_THROWS_NOTHING(alg.setProperty("InputWorkspace", workspace));
TS_ASSERT_THROWS_NOTHING(alg.setProperty("IgnoreInvalidData", ignoreInvalidData));
TS_ASSERT_THROWS_NOTHING(alg.setProperty("DivideByDOF", divideByDOF));
if (dynamic_cast<IFunction1D *>(function.get())) {
TS_ASSERT_THROWS_NOTHING(alg.setProperty("WorkspaceIndex", workspaceIndex));
TS_ASSERT_THROWS_NOTHING(alg.setProperty("StartX", StartX));
Expand Down Expand Up @@ -236,11 +278,22 @@ class CalculateChiSquaredTest : public CxxTest::TestSuite {
ignoreInvalidData = true;
}

void setDivideByDOF() {
divideByDOF = true;
}

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);
if (nParams < function->nParams()) {
const size_t dn = function->nParams() - nParams;
for(size_t i = 0; i < dn; ++i) {
function->fix(i);
}
} else {
TS_ASSERT_EQUALS(function->nParams(), nParams);
}
}

void set1DSpectrumEmpty() {
Expand Down Expand Up @@ -290,10 +343,15 @@ class CalculateChiSquaredTest : public CxxTest::TestSuite {
FunctionValues y(x);
function->function(x, y);
const double tmp = yValues[i] - y[0];
//std::cerr << "test " << yValues[i] << ' ' << y[0] << std::endl;
//std::cerr << "test " << xValue << ' ' << yValues[i] << ' ' << y[0] << std::endl;
sum2 += tmp * tmp;
}
}
if (divideByDOF) {
double dof = double(nData) - double(nParams);
if (dof <= 0.0) dof = 1.0;
sum2 /= dof;
}
TS_ASSERT_DIFFERS(sum2, 0);
TS_ASSERT_DELTA(sum2, chiSquared, 1e-10);
}
Expand Down

0 comments on commit 8f220b0

Please sign in to comment.