Mantid
Loading...
Searching...
No Matches
LoadBBY.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 +
7#include <cmath>
8#include <cstdio>
9
11#include "MantidAPI/Axis.h"
15#include "MantidAPI/Run.h"
25
26#include <boost/algorithm/string.hpp>
27#include <boost/algorithm/string/trim.hpp>
28#include <boost/math/special_functions/round.hpp>
29
30#include <Poco/AutoPtr.h>
31#include <Poco/DOM/AutoPtr.h>
32#include <Poco/DOM/DOMParser.h>
33#include <Poco/DOM/Document.h>
34#include <Poco/DOM/Element.h>
35#include <Poco/DOM/NodeFilter.h>
36#include <Poco/DOM/NodeIterator.h>
37#include <Poco/DOM/NodeList.h>
38#include <Poco/TemporaryFile.h>
39#include <Poco/Util/PropertyFileConfiguration.h>
40
41namespace Mantid::DataHandling {
42
43// register the algorithm into the AlgorithmFactory
45
46// consts
47static const size_t HISTO_BINS_X = 240;
48static const size_t HISTO_BINS_Y = 256;
49// 100 = 40 + 20 + 40
50static const size_t Progress_LoadBinFile = 48;
51static const size_t Progress_ReserveMemory = 4;
53
54static char const *const FilenameStr = "Filename";
55static char const *const MaskStr = "Mask";
56
57static char const *const FilterByTofMinStr = "FilterByTofMin";
58static char const *const FilterByTofMaxStr = "FilterByTofMax";
59
60static char const *const FilterByTimeStartStr = "FilterByTimeStart";
61static char const *const FilterByTimeStopStr = "FilterByTimeStop";
62
63using ANSTO::EventVector_pt;
64
65template <typename TYPE>
66void AddSinglePointTimeSeriesProperty(API::LogManager &logManager, const std::string &time, const std::string &name,
67 const TYPE value) {
68 // create time series property and add single value
69 auto p = new Kernel::TimeSeriesProperty<TYPE>(name);
70 p->addValue(time, value);
71
72 // add to log manager
73 logManager.addProperty(p);
74}
75
83 if (descriptor.extension() != ".tar")
84 return 0;
85
86 ANSTO::Tar::File file(descriptor.filename());
87 if (!file.good())
88 return 0;
89
90 size_t hdfFiles = 0;
91 size_t binFiles = 0;
92 const std::vector<std::string> &subFiles = file.files();
93 for (const auto &subFile : subFiles) {
94 auto len = subFile.length();
95 if ((len > 4) && (subFile.find_first_of("\\/", 0, 2) == std::string::npos)) {
96 if ((subFile.rfind(".hdf") == len - 4) && (subFile.compare(0, 3, "BBY") == 0))
97 hdfFiles++;
98 else if (subFile.rfind(".bin") == len - 4)
99 binFiles++;
100 }
101 }
102
103 return (hdfFiles == 1) && (binFiles == 1) ? 50 : 0;
104}
111 // Specify file extensions which can be associated with a specific file.
112 std::vector<std::string> exts;
113
114 // Declare the Filename algorithm property. Mandatory. Sets the path to the
115 // file to load.
116 exts.clear();
117 exts.emplace_back(".tar");
118 declareProperty(std::make_unique<API::FileProperty>(FilenameStr, "", API::FileProperty::Load, exts),
119 "The input filename of the stored data");
120
121 // mask
122 exts.clear();
123 exts.emplace_back(".xml");
124 declareProperty(std::make_unique<API::FileProperty>(MaskStr, "", API::FileProperty::OptionalLoad, exts),
125 "The input filename of the mask data");
126
127 // OutputWorkspace
129 std::make_unique<API::WorkspaceProperty<API::IEventWorkspace>>("OutputWorkspace", "", Kernel::Direction::Output));
130
131 // FilterByTofMin
133 "Optional: To exclude events that do not fall within a range "
134 "of times-of-flight. "
135 "This is the minimum accepted value in microseconds. Keep "
136 "blank to load all events.");
137
138 // FilterByTofMax
141 "Optional: To exclude events that do not fall within a range "
142 "of times-of-flight. "
143 "This is the maximum accepted value in microseconds. Keep "
144 "blank to load all events.");
145
146 // FilterByTimeStart
149 "Optional: To only include events after the provided start time, in "
150 "seconds (relative to the start of the run).");
151
152 // FilterByTimeStop
155 "Optional: To only include events before the provided stop time, in "
156 "seconds (relative to the start of the run).");
157
158 std::string grpOptional = "Filters";
163}
168 // Delete the output workspace name if it existed
169 std::string outName = getPropertyValue("OutputWorkspace");
170 if (API::AnalysisDataService::Instance().doesExist(outName))
171 API::AnalysisDataService::Instance().remove(outName);
172
173 // Get the name of the data file.
174 std::string filename = getPropertyValue(FilenameStr);
175 ANSTO::Tar::File tarFile(filename);
176 if (!tarFile.good())
177 throw std::invalid_argument("invalid BBY file");
178
179 // region of intreset
180 std::vector<bool> roi = createRoiVector(getPropertyValue(MaskStr));
181
182 double tofMinBoundary = getProperty(FilterByTofMinStr);
183 double tofMaxBoundary = getProperty(FilterByTofMaxStr);
184
185 double timeMinBoundary = getProperty(FilterByTimeStartStr);
186 double timeMaxBoundary = getProperty(FilterByTimeStopStr);
187
188 if (isEmpty(tofMaxBoundary))
189 tofMaxBoundary = std::numeric_limits<double>::infinity();
190 if (isEmpty(timeMaxBoundary))
191 timeMaxBoundary = std::numeric_limits<double>::infinity();
192
193 API::Progress prog(this, 0.0, 1.0, Progress_Total);
194 prog.doReport("creating instrument");
195
196 // create workspace
197 DataObjects::EventWorkspace_sptr eventWS = std::make_shared<DataObjects::EventWorkspace>();
198
199 eventWS->initialize(HISTO_BINS_Y * HISTO_BINS_X,
200 2, // number of TOF bin boundaries
201 1);
202
203 // create instrument
204 InstrumentInfo instrumentInfo;
205 std::map<std::string, double> logParams;
206 std::map<std::string, std::string> logStrings;
207 std::map<std::string, std::string> allParams;
208 createInstrument(tarFile, instrumentInfo, logParams, logStrings, allParams);
209
210 // set the units
211 if (instrumentInfo.is_tof)
212 eventWS->getAxis(0)->unit() = Kernel::UnitFactory::Instance().create("TOF");
213 else
214 eventWS->getAxis(0)->unit() = Kernel::UnitFactory::Instance().create("Wavelength");
215
216 eventWS->setYUnit("Counts");
217
218 // set title
219 const std::vector<std::string> &subFiles = tarFile.files();
220 for (const auto &subFile : subFiles)
221 if (subFile.compare(0, 3, "BBY") == 0) {
222 std::string title = subFile;
223
224 if (title.rfind(".hdf") == title.length() - 4)
225 title.resize(title.length() - 4);
226
227 if (title.rfind(".nx") == title.length() - 3)
228 title.resize(title.length() - 3);
229
230 eventWS->setTitle(title);
231 break;
232 }
233
234 // load events
235 size_t numberHistograms = eventWS->getNumberHistograms();
236
237 std::vector<EventVector_pt> eventVectors(numberHistograms, nullptr);
238 std::vector<size_t> eventCounts(numberHistograms, 0);
239
240 // phase correction
241
242 double periodMaster = instrumentInfo.period_master;
243 double periodSlave = instrumentInfo.period_slave;
244 double phaseSlave = instrumentInfo.phase_slave;
245
246 double period = periodSlave;
247 double shift = -1.0 / 6.0 * periodMaster - periodSlave * phaseSlave / 360.0;
248
249 // get the start time from the file
250 Types::Core::DateAndTime startTime(instrumentInfo.start_time);
251 auto startInNanosec = startTime.totalNanoseconds();
252
253 // count total events per pixel to reserve necessary memory
254 ANSTO::EventCounter eventCounter(roi, HISTO_BINS_Y, period, shift, startInNanosec, tofMinBoundary, tofMaxBoundary,
255 timeMinBoundary, timeMaxBoundary, eventCounts);
256
257 loadEvents(prog, "loading neutron counts", tarFile, eventCounter);
258
259 // prepare event storage
260 ANSTO::ProgressTracker progTracker(prog, "creating neutron event lists", numberHistograms, Progress_ReserveMemory);
261
262 for (size_t i = 0; i != numberHistograms; ++i) {
263 DataObjects::EventList &eventList = eventWS->getSpectrum(i);
264
266 eventList.reserve(eventCounts[i]);
267
268 eventList.setDetectorID(static_cast<detid_t>(i));
269 eventList.setSpectrumNo(static_cast<detid_t>(i));
270
271 DataObjects::getEventsFrom(eventList, eventVectors[i]);
272
273 progTracker.update(i);
274 }
275 progTracker.complete();
276
277 if (instrumentInfo.is_tof) {
278 ANSTO::EventAssigner eventAssigner(roi, HISTO_BINS_Y, period, shift, startInNanosec, tofMinBoundary, tofMaxBoundary,
279 timeMinBoundary, timeMaxBoundary, eventVectors);
280
281 loadEvents(prog, "loading neutron events (TOF)", tarFile, eventAssigner);
282 } else {
283 ANSTO::EventAssignerFixedWavelength eventAssigner(roi, HISTO_BINS_Y, instrumentInfo.wavelength, period, shift,
284 startInNanosec, tofMinBoundary, tofMaxBoundary, timeMinBoundary,
285 timeMaxBoundary, eventVectors);
286
287 loadEvents(prog, "loading neutron events (Wavelength)", tarFile, eventAssigner);
288 }
289
290 auto getParam = [&allParams](const std::string &tag, double defValue) {
291 try {
292 return std::stod(allParams[tag]);
293 } catch (const std::invalid_argument &) {
294 return defValue;
295 }
296 };
297 if (instrumentInfo.is_tof) {
298 // just to make sure the bins hold it all
299 eventWS->setAllX(HistogramData::BinEdges{std::max(0.0, floor(eventCounter.tofMin())), eventCounter.tofMax() + 1});
300 } else {
301 double lof = getParam("wavelength_extn_lo", 0.95);
302 double hif = getParam("wavelength_extn_hi", 1.05);
303 eventWS->setAllX(HistogramData::BinEdges{instrumentInfo.wavelength * lof, instrumentInfo.wavelength * hif});
304 }
305
306 // count total number of masked bins
307 size_t maskedBins = 0;
308 for (size_t i = 0; i != roi.size(); i++)
309 if (!roi[i])
310 maskedBins++;
311
312 if (maskedBins > 0) {
313 // create list of masked bins
314 std::vector<size_t> maskIndexList(maskedBins);
315 size_t maskIndex = 0;
316
317 for (size_t i = 0; i != roi.size(); i++)
318 if (!roi[i])
319 maskIndexList[maskIndex++] = i;
320
321 auto maskingAlg = createChildAlgorithm("MaskDetectors");
322 maskingAlg->setProperty("Workspace", eventWS);
323 maskingAlg->setProperty("WorkspaceIndexList", maskIndexList);
324 maskingAlg->executeAsChildAlg();
325 }
326
327 // set log values
328 API::LogManager &logManager = eventWS->mutableRun();
329
330 auto frame_count = static_cast<int>(eventCounter.numFrames());
331
332 logManager.addProperty("filename", filename);
333 logManager.addProperty("att_pos", static_cast<int>(instrumentInfo.att_pos));
334 logManager.addProperty("frame_count", frame_count);
335 logManager.addProperty("period", period);
336
337 // currently beam monitor counts are not available, instead number of frames
338 // times period is used
339 logManager.addProperty("bm_counts", static_cast<double>(frame_count) * period /
340 1.0e6); // static_cast<double>(instrumentInfo.bm_counts)
341
342 Types::Core::time_duration duration =
343 boost::posix_time::microseconds(static_cast<boost::int64_t>(static_cast<double>(frame_count) * period));
344
345 Types::Core::DateAndTime start_time(instrumentInfo.start_time);
346 Types::Core::DateAndTime end_time(start_time + duration);
347
348 logManager.addProperty("start_time", start_time.toISO8601String());
349 logManager.addProperty("run_start", start_time.toISO8601String());
350 logManager.addProperty("end_time", end_time.toISO8601String());
351 logManager.addProperty("is_tof", instrumentInfo.is_tof);
352
353 std::string time_str = start_time.toISO8601String();
354
355 logManager.addProperty("sample_name", instrumentInfo.sample_name);
356 logManager.addProperty("sample_description", instrumentInfo.sample_description);
357 AddSinglePointTimeSeriesProperty(logManager, time_str, "wavelength", instrumentInfo.wavelength);
358 AddSinglePointTimeSeriesProperty(logManager, time_str, "master1_chopper_id", instrumentInfo.master1_chopper_id);
359 AddSinglePointTimeSeriesProperty(logManager, time_str, "master2_chopper_id", instrumentInfo.master2_chopper_id);
360
361 for (auto &x : logStrings) {
362 logManager.addProperty(x.first, x.second);
363 }
364 for (auto &x : logParams) {
365 AddSinglePointTimeSeriesProperty(logManager, time_str, x.first, x.second);
366 }
367
368 auto loadInstrumentAlg = createChildAlgorithm("LoadInstrument");
369 loadInstrumentAlg->setProperty("Workspace", eventWS);
370 loadInstrumentAlg->setPropertyValue("InstrumentName", "BILBY");
371 loadInstrumentAlg->setProperty("RewriteSpectraMap", Mantid::Kernel::OptionalBool(false));
372 loadInstrumentAlg->executeAsChildAlg();
373
374 setProperty("OutputWorkspace", eventWS);
375}
376
377// region of intreset
378std::vector<bool> LoadBBY::createRoiVector(const std::string &maskfile) {
379 std::vector<bool> result(HISTO_BINS_Y * HISTO_BINS_X, true);
380
381 if (maskfile.length() == 0)
382 return result;
383
384 std::ifstream input(maskfile.c_str());
385 if (!input.good())
386 throw std::invalid_argument("invalid mask file");
387
388 std::string line;
389 while (std::getline(input, line)) {
390 auto i0 = line.find("<detids>");
391 auto iN = line.find("</detids>");
392
393 if ((i0 != std::string::npos) && (iN != std::string::npos) && (i0 < iN)) {
394 line = line.substr(i0 + 8, iN - i0 - 8); // 8 = len("<detids>")
395 std::stringstream ss(line);
396
397 std::string item;
398 while (std::getline(ss, item, ',')) {
399 auto k = item.find('-');
400
401 size_t p0, p1;
402 if (k != std::string::npos) {
403 p0 = boost::lexical_cast<size_t>(item.substr(0, k));
404 p1 = boost::lexical_cast<size_t>(item.substr(k + 1, item.size() - k - 1));
405
406 if (p0 > p1)
407 std::swap(p0, p1);
408 } else {
409 p0 = boost::lexical_cast<size_t>(item);
410 p1 = p0;
411 }
412
413 if (p0 < result.size()) {
414 if (p1 >= result.size())
415 p1 = result.size() - 1;
416
417 while (p0 <= p1)
418 result[p0++] = false;
419 }
420 }
421 }
422 }
423
424 return result;
425}
426
427// loading instrument parameters
428void LoadBBY::loadInstrumentParameters(const NeXus::NXEntry &entry, std::map<std::string, double> &logParams,
429 std::map<std::string, std::string> &logStrings,
430 std::map<std::string, std::string> &allParams) {
431 using namespace Poco::XML;
432 std::string idfDirectory = Mantid::Kernel::ConfigService::Instance().getString("instrumentDefinition.directory");
433
434 try {
435 std::string parameterFilename = idfDirectory + "BILBY_Parameters.xml";
436 // Set up the DOM parser and parse xml file
437 DOMParser pParser;
438 Poco::AutoPtr<Poco::XML::Document> pDoc;
439 try {
440 pDoc = pParser.parse(parameterFilename);
441 } catch (...) {
442 throw Kernel::Exception::FileError("Unable to parse File:", parameterFilename);
443 }
444 NodeIterator it(pDoc, Poco::XML::NodeFilter::SHOW_ELEMENT);
445 Node *pNode = it.nextNode();
446 while (pNode) {
447 if (pNode->nodeName() == "parameter") {
448 auto pElem = dynamic_cast<Element *>(pNode);
449 std::string name = pElem->getAttribute("name");
450 Poco::AutoPtr<NodeList> nodeList = pElem->childNodes();
451 for (unsigned long i = 0; i < nodeList->length(); i++) {
452 auto cNode = nodeList->item(i);
453 if (cNode->nodeName() == "value") {
454 auto cElem = dynamic_cast<Poco::XML::Element *>(cNode);
455 std::string value = cElem->getAttribute("val");
456 allParams[name] = value;
457 }
458 }
459 }
460 pNode = it.nextNode();
461 }
462
463 auto isNumeric = [](const std::string &tag) {
464 try {
465 auto stag = boost::algorithm::trim_copy(tag);
466 size_t sz = 0;
467 auto value = std::stod(stag, &sz);
468 return sz > 0 && stag.size() == sz && std::isfinite(value);
469 } catch (const std::invalid_argument &) {
470 return false;
471 }
472 };
473
474 std::string tmpString;
475 float tmpFloat = 0.0f;
476 for (auto &x : allParams) {
477 if (x.first.find("log_") == 0) {
478 auto logTag = boost::algorithm::trim_copy(x.first.substr(4));
479 auto line = x.second;
480
481 // comma separated details
482 std::vector<std::string> details;
483 boost::split(details, line, boost::is_any_of(","));
484 if (details.size() < 3) {
485 g_log.warning() << "Invalid format for BILBY parameter " << x.first << std::endl;
486 continue;
487 }
488 auto hdfTag = boost::algorithm::trim_copy(details[0]);
489 try {
490 // extract the parameter and add it to the parameter dictionary,
491 // check the scale factor for numeric and string
492 auto updateOk = false;
493 if (!hdfTag.empty()) {
494 if (isNumeric(details[1])) {
495 if (loadNXDataSet(entry, hdfTag, tmpFloat)) {
496 auto factor = std::stod(details[1]);
497 logParams[logTag] = factor * tmpFloat;
498 updateOk = true;
499 }
500 } else if (loadNXString(entry, hdfTag, tmpString)) {
501 logStrings[logTag] = tmpString;
502 updateOk = true;
503 }
504 }
505 if (!updateOk) {
506 // if the hdf is missing the tag then add the default if
507 // it is provided
508 auto defValue = boost::algorithm::trim_copy(details[2]);
509 if (!defValue.empty()) {
510 if (isNumeric(defValue))
511 logParams[logTag] = std::stod(defValue);
512 else
513 logStrings[logTag] = defValue;
514 if (!hdfTag.empty())
515 g_log.warning() << "Cannot find hdf parameter " << hdfTag << ", using default.\n";
516 }
517 }
518 } catch (const std::invalid_argument &) {
519 g_log.warning() << "Invalid format for BILBY parameter " << x.first << std::endl;
520 }
521 }
522 }
523 } catch (std::exception &ex) {
524 g_log.warning() << "Failed to load instrument with error: " << ex.what()
525 << ". The current facility may not be fully "
526 "supported.\n";
527 }
528}
529
530// instrument creation
532 std::map<std::string, double> &logParams, std::map<std::string, std::string> &logStrings,
533 std::map<std::string, std::string> &allParams) {
534
535 const double toMeters = 1.0 / 1000;
536
537 instrumentInfo.sample_name = "UNKNOWN";
538 instrumentInfo.sample_description = "UNKNOWN";
539 instrumentInfo.start_time = "2000-01-01T00:00:00";
540
541 instrumentInfo.bm_counts = 0;
542 instrumentInfo.att_pos = 0;
543 instrumentInfo.master1_chopper_id = -1;
544 instrumentInfo.master2_chopper_id = -1;
545
546 instrumentInfo.is_tof = true;
547 instrumentInfo.wavelength = 0.0;
548
549 instrumentInfo.period_master = 0.0;
550 instrumentInfo.period_slave = (1.0 / 50.0) * 1.0e6;
551 instrumentInfo.phase_slave = 0.0;
552
553 // extract log and hdf file
554 const std::vector<std::string> &files = tarFile.files();
555 auto file_it = std::find_if(files.cbegin(), files.cend(),
556 [](const std::string &file) { return file.rfind(".hdf") == file.length() - 4; });
557 if (file_it != files.end()) {
558 tarFile.select(file_it->c_str());
559 // extract hdf file into tmp file
560 Poco::TemporaryFile hdfFile;
561 std::shared_ptr<FILE> handle(fopen(hdfFile.path().c_str(), "wb"), fclose);
562 if (handle) {
563 // copy content
564 char buffer[4096];
565 size_t bytesRead;
566 while (0 != (bytesRead = tarFile.read(buffer, sizeof(buffer))))
567 fwrite(buffer, bytesRead, 1, handle.get());
568 handle.reset();
569
570 NeXus::NXRoot root(hdfFile.path());
571 NeXus::NXEntry entry = root.openFirstEntry();
572
573 float tmp_float = 0.0f;
574 int32_t tmp_int32 = 0;
575 std::string tmp_str;
576
577 if (loadNXDataSet(entry, "monitor/bm1_counts", tmp_int32))
578 instrumentInfo.bm_counts = tmp_int32;
579 if (loadNXDataSet(entry, "instrument/att_pos", tmp_float))
580 instrumentInfo.att_pos = boost::math::iround(tmp_float); // [1.0, 2.0, ..., 5.0]
581
582 if (loadNXString(entry, "sample/name", tmp_str))
583 instrumentInfo.sample_name = tmp_str;
584 if (loadNXString(entry, "sample/description", tmp_str))
585 instrumentInfo.sample_description = tmp_str;
586 if (loadNXString(entry, "start_time", tmp_str))
587 instrumentInfo.start_time = tmp_str;
588
589 if (loadNXDataSet(entry, "instrument/master1_chopper_id", tmp_int32))
590 instrumentInfo.master1_chopper_id = tmp_int32;
591 if (loadNXDataSet(entry, "instrument/master2_chopper_id", tmp_int32))
592 instrumentInfo.master2_chopper_id = tmp_int32;
593
594 if (loadNXString(entry, "instrument/detector/frame_source", tmp_str))
595 instrumentInfo.is_tof = tmp_str == "EXTERNAL";
596 if (loadNXDataSet(entry, "instrument/nvs067/lambda", tmp_float))
597 instrumentInfo.wavelength = tmp_float;
598
599 if (loadNXDataSet(entry, "instrument/master_chopper_freq", tmp_float) && (tmp_float > 0.0f))
600 instrumentInfo.period_master = 1.0 / tmp_float * 1.0e6;
601 if (loadNXDataSet(entry, "instrument/t0_chopper_freq", tmp_float) && (tmp_float > 0.0f))
602 instrumentInfo.period_slave = 1.0 / tmp_float * 1.0e6;
603 if (loadNXDataSet(entry, "instrument/t0_chopper_phase", tmp_float))
604 instrumentInfo.phase_slave = tmp_float < 999.0 ? tmp_float : 0.0;
605
606 loadInstrumentParameters(entry, logParams, logStrings, allParams);
607
608 // Ltof_det_value is not present for monochromatic data so check
609 // and replace with default
610 auto findLtof = logParams.find("Ltof_det_value");
611 if (findLtof != logParams.end()) {
612 logParams["L1_chopper_value"] = logParams["Ltof_det_value"] - logParams["L2_det_value"];
613 } else {
614 logParams["L1_chopper_value"] = 18.4726;
615 g_log.warning() << "Cannot recover parameter 'L1_chopper_value'"
616 << ", using default.\n";
617 }
618 }
619 }
620
621 // patching
622 file_it = std::find(files.cbegin(), files.cend(), "History.log");
623 if (file_it != files.cend()) {
624 tarFile.select(file_it->c_str());
625 std::string logContent;
626 logContent.resize(tarFile.selected_size());
627 tarFile.read(&logContent[0], logContent.size());
628 std::istringstream data(logContent);
629 Poco::AutoPtr<Poco::Util::PropertyFileConfiguration> conf(new Poco::Util::PropertyFileConfiguration(data));
630
631 if (conf->hasProperty("Bm1Counts"))
632 instrumentInfo.bm_counts = conf->getInt("Bm1Counts");
633 if (conf->hasProperty("AttPos"))
634 instrumentInfo.att_pos = boost::math::iround(conf->getDouble("AttPos"));
635
636 if (conf->hasProperty("SampleName"))
637 instrumentInfo.sample_name = conf->getString("SampleName");
638
639 if (conf->hasProperty("MasterChopperFreq")) {
640 auto tmp = conf->getDouble("MasterChopperFreq");
641 if (tmp > 0.0f)
642 instrumentInfo.period_master = 1.0 / tmp * 1.0e6;
643 }
644 if (conf->hasProperty("T0ChopperFreq")) {
645 auto tmp = conf->getDouble("T0ChopperFreq");
646 if (tmp > 0.0f)
647 instrumentInfo.period_slave = 1.0 / tmp * 1.0e6;
648 }
649 if (conf->hasProperty("T0ChopperPhase")) {
650 auto tmp = conf->getDouble("T0ChopperPhase");
651 instrumentInfo.phase_slave = tmp < 999.0 ? tmp : 0.0;
652 }
653
654 if (conf->hasProperty("FrameSource"))
655 instrumentInfo.is_tof = conf->getString("FrameSource") == "EXTERNAL";
656 if (conf->hasProperty("Wavelength"))
657 instrumentInfo.wavelength = conf->getDouble("Wavelength");
658
659 if (conf->hasProperty("SampleAperture"))
660 logParams["sample_aperture"] = conf->getDouble("SampleAperture");
661 if (conf->hasProperty("SourceAperture"))
662 logParams["source_aperture"] = conf->getDouble("SourceAperture");
663 if (conf->hasProperty("L1"))
664 logParams["L1_source_value"] = conf->getDouble("L1") * toMeters;
665 if (conf->hasProperty("LTofDet"))
666 logParams["L1_chopper_value"] = conf->getDouble("LTofDet") * toMeters - logParams["L2_det_value"];
667 if (conf->hasProperty("L2Det"))
668 logParams["L2_det_value"] = conf->getDouble("L2Det") * toMeters;
669
670 if (conf->hasProperty("L2CurtainL"))
671 logParams["L2_curtainl_value"] = conf->getDouble("L2CurtainL") * toMeters;
672 if (conf->hasProperty("L2CurtainR"))
673 logParams["L2_curtainr_value"] = conf->getDouble("L2CurtainR") * toMeters;
674 if (conf->hasProperty("L2CurtainU"))
675 logParams["L2_curtainu_value"] = conf->getDouble("L2CurtainU") * toMeters;
676 if (conf->hasProperty("L2CurtainD"))
677 logParams["L2_curtaind_value"] = conf->getDouble("L2CurtainD") * toMeters;
678
679 if (conf->hasProperty("CurtainL"))
680 logParams["D_curtainl_value"] = conf->getDouble("CurtainL") * toMeters;
681 if (conf->hasProperty("CurtainR"))
682 logParams["D_curtainr_value"] = conf->getDouble("CurtainR") * toMeters;
683 if (conf->hasProperty("CurtainU"))
684 logParams["D_curtainu_value"] = conf->getDouble("CurtainU") * toMeters;
685 if (conf->hasProperty("CurtainD"))
686 logParams["D_curtaind_value"] = conf->getDouble("CurtainD") * toMeters;
687 }
688}
689
690// load nx dataset
691template <class T> bool LoadBBY::loadNXDataSet(const NeXus::NXEntry &entry, const std::string &path, T &value) {
692 try {
693 NeXus::NXDataSetTyped<T> dataSet = entry.openNXDataSet<T>(path);
694 dataSet.load();
695
696 value = *dataSet();
697 return true;
698 } catch (std::runtime_error &) {
699 return false;
700 }
701}
702bool LoadBBY::loadNXString(const NeXus::NXEntry &entry, const std::string &path, std::string &value) {
703 try {
704 NeXus::NXChar dataSet = entry.openNXChar(path);
705 dataSet.load();
706
707 value = std::string(dataSet(), dataSet.dim0());
708 return true;
709 } catch (std::runtime_error &) {
710 return false;
711 }
712}
713
714// read counts/events from binary file
715template <class EventProcessor>
716void LoadBBY::loadEvents(API::Progress &prog, const char *progMsg, ANSTO::Tar::File &tarFile,
717 EventProcessor &eventProcessor) {
718 prog.doReport(progMsg);
719
720 // select bin file
721 int64_t fileSize = 0;
722 const std::vector<std::string> &files = tarFile.files();
723 const auto found = std::find_if(files.cbegin(), files.cend(),
724 [](const auto &file) { return file.rfind(".bin") == file.length() - 4; });
725 if (found != files.cend()) {
726 tarFile.select(found->c_str());
727 fileSize = tarFile.selected_size();
728 }
729
730 // for progress notifications
731 ANSTO::ProgressTracker progTracker(prog, progMsg, fileSize, Progress_LoadBinFile);
732
733 uint32_t x = 0; // 9 bits [0-239] tube number
734 uint32_t y = 0; // 8 bits [0-255] position along tube
735
736 // uint v = 0; // 0 bits [ ]
737 // uint w = 0; // 0 bits [ ] energy
738 uint32_t dt = 0;
739 double tof = 0.0;
740
741 if ((fileSize == 0) || !tarFile.skip(128))
742 return;
743
744 int state = 0;
745 int invalidEvents = 0;
746 uint32_t c;
747 while ((c = static_cast<uint32_t>(tarFile.read_byte())) != static_cast<uint32_t>(-1)) {
748
749 bool event_ended = false;
750 switch (state) {
751 case 0:
752 x = (c & 0xFF) >> 0; // set bit 1-8
753 break;
754
755 case 1:
756 x |= (c & 0x01) << 8; // set bit 9
757 y = (c & 0xFE) >> 1; // set bit 1-7
758 break;
759
760 case 2:
761 event_ended = (c & 0xC0) != 0xC0;
762 if (!event_ended)
763 c &= 0x3F;
764
765 y |= (c & 0x01) << 7; // set bit 8
766 dt = (c & 0xFE) >> 1; // set bit 1-5(7)
767 break;
768
769 case 3:
770 case 4:
771 case 5:
772 case 6:
773 case 7:
774 event_ended = (c & 0xC0) != 0xC0;
775 if (!event_ended)
776 c &= 0x3F;
777
778 // state is either 3, 4, 5, 6 or 7
779 dt |= (c & 0xFF) << (5 + 6 * (state - 3)); // set bit 6...
780 break;
781 }
782 state++;
783
784 if (event_ended || (state == 8)) {
785 state = 0;
786
787 if ((x == 0) && (y == 0) && (dt == 0xFFFFFFFF)) {
788 tof = 0.0;
789 eventProcessor.newFrame();
790 } else if ((x >= HISTO_BINS_X) || (y >= HISTO_BINS_Y)) {
791 // cannot ignore the dt contrbition even if the event
792 // is out of bounds as it is used in the encoding process
793 tof += static_cast<int>(dt) * 0.1;
794 invalidEvents++;
795 } else {
796 // conversion from 100 nanoseconds to 1 microsecond
797 tof += static_cast<int>(dt) * 0.1;
798
799 eventProcessor.addEvent(x, y, tof);
800 }
801
802 progTracker.update(tarFile.selected_position());
803 }
804 }
805 if (invalidEvents > 0) {
806 g_log.warning() << "BILBY loader dropped " << invalidEvents << " invalid event(s)" << std::endl;
807 }
808}
809
810} // namespace Mantid::DataHandling
gsl_vector * tmp
double value
The value of the point.
Definition: FitMW.cpp:51
#define DECLARE_FILELOADER_ALGORITHM(classname)
DECLARE_FILELOADER_ALGORITHM should be used in place of the standard DECLARE_ALGORITHM macro when wri...
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
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 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
Definition: FileProperty.h:53
@ Load
allowed here which will be passed to the algorithm
Definition: FileProperty.h:52
void setDetectorID(const detid_t detID)
Clear the list of detector IDs, then add one.
Definition: ISpectrum.cpp:84
void setSpectrumNo(specnum_t num)
Sets the the spectrum number of this spectrum.
Definition: ISpectrum.cpp:127
This class contains the information about the log entries.
Definition: LogManager.h:44
void addProperty(Kernel::Property *prop, bool overwrite=false)
Add data to the object in the form of a property.
Definition: LogManager.h:79
Helper class for reporting progress from algorithms.
Definition: Progress.h:25
void doReport(const std::string &msg="") override
Actually do the reporting, without changing the loop counter.
Definition: Progress.cpp:70
A property class for workspaces.
helper class to keep track of progress
size_t read(void *dst, size_t size)
const std::vector< std::string > & files() const
static bool loadNXDataSet(const NeXus::NXEntry &entry, const std::string &path, T &value)
Definition: LoadBBY.cpp:691
static std::vector< bool > createRoiVector(const std::string &maskfile)
Definition: LoadBBY.cpp:378
void init() override
Initialise the algorithm.
Definition: LoadBBY.cpp:110
void createInstrument(ANSTO::Tar::File &tarFile, InstrumentInfo &instrumentInfo, std::map< std::string, double > &logParams, std::map< std::string, std::string > &logStrings, std::map< std::string, std::string > &allParams)
Definition: LoadBBY.cpp:531
const std::string name() const override
function to return a name of the algorithm, must be overridden in all algorithms
Definition: LoadBBY.h:54
void loadEvents(API::Progress &prog, const char *progMsg, ANSTO::Tar::File &tarFile, EventProcessor &eventProcessor)
Definition: LoadBBY.cpp:716
int confidence(Kernel::FileDescriptor &descriptor) const override
Return the confidence value that this algorithm can load the file.
Definition: LoadBBY.cpp:82
void loadInstrumentParameters(const NeXus::NXEntry &entry, std::map< std::string, double > &logParams, std::map< std::string, std::string > &logStrings, std::map< std::string, std::string > &allParams)
Definition: LoadBBY.cpp:428
void exec() override
Execute the algorithm.
Definition: LoadBBY.cpp:167
bool loadNXString(const NeXus::NXEntry &entry, const std::string &path, std::string &value)
Definition: LoadBBY.cpp:702
A class for holding :
Definition: EventList.h:56
void setSortOrder(const EventSortType order) const
Manually set the event list sort order value.
Definition: EventList.cpp:943
void reserve(size_t num) override
Reserve a certain number of entries in event list of the specified eventType.
Definition: EventList.cpp:895
Records the filename and the description of failure.
Definition: Exception.h:98
Defines a wrapper around an open file.
const std::string & filename() const
Access the filename.
const std::string & extension() const
Access the file extension.
IPropertyManager * setProperty(const std::string &name, const T &value)
Templated method to set the value of a PropertyWithValue.
void setPropertyGroup(const std::string &name, const std::string &group)
Set the group for a given property.
void warning(const std::string &msg)
Logs at warning level.
Definition: Logger.cpp:86
OptionalBool : Tri-state bool.
Definition: OptionalBool.h:25
The concrete, templated 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.
NXDataSetTyped< T > openNXDataSet(const std::string &name) const
Templated method for creating datasets.
Definition: NexusClasses.h:536
NXChar openNXChar(const std::string &name) const
Creates and opens a char dataset.
Definition: NexusClasses.h:561
Templated class implementation of NXDataSet.
Definition: NexusClasses.h:203
void load(const int blocksize=1, int i=-1, int j=-1, int k=-1, int l=-1) override
Implementation of the virtual NXDataSet::load(...) method.
Definition: NexusClasses.h:289
int dim0() const
Returns the number of elements along the first dimension.
Implements NXentry Nexus class.
Definition: NexusClasses.h:898
Implements NXroot Nexus class.
Definition: NexusClasses.h:922
NXEntry openFirstEntry()
Open the first NXentry in the file.
Kernel::Logger g_log("ExperimentInfo")
static logger object
static char const *const FilenameStr
Definition: LoadBBY.cpp:54
static char const *const FilterByTofMaxStr
Definition: LoadBBY.cpp:58
static const size_t HISTO_BINS_X
Definition: LoadBBY.cpp:47
static const size_t HISTO_BINS_Y
Definition: LoadBBY.cpp:48
static const size_t Progress_LoadBinFile
Definition: LoadBBY.cpp:50
static char const *const FilterByTofMinStr
Definition: LoadBBY.cpp:57
void AddSinglePointTimeSeriesProperty(API::LogManager &logManager, const std::string &time, const std::string &name, const TYPE value)
Definition: LoadBBY.cpp:66
static char const *const MaskStr
Definition: LoadBBY.cpp:55
static char const *const FilterByTimeStartStr
Definition: LoadBBY.cpp:60
static const size_t Progress_Total
Definition: LoadBBY.cpp:52
static char const *const FilterByTimeStopStr
Definition: LoadBBY.cpp:61
static const size_t Progress_ReserveMemory
Definition: LoadBBY.cpp:51
DLLExport void getEventsFrom(EventList &el, std::vector< Types::Event::TofEvent > *&events)
std::shared_ptr< EventWorkspace > EventWorkspace_sptr
shared pointer to the EventWorkspace class
int32_t detid_t
Typedef for a detector ID.
Definition: SpectrumInfo.h:21
constexpr double EMPTY_DBL() noexcept
Returns what we consider an "empty" double within a property.
Definition: EmptyValues.h:43
STL namespace.
@ Input
An input workspace.
Definition: Property.h:53
@ Output
An output workspace.
Definition: Property.h:54