Mantid
Loading...
Searching...
No Matches
HFIRDarkCurrentSubtraction.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 +
11#include "MantidAPI/Run.h"
15#include "Poco/Path.h"
16#include "Poco/String.h"
17
19
20// Register the algorithm into the AlgorithmFactory
21DECLARE_ALGORITHM(HFIRDarkCurrentSubtraction)
22
23using namespace Kernel;
24using namespace API;
25using namespace Geometry;
26
28 auto wsValidator = std::make_shared<WorkspaceUnitValidator>("Wavelength");
29 declareProperty(std::make_unique<WorkspaceProperty<>>("InputWorkspace", "", Direction::Input, wsValidator));
30
31 declareProperty(std::make_unique<API::FileProperty>("Filename", "", API::FileProperty::Load, ".xml"),
32 "The name of the input event Nexus file to load as dark current.");
33
34 declareProperty(std::make_unique<WorkspaceProperty<>>("OutputWorkspace", "", Direction::Output));
35 declareProperty("PersistentCorrection", true,
36 "If true, the algorithm will be persistent and re-used when "
37 "other data sets are processed");
38 declareProperty("ReductionProperties", "__sans_reduction_properties", Direction::Input);
39 declareProperty(std::make_unique<WorkspaceProperty<MatrixWorkspace>>("OutputDarkCurrentWorkspace", "",
41 declareProperty("OutputMessage", "", Direction::Output);
42}
43
45 std::string output_message;
46 // Reduction property manager
47 const std::string reductionManagerName = getProperty("ReductionProperties");
48 std::shared_ptr<PropertyManager> reductionManager;
49 if (PropertyManagerDataService::Instance().doesExist(reductionManagerName)) {
50 reductionManager = PropertyManagerDataService::Instance().retrieve(reductionManagerName);
51 } else {
52 reductionManager = std::make_shared<PropertyManager>();
53 PropertyManagerDataService::Instance().addOrReplace(reductionManagerName, reductionManager);
54 }
55
56 // If the load algorithm isn't in the reduction properties, add it
57 const bool persistent = getProperty("PersistentCorrection");
58 if (!reductionManager->existsProperty("DarkCurrentAlgorithm") && persistent) {
59 auto algProp = std::make_unique<AlgorithmProperty>("DarkCurrentAlgorithm");
60 algProp->setValue(toString());
61 reductionManager->declareProperty(std::move(algProp));
62 }
63
64 Progress progress(this, 0.0, 1.0, 10);
65
66 MatrixWorkspace_sptr inputWS = getProperty("InputWorkspace");
67 const std::string fileName = getPropertyValue("Filename");
69 std::string darkWSName = getPropertyValue("OutputDarkCurrentWorkspace");
70
71 progress.report("Subtracting dark current");
72
73 // Look for an entry for the dark current in the reduction table
74 Poco::Path path(fileName);
75 const std::string entryName = "DarkCurrent" + path.getBaseName();
76
77 if (reductionManager->existsProperty(entryName)) {
78 darkWS = reductionManager->getProperty(entryName);
79 darkWSName = reductionManager->getPropertyValue(entryName);
80 output_message += darkWSName + '\n';
81 } else {
82 // Load the dark current if we don't have it already
83 if (darkWSName.empty()) {
84 darkWSName = "__dark_current_" + path.getBaseName();
85 setPropertyValue("OutputDarkCurrentWorkspace", darkWSName);
86 }
87
88 IAlgorithm_sptr loadAlg;
89 if (!reductionManager->existsProperty("LoadAlgorithm")) {
90 loadAlg = createChildAlgorithm("HFIRLoad", 0.1, 0.3);
91 loadAlg->setProperty("Filename", fileName);
92 loadAlg->setProperty("ReductionProperties", reductionManagerName);
93 loadAlg->executeAsChildAlg();
94 } else {
95 IAlgorithm_sptr loadAlg0 = reductionManager->getProperty("LoadAlgorithm");
96 const std::string loadString = loadAlg0->toString();
97 loadAlg = Algorithm::fromString(loadString);
98 loadAlg->setChild(true);
99 loadAlg->setProperty("Filename", fileName);
100 loadAlg->setProperty("ReductionProperties", reductionManagerName);
101 loadAlg->setPropertyValue("OutputWorkspace", darkWSName);
102 loadAlg->execute();
103 }
104 darkWS = loadAlg->getProperty("OutputWorkspace");
105 output_message += "\n Loaded " + fileName + "\n";
106 if (loadAlg->existsProperty("OutputMessage")) {
107 std::string msg = loadAlg->getPropertyValue("OutputMessage");
108 output_message += " |" + Poco::replace(msg, "\n", "\n |") + "\n";
109 }
110
111 setProperty("OutputDarkCurrentWorkspace", darkWS);
112 reductionManager->declareProperty(std::make_unique<WorkspaceProperty<>>(entryName, "", Direction::Output));
113 reductionManager->setPropertyValue(entryName, darkWSName);
114 reductionManager->setProperty(entryName, darkWS);
115 }
116 progress.report(3, "Loaded dark current");
117
118 // Perform subtraction
119 double darkTimer = getCountingTime(darkWS);
120 double dataTimer = getCountingTime(inputWS);
121 auto scaleAlg = createChildAlgorithm("Scale", 0.3, 0.5);
122 scaleAlg->setProperty("InputWorkspace", darkWS);
123 scaleAlg->setProperty("Factor", dataTimer / darkTimer);
124 scaleAlg->setProperty("Operation", "Multiply");
125 scaleAlg->executeAsChildAlg();
126 MatrixWorkspace_sptr scaledDarkWS = scaleAlg->getProperty("OutputWorkspace");
127
128 // Zero out timer and monitor so that we don't subtract them out
129 for (size_t i = 0; i < scaledDarkWS->dataY(0).size(); i++) {
130 scaledDarkWS->dataY(DEFAULT_TIMER_ID)[i] = 0.0;
131 scaledDarkWS->dataE(DEFAULT_TIMER_ID)[i] = 0.0;
132 scaledDarkWS->dataY(DEFAULT_MONITOR_ID)[i] = 0.0;
133 scaledDarkWS->dataE(DEFAULT_MONITOR_ID)[i] = 0.0;
134 }
135 auto minusAlg = createChildAlgorithm("Minus", 0.5, 0.7);
136 minusAlg->setProperty("LHSWorkspace", inputWS);
137 minusAlg->setProperty("RHSWorkspace", scaledDarkWS);
138 MatrixWorkspace_sptr outputWS = getProperty("OutputWorkspace");
139 minusAlg->setProperty("OutputWorkspace", outputWS);
140 minusAlg->executeAsChildAlg();
141 MatrixWorkspace_sptr correctedWS = minusAlg->getProperty("OutputWorkspace");
142 setProperty("OutputWorkspace", correctedWS);
143 setProperty("OutputMessage", "Dark current subtracted: " + output_message);
144
145 progress.report("Subtracted dark current");
146}
147
151 // First, look whether we have the information in the log
152 if (inputWS->run().hasProperty("timer")) {
153 return inputWS->run().getPropertyValueAsType<double>("timer");
154 } else {
155 // If we don't have the information in the log, use the default timer
156 // spectrum
157 MantidVec &timer = inputWS->dataY(DEFAULT_TIMER_ID);
158 return timer[0];
159 }
160}
161
162} // namespace Mantid::WorkflowAlgorithms
#define DECLARE_ALGORITHM(classname)
Definition: Algorithm.h:576
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
std::string getPropertyValue(const std::string &name) const override
Get the value of a property as a string.
Definition: Algorithm.cpp:2026
TypedValue getProperty(const std::string &name) const override
Get the value of a property.
Definition: Algorithm.cpp:2076
std::string toString() const override
Serialize an object to a string.
Definition: Algorithm.cpp:905
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.
Definition: Algorithm.cpp:842
static IAlgorithm_sptr fromString(const std::string &input)
De-serialize an object from a string.
Definition: Algorithm.cpp:967
void progress(double p, const std::string &msg="", double estimatedTime=0.0, int progressPrecision=0)
Sends ProgressNotification.
Definition: Algorithm.cpp:231
void setPropertyValue(const std::string &name, const std::string &value) override
Set the value of a property by string N.B.
Definition: Algorithm.cpp:1975
@ Load
allowed here which will be passed to the algorithm
Definition: FileProperty.h:52
Helper class for reporting progress from algorithms.
Definition: Progress.h:25
A property class for workspaces.
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...
double getCountingTime(const API::MatrixWorkspace_sptr &inputWS)
Get the counting time from a workspace.
std::shared_ptr< IAlgorithm > IAlgorithm_sptr
shared pointer to Mantid::API::IAlgorithm
std::shared_ptr< MatrixWorkspace > MatrixWorkspace_sptr
shared pointer to the matrix workspace base class
std::vector< double > MantidVec
typedef for the data storage used in Mantid matrix workspaces
Definition: cow_ptr.h:172
@ Input
An input workspace.
Definition: Property.h:53
@ Output
An output workspace.
Definition: Property.h:54