ExperimentInfo.cpp

#include "MantidAPI/ExperimentInfo.h"
#include "MantidAPI/ChopperModel.h"
#include "MantidAPI/ComponentInfo.h"
#include "MantidAPI/DetectorInfo.h"
#include "MantidAPI/InfoComponentVisitor.h"
#include "MantidAPI/InstrumentDataService.h"
#include "MantidAPI/ModeratorModel.h"
#include "MantidAPI/ResizeRectangularDetectorHelper.h"
#include "MantidAPI/Run.h"
#include "MantidAPI/Sample.h"
#include "MantidAPI/SpectrumInfo.h"

#include "MantidGeometry/IComponent.h"
#include "MantidGeometry/ICompAssembly.h"
#include "MantidGeometry/IDetector.h"
#include "MantidGeometry/Instrument/InstrumentDefinitionParser.h"
#include "MantidGeometry/Crystal/OrientedLattice.h"
#include "MantidGeometry/Instrument/ComponentVisitor.h"
#include "MantidGeometry/Instrument/Detector.h"
#include "MantidGeometry/Instrument/ParameterFactory.h"
#include "MantidGeometry/Instrument/ParameterMap.h"
#include "MantidGeometry/Instrument/ParComponentFactory.h"
#include "MantidGeometry/Instrument/RectangularDetector.h"
#include "MantidGeometry/Instrument/XMLInstrumentParameter.h"

#include "MantidBeamline/ComponentInfo.h"
#include "MantidBeamline/DetectorInfo.h"
#include "MantidBeamline/SpectrumInfo.h"

#include "MantidKernel/ConfigService.h"
#include "MantidKernel/DateAndTime.h"
#include "MantidKernel/EigenConversionHelpers.h"
#include "MantidKernel/InstrumentInfo.h"
#include "MantidKernel/Property.h"
#include "MantidKernel/Strings.h"
#include "MantidKernel/StringTokenizer.h"
#include "MantidKernel/make_unique.h"

#include "MantidTypes/SpectrumDefinition.h"

#include <boost/algorithm/string.hpp>
#include <boost/lexical_cast.hpp>
#include <boost/make_shared.hpp>
#include <boost/regex.hpp>

#include <Poco/DirectoryIterator.h>
#include <Poco/Path.h>
#include <Poco/SAX/Attributes.h>
#include <Poco/SAX/ContentHandler.h>
#include <Poco/SAX/SAXParser.h>
#include <nexus/NeXusException.hpp>

using namespace Mantid::Geometry;
using namespace Mantid::Kernel;
using namespace Poco::XML;

namespace Mantid {

namespace API {
namespace {
/// static logger object
Kernel::Logger g_log("ExperimentInfo");
}

/** Constructor
 */
ExperimentInfo::ExperimentInfo()
    : m_moderatorModel(), m_choppers(), m_sample(new Sample()),
      m_run(new Run()), m_parmap(new ParameterMap()),
      sptr_instrument(new Instrument()),
      m_detectorInfo(boost::make_shared<Beamline::DetectorInfo>()),
      m_componentInfo(boost::make_shared<Beamline::ComponentInfo>()) {
  m_parmap->setDetectorInfo(m_detectorInfo);
  m_parmap->setComponentInfo(m_componentInfo, std::vector<ComponentID>{});
  m_detectorInfoWrapper = Kernel::make_unique<DetectorInfo>(
      *m_detectorInfo, getInstrument(), m_parmap.get());
}

/**
 * Constructs the object from a copy if the input. This leaves the new mutex
 * unlocked.
 * @param source The source object from which to initialize
 */
ExperimentInfo::ExperimentInfo(const ExperimentInfo &source) {
  this->copyExperimentInfoFrom(&source);
  setSpectrumDefinitions(source.spectrumInfo().sharedSpectrumDefinitions());
}

// Defined as default in source for forward declaration with std::unique_ptr.
ExperimentInfo::~ExperimentInfo() = default;

/** Copy the experiment info data from another ExperimentInfo instance,
 * e.g. a MatrixWorkspace.
 * @param other :: the source from which to copy ExperimentInfo
 */
void ExperimentInfo::copyExperimentInfoFrom(const ExperimentInfo *other) {
  m_sample = other->m_sample;
  m_run = other->m_run->clone();
  this->setInstrument(other->getInstrument());
  if (other->m_moderatorModel)
    m_moderatorModel = other->m_moderatorModel->clone();
  m_choppers.clear();
  for (const auto &chopper : other->m_choppers) {
    m_choppers.push_back(chopper->clone());
  }
  // We do not copy Beamline::SpectrumInfo (which contains detector grouping
  // information) for now:
  // - For MatrixWorkspace, grouping information is still stored in ISpectrum
  //   and should not be overridden (copy is done in ExperimentInfo ctor, but
  //   not here since we just copy the experiment data).
  // - For cached groupings (for MDWorkspaces), grouping was not copied in the
  //   old implementation either.
}

/** Clone this ExperimentInfo class into a new one
 */
ExperimentInfo *ExperimentInfo::cloneExperimentInfo() const {
  return new ExperimentInfo(*this);
}

/// @returns A human-readable description of the object
const std::string ExperimentInfo::toString() const {
  try {
    populateIfNotLoaded();
  } catch (std::exception &) {
    // Catch any errors so that the string returned has as much information
    // as possible
  }

  std::ostringstream out;

  Geometry::Instrument_const_sptr inst = this->getInstrument();
  const auto instName = inst->getName();
  out << "Instrument: ";
  if (!instName.empty()) {
    out << instName << " ("
        << inst->getValidFromDate().toFormattedString("%Y-%b-%d") << " to "
        << inst->getValidToDate().toFormattedString("%Y-%b-%d") << ")";
    const auto instFilename = inst->getFilename();
    if (!instFilename.empty()) {
      out << "Instrument from: " << instFilename;
      out << "\n";
    }
  } else {
    out << "None";
  }
  out << "\n";

  // parameter files loaded
  auto paramFileVector =
      this->constInstrumentParameters().getParameterFilenames();
  for (auto &itFilename : paramFileVector) {
    out << "Parameters from: " << itFilename;
    out << "\n";
  }

  std::string runStart = getAvailableWorkspaceStartDate();
  std::string runEnd = getAvailableWorkspaceEndDate();
  std::string msgNA = "not available";
  if (runStart.empty())
    runStart = msgNA;
  if (runEnd.empty())
    runEnd = msgNA;
  out << "Run start: " << runStart << "\n";
  out << "Run end:  " << runEnd
      << "\n"; // note extra space for pseudo/approx-alignment

  if (this->sample().hasOrientedLattice()) {
    const Geometry::OrientedLattice &latt = this->sample().getOrientedLattice();
    out << "Sample: a " << std::fixed << std::setprecision(1) << latt.a()
        << ", b " << latt.b() << ", c " << latt.c();
    out << "; alpha " << std::fixed << std::setprecision(0) << latt.alpha()
        << ", beta " << latt.beta() << ", gamma " << latt.gamma();
    out << "\n";
  }
  return out.str();
}

// Helpers for setInstrument and getInstrument
namespace {
void checkDetectorInfoSize(const Instrument &instr,
                           const Beamline::DetectorInfo &detInfo) {
  const auto numDets = instr.getNumberDetectors();
  if (numDets != detInfo.size())
    throw std::runtime_error("ExperimentInfo: size mismatch between "
                             "DetectorInfo and number of detectors in "
                             "instrument");
}

boost::shared_ptr<Beamline::ComponentInfo>
makeComponentInfo(const Instrument &instrument,
                  const API::DetectorInfo &detectorInfo,
                  std::vector<Geometry::ComponentID> &componentIds) {

  if (instrument.hasComponentInfo()) {
    const auto &componentInfo = instrument.componentInfo();
    componentIds = instrument.componentIds();
    return boost::make_shared<Beamline::ComponentInfo>(componentInfo);
  } else {
    InfoComponentVisitor visitor(
        detectorInfo.size(), std::bind(&DetectorInfo::indexOf, &detectorInfo,
                                       std::placeholders::_1));

    if (instrument.isParametrized()) {
      // Register everything via visitor
      instrument.baseInstrument()->registerContents(visitor);
    } else {
      instrument.registerContents(visitor);
    }

    // Extract component ids. We need this for the ComponentInfo wrapper.
    componentIds = visitor.componentIds();

    return boost::make_shared<Mantid::Beamline::ComponentInfo>(
        visitor.assemblySortedDetectorIndices(),
        visitor.componentDetectorRanges());
  }
}

void clearPositionAndRotationsParameters(ParameterMap &pmap,
                                         const IDetector &det) {
  pmap.clearParametersByName(ParameterMap::pos(), &det);
  pmap.clearParametersByName(ParameterMap::posx(), &det);
  pmap.clearParametersByName(ParameterMap::posy(), &det);
  pmap.clearParametersByName(ParameterMap::posz(), &det);
  pmap.clearParametersByName(ParameterMap::rot(), &det);
  pmap.clearParametersByName(ParameterMap::rotx(), &det);
  pmap.clearParametersByName(ParameterMap::roty(), &det);
  pmap.clearParametersByName(ParameterMap::rotz(), &det);
}

std::unique_ptr<Beamline::DetectorInfo>
makeDetectorInfo(const Instrument &oldInstr, const Instrument &newInstr) {
  if (newInstr.hasDetectorInfo()) {
    // We allocate a new DetectorInfo in case there is an Instrument holding a
    // reference to our current DetectorInfo.
    const auto &detInfo = newInstr.detectorInfo();
    checkDetectorInfoSize(oldInstr, detInfo);
    return Kernel::make_unique<Beamline::DetectorInfo>(detInfo);
  } else {
    // If there is no DetectorInfo in the instrument we create a default one.
    const auto numDets = oldInstr.getNumberDetectors();
    const auto pmap = oldInstr.getParameterMap();
    // Currently monitors flags are stored in the detector cache of the base
    // instrument. The copy being made here is strictly speaking duplicating
    // that data, but with future refactoring this will no longer be the case.
    // Note that monitors will not change after creating a workspace.
    // Instrument::markAsMonitor works only for the base instrument and it is
    // not possible to obtain a non-const reference to the base instrument in a
    // workspace. Thus we do not need to worry about the two copies of monitor
    // flags running out of sync.
    std::vector<size_t> monitors;
    for (size_t i = 0; i < numDets; ++i)
      if (newInstr.isMonitorViaIndex(i))
        monitors.push_back(i);
    std::vector<Eigen::Vector3d> positions;
    std::vector<Eigen::Quaterniond> rotations;
    const auto &detIDs = newInstr.getDetectorIDs();
    for (const auto &id : detIDs) {
      const auto &det = newInstr.getDetector(id);
      // In the case of RectangularDetectorPixel the position is also affected
      // by the parameters scalex and scaly, but `getPos()` takes that into
      // account (if no DetectorInfo is set in the ParameterMap).
      positions.emplace_back(toVector3d(det->getPos()));
      rotations.emplace_back(toQuaterniond(det->getRotation()));
      clearPositionAndRotationsParameters(*pmap, *det);
      // Note that scalex and scaley also affect positions when set for a
      // RectangularDetector, but those are not parameters of the detector
      // itself so they are not cleared.
    }

    return Kernel::make_unique<Beamline::DetectorInfo>(
        std::move(positions), std::move(rotations), monitors);
  }
}
}
/** Set the instrument
* @param instr :: Shared pointer to an instrument.
*/
void ExperimentInfo::setInstrument(const Instrument_const_sptr &instr) {
  m_spectrumInfoWrapper = nullptr;
  if (instr->isParametrized()) {
    sptr_instrument = instr->baseInstrument();
    // We take a *copy* of the ParameterMap since we are modifying it by setting
    // a pointer to our DetectorInfo, and in case it contains legacy parameters
    // such as positions or rotations.
    m_parmap = boost::make_shared<ParameterMap>(*instr->getParameterMap());
  } else {
    sptr_instrument = instr;
    m_parmap = boost::make_shared<ParameterMap>();
  }
  const auto parInstrument = Geometry::ParComponentFactory::createInstrument(
      sptr_instrument, m_parmap);
  m_detectorInfo = makeDetectorInfo(*parInstrument, *instr);
  m_parmap->setDetectorInfo(m_detectorInfo);
  m_detectorInfoWrapper = Kernel::make_unique<DetectorInfo>(
      *m_detectorInfo, getInstrument(), m_parmap.get());

  m_componentInfo = makeComponentInfo(*instr, detectorInfo(), m_componentIds);
  m_parmap->setComponentInfo(m_componentInfo, m_componentIds);
  m_componentInfoWrapper =
      Kernel::make_unique<ComponentInfo>(*m_componentInfo, m_componentIds);
  // Detector IDs that were previously dropped because they were not part of the
  // instrument may now suddenly be valid, so we have to reinitialize the
  // detector grouping. Also the index corresponding to specific IDs may have
  // changed.
  invalidateAllSpectrumDefinitions();
}

/** Get a shared pointer to the parametrized instrument associated with this
*workspace
*
*  @return The instrument class
*/
Instrument_const_sptr ExperimentInfo::getInstrument() const {
  populateIfNotLoaded();
  checkDetectorInfoSize(*sptr_instrument, *m_detectorInfo);
  auto instrument = Geometry::ParComponentFactory::createInstrument(
      sptr_instrument, m_parmap);
  instrument->setDetectorInfo(m_detectorInfo);
  instrument->setComponentInfo(m_componentInfo, m_componentIds);
  return instrument;
}

/**  Returns a new copy of the instrument parameters
*    @return a (new) copy of the instruments parameter map
*/
Geometry::ParameterMap &ExperimentInfo::instrumentParameters() {
  populateIfNotLoaded();
  return *m_parmap;
}

/**  Returns a const reference to the instrument parameters.
*    @return a const reference to the instrument ParameterMap.
*/
const Geometry::ParameterMap &ExperimentInfo::instrumentParameters() const {
  populateIfNotLoaded();
  return *m_parmap;
}

/**  Returns a const reference to the instrument parameters.
*    @return a const reference to the instrument ParameterMap.
*/
const Geometry::ParameterMap &
ExperimentInfo::constInstrumentParameters() const {
  populateIfNotLoaded();
  return *m_parmap;
}

namespace {
///@cond

/// Used for storing info about "r-position", "t-position" and "p-position"
/// parameters
/// These are translated to X,Y,Z and so must all be processed together
struct RTP {
  RTP() : radius(0.0), haveRadius(false), theta(0.0), phi(0.0) {}
  double radius;
  bool haveRadius;
  double theta;
  double phi;
};

struct ParameterValue {
  ParameterValue(const Geometry::XMLInstrumentParameter &paramInfo,
                 const API::Run &run)
      : info(paramInfo), runData(run) {}

  operator double() {
    if (info.m_logfileID.empty())
      return boost::lexical_cast<double>(info.m_value);
    else
      return info.createParamValue(
          runData.getTimeSeriesProperty<double>(info.m_logfileID));
  }
  operator int() { return boost::lexical_cast<int>(info.m_value); }
  operator bool() {
    if (boost::iequals(info.m_value, "true"))
      return true;
    else if (boost::iequals(info.m_value, "yes"))
      return true;
    else
      return false;
  }
  const Geometry::XMLInstrumentParameter &info;
  const Run &runData;
};
///@endcond
}

namespace {
bool isPositionParameter(const std::string &name) {
  return ParameterMap::pos() == name;
}

bool isRotationParameter(const std::string &name) {
  return ParameterMap::rot() == name;
}

bool isScaleParameter(const std::string &name) {
  return (name == "scalex" || name == "scaley");
}

bool isRedundantPosOrRot(const std::string &name) {
  // Check size first as a small optimization.
  return (name.size() == 4) &&
         (name == "posx" || name == "posy" || name == "posz" ||
          name == "rotx" || name == "roty" || name == "rotz");
}

template <class T>
T getParam(const std::string &paramType, const std::string &paramValue) {
  const std::string name = "dummy";
  auto param = ParameterFactory::create(paramType, name);
  param->fromString(paramValue);
  return param->value<T>();
}

void updatePosition(DetectorInfo &detectorInfo, const Instrument &instrument,
                    const IComponent *component, const V3D &newRelPos) {
  // Important: Get component WITH ParameterMap so we see correct parent
  // positions and rotations and DetectorInfo updates correctly.
  if (!instrument.isParametrized())
    throw std::runtime_error(
        "Need parametrized instrument for updating positions");
  const auto &parComp = instrument.getComponentByID(component);
  auto position = newRelPos;
  if (auto parent = parComp->getParent()) {
    parent->getRotation().rotate(position);
    position += parent->getPos();
  }
  detectorInfo.setPosition(*parComp, position);
}

void updateRotation(DetectorInfo &detectorInfo, const Instrument &instrument,
                    const IComponent *component, const Quat &newRelRot) {
  // Important: Get component WITH ParameterMap so we see correct parent
  // positions and rotations and DetectorInfo updates correctly.
  if (!instrument.isParametrized())
    throw std::runtime_error(
        "Need parametrized instrument for updating rotations");
  const auto &parComp = instrument.getComponentByID(component);
  auto rotation = newRelRot;
  if (auto parent = parComp->getParent()) {
    // Note the unusual order. This is what Component::getRotation does.
    rotation = parent->getRotation() * newRelRot;
  }
  detectorInfo.setRotation(*parComp, rotation);
}

void adjustPositionsFromScaleFactor(DetectorInfo &detectorInfo,
                                    const Instrument &instrument,
                                    const IComponent *component,
                                    const std::string &paramName,
                                    double factor) {
  // Important: Get component WITH ParameterMap so we see correct parent
  // positions and rotations and DetectorInfo updates correctly.
  if (!instrument.isParametrized())
    throw std::runtime_error("Need parametrized instrument for updating "
                             "positions from scale factors");
  const auto &parComp = instrument.getComponentByID(component);
  const auto &det = dynamic_cast<const RectangularDetector &>(*parComp);
  double ScaleX = 1.0;
  double ScaleY = 1.0;
  if (paramName == "scalex")
    ScaleX = factor;
  else
    ScaleY = factor;
  applyRectangularDetectorScaleToDetectorInfo(detectorInfo, det, ScaleX,
                                              ScaleY);
}
}

/** Add parameters to the instrument parameter map that are defined in
* instrument
*   definition file or parameter file, which may contain parameters that require
*   logfile data to be available. Logs must be loaded before running this
* method.
*/
void ExperimentInfo::populateInstrumentParameters() {
  populateIfNotLoaded();

  // Reference to the run
  const auto &runData = run();

  // Get pointer to parameter map that we may add parameters to and information
  // about
  // the parameters that my be specified in the instrument definition file (IDF)
  Geometry::ParameterMap &paramMap = instrumentParameters();
  Geometry::ParameterMap paramMapForPosAndRot;

  // Get instrument and sample
  auto &detectorInfo = mutableDetectorInfo();
  const auto parInstrument = getInstrument();
  const auto instrument = parInstrument->baseInstrument();
  const auto &paramInfoFromIDF = instrument->getLogfileCache();

  std::map<const IComponent *, RTP> rtpParams;

  // In this loop position and rotation parameters are inserted into the
  // temporary map paramMapForPosAndRot. In the subsequent loop, after all
  // parameters have been parsed, we update positions and rotations in
  // DetectorInfo and the temporary map goes out of scope. The main reason for
  // this is that ParameterMap will then take care of assembling parameters for
  // individual position or rotation components into a vector or quaternion. In
  // particular, we cannot directly change DetectorInfo since the order of
  // rotation components is not guaranteed.
  for (const auto &item : paramInfoFromIDF) {
    const auto &nameComp = item.first;
    const auto &paramInfo = item.second;
    const std::string &paramN = nameComp.first;

    try {
      // Special case where user has specified r-position,t-position, and/or
      // p-position.
      // We need to know all three first to calculate a set of X,Y,Z
      if (paramN.compare(1, 9, "-position") == 0) {
        auto &rtpValues = rtpParams[paramInfo->m_component]; // If not found,
                                                             // constructs
                                                             // default
        double value = ParameterValue(*paramInfo, runData);
        if (paramN.compare(0, 1, "r") == 0) {
          rtpValues.radius = value;
          rtpValues.haveRadius = true;
        } else if (paramN.compare(0, 1, "t") == 0)
          rtpValues.theta = value;
        else if (paramN.compare(0, 1, "p") == 0)
          rtpValues.phi = value;
        if (rtpValues.haveRadius) {
          V3D pos;
          pos.spherical(rtpValues.radius, rtpValues.theta, rtpValues.phi);
          paramMapForPosAndRot.addV3D(paramInfo->m_component,
                                      ParameterMap::pos(), pos);
        }
      } else {
        populateWithParameter(paramMap, paramMapForPosAndRot, paramN,
                              *paramInfo, runData);
      }
    } catch (std::exception &exc) {
      g_log.information() << "Unable to add component parameter '"
                          << nameComp.first << "'. Error: " << exc.what();
      continue;
    }
  }
  for (const auto &item : paramMapForPosAndRot) {
    if (isPositionParameter(item.second->name())) {
      const auto newRelPos = item.second->value<V3D>();
      updatePosition(detectorInfo, *parInstrument, item.first, newRelPos);
    } else if (isRotationParameter(item.second->name())) {
      const auto newRelRot = item.second->value<Quat>();
      updateRotation(detectorInfo, *parInstrument, item.first, newRelRot);
    }
    // Parameters for individual components (x,y,z) are ignored. ParameterMap
    // did compute the correct compound positions and rotations internally.
  }
  // Special case RectangularDetector: Parameters scalex and scaley affect pixel
  // positions.
  for (const auto &item : paramMap) {
    if (isScaleParameter(item.second->name()))
      adjustPositionsFromScaleFactor(detectorInfo, *parInstrument, item.first,
                                     item.second->name(),
                                     item.second->value<double>());
  }
  // paramMapForPosAndRot goes out of scope, dropping all position and rotation
  // parameters of detectors (parameters for non-detector components have been
  // inserted into paramMap via DetectorInfo::setPosition(IComponent *)).
}

/**
 * Caches a lookup for the detector IDs of the members that are part of the same
 * group
 * @param mapping :: A map between a detector ID and the other IDs that are part
 * of the same
 * group.
 */
void ExperimentInfo::cacheDetectorGroupings(const det2group_map &mapping) {
  populateIfNotLoaded();
  if (mapping.empty()) {
    cacheDefaultDetectorGrouping();
    return;
  }
  setNumberOfDetectorGroups(mapping.size());
  size_t specIndex = 0;
  for (const auto &item : mapping) {
    m_det2group[item.first] = specIndex;
    setDetectorGrouping(specIndex, item.second);
    specIndex++;
  }
}

/** Sets the number of detector groups.
 *
 * This method should not need to be called explicitly. The number of detector
 * groups will be set either when initializing a MatrixWorkspace, or by calling
 * `cacheDetectorGroupings` for an ExperimentInfo stored in an MDWorkspace. */
void ExperimentInfo::setNumberOfDetectorGroups(const size_t count) const {
  populateIfNotLoaded();
  if (m_spectrumInfo)
    m_spectrumDefinitionNeedsUpdate.clear();
  m_spectrumDefinitionNeedsUpdate.resize(count, 1);
  m_spectrumInfo = Kernel::make_unique<Beamline::SpectrumInfo>(count);
  m_spectrumInfoWrapper = nullptr;
}

/** Sets the detector grouping for the spectrum with the given `index`.
 *
 * This method should not need to be called explicitly. Groupings are updated
 * automatically when modifying detector IDs in a workspace (via ISpectrum). */
void ExperimentInfo::setDetectorGrouping(
    const size_t index, const std::set<detid_t> &detIDs) const {
  SpectrumDefinition specDef;
  // Wrap translation in check for detector count as an optimization of
  // otherwise slow failures via exceptions.
  if (detectorInfo().size() > 0) {
    for (const auto detID : detIDs) {
      try {
        const size_t detIndex = detectorInfo().indexOf(detID);
        specDef.add(detIndex);
      } catch (std::out_of_range &) {
        // Silently strip bad detector IDs
      }
    }
  }
  m_spectrumInfo->setSpectrumDefinition(index, std::move(specDef));
  m_spectrumDefinitionNeedsUpdate.at(index) = 0;
}

/** Update detector grouping for spectrum with given index.
 *
 * This method is called when the detector grouping stored in SpectrumDefinition
 * at `index` in Beamline::SpectrumInfo is not initialized or outdated. The
 * implementation throws, since no grouping information for update is available
 * when grouping comes from a call to `cacheDetectorGroupings`. This method is
 * overridden in MatrixWorkspace. */
void ExperimentInfo::updateCachedDetectorGrouping(const size_t) const {
  throw std::runtime_error("ExperimentInfo::updateCachedDetectorGrouping: "
                           "Cannot update -- grouping information not "
                           "available");
}

/**
 * Set an object describing the moderator properties and take ownership
 * @param source :: A pointer to an object describing the source. Ownership is
 * transferred to this object
 */
void ExperimentInfo::setModeratorModel(ModeratorModel *source) {
  populateIfNotLoaded();
  if (!source) {
    throw std::invalid_argument(
        "ExperimentInfo::setModeratorModel - NULL source object found.");
  }
  m_moderatorModel = boost::shared_ptr<ModeratorModel>(source);
}

/// Returns a reference to the source properties object
ModeratorModel &ExperimentInfo::moderatorModel() const {
  populateIfNotLoaded();
  if (!m_moderatorModel) {
    throw std::runtime_error("ExperimentInfo::moderatorModel - No source "
                             "desciption has been defined");
  }
  return *m_moderatorModel;
}

/**
 * Sets a new chopper description at a given point. The point is given by index
 * where 0 is
 * closest to the source
 * @param chopper :: A pointer to a new chopper object, this class takes
 * ownership of the pointer
 * @param index :: An optional index that specifies which chopper point the
 * chopper belongs to (default=0)
 */
void ExperimentInfo::setChopperModel(ChopperModel *chopper,
                                     const size_t index) {
  populateIfNotLoaded();
  if (!chopper) {
    throw std::invalid_argument(
        "ExperimentInfo::setChopper - NULL chopper object found.");
  }
  auto iter = m_choppers.begin();
  std::advance(iter, index);
  if (index < m_choppers.size()) // Replacement
  {
    (*iter) = boost::shared_ptr<ChopperModel>(chopper);
  } else // Insert it
  {
    m_choppers.insert(iter, boost::shared_ptr<ChopperModel>(chopper));
  }
}

/**
 * Returns a const reference to a chopper description
 * @param index :: An optional index giving the point within the instrument this
 * chopper describes (default=0)
 * @return A reference to a const chopper object
 */
ChopperModel &ExperimentInfo::chopperModel(const size_t index) const {
  populateIfNotLoaded();
  if (index < m_choppers.size()) {
    auto iter = m_choppers.begin();
    std::advance(iter, index);
    return **iter;
  } else {
    std::ostringstream os;
    os << "ExperimentInfo::chopper - Invalid index=" << index << ". "
       << m_choppers.size() << " chopper descriptions have been set.";
    throw std::invalid_argument(os.str());
  }
}

/** Get a constant reference to the Sample associated with this workspace.
* @return const reference to Sample object
*/
const Sample &ExperimentInfo::sample() const {
  populateIfNotLoaded();
  std::lock_guard<std::recursive_mutex> lock(m_mutex);
  return *m_sample;
}

/** Get a reference to the Sample associated with this workspace.
*  This non-const method will copy the sample if it is shared between
*  more than one workspace, and the reference returned will be to the copy.
*  Can ONLY be taken by reference!
* @return reference to sample object
*/
Sample &ExperimentInfo::mutableSample() {
  populateIfNotLoaded();
  // Use a double-check for sharing so that we only
  // enter the critical region if absolutely necessary
  if (!m_sample.unique()) {
    std::lock_guard<std::recursive_mutex> lock(m_mutex);
    // Check again because another thread may have taken copy
    // and dropped reference count since previous check
    if (!m_sample.unique()) {
      boost::shared_ptr<Sample> oldData = m_sample;
      m_sample = boost::make_shared<Sample>(*oldData);
    }
  }
  return *m_sample;
}

/** Get a constant reference to the Run object associated with this workspace.
* @return const reference to run object
*/
const Run &ExperimentInfo::run() const {
  populateIfNotLoaded();
  std::lock_guard<std::recursive_mutex> lock(m_mutex);
  return *m_run;
}

/** Get a reference to the Run object associated with this workspace.
*  This non-const method will copy the Run object if it is shared between
*  more than one workspace, and the reference returned will be to the copy.
*  Can ONLY be taken by reference!
* @return reference to Run object
*/
Run &ExperimentInfo::mutableRun() {
  populateIfNotLoaded();
  // Use a double-check for sharing so that we only
  // enter the critical region if absolutely necessary
  if (!m_run.unique()) {
    std::lock_guard<std::recursive_mutex> lock(m_mutex);
    // Check again because another thread may have taken copy
    // and dropped reference count since previous check
    if (!m_run.unique()) {
      boost::shared_ptr<Run> oldData = m_run;
      m_run = boost::make_shared<Run>(*oldData);
    }
  }
  return *m_run;
}

/**
 * Get an experimental log either by log name or by type, e.g.
 *   - temperature_log
 *   - chopper_speed_log
 * The logs are first checked for one matching the given string and if that
 * fails then the instrument is checked for a parameter of the same name
 * and if this exists then its value is assume to be the actual log required
 * @param log :: A string giving either a specific log name or instrument
 * parameter whose
 * value is to be retrieved
 * @return A pointer to the property
 */
Kernel::Property *ExperimentInfo::getLog(const std::string &log) const {
  populateIfNotLoaded();
  try {
    return run().getProperty(log);
  } catch (Kernel::Exception::NotFoundError &) {
    // No log with that name
  }
  // If the instrument has a parameter with that name then take the value as a
  // log name
  const std::string logName =
      constInstrumentParameters().getString(sptr_instrument.get(), log);
  if (logName.empty()) {
    throw std::invalid_argument(
        "ExperimentInfo::getLog - No instrument parameter named \"" + log +
        "\". Cannot access full log name");
  }
  return run().getProperty(logName);
}

/**
 * Get an experimental log as a single value either by log name or by type. @see
 * getLog
 * @param log :: A string giving either a specific log name or instrument
 * parameter whose
 * value is to be retrieved
 * @return A pointer to the property
 */
double ExperimentInfo::getLogAsSingleValue(const std::string &log) const {
  populateIfNotLoaded();
  try {
    return run().getPropertyAsSingleValue(log);
  } catch (Kernel::Exception::NotFoundError &) {
    // No log with that name
  }
  // If the instrument has a parameter with that name then take the value as a
  // log name
  const std::string logName =
      constInstrumentParameters().getString(sptr_instrument.get(), log);
  if (logName.empty()) {
    throw std::invalid_argument(
        "ExperimentInfo::getLog - No instrument parameter named \"" + log +
        "\". Cannot access full log name");
  }
  return run().getPropertyAsSingleValue(logName);
}

/** Utility method to get the run number
 *
 * @return the run number (int) or 0 if not found.
 */
int ExperimentInfo::getRunNumber() const {
  populateIfNotLoaded();
  const Run &thisRun = run();
  if (!thisRun.hasProperty("run_number")) {
    // No run_number property, default to 0
    return 0;
  } else {
    Property *prop = m_run->getProperty("run_number");
    if (prop) {
      // Use the string representation. That way both a string and a number
      // property will work.
      int val;
      if (Strings::convert(prop->value(), val))
        return val;
      else
        return 0;
    }
  }
  return 0;
}

/**
 * Returns the emode for this run. It first searchs the run logs for a
 * "deltaE-mode" log and falls back to
 * the instrument if one is not found. If neither exist then the run is
 * considered Elastic.
 * @return The emode enum for the energy transfer mode of this run. Currently
 * checks the sample log & instrument in this order
 */
Kernel::DeltaEMode::Type ExperimentInfo::getEMode() const {
  populateIfNotLoaded();
  static const char *emodeTag = "deltaE-mode";
  std::string emodeStr;
  if (run().hasProperty(emodeTag)) {
    emodeStr = run().getPropertyValueAsType<std::string>(emodeTag);
  } else if (sptr_instrument &&
             constInstrumentParameters().contains(sptr_instrument.get(),
                                                  emodeTag)) {
    Geometry::Parameter_sptr param =
        constInstrumentParameters().get(sptr_instrument.get(), emodeTag);
    emodeStr = param->asString();
  } else {
    return Kernel::DeltaEMode::Elastic;
  }
  return Kernel::DeltaEMode::fromString(emodeStr);
}

/**
 * Easy access to the efixed value for this run & detector ID
 * @param detID :: The detector ID to ask for the efixed mode (ignored in Direct
 * & Elastic mode). The
 * detector with ID matching that given is pulled from the instrument with this
 * method and it will
 * throw a Exception::NotFoundError if the ID is unknown.
 * @return The current EFixed value
 */
double ExperimentInfo::getEFixed(const detid_t detID) const {
  populateIfNotLoaded();
  IDetector_const_sptr det = getInstrument()->getDetector(detID);
  return getEFixed(det);
}

/**
 * Easy access to the efixed value for this run & detector
 * @param detector :: The detector object to ask for the efixed mode. Only
 * required for Indirect mode
 * @return The current efixed value
 */
double
ExperimentInfo::getEFixed(const Geometry::IDetector_const_sptr detector) const {
  populateIfNotLoaded();
  Kernel::DeltaEMode::Type emode = getEMode();
  if (emode == Kernel::DeltaEMode::Direct) {
    try {
      return this->run().getPropertyValueAsType<double>("Ei");
    } catch (Kernel::Exception::NotFoundError &) {
      throw std::runtime_error(
          "Experiment logs do not contain an Ei value. Have you run GetEi?");
    }
  } else if (emode == Kernel::DeltaEMode::Indirect) {
    if (!detector)
      throw std::runtime_error("ExperimentInfo::getEFixed - Indirect mode "
                               "efixed requested without a valid detector.");
    Parameter_sptr par =
        constInstrumentParameters().getRecursive(detector.get(), "Efixed");
    if (par) {
      return par->value<double>();
    } else {
      std::vector<double> efixedVec = detector->getNumberParameter("Efixed");
      if (efixedVec.empty()) {
        int detid = detector->getID();
        IDetector_const_sptr detectorSingle =
            getInstrument()->getDetector(detid);
        efixedVec = detectorSingle->getNumberParameter("Efixed");
      }
      if (!efixedVec.empty()) {
        return efixedVec.at(0);
      } else {
        std::ostringstream os;
        os << "ExperimentInfo::getEFixed - Indirect mode efixed requested but "
              "detector has no Efixed parameter attached. ID="
           << detector->getID();
        throw std::runtime_error(os.str());
      }
    }
  } else {
    throw std::runtime_error("ExperimentInfo::getEFixed - EFixed requested for "
                             "elastic mode, don't know what to do!");
  }
}

void ExperimentInfo::setEFixed(const detid_t detID, const double value) {
  populateIfNotLoaded();
  IDetector_const_sptr det = getInstrument()->getDetector(detID);
  Geometry::ParameterMap &pmap = instrumentParameters();
  pmap.addDouble(det.get(), "Efixed", value);
}

// used to terminate SAX process
class DummyException {
public:
  std::string m_validFrom;
  std::string m_validTo;
  DummyException(const std::string &validFrom, const std::string &validTo)
      : m_validFrom(validFrom), m_validTo(validTo) {}
};

// SAX content handler for grapping stuff quickly from IDF
class myContentHandler : public Poco::XML::ContentHandler {
  void startElement(const XMLString &, const XMLString &localName,
                    const XMLString &, const Attributes &attrList) override {
    if (localName == "instrument") {
      throw DummyException(
          static_cast<std::string>(attrList.getValue("", "valid-from")),
          static_cast<std::string>(attrList.getValue("", "valid-to")));
    }
  }
  void endElement(const XMLString &, const XMLString &,
                  const XMLString &) override {}
  void startDocument() override {}
  void endDocument() override {}
  void characters(const XMLChar[], int, int) override {}
  void endPrefixMapping(const XMLString &) override {}
  void ignorableWhitespace(const XMLChar[], int, int) override {}
  void processingInstruction(const XMLString &, const XMLString &) override {}
  void setDocumentLocator(const Locator *) override {}
  void skippedEntity(const XMLString &) override {}
  void startPrefixMapping(const XMLString &, const XMLString &) override {}
};

/** Return from an IDF the values of the valid-from and valid-to attributes
*
*  @param IDFfilename :: Full path of an IDF
*  @param[out] outValidFrom :: Used to return valid-from date
*  @param[out] outValidTo :: Used to return valid-to date
*/
void ExperimentInfo::getValidFromTo(const std::string &IDFfilename,
                                    std::string &outValidFrom,
                                    std::string &outValidTo) {
  SAXParser pParser;
  // Create on stack to ensure deletion. Relies on pParser also being local
  // variable.
  myContentHandler conHand;
  pParser.setContentHandler(&conHand);