Mantid
Loading...
Searching...
No Matches
SaveSESANS.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 +
8
11#include "MantidAPI/Sample.h"
15
16#include <algorithm>
17#include <cmath>
18#include <fstream>
19#include <iomanip>
20#include <map>
21#include <set>
22#include <string>
23
24namespace Mantid::DataHandling {
25
26// Register the algorithm with the AlgorithmFactory
27DECLARE_ALGORITHM(SaveSESANS)
28
29
30const std::string SaveSESANS::name() const { return "SaveSESANS"; }
31
33const std::string SaveSESANS::summary() const { return "Save a file using the SESANS format"; }
34
36int SaveSESANS::version() const { return 1; }
37
39const std::string SaveSESANS::category() const { return "DataHandling\\Text"; }
40
41std::map<std::string, std::string> SaveSESANS::validateInputs() {
42 std::map<std::string, std::string> invalidInputs;
43
44 for (const auto &propertyName : mandatoryDoubleProperties) {
45 double value = getProperty(propertyName);
46 if (value == EMPTY_DBL()) {
47 invalidInputs[propertyName] = propertyName + " must be set";
48 }
49 }
50
51 if (getPropertyValue("Sample").empty()) {
52 invalidInputs["Sample"] = "Sample must be set";
53 }
54
55 m_sampleThickness = static_cast<double>(getProperty("OverrideSampleThickness"));
57 invalidInputs["OverrideSampleThickness"] = "OverrideSampleThickness value must be greater than 0";
58 }
59
60 return invalidInputs;
61}
62
67 auto validOrientation = std::make_shared<Kernel::StringListValidator>(std::set<std::string>{"X", "Y", "Z"});
68
69 declareProperty(std::make_unique<API::WorkspaceProperty<>>("InputWorkspace", "", Kernel::Direction::Input),
70 "The name of the workspace to save");
71 declareProperty(std::make_unique<API::FileProperty>("Filename", "", API::FileProperty::Save, fileExtensions),
72 "The name to use when saving the file");
73
74 declareProperty("ThetaZMax", EMPTY_DBL(), "The angular acceptance in the encoding direction",
76 declareProperty("ThetaZMaxUnit", "radians", Kernel::Direction::Input);
77 declareProperty("ThetaYMax", EMPTY_DBL(), "The angular acceptance in the non-encoding direction",
79 declareProperty("ThetaYMaxUnit", "radians", Kernel::Direction::Input);
80 declareProperty("EchoConstant", EMPTY_DBL(), "The spin echo length, in nanometers, probed by a 1A neutron",
82
83 declareProperty<std::string>("Sample", "", "Sample name", Kernel::Direction::Input);
84
85 declareProperty<std::string>("Orientation", "Z", validOrientation, "Orientation of the instrument");
86
87 declareProperty("OverrideSampleThickness", EMPTY_DBL(),
88 "The sample thickness in mm. If set, this value will be used instead of the thickness from the input "
89 "workspace sample",
91}
92
97 API::MatrixWorkspace_const_sptr ws = getProperty("InputWorkspace");
98
99 // Check workspace has only one spectrum
100 if (ws->getNumberHistograms() != 1) {
101 g_log.error("This algorithm expects a workspace with exactly 1 spectrum");
102 throw std::runtime_error("SaveSESANS passed workspace with incorrect "
103 "number of spectra, expected 1");
104 }
105
106 if (m_sampleThickness == EMPTY_DBL()) {
107 m_sampleThickness = ws->sample().getThickness();
109 g_log.error("The workspace passed in does not provide the sample thickness. Please use the "
110 "OverrideSampleThickness property to "
111 "provide this value instead");
112 throw std::runtime_error("SaveSESANS passed workspace with sample thickness less than "
113 "or equal to 0 and no value provided for OverrideSampleThickness property");
114 }
115 }
116
117 auto filename = getPropertyValue("Filename");
118 std::ofstream outfile(filename, std::ofstream::trunc);
119 if (outfile.fail()) {
120 const std::string error = strerror(errno);
121 g_log.error("Failed to open file. Error was: " + error);
122 throw std::runtime_error("Could not open file at the following path: " + filename);
123 }
124
125 writeHeaders(outfile, ws);
126 outfile << "\n"
127 << "BEGIN_DATA"
128 << "\n";
129
130 const auto &wavelength = ws->points(0);
131 const auto &yValues = ws->y(0);
132 const auto &eValues = ws->e(0);
133
134 const auto spinEchoLength = calculateSpinEchoLength(wavelength);
135 const auto depolarisation = calculateDepolarisation(yValues, wavelength);
136 const auto error = calculateError(eValues, yValues, wavelength);
137
138 outfile << "SpinEchoLength Depolarisation Depolarisation_error Wavelength\n";
139
140 for (size_t i = 0; i < spinEchoLength.size(); ++i) {
141 outfile << spinEchoLength[i] << " ";
142 outfile << depolarisation[i] << " ";
143 outfile << error[i] << " ";
144 outfile << wavelength[i] << "\n";
145 }
146
147 outfile.close();
148}
149
156 writeHeader(outfile, "FileFormatVersion", "1.0");
157 writeHeader(outfile, "DataFileTitle", ws->getTitle());
158 writeHeader(outfile, "Sample", getPropertyValue("Sample"));
159 writeHeader(outfile, "Thickness", std::to_string(m_sampleThickness));
160 writeHeader(outfile, "Thickness_unit", "mm");
161 writeHeader(outfile, "Theta_zmax", getPropertyValue("ThetaZMax"));
162 writeHeader(outfile, "Theta_zmax_unit", getPropertyValue("ThetaZMaxUnit"));
163 writeHeader(outfile, "Theta_ymax", getPropertyValue("ThetaYMax"));
164 writeHeader(outfile, "Theta_ymax_unit", getPropertyValue("ThetaYMaxUnit"));
165 writeHeader(outfile, "Orientation", "Z");
166 writeHeader(outfile, "SpinEchoLength_unit", "A");
167 writeHeader(outfile, "Depolarisation_unit", "A-2 cm-1");
168 writeHeader(outfile, "Wavelength_unit", "A");
169 writeHeader(outfile, "Echo_constant", getPropertyValue("EchoConstant"));
170}
171
178void SaveSESANS::writeHeader(std::ofstream &outfile, const std::string &name, const std::string &value) {
179 outfile << std::setfill(' ') << std::setw(MAX_HDR_LENGTH) << std::left << name << value << "\n";
180}
181
188std::vector<double> SaveSESANS::calculateSpinEchoLength(const HistogramData::Points &wavelength) const {
189 std::vector<double> spinEchoLength;
190 const double echoConstant = getProperty("EchoConstant");
191
192 // SEL is calculated as wavelength^2 * echoConstant
193 transform(wavelength.begin(), wavelength.end(), back_inserter(spinEchoLength),
194 [&](double w) { return w * w * echoConstant; });
195 return spinEchoLength;
196}
197
206std::vector<double> SaveSESANS::calculateDepolarisation(const HistogramData::HistogramY &yValues,
207 const HistogramData::Points &wavelength) const {
208 Mantid::MantidVec depolarisation;
209 // The sample thickness is provided in mm but we need to use cm here
210 double thickness = m_sampleThickness * 0.1;
211
212 // Depol is calculated as ln(y) / wavelength^2 / sample thickness in cm
213 transform(yValues.begin(), yValues.end(), wavelength.begin(), back_inserter(depolarisation),
214 [&](double y, double w) { return log(y) / (w * w) / thickness; });
215 return depolarisation;
216}
217
226Mantid::MantidVec SaveSESANS::calculateError(const HistogramData::HistogramE &eValues,
227 const HistogramData::HistogramY &yValues,
228 const HistogramData::Points &wavelength) const {
230
231 // The sample thickness is provided in mm but we need to use cm here
232 double thickness = m_sampleThickness * 0.1;
233
234 // Error is calculated as e / (y * wavelength^2) / sample thickness in cm
235 for (size_t i = 0; i < eValues.size(); i++) {
236 error.emplace_back(eValues[i] / (yValues[i] * wavelength[i] * wavelength[i]) / thickness);
237 }
238 return error;
239}
240
241} // namespace Mantid::DataHandling
#define DECLARE_ALGORITHM(classname)
Definition: Algorithm.h:576
double value
The value of the point.
Definition: FitMW.cpp:51
double error
Definition: IndexPeaks.cpp:133
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
Kernel::Logger & g_log
Definition: Algorithm.h:451
@ Save
to specify a file to write to, the file may or may not exist
Definition: FileProperty.h:49
A property class for workspaces.
SaveSESANS : Save a workspace in the SESANS file format.
Definition: SaveSESANS.h:36
const std::string category() const override
Get the algorithm's category.
Definition: SaveSESANS.cpp:39
std::vector< double > calculateSpinEchoLength(const HistogramData::Points &wavelength) const
Calculate SpinEchoLength column from wavelength ( SEL = wavelength * wavelength * echoConstant)
Definition: SaveSESANS.cpp:188
Mantid::MantidVec calculateError(const HistogramData::HistogramE &eValues, const HistogramData::HistogramY &yValues, const HistogramData::Points &wavelength) const
Calculate the error column from the workspace values (error = e / (y * wavelength^2)) / sample thickn...
Definition: SaveSESANS.cpp:226
std::vector< double > calculateDepolarisation(const HistogramData::HistogramY &yValues, const HistogramData::Points &wavelength) const
Calculate the depolarisation column from wavelength and Y values in the workspace (depol = ln(y) / wa...
Definition: SaveSESANS.cpp:206
int version() const override
Get the version number of the algorithm.
Definition: SaveSESANS.cpp:36
void init() override
Initialise the algorithm.
Definition: SaveSESANS.cpp:66
const std::vector< std::string > fileExtensions
Definition: SaveSESANS.h:51
std::map< std::string, std::string > validateInputs() override
Perform validation of ALL the input properties of the algorithm.
Definition: SaveSESANS.cpp:41
const std::string name() const override
Get the algorithm's name.
Definition: SaveSESANS.cpp:30
void writeHeader(std::ofstream &outfile, const std::string &name, const std::string &value)
Write a single header to the output file.
Definition: SaveSESANS.cpp:178
void writeHeaders(std::ofstream &outfile, API::MatrixWorkspace_const_sptr &ws)
Write header values to the output file.
Definition: SaveSESANS.cpp:155
void exec() override
Execute the algorithm.
Definition: SaveSESANS.cpp:96
const std::vector< std::string > mandatoryDoubleProperties
Definition: SaveSESANS.h:52
const std::string summary() const override
Get a summary of the algorithm.
Definition: SaveSESANS.cpp:33
void error(const std::string &msg)
Logs at error level.
Definition: Logger.cpp:77
std::shared_ptr< const MatrixWorkspace > MatrixWorkspace_const_sptr
shared pointer to the matrix workspace base class (const version)
std::vector< double > MantidVec
typedef for the data storage used in Mantid matrix workspaces
Definition: cow_ptr.h:172
constexpr double EMPTY_DBL() noexcept
Returns what we consider an "empty" double within a property.
Definition: EmptyValues.h:43
STL namespace.
std::string to_string(const wide_integer< Bits, Signed > &n)
@ Input
An input workspace.
Definition: Property.h:53