Initial commit of LoadILLTOF2, version 2 of ILL's TOF loader.

Re #18091

Framework/DataHandling/CMakeLists.txt

@@ -56,6 +56,7 @@ set ( SRC_FILES
 	src/LoadILL.cpp
 	src/LoadILLIndirect.cpp
 	src/LoadILLTOF.cpp
+	src/LoadILLTOF2.cpp
 	src/LoadILLReflectometry.cpp
 	src/LoadILLSANS.cpp
 	src/LoadISISNexus2.cpp
@@ -234,6 +235,7 @@ set ( INC_FILES
 	inc/MantidDataHandling/LoadILL.h
 	inc/MantidDataHandling/LoadILLIndirect.h
 	inc/MantidDataHandling/LoadILLTOF.h
+	inc/MantidDataHandling/LoadILLTOF2.h
 	inc/MantidDataHandling/LoadILLReflectometry.h
 	inc/MantidDataHandling/LoadILLSANS.h
 	inc/MantidDataHandling/LoadISISNexus2.h
@@ -405,6 +407,7 @@ set ( TEST_FILES
 	LoadIDFFromNexusTest.h
 	LoadILLIndirectTest.h
 	LoadILLTOFTest.h
+	LoadILLTOF2Test.h
 	LoadILLReflectometryTest.h
 	LoadILLSANSTest.h
 	LoadISISNexusTest.h

Framework/DataHandling/inc/MantidDataHandling/LoadILLTOF2.h

+#ifndef MANTID_DATAHANDLING_LOADILLTOF2_H_
+#define MANTID_DATAHANDLING_LOADILLTOF2_H_
+
+//---------------------------------------------------
+// Includes
+//---------------------------------------------------
+#include "MantidAPI/IFileLoader.h"
+#include "MantidNexus/NexusClasses.h"
+#include "MantidDataHandling/LoadHelper.h"
+#include "MantidGeometry/IDTypes.h"
+#include "MantidAPI/Progress.h"
+
+namespace Mantid {
+namespace DataHandling {
+/**
+ Loads an ILL IN4/5/6 nexus file into a Mantid workspace.
+
+ Copyright © 2016 ISIS Rutherford Appleton Laboratory, NScD Oak Ridge
+ National Laboratory & European Spallation Source
+
+ This file is part of Mantid.
+
+ Mantid is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 3 of the License, or
+ (at your option) any later version.
+
+ Mantid is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ GNU General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with this program.  If not, see <http://www.gnu.org/licenses/>.
+
+ File change history is stored at: <https://github.com/mantidproject/mantid>
+ Code Documentation is available at: <http://doxygen.mantidproject.org>
+ */
+class DLLExport LoadILLTOF2 : public API::IFileLoader<Kernel::NexusDescriptor> {
+public:
+  /// Constructor
+  LoadILLTOF2();
+  /// Algorithm's name
+  const std::string name() const override { return "LoadILLTOF"; }
+  /// Summary of algorithms purpose
+  const std::string summary() const override {
+    return "Loads an ILL TOF NeXus file.";
+  }
+
+  /// Algorithm's version
+  int version() const override { return 2; }
+  /// Algorithm's category for identification
+  const std::string category() const override { return "DataHandling\\Nexus"; }
+
+  /// Returns a confidence value that this algorithm can load a file
+  int confidence(Kernel::NexusDescriptor &descriptor) const override;
+
+private:
+  const static std::vector<std::string> SUPPORTED_INSTRUMENTS;
+
+  // Initialisation code
+  void init() override;
+  // Execution code
+  void exec() override;
+
+  void loadInstrumentDetails(NeXus::NXEntry &);
+  std::vector<std::vector<int>> getMonitorInfo(NeXus::NXEntry &firstEntry);
+  void initWorkSpace(NeXus::NXEntry &entry,
+                     const std::vector<std::vector<int>> &);
+  void initInstrumentSpecific();
+  void addAllNexusFieldsAsProperties(std::string filename);
+  void addEnergyToRun();
+  void addFacility();
+  void addPulseInterval();
+
+  void loadTimeDetails(NeXus::NXEntry &entry);
+  void
+  loadDataIntoTheWorkSpace(NeXus::NXEntry &entry,
+                           const std::vector<std::vector<int>> &);
+  void loadSpectra(size_t &spec, size_t numberOfMonitors, size_t numberOfTubes,
+                   std::vector<Mantid::detid_t> &detectorIDs, NeXus::NXInt data,
+                   Mantid::API::Progress progress);
+
+  void runLoadInstrument();
+
+  /// Calculate error for y
+  static double calculateError(double in) { return sqrt(in); }
+
+  API::MatrixWorkspace_sptr m_localWorkspace;
+
+  std::string m_instrumentName = ""; ///< Name of the instrument
+  std::string m_instrumentPath = ""; ///< Name of the instrument path
+
+  // Variables describing the data in the detector
+  size_t m_numberOfTubes = 0;         // number of tubes - X
+  size_t m_numberOfPixelsPerTube = 0; // number of pixels per tube - Y
+  size_t m_numberOfChannels = 0;      // time channels - Z
+  size_t m_numberOfHistograms = 0;
+
+  /* Values parsed from the nexus file */
+  double m_wavelength = 0;
+  double m_channelWidth = 0;
+  double m_timeOfFlightDelay = 0;
+
+  LoadHelper m_loader;
+};
+
+} // namespace DataHandling
+} // namespace Mantid
+
+#endif /*MANTID_DATAHANDLING_LOADILLTOF2_H_*/

Framework/DataHandling/src/LoadILLTOF.cpp

@@ -32,22 +32,9 @@ DECLARE_NEXUS_FILELOADER_ALGORITHM(LoadILLTOF)
  * be used
  */
 int LoadILLTOF::confidence(Kernel::NexusDescriptor &descriptor) const {
-
-  // fields existent only at the ILL
-  if (descriptor.pathExists("/entry0/wavelength") &&
-      descriptor.pathExists("/entry0/experiment_identifier") &&
-      descriptor.pathExists("/entry0/mode") &&
-      !descriptor.pathExists(
-          "/entry0/dataSD") // This one is for LoadILLIndirect
-      &&
-      !descriptor.pathExists(
-          "/entry0/instrument/VirtualChopper") // This one is for
-                                               // LoadILLReflectometry
-      ) {
-    return 80;
-  } else {
-    return 0;
-  }
+  UNUSED_ARG( descriptor )
+  // This loader is deprecated.
+  return 0;
 }
 
 LoadILLTOF::LoadILLTOF() : API::IFileLoader<Kernel::NexusDescriptor>() {

Framework/DataHandling/src/LoadILLTOF2.cpp

+#include "MantidDataHandling/LoadILLTOF2.h"
+
+#include "MantidAPI/Axis.h"
+#include "MantidAPI/FileProperty.h"
+#include "MantidAPI/MatrixWorkspace.h"
+#include "MantidAPI/RegisterFileLoader.h"
+#include "MantidAPI/WorkspaceFactory.h"
+#include "MantidGeometry/Instrument.h"
+#include "MantidKernel/UnitFactory.h"
+
+namespace Mantid {
+namespace DataHandling {
+
+using namespace Kernel;
+using namespace API;
+using namespace NeXus;
+using namespace HistogramData;
+
+DECLARE_NEXUS_FILELOADER_ALGORITHM(LoadILLTOF2)
+
+/// A vector containing the supported instrument names
+const std::vector<std::string> LoadILLTOF2::SUPPORTED_INSTRUMENTS = {"IN4", "IN5", "IN6"};
+
+/**
+ * Return the confidence with with this algorithm can load the file
+ *
+ * @param descriptor A descriptor for the file
+ *
+ * @return An integer specifying the confidence level. 0 indicates it will not
+ * be used
+ */
+int LoadILLTOF2::confidence(Kernel::NexusDescriptor &descriptor) const {
+
+  // fields existent only at the ILL
+  if (descriptor.pathExists("/entry0/wavelength") &&
+      descriptor.pathExists("/entry0/experiment_identifier") &&
+      descriptor.pathExists("/entry0/mode") &&
+      !descriptor.pathExists(
+          "/entry0/dataSD") // This one is for LoadILLIndirect
+      &&
+      !descriptor.pathExists(
+          "/entry0/instrument/VirtualChopper") // This one is for
+                                               // LoadILLReflectometry
+      ) {
+    return 80;
+  } else {
+    return 0;
+  }
+}
+
+LoadILLTOF2::LoadILLTOF2() : API::IFileLoader<Kernel::NexusDescriptor>() {
+}
+
+/**
+ * Initialises the algorithm
+ */
+void LoadILLTOF2::init() {
+  declareProperty(
+      make_unique<FileProperty>("Filename", "", FileProperty::Load, ".nxs"),
+      "File path of the Data file to load");
+
+  declareProperty(make_unique<WorkspaceProperty<>>("OutputWorkspace", "",
+                                                   Direction::Output),
+                  "The name to use for the output workspace");
+}
+
+/**
+ * Executes the algorithm
+ */
+void LoadILLTOF2::exec() {
+  // Retrieve filename
+  std::string filenameData = getPropertyValue("Filename");
+
+  // open the root node
+  NeXus::NXRoot dataRoot(filenameData);
+  NXEntry dataFirstEntry = dataRoot.openFirstEntry();
+
+  loadInstrumentDetails(dataFirstEntry);
+  loadTimeDetails(dataFirstEntry);
+
+  std::vector<std::vector<int>> monitors = getMonitorInfo(dataFirstEntry);
+
+  initWorkSpace(dataFirstEntry, monitors);
+
+  addAllNexusFieldsAsProperties(filenameData);
+  addFacility();
+
+  runLoadInstrument(); // just to get IDF contents
+
+  loadDataIntoTheWorkSpace(dataFirstEntry, monitors);
+
+  addEnergyToRun();
+  addPulseInterval();
+
+  // load the instrument from the IDF if it exists
+  runLoadInstrument();
+
+  // Set the output workspace property
+  setProperty("OutputWorkspace", m_localWorkspace);
+}
+
+/**
+ * Loads Monitor data into an vector of monitor data
+ *
+ * @param firstEntry The NeXus entry
+ *
+ * @return List of monitor data
+ */
+std::vector<std::vector<int>>
+LoadILLTOF2::getMonitorInfo(NeXus::NXEntry &firstEntry) {
+
+  std::vector<std::vector<int>> monitorList;
+
+  for (std::vector<NXClassInfo>::const_iterator it =
+           firstEntry.groups().begin();
+       it != firstEntry.groups().end(); ++it) {
+
+    if (it->nxclass == "NXmonitor" ||
+        boost::starts_with(it->nxname, "monitor")) {
+
+      g_log.debug() << "Load monitor data from " + it->nxname;
+
+      NXData dataGroup = firstEntry.openNXData(it->nxname + "/data");
+      NXInt data = dataGroup.openIntData();
+      // load the counts from the file into memory
+      data.load();
+
+      std::vector<int> thisMonitor(data(), data() + data.size());
+      monitorList.push_back(thisMonitor);
+    }
+  }
+  return monitorList;
+}
+
+/**
+ * Sets the instrument name along with its path in the nexus file
+ *
+ * @param firstEntry The NeXus entry
+ */
+void LoadILLTOF2::loadInstrumentDetails(NeXus::NXEntry &firstEntry) {
+
+  m_instrumentPath = m_loader.findInstrumentNexusPath(firstEntry);
+
+  if (m_instrumentPath == "") {
+    throw std::runtime_error(
+        "Cannot set the instrument name from the Nexus file!");
+  }
+
+  m_instrumentName =
+      m_loader.getStringFromNexusPath(firstEntry, m_instrumentPath + "/name");
+
+  if (std::find(SUPPORTED_INSTRUMENTS.begin(), SUPPORTED_INSTRUMENTS.end(),
+                m_instrumentName) == SUPPORTED_INSTRUMENTS.end()) {
+    std::string message =
+        "The instrument " + m_instrumentName + " is not valid for this loader!";
+    throw std::runtime_error(message);
+  }
+
+  g_log.debug() << "Instrument name set to: " + m_instrumentName << '\n';
+}
+
+/**
+ * Creates the workspace and initialises member variables with
+ * the corresponding values
+ *
+ * @param entry The NeXus entry
+ * @param monitors List of monitor data
+ */
+void LoadILLTOF2::initWorkSpace(NeXus::NXEntry &entry,
+                               const std::vector<std::vector<int>> &monitors) {
+
+  // read in the data
+  NXData dataGroup = entry.openNXData("data");
+  NXInt data = dataGroup.openIntData();
+
+  m_numberOfTubes = static_cast<size_t>(data.dim0());
+  m_numberOfPixelsPerTube = static_cast<size_t>(data.dim1());
+  m_numberOfChannels = static_cast<size_t>(data.dim2());
+  size_t numberOfMonitors = monitors.size();
+
+  /**
+   * IN4 : Rosace detector is in a different field.
+   */
+  size_t numberOfTubesInRosace = 0;
+  if (m_instrumentName == "IN4") {
+    NXData dataGroupRosace =
+        entry.openNXData("instrument/Detector_Rosace/data");
+    NXInt dataRosace = dataGroupRosace.openIntData();
+    numberOfTubesInRosace += static_cast<size_t>(dataRosace.dim0());
+  }
+
+  // dim0 * m_numberOfPixelsPerTube is the total number of detectors
+  m_numberOfHistograms =
+      (m_numberOfTubes + numberOfTubesInRosace) * m_numberOfPixelsPerTube;
+
+  g_log.debug() << "NumberOfTubes: " << m_numberOfTubes << '\n';
+  g_log.debug() << "NumberOfPixelsPerTube: " << m_numberOfPixelsPerTube << '\n';
+  g_log.debug() << "NumberOfChannels: " << m_numberOfChannels << '\n';
+
+  // Now create the output workspace
+  // total number of spectra + number of monitors,
+  // bin boundaries = m_numberOfChannels + 1
+  // Z/time dimension
+  m_localWorkspace = WorkspaceFactory::Instance().create(
+      "Workspace2D", m_numberOfHistograms + numberOfMonitors,
+      m_numberOfChannels + 1, m_numberOfChannels);
+  m_localWorkspace->getAxis(0)->unit() = UnitFactory::Instance().create("TOF");
+  m_localWorkspace->setYUnitLabel("Counts");
+}
+
+/**
+ * Load the time details from the nexus file.
+ *
+ * @param entry :: The Nexus entry
+ */
+void LoadILLTOF2::loadTimeDetails(NeXus::NXEntry &entry) {
+
+  m_wavelength = entry.getFloat("wavelength");
+
+  // Monitor can be monitor (IN5) or monitor1 (IN6)
+  std::string monitorName;
+  if (entry.containsGroup("monitor"))
+    monitorName = "monitor";
+  else if (entry.containsGroup("monitor1"))
+    monitorName = "monitor1";
+  else {
+    std::string message("Cannot find monitor[1] in the Nexus file!");
+    g_log.error(message);
+    throw std::runtime_error(message);
+  }
+
+  NXFloat time_of_flight_data =
+      entry.openNXFloat(monitorName + "/time_of_flight");
+  time_of_flight_data.load();
+
+  // The entry "monitor/time_of_flight", has 3 fields:
+  // channel width , number of channels, Time of flight delay
+  m_channelWidth = time_of_flight_data[0];
+  m_timeOfFlightDelay = time_of_flight_data[2];
+
+  g_log.debug("Nexus Data:");
+  g_log.debug() << " ChannelWidth: " << m_channelWidth << '\n';
+  g_log.debug() << " TimeOfFlightDealy: " << m_timeOfFlightDelay << '\n';
+  g_log.debug() << " Wavelength: " << m_wavelength << '\n';
+}
+
+/**
+ * Goes through all the fields of the NeXus file and adds them
+ * as parameters in the workspace
+ *
+ * @param filename The NeXus file
+ */
+void LoadILLTOF2::addAllNexusFieldsAsProperties(std::string filename) {
+
+  API::Run &runDetails = m_localWorkspace->mutableRun();
+
+  // Open NeXus file
+  NXhandle nxfileID;
+  NXstatus stat = NXopen(filename.c_str(), NXACC_READ, &nxfileID);
+
+  g_log.debug() << "Starting parsing properties from : " << filename << '\n';
+  if (stat == NX_ERROR) {
+    g_log.debug() << "convertNexusToProperties: Error loading " << filename;
+    throw Kernel::Exception::FileError("Unable to open File:", filename);
+  }
+  m_loader.addNexusFieldsToWsRun(nxfileID, runDetails);
+
+  g_log.debug() << "End parsing properties from : " << filename << '\n';
+}
+
+/**
+ * Calculates the incident energy from the wavelength and adds
+ * it as sample log 'Ei'
+ */
+void LoadILLTOF2::addEnergyToRun() {
+
+  API::Run &runDetails = m_localWorkspace->mutableRun();
+  double ei = m_loader.calculateEnergy(m_wavelength);
+  runDetails.addProperty<double>("Ei", ei, true); // overwrite
+}
+
+/**
+ * Adds facility info to the sample logs
+ */
+void LoadILLTOF2::addFacility() {
+  API::Run &runDetails = m_localWorkspace->mutableRun();
+  runDetails.addProperty("Facility", std::string("ILL"));
+}
+
+/**
+ * Calculates and adds the pulse intervals for the run
+ */
+void LoadILLTOF2::addPulseInterval() {
+  API::Run &runDetails = m_localWorkspace->mutableRun();
+  double pulseInterval;
+  double n_pulses;
+  double fermiChopperSpeed;
+
+  if (m_instrumentName == "IN4") {
+    fermiChopperSpeed =
+        runDetails.getPropertyAsSingleValue("FC.rotation_speed");
+    double bkgChopper1Speed =
+        runDetails.getPropertyAsSingleValue("BC1.rotation_speed");
+    double bkgChopper2Speed =
+        runDetails.getPropertyAsSingleValue("BC2.rotation_speed");
+
+    if (std::abs(bkgChopper1Speed - bkgChopper2Speed) > 1) {
+      throw std::invalid_argument(
+          "Background choppers 1 and 2 have different speeds");
+    }
+
+    n_pulses = fermiChopperSpeed / bkgChopper1Speed / 4;
+  } else if (m_instrumentName == "IN6") {
+    fermiChopperSpeed =
+        runDetails.getPropertyAsSingleValue("Fermi.rotation_speed");
+    double suppressorSpeed =
+        runDetails.getPropertyAsSingleValue("Suppressor.rotation_speed");
+
+    n_pulses = fermiChopperSpeed / suppressorSpeed;
+  } else {
+    return;
+  }
+
+  pulseInterval = 60.0 / (2 * fermiChopperSpeed) * n_pulses;
+  runDetails.addProperty<double>("pulse_interval", pulseInterval);
+}
+
+/**
+ * Loads all the spectra into the workspace, including that from the monitor
+ *
+ * @param entry The Nexus entry
+ * @param monitors List of monitor data
+ */
+void LoadILLTOF2::loadDataIntoTheWorkSpace(
+    NeXus::NXEntry &entry, const std::vector<std::vector<int>> &monitors) {
+
+  g_log.debug() << "Loading data into the workspace...\n";
+  // read in the data
+  NXData dataGroup = entry.openNXData("data");
+  NXInt data = dataGroup.openIntData();
+  // load the counts from the file into memory
+  data.load();
+
+  // Put tof in an array
+  auto &X0 = m_localWorkspace->mutableX(0);
+  for (size_t i = 0; i < m_numberOfChannels + 1; ++i) {
+    X0[i] = m_timeOfFlightDelay +
+            m_channelWidth *
+                static_cast<double>(i) -
+            m_channelWidth /
+                2; // to make sure the bin centre is correct
+  }
+
+  // The binning for monitors is considered the same as for detectors
+  size_t spec = 0;
+
+  auto const &instrument = m_localWorkspace->getInstrument();
+
+  std::vector<detid_t> monitorIDs = instrument->getMonitors();
+
+  for (const auto &monitor : monitors) {
+    m_localWorkspace->setHistogram(spec, m_localWorkspace->binEdges(0),
+                                   Counts(monitor.begin(), monitor.end()));
+    m_localWorkspace->getSpectrum(spec).setDetectorID(monitorIDs[spec]);
+    spec++;
+  }
+
+  std::vector<detid_t> detectorIDs = instrument->getDetectorIDs(true);
+  size_t numberOfMonitors = monitors.size();
+
+  Progress progress(this, 0, 1, m_numberOfTubes * m_numberOfPixelsPerTube);
+
+  loadSpectra(spec, numberOfMonitors, m_numberOfTubes, detectorIDs, data,
+              progress);
+
+  g_log.debug() << "Loading data into the workspace: DONE!\n";
+
+  /**
+   * IN4 Rosace detectors are in a different NeXus entry
+   */
+  if (m_instrumentName == "IN4") {
+    g_log.debug() << "Loading data into the workspace: IN4 Rosace!\n";
+    // read in the data
+    NXData dataGroupRosace =
+        entry.openNXData("instrument/Detector_Rosace/data");
+    NXInt dataRosace = dataGroupRosace.openIntData();
+    auto numberOfTubes = static_cast<size_t>(dataRosace.dim0());
+    // load the counts from the file into memory
+    dataRosace.load();
+
+    Progress progressRosace(this, 0, 1,
+                            numberOfTubes * m_numberOfPixelsPerTube);
+
+    loadSpectra(spec, numberOfMonitors, numberOfTubes, detectorIDs, dataRosace,
+                progressRosace);
+  }
+}
+
+/**
+ * Loops over all the pixels and loads the correct spectra. Called for each set
+ * of detector types in the workspace.
+ *
+ * @param spec The current spectrum id
+ * @param numberOfMonitors The number of monitors in the workspace
+ * @param numberOfTubes The number of detector tubes in the workspace
+ * @param detectorIDs A list of all of the detector IDs
+ * @param data The NeXus data to load into the workspace
+ * @param progress The progress monitor
+ */
+void LoadILLTOF2::loadSpectra(size_t &spec, size_t numberOfMonitors,
+                             size_t numberOfTubes,
+                             std::vector<detid_t> &detectorIDs, NXInt data,
+                             Progress progress) {
+  for (size_t i = 0; i < numberOfTubes; ++i) {
+    for (size_t j = 0; j < m_numberOfPixelsPerTube; ++j) {
+      int *data_p = &data(static_cast<int>(i), static_cast<int>(j), 0);
+      m_localWorkspace->setHistogram(
+          spec, m_localWorkspace->binEdges(0),
+          Counts(data_p, data_p + m_numberOfChannels));
+      m_localWorkspace->getSpectrum(spec)
+          .setDetectorID(detectorIDs[spec - numberOfMonitors]);
+      spec++;
+      progress.report();
+    }
+  }
+}
+
+/**
+ * Runs LoadInstrument to attach the instrument to the workspace
+ */
+void LoadILLTOF2::runLoadInstrument() {
+
+  IAlgorithm_sptr loadInst = createChildAlgorithm("LoadInstrument");
+
+  // Execute the child algorithm. Catch and log any error, but don't stop.
+  try {
+    loadInst->setPropertyValue("InstrumentName", m_instrumentName);
+    loadInst->setProperty<MatrixWorkspace_sptr>("Workspace", m_localWorkspace);
+    loadInst->setProperty("RewriteSpectraMap",
+                          Mantid::Kernel::OptionalBool(false));
+    loadInst->execute();
+  } catch (...) {
+    g_log.information("Cannot load the instrument definition.");
+  }
+}