Skip to content
Snippets Groups Projects
SetSample.cpp 16.9 KiB
Newer Older
Martyn Gigg's avatar
Martyn Gigg committed
#include "MantidDataHandling/SetSample.h"

#include "MantidAPI/MatrixWorkspace.h"
#include "MantidAPI/Sample.h"
Martyn Gigg's avatar
Martyn Gigg committed
#include "MantidGeometry/Instrument.h"
#include "MantidGeometry/Instrument/Goniometer.h"
#include "MantidGeometry/Instrument/ReferenceFrame.h"
Martyn Gigg's avatar
Martyn Gigg committed
#include "MantidGeometry/Instrument/SampleEnvironmentFactory.h"
#include "MantidKernel/ConfigService.h"
#include "MantidKernel/FacilityInfo.h"
#include "MantidKernel/InstrumentInfo.h"
#include "MantidKernel/Logger.h"
#include "MantidKernel/Material.h"
#include "MantidKernel/Matrix.h"
#include "MantidKernel/PropertyManager.h"
Martyn Gigg's avatar
Martyn Gigg committed
#include "MantidKernel/PropertyManagerProperty.h"

#include <boost/algorithm/string/case_conv.hpp>
#include <boost/algorithm/string/predicate.hpp>
Martyn Gigg's avatar
Martyn Gigg committed
#include <Poco/Path.h>

namespace Mantid {
namespace DataHandling {

using Geometry::Goniometer;
using Geometry::ReferenceFrame;
using Geometry::SampleEnvironment;
using Kernel::Logger;
using Kernel::PropertyWithValue;
using Kernel::Quat;
/**
  * Return the centre coordinates of the base of a cylinder given the
  * coordinates of the centre of the cylinder
  * @param cylCentre Coordinates of centre of the cylinder (X,Y,Z) (in metres)
  * @param height Height of the cylinder (in metres)
  * @param axis The index of the height-axis of the cylinder
  */
V3D cylBaseCentre(const std::vector<double> &cylCentre, double height,
  const V3D halfHeight = [&]() {
    switch (axisIdx) {
    case 0:
      return V3D(0.5 * height, 0, 0);
    case 1:
      return V3D(0, 0.5 * height, 0);
    case 2:
      return V3D(0, 0, 0.5 * height);
    default:
      return V3D();
    }
  }();
  return V3D(cylCentre[0], cylCentre[1], cylCentre[2]) - halfHeight;
/**
 * Create the xml tag require for a given axis index
 * @param axisIdx Index 0,1,2 for the axis of a cylinder
 * @return A string containing the axis tag for this index
 */
std::string axisXML(unsigned axisIdx) {
  switch (axisIdx) {
  case 0:
    return "<axis x=\"1\" y=\"0\" z=\"0\" />";
  case 1:
    return "<axis x=\"0\" y=\"1\" z=\"0\" />";
  case 2:
    return "<axis x=\"0\" y=\"0\" z=\"1\" />";
  default:
    return "";
  }
}
}
Martyn Gigg's avatar
Martyn Gigg committed

// Register the algorithm into the AlgorithmFactory
DECLARE_ALGORITHM(SetSample)

/// Algorithms name for identification. @see Algorithm::name
const std::string SetSample::name() const { return "SetSample"; }

/// Algorithm's version for identification. @see Algorithm::version
int SetSample::version() const { return 1; }

/// Algorithm's category for identification. @see Algorithm::category
const std::string SetSample::category() const { return "Sample"; }

/// Algorithm's summary for use in the GUI and help. @see Algorithm::summary
const std::string SetSample::summary() const {
  return "Set properties of the sample and its environment for a workspace";
}

/// Validate the inputs against each other @see Algorithm::validateInputs
std::map<std::string, std::string> SetSample::validateInputs() {
  using Kernel::PropertyManager_sptr;
  std::map<std::string, std::string> errors;

  // Validate Environment
  PropertyManager_sptr environArgs = getProperty("Environment");
  if (environArgs) {
    if (!environArgs->existsProperty("Name")) {
      errors["Environment"] = "Environment flags must contain a 'Name' entry.";
    } else {
      std::string name = environArgs->getPropertyValue("Name");
      if (name.empty()) {
        errors["Environment"] = "Environment 'Name' flag is an empty string!";
      }
    }

Martyn Gigg's avatar
Martyn Gigg committed
    if (!environArgs->existsProperty("Container")) {
      errors["Environment"] =
          "Environment flags must contain a 'Container' entry.";
Martyn Gigg's avatar
Martyn Gigg committed
    } else {
Martyn Gigg's avatar
Martyn Gigg committed
      std::string name = environArgs->getPropertyValue("Container");
Martyn Gigg's avatar
Martyn Gigg committed
      if (name.empty()) {
        errors["Environment"] = "Environment 'Can' flag is an empty string!";
      }
    }
  }

  return errors;
}

/**
 * Initialize the algorithm's properties.
 */
void SetSample::init() {
  using API::WorkspaceProperty;
  using Kernel::Direction;
  using Kernel::PropertyManagerProperty;

  // Inputs
  declareProperty(Kernel::make_unique<WorkspaceProperty<>>("InputWorkspace", "",
                                                           Direction::InOut),
                  "A workspace whose sample properties will be updated");
  declareProperty(Kernel::make_unique<PropertyManagerProperty>(
                      "Geometry", Direction::Input),
                  "A dictionary of geometry parameters for the sample.");
  declareProperty(Kernel::make_unique<PropertyManagerProperty>(
                      "Material", Direction::Input),
                  "A dictionary of material parameters for the sample. See "
                  "SetSampleMaterial for all accepted parameters");
  declareProperty(
      Kernel::make_unique<PropertyManagerProperty>("Environment",
                                                   Direction::Input),
      "A dictionary of parameters to configure the sample environment");
}

/**
 * Execute the algorithm.
 */
void SetSample::exec() {
  using API::MatrixWorkspace_sptr;
  using Kernel::PropertyManager_sptr;

  MatrixWorkspace_sptr workspace = getProperty("InputWorkspace");
  PropertyManager_sptr environArgs = getProperty("Environment");
  PropertyManager_sptr geometryArgs = getProperty("Geometry");
  PropertyManager_sptr materialArgs = getProperty("Material");

  // The order here is important. Se the environment first. If this
  // defines a sample geometry then we can process the Geometry flags
  // combined with this
  const SampleEnvironment *sampleEnviron(nullptr);
  if (environArgs) {
    sampleEnviron = setSampleEnvironment(workspace, environArgs);
Martyn Gigg's avatar
Martyn Gigg committed
  }

  if (geometryArgs || sampleEnviron) {
    setSampleShape(workspace, geometryArgs, sampleEnviron);
  }

  // Finally the material arguments
  if (materialArgs) {
    runChildAlgorithm("SetSampleMaterial", workspace, *materialArgs);
  }
}

/**
 * Set the requested sample environment on the workspace
 * @param workspace A pointer to the workspace to be affected
 * @param args The dictionary of flags for the environment
 * @return A pointer to the new sample environment
 */
const Geometry::SampleEnvironment *
SetSample::setSampleEnvironment(API::MatrixWorkspace_sptr &workspace,
                                const Kernel::PropertyManager_sptr &args) {
Martyn Gigg's avatar
Martyn Gigg committed
  using Geometry::SampleEnvironmentSpecFileFinder;
  using Geometry::SampleEnvironmentFactory;
  using Kernel::ConfigService;

  const std::string envName = args->getPropertyValue("Name");
  const std::string canName = args->getPropertyValue("Container");
Martyn Gigg's avatar
Martyn Gigg committed
  // The specifications need to be qualified by the facility and instrument.
  // Check instrument for name and then lookup facility if facility
  // is unknown then set to default facility & instrument.
Martyn Gigg's avatar
Martyn Gigg committed
  auto instrument = workspace->getInstrument();
  const auto &instOnWS = instrument->getName();
Martyn Gigg's avatar
Martyn Gigg committed
  const auto &config = ConfigService::Instance();
  std::string facilityName, instrumentName;
  try {
    const auto &instInfo = config.getInstrument(instOnWS);
    instrumentName = instInfo.name();
    facilityName = instInfo.facility().name();
  } catch (std::runtime_error &) {
    // use default facility/instrument
    facilityName = config.getFacility().name();
    instrumentName = config.getInstrument().name();
  }
Martyn Gigg's avatar
Martyn Gigg committed

  const auto &instDirs = config.getInstrumentDirectories();
  std::vector<std::string> environDirs(instDirs);
  for (auto &direc : environDirs) {
    direc = Poco::Path(direc).append("sampleenvironments").toString();
  }
  auto finder =
      Kernel::make_unique<SampleEnvironmentSpecFileFinder>(environDirs);
  SampleEnvironmentFactory factory(std::move(finder));
  auto sampleEnviron =
      factory.create(facilityName, instrumentName, envName, canName);
Martyn Gigg's avatar
Martyn Gigg committed
  workspace->mutableSample().setEnvironment(sampleEnviron.release());

  return &(workspace->sample().getEnvironment());
}

/**
 * @param workspace A pointer to the workspace to be affected
 * @param args The user-supplied dictionary of flags
 * @param sampleEnv A pointer to the sample environment if one exists, otherwise
 * null
 * @return A string containing the XML definition of the shape
 */
void SetSample::setSampleShape(API::MatrixWorkspace_sptr &workspace,
                               Kernel::PropertyManager_sptr &args,
Martyn Gigg's avatar
Martyn Gigg committed
                               const Geometry::SampleEnvironment *sampleEnv) {
Martyn Gigg's avatar
Martyn Gigg committed
  using Geometry::Container;
Martyn Gigg's avatar
Martyn Gigg committed
  /* The sample geometry can be specified in two ways:
     - a known set of primitive shapes with values or CSG string
     - or a <samplegeometry> field sample environment can, with values possible
       overridden by the Geometry flags
  */

  // Try known shapes or CSG first if supplied
  if (args) {
    auto refFrame = workspace->getInstrument()->getReferenceFrame();
    const auto xml = tryCreateXMLFromArgsOnly(*args, *refFrame);
    if (!xml.empty()) {
      runSetSampleShape(workspace, xml);
      return;
    }
Martyn Gigg's avatar
Martyn Gigg committed
  }
  // Any arguments in the args dict are assumed to be values that should
  // override the default set by the sampleEnv samplegeometry if it exists
  if (sampleEnv) {
Martyn Gigg's avatar
Martyn Gigg committed
    if (sampleEnv->container()->hasSampleShape()) {
      const auto &can = sampleEnv->container();
      Container::ShapeArgs shapeArgs;
Martyn Gigg's avatar
Martyn Gigg committed
      if (args) {
        const auto &props = args->getProperties();
        for (const auto &prop : props) {
          // assume in cm
Martyn Gigg's avatar
Martyn Gigg committed
          const double val = args->getProperty(prop->name());
          shapeArgs.emplace(boost::algorithm::to_lower_copy(prop->name()),
                            val * 0.01);
Martyn Gigg's avatar
Martyn Gigg committed
        }
      }
      auto shapeObject = can->createSampleShape(shapeArgs);
      // Set the object directly on the sample ensuring we preserve the
      // material
      const auto mat = workspace->sample().getMaterial();
Martyn Gigg's avatar
Martyn Gigg committed
      shapeObject->setMaterial(mat);
      workspace->mutableSample().setShape(*shapeObject);
    } else {
      throw std::runtime_error("The can does not define the sample shape. "
                               "Please either provide a 'Shape' argument "
                               "or update the environment definition with "
                               "this information.");
    }
  } else {
    throw std::runtime_error("No sample environment defined, please provide "
                             "a 'Shape' argument to define the sample "
                             "shape.");
  }
}

/**
 * Create the required XML for a given shape type plus its arguments
 * @param args A dict of flags defining the shape
 * @param refFrame Defines the reference frame for the shape
Martyn Gigg's avatar
Martyn Gigg committed
 * @return A string containing the XML if possible or an empty string
 */
std::string
SetSample::tryCreateXMLFromArgsOnly(const Kernel::PropertyManager &args,
                                    const Geometry::ReferenceFrame &refFrame) {
Martyn Gigg's avatar
Martyn Gigg committed
  std::string result;
  if (!args.existsProperty("Shape")) {
Martyn Gigg's avatar
Martyn Gigg committed
    return result;
  }

  const auto shape = args.getPropertyValue("Shape");
Martyn Gigg's avatar
Martyn Gigg committed
  if (shape == "CSG") {
    result = args.getPropertyValue("Value");
Martyn Gigg's avatar
Martyn Gigg committed
  } else if (shape == "FlatPlate") {
    result = createFlatPlateXML(args, refFrame);
  } else if (boost::algorithm::ends_with(shape, "Cylinder")) {
    result = createCylinderLikeXML(
        args, refFrame, boost::algorithm::starts_with(shape, "Hollow"));
Martyn Gigg's avatar
Martyn Gigg committed
  } else {
    throw std::invalid_argument(
        "Unknown 'Shape' argument provided in "
        "'Geometry'. Allowed "
        "values=FlatPlate,CSG,Cylinder,HollowCylinder.");
Martyn Gigg's avatar
Martyn Gigg committed
  }
  if (g_log.is(Logger::Priority::PRIO_DEBUG)) {
    g_log.debug("XML shape definition:\n" + result + '\n');
  }
Martyn Gigg's avatar
Martyn Gigg committed
  return result;
}

/**
 * Create the XML required to define a flat plate from the given args
 * @param args A user-supplied dict of args
 * @param refFrame Defines the reference frame for the shape
Martyn Gigg's avatar
Martyn Gigg committed
 * @return The XML definition string
 */
std::string
SetSample::createFlatPlateXML(const Kernel::PropertyManager &args,
                              const Geometry::ReferenceFrame &refFrame) const {
  // Helper to take 3 coordinates and turn them to a V3D respecting the
  // current reference frame
  auto makeV3D = [&refFrame](double x, double y, double z) {
    V3D v;
    v[refFrame.pointingHorizontal()] = x;
    v[refFrame.pointingUp()] = y;
    v[refFrame.pointingAlongBeam()] = z;
    return v;
  };
  const double widthInCM = args.getProperty("Width");
  const double heightInCM = args.getProperty("Height");
  const double thickInCM = args.getProperty("Thick");
  // Convert to half-"width" in metres
  const double szX = (widthInCM * 5e-3);
  const double szY = (heightInCM * 5e-3);
  const double szZ = (thickInCM * 5e-3);
  // Contruct cuboid corners. Define points about origin, rotate and then
  // translate to final center position
  auto lfb = makeV3D(szX, -szY, -szZ);
  auto lft = makeV3D(szX, szY, -szZ);
  auto lbb = makeV3D(szX, -szY, szZ);
  auto rfb = makeV3D(-szX, -szY, -szZ);
  // optional rotation about the center of object
  if (args.existsProperty("Angle")) {
    Goniometer gr;
    const auto upAxis = makeV3D(0, 1, 0);
    gr.pushAxis("up", upAxis.X(), upAxis.Y(), upAxis.Z(),
                args.getProperty("Angle"), Geometry::CCW, Geometry::angDegrees);
    auto &rotation = gr.getR();
    lfb.rotate(rotation);
    lft.rotate(rotation);
    lbb.rotate(rotation);
    rfb.rotate(rotation);
  }
  std::vector<double> center = args.getProperty("Center");
  const V3D centrePos(center[0] * 0.01, center[1] * 0.01, center[2] * 0.01);
  // translate to true center after rotation
  lfb += centrePos;
  lft += centrePos;
  lbb += centrePos;
  rfb += centrePos;

  std::ostringstream xmlShapeStream;
  xmlShapeStream << " <cuboid id=\"sample-shape\"> "
                 << "<left-front-bottom-point x=\"" << lfb.X() << "\" y=\""
                 << lfb.Y() << "\" z=\"" << lfb.Z() << "\"  /> "
                 << "<left-front-top-point  x=\"" << lft.X() << "\" y=\""
                 << lft.Y() << "\" z=\"" << lft.Z() << "\"  /> "
                 << "<left-back-bottom-point  x=\"" << lbb.X() << "\" y=\""
                 << lbb.Y() << "\" z=\"" << lbb.Z() << "\"  /> "
                 << "<right-front-bottom-point  x=\"" << rfb.X() << "\" y =\""
                 << rfb.Y() << "\" z=\"" << rfb.Z() << "\"  /> "
                 << "</cuboid>";
  return xmlShapeStream.str();
Martyn Gigg's avatar
Martyn Gigg committed
}

/**
 * Create the XML required to define a cylinder from the given args
 * @param args A user-supplied dict of args
 * @param hollow True if an annulus is to be created
Martyn Gigg's avatar
Martyn Gigg committed
 * @return The XML definition string
 */
std::string
SetSample::createCylinderLikeXML(const Kernel::PropertyManager &args,
                                 const Geometry::ReferenceFrame &refFrame,
                                 bool hollow) const {
  const std::string tag = hollow ? "hollow-cylinder" : "cylinder";
  double height = args.getProperty("Height");
  double innerRadius = hollow ? args.getProperty("InnerRadius") : 0.0;
  double outerRadius =
      hollow ? args.getProperty("OuterRadius") : args.getProperty("Radius");
  std::vector<double> centre = args.getProperty("Center");
  // convert to metres
  height *= 0.01;
  innerRadius *= 0.01;
  outerRadius *= 0.01;
  std::transform(centre.begin(), centre.end(), centre.begin(),
                 [](double val) { return val *= 0.01; });
  // XML needs center position of bottom base but user specifies center of
  // cylinder
  const unsigned axisIdx = static_cast<unsigned>(refFrame.pointingUp());
  const V3D baseCentre = cylBaseCentre(centre, height, axisIdx);

  std::ostringstream xmlShapeStream;
  xmlShapeStream << "<" << tag << " id=\"sample-shape\"> "
                 << "<centre-of-bottom-base x=\"" << baseCentre.X() << "\" y=\""
                 << baseCentre.Y() << "\" z=\"" << baseCentre.Z() << "\" /> "
                 << axisXML(axisIdx) << "<height val=\"" << height << "\" /> ";
  if (hollow) {
    xmlShapeStream << "<inner-radius val=\"" << innerRadius << "\"/>"
                   << "<outer-radius val=\"" << outerRadius << "\"/>";
  } else {
    xmlShapeStream << "<radius val=\"" << outerRadius << "\"/>";
  }
  xmlShapeStream << "</" << tag << ">";
  return xmlShapeStream.str();
Martyn Gigg's avatar
Martyn Gigg committed
}

/**
 * Run SetSampleShape as an algorithm to set the shape of the sample
 * @param workspace A reference to the workspace
 * @param xml A string containing the XML definition
 */
void SetSample::runSetSampleShape(API::MatrixWorkspace_sptr &workspace,
                                  const std::string &xml) {
  auto alg = createChildAlgorithm("CreateSampleShape");
  alg->setProperty("InputWorkspace", workspace);
  alg->setProperty("ShapeXML", xml);
  alg->executeAsChildAlg();
}

/**
 * Run the named child algorithm on the given workspace. It assumes an in/out
 * workspace property called InputWorkspace
 * @param name The name of the algorithm to run
 * @param workspace A reference to the workspace
 * @param args A PropertyManager specifying the required arguments
 */
void SetSample::runChildAlgorithm(const std::string &name,
                                  API::MatrixWorkspace_sptr &workspace,
                                  const Kernel::PropertyManager &args) {
  auto alg = createChildAlgorithm(name);
  alg->setProperty("InputWorkspace", workspace);
  alg->updatePropertyValues(args);
  alg->executeAsChildAlg();
}

} // namespace DataHandling
} // namespace Mantid