Mantid
Loading...
Searching...
No Matches
ConvolveWorkspaces.cpp
Go to the documentation of this file.
1// Mantid Repository : https://github.com/mantidproject/mantid
2//
3// Copyright © 2018 ISIS Rutherford Appleton Laboratory UKRI,
4// NScD Oak Ridge National Laboratory, European Spallation Source,
5// Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
6// SPDX - License - Identifier: GPL - 3.0 +
7//----------------------------------------------------------------------
8// Includes
9//----------------------------------------------------------------------
14
15#include <gsl/gsl_errno.h>
16#include <gsl/gsl_fft_halfcomplex.h>
17#include <gsl/gsl_fft_real.h>
18
19#include <sstream>
20
22
23// Register the algorithm into the AlgorithmFactory
24DECLARE_ALGORITHM(ConvolveWorkspaces)
25
26using namespace Kernel;
27using namespace API;
28using namespace DataObjects;
29using namespace Geometry;
30using namespace Functions;
31
32void ConvolveWorkspaces::init() {
33 declareProperty(std::make_unique<WorkspaceProperty<Workspace2D>>("Workspace1", "", Direction::Input),
34 "The name of the resolution workspace");
35 declareProperty(std::make_unique<WorkspaceProperty<Workspace2D>>("Workspace2", "", Direction::Input),
36 "The name of the data workspace.");
37
38 declareProperty(std::make_unique<WorkspaceProperty<Workspace2D>>("OutputWorkspace", "", Direction::Output),
39 "The name of the output workspace.");
40}
41
42void ConvolveWorkspaces::exec() {
43 std::string ws1name = getProperty("Workspace1");
44 std::string ws2name = getProperty("Workspace2");
45 Workspace2D_sptr ws1 = getProperty("Workspace1");
46 Workspace2D_sptr ws2 = getProperty("Workspace2");
47
48 // Cache a few things for later use
49 const size_t numHists = ws1->getNumberHistograms();
50 Workspace2D_sptr outputWS = std::dynamic_pointer_cast<Workspace2D>(WorkspaceFactory::Instance().create(ws2));
51
52 // First check that the workspace are the same size
53 if (numHists != ws2->getNumberHistograms()) {
54 throw std::runtime_error("Size mismatch");
55 }
56
57 m_progress = std::make_unique<Progress>(this, 0.0, 1.0, numHists);
58 // Now convolve the histograms
59 PARALLEL_FOR_IF(Kernel::threadSafe(*ws1, *ws2, *outputWS))
60 for (int i = 0; i < static_cast<int>(numHists); ++i) {
62 m_progress->report();
63 outputWS->setSharedX(i, ws2->sharedX(i));
64 auto &Yout = outputWS->mutableY(i);
65 Convolution conv;
66
67 auto res = std::make_shared<TabulatedFunction>();
68 res->setAttributeValue("Workspace", ws1name);
69 res->setAttributeValue("WorkspaceIndex", i);
70
71 conv.addFunction(res);
72
73 auto fun = std::make_shared<TabulatedFunction>();
74 fun->setAttributeValue("Workspace", ws2name);
75 fun->setAttributeValue("WorkspaceIndex", i);
76
77 conv.addFunction(fun);
78 size_t N = Yout.size();
79 size_t Nx = outputWS->mutableX(i).size();
80 const double *firstX = &outputWS->mutableX(i)[0];
81 FunctionDomain1DView xView(firstX, Nx);
82 FunctionValues out(xView);
83 conv.function(xView, out);
84
85 for (size_t j = 0; j < N; j++) {
86 Yout[j] = out.getCalculated(j);
87 }
89 }
91 // Assign it to the output workspace property
92 setProperty("OutputWorkspace", outputWS);
93}
94
95} // namespace Mantid::CurveFitting::Algorithms
#define DECLARE_ALGORITHM(classname)
Definition: Algorithm.h:576
#define PARALLEL_START_INTERRUPT_REGION
Begins a block to skip processing is the algorithm has been interupted Note the end of the block if n...
Definition: MultiThreaded.h:94
#define PARALLEL_END_INTERRUPT_REGION
Ends a block to skip processing is the algorithm has been interupted Note the start of the block if n...
#define PARALLEL_FOR_IF(condition)
Empty definitions - to enable set your complier to enable openMP.
#define PARALLEL_CHECK_INTERRUPT_REGION
Adds a check after a Parallel region to see if it was interupted.
void declareProperty(std::unique_ptr< Kernel::Property > p, const std::string &doc="") override
Add a property to the list of managed properties.
Definition: Algorithm.cpp:1913
TypedValue getProperty(const std::string &name) const override
Get the value of a property.
Definition: Algorithm.cpp:2076
1D domain - a wrapper around an array of doubles.
A class to store values calculated by a function.
double getCalculated(size_t i) const
Get i-th calculated value.
A property class for workspaces.
Performes convolution of two functions.
Definition: Convolution.h:28
size_t addFunction(API::IFunction_sptr f) override
Add a function.
void function(const API::FunctionDomain &domain, API::FunctionValues &values) const override
Function you want to fit to.
IPropertyManager * setProperty(const std::string &name, const T &value)
Templated method to set the value of a PropertyWithValue.
static T & Instance()
Return a reference to the Singleton instance, creating it if it does not already exist Creation is do...
std::shared_ptr< Workspace2D > Workspace2D_sptr
shared pointer to Mantid::DataObjects::Workspace2D
std::unique_ptr< T > create(const P &parent, const IndexArg &indexArg, const HistArg &histArg)
This is the create() method that all the other create() methods call.
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.
Definition: MultiThreaded.h:22
@ Input
An input workspace.
Definition: Property.h:53
@ Output
An output workspace.
Definition: Property.h:54