20#include <boost/regex.hpp>
22#include "Poco/DirectoryIterator.h"
23#include "Poco/NumberFormatter.h"
24#include "Poco/NumberParser.h"
25#include "Poco/String.h"
34using namespace Kernel;
36using namespace Geometry;
37using namespace DataObjects;
41 "The name of the input event Nexus file to load");
43 auto wsValidator = std::make_shared<WorkspaceUnitValidator>(
"TOF");
46 "Input event workspace. Assumed to be unmodified events "
47 "straight from LoadEventNexus");
50 "Then name of the output EventWorkspace");
51 declareProperty(
"NoBeamCenter",
false,
"If true, the detector will not be moved according to the beam center");
53 "If true, the beam center defined in "
54 "the configuration file will be "
57 "Beam position in X pixel "
58 "coordinates (used only if "
59 "UseConfigBeam is false)");
61 "Beam position in Y pixel "
62 "coordinates (used only if "
63 "UseConfigBeam is false)");
65 "If true, the edges of the TOF distribution will be cut "
66 "according to the configuration file");
68 "TOF value below which events will not be "
69 "loaded into the workspace at load-time");
71 "TOF value above which events will not be "
72 "loaded into the workspace at load-time");
73 declareProperty(
"SkipTOFCorrection",
false,
"IF true, the EQSANS TOF correction will be skipped");
75 "Wavelength steps to be used when "
76 "rebinning the data before performing "
79 "If true, the masking information "
80 "found in the configuration file "
82 declareProperty(
"UseConfig",
true,
"If true, the best configuration file found will be used");
84 "If true, the TOF will be modified for the true flight path "
85 "from the sample to the detector pixel");
86 declareProperty(
"PreserveEvents",
true,
"If true, the output workspace will be an event workspace");
88 "Sample to detector distance to use (overrides meta data), in mm");
90 "Offset to the sample to detector distance (use only when "
91 "using the distance found in the meta data), in mm");
93 "Offset to be applied to the sample position (use only when "
94 "using the detector distance found in the meta data), in mm");
96 "Offset to be applied to the detector position (use only when "
97 "using the distance found in the meta data), in mm");
98 declareProperty(
"LoadMonitors",
true,
"If true, the monitor workspace will be loaded");
100 "Reads the embedded Instrument XML from the NeXus file "
101 "(optional, default True). ");
113 throw std::runtime_error(
"Could not cast (interpret) the property " + pname +
114 " as a floating point numeric value.");
123 std::string sns_folder =
"/SNS/EQSANS/shared/instrument_configuration";
124 if (Poco::File(sns_folder).
exists())
129 int max_run_number = 0;
130 std::string config_file;
131 static boost::regex re1(
"eqsans_configuration\\.([0-9]+)$");
132 boost::smatch matches;
133 for (
const auto &searchPath : searchPaths) {
134 if (Poco::File(searchPath).
exists()) {
135 Poco::DirectoryIterator file_it(searchPath);
136 Poco::DirectoryIterator end;
137 for (; file_it != end; ++file_it) {
138 if (boost::regex_search(file_it.name(), matches, re1)) {
139 std::string s = matches[1];
141 Poco::NumberParser::tryParse(s, run_number);
142 if (run_number > max_run_number && run_number <= run) {
143 max_run_number = run_number;
144 config_file = file_it.path().toString();
158 boost::regex re_key(
"rectangular mask", boost::regex::icase);
159 boost::regex re_key_alt(
"elliptical mask", boost::regex::icase);
160 if (boost::regex_search(line, re_key) || boost::regex_search(line, re_key_alt)) {
161 boost::regex re_sig(
"=[ ]*([0-9]+)[ ]*[ ,][ ]*([0-9]+)[ ]*[ ;,][ "
162 "]*([0-9]+)[ ]*[ ,][ ]*([0-9]+)");
163 boost::smatch posVec;
164 if (boost::regex_search(line, posVec, re_sig)) {
165 if (posVec.size() == 5) {
166 for (
int i = 0; i < 4; i++) {
167 std::string num_str = posVec[i + 1];
179 boost::regex re_key(
"tof edge discard", boost::regex::icase);
180 if (boost::regex_search(line, re_key)) {
181 boost::regex re_sig(
"=[ ]*([0-9]+)[ ]*[ ,][ ]*([0-9]+)");
182 boost::smatch posVec;
183 if (boost::regex_search(line, posVec, re_sig)) {
184 if (posVec.size() == 3) {
185 std::string num_str = posVec[1];
197 boost::regex re_key(
"spectrum center", boost::regex::icase);
198 if (boost::regex_search(line, re_key)) {
199 boost::regex re_sig(
"=[ ]*([0-9]+.[0-9]*)[ ]*[ ,][ ]*([0-9]+.[0-9]+)");
200 boost::smatch posVec;
201 if (boost::regex_search(line, posVec, re_sig)) {
202 if (posVec.size() == 3) {
203 std::string num_str = posVec[1];
204 Poco::NumberParser::tryParseFloat(num_str,
m_center_x);
206 Poco::NumberParser::tryParseFloat(num_str,
m_center_y);
215 boost::regex re_key(
"sample location", boost::regex::icase);
216 if (boost::regex_search(line, re_key)) {
217 boost::regex re_sig(
"=[ ]*([0-9]+)");
218 boost::smatch posVec;
219 if (boost::regex_search(line, posVec, re_sig)) {
220 if (posVec.size() == 2) {
221 std::string num_str = posVec[1];
232 boost::regex re_key(
"wheel", boost::regex::icase);
233 if (boost::regex_search(line, re_key)) {
234 boost::regex re_sig(R
"(([1-8]) wheel[ ]*([1-3])[ \t]*=[ \t]*(\w+))");
235 boost::smatch posVec;
236 if (boost::regex_search(line, posVec, re_sig)) {
237 if (posVec.size() == 4) {
238 std::string num_str = posVec[1];
240 Poco::NumberParser::tryParse(num_str, slit_number);
244 int wheel_number = 0;
245 Poco::NumberParser::tryParse(num_str, wheel_number);
249 boost::regex re_size(
"\\w*?([0-9]+)mm");
251 if (boost::regex_search(num_str, posVec, re_size)) {
252 if (posVec.size() == 2) {
254 Poco::NumberParser::tryParse(num_str, slit_size);
265 if (!
dataWS->run().hasProperty(
"vBeamSlit")) {
271 const std::string slit1Name =
"vBeamSlit";
277 slit1 =
static_cast<int>(dp->getStatistics().mean);
281 throw std::runtime_error(
"Could not cast (interpret) the property " + slit1Name +
282 " as a time series property with "
283 "int or floating point values.");
285 const std::string slit2Name =
"vBeamSlit2";
286 prop =
dataWS->run().getProperty(slit2Name);
291 slit2 =
static_cast<int>(dp->getStatistics().mean);
295 throw std::runtime_error(
"Could not cast (interpret) the property " + slit2Name +
296 " as a time series property with "
297 "int or floating point values.");
299 const std::string slit3Name =
"vBeamSlit3";
300 prop =
dataWS->run().getProperty(slit3Name);
305 slit3 =
static_cast<int>(dp->getStatistics().mean);
309 throw std::runtime_error(
"Could not cast (interpret) the property " + slit3Name +
310 " as a time series property with "
311 "int or floating point values.");
313 if (slit1 < 0 && slit2 < 0 && slit3 < 0) {
321 const double ssd =
fabs(
dataWS->getInstrument()->getSource()->getPos().Z()) * 1000.0;
322 int slits[3] = {slit1, slit2, slit3};
323 for (
int i = 0; i < 3; i++) {
324 int m = slits[i] - 1;
325 if (
m >= 0 &&
m < 6) {
328 if (L1 < 0 ||
x /
y < S1 / L1) {
334 dataWS->mutableRun().addProperty(
"source-aperture-diameter", S1,
"mm",
true);
351 int nx_pixels =
static_cast<int>(
dataWS->getInstrument()->getNumberParameter(
"number-of-x-pixels")[0]);
352 int ny_pixels =
static_cast<int>(
dataWS->getInstrument()->getNumberParameter(
"number-of-y-pixels")[0]);
353 const auto &detInfo =
dataWS->detectorInfo();
354 const V3D pixel_first = detInfo.position(detInfo.indexOf(0));
358 const V3D pixel_last_x = detInfo.position(detInfo.indexOf(detIDx));
359 const V3D pixel_last_y = detInfo.position(detInfo.indexOf(detIDy));
360 double x_offset = (pixel_first.
X() + pixel_last_x.
X()) / 2.0;
361 double y_offset = (pixel_first.
Y() + pixel_last_y.
Y()) / 2.0;
362 double beam_ctr_x = 0.0;
363 double beam_ctr_y = 0.0;
368 mvAlg->setProperty(
"ComponentName",
"detector1");
369 mvAlg->setProperty(
"X", -x_offset - beam_ctr_x);
370 mvAlg->setProperty(
"Y", -y_offset - beam_ctr_y);
371 mvAlg->setProperty(
"RelativePosition",
true);
372 mvAlg->executeAsChildAlg();
373 m_output_message +=
" Beam center offset: " + Poco::NumberFormatter::format(x_offset) +
", " +
374 Poco::NumberFormatter::format(y_offset) +
" m\n";
380 dataWS->mutableRun().addProperty(
"beam_center_x",
m_center_x,
"pixel",
true);
381 dataWS->mutableRun().addProperty(
"beam_center_y",
m_center_y,
"pixel",
true);
383 Poco::NumberFormatter::format(
m_center_y) +
"\n";
394 bool use_config_mask =
getProperty(
"UseConfigMask");
395 bool use_config_cutoff =
getProperty(
"UseConfigTOFCuts");
396 bool use_config_center =
getProperty(
"UseConfigBeam");
398 std::ifstream file(filePath.c_str());
400 g_log.
error() <<
"Unable to open file: " << filePath <<
'\n';
407 while (getline(file, line)) {
409 std::string comment = line.substr(0, 1);
410 if (Poco::icompare(comment,
"#") == 0)
414 if (use_config_cutoff)
416 if (use_config_center)
440 if (fileName.empty() && !inputEventWS) {
441 g_log.
error() <<
"EQSANSLoad input error: Either a valid file path or an "
442 "input workspace must be provided\n";
443 throw std::runtime_error(
"EQSANSLoad input error: Either a valid file path "
444 "or an input workspace must be provided");
445 }
else if (!fileName.empty() && inputEventWS) {
446 g_log.
error() <<
"EQSANSLoad input error: Either a valid file path or an "
447 "input workspace must be provided, but not both\n";
448 throw std::runtime_error(
"EQSANSLoad input error: Either a valid file path "
449 "or an input workspace must be provided, but not "
454 const bool skipTOFCorrection =
getProperty(
"SkipTOFCorrection");
461 const bool noBeamCenter =
getProperty(
"NoBeamCenter");
464 const std::string reductionManagerName =
getProperty(
"ReductionProperties");
465 std::shared_ptr<PropertyManager> reductionManager;
469 reductionManager = std::make_shared<PropertyManager>();
473 if (!reductionManager->existsProperty(
"LoadAlgorithm")) {
474 auto loadProp = std::make_unique<AlgorithmProperty>(
"LoadAlgorithm");
478 reductionManager->declareProperty(std::move(loadProp));
481 if (!reductionManager->existsProperty(
"InstrumentName"))
489 const bool loadMonitors =
getProperty(
"LoadMonitors");
490 const bool loadNexusInstrumentXML =
getProperty(
"LoadNexusInstrumentXML");
492 loadAlg->setProperty(
"LoadMonitors", loadMonitors);
493 loadAlg->setProperty(
"Filename", fileName);
494 loadAlg->setProperty(
"LoadNexusInstrumentXML", loadNexusInstrumentXML);
495 if (skipTOFCorrection) {
502 Workspace_sptr dataWS_asWks = loadAlg->getProperty(
"OutputWorkspace");
503 dataWS = std::dynamic_pointer_cast<MatrixWorkspace>(dataWS_asWks);
507 if (loadMonitors && loadAlg->existsProperty(
"MonitorWorkspace")) {
508 Workspace_sptr monWSOutput = loadAlg->getProperty(
"MonitorWorkspace");
510 if ((monWSOutput) && (!monWS))
514 "data, support for this is not "
515 "implemented in EQSANSLoad yet");
517 "Monitors from the Event NeXus file");
523 if (inputEventWS != outputEventWS) {
525 copyAlg->setProperty(
"InputWorkspace", inputEventWS);
526 copyAlg->executeAsChildAlg();
527 Workspace_sptr dataWS_asWks = copyAlg->getProperty(
"OutputWorkspace");
528 dataWS = std::dynamic_pointer_cast<MatrixWorkspace>(dataWS_asWks);
530 dataWS = std::dynamic_pointer_cast<MatrixWorkspace>(inputEventWS);
538 const double sampleflange_det_dist =
getProperty(
"SampleDetectorDistance");
539 if (!
isEmpty(sampleflange_det_dist)) {
540 sfdd = sampleflange_det_dist;
542 if (!
dataWS->run().hasProperty(
"detectorZ")) {
543 g_log.
error() <<
"Could not determine Z position: the "
544 "SampleDetectorDistance property was not set "
545 "and the run logs do not contain the detectorZ property\n";
546 throw std::invalid_argument(
"Could not determine Z position: stopping execution");
549 const std::string dzName =
"detectorZ";
553 throw std::runtime_error(
"Could not cast (interpret) the property " + dzName +
554 " as a time series property value.");
555 sfdd = dp->getStatistics().mean;
558 const double sampleflange_det_offset =
getProperty(
"DetectorOffset");
559 if (!
isEmpty(sampleflange_det_offset))
560 sfdd += sampleflange_det_offset;
564 const double sample_det_offset =
getProperty(
"SampleDetectorDistanceOffset");
565 if (!
isEmpty(sample_det_offset))
566 sfdd += sample_det_offset;
567 if (!
isEmpty(sample_det_offset) && !
isEmpty(sampleflange_det_offset))
568 g_log.
error() <<
"Both DetectorOffset and SampleDetectorDistanceOffset "
569 "are set. Only one should be used.\n";
574 const double sampleflange_sample_offset =
getProperty(
"SampleOffset");
575 if (!
isEmpty(sampleflange_sample_offset)) {
576 s2d -= sampleflange_sample_offset;
581 mvAlg->setProperty(
"ComponentName",
"sample-position");
582 mvAlg->setProperty(
"Z", sampleflange_sample_offset / 1000.0);
583 mvAlg->setProperty(
"RelativePosition",
false);
584 mvAlg->executeAsChildAlg();
585 g_log.
information() <<
"Moving sample to " << sampleflange_sample_offset / 1000.0 <<
" meters\n";
587 " Sample position: " + Poco::NumberFormatter::format(sampleflange_sample_offset / 1000.0, 3) +
" m\n";
590 dataWS->mutableRun().addProperty(
"sampleflange_detector_distance", sfdd,
"mm",
true);
591 dataWS->mutableRun().addProperty(
"sample_detector_distance", s2d,
"mm",
true);
596 mvAlg->setProperty(
"ComponentName",
"detector1");
597 mvAlg->setProperty(
"Z", sfdd / 1000.0);
598 mvAlg->setProperty(
"RelativePosition",
false);
599 mvAlg->executeAsChildAlg();
600 g_log.
information() <<
"Moving detector to " << sfdd / 1000.0 <<
" meters\n";
601 m_output_message +=
" Detector position: " + Poco::NumberFormatter::format(sfdd / 1000.0, 3) +
" m\n";
605 std::string config_file;
606 if (
dataWS->run().hasProperty(
"run_number")) {
607 const std::string run_str =
dataWS->run().getPropertyValueAsType<std::string>(
"run_number");
608 Poco::NumberParser::tryParse(run_str, run_number);
618 if (use_config && !config_file.empty()) {
626 }
else if (use_config) {
630 " Could not find configuration file for run " + Poco::NumberFormatter::format(run_number) +
"\n";
636 g_log.
error() <<
"Moderator position seems close to the sample, please check\n";
641 mvAlg->setProperty(
"ComponentName",
"moderator");
643 mvAlg->setProperty(
"RelativePosition",
false);
644 mvAlg->executeAsChildAlg();
653 if (reductionManager->existsProperty(
"LatestBeamCenterX") &&
654 reductionManager->existsProperty(
"LatestBeamCenterY")) {
655 m_center_x = reductionManager->getProperty(
"LatestBeamCenterX");
656 m_center_y = reductionManager->getProperty(
"LatestBeamCenterY");
664 if (!reductionManager->existsProperty(
"LatestBeamCenterX"))
667 reductionManager->setProperty(
"LatestBeamCenterX",
m_center_x);
668 if (!reductionManager->existsProperty(
"LatestBeamCenterY"))
671 reductionManager->setProperty(
"LatestBeamCenterY",
m_center_y);
675 const bool correct_for_flight_path =
getProperty(
"CorrectForFlightPath");
678 double wl_combined_max = 0.0;
679 if (skipTOFCorrection) {
680 m_output_message +=
" Skipping EQSANS TOF correction: assuming a single frame\n";
681 dataWS->mutableRun().addProperty(
"is_frame_skipping", 0,
true);
682 if (correct_for_flight_path) {
683 g_log.
error() <<
"CorrectForFlightPath and SkipTOFCorrection can't be "
684 "set to true at the same time\n";
689 if (!correct_for_flight_path)
697 tofAlg->setProperty(
"FlightPathCorrection", correct_for_flight_path);
698 tofAlg->executeAsChildAlg();
699 wl_min = tofAlg->getProperty(
"WavelengthMin");
700 wl_max = tofAlg->getProperty(
"WavelengthMax");
701 if (wl_min >= wl_max) {
706 const bool frame_skipping = tofAlg->getProperty(
"FrameSkipping");
707 dataWS->mutableRun().addProperty(
"wavelength_min", wl_min,
"Angstrom",
true);
708 dataWS->mutableRun().addProperty(
"wavelength_max", wl_max,
"Angstrom",
true);
709 dataWS->mutableRun().addProperty(
"is_frame_skipping",
int(frame_skipping),
true);
710 wl_combined_max = wl_max;
712 " Wavelength range: " + Poco::NumberFormatter::format(wl_min) +
" - " + Poco::NumberFormatter::format(wl_max);
713 if (frame_skipping) {
714 const double wl_min2 = tofAlg->getProperty(
"WavelengthMinFrame2");
715 const double wl_max2 = tofAlg->getProperty(
"WavelengthMaxFrame2");
716 wl_combined_max = wl_max2;
717 dataWS->mutableRun().addProperty(
"wavelength_min_frame2", wl_min2,
"Angstrom",
true);
718 dataWS->mutableRun().addProperty(
"wavelength_max_frame2", wl_max2,
"Angstrom",
true);
719 m_output_message +=
" and " + Poco::NumberFormatter::format(wl_min2) +
" - " +
720 Poco::NumberFormatter::format(wl_max2) +
" Angstrom\n";
728 const double ssd =
fabs(
dataWS->getInstrument()->getSource()->getPos().Z()) * 1000.0;
729 const double conversion_factor = 3.9560346 / (sfdd + ssd);
731 " TOF to wavelength conversion factor: " + Poco::NumberFormatter::format(conversion_factor) +
"\n";
733 if (skipTOFCorrection) {
735 if (dataWS_evt->getNumberEvents() == 0)
736 throw std::invalid_argument(
"No event to process: check your TOF cuts");
737 wl_min = dataWS_evt->getTofMin() * conversion_factor;
738 wl_max = dataWS_evt->getTofMax() * conversion_factor;
739 wl_combined_max = wl_max;
740 g_log.
information() <<
"Wavelength range: " << wl_min <<
" to " << wl_max <<
'\n';
741 dataWS->mutableRun().addProperty(
"wavelength_min", wl_min,
"Angstrom",
true);
742 dataWS->mutableRun().addProperty(
"wavelength_max", wl_max,
"Angstrom",
true);
748 scAlg->setProperty(
"Factor", conversion_factor);
749 scAlg->executeAsChildAlg();
750 dataWS->getAxis(0)->setUnit(
"Wavelength");
753 const bool preserveEvents =
getProperty(
"PreserveEvents");
754 const double wl_step =
getProperty(
"WavelengthStep");
756 const double wl_min_rounded = round(wl_min * 100.0) / 100.0;
757 const double wl_max_rounded = round(wl_combined_max * 100.0) / 100.0;
758 std::string params = Poco::NumberFormatter::format(wl_min_rounded, 2) +
"," + Poco::NumberFormatter::format(wl_step) +
759 "," + Poco::NumberFormatter::format(wl_max_rounded, 2);
765 rebinAlg->setPropertyValue(
"Params", params);
766 rebinAlg->setProperty(
"PreserveEvents", preserveEvents);
767 rebinAlg->executeAsChildAlg();
770 dataWS = rebinAlg->getProperty(
"OutputWorkspace");
773 setProperty<MatrixWorkspace_sptr>(
"OutputWorkspace", std::dynamic_pointer_cast<MatrixWorkspace>(
dataWS));
#define DECLARE_ALGORITHM(classname)
void declareProperty(std::unique_ptr< Kernel::Property > p, const std::string &doc="") override
Add a property to the list of managed properties.
std::string getPropertyValue(const std::string &name) const override
Get the value of a property as a string.
TypedValue getProperty(const std::string &name) const override
Get the value of a property.
std::string toString() const override
Serialize an object to a string.
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.
void setPropertyValue(const std::string &name, const std::string &value) override
Set the value of a property by string N.B.
static bool isEmpty(const NumT toCheck)
checks that the value was not set by users, uses the value in empty double/int.
@ OptionalLoad
to specify a file to read but the file doesn't have to exist
A property class for workspaces.
Records the filename and the description of failure.
Marks code as not implemented yet.
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.
void information(const std::string &msg)
Logs at information level.
The concrete, templated class for properties.
Base class for properties.
static T & Instance()
Return a reference to the Singleton instance, creating it if it does not already exist Creation is do...
A specialised Property class for holding a series of time-value pairs.
TimeSeriesPropertyStatistics getStatistics() const
Return a TimeSeriesPropertyStatistics object.
constexpr double X() const noexcept
Get x.
constexpr double Y() const noexcept
Get y.
void readBeamCenter(const std::string &line)
Read the beam center from a config file string.
double m_moderator_position
void readTOFcuts(const std::string &line)
Read the TOF cuts from a config file string.
void exec() override
Execution code.
API::MatrixWorkspace_sptr dataWS
void readModeratorPosition(const std::string &line)
Read the moderator position from a config file string.
std::string m_mask_as_string
double m_slit_positions[3][8]
void init() override
Initialisation code.
std::string m_output_message
void getSourceSlitSize()
Get the source slit size from the slit information of the run properties.
void readSourceSlitSize(const std::string &line)
Read the source slit sizes from a config file string.
void moveToBeamCenter()
Move the detector according to the beam center.
void readRectangularMasks(const std::string &line)
Read rectangular masks from a config file string.
void readConfigFile(const std::string &filePath)
Read a config file.
std::string findConfigFile(const int &run)
Find the most appropriate configuration file for a given run.
std::shared_ptr< Workspace > Workspace_sptr
shared pointer to Mantid::API::Workspace
std::shared_ptr< MatrixWorkspace > MatrixWorkspace_sptr
shared pointer to the matrix workspace base class
bool exists(::NeXus::File &file, const std::string &name)
Based on the current group in the file, does the named sub-entry exist?
std::shared_ptr< EventWorkspace > EventWorkspace_sptr
shared pointer to the EventWorkspace class
void getDefaultBeamCenter(const API::MatrixWorkspace_sptr &dataWS, double &pixel_x, double &pixel_y)
int getDetectorFromPixel(const int &pixel_x, const int &pixel_y, const API::MatrixWorkspace_sptr &dataWS)
void getCoordinateFromPixel(const double &pixel_x, const double &pixel_y, const API::MatrixWorkspace_sptr &dataWS, double &x, double &y)
double getRunPropertyDbl(const MatrixWorkspace_sptr &inputWS, const std::string &pname)
Returns the value of a run property from a given workspace.
constexpr double EMPTY_DBL() noexcept
Returns what we consider an "empty" double within a property.
@ Input
An input workspace.
@ Output
An output workspace.