diff --git a/Code/Mantid/Algorithms/inc/MantidAlgorithms/InterpolatingRebin.h b/Code/Mantid/Algorithms/inc/MantidAlgorithms/InterpolatingRebin.h
index b2bfab663149f7cd8a41e0f1396a5302c78b2897..72eaacd1e9de9c5c0d326bec19e3a2fca549bc9e 100644
--- a/Code/Mantid/Algorithms/inc/MantidAlgorithms/InterpolatingRebin.h
+++ b/Code/Mantid/Algorithms/inc/MantidAlgorithms/InterpolatingRebin.h
@@ -80,7 +80,8 @@ protected:
void outputYandEValues(API::MatrixWorkspace_const_sptr inputW, const DataObjects::Histogram1D::RCtype &XValues_new, API::MatrixWorkspace_sptr outputW);
void cubicInterpolation(const MantidVec &xOld, const MantidVec &yOld, const MantidVec &eOld, const MantidVec& xNew, MantidVec &yNew, MantidVec &eNew) const;
- double estimateError(const std::vector<double>& xsOld, const std::vector<double>& esOld, const double xNew) const;
+ void noInterpolation(const MantidVec &xOld, const double yOld, const MantidVec &eOld, const MantidVec& xNew, MantidVec &yNew, MantidVec &eNew) const;
+ double estimateError(const MantidVec &xsOld, const MantidVec &esOld, const double xNew) const;
};
} // namespace Algorithms
diff --git a/Code/Mantid/Algorithms/src/InterpolatingRebin.cpp b/Code/Mantid/Algorithms/src/InterpolatingRebin.cpp
index 2210d1a87e57879c8dc61ec7fb621d9ba016308a..9505b97d42d868a2e0e840482feae3d63aa25969 100644
--- a/Code/Mantid/Algorithms/src/InterpolatingRebin.cpp
+++ b/Code/Mantid/Algorithms/src/InterpolatingRebin.cpp
@@ -118,10 +118,8 @@ namespace Mantid
const int histnumber = inputW->getNumberHistograms();
Progress prog(this,0.0,1.0,histnumber);
-// PARALLEL_FOR2(inputW,outputW)
for (int hist=0; hist < histnumber;++hist)
{
-// PARALLEL_START_INTERUPT_REGION
// get const references to input Workspace arrays (no copying)
const MantidVec& XValues = inputW->readX(hist);
const MantidVec& YValues = inputW->readY(hist);
@@ -147,9 +145,7 @@ namespace Mantid
outputW->setX(hist, XValues_new);
prog.report();
-// PARALLEL_END_INTERUPT_REGION
}
-// PARALLEL_CHECK_INTERUPT_REGION
gsl_set_error_handler(old_handler);
}
@@ -181,11 +177,13 @@ namespace Mantid
{
// Make sure y and e vectors are of correct sizes
const size_t size_old = yOld.size();
+ if (size_old == 0 )
+ throw std::runtime_error("Empty spectrum found aborting");
if (size_old != (xOld.size() - 1) || size_old != eOld.size() )
- throw std::runtime_error("rebin: y and error vectors should be of same size & 1 shorter than x");
+ throw std::runtime_error("y and error vectors must be of same size & 1 shorter than x");
const size_t size_new = yNew.size();
if (size_new != (xNew.size() - 1) || size_new != eNew.size() )
- throw std::runtime_error("rebin: y and error vectors should be of same size & 1 shorter than x");
+ throw std::runtime_error("y and error vectors must be of same size & 1 shorter than x");
// get the bin centres of the input data
std::vector<double> xCensOld(size_new);
@@ -204,34 +202,81 @@ namespace Mantid
//check that the intepolation points fit well enough within the data for reliable intepolation to be done
if ( oldIn1<2 || oldIn2>=size_old-1 || oldIn1>oldIn2 )
{
- throw std::invalid_argument("Problem with the requested x-values to intepolate to: There must be at\nleast two input data points below the range of intepolation points and\ntwo higher. Also the intepolation points must have monatomically increasing x-values.");
+ double constantVal = -1;
+ if (VectorHelper::isConstantValue(yOld, constantVal))
+ {
+ //this copies the single y-value into the output array, errors are still calculated from the nearest input data points
+ noInterpolation(xOld, constantVal, eOld, xNew, yNew, eNew);
+ //this is as much as we need to do in this (trival) case
+ return;
+ }
+ else
+ {
+ throw std::invalid_argument("Problem with the requested x-values to intepolate to: There must be at\nleast two input data points below the range of intepolation points and\ntwo higher. Also the intepolation points must have monatomically increasing x-values.");
+ }
}
//bring one point before and one point after into the inpolation to reduce any possible errors from running out of data
oldIn1 -= 2;
oldIn2 ++;
- //get the GSL to allocate the memory, if this wasn't already done
- gsl_interp_accel *acc = gsl_interp_accel_alloc();
- const size_t nPoints = oldIn2 - oldIn1 + 1;
- gsl_spline *spline = gsl_spline_alloc(gsl_interp_cspline, nPoints);
- if ( ! acc || ! spline ||
- //GSL calculates the splines
- gsl_spline_init(spline, &xCensOld[oldIn1], &yOld[oldIn1], nPoints) )
+ //get the GSL to allocate the memory
+ gsl_interp_accel *acc = NULL;
+ gsl_spline *spline = NULL;
+ try
{
- throw std::runtime_error("Error setting up GSL spline functions");
- }
+ acc = gsl_interp_accel_alloc();
+ const size_t nPoints = oldIn2 - oldIn1 + 1;
+ spline = gsl_spline_alloc(gsl_interp_cspline, nPoints);
- for ( size_t i = 0; i < size_new; ++i )
- {
- yNew[i] = gsl_spline_eval(spline, xCensNew[i], acc);
- //(basic) error estimate the based on a weighted mean of the errors of the surrounding input data points
- eNew[i] = estimateError(xCensOld, eOld, xCensNew[i]);
- }
+ //test the allocation
+ if ( ! acc || ! spline ||
+ //calculate those splines, GSL uses pointers to the vector array (which is always contiguous)
+ gsl_spline_init(spline, &xCensOld[oldIn1], &yOld[oldIn1], nPoints) )
+ {
+ throw std::runtime_error("Error setting up GSL spline functions");
+ }
+ for ( size_t i = 0; i < size_new; ++i )
+ {
+ yNew[i] = gsl_spline_eval(spline, xCensNew[i], acc);
+ //(basic) error estimate the based on a weighted mean of the errors of the surrounding input data points
+ eNew[i] = estimateError(xCensOld, eOld, xCensNew[i]);
+ }
+ }
//for GSL to clear up its memory use
- gsl_spline_free (spline);
- gsl_interp_accel_free (acc);
+ catch (std::exception &)
+ {
+ if (acc)
+ {
+ if (spline)
+ {
+ gsl_spline_free(spline);
+ }
+ gsl_interp_accel_free(acc);
+ }
+ throw;
+ }
+ gsl_spline_free(spline);
+ gsl_interp_accel_free(acc);
+ }
+ /** This can be used whenever the original spectrum is filled with only one value. It is intended allow
+ * some spectra with null like values, for example all zeros
+ * @param[in] xOld the x-values of the input data
+ * @param[in] yOld the value of y that will be copied to the output array
+ * @param[in] eOld the error on each y-value, must be same length as yOld.
+ * @param[in] xNew x-values to rebin to, must be monotonically increasing
+ * @param[out] yNew is overwritten with the value repeated for as many times as there are x-values
+ * @param[out] eNew is overwritten with errors from the errors on the nearest input data points
+ */
+ void InterpolatingRebin::noInterpolation(const MantidVec &xOld, const double yOld, const MantidVec &eOld,
+ const MantidVec& xNew, MantidVec &yNew, MantidVec &eNew) const
+ {
+ yNew.assign(yNew.size(), yOld);
+ for ( MantidVec::size_type i = 0; i < eNew.size(); ++i )
+ {
+ eNew[i] = estimateError(xOld, eOld, xNew[i]);
+ }
}
/**Estimates the error on each interpolated point by assuming it is similar to the errors in
* near by input data points. Output points with the same x-value as an input point have the
@@ -242,13 +287,23 @@ namespace Mantid
* @param[in] xNew the value of x for at the point of interest
* @return the estimated error at that point
*/
- double InterpolatingRebin::estimateError(const std::vector<double>& xsOld, const std::vector<double>& esOld,
- const double xNew) const
+ double InterpolatingRebin::estimateError(const MantidVec &xsOld, const MantidVec &esOld, const double xNew) const
{
- //get the index of the first point that is higher in x, we'll base some of the error estimate on the error on this point
+ //find the first point in the array that has a higher value of x, we'll base some of the error estimate on the error on this point
const size_t indAbove =
std::lower_bound(xsOld.begin(), xsOld.end(), xNew) - xsOld.begin();
-
+
+ //if the point's x-value is out of the range covered by the x-values in the input data return the error value at the end of the range
+ if ( indAbove == 0 )
+ {
+ return esOld.front();
+ }
+ //xsOld is 1 longer than xsOld
+ if ( indAbove >= esOld.size() )
+ {//cubicInterpolation() checks that that there are no empty histograms
+ return esOld.back();
+ }
+
const double error1 = esOld[indAbove];
// ratio of weightings will be inversely proportional to the distance between the points
double weight1 = xsOld[indAbove] - xNew;
diff --git a/Code/Mantid/Algorithms/test/InterpolatingRebinTest.h b/Code/Mantid/Algorithms/test/InterpolatingRebinTest.h
index b9cbba4a0584aaf639adf372e448b1742bd6cfac..4964ec672d7c6cd3f4793944999db276a199d7e8 100644
--- a/Code/Mantid/Algorithms/test/InterpolatingRebinTest.h
+++ b/Code/Mantid/Algorithms/test/InterpolatingRebinTest.h
@@ -7,19 +7,20 @@
#include "MantidAPI/AnalysisDataService.h"
#include "MantidAlgorithms/InterpolatingRebin.h"
#include "MantidAPI/WorkspaceProperty.h"
+#include <cmath>
+#include <limits>
using namespace Mantid::Kernel;
using namespace Mantid::DataObjects;
using namespace Mantid::API;
using namespace Mantid::Algorithms;
-
class InterpolatingRebinTest : public CxxTest::TestSuite
{
public:
- void testworkspace_dist()
+ void testWorkspace_dist()
{
- Workspace2D_sptr test_in1D = Create1DWorkspace(50);
+ Workspace2D_sptr test_in1D = Create1DData();
test_in1D->isDistribution(true);
AnalysisDataService::Instance().add("InterpolatingRebinTest_indist", test_in1D);
@@ -43,14 +44,20 @@ public:
// set the new bins to be less than half the size of the old, one in every 2 old bins and one in every 5 old will coinside
rebin.setPropertyValue("Params", "2.225,0.2,15");
- TS_ASSERT(rebin.execute())
+ TS_ASSERT_THROWS_NOTHING(rebin.execute())
TS_ASSERT(rebin.isExecuted())
+ //get the output workspace and test it
MatrixWorkspace_sptr rebindata =
boost::dynamic_pointer_cast<MatrixWorkspace>(AnalysisDataService::Instance().retrieve("InterpolatingRebinTest_outdist"));
+ TS_ASSERT_EQUALS(rebindata->getNumberHistograms(), 1);
+
const Mantid::MantidVec outX=rebindata->dataX(0);
const Mantid::MantidVec outY=rebindata->dataY(0);
const Mantid::MantidVec outE=rebindata->dataE(0);
+ TS_ASSERT_EQUALS(outX.size(), ceil((15-2.225)/0.2)+1);
+ TS_ASSERT_EQUALS(outY.size(), ceil((15-2.225)/0.2));
+ TS_ASSERT_EQUALS(outE.size(), ceil((15-2.225)/0.2));
//this intepolated data was found by running the debugger on this test
TS_ASSERT_DELTA(outX[0], 2.225, 0.00001);
@@ -84,10 +91,10 @@ public:
AnalysisDataService::Instance().remove("InterpolatingRebinTest_outdist");
}
- void testworkspace_nondist()
+ void testWorkspace_nondist()
{
- Workspace2D_sptr test_in1D = Create1DWorkspace(50);
+ Workspace2D_sptr test_in1D = Create1DData();
test_in1D->isDistribution(false);
AnalysisDataService::Instance().add("InterpolatingRebinTest_in_nondist", test_in1D);
@@ -98,9 +105,10 @@ public:
// set the new bins to be less than half the size of the old, one in every 2 old bins and one in every 5 old will coinside
rebin.setPropertyValue("Params", "2.225,0.2,15");
- TS_ASSERT(rebin.execute())
+ TS_ASSERT_THROWS_NOTHING(rebin.execute())
TS_ASSERT(rebin.isExecuted())
+ //get the output workspace and test it
MatrixWorkspace_sptr rebindata =
boost::dynamic_pointer_cast<MatrixWorkspace>(AnalysisDataService::Instance().retrieve("InterpolatingRebinTest_out_nondist"));
const Mantid::MantidVec outX=rebindata->dataX(0);
@@ -141,14 +149,74 @@ public:
AnalysisDataService::Instance().remove("InterpolatingRebinTest_out_nondist");
}
+ void testNullDataHandling()
+ {
+
+ Workspace2D_sptr test_in1D = badData();
+ test_in1D->isDistribution(true);
+ AnalysisDataService::Instance().add("InterpolatingRebinTest_in_nulldata", test_in1D);
+
+ InterpolatingRebin rebin;
+ rebin.initialize();
+ rebin.setPropertyValue("InputWorkspace","InterpolatingRebinTest_in_nulldata");
+ rebin.setPropertyValue("OutputWorkspace","InterpolatingRebinTest_out_nulldata");
+
+ // set the new bins to be less than half the size of the old, one in every 2 old bins and one in every 5 old will coinside
+ rebin.setPropertyValue("Params", "2,0.2,11");
+ TS_ASSERT_THROWS_NOTHING(rebin.execute())
+ TS_ASSERT(rebin.isExecuted())
+
+ //get the output workspace and test it
+ MatrixWorkspace_sptr rebindata =
+ boost::dynamic_pointer_cast<MatrixWorkspace>(AnalysisDataService::Instance().retrieve("InterpolatingRebinTest_out_nulldata"));
+ TS_ASSERT_EQUALS(rebindata->getNumberHistograms(), 2);
+
+ Mantid::MantidVec outX=rebindata->dataX(0);
+ Mantid::MantidVec outY=rebindata->dataY(0);
+ Mantid::MantidVec outE=rebindata->dataE(0);
+ TS_ASSERT_EQUALS(outX.size(), ceil((11-2.0)/0.2)+1);
+ TS_ASSERT_EQUALS(outY.size(), ceil((11-2.0)/0.2));
+ TS_ASSERT_EQUALS(outE.size(), ceil((11-2.0)/0.2));
+
+ //the first spectrum should be only zeros test the first spectrum
+ TS_ASSERT_DELTA(outX[0], 2, 0.00001);
+ TS_ASSERT_DELTA(outY[0], 0, 0.0001);
+ TS_ASSERT_DELTA(outE[0], 0, 0.0001);
+
+ //test a random location
+ TS_ASSERT_DELTA(outX[2], 2.4, 0.00001);
+ TS_ASSERT_DELTA(outY[2], 0, 0.0001);
+ TS_ASSERT_DELTA(outE[2], 0, 0.0001);
+
+ //check the last point
+ TS_ASSERT_DELTA(outX[45], 11, 0.00001);
+ TS_ASSERT_DELTA(outY[44], 0, 0.0001);
+ TS_ASSERT_DELTA(outE[44], 0,0.0001);
+
+ //the second spectrum is NAN
+ outX=rebindata->dataX(1);
+ outY=rebindata->dataY(1);
+ outE=rebindata->dataE(1);
+ //test a random one
+ TS_ASSERT_DELTA(outX[7], 3.4, 0.00001);
+ //check for numeric_limits<double>::quiet_NaN()
+ TS_ASSERT(outY[7] != outY[7]);
+ TS_ASSERT_DELTA(outE[7], 2, 0.00001);
+
+ AnalysisDataService::Instance().remove("InterpolatingRebinTest_in_nulldata");
+ AnalysisDataService::Instance().remove("InterpolatingRebinTest_out_nulldata");
+ }
+
private:
- Workspace2D_sptr Create1DWorkspace(int size)
+ Workspace2D_sptr Create1DData()
{
+ static const int nBins = 50;
Workspace2D_sptr retVal(new Workspace2D);
- retVal->initialize(1,size,size-1);
+ retVal->initialize(1, nBins+1, nBins);
+
double j=1.0;
int i = 0;
- for ( ; i<size-1; i++)
+ for ( ; i<nBins; i++)
{
retVal->dataX(0)[i]=j*0.5;
retVal->dataY(0)[i] = j;
@@ -156,6 +224,31 @@ private:
j+=1.5;
}
retVal->dataX(0)[i]=j*0.5;
+
+ return retVal;
+ }
+
+ Workspace2D_sptr badData()
+ {
+ static const int nBins = 24;
+ Workspace2D_sptr retVal(new Workspace2D);
+ retVal->initialize(2, nBins+1, nBins);
+
+ //the first histogram has all zeros
+ for ( int i = 0; i<nBins-1; i++)
+ {
+ retVal->dataX(0)[i] = 0;
+ retVal->dataY(0)[i] = 0;
+ retVal->dataE(0)[i] = 0;
+ }
+ //the second has NAN values
+ for ( int i = 0; i<nBins-1; i++)
+ {
+ retVal->dataX(1)[i] = i;
+ retVal->dataY(1)[i] = std::numeric_limits<double>::quiet_NaN();
+ retVal->dataE(1)[i] = 2;
+ }
+
return retVal;
}
diff --git a/Code/Mantid/Kernel/inc/MantidKernel/VectorHelper.h b/Code/Mantid/Kernel/inc/MantidKernel/VectorHelper.h
index 324b7cfb6e134ea42726c4b2eebf200669d1c9f6..1338a44bfc8854f31df94aadf28f2da2fc62f07a 100644
--- a/Code/Mantid/Kernel/inc/MantidKernel/VectorHelper.h
+++ b/Code/Mantid/Kernel/inc/MantidKernel/VectorHelper.h
@@ -54,6 +54,8 @@ namespace VectorHelper
/// Convert an array of bin boundaries to bin centre values.
void DLLExport convertToBinCentre(const std::vector<double> & bin_edges, std::vector<double> & bin_centres);
+ bool DLLExport isConstantValue(const std::vector<double> &arra, double &val);
+
//! Functor used for computing the sum of the square values of a vector, using the accumulate algorithm
template <class T> struct SumGaussError: public std::binary_function<T,T,T>
{
diff --git a/Code/Mantid/Kernel/src/VectorHelper.cpp b/Code/Mantid/Kernel/src/VectorHelper.cpp
index 8459d8fdbd59fee3066a5b0b54c3fcbf6f2ffbeb..5ecf06a40eaf6cdda0b2153a15cc3070868a1a55 100644
--- a/Code/Mantid/Kernel/src/VectorHelper.cpp
+++ b/Code/Mantid/Kernel/src/VectorHelper.cpp
@@ -4,6 +4,7 @@
#include "MantidKernel/VectorHelper.h"
#include <algorithm>
#include <numeric>
+#include <limits>
#include <iostream>
namespace Mantid
@@ -305,6 +306,68 @@ void convertToBinCentre(const std::vector<double> & bin_edges, std::vector<doubl
bin_centres.erase(bin_centres.begin());
}
+/** Assess if all the values in the vector are equal or if there are some different values
+
+* @param[in] arra the vector to examine
+
+* @param[out] val if there is only one value this variable takes that value, otherwise its contents are undefined
+
+*/
+
+bool isConstantValue(const std::vector<double> &arra, double &val)
+
+{
+
+ //make comparisons with the first value
+
+ std::vector<double>::const_iterator i = arra.begin();
+
+
+ if ( i == arra.end() )
+
+ {//empty array
+ return true;
+
+ }
+
+
+
+ //this loop can be entered! NAN values make comparisons difficult because nan != nan, deal with these first
+
+ for ( val = *i; val != val ; )
+
+ {
+
+ ++i;
+
+ if ( i == arra.end() )
+
+ {
+
+ //all values are contant (NAN)
+
+ return true;
+
+ }
+
+ val = *i;
+
+ }
+
+
+
+ for ( ; i != arra.end() ; ++i )
+ {
+ if ( *i != val )
+ {
+ return false;
+ }
+ }
+ //no different value was found and so every must be equal to c
+ return true;
+}
+
+
} // End namespace VectorHelper
} // End namespace Kernel
} // End namespace Mantid
diff --git a/Code/Mantid/PythonAPI/scripts/SANS/SANSReduction.py b/Code/Mantid/PythonAPI/scripts/SANS/SANSReduction.py
index 03fd025e7349a4fcd9ffdcf6a1b13fe867e7ba8e..1c903074e1c95cbca413c726f6d5504b38bde002 100644
--- a/Code/Mantid/PythonAPI/scripts/SANS/SANSReduction.py
+++ b/Code/Mantid/PythonAPI/scripts/SANS/SANSReduction.py
@@ -104,11 +104,12 @@ BACKMON_END = None
# S2D: front-detector or rear-detector
DETBANK = None
-# The monitor spectrum taken from the GUI. Is this still necessary?? or can I just deduce
-# it from the instrument name
+# The monitor spectrum taken from the GUI
MONITORSPECTRUM = None
# agruments after MON/LENGTH need to take precendence over those after MON/SPECTRUM and this variable ensures that
MONITORSPECLOCKED = False
+# if this is set InterpolationRebin will be used on the monitor spectrum used to normalise the sample
+SAMP_INTERPOLATE = False
# Detector position information for SANS2D
FRONT_DET_RADIUS = 306.0
@@ -155,6 +156,8 @@ TRANS_WAV2_FULL = None
# Mon/Det for SANS2D
TRANS_UDET_MON = 2
TRANS_UDET_DET = 3
+# this is if to use InterpolatingRebin on the monitor spectrum used to normalise the transmission
+TRANS_INTERPOLATE = False
###################################################################################################################
#
@@ -229,9 +232,8 @@ def printParameter(var):
###########################
def SANS2D():
_printMessage('SANS2D()')
- global INSTR_NAME, MONITORSPECTRUM, TRANS_WAV1, TRANS_WAV2, TRANS_WAV1_FULL, TRANS_WAV2_FULL
+ global INSTR_NAME, TRANS_WAV1, TRANS_WAV2, TRANS_WAV1_FULL, TRANS_WAV2_FULL
INSTR_NAME = 'SANS2D'
- MONITORSPECTRUM = 2
TRANS_WAV1_FULL = TRANS_WAV1 = 2.0
TRANS_WAV2_FULL = TRANS_WAV2 = 14.0
if DETBANK != 'rear-detector':
@@ -741,6 +743,15 @@ def clearCurrentMaskDefaults():
global BACKMON_START, BACKMON_END
BACKMON_START = BACKMON_END = None
+
+ global MONITORSPECTRUM, MONITORSPECLOCKED, SAMP_INTERPOLATE, TRANS_UDET_MON, TRANS_UDET_DET, TRANS_INTERPOLATE
+ MONITORSPECTRUM = None
+ MONITORSPECLOCKED = False
+ SAMP_INTERPOLATE = False
+
+ TRANS_UDET_MON = 2
+ TRANS_UDET_DET = 3
+ TRANS_INTERPOLATE = False
####################################
# Add a mask to the correct string
@@ -830,56 +841,25 @@ def MaskFile(filename):
upper_line = line.upper()
if upper_line.startswith('L/'):
_readLimitValues(line)
+
elif upper_line.startswith('MON/'):
- details = line[4:]
- if details.upper().startswith('LENGTH'):
- SuggestMonitorSpectrum(int(details.split()[1]))
- elif details.upper().startswith('SPECTRUM'):
- SetMonitorSpectrum(int(details.split('=')[1]))
- elif 'DIRECT' in details.upper():
- parts = details.split("=")
- if len(parts) == 2:
- filepath = parts[1].rstrip()
- if '[' in filepath:
- idx = filepath.rfind(']')
- filepath = filepath[idx + 1:]
- if not os.path.isabs(filepath):
- filepath = os.path.join(USER_PATH, filepath)
- type = parts[0]
- parts = type.split("/")
- if len(parts) == 1:
- if parts[0] == 'DIRECT':
- SetRearEfficiencyFile(filepath)
- SetFrontEfficiencyFile(filepath)
- elif parts[0] == 'HAB':
- SetFrontEfficiencyFile(filepath)
- else:
- pass
- elif len(parts) == 2:
- detname = parts[1]
- if detname == 'REAR':
- SetRearEfficiencyFile(filepath)
- elif detname == 'FRONT' or detname == 'HAB':
- SetFrontEfficiencyFile(filepath)
- else:
- _issueWarning('Incorrect detector specified for efficiency file "' + line + '"')
- else:
- _issueWarning('Unable to parse monitor line "' + line + '"')
- else:
- _issueWarning('Unable to parse monitor line "' + line + '"')
- else:
- continue
+ _readMONValues(line)
+
elif upper_line.startswith('MASK'):
Mask(upper_line)
+
elif upper_line.startswith('SET CENTRE'):
values = upper_line.split()
SetCentre(float(values[2]), float(values[3]))
+
elif upper_line.startswith('SET SCALES'):
values = upper_line.split()
_SetScales(float(values[2]))
+
elif upper_line.startswith('SAMPLE/OFFSET'):
values = upper_line.split()
SetSampleOffset(values[1])
+
elif upper_line.startswith('DET/'):
type = upper_line[4:]
if type.startswith('CORR'):
@@ -887,6 +867,7 @@ def MaskFile(filename):
else:
# This checks whether the type is correct and issues warnings if it is not
Detector(type)
+
elif upper_line.startswith('GRAVITY'):
flag = upper_line[8:]
if flag == 'ON':
@@ -896,6 +877,7 @@ def MaskFile(filename):
else:
_issueWarning("Gravity flag incorrectly specified, disabling gravity correction")
Gravity(False)
+
elif upper_line.startswith('BACK/MON/TIMES'):
tokens = upper_line.split()
global BACKMON_START, BACKMON_END
@@ -906,6 +888,7 @@ def MaskFile(filename):
_issueWarning('Incorrectly formatted BACK/MON/TIMES line, not running FlatBackground.')
BACKMON_START = None
BACKMON_END = None
+
elif upper_line.startswith("FIT/TRANS/"):
params = upper_line[10:].split()
if len(params) == 3:
@@ -914,6 +897,7 @@ def MaskFile(filename):
else:
_issueWarning('Incorrectly formatted FIT/TRANS line, setting defaults to LOG and full range')
TransFit(TRANS_FIT_DEF)
+
else:
continue
@@ -996,6 +980,59 @@ def _readLimitValues(limit_line):
else:
pass
+def _readMONValues(line):
+ details = line[4:].upper()
+
+ #MON/LENTH, MON/SPECTRUM and MON/TRANS all accept the INTERPOLATE option
+ interpolate = False
+ if details.endswith('/INTERPOLATE') :
+ interpolate = True
+ details = details.split('/INTERPOLATE')[0]
+
+ if details.startswith('LENGTH'):
+ SuggestMonitorSpectrum(int(details.split()[1]), interpolate)
+
+ elif details.startswith('SPECTRUM'):
+ SetMonitorSpectrum(int(details.split('=')[1]), interpolate)
+
+ elif details.startswith('TRANS'):
+ parts = details.split('=')
+ if len(parts) < 2 or parts[0] != 'TRANS/SPECTRUM' :
+ _issueWarning('Unable to parse MON/TRANS line, needs MON/TRANS/SPECTRUM=')
+ SetTransSpectrum(int(parts[1]), interpolate)
+
+ elif 'DIRECT' in details:
+ parts = details.split("=")
+ if len(parts) == 2:
+ filepath = parts[1].rstrip()
+ if '[' in filepath:
+ idx = filepath.rfind(']')
+ filepath = filepath[idx + 1:]
+ if not os.path.isabs(filepath):
+ filepath = os.path.join(USER_PATH, filepath)
+ type = parts[0]
+ parts = type.split("/")
+ if len(parts) == 1:
+ if parts[0] == 'DIRECT':
+ SetRearEfficiencyFile(filepath)
+ SetFrontEfficiencyFile(filepath)
+ elif parts[0] == 'HAB':
+ SetFrontEfficiencyFile(filepath)
+ else:
+ pass
+ elif len(parts) == 2:
+ detname = parts[1]
+ if detname == 'REAR':
+ SetRearEfficiencyFile(filepath)
+ elif detname == 'FRONT' or detname == 'HAB':
+ SetFrontEfficiencyFile(filepath)
+ else:
+ _issueWarning('Incorrect detector specified for efficiency file "' + line + '"')
+ else:
+ _issueWarning('Unable to parse monitor line "' + line + '"')
+ else:
+ _issueWarning('Unable to parse monitor line "' + line + '"')
+
def _readDetectorCorrections(details):
values = details.split()
det_name = values[0]
@@ -1031,17 +1068,33 @@ def SetSampleOffset(value):
global SAMPLE_Z_CORR
SAMPLE_Z_CORR = float(value)/1000.
-def SetMonitorSpectrum(spec):
+def SetMonitorSpectrum(specNum, interp=False):
global MONITORSPECTRUM
- MONITORSPECTRUM = spec
+ MONITORSPECTRUM = int(specNum)
+
+ global SAMP_INTERPOLATE
+ SAMP_INTERPOLATE = bool(interp)
+
global MONITORSPECLOCKED
MONITORSPECLOCKED = True
-def SuggestMonitorSpectrum(spec):
+def SuggestMonitorSpectrum(specNum, interp=False):
global MONITORSPECLOCKED
- if not MONITORSPECLOCKED :
- global MONITORSPECTRUM
- MONITORSPECTRUM = spec
+ if MONITORSPECLOCKED :
+ return
+
+ global SAMP_INTERPOLATE
+ SAMP_INTERPOLATE = bool(interp)
+
+ global MONITORSPECTRUM
+ MONITORSPECTRUM = int(specNum)
+
+def SetTransSpectrum(specNum, interp=False):
+ global TRANS_UDET_MON
+ TRANS_UDET_MON = int(specNum)
+
+ global TRANS_INTERPOLATE
+ TRANS_INTERPOLATE = bool(interp)
def SetRearEfficiencyFile(filename):
global DIRECT_BEAM_FILE_R
@@ -1225,17 +1278,19 @@ def CalculateTransmissionCorrection(run_setup, lambdamin, lambdamax, use_def_tra
fit_type = 'Linear'
else:
fit_type = TRANS_FIT
+ #retrieve the user setting that tells us whether Rebin or InterpolatingRebin will be used during the normalisation
+ global TRANS_INTERPOLATE
if INSTR_NAME == 'LOQ':
# Change the instrument definition to the correct one in the LOQ case
LoadInstrument(trans_raw, INSTR_DIR + "/LOQ_trans_Definition.xml")
LoadInstrument(direct_raw, INSTR_DIR + "/LOQ_trans_Definition.xml")
- trans_tmp_out = SANSUtility.SetupTransmissionWorkspace(trans_raw, '1,2', BACKMON_START, BACKMON_END, wavbin, True)
- direct_tmp_out = SANSUtility.SetupTransmissionWorkspace(direct_raw, '1,2', BACKMON_START, BACKMON_END, wavbin, True)
+ trans_tmp_out = SANSUtility.SetupTransmissionWorkspace(trans_raw, '1,2', BACKMON_START, BACKMON_END, wavbin, TRANS_INTERPOLATE, True)
+ direct_tmp_out = SANSUtility.SetupTransmissionWorkspace(direct_raw, '1,2', BACKMON_START, BACKMON_END, wavbin, TRANS_INTERPOLATE, True)
CalculateTransmission(trans_tmp_out,direct_tmp_out, fittedtransws, MinWavelength = translambda_min, MaxWavelength = translambda_max, \
FitMethod = fit_type, OutputUnfittedData=True)
else:
- trans_tmp_out = SANSUtility.SetupTransmissionWorkspace(trans_raw, '1,2', BACKMON_START, BACKMON_END, wavbin, False)
- direct_tmp_out = SANSUtility.SetupTransmissionWorkspace(direct_raw, '1,2', BACKMON_START, BACKMON_END, wavbin, False)
+ trans_tmp_out = SANSUtility.SetupTransmissionWorkspace(trans_raw, '1,2', BACKMON_START, BACKMON_END, wavbin, TRANS_INTERPOLATE, False)
+ direct_tmp_out = SANSUtility.SetupTransmissionWorkspace(direct_raw, '1,2', BACKMON_START, BACKMON_END, wavbin, TRANS_INTERPOLATE, False)
CalculateTransmission(trans_tmp_out,direct_tmp_out, fittedtransws, TRANS_UDET_MON, TRANS_UDET_DET, MinWavelength = translambda_min, \
MaxWavelength = translambda_max, FitMethod = fit_type, OutputUnfittedData=True)
# Remove temporaries
@@ -1388,7 +1443,11 @@ def Correct(run_setup, wav_start, wav_end, use_def_trans, finding_centre = False
# ConvertUnits does have a rebin option, but it's crude. In particular it rebins on linear scale.
ConvertUnits(monitorWS, monitorWS, "Wavelength")
wavbin = str(wav_start) + "," + str(DWAV) + "," + str(wav_end)
- Rebin(monitorWS, monitorWS,wavbin)
+ if SAMP_INTERPOLATE :
+ InterpolatingRebin(monitorWS, monitorWS,wavbin)
+ else :
+ Rebin(monitorWS, monitorWS,wavbin)
+
ConvertUnits(final_result,final_result,"Wavelength")
Rebin(final_result,final_result,wavbin)
####################################################################################
diff --git a/Code/Mantid/PythonAPI/scripts/SANS/SANSUtility.py b/Code/Mantid/PythonAPI/scripts/SANS/SANSUtility.py
index c6c5cbe75db4ea69983b27a5e7db25a57a2f38a5..10f263e5d64b655529ffa49a5eb3ff8889159a13 100644
--- a/Code/Mantid/PythonAPI/scripts/SANS/SANSUtility.py
+++ b/Code/Mantid/PythonAPI/scripts/SANS/SANSUtility.py
@@ -234,7 +234,7 @@ def MaskByBinRange(workspace, timemask):
MaskBins(workspace, workspace, XMin= limits[0] ,XMax=limits[1])
# Setup the transmission workspace
-def SetupTransmissionWorkspace(inputWS, spec_list, backmon_start, backmon_end, wavbining, loqremovebins):
+def SetupTransmissionWorkspace(inputWS, spec_list, backmon_start, backmon_end, wavbining, interpolate, loqremovebins):
tmpWS = inputWS + '_tmp'
if loqremovebins == True:
RemoveBins(inputWS,tmpWS, 19900, 20500, Interpolation='Linear')
@@ -244,7 +244,12 @@ def SetupTransmissionWorkspace(inputWS, spec_list, backmon_start, backmon_end, w
inputWS = tmpWS
# Convert and rebin
ConvertUnits(inputWS,tmpWS,"Wavelength")
- Rebin(tmpWS, tmpWS, wavbining)
+
+ if interpolate :
+ InterpolatingRebin(tmpWS, tmpWS, wavbining)
+ else :
+ Rebin(tmpWS, tmpWS, wavbining)
+
return tmpWS
# Correct of for the volume of the sample/can. Dimensions should be in order: width, height, thickness
diff --git a/Code/qtiplot/MantidQt/CustomInterfaces/inc/SANSRunWindow.ui b/Code/qtiplot/MantidQt/CustomInterfaces/inc/SANSRunWindow.ui
index f4589eb4d83066622503f0e6d8fa90e145f15a4f..8c03506835e99f92521f395d5d8a316c4bfd0615 100644
--- a/Code/qtiplot/MantidQt/CustomInterfaces/inc/SANSRunWindow.ui
+++ b/Code/qtiplot/MantidQt/CustomInterfaces/inc/SANSRunWindow.ui
@@ -1687,14 +1687,14 @@ p, li { white-space: pre-wrap; }
</property>
</widget>
</item>
- <item row="1" column="0" >
+ <item row="2" column="0" >
<widget class="QLabel" name="label_43" >
<property name="text" >
<string>Detector</string>
</property>
</widget>
</item>
- <item row="1" column="1" colspan="2" >
+ <item row="2" column="1" colspan="3" >
<widget class="QComboBox" name="detbank_sel" >
<property name="sizePolicy" >
<sizepolicy vsizetype="Maximum" hsizetype="Preferred" >
@@ -1714,7 +1714,7 @@ p, li { white-space: pre-wrap; }
</item>
</widget>
</item>
- <item row="0" column="2" >
+ <item row="0" column="3" >
<spacer name="horizontalSpacer_6" >
<property name="orientation" >
<enum>Qt::Horizontal</enum>
@@ -1727,6 +1727,43 @@ p, li { white-space: pre-wrap; }
</property>
</spacer>
</item>
+ <item row="0" column="2" >
+ <widget class="QCheckBox" name="monitor_interp" >
+ <property name="text" >
+ <string>Interpolate</string>
+ </property>
+ </widget>
+ </item>
+ <item row="1" column="0" >
+ <widget class="QLabel" name="label_39" >
+ <property name="text" >
+ <string>Transmission monitor</string>
+ </property>
+ </widget>
+ </item>
+ <item row="1" column="1" >
+ <widget class="QLineEdit" name="trans_monitor" >
+ <property name="sizePolicy" >
+ <sizepolicy vsizetype="Fixed" hsizetype="Maximum" >
+ <horstretch>0</horstretch>
+ <verstretch>0</verstretch>
+ </sizepolicy>
+ </property>
+ <property name="maximumSize" >
+ <size>
+ <width>40</width>
+ <height>16777215</height>
+ </size>
+ </property>
+ </widget>
+ </item>
+ <item row="1" column="2" >
+ <widget class="QCheckBox" name="trans_interp" >
+ <property name="text" >
+ <string>Interpolate</string>
+ </property>
+ </widget>
+ </item>
</layout>
</widget>
</item>
diff --git a/Code/qtiplot/MantidQt/CustomInterfaces/src/SANSRunWindow.cpp b/Code/qtiplot/MantidQt/CustomInterfaces/src/SANSRunWindow.cpp
index 961601f53a52b3990da29a55d8eba1e6ce84dc99..a473d35eb562f4a60ac29278f12926707852c733 100644
--- a/Code/qtiplot/MantidQt/CustomInterfaces/src/SANSRunWindow.cpp
+++ b/Code/qtiplot/MantidQt/CustomInterfaces/src/SANSRunWindow.cpp
@@ -454,8 +454,13 @@ bool SANSRunWindow::loadUserFile()
m_uiForm.trans_opt->setCurrentIndex(index);
}
- //Monitor spectrum
+ //Monitor spectra
m_uiForm.monitor_spec->setText(runReduceScriptFunction("printParameter('MONITORSPECTRUM'),"));
+ m_uiForm.trans_monitor->setText(runReduceScriptFunction("printParameter('TRANS_UDET_MON'),"));
+ m_uiForm.monitor_interp->setChecked(
+ runReduceScriptFunction("printParameter('SAMP_INTERPOLATE'),") == "True");
+ m_uiForm.trans_interp->setChecked(
+ runReduceScriptFunction("printParameter('TRANS_INTERPOLATE'),") == "True");
//Direct efficiency correction
m_uiForm.direct_file->setText(runReduceScriptFunction("printParameter('DIRECT_BEAM_FILE_R'),"));
@@ -1458,7 +1463,13 @@ QString SANSRunWindow::createAnalysisDetailsScript(const QString & type)
//Sample offset
exec_reduce += "SetSampleOffset(" + m_uiForm.smpl_offset->text() + ")\n";
//Monitor spectrum
- exec_reduce += "SetMonitorSpectrum(" + m_uiForm.monitor_spec->text() + ")\n";
+ exec_reduce += "SetMonitorSpectrum(" + m_uiForm.monitor_spec->text() + ", ";
+ exec_reduce += m_uiForm.monitor_interp->isChecked() ? "True" : "False";
+ exec_reduce += ")\n";
+ //the monitor to normalise the tranmission spectrum against
+ exec_reduce += "SetTransSpectrum(" + m_uiForm.trans_monitor->text() + ", ";
+ exec_reduce += m_uiForm.trans_interp->isChecked() ? "True" : "False";
+ exec_reduce += ")\n";
//Extra mask information
addUserMaskStrings(exec_reduce);