Mantid
Loading...
Searching...
No Matches
EQSANSDarkCurrentSubtraction.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"
17#include "Poco/Path.h"
18#include "Poco/String.h"
19
21
22// Register the algorithm into the AlgorithmFactory
23DECLARE_ALGORITHM(EQSANSDarkCurrentSubtraction)
24
25using namespace Kernel;
26using namespace API;
27using namespace Geometry;
28using namespace DataObjects;
29
31 auto wsValidator = std::make_shared<WorkspaceUnitValidator>("Wavelength");
32 declareProperty(std::make_unique<WorkspaceProperty<>>("InputWorkspace", "", Direction::Input, wsValidator));
33
34 declareProperty(std::make_unique<API::FileProperty>("Filename", "", API::FileProperty::Load, "_event.nxs"),
35 "The name of the input event Nexus file to load as dark current.");
36
37 declareProperty(std::make_unique<WorkspaceProperty<>>("OutputWorkspace", "", Direction::Output));
38 declareProperty("PersistentCorrection", true,
39 "If true, the algorithm will be persistent and re-used when "
40 "other data sets are processed");
41 declareProperty("ReductionProperties", "__sans_reduction_properties", Direction::Input);
42 declareProperty(std::make_unique<WorkspaceProperty<MatrixWorkspace>>("OutputDarkCurrentWorkspace", "",
44 declareProperty("OutputMessage", "", Direction::Output);
45}
46
48 std::string output_message;
49 // Reduction property manager
50 const std::string reductionManagerName = getProperty("ReductionProperties");
51 std::shared_ptr<PropertyManager> reductionManager;
52 if (PropertyManagerDataService::Instance().doesExist(reductionManagerName)) {
53 reductionManager = PropertyManagerDataService::Instance().retrieve(reductionManagerName);
54 } else {
55 reductionManager = std::make_shared<PropertyManager>();
56 PropertyManagerDataService::Instance().addOrReplace(reductionManagerName, reductionManager);
57 }
58
59 // If the load algorithm isn't in the reduction properties, add it
60 const bool persistent = getProperty("PersistentCorrection");
61 if (!reductionManager->existsProperty("DarkCurrentAlgorithm") && persistent) {
62 auto algProp = std::make_unique<AlgorithmProperty>("DarkCurrentAlgorithm");
63 algProp->setValue(toString());
64 reductionManager->declareProperty(std::move(algProp));
65 }
66
67 Progress progress(this, 0.0, 1.0, 10);
68
69 MatrixWorkspace_sptr inputWS = getProperty("InputWorkspace");
70
71 const std::string fileName = getPropertyValue("Filename");
73
74 progress.report("Subtracting dark current");
75
76 // Look for an entry for the dark current in the reduction table
77 Poco::Path path(fileName);
78 const std::string entryName = "DarkCurrent" + path.getBaseName();
79 std::string darkWSName = "__dark_current_" + path.getBaseName();
80
81 if (reductionManager->existsProperty(entryName)) {
82 darkWS = reductionManager->getProperty(entryName);
83 darkWSName = reductionManager->getPropertyValue(entryName);
84 output_message += darkWSName + '\n';
85 } else {
86 // Load the dark current if we don't have it already
87 IAlgorithm_sptr loadAlg;
88 if (!reductionManager->existsProperty("LoadAlgorithm")) {
89 loadAlg = createChildAlgorithm("EQSANSLoad", 0.1, 0.3);
90 loadAlg->setProperty("Filename", fileName);
91 if (loadAlg->existsProperty("LoadMonitors"))
92 loadAlg->setProperty("LoadMonitors", false);
93 loadAlg->executeAsChildAlg();
94 darkWS = loadAlg->getProperty("OutputWorkspace");
95 } else {
96 // Get load algorithm as a string so that we can create a completely
97 // new proxy and ensure that we don't overwrite existing properties
98 IAlgorithm_sptr loadAlg0 = reductionManager->getProperty("LoadAlgorithm");
99 const std::string loadString = loadAlg0->toString();
100 loadAlg = Algorithm::fromString(loadString);
101 loadAlg->setChild(true);
102 loadAlg->setProperty("Filename", fileName);
103 if (loadAlg->existsProperty("LoadMonitors"))
104 loadAlg->setProperty("LoadMonitors", false);
105 loadAlg->setPropertyValue("OutputWorkspace", darkWSName);
106 loadAlg->execute();
107 darkWS = loadAlg->getProperty("OutputWorkspace");
108 }
109
110 output_message += "\n Loaded " + fileName + "\n";
111 if (loadAlg->existsProperty("OutputMessage")) {
112 std::string msg = loadAlg->getPropertyValue("OutputMessage");
113 output_message += " |" + Poco::replace(msg, "\n", "\n |") + "\n";
114 }
115
116 std::string darkWSOutputName = getPropertyValue("OutputDarkCurrentWorkspace");
117 if (!darkWSOutputName.empty())
118 setProperty("OutputDarkCurrentWorkspace", darkWS);
119 AnalysisDataService::Instance().addOrReplace(darkWSName, darkWS);
120 reductionManager->declareProperty(std::make_unique<WorkspaceProperty<>>(entryName, "", Direction::Output));
121 reductionManager->setPropertyValue(entryName, darkWSName);
122 reductionManager->setProperty(entryName, darkWS);
123 }
124 progress.report(3, "Loaded dark current");
125
126 // Normalize the dark current and data to counting time
127 double scaling_factor = 1.0;
128 if (inputWS->run().hasProperty("proton_charge")) {
129 auto dp = inputWS->run().getTimeSeriesProperty<double>("proton_charge");
130 double duration = dp->getStatistics().duration;
131
132 dp = darkWS->run().getTimeSeriesProperty<double>("proton_charge");
133 double dark_duration = dp->getStatistics().duration;
134 scaling_factor = duration / dark_duration;
135 } else if (inputWS->run().hasProperty("timer")) {
136 auto duration = inputWS->run().getPropertyValueAsType<double>("timer");
137 auto dark_duration = darkWS->run().getPropertyValueAsType<double>("timer");
138 ;
139 scaling_factor = duration / dark_duration;
140 } else {
141 output_message += "\n Could not find proton charge or duration in sample logs";
142 g_log.error() << "ERROR: Could not find proton charge or duration in sample logs\n";
143 };
144
145 progress.report("Scaling dark current");
146
147 // Scale the stored dark current by the counting time
148 auto rebinAlg = createChildAlgorithm("RebinToWorkspace", 0.4, 0.5);
149 rebinAlg->setProperty("WorkspaceToRebin", darkWS);
150 rebinAlg->setProperty("WorkspaceToMatch", inputWS);
151 rebinAlg->setProperty("OutputWorkspace", darkWS);
152 rebinAlg->executeAsChildAlg();
153 MatrixWorkspace_sptr scaledDarkWS = rebinAlg->getProperty("OutputWorkspace");
154
155 // Perform subtraction
156 auto scaleAlg = createChildAlgorithm("Scale", 0.5, 0.6);
157 scaleAlg->setProperty("InputWorkspace", scaledDarkWS);
158 scaleAlg->setProperty("Factor", scaling_factor);
159 scaleAlg->setProperty("OutputWorkspace", scaledDarkWS);
160 scaleAlg->setProperty("Operation", "Multiply");
161 scaleAlg->executeAsChildAlg();
162 scaledDarkWS = rebinAlg->getProperty("OutputWorkspace");
163
164 auto minusAlg = createChildAlgorithm("Minus", 0.6, 0.7);
165 minusAlg->setProperty("LHSWorkspace", inputWS);
166 minusAlg->setProperty("RHSWorkspace", scaledDarkWS);
167 const std::string outputWSname = getPropertyValue("OutputWorkspace");
168 minusAlg->setPropertyValue("OutputWorkspace", outputWSname);
169 minusAlg->executeAsChildAlg();
170 MatrixWorkspace_sptr outputWS = minusAlg->getProperty("OutputWorkspace");
171
172 setProperty("OutputWorkspace", outputWS);
173 setProperty("OutputMessage", "Dark current subtracted: " + output_message);
174
175 progress.report("Subtracted dark current");
176}
177
178} // 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
Kernel::Logger & g_log
Definition: Algorithm.h:451
void progress(double p, const std::string &msg="", double estimatedTime=0.0, int progressPrecision=0)
Sends ProgressNotification.
Definition: Algorithm.cpp:231
@ 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.
void error(const std::string &msg)
Logs at error level.
Definition: Logger.cpp:77
static T & Instance()
Return a reference to the Singleton instance, creating it if it does not already exist Creation is do...
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
@ Input
An input workspace.
Definition: Property.h:53
@ Output
An output workspace.
Definition: Property.h:54