Mantid
Loading...
Searching...
No Matches
ConvertToConstantL2.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 +
9#include "MantidAPI/Run.h"
19#include "MantidTypes/SpectrumDefinition.h"
20
21#include <boost/format.hpp>
22#include <cmath>
23
24namespace Mantid::Algorithms {
25
26using namespace Kernel;
27using namespace API;
28using namespace Geometry;
29using namespace DataObjects;
30
31// Register the class into the algorithm factory
32DECLARE_ALGORITHM(ConvertToConstantL2)
33
34// Constructor
36 : API::Algorithm(), m_inputWS(), m_outputWS(), m_instrument(), m_l2(0.), m_wavelength(0.) {}
37
42 auto wsValidator = std::make_shared<CompositeValidator>();
43 wsValidator->add<WorkspaceUnitValidator>("TOF");
44 wsValidator->add<HistogramValidator>();
46 std::make_unique<WorkspaceProperty<API::MatrixWorkspace>>("InputWorkspace", "", Direction::Input, wsValidator),
47 "Name of the input workspace");
48 declareProperty(std::make_unique<WorkspaceProperty<API::MatrixWorkspace>>("OutputWorkspace", "", Direction::Output),
49 "Name of the output workspace, can be the same as the input");
50}
51
57 // Get the workspaces
58 m_inputWS = this->getProperty("InputWorkspace");
59 m_outputWS = this->getProperty("OutputWorkspace");
60 m_instrument = m_inputWS->getInstrument();
61 // If input and output workspaces are not the same, create a new workspace for
62 // the output
63 if (m_outputWS != this->m_inputWS) {
64 m_outputWS = create<MatrixWorkspace>(*m_inputWS);
65 }
66
67 m_wavelength = getRunProperty("wavelength");
68 g_log.debug() << "Wavelength = " << m_wavelength;
70 g_log.debug() << " L2 = " << m_l2 << '\n';
71}
72
78
80
81 // Calculate the number of spectra in this workspace
82 const size_t numberOfSpectra = m_inputWS->getNumberHistograms();
83 API::Progress prog(this, 0.0, 1.0, numberOfSpectra);
84
85 auto numberOfSpectra_i = static_cast<int64_t>(numberOfSpectra); // cast to make openmp happy
86
87 const auto &inputSpecInfo = m_inputWS->spectrumInfo();
88 auto &outputDetInfo = m_outputWS->mutableDetectorInfo();
89
90 // Loop over the histograms (detector spectra)
92 for (int64_t i = 0; i < numberOfSpectra_i; ++i) {
94 m_outputWS->setHistogram(i, m_inputWS->histogram(i));
95
96 // Should not move the monitors
97 if (inputSpecInfo.isMonitor(i))
98 continue;
99
100 // Throw if detector doesn't exist or is a group
101 if (!inputSpecInfo.hasUniqueDetector(i)) {
102 const auto errorMsg = boost::format("The detector for spectrum number %d was either not "
103 "found, or is a group.") %
104 i;
105 throw std::runtime_error(errorMsg.str());
106 }
107
108 // subract the diference in l2
109 double thisDetL2 = inputSpecInfo.l2(i);
110 double deltaL2 = std::abs(thisDetL2 - m_l2);
111 double deltaTOF = calculateTOF(deltaL2);
112 deltaTOF *= 1e6; // micro sec
113
114 // position - set all detector distance to constant l2
115 double r, theta, phi;
116 V3D oldPos = inputSpecInfo.position(i);
117 oldPos.getSpherical(r, theta, phi);
118 V3D newPos;
119 newPos.spherical(m_l2, theta, phi);
120
121 const auto detIndex = inputSpecInfo.spectrumDefinition(i)[0];
122 outputDetInfo.setPosition(detIndex, newPos);
123
124 m_outputWS->mutableX(i) -= deltaTOF;
125
126 prog.report("Aligning elastic line...");
128 } // end for i
130
131 this->setProperty("OutputWorkspace", this->m_outputWS);
132}
133
134/*
135 * Get run property as double
136 * @s - input property name
137 *
138 */
139double ConvertToConstantL2::getRunProperty(const std::string &s) {
140 const auto &run = m_inputWS->run();
141 if (!run.hasProperty(s)) {
142 throw Exception::NotFoundError("Sample log property not found", s);
143 }
144 Mantid::Kernel::Property *prop = run.getProperty(s);
145 double val;
146 if (!Strings::convert(prop->value(), val)) {
147 const std::string mesg = "Cannot convert sample log '" + s + "' to a number.";
148 throw std::runtime_error(mesg);
149 }
150 return val;
151}
152/*
153 * Get instrument property as double
154 * @s - input property name
155 *
156 */
157double ConvertToConstantL2::getInstrumentProperty(const std::string &s) {
158 std::vector<std::string> prop = m_instrument->getStringParameter(s);
159 if (prop.empty()) {
160 const std::string mesg = "Property <" + s + "> doesn't exist!";
161 throw std::runtime_error(mesg);
162 }
163 g_log.debug() << "prop[0] = " << prop[0] << '\n';
164 return boost::lexical_cast<double>(prop[0]);
165}
166
167/*
168 * Returns the neutron TOF
169 * @distance - Distance in meters
170 */
171double ConvertToConstantL2::calculateTOF(double distance) {
172 double velocity = PhysicalConstants::h / (PhysicalConstants::NeutronMass * m_wavelength * 1e-10); // m/s
173
174 return distance / velocity;
175}
176
177} // 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.
Base class from which all concrete algorithm classes should be derived.
Definition: Algorithm.h:85
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
Kernel::Logger & g_log
Definition: Algorithm.h:451
A validator which checks that a workspace contains histogram data (the default) or point data as requ...
Helper class for reporting progress from algorithms.
Definition: Progress.h:25
A property class for workspaces.
A validator which checks that the unit of the workspace referred to by a WorkspaceProperty is the exp...
Convert workspace to have a constant L2.
API::MatrixWorkspace_sptr m_outputWS
The output workspace, maybe the same as the input one.
API::MatrixWorkspace_const_sptr m_inputWS
The user selected (input) workspace.
void exec() override
Executes the algorithm.
Geometry::Instrument_const_sptr m_instrument
double getInstrumentProperty(const std::string &)
void init() override
Initialisation method.
void initWorkspaces()
Initialises input and output workspaces.
Exception for when an item is not found in a collection.
Definition: Exception.h:145
IPropertyManager * setProperty(const std::string &name, const T &value)
Templated method to set the value of a PropertyWithValue.
void debug(const std::string &msg)
Logs at debug level.
Definition: Logger.cpp:114
void report()
Increments the loop counter by 1, then sends the progress notification on behalf of its algorithm.
Definition: ProgressBase.h:51
Base class for properties.
Definition: Property.h:94
virtual std::string value() const =0
Returns the value of the property as a string.
Class for 3D vectors.
Definition: V3D.h:34
void spherical(const double R, const double theta, const double phi) noexcept
Sets the vector position based on spherical coordinates.
Definition: V3D.cpp:57
void getSpherical(double &R, double &theta, double &phi) const noexcept
Return the vector's position in spherical coordinates.
Definition: V3D.cpp:117
int convert(const std::string &A, T &out)
Convert a string into a number.
Definition: Strings.cpp:665
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
static constexpr double NeutronMass
Mass of the neutron in kg.
static constexpr double h
Planck constant in J*s.
@ Input
An input workspace.
Definition: Property.h:53
@ Output
An output workspace.
Definition: Property.h:54