Commit d12a0b4f authored by Mathieu Tillet's avatar Mathieu Tillet
Browse files

Change computation method for revlog in rebin

parent cd464a21
......@@ -507,14 +507,13 @@ public:
MatrixWorkspace_sptr out = AnalysisDataService::Instance().retrieveWS<MatrixWorkspace>("test_Rebin_revLog");
auto &outX = out->x(0);
TS_ASSERT_EQUALS(outX.size(), 7);
TS_ASSERT_EQUALS(outX[0], 1);
TS_ASSERT_EQUALS(outX[1], 6);
TS_ASSERT_EQUALS(outX[2], 22);
TS_ASSERT_EQUALS(outX[3], 30);
TS_ASSERT_EQUALS(outX[4], 34);
TS_ASSERT_EQUALS(outX[5], 36);
TS_ASSERT_EQUALS(outX[6], 37);
TS_ASSERT_EQUALS(outX.size(), 6);
TS_ASSERT_DELTA(outX[0], 1, 1e-5);
TS_ASSERT_DELTA(outX[1], 22, 1e-5);
TS_ASSERT_DELTA(outX[2], 30, 1e-5);
TS_ASSERT_DELTA(outX[3], 34, 1e-5);
TS_ASSERT_DELTA(outX[4], 36, 1e-5);
TS_ASSERT_DELTA(outX[5], 37, 1e-5);
AnalysisDataService::Instance().remove("test_Rebin_revLog");
}
......@@ -528,24 +527,24 @@ public:
rebin.initialize();
rebin.setPropertyValue("InputWorkspace", "test_Rebin_revLog");
rebin.setPropertyValue("OutputWorkspace", "test_Rebin_revLog");
rebin.setPropertyValue("Params", "1, -0.5, 37");
rebin.setPropertyValue("Params", "1, -2, 42");
rebin.setProperty("UseReverseLogarithmic", true);
TS_ASSERT_THROWS_NOTHING(rebin.execute());
MatrixWorkspace_sptr out = AnalysisDataService::Instance().retrieveWS<MatrixWorkspace>("test_Rebin_revLog");
auto &outX = out->x(0);
TS_ASSERT_EQUALS(outX.size(), 10);
TS_ASSERT_EQUALS(outX[0], 1);
TS_ASSERT_EQUALS(outX[7], 35.75);
TS_ASSERT_EQUALS(outX[8], 36.5);
TS_ASSERT_EQUALS(outX[9], 37);
TS_ASSERT_EQUALS(outX.size(), 4);
TS_ASSERT_DELTA(outX[0], 1, 1e-5);
TS_ASSERT_DELTA(outX[1], 34, 1e-5);
TS_ASSERT_DELTA(outX[2], 40, 1e-5);
TS_ASSERT_DELTA(outX[3], 42, 1e-5);
AnalysisDataService::Instance().remove("test_Rebin_revLog");
}
void test_reverseLogEdgeCase() {
// Check that if the parameters given are so that the edges of the bins fall perfectly, there is no offset by one in
// the split
// Check the case where the parameters given are so that the edges of the bins fall perfectly, and so no padding is
// needed
Workspace2D_sptr test_1D = Create1DWorkspace(51);
test_1D->setDistribution(false);
AnalysisDataService::Instance().add("test_Rebin_revLog", test_1D);
......@@ -554,18 +553,19 @@ public:
rebin.initialize();
rebin.setPropertyValue("InputWorkspace", "test_Rebin_revLog");
rebin.setPropertyValue("OutputWorkspace", "test_Rebin_revLog");
rebin.setPropertyValue("Params", "30, -1, 37");
rebin.setPropertyValue("Params", "1, -1, 16");
rebin.setProperty("UseReverseLogarithmic", true);
TS_ASSERT_THROWS_NOTHING(rebin.execute());
MatrixWorkspace_sptr out = AnalysisDataService::Instance().retrieveWS<MatrixWorkspace>("test_Rebin_revLog");
auto &outX = out->x(0);
TS_ASSERT_EQUALS(outX.size(), 4);
TS_ASSERT_EQUALS(outX[0], 30);
TS_ASSERT_EQUALS(outX[1], 34);
TS_ASSERT_EQUALS(outX[2], 36);
TS_ASSERT_EQUALS(outX[3], 37);
TS_ASSERT_EQUALS(outX.size(), 5);
TS_ASSERT_DELTA(outX[0], 1, 1e-5)
TS_ASSERT_DELTA(outX[1], 9, 1e-5);
TS_ASSERT_DELTA(outX[2], 13, 1e-5);
TS_ASSERT_DELTA(outX[3], 15, 1e-5);
TS_ASSERT_DELTA(outX[4], 16, 1e-5);
AnalysisDataService::Instance().remove("test_Rebin_revLog");
}
......@@ -579,22 +579,21 @@ public:
rebin.initialize();
rebin.setPropertyValue("InputWorkspace", "test_Rebin_revLog");
rebin.setPropertyValue("OutputWorkspace", "test_Rebin_revLog");
rebin.setPropertyValue("Params", "1, 2, 15, -1, 37");
rebin.setPropertyValue("Params", "1, 2, 5, -1, 100");
rebin.setProperty("UseReverseLogarithmic", true);
TS_ASSERT_THROWS_NOTHING(rebin.execute());
MatrixWorkspace_sptr out = AnalysisDataService::Instance().retrieveWS<MatrixWorkspace>("test_Rebin_revLog");
auto &outX = out->x(0);
TS_ASSERT_EQUALS(outX.size(), 13);
TS_ASSERT_EQUALS(outX[0], 1);
TS_ASSERT_EQUALS(outX[6], 13);
TS_ASSERT_EQUALS(outX[7], 15);
TS_ASSERT_EQUALS(outX[8], 22);
TS_ASSERT_EQUALS(outX[9], 30);
TS_ASSERT_EQUALS(outX[10], 34);
TS_ASSERT_EQUALS(outX[11], 36);
TS_ASSERT_EQUALS(outX[12], 37);
TS_ASSERT_EQUALS(outX.size(), 7); // 2 lin + 4 log
TS_ASSERT_DELTA(outX[0], 1, 1e-5);
TS_ASSERT_DELTA(outX[1], 3, 1e-5);
TS_ASSERT_DELTA(outX[2], 5, 1e-5);
TS_ASSERT_DELTA(outX[3], 65, 1e-5);
TS_ASSERT_DELTA(outX[4], 85, 1e-5);
TS_ASSERT_DELTA(outX[5], 95, 1e-5);
TS_ASSERT_DELTA(outX[6], 100, 1e-5);
AnalysisDataService::Instance().remove("test_Rebin_revLog");
}
......@@ -609,27 +608,30 @@ public:
rebin.initialize();
rebin.setPropertyValue("InputWorkspace", "test_Rebin_revLog");
rebin.setPropertyValue("OutputWorkspace", "test_Rebin_revLog");
rebin.setPropertyValue("Params", "1, 2, 15, -1, 36, 0.5, 38");
rebin.setPropertyValue("Params", "1, 2, 5, -1, 100, 2, 110");
rebin.setProperty("UseReverseLogarithmic", true);
TS_ASSERT_THROWS_NOTHING(rebin.execute());
MatrixWorkspace_sptr out = AnalysisDataService::Instance().retrieveWS<MatrixWorkspace>("test_Rebin_revLog");
auto &outX = out->x(0);
TS_ASSERT_EQUALS(outX.size(), 16);
TS_ASSERT_EQUALS(outX[0], 1);
TS_ASSERT_EQUALS(outX[6], 13);
TS_ASSERT_EQUALS(outX[7], 15);
TS_ASSERT_EQUALS(outX[8], 22);
TS_ASSERT_EQUALS(outX[9], 30);
TS_ASSERT_EQUALS(outX[10], 34);
TS_ASSERT_EQUALS(outX[11], 36);
TS_ASSERT_EQUALS(outX[12], 36.5);
TS_ASSERT_EQUALS(outX.size(), 12);
TS_ASSERT_DELTA(outX[0], 1, 1e-5);
TS_ASSERT_DELTA(outX[1], 3, 1e-5);
TS_ASSERT_DELTA(outX[2], 5, 1e-5);
TS_ASSERT_DELTA(outX[3], 65, 1e-5);
TS_ASSERT_DELTA(outX[4], 85, 1e-5);
TS_ASSERT_DELTA(outX[5], 95, 1e-5);
TS_ASSERT_DELTA(outX[6], 100, 1e-5);
TS_ASSERT_DELTA(outX[7], 102, 1e-5);
TS_ASSERT_DELTA(outX[11], 110, 1e-5);
AnalysisDataService::Instance().remove("test_Rebin_revLog");
}
void test_reverseLogFullBinsOnly() {
// Test UseReverseLogarithmic with the FullBinsOnly option checked. It should change absolutely nothing
// Test UseReverseLogarithmic with the FullBinsOnly option checked. It should not change anything from the non
// checked version.
Workspace2D_sptr test_1D = Create1DWorkspace(51);
test_1D->setDistribution(false);
AnalysisDataService::Instance().add("test_Rebin_revLog", test_1D);
......@@ -646,14 +648,44 @@ public:
MatrixWorkspace_sptr out = AnalysisDataService::Instance().retrieveWS<MatrixWorkspace>("test_Rebin_revLog");
auto &outX = out->x(0);
TS_ASSERT_EQUALS(outX.size(), 6);
TS_ASSERT_DELTA(outX[0], 1, 1e-5);
TS_ASSERT_DELTA(outX[1], 22, 1e-5);
TS_ASSERT_DELTA(outX[2], 30, 1e-5);
TS_ASSERT_DELTA(outX[3], 34, 1e-5);
TS_ASSERT_DELTA(outX[4], 36, 1e-5);
TS_ASSERT_DELTA(outX[5], 37, 1e-5);
AnalysisDataService::Instance().remove("test_Rebin_revLog");
}
void test_reverseLogIncompleteFirstBin() {
// Test UseReverseLogarithmic with a first bin that is incomplete, but still bigger than the next one so not merged.
Workspace2D_sptr test_1D = Create1DWorkspace(51);
test_1D->setDistribution(false);
AnalysisDataService::Instance().add("test_Rebin_revLog", test_1D);
Rebin rebin;
rebin.initialize();
rebin.setPropertyValue("InputWorkspace", "test_Rebin_revLog");
rebin.setPropertyValue("OutputWorkspace", "test_Rebin_revLog");
rebin.setPropertyValue("Params", "1, -1, 60");
rebin.setProperty("UseReverseLogarithmic", true);
rebin.setProperty("FullBinsOnly", true);
TS_ASSERT_THROWS_NOTHING(rebin.execute());
MatrixWorkspace_sptr out = AnalysisDataService::Instance().retrieveWS<MatrixWorkspace>("test_Rebin_revLog");
auto &outX = out->x(0);
TS_ASSERT_EQUALS(outX.size(), 7);
TS_ASSERT_EQUALS(outX[0], 1);
TS_ASSERT_EQUALS(outX[1], 6);
TS_ASSERT_EQUALS(outX[2], 22);
TS_ASSERT_EQUALS(outX[3], 30);
TS_ASSERT_EQUALS(outX[4], 34);
TS_ASSERT_EQUALS(outX[5], 36);
TS_ASSERT_EQUALS(outX[6], 37);
TS_ASSERT_DELTA(outX[0], 1, 1e-5);
TS_ASSERT_DELTA(outX[1], 29, 1e-5);
TS_ASSERT_DELTA(outX[2], 45, 1e-5);
TS_ASSERT_DELTA(outX[3], 53, 1e-5);
TS_ASSERT_DELTA(outX[4], 57, 1e-5);
TS_ASSERT_DELTA(outX[5], 59, 1e-5);
TS_ASSERT_DELTA(outX[6], 60, 1e-5);
AnalysisDataService::Instance().remove("test_Rebin_revLog");
}
......
......@@ -70,30 +70,27 @@ int DLLExport createAxisFromRebinParams(const std::vector<double> &params, std::
if (resize_xnew)
xnew.emplace_back(xcurr);
int n0 = -1;
while ((ibound <= ibounds) && (istep <= isteps)) {
// if step is negative then it is logarithmic step
bool isLogBin = (fullParams[istep] < 0.0);
bool isReverseLogBin = isLogBin && useReverseLogarithmic;
double alpha = std::fabs(fullParams[istep]);
if (isReverseLogBin && xcurr == fullParams[ibound - 2]) {
// we are starting a new bin, but since it is a rev log, xcurr needs to be at its end
xcurr = fullParams[ibound];
}
if (!isLogBin)
xs = fullParams[istep];
else {
if (useReverseLogarithmic) {
double alpha = std::fabs(fullParams[istep]);
double x0 = fullParams[ibound];
// we go through a reverse log bin by starting from its end, and working our way back to the beginning
// this way we can define the bins in a reccuring way, and with a more obvious closeness with the usual log.
double x0 = fullParams[ibound - 2];
double step = x0 + fullParams[ibound] - xcurr;
if (n0 < 0) {
// if n0 is negative, that means we are starting a new reverse log interval split. In this case, we compute
// n0, which is the number of bins this interval still has to be split in.
// this is the formula for n0
n0 = static_cast<int>(std::ceil(std::log((xcurr - xMaxHint) / (x0 - xMaxHint)) / std::log1p(alpha) - 1));
}
xs = std::pow(1 + alpha, n0) * (x0 - xMaxHint) + xMaxHint - xcurr;
n0 -= 1;
xs = -step * alpha;
} else
xs = xcurr * fabs(fullParams[istep]);
......@@ -106,10 +103,9 @@ int DLLExport createAxisFromRebinParams(const std::vector<double> &params, std::
throw std::runtime_error("An infinite or NaN value was found in the binning parameters.");
}
if ((xcurr + xs * (1.0 + lastBinCoef) <= fullParams[ibound]) ||
(isReverseLogBin && (xcurr + xs < fullParams[ibound]))) {
if ((!isReverseLogBin && xcurr + xs * (1.0 + lastBinCoef) <= fullParams[ibound]) ||
(isReverseLogBin && xcurr + 2 * xs >= fullParams[ibound - 2])) {
// If we can still fit current bin _plus_ specified portion of a last bin, continue
// For reverse log, as long as the next xcurr is within the bound, continue.
xcurr += xs;
} else {
......@@ -123,9 +119,10 @@ int DLLExport createAxisFromRebinParams(const std::vector<double> &params, std::
// For non full_bins_only, finish by adding as much as is left from the range
xcurr = fullParams[ibound];
} else {
// we have finished this range, because its starting time has already been added, so we jump back to the last
// value of the bin and resume normal behaviour
xcurr = fullParams[ibound];
}
istep += 2;
ibound += 2;
}
......@@ -133,6 +130,7 @@ int DLLExport createAxisFromRebinParams(const std::vector<double> &params, std::
xnew.emplace_back(xcurr);
inew++;
}
std::sort(xnew.begin(), xnew.end());
return inew;
}
......
......@@ -17,8 +17,8 @@ create logarithmic binning using the formula
If UseReverseLogarithmic is set to True, the negative values create reverse logarithmic binning, i.e. a binning whose
bins are first big and then get exponentially smaller as they become close to the end of the specified interval.
The formula used is :math:`x(j) = (1+|\Delta x_i\|)^j * (x_{i+1} - f) + f` where :math:`f` is the last time of the
workspace.
The bins are the same as those generated by a normal log with these parameters, but arranged in the reversed order.
The first bin will be merged with the second one if the former is smaller than the latter.
This algorithms is useful both in data reduction, but also in remapping
:ref:`ragged workspace <Ragged_Workspace>` to a regular set of bin
......@@ -92,8 +92,7 @@ following will happen:
Hence the actual *Param* string used is "0, 2, 4, 3, 10".
This flag is ignored by the UseReverseLogarithmic flag on logarithmic bin, ie reverse logarithmic splits will not
necessarily begin with a full bin, no matter this flag. Other splits still follow the normal behavior though.
This flag is ignored if UseReverseLogarithm is checked.
.. _rebin-usage:
......@@ -151,7 +150,7 @@ Output:
dataY = [1,2,3,4,5,6,7,8,9] # or use dataY=range(1,10)
ws = CreateWorkspace(dataX, dataY)
# rebin from min to max - 1 with reverse logarithmic bins of 1
# rebin from min to max - 1 with reverse logarithmic bins of growth factor 1
ws = Rebin(ws, "1, -1, 9", UseReverseLogarithmic=True)
print("The rebinned X values are: {}".format(ws.readX(0)))
......@@ -160,7 +159,7 @@ Output:
.. testoutput:: ExHistRevLog
The rebinned X values are: [ 1. 2. 6. 8. 9.]
The rebinned X values are: [ 1. 6. 8. 9.]
**Example - custom two regions rebinning:**
......
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