-
Notifications
You must be signed in to change notification settings - Fork 128
Fix relative comparison in CompareWorkspaces
#37889
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Changes from all commits
5315c8d
95fdcdc
eb5e0ae
b721688
0701700
2c83bb8
eaf3985
File filter
Filter by extension
Conversations
Jump to
Diff view
Diff view
There are no files selected for viewing
Original file line number | Diff line number | Diff line change |
---|---|---|
|
@@ -78,13 +78,13 @@ int compareEventLists(Kernel::Logger &logger, const EventList &el1, const EventL | |
diffpulse = true; | ||
++numdiffpulse; | ||
} | ||
if (fabs(e1.tof() - e2.tof()) > tolTof) { | ||
if (std::fabs(e1.tof() - e2.tof()) > tolTof) { | ||
difftof = true; | ||
++numdifftof; | ||
} | ||
if (diffpulse && difftof) | ||
++numdiffboth; | ||
if (fabs(e1.weight() - e2.weight()) > tolWeight) { | ||
if (std::fabs(e1.weight() - e2.weight()) > tolWeight) { | ||
diffweight = true; | ||
++numdiffweight; | ||
} | ||
|
@@ -265,7 +265,7 @@ void CompareWorkspaces::processGroups(const std::shared_ptr<const API::Workspace | |
checker->setPropertyValue("Workspace1", namesOne[i]); | ||
checker->setPropertyValue("Workspace2", namesTwo[i]); | ||
for (size_t j = 0; j < numNonDefault; ++j) { | ||
Property *p = nonDefaultProps[j]; | ||
Property const *p = nonDefaultProps[j]; | ||
checker->setPropertyValue(p->name(), p->value()); | ||
} | ||
checker->execute(); | ||
|
@@ -592,15 +592,26 @@ bool CompareWorkspaces::checkData(const API::MatrixWorkspace_const_sptr &ws1, | |
// Cache a few things for later use | ||
const size_t numHists = ws1->getNumberHistograms(); | ||
bool raggedWorkspace{false}; | ||
size_t numBins; | ||
size_t numBins(0UL); | ||
try { | ||
numBins = ws1->blocksize(); | ||
} catch (std::length_error &) { | ||
raggedWorkspace = true; | ||
} | ||
const bool histogram = ws1->isHistogramData(); | ||
const bool checkAllData = getProperty("CheckAllData"); | ||
const bool RelErr = getProperty("ToleranceRelErr"); | ||
const bool isRelErr = getProperty("ToleranceRelErr"); | ||
const double tolerance = getProperty("Tolerance"); | ||
std::function<bool(double const, double const)> compare; | ||
if (isRelErr) { | ||
compare = [tolerance](double const x1, double const x2) -> bool { | ||
return CompareWorkspaces::withinRelativeTolerance(x1, x2, tolerance); | ||
}; | ||
} else { | ||
compare = [tolerance](double const x1, double const x2) -> bool { | ||
return CompareWorkspaces::withinAbsoluteTolerance(x1, x2, tolerance); | ||
}; | ||
} | ||
|
||
// First check that the workspace are the same size | ||
if (numHists != ws2->getNumberHistograms() || | ||
|
@@ -615,8 +626,8 @@ bool CompareWorkspaces::checkData(const API::MatrixWorkspace_const_sptr &ws1, | |
return false; | ||
} | ||
|
||
const double tolerance = getProperty("Tolerance"); | ||
bool resultBool = true; | ||
bool logDebug = g_log.is(Logger::Priority::PRIO_DEBUG); | ||
|
||
// Now check the data itself | ||
PARALLEL_FOR_IF(m_parallelComparison && ws1->threadSafe() && ws2->threadSafe()) | ||
|
@@ -642,29 +653,26 @@ bool CompareWorkspaces::checkData(const API::MatrixWorkspace_const_sptr &ws1, | |
} else { | ||
|
||
for (int j = 0; j < static_cast<int>(Y1.size()); ++j) { | ||
bool err; | ||
if (RelErr) { | ||
err = | ||
(relErr(X1[j], X2[j], tolerance) || relErr(Y1[j], Y2[j], tolerance) || relErr(E1[j], E2[j], tolerance)); | ||
} else | ||
err = (std::fabs(X1[j] - X2[j]) > tolerance || std::fabs(Y1[j] - Y2[j]) > tolerance || | ||
std::fabs(E1[j] - E2[j]) > tolerance); | ||
|
||
bool err = (!compare(X1[j], X2[j]) || !compare(Y1[j], Y2[j]) || !compare(E1[j], E2[j])); | ||
if (err) { | ||
g_log.debug() << "Data mismatch at cell (hist#,bin#): (" << i << "," << j << ")\n"; | ||
g_log.debug() << " Dataset #1 (X,Y,E) = (" << X1[j] << "," << Y1[j] << "," << E1[j] << ")\n"; | ||
g_log.debug() << " Dataset #2 (X,Y,E) = (" << X2[j] << "," << Y2[j] << "," << E2[j] << ")\n"; | ||
g_log.debug() << " Difference (X,Y,E) = (" << std::fabs(X1[j] - X2[j]) << "," << std::fabs(Y1[j] - Y2[j]) | ||
<< "," << std::fabs(E1[j] - E2[j]) << ")\n"; | ||
if (logDebug) { | ||
g_log.debug() << "Data mismatch at cell (hist#,bin#): (" << i << "," << j << ")\n"; | ||
g_log.debug() << " Dataset #1 (X,Y,E) = (" << X1[j] << "," << Y1[j] << "," << E1[j] << ")\n"; | ||
g_log.debug() << " Dataset #2 (X,Y,E) = (" << X2[j] << "," << Y2[j] << "," << E2[j] << ")\n"; | ||
g_log.debug() << " Difference (X,Y,E) = (" << std::fabs(X1[j] - X2[j]) << "," << std::fabs(Y1[j] - Y2[j]) | ||
<< "," << std::fabs(E1[j] - E2[j]) << ")\n"; | ||
} | ||
PARALLEL_CRITICAL(resultBool) | ||
resultBool = false; | ||
} | ||
} | ||
|
||
// Extra one for histogram data | ||
if (histogram && std::fabs(X1.back() - X2.back()) > tolerance) { | ||
g_log.debug() << " Data ranges mismatch for spectra N: (" << i << ")\n"; | ||
g_log.debug() << " Last bin ranges (X1_end vs X2_end) = (" << X1.back() << "," << X2.back() << ")\n"; | ||
if (histogram && !compare(X1.back(), X2.back())) { | ||
if (logDebug) { | ||
g_log.debug() << " Data ranges mismatch for spectra N: (" << i << ")\n"; | ||
g_log.debug() << " Last bin ranges (X1_end vs X2_end) = (" << X1.back() << "," << X2.back() << ")\n"; | ||
} | ||
PARALLEL_CRITICAL(resultBool) | ||
resultBool = false; | ||
} | ||
|
@@ -676,7 +684,7 @@ bool CompareWorkspaces::checkData(const API::MatrixWorkspace_const_sptr &ws1, | |
|
||
if (!resultBool) | ||
recordMismatch("Data mismatch"); | ||
// If all is well, return true | ||
// return result | ||
return resultBool; | ||
} | ||
|
||
|
@@ -944,7 +952,7 @@ bool CompareWorkspaces::checkRunProperties(const API::Run &run1, const API::Run | |
return false; | ||
} else { | ||
// Sort logs by name before one-by-one comparison | ||
auto compareNames = [](Kernel::Property *p1, Kernel::Property *p2) { return p1->name() < p2->name(); }; | ||
auto compareNames = [](Kernel::Property const *p1, Kernel::Property const *p2) { return p1->name() < p2->name(); }; | ||
std::sort(ws1logs.begin(), ws1logs.end(), compareNames); | ||
std::sort(ws2logs.begin(), ws2logs.end(), compareNames); | ||
for (size_t i = 0; i < ws1logs.size(); ++i) { | ||
|
@@ -1097,17 +1105,13 @@ void CompareWorkspaces::doPeaksComparison(PeaksWorkspace_sptr tws1, PeaksWorkspa | |
} | ||
bool mismatch = false; | ||
if (isRelErr) { | ||
if (relErr(s1, s2, tolerance)) { | ||
mismatch = true; | ||
} | ||
} else if (std::fabs(s1 - s2) > tolerance) { | ||
mismatch = true; | ||
} else if (std::fabs(v1[0] - v2[0]) > tolerance) { | ||
mismatch = true; | ||
} else if (std::fabs(v1[1] - v2[1]) > tolerance) { | ||
mismatch = true; | ||
} else if (std::fabs(v1[2] - v2[2]) > tolerance) { | ||
mismatch = true; | ||
mismatch = !withinRelativeTolerance(s1, s2, tolerance); | ||
// Q: whould we not also compare the vectors? | ||
} else { | ||
mismatch = !withinAbsoluteTolerance(s1, s2, tolerance) || // | ||
!withinAbsoluteTolerance(v1[0], v2[0], tolerance) || // | ||
!withinAbsoluteTolerance(v1[1], v2[1], tolerance) || // | ||
!withinAbsoluteTolerance(v1[2], v2[2], tolerance); // | ||
} | ||
if (mismatch) { | ||
g_log.notice(name); | ||
|
@@ -1213,11 +1217,12 @@ void CompareWorkspaces::doLeanElasticPeaksComparison(const LeanElasticPeaksWorks | |
g_log.information() << "Column " << name << " is not compared\n"; | ||
} | ||
bool mismatch = false; | ||
// Q: why does it not perform the user-specified operation for QLab and QSample? | ||
if (isRelErr && name != "QLab" && name != "QSample") { | ||
if (relErr(s1, s2, tolerance)) { | ||
if (!withinRelativeTolerance(s1, s2, tolerance)) { | ||
mismatch = true; | ||
} | ||
} else if (std::fabs(s1 - s2) > tolerance) { | ||
} else if (!withinAbsoluteTolerance(s1, s2, tolerance)) { | ||
mismatch = true; | ||
} | ||
if (mismatch) { | ||
|
@@ -1275,19 +1280,13 @@ void CompareWorkspaces::doTableComparison(const API::ITableWorkspace_const_sptr | |
const auto c2 = tws2->getColumn(i); | ||
|
||
if (isRelErr) { | ||
if (!c1->equalsRelErr(*c2, tolerance)) { | ||
mismatch = true; | ||
} | ||
mismatch = !c1->equalsRelErr(*c2, tolerance); | ||
} else { | ||
|
||
if (!c1->equals(*c2, tolerance)) { | ||
mismatch = true; | ||
} | ||
mismatch = !c1->equals(*c2, tolerance); | ||
} | ||
if (mismatch) { | ||
g_log.debug() << "Table data mismatch at column " << i << "\n"; | ||
recordMismatch("Table data mismatch"); | ||
mismatch = false; | ||
if (!checkAllData) { | ||
return; | ||
} | ||
|
@@ -1340,27 +1339,41 @@ void CompareWorkspaces::recordMismatch(const std::string &msg, std::string ws1, | |
m_result = false; | ||
} | ||
|
||
//------------------------------------------------------------------------------------------------ | ||
/** Function which calculates absolute error between two values and analyses if | ||
this error is within the limits requested. | ||
|
||
@param x1 -- first value to check difference | ||
@param x2 -- second value to check difference | ||
@param atol -- the tolerance of the comparison. Must be nonnegative | ||
|
||
@returns true if absolute difference is within the tolerance; false otherwise | ||
*/ | ||
bool CompareWorkspaces::withinAbsoluteTolerance(double const x1, double const x2, double const atol) { | ||
// NOTE !(|x1-x2| > atol) is not the same as |x1-x2| <= atol | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. very good comment There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Thanks. For anyone from the future coming to this PR and wondering about this comment. Mathematically they are identical, but they do not evaluate the same in some of the tests. The checks throughout the algo were originally written to expect It seemed to make more sense to return Two possibilities are values such as There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. I was curious why. |
||
return !(std::fabs(x1 - x2) > atol); | ||
} | ||
|
||
//------------------------------------------------------------------------------------------------ | ||
/** Function which calculates relative error between two values and analyses if | ||
this error is within the limits | ||
* requested. When the absolute value of the difference is smaller then the value | ||
of the error requested, | ||
* absolute error is used instead of relative error. | ||
this error is within the limits requested. | ||
|
||
@param x1 -- first value to check difference | ||
@param x2 -- second value to check difference | ||
@param errorVal -- the value of the error, to check against. Should be large | ||
then 0 | ||
@param x1 -- first value to check difference | ||
@param x2 -- second value to check difference | ||
@param rtol -- the tolerance of the comparison. Must be nonnegative | ||
|
||
@returns true if error or false if the value is within the limits requested | ||
@returns true if relative difference is within the tolerance; false otherwise | ||
@returns true if error or false if the relative value is within the limits requested | ||
*/ | ||
bool CompareWorkspaces::relErr(double x1, double x2, double errorVal) const { | ||
double num = std::fabs(x1 - x2); | ||
// how to treat x1<0 and x2 > 0 ? probably this way | ||
double den = 0.5 * (std::fabs(x1) + std::fabs(x2)); | ||
if (den < errorVal) | ||
return (num > errorVal); | ||
|
||
return (num / den > errorVal); | ||
bool CompareWorkspaces::withinRelativeTolerance(double const x1, double const x2, double const rtol) { | ||
// calculate difference | ||
double const num = std::fabs(x1 - x2); | ||
// return early if the values are equal | ||
if (num == 0.0) | ||
return true; | ||
// compare the difference to the midpoint value -- could lead to issues for negative values | ||
double const den = 0.5 * std::fabs(x1 + x2); | ||
// NOTE !(num > rtol*den) is not the same as (num <= rtol*den) | ||
return !(num > (rtol * den)); | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. please add the comment about this not being the same as There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Will do. People from the future, please see this note above. |
||
} | ||
} // namespace Mantid::Algorithms |
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1 @@ | ||
- fix :ref:`CompareWorkspaces <algm-CompareWorkspaces-v1>` for relative differences of small values. |
Uh oh!
There was an error while loading. Please reload this page.