Mantid
Loading...
Searching...
No Matches
CarpenterSampleCorrection.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 +
12
13#include <stdexcept>
14
15namespace Mantid::Algorithms {
16DECLARE_ALGORITHM(CarpenterSampleCorrection) // Register the class
17 // into the algorithm
18 // factory
19
20using namespace Kernel;
21using namespace API;
24using namespace Geometry;
25
26const std::string CarpenterSampleCorrection::name() const { return "CarpenterSampleCorrection"; }
27
28int CarpenterSampleCorrection::version() const { return 1; }
29
30const std::string CarpenterSampleCorrection::category() const { return "CorrectionFunctions\\AbsorptionCorrections"; }
31
36 // The input workspace must have an instrument and units of wavelength
37 auto wsValidator = std::make_shared<CompositeValidator>();
38 wsValidator->add<WorkspaceUnitValidator>("Wavelength");
39 wsValidator->add<InstrumentValidator>();
40
41 auto algCalcCarpenter = AlgorithmManager::Instance().createUnmanaged("CalculateCarpenterSampleCorrection");
42 algCalcCarpenter->initialize();
43
44 declareProperty(
45 std::make_unique<WorkspaceProperty<API::MatrixWorkspace>>("InputWorkspace", "", Direction::Input, wsValidator),
46 "The name of the input workspace.");
47 declareProperty(std::make_unique<WorkspaceProperty<API::MatrixWorkspace>>("OutputWorkspace", "", Direction::Output),
48 "The name of the output workspace.");
49
50 copyProperty(algCalcCarpenter, "AttenuationXSection");
51 copyProperty(algCalcCarpenter, "ScatteringXSection");
52 copyProperty(algCalcCarpenter, "SampleNumberDensity");
53 copyProperty(algCalcCarpenter, "CylinderSampleRadius");
54}
55
60 // common information
61 MatrixWorkspace_sptr inputWksp = getProperty("InputWorkspace");
62 double radius = getProperty("CylinderSampleRadius");
63 double coeff1 = getProperty("AttenuationXSection");
64 double coeff2 = getProperty("SampleNumberDensity");
65 double coeff3 = getProperty("ScatteringXSection");
66
67 // Calculate the absorption and multiple scattering corrections
68 WorkspaceGroup_sptr calcOutput = calculateCorrection(inputWksp, radius, coeff1, coeff2, coeff3, true, true);
69 Workspace_sptr absPtr = calcOutput->getItem(0);
70 Workspace_sptr msPtr = calcOutput->getItem(1);
71 auto absWksp = std::dynamic_pointer_cast<MatrixWorkspace>(absPtr);
72 auto msWksp = std::dynamic_pointer_cast<MatrixWorkspace>(msPtr);
73
74 EventWorkspace_sptr inputWkspEvent = std::dynamic_pointer_cast<EventWorkspace>(inputWksp);
75
76 // Inverse the absorption correction ( 1/A)
77 const auto NUM_HIST = static_cast<int64_t>(inputWksp->getNumberHistograms());
79 for (int64_t i = 0; i < NUM_HIST; ++i) {
81 auto &y = absWksp->mutableY(i);
82 for (size_t j = 0; j < y.size(); j++) {
83 y[j] = 1.0 / y[j];
84 }
86 }
88
89 // Compute the overall correction (= 1/A - MS ) to multiply by
90 auto correctionWksp = minus(absWksp, msWksp);
91
92 // Apply the correction to the sample workspace
93 // = (1/A - MS) * wksp
94 // = wksp/A - MS * wksp
95 MatrixWorkspace_sptr outputWksp = multiply(inputWksp, correctionWksp);
96
97 // Output workspace
98 if (inputWkspEvent) {
99 auto outputWkspEvent = std::dynamic_pointer_cast<EventWorkspace>(outputWksp);
100 setProperty("OutputWorkspace", outputWkspEvent);
101 }
102 setProperty("OutputWorkspace", outputWksp);
103}
104
106 double coeff1, double coeff2, double coeff3,
107 bool doAbs, bool doMS) {
108 auto calculate = this->createChildAlgorithm("CalculateCarpenterSampleCorrection", 0.0, 0.25);
109 calculate->setProperty("InputWorkspace", inputWksp);
110 calculate->setProperty("CylinderSampleRadius", radius);
111 calculate->setProperty("AttenuationXSection", coeff1);
112 calculate->setProperty("SampleNumberDensity", coeff2);
113 calculate->setProperty("ScatteringXSection", coeff3);
114 calculate->setProperty("Absorption", doAbs);
115 calculate->setProperty("MultipleScattering", doMS);
116 calculate->execute();
117 WorkspaceGroup_sptr calcOutput = calculate->getProperty("OutputWorkspaceBaseName");
118 return calcOutput;
119}
120
122 const MatrixWorkspace_sptr &rhsWS) {
123 auto minus = this->createChildAlgorithm("Minus", 0.5, 0.75);
124 minus->setProperty("LHSWorkspace", lhsWS);
125 minus->setProperty("RHSWorkspace", rhsWS);
126 minus->execute();
127 MatrixWorkspace_sptr outWS = minus->getProperty("OutputWorkspace");
128 return outWS;
129}
130
132 const MatrixWorkspace_sptr &rhsWS) {
133 auto multiply = this->createChildAlgorithm("Multiply", 0.75, 1.0);
134 multiply->setProperty("LHSWorkspace", lhsWS);
135 multiply->setProperty("RHSWorkspace", rhsWS);
136 multiply->execute();
137 MatrixWorkspace_sptr outWS = multiply->getProperty("OutputWorkspace");
138 return outWS;
139}
140
141} // namespace Mantid::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.
double radius
Definition: Rasterize.cpp:31
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) override
Create a Child Algorithm.
void copyProperty(const API::Algorithm_sptr &alg, const std::string &name)
Copy a property from an existing algorithm.
Kernel::IPropertyManager::TypedValue getProperty(const std::string &name) const override
Get the property held by this object.
A validator which checks that a workspace has a valid instrument.
A property class for workspaces.
A validator which checks that the unit of the workspace referred to by a WorkspaceProperty is the exp...
API::MatrixWorkspace_sptr multiply(const API::MatrixWorkspace_sptr &lhsWS, const API::MatrixWorkspace_sptr &rhsWS)
const std::string category() const override
Algorithm's category for identification overriding a virtual method.
const std::string name() const override
Algorithm's name for identification overriding a virtual method.
int version() const override
Algorithm's version for identification overriding a virtual method.
void init() override
Initialize the properties to default values.
API::MatrixWorkspace_sptr minus(const API::MatrixWorkspace_sptr &lhsWS, const API::MatrixWorkspace_sptr &rhsWS)
API::WorkspaceGroup_sptr calculateCorrection(API::MatrixWorkspace_sptr &inputWksp, double radius, double coeff1, double coeff2, double coeff3, bool doAbs, bool doMS)
This class is intended to fulfill the design specified in <https://github.com/mantidproject/documents...
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
Definition: Workspace_fwd.h:20
std::shared_ptr< MatrixWorkspace > MatrixWorkspace_sptr
shared pointer to the matrix workspace base class
std::shared_ptr< EventWorkspace > EventWorkspace_sptr
shared pointer to the EventWorkspace class
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