Skip to content

Commit

Permalink
Allow ScaleX to take its factor from an instrument parameter.
Browse files Browse the repository at this point in the history 
Improved unit test coverage at the same time as there were no tests for
the EventWorkspace case.
Refs #8057
  • Loading branch information
@martyngigg
martyngigg committed Oct 2, 2013
1 parent 62ffe7b commit a53ac3a
Show file tree
Hide file tree
Showing 3 changed files with 316 additions and 98 deletions.
14 changes: 8 additions & 6 deletions Code/Mantid/Framework/Algorithms/inc/MantidAlgorithms/ScaleX.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -63,23 +63,25 @@ namespace Mantid
// Overridden Algorithm methods
void init();
void exec();

/// Execute algorithm for EventWorkspaces
void execEvent();

/// Create output workspace
API::MatrixWorkspace_sptr createOutputWS(API::MatrixWorkspace_sptr input);
API::MatrixWorkspace_sptr createOutputWS(const API::MatrixWorkspace_sptr & input);
/// Get the scale factor for the given spectrum
double getScaleFactor(const API::MatrixWorkspace_const_sptr & inputWS, const size_t index);

/// The progress reporting object
API::Progress *m_progress;

/// Scaling factor
double factor;
double m_factor;
/// instrument parameter name
std::string m_parname;
/// Start workspace index
int wi_min;
int m_wi_min;
/// Stop workspace index
int wi_max;

int m_wi_max;
};

} // namespace Algorithm
Expand Down
151 changes: 105 additions & 46 deletions Code/Mantid/Framework/Algorithms/src/ScaleX.cpp
View file
Original file line number Diff line number Diff line change
@@ -1,7 +1,8 @@
/*WIKI*
Scales the X axis of the input workspace by the amount requested.
Scales the X axis of the input workspace by the amount requested. The amount can be specified either as:
* an absolute numerical value via the "Factor" argument or
* an detector parameter name whose value is retrieved from the instrument.
*WIKI*/
//----------------------------------------------------------------------
Expand All @@ -12,6 +13,8 @@ Scales the X axis of the input workspace by the amount requested.
#include "MantidKernel/BoundedValidator.h"
#include "MantidKernel/ListValidator.h"

#include <boost/function.hpp>

namespace Mantid
{
namespace Algorithms
Expand Down Expand Up @@ -63,6 +66,8 @@ void ScaleX::init()
mustBePositive->setLower(0);
declareProperty("IndexMin", 0, mustBePositive, "The workspace index of the first spectrum to scale. Only used if IndexMax is set.");
declareProperty("IndexMax", Mantid::EMPTY_INT(), mustBePositive, "The workspace index of the last spectrum to scale. Only used if explicitly set.");
//Add InstrumentParameter property here so as not to mess with the parameter order for current scripts
declareProperty("InstrumentParameter", "", "The name of an instrument parameter whose value is used to scale as the input factor");
}

/**
Expand All @@ -72,24 +77,26 @@ void ScaleX::exec()
{
//Get input workspace and offset
const MatrixWorkspace_sptr inputW = getProperty("InputWorkspace");
factor = getProperty("Factor");
m_factor = getProperty("Factor");
m_parname = getPropertyValue("InstrumentParameter");

const std::string op = getPropertyValue("Operation");
API::MatrixWorkspace_sptr outputW = createOutputWS(inputW);
//Get number of histograms
int histnumber = static_cast<int>(inputW->getNumberHistograms());
m_progress = new API::Progress(this, 0.0, 1.0, histnumber+1);
m_progress->report("Scaling X");
wi_min = 0;
wi_max = histnumber-1;
m_wi_min = 0;
m_wi_max = histnumber-1;
//check if workspace indexes have been set
int tempwi_min = getProperty("IndexMin");
int tempwi_max = getProperty("IndexMax");
if ( tempwi_max != Mantid::EMPTY_INT() )
{
if ((wi_min <= tempwi_min) && (tempwi_min <= tempwi_max) && (tempwi_max <= wi_max))
if ((m_wi_min <= tempwi_min) && (tempwi_min <= tempwi_max) && (tempwi_max <= m_wi_max))
{
wi_min = tempwi_min;
wi_max = tempwi_max;
m_wi_min = tempwi_min;
m_wi_max = tempwi_max;
}
else
{
Expand All @@ -104,42 +111,50 @@ void ScaleX::exec()
this->execEvent();
return;
}

const bool multiply = (op=="Multiply");
boost::function<double (double x,double factor)> scale;
if(multiply) scale = std::multiplies<double>();
else scale = std::plus<double>();

// do the shift in X
PARALLEL_FOR2(inputW, outputW)
for (int i=0; i < histnumber; ++i)
for (int i = 0; i < histnumber; ++i)
{
PARALLEL_START_INTERUPT_REGION
//Do the offsetting
for (int j=0; j < static_cast<int>(inputW->readX(i).size()); ++j)

//Copy y and e data
auto & outY = outputW->dataY(i);
outY = inputW->dataY(i);
auto & outE = outputW->dataE(i);
outE = inputW->dataE(i);

auto & outX = outputW->dataX(i);
const auto & inX = inputW->readX(i);
//Change bin value by offset
if ((i >= m_wi_min) && (i <= m_wi_max))
{
//Change bin value by offset
if ((i >= wi_min) && (i <= wi_max))
double factor = getScaleFactor(inputW, i);
// Do the offsetting
std::transform(inX.begin(), inX.end(), outX.begin(), std::bind2nd(scale, factor));
// reverse the vector if multiplicative factor was negative
if(multiply && m_factor < 0.0)
{
if(op=="Multiply")
{
outputW->dataX(i)[j] = inputW->readX(i)[j] * factor;
}
else if(op=="Add")
{
outputW->dataX(i)[j] = inputW->readX(i)[j] + factor;
}
std::reverse( outX.begin(), outX.end() );
std::reverse( outY.begin(), outY.end() );
std::reverse( outE.begin(), outE.end() );
}
else outputW->dataX(i)[j] = inputW->readX(i)[j];
}
//Copy y and e data
outputW->dataY(i) = inputW->dataY(i);
outputW->dataE(i) = inputW->dataE(i);
// reverse the vector if multiplicative factor was negative
if( (i >= wi_min) && (i <= wi_max) && op=="Multiply" && factor<0 )
else
{
std::reverse( outputW->dataX(i).begin(), outputW->dataX(i).end() );
std::reverse( outputW->dataY(i).begin(), outputW->dataY(i).end() );
std::reverse( outputW->dataE(i).begin(), outputW->dataE(i).end() );
outX = inX; //copy
}
m_progress->report("Scaling X");

PARALLEL_END_INTERUPT_REGION
}
PARALLEL_CHECK_INTERUPT_REGION

// Copy units
if (outputW->getAxis(0)->unit().get())
outputW->getAxis(0)->unit() = inputW->getAxis(0)->unit();
Expand All @@ -156,18 +171,6 @@ void ScaleX::exec()
setProperty("OutputWorkspace",outputW);
}

API::MatrixWorkspace_sptr ScaleX::createOutputWS(API::MatrixWorkspace_sptr input)
{
//Check whether input = output to see whether a new workspace is required.
MatrixWorkspace_sptr output = getProperty("OutputWorkspace");
if ( input != output )
{
//Create new workspace for output from old
output = API::WorkspaceFactory::Instance().create(input);
}
return output;
}

void ScaleX::execEvent()
{
g_log.information("Processing event workspace");
Expand Down Expand Up @@ -197,19 +200,19 @@ void ScaleX::execEvent()
{
PARALLEL_START_INTERUPT_REGION
//Do the offsetting
if ((i >= wi_min) && (i <= wi_max))
if ((i >= m_wi_min) && (i <= m_wi_max))
{
if(op=="Multiply")
{
outputWS->getEventList(i).scaleTof(factor);
if( factor < 0 )
outputWS->getEventList(i).scaleTof(getScaleFactor(inputWS, i));
if( m_factor < 0 )
{
outputWS->getEventList(i).reverse();
}
}
else if(op=="Add")
{
outputWS->getEventList(i).addTof(factor);
outputWS->getEventList(i).addTof(getScaleFactor(inputWS, i));
}
}
m_progress->report("Scaling X");
Expand All @@ -219,6 +222,62 @@ void ScaleX::execEvent()
outputWS->clearMRU();
}

API::MatrixWorkspace_sptr ScaleX::createOutputWS(const API::MatrixWorkspace_sptr & input)
{
//Check whether input = output to see whether a new workspace is required.
MatrixWorkspace_sptr output = getProperty("OutputWorkspace");
if ( input != output )
{
//Create new workspace for output from old
output = API::WorkspaceFactory::Instance().create(input);
}
return output;
}


/**
* If the InstrumentParameter property is set then it attempts to retrieve the parameter
* from the component, else it returns the value of the Factor property
* @param inputWS A pointer to the input workspace
* @param index The current index to inspect
* @return Value for the scale factor
*/
double ScaleX::getScaleFactor(const API::MatrixWorkspace_const_sptr & inputWS, const size_t index)
{
if(m_parname.empty()) return m_factor;

// Try and get factor from component. If we see a DetectorGroup use this will use the first component
Geometry::IDetector_const_sptr det;
auto inst = inputWS->getInstrument();

auto *spec = inputWS->getSpectrum(index);
const auto & ids = spec->getDetectorIDs();
const size_t ndets(ids.size());
if(ndets > 0)
{
try
{
det = inst->getDetector(*ids.begin());
}
catch(Exception::NotFoundError&)
{
return 0.0;
}
}
else return 0.0;

const auto & pmap = inputWS->constInstrumentParameters();
auto par = pmap.getRecursive(det->getComponentID(), m_parname);
if(par) return par->value<double>();
else
{
std::ostringstream os;
os << "Spectrum at index '" << index << "' has no parameter named '" << m_parname << "'\n";
throw std::runtime_error(os.str());
}
}


} // namespace Algorithm
} // namespace Mantid

Expand Down

0 comments on commit a53ac3a

Please sign in to comment.