22#include "MantidHistogramData/LinearGenerator.h"
51 "An input workspace.");
53 "An output workspace.");
55 auto twoOrThreeElements = std::make_shared<ArrayLengthValidator<double>>(2, 3);
56 std::vector<double> myInput(3, 0);
58 "Coordinate of the centre of the ring");
60 auto nonNegative = std::make_shared<BoundedValidator<double>>();
61 nonNegative->setLower(0);
62 declareProperty(
"MinRadius", 0.0, nonNegative->clone(),
"Length of the inner ring. Default=0");
64 std::move(nonNegative)),
65 "Length of the outer ring. Default=ImageSize.");
67 auto nonNegativeInt = std::make_shared<BoundedValidator<int>>();
68 nonNegativeInt->setLower(1);
69 declareProperty(
"NumBins", 100, std::move(nonNegativeInt),
"Number of slice bins for the output. Default=100");
71 const char *normBy =
"NormalizeByRadius";
72 const char *normOrder =
"NormalizationOrder";
75 "Divide the sum of each ring by the radius "
76 "powered by Normalization Order");
78 "If 2, the normalization will be divided by "
79 "the quadratic value of the ring for each "
83 const char *groupNorm =
"Normalization";
95 std::vector<double> output;
108 g_log.
debug() <<
"Process Instrument related image\n";
110 std::vector<double> accumulator(
num_bins, 0);
112 g_log.
debug() <<
"For every detector in the image get its position "
113 <<
" and sum up all the counts inside the related spectrum\n";
115 const auto &spectrumInfo =
inputWS->spectrumInfo();
116 for (
size_t i = 0; i <
inputWS->getNumberHistograms(); i++) {
117 if (!spectrumInfo.hasDetectors(i)) {
123 if (spectrumInfo.isMonitor(i))
131 auto &refY =
inputWS->getSpectrum(i).y();
132 accumulator[bin_n] += std::accumulate(refY.cbegin(), refY.cend(), 0.0);
140 std::vector<double> accumulator(
num_bins, 0);
153 g_log.
debug() <<
"Define the X positions of the pixels\n";
154 auto x_pos =
inputWS->points(0);
156 g_log.
debug() <<
"For every pixel define its bin position and sum them up\n";
158 for (
size_t i = 0; i <
inputWS->getNumberHistograms(); i++) {
161 auto &refY =
inputWS->getSpectrum(i).y();
164 for (
size_t j = 0; j < refY.size(); j++) {
168 V3D pixelPos(x_pos[j], (*verticalAxis)(i), 0);
175 accumulator[bin_pos] += refY[j];
183 double radio, theta, phi;
189 double effect_distance = radio * sin(theta * M_PI / 180);
191 if (effect_distance < min_radius || effect_distance >
max_radius)
198 g_log.
debug() <<
"Copy the input property values\n";
201 g_log.
debug() <<
"Extract the center and make it a V3D object\n";
202 std::vector<double> centre_aux =
getProperty(
"Centre");
203 if (centre_aux.size() == 2)
204 centre =
V3D(centre_aux[0], centre_aux[1], 0);
206 centre =
V3D(centre_aux[0], centre_aux[1], centre_aux[2]);
208 g_log.
debug() <<
"Copy the remaning properties: MinRadius, MaxRadius and NumBins\n";
217 g_log.
debug() <<
"Check MinRadius < MaxRadius\n";
220 s <<
"Wrong definition of the radius min and max. The minimum radius can "
221 "not be bigger than maximum. "
223 throw std::invalid_argument(s.str());
228 std::copy(boundary_limits.begin(), std::prev(boundary_limits.end()), std::ostream_iterator<double>(s,
", "));
229 s << boundary_limits.back();
232 g_log.
debug() <<
"Check: centre is defined inside the region defined by the "
233 "image or instrument\n";
236 g_log.
debug() <<
"Recalculate MaxRadius if default value is given\n";
237 if (
max_radius > 0.9 * std::numeric_limits<double>::max()) {
242 g_log.
debug() <<
"Check number of bins to alert user if many bins will end up empty\n";
260 return inWS->getAxis(1)->isSpectra();
293 const auto &refX = inWS->x(inWS->getNumberHistograms() / 2);
297 const double &first_x = refX.front();
298 const double &last_x = refX.back();
299 if (first_x < last_x) {
310 throw std::invalid_argument(
"Vertical axis is not a numeric axis. Can not "
311 "find the limits of the image.");
314 min_y = verticalAxis->
getMin();
315 max_y = verticalAxis->
getMax();
320 throw std::logic_error(
"Failure to get the boundaries of this image. "
321 "Internal logic error. Please, inform MantidHelp");
324 std::vector<double> output(4);
378 double first_x, first_y, first_z;
380 const auto &spectrumInfo = inWS->spectrumInfo();
383 if (i >= inWS->getNumberHistograms())
384 throw std::invalid_argument(
"Did not find any non monitor detector. "
385 "Failed to identify the boundaries of this "
388 if (spectrumInfo.isMonitor(i))
391 const auto pos = spectrumInfo.position(i);
399 double last_x, last_y, last_z;
400 i = inWS->getNumberHistograms() - 1;
404 throw std::invalid_argument(
"There is no region defined for the "
405 "instrument of this workspace. Failed to "
406 "identify the boundaries of this instrument");
408 if (spectrumInfo.isMonitor(i))
411 const auto pos = spectrumInfo.position(i);
420 double xMax, yMax, zMax;
421 double xMin, yMin, zMin;
422 xMax = std::max(first_x, last_x);
423 yMax = std::max(first_y, last_y);
424 zMax = std::max(first_z, last_z);
425 xMin = std::min(first_x, last_x);
426 yMin = std::min(first_y, last_y);
427 zMin = std::min(first_z, last_z);
429 std::vector<double> output(6);
464 if ((2 *
centre.size()) != boundaries.size())
465 throw std::invalid_argument(
"CenterIsInsideLimits: The centre and boundaries were not given in the "
466 "correct as {x1,x2,...} and {x1_min, x1_max, x2_min, x2_max, ... }");
468 bool check_passed =
true;
471 for (
size_t i = 0; i < 2; i++) {
472 if (
centre[i] < boundaries[2 * i] ||
centre[i] > boundaries[2 * i + 1]) {
473 check_passed =
false;
474 s <<
"The position for axis " << i + 1 <<
" (" <<
centre[i] <<
") is outside the limits [" << boundaries[2 * i]
475 <<
", " << boundaries[2 * i + 1] <<
"]. \n";
480 throw std::invalid_argument(s.str());
492 if (size_bins < min_bin_size)
493 g_log.
warning() <<
"The number of bins provided is too big. "
494 <<
"It corresponds to a separation smaller than the image "
495 "resolution (detector size). "
496 <<
"A resonable number is smaller than "
505 const auto &spectrumInfo = inWS->spectrumInfo();
506 for (
size_t i = 0; i < inWS->getNumberHistograms(); ++i) {
507 if (spectrumInfo.isMonitor(i))
510 spectrumInfo.detector(i).getBoundingBox(bbox);
514 throw std::invalid_argument(
"Did not find any non monitor detector position");
524 const auto &refX = inWS->x(
inputWS->getNumberHistograms() / 2);
525 auto nX =
static_cast<int>(refX.size());
526 int nY =
static_cast<int>(inWS->getAxis(1)->length());
529 return std::min(((boundaries[1] - boundaries[0]) / nX), ((boundaries[3] - boundaries[2]) / nY));
533 g_log.
debug() <<
"Normalization of the output in relation to the 'radius' (distance)\n";
538 double first_radius =
min_radius + bin_size / 2;
540 g_log.
debug() <<
"Calculate Output[i] = Counts[i] / (Radius[i] ^ " << exp_power <<
") << \n";
541 if (exp_power > 1.00001 || exp_power < 0.99999) {
542 for (
int i = 0; i < static_cast<int>(values.size()); i++) {
543 values[i] = values[i] / pow(first_radius + i * bin_size, exp_power);
546 for (
int i = 0; i < static_cast<int>(values.size()); i++) {
547 values[i] = values[i] / (first_radius + i * bin_size);
554 std::array<double, 2> Xs = {{boundary_limits[0], boundary_limits[1]}};
555 std::array<double, 2> Ys = {{boundary_limits[2], boundary_limits[3]}};
556 std::array<double, 2> Zs = {{0., 0.}};
558 if (boundary_limits.size() == 6) {
559 Zs[0] = boundary_limits[4];
560 Zs[1] = boundary_limits[5];
563 double max_distance = 0;
571 if (curr_distance > max_distance)
572 max_distance = curr_distance;
579 g_log.
debug() <<
"Output calculated, setting up the output workspace\n";
585 outputWS->mutableY(0) = values;
589 auto xSize = outputWS->mutableX(0).size();
591 outputWS->setBinEdges(0, xSize, HistogramData::LinearGenerator(
min_radius, bin_size));
600 auto horizontal = std::make_unique<API::NumericAxis>(xSize);
602 std::dynamic_pointer_cast<Units::Label>(labelX)->setLabel(
"Radius");
603 horizontal->unit() = labelX;
604 outputWS->replaceAxis(0, std::move(horizontal));
#define DECLARE_ALGORITHM(classname)
void declareProperty(std::unique_ptr< Kernel::Property > p, const std::string &doc="") override
Add a property to the list of managed properties.
TypedValue getProperty(const std::string &name) const override
Get the value of a property.
Class to represent a numeric axis of a workspace.
double getMax() const override
returns max value defined on axis
double getMin() const override
returns min value defined on axis
A property class for workspaces.
int fromDistanceToBin(double distance)
Return the bin position for a given distance From the input, it is defined the limits of distances as...
double getMaxDistance(const Kernel::V3D ¢re, const std::vector< double > &boundary_limits)
static void centerIsInsideLimits(const std::vector< double > ¢re, const std::vector< double > &boundaries)
Check if a given position is inside the limits defined by the boundaries.
std::vector< double > processNumericImageRadiusSum()
void normalizeOutputByRadius(std::vector< double > &values, double exp_power)
static bool inputWorkspaceHasInstrumentAssociated(const API::MatrixWorkspace_sptr &)
Differentiate between Instrument related image (where the position of the pixels depend on the instru...
void init() override
Initialize the algorithm's properties.
double getMinBinSizeForInstrument(const API::MatrixWorkspace_sptr &)
void cacheInputPropertyValues()
API::MatrixWorkspace_sptr inputWS
double getMinBinSizeForNumericImage(const API::MatrixWorkspace_sptr &)
const std::string category() const override
Algorithm's category for identification.
int version() const override
Algorithm's version for identification.
void numBinsIsReasonable()
static std::vector< double > getBoundariesOfNumericImage(const API::MatrixWorkspace_sptr &)
Assuming that the input workspace is a Numeric Image where the pixel positions depend on their relati...
void setUpOutputWorkspace(const std::vector< double > &values)
std::vector< double > processInstrumentRadiusSum()
void exec() override
Execute the algorithm.
int getBinForPixelPos(const Kernel::V3D &pos)
static std::vector< double > getBoundariesOfInstrument(const API::MatrixWorkspace_sptr &)
Assuming that the workspace has an instrument associated with it from which the pixel positions has t...
std::vector< double > getBoundariesOfInputWorkspace()
void inputValidationSanityCheck()
A simple structure that defines an axis-aligned cuboid shaped bounding box for a geometrical object.
Kernel::V3D width() const
Returns the width of the box.
Support for a property that holds an array of values.
IPropertyManager * setProperty(const std::string &name, const T &value)
Templated method to set the value of a PropertyWithValue.
void setPropertySettings(const std::string &name, std::unique_ptr< IPropertySettings > settings)
void setPropertyGroup(const std::string &name, const std::string &group)
Set the group for a given property.
void debug(const std::string &msg)
Logs at debug level.
void notice(const std::string &msg)
Logs at notice level.
void warning(const std::string &msg)
Logs at warning level.
void information(const std::string &msg)
Logs at information level.
The concrete, templated class for properties.
static T & Instance()
Return a reference to the Singleton instance, creating it if it does not already exist Creation is do...
double distance(const V3D &v) const noexcept
Calculates the distance between two vectors.
double norm() const noexcept
void getSpherical(double &R, double &theta, double &phi) const noexcept
Return the vector's position in spherical coordinates.
std::shared_ptr< MatrixWorkspace > MatrixWorkspace_sptr
shared pointer to the matrix workspace base class
@ Input
An input workspace.
@ Output
An output workspace.