Skip to content
Snippets Groups Projects
GroupDetectors2.cpp 51.4 KiB
Newer Older
    Mantid::API::IAlgorithm_sptr divide = createChildAlgorithm("Divide");
    divide->initialize();
    divide->setProperty<API::MatrixWorkspace_sptr>("LHSWorkspace", outputWS);
    divide->setProperty<API::MatrixWorkspace_sptr>("RHSWorkspace", beh);
    divide->setProperty<API::MatrixWorkspace_sptr>("OutputWorkspace", outputWS);
    divide->execute();
  }

  g_log.debug() << name() << " created " << outIndex
                << " new grouped spectra\n";
*  Move the user selected spectra in the input workspace into groups in the
* output workspace
*  @param inputWS :: user selected input workspace for the algorithm
*  @param outputWS :: user selected output workspace for the algorithm
*  @param prog4Copy :: the amount of algorithm progress to attribute to moving a
* single spectra
*  @return number of new grouped spectra
*/
size_t
GroupDetectors2::formGroupsEvent(DataObjects::EventWorkspace_const_sptr inputWS,
                                 DataObjects::EventWorkspace_sptr outputWS,
                                 const double prog4Copy) {
  // get "Behaviour" string
  const std::string behaviour = getProperty("Behaviour");
  int bhv = 0;
  if (behaviour == "Average")
    bhv = 1;

  API::MatrixWorkspace_sptr beh = API::WorkspaceFactory::Instance().create(
      "Workspace2D", static_cast<int>(m_GroupWsInds.size()), 1, 1);
  g_log.debug() << name() << ": Preparing to group spectra into "
                << m_GroupWsInds.size() << " groups\n";
  // where we are copying spectra to, we start copying to the start of the
  // output workspace
  size_t outIndex = 0;
  // Only used for averaging behaviour. We may have a 1:1 map where a Divide
  // would be waste as it would be just dividing by 1
  bool requireDivide(false);
  for (storage_map::const_iterator it = m_GroupWsInds.begin();
       it != m_GroupWsInds.end(); ++it) {
    // This is the grouped spectrum
    EventList &outEL = outputWS->getSpectrum(outIndex);

    // The spectrum number of the group is the key
    outEL.setSpectrumNo(it->first);
    // Start fresh with no detector IDs
    outEL.clearDetectorIDs();

    // the Y values and errors from spectra being grouped are combined in the
    // output spectrum
    // Keep track of number of detectors required for masking
    size_t nonMaskedSpectra(0);
    beh->dataX(outIndex)[0] = 0.0;
    beh->dataE(outIndex)[0] = 0.0;
    for (auto originalWI : it->second) {
      const EventList &fromEL = inputWS->getSpectrum(originalWI);
      // Add the event lists with the operator
      outEL += fromEL;

      // detectors to add to the output spectrum
      outEL.addDetectorIDs(fromEL.getDetectorIDs());
      try {
        Geometry::IDetector_const_sptr det = inputWS->getDetector(originalWI);
        if (!det->isMasked())
          ++nonMaskedSpectra;
      } catch (Exception::NotFoundError &) {
        // If a detector cannot be found, it cannot be masked
        ++nonMaskedSpectra;
      }
    }
    if (nonMaskedSpectra == 0)
      ++nonMaskedSpectra; // Avoid possible divide by zero
    if (!requireDivide)
      requireDivide = (nonMaskedSpectra > 1);
    beh->dataY(outIndex)[0] = static_cast<double>(nonMaskedSpectra);

    // make regular progress reports and check for cancelling the algorithm
    if (outIndex % INTERVAL == 0) {
      m_FracCompl += INTERVAL * prog4Copy;
      if (m_FracCompl > 1.0)
        m_FracCompl = 1.0;
      progress(m_FracCompl);
      interruption_point();
    }
  if (bhv == 1 && requireDivide) {
    g_log.debug() << "Running Divide algorithm to perform averaging.\n";
    Mantid::API::IAlgorithm_sptr divide = createChildAlgorithm("Divide");
    divide->initialize();
    divide->setProperty<API::MatrixWorkspace_sptr>("LHSWorkspace", outputWS);
    divide->setProperty<API::MatrixWorkspace_sptr>("RHSWorkspace", beh);
    divide->setProperty<API::MatrixWorkspace_sptr>("OutputWorkspace", outputWS);
    divide->execute();
  }

  g_log.debug() << name() << " created " << outIndex
                << " new grouped spectra\n";
// RangeHelper
/** Expands any ranges in the input string of non-negative integers, eg. "1 3-5
* 4" -> "1 3 4 5 4"
*  @param line :: a line of input that is interpreted and expanded
*  @param outList :: all integers specified both as ranges and individually in
* order
*  @throw invalid_argument if a character is found that is not an integer or
* hypehn and when a hyphen occurs at the start or the end of the line
void GroupDetectors2::RangeHelper::getList(const std::string &line,
                                           std::vector<size_t> &outList) {
  if (line.empty()) { // it is not an error to have an empty line but it would
                      // cause problems with an error check a the end of this
                      // function
  Mantid::Kernel::StringTokenizer ranges(line, "-");
Nick Draper's avatar
Nick Draper committed
    size_t loop = 0;
      Mantid::Kernel::StringTokenizer beforeHyphen(ranges[loop], " ",
                                                   IGNORE_SPACES);
      auto readPostion = beforeHyphen.begin();
      if (readPostion == beforeHyphen.end()) {
        throw std::invalid_argument("'-' found at the start of a list, can't "
                                    "interpret range specification");
      for (; readPostion != beforeHyphen.end(); ++readPostion) {
        outList.push_back(boost::lexical_cast<size_t>(*readPostion));
      // this will be the start of a range if it was followed by a - i.e.
      // another token was captured
      const size_t rangeStart = outList.back();
      if (loop + 1 == ranges.count()) { // there is no more input
      Mantid::Kernel::StringTokenizer afterHyphen(ranges[loop + 1], " ",
                                                  IGNORE_SPACES);
      readPostion = afterHyphen.begin();
      if (readPostion == afterHyphen.end()) {
        throw std::invalid_argument("A '-' follows straight after another '-', "
                                    "can't interpret range specification");
      }

      // the tokenizer will always return at least on string
      const size_t rangeEnd = boost::lexical_cast<size_t>(*readPostion);
      // this is unanticipated and marked as an error, it would be easy to
      // change this to count down however
      if (rangeStart > rangeEnd) {
        throw std::invalid_argument("A range where the first integer is larger "
                                    "than the second is not allowed");
      // expand the range
      for (size_t j = rangeStart + 1; j < rangeEnd; j++) {
        outList.push_back(j);
      }
      loop++;
    } while (loop < ranges.count());
  } catch (boost::bad_lexical_cast &e) {
    throw std::invalid_argument(
        std::string("Expected list of integers, exception thrown: ") +
        e.what());
  if (*(line.end() - 1) == '-') {
    throw std::invalid_argument(
        "'-' found at the end of a list, can't interpret range specification");
namespace {

/* The following functions are used to translate single operators into
 * groups, just like the ones this algorithm loads from .map files.
 *
 * Each function takes a string, such as "3+4", or "6:10" and then adds
 * the resulting groups of spectra to outGroups.
 */

// An add operation, i.e. "3+4" -> [3+4]
void translateAdd(const std::string &instructions,
                  std::vector<std::vector<int>> &outGroups) {
  std::vector<std::string> spectra;
  boost::split(spectra, instructions, boost::is_any_of("+"));

  std::vector<int> outSpectra;
  for (auto spectrum : spectra) {
    // remove leading/trailing whitespace
    boost::trim(spectrum);
    // add this spectrum to the group we're about to add
    outSpectra.push_back(boost::lexical_cast<int>(spectrum));
  }
  outGroups.push_back(outSpectra);
}

// A range summation, i.e. "3-6" -> [3+4+5+6]
void translateSumRange(const std::string &instructions,
                       std::vector<std::vector<int>> &outGroups) {
  // add a group with the sum of the spectra in the range
  std::vector<std::string> spectra;
  boost::split(spectra, instructions, boost::is_any_of("-"));
  if (spectra.size() != 2)
    throw std::runtime_error("Malformed range (-) operation.");
  // fetch the start and stop spectra
  int first = boost::lexical_cast<int>(spectra[0]);
  int last = boost::lexical_cast<int>(spectra[1]);
  // swap if they're back to front
  if (first > last)
    std::swap(first, last);

  // add all the spectra in the range to the output group
  std::vector<int> outSpectra;
  for (int i = first; i <= last; ++i)
    outSpectra.push_back(i);
  if (!outSpectra.empty())
    outGroups.push_back(outSpectra);
}

// A range insertion, i.e. "3:6" -> [3,4,5,6]
void translateRange(const std::string &instructions,
                    std::vector<std::vector<int>> &outGroups) {
  // add a group per spectra
  std::vector<std::string> spectra;
  boost::split(spectra, instructions, boost::is_any_of(":"));
  if (spectra.size() != 2)
    throw std::runtime_error("Malformed range (:) operation.");
  // fetch the start and stop spectra
  int first = boost::lexical_cast<int>(spectra[0]);
  int last = boost::lexical_cast<int>(spectra[1]);
  // swap if they're back to front
  if (first > last)
    std::swap(first, last);
  // add all the spectra in the range to separate output groups
  for (int i = first; i <= last; ++i) {
    // create group of size 1 with the spectrum in it
    std::vector<int> newGroup(1, i);
    // and add it to output
    outGroups.push_back(newGroup);
  }
}
} // anonymous namespace

/**
 * Used to validate the inputs for GroupDetectors2
 *
 * @returns : A map of the invalid property names to what the problem is.
 */
std::map<std::string, std::string> GroupDetectors2::validateInputs() {
  std::map<std::string, std::string> errors;

  const std::string pattern = getPropertyValue("GroupingPattern");

  boost::regex re(
      "^\\s*[0-9]+\\s*$|^(\\s*,*[0-9]+(\\s*(,|:|\\+|\\-)\\s*)*[0-9]*)*$");
  if (!pattern.empty() && !boost::regex_match(pattern, re)) {
    errors["GroupingPattern"] =
        "GroupingPattern is not well formed: " + pattern;
  }

  return errors;
}

/**
 * Translate the PerformIndexOperations processing instructions into a format
 * usable by GroupDetectors.
 *
 * @param instructions : Instructions to translate
 * @param commands : A stringstream to be filled
 */
void GroupDetectors2::translateInstructions(const std::string &instructions,
                                            std::stringstream &commands) {
  // vector of groups, each group being a vector of its spectra
  std::vector<std::vector<int>> outGroups;

  // split into comma separated groups, each group potentially containing
  // an operation (+-:) that produces even more groups.
  std::vector<std::string> groups;
  boost::split(groups, instructions, boost::is_any_of(","));

  for (auto groupStr : groups) {
    // remove leading/trailing whitespace
    boost::trim(groupStr);

    // Look for the various operators in the string. If one is found then
    // do the necessary translation into groupings.
    if (groupStr.find('+') != std::string::npos) {
      // add a group with the given spectra
      translateAdd(groupStr, outGroups);
    } else if (groupStr.find('-') != std::string::npos) {
      translateSumRange(groupStr, outGroups);
    } else if (groupStr.find(':') != std::string::npos) {
      translateRange(groupStr, outGroups);
    } else if (!groupStr.empty()) {
      // contains no instructions, just add this spectrum as a new group
      // create group of size 1 with the spectrum in it
      std::vector<int> newGroup(1, boost::lexical_cast<int>(groupStr));
      // and add it to output
      outGroups.push_back(newGroup);
    }
  }

  // We now have the groups as a vector of a vector of ints. Turn this into a
  // string, just like the contents of a map file.
  commands << outGroups.size() << "\n";
  for (auto &outGroup : outGroups) {
    const int groupId = outGroup[0] + 1;
    const int groupSize = static_cast<int>(outGroup.size());

    // Comment the output for readability
    commands << "# Group " << groupId;
    commands << ", contains " << groupSize << " spectra.\n";

    commands << groupId << "\n";
    commands << groupSize << "\n";

    // Group members
    // So far we've been using 0-indexed ids, but the mapfile syntax expects
    // 1-indexed ids, so we add 1 to the spectra ids here.
    for (size_t j = 0; j < outGroup.size(); ++j)
      commands << (j > 0 ? " " : "") << outGroup[j] + 1;