21using namespace Kernel;
23using namespace Geometry;
24using namespace DataObjects;
36static void applyBadPixelThreshold(
MatrixWorkspace &outputWS,
double minThreshold,
double maxThreshold) {
42 for (
size_t i = 0; i < numberOfSpectra; i++) {
47 if (spectrumInfo.isMonitor(i)) {
51 }
else if (spectrumInfo.isMasked(i)) {
57 const auto y = YOut.front();
58 if (y < minThreshold || y > maxThreshold) {
71 "The workspace containing the flood data");
74 "The name of the workspace to be created as the output of the algorithm");
76 auto positiveDouble = std::make_shared<BoundedValidator<double>>();
77 positiveDouble->setLower(0);
80 "Maximum threshold for a pixel to be considered");
83 "Whether to merge entries when WorkspaceGroup is specified as input.");
87 std::map<std::string, std::string> result;
92 auto oneBinMsg =
"Input workspace must have only one bin. Consider "
93 "integrating the input over all the bins.";
96 if (inputWS ==
nullptr) {
98 if (inputGroup !=
nullptr) {
99 for (
auto entryNo = 0; entryNo < inputGroup->getNumberOfEntries(); ++entryNo) {
100 auto const entry = std::static_pointer_cast<const MatrixWorkspace>(inputGroup->getItem(entryNo));
101 if (entry->blocksize() > 1) {
108 "WorkspaceGroup containing MatrixWorkspaces";
110 }
else if (inputWS->blocksize() > 1) {
122 const std::string msg{
"MinThreshold must be less than MaxThreshold"};
135 auto inputWS = std::static_pointer_cast<MatrixWorkspace>(ws1);
142 double startProgress,
double stepProgress) {
145 if (std::dynamic_pointer_cast<EventWorkspace>(inputWorkspace)) {
148 childAlg->setProperty(
"WorkspaceToRebin", inputWorkspace);
149 childAlg->setProperty(
"WorkspaceToMatch", inputWorkspace);
151 childAlg->setProperty(
"PreserveEvents",
false);
152 childAlg->executeAsChildAlg();
153 outputWS = childAlg->getProperty(
"OutputWorkspace");
156 outputWS = inputWorkspace->clone();
162 progress(startProgress + 0.1 * stepProgress,
"Computing the counts.");
164 if (counts.nPixels == 0) {
165 throw std::runtime_error(
"No pixels being used for calculation");
168 progress(startProgress + 0.3 * stepProgress,
"Normalising the detectors.");
171 progress(startProgress + 0.5 * stepProgress,
"Applying bad pixel threshold.");
175 progress(startProgress + 0.7 * stepProgress,
"Computing the counts.");
177 if (counts.nPixels == 0) {
178 throw std::runtime_error(
"All pixels are outside of the threshold values");
181 progress(startProgress + 0.9 * stepProgress,
"Normalising the detectors.");
194 for (
const auto &result : results) {
195 throw std::runtime_error(
"Issue in " + result.first +
" property: " + result.second);
219 auto outputGroup = std::make_shared<WorkspaceGroup>();
220 auto const nEntries = inputWS->getNumberOfEntries();
221 auto const stepProgress = 1.0 / nEntries;
222 for (
auto entryNo = 0; entryNo < nEntries; ++entryNo) {
223 auto entryWS = std::static_pointer_cast<API::MatrixWorkspace>(inputWS->getItem(entryNo));
224 auto const startProgress =
static_cast<double>(entryNo) / nEntries;
226 outputGroup->addWorkspace(outputWS);
243 const size_t numberOfSpectra =
workspace.getNumberHistograms();
244 SummedResults results;
246 const auto &spectrumInfo =
workspace.spectrumInfo();
247 for (
size_t i = 0; i < numberOfSpectra; i++) {
255 if (spectrumInfo.isMonitor(i) || spectrumInfo.isMasked(i) ||
isEmpty(YValues.front()))
258 results.sum += YValues.front();
259 results.error += YErrors.front() * YErrors.front();
262 results.error = std::sqrt(results.error);
264 g_log.
debug() <<
"Total of unmasked/dead pixels = " << results.nPixels <<
" from a total of " << numberOfSpectra
273 const size_t numberOfSpectra =
workspace.getNumberHistograms();
274 const auto &spectrumInfo =
workspace.spectrumInfo();
276 const double averageY = counts.sum /
static_cast<double>(counts.nPixels);
277 const double averageE = counts.error /
static_cast<double>(counts.nPixels);
279 for (
size_t i = 0; i < numberOfSpectra; i++) {
283 const auto yOriginal =
y.front();
287 if (spectrumInfo.isMasked(i) || spectrumInfo.isMonitor(i) ||
isEmpty(yOriginal))
290 const auto eOriginal = e.front();
293 y.front() = yOriginal / averageY;
294 const double signalToNoiseOrig = eOriginal / yOriginal;
295 const double signalToNoiseAvg = averageE / averageY;
296 e.front() =
y.front() * std::sqrt((signalToNoiseOrig * signalToNoiseOrig) + (signalToNoiseAvg * signalToNoiseAvg));
299 g_log.
debug() <<
"Averages :: counts = " << averageY <<
"; error = " << averageE <<
"\n";
314 mergeRuns->setProperty(
"InputWorkspaces", input.
getName());
315 mergeRuns->executeAsChildAlg();
316 Workspace_sptr mergedWs = mergeRuns->getProperty(
"OutputWorkspace");
319 createSingleAlg->setProperty(
"DataValue",
static_cast<double>(nEntries));
320 createSingleAlg->executeAsChildAlg();
324 divideAlg->setProperty(
"LHSWorkspace", mergedWs);
325 divideAlg->setProperty(
"RHSWorkspace", normWs);
326 divideAlg->executeAsChildAlg();
329 auto spectrumInfo = mergedNormalisedWs->spectrumInfo();
331 for (
auto spectrumNo = 0; spectrumNo < static_cast<int>(mergedNormalisedWs->getNumberHistograms()); ++spectrumNo) {
332 if (spectrumInfo.isMasked(spectrumNo)) {
333 auto &detDataY = mergedNormalisedWs->mutableY(spectrumNo);
334 auto &detDataErr = mergedNormalisedWs->mutableE(spectrumNo);
337 auto nonMaskedEntries = 0;
338 for (
auto entryNo = 0; entryNo < nEntries; ++entryNo) {
340 auto spectrumInfoEntry = entry->spectrumInfo();
341 if (!spectrumInfoEntry.isMasked(spectrumNo)) {
342 dataY += entry->readY(spectrumNo)[0];
343 dataE += pow(entry->readE(spectrumNo)[0], 2);
347 if (nonMaskedEntries != 0) {
349 spectrumInfo.setMasked(spectrumNo,
false);
350 detDataY.front() = dataY / nonMaskedEntries;
351 detDataErr.front() = sqrt(dataE) / nonMaskedEntries;
355 return mergedNormalisedWs;
#define DECLARE_ALGORITHM(classname)
IPeaksWorkspace_sptr workspace
#define PARALLEL_FOR_IF(condition)
Empty definitions - to enable set your complier to enable openMP.
void declareProperty(std::unique_ptr< Kernel::Property > p, const std::string &doc="") override
Add a property to the list of managed properties.
std::string getPropertyValue(const std::string &name) const override
Get the value of a property as a string.
TypedValue getProperty(const std::string &name) const override
Get the value of a property.
virtual std::shared_ptr< Algorithm > createChildAlgorithm(const std::string &name, const double startProgress=-1., const double endProgress=-1., const bool enableLogging=true, const int &version=-1)
Create a Child Algorithm.
void progress(double p, const std::string &msg="", double estimatedTime=0.0, int progressPrecision=0)
Sends ProgressNotification.
static bool isEmpty(const NumT toCheck)
checks that the value was not set by users, uses the value in empty double/int.
const SpectrumInfo & spectrumInfo() const
Return a reference to the SpectrumInfo object.
Base MatrixWorkspace Abstract Class.
virtual std::size_t getNumberHistograms() const =0
Returns the number of histograms in the workspace.
HistogramData::HistogramE & mutableE(const size_t index) &
HistogramData::HistogramY & mutableY(const size_t index) &
Class to hold a set of workspaces.
int getNumberOfEntries() const
Return the number of entries within the group.
Workspace_sptr getItem(const size_t index) const
Return the ith workspace.
A property class for workspaces.
const std::string & getName() const override
Get the workspace name.
void averageAndNormalizePixels(API::MatrixWorkspace &workspace, const SummedResults &results)
void exec() override
Executes the algorithm.
double m_maxThreshold
Maximum efficiency. Pixels with higher efficiency will be masked.
std::map< std::string, std::string > validateInputs() override
Method checking errors on ALL the inputs, before execution.
API::MatrixWorkspace_sptr mergeGroup(API::WorkspaceGroup &)
Merges the input workspace group and tries to remove masked spectra and replace the masked value with...
bool processGroups() override
Process groups and merge entries if MergeGroup property is set.
double m_minThreshold
Minimum efficiency. Pixels with lower efficiency will be masked.
SummedResults sumUnmaskedAndDeadPixels(const API::MatrixWorkspace &workspace)
Sum all detectors, excluding monitors and masked detectors.
void validateGroupInput()
Explicitly calls validateInputs and throws runtime error in case of issues in the input properties.
API::MatrixWorkspace_sptr calculateEfficiency(const API::MatrixWorkspace_sptr &, double startProgress=0.0, double stepProgress=1.0)
void init() override
Initialization method.
IPropertyManager * setProperty(const std::string &name, const T &value)
Templated method to set the value of a PropertyWithValue.
void debug(const std::string &msg)
Logs at debug level.
static T & Instance()
Return a reference to the Singleton instance, creating it if it does not already exist Creation is do...
std::shared_ptr< WorkspaceGroup > WorkspaceGroup_sptr
shared pointer to Mantid::API::WorkspaceGroup
std::shared_ptr< Workspace > Workspace_sptr
shared pointer to Mantid::API::Workspace
std::shared_ptr< const Workspace > Workspace_const_sptr
shared pointer to Mantid::API::Workspace (const version)
std::shared_ptr< const MatrixWorkspace > MatrixWorkspace_const_sptr
shared pointer to the matrix workspace base class (const version)
std::shared_ptr< const WorkspaceGroup > WorkspaceGroup_const_sptr
shared pointer to Mantid::API::WorkspaceGroup, pointer to const version
std::shared_ptr< MatrixWorkspace > MatrixWorkspace_sptr
shared pointer to the matrix workspace base class
static const std::string OUTPUT_WORKSPACE
static const std::string MERGE_GROUP
static const std::string INPUT_WORKSPACE
static const std::string MIN_THRESHOLD
static const std::string MAX_THRESHOLD
std::enable_if< std::is_pointer< Arg >::value, bool >::type threadSafe(Arg workspace)
Thread-safety check Checks the workspace to ensure it is suitable for multithreaded access.
constexpr double EMPTY_DBL() noexcept
Returns what we consider an "empty" double within a property.
@ Input
An input workspace.
@ Output
An output workspace.