C++ MatrixWorkspace_const_sptr::blocksize方法代码示例

/**
 * Return true if the two workspaces are compatible for this operation
 * Virtual: will be overridden as needed.
 * @param lhs :: left-hand workspace to check
 * @param rhs :: right-hand workspace to check
 * @return flag for the compatibility to the two workspaces
 */
bool BinaryOperation::checkCompatibility(
    const API::MatrixWorkspace_const_sptr lhs,
    const API::MatrixWorkspace_const_sptr rhs) const {
  Unit_const_sptr lhs_unit;
  Unit_const_sptr rhs_unit;
  if (lhs->axes() && rhs->axes()) // If one of these is a WorkspaceSingleValue
                                  // then we don't want to check units match
  {
    lhs_unit = lhs->getAxis(0)->unit();
    rhs_unit = rhs->getAxis(0)->unit();
  }

  const std::string lhs_unitID = (lhs_unit ? lhs_unit->unitID() : "");
  const std::string rhs_unitID = (rhs_unit ? rhs_unit->unitID() : "");

  // Check the workspaces have the same units and distribution flag
  if (lhs_unitID != rhs_unitID && lhs->blocksize() > 1 &&
      rhs->blocksize() > 1) {
    g_log.error("The two workspace are not compatible because they have "
                "different units on the X axis.");
    return false;
  }

  // Check the size compatibility
  std::string checkSizeCompatibilityResult = checkSizeCompatibility(lhs, rhs);
  if (!checkSizeCompatibilityResult.empty()) {
    throw std::invalid_argument(checkSizeCompatibilityResult);
  }

  return true;
}

    /**
     * Return true if the two workspaces are compatible for this operation
     * Virtual: will be overridden as needed.
     * @param lhs :: left-hand workspace to check
     * @param rhs :: right-hand workspace to check
     * @return flag for the compatibility to the two workspaces
     */
    bool BinaryOperation::checkCompatibility(const API::MatrixWorkspace_const_sptr lhs,const API::MatrixWorkspace_const_sptr rhs) const
    {
      Unit_const_sptr lhs_unit;
      Unit_const_sptr rhs_unit;
      if ( lhs->axes() && rhs->axes() ) // If one of these is a WorkspaceSingleValue then we don't want to check units match
      {
        lhs_unit = lhs->getAxis(0)->unit();
        rhs_unit = rhs->getAxis(0)->unit();
      }

      const std::string lhs_unitID = ( lhs_unit ? lhs_unit->unitID() : "" );
      const std::string rhs_unitID = ( rhs_unit ? rhs_unit->unitID() : "" );

      // Check the workspaces have the same units and distribution flag
      if ( lhs_unitID != rhs_unitID && lhs->blocksize() > 1 && rhs->blocksize() > 1 )
      {
        g_log.error("The two workspace are not compatible because they have different units on the X axis.");
        return false;
      }

      // Check the size compatibility
      if (!checkSizeCompatibility(lhs,rhs))
      {
        std::ostringstream ostr;
        ostr<<"The sizes of the two workspaces " <<
            "(" << lhs->getName() << ": " << lhs->getNumberHistograms() << " spectra, blocksize " << lhs->blocksize() << ")"
            << " and " <<
            "(" << rhs->getName() << ": " << rhs->getNumberHistograms() << " spectra, blocksize " << rhs->blocksize() << ")"
            << " are not compatible for algorithm "<<this->name();
        g_log.error() << ostr.str() << std::endl;
        throw std::invalid_argument( ostr.str() );
      }

      return true;
    }

/** Performs a simple check to see if the sizes of two workspaces are compatible
 * for a binary operation
 *  In order to be size compatible then the larger workspace
 *  must divide be the size of the smaller workspace leaving no remainder
 *  @param lhs :: the first workspace to compare
 *  @param rhs :: the second workspace to compare
 *  @retval "" The two workspaces are size compatible
 *  @retval "<reason why not compatible>" The two workspaces are NOT size
 * compatible
 */
std::string BinaryOperation::checkSizeCompatibility(
    const API::MatrixWorkspace_const_sptr lhs,
    const API::MatrixWorkspace_const_sptr rhs) const {
  const size_t lhsSize = lhs->size();
  const size_t rhsSize = rhs->size();
  // A SingleValueWorkspace on the right matches anything
  if (rhsSize == 1)
    return "";
  // The lhs must not be smaller than the rhs
  if (lhsSize < rhsSize)
    return "Left hand side smaller than right hand side.";

  // Did checkRequirements() tell us that the X histogram size did not matter?
  if (!m_matchXSize) {
    // If so, only the vertical # needs to match

    if (lhs->getNumberHistograms() == rhs->getNumberHistograms()) {
      return "";
    } else {
      return "Number of histograms not identical.";
    }
  }
  // Otherwise they must match both ways, or horizontally or vertically with the
  // other rhs dimension=1
  if (rhs->blocksize() == 1 &&
      lhs->getNumberHistograms() == rhs->getNumberHistograms())
    return "";
  // Past this point, we require the X arrays to match. Note this only checks
  // the first spectrum
  if (!WorkspaceHelpers::matchingBins(*lhs, *rhs, true)) {
    return "X arrays must match when performing this operation on a 2D "
           "workspaces.";
  }

  const size_t rhsSpec = rhs->getNumberHistograms();

  if (lhs->blocksize() == rhs->blocksize()) {
    if (rhsSpec == 1 || lhs->getNumberHistograms() == rhsSpec) {
      return "";
    } else {
      // can't be more specific as if this is reached both failed and only one
      // or both are needed
      return "Left and right sides should contain the same amount of spectra "
             "or the right side should contian only one spectra.";
    }
  } else {
    // blocksize check failed, but still check the number of spectra to see if
    // that was wrong too
    if (rhsSpec == 1 || lhs->getNumberHistograms() == rhsSpec) {
      return "Number of y values not equal on left and right sides.";
    } else {
      // can't be more specific as if this is reached both failed and only one
      // or both are needed
      return "Number of y values not equal on left and right sides and the "
             "right side contained neither only one spectra or the same amount "
             "of spectra as the left.";
    }
  }
}

/** Calculates rebin parameters: the min and max bin boundaries and the
   logarithmic step. The aim is to have approx.
    the same number of bins as in the input workspace.
    @param workspace :: The workspace being rebinned
    @param min ::       (return) The calculated frame starting point
    @param max ::       (return) The calculated frame ending point
    @param step ::      (return) The calculated bin width
 */
void DiffractionFocussing::calculateRebinParams(
    const API::MatrixWorkspace_const_sptr &workspace, double &min, double &max,
    double &step) {

  min = std::numeric_limits<double>::max();
  // for min and max we need to iterate over the data block and investigate each
  // one
  int64_t length = workspace->getNumberHistograms();
  for (int64_t i = 0; i < length; i++) {
    const MantidVec &xVec = workspace->readX(i);
    const double &localMin = xVec[0];
    const double &localMax = xVec[xVec.size() - 1];
    if (localMin != std::numeric_limits<double>::infinity() &&
        localMax != std::numeric_limits<double>::infinity()) {
      if (localMin < min)
        min = localMin;
      if (localMax > max)
        max = localMax;
    }
  }

  if (min <= 0.)
    min = 1e-6;

  // step is easy
  double n = static_cast<double>(workspace->blocksize());
  step = (log(max) - log(min)) / n;
}

/** Do the initial copy of the data from the input to the output workspace for
 * histogram workspaces.
 *  Takes out the bin width if necessary.
 *  @param inputWS  The input workspace
 *  @param outputWS The output workspace
 */
void ConvertUnitsUsingDetectorTable::fillOutputHist(
    const API::MatrixWorkspace_const_sptr inputWS,
    const API::MatrixWorkspace_sptr outputWS) {
  const int size = static_cast<int>(inputWS->blocksize());

  // Loop over the histograms (detector spectra)
  Progress prog(this, 0.0, 0.2, m_numberOfSpectra);
  int64_t numberOfSpectra_i =
      static_cast<int64_t>(m_numberOfSpectra); // cast to make openmp happy
  PARALLEL_FOR2(inputWS, outputWS)
  for (int64_t i = 0; i < numberOfSpectra_i; ++i) {
    PARALLEL_START_INTERUPT_REGION
    // Take the bin width dependency out of the Y & E data
    if (m_distribution) {
      for (int j = 0; j < size; ++j) {
        const double width =
            std::abs(inputWS->dataX(i)[j + 1] - inputWS->dataX(i)[j]);
        outputWS->dataY(i)[j] = inputWS->dataY(i)[j] * width;
        outputWS->dataE(i)[j] = inputWS->dataE(i)[j] * width;
      }
    } else {
      // Just copy over
      outputWS->dataY(i) = inputWS->readY(i);
      outputWS->dataE(i) = inputWS->readE(i);
    }
    // Copy over the X data
    outputWS->setX(i, inputWS->refX(i));

    prog.report("Convert to " + m_outputUnit->unitID());
    PARALLEL_END_INTERUPT_REGION
  }
  PARALLEL_CHECK_INTERUPT_REGION
}

/** Calculates rebin parameters: the min and max bin boundaries and the
   logarithmic step. The aim is to have approx.
    the same number of bins as in the input workspace.
    @param workspace :: The workspace being rebinned
    @param min ::       (return) The calculated frame starting point
    @param max ::       (return) The calculated frame ending point
    @param step ::      (return) The calculated bin width
 */
void DiffractionFocussing::calculateRebinParams(
    const API::MatrixWorkspace_const_sptr &workspace, double &min, double &max,
    double &step) {

  min = std::numeric_limits<double>::max();
  // for min and max we need to iterate over the data block and investigate each
  // one
  int64_t length = workspace->getNumberHistograms();
  for (int64_t i = 0; i < length; i++) {
    auto &xVec = workspace->x(i);
    const double &localMin = xVec.front();
    const double &localMax = xVec.back();
    if (std::isfinite(localMin) && std::isfinite(localMax)) {
      min = std::min(min, localMin);
      max = std::max(max, localMax);
    }
  }

  if (min <= 0.)
    min = 1e-6;

  // step is easy
  double n = static_cast<double>(workspace->blocksize());
  step = (log(max) - log(min)) / n;
}

/**
 * Creates the output workspace for this algorithm
 * @param inputWorkspace A parent workspace to initialize from.
 * @return A pointer to the output workspace.
 */
API::MatrixWorkspace_sptr Transpose::createOutputWorkspace(
    API::MatrixWorkspace_const_sptr inputWorkspace) {
  Mantid::API::Axis *yAxis = getVerticalAxis(inputWorkspace);
  const size_t oldNhist = inputWorkspace->getNumberHistograms();
  const auto &inX = inputWorkspace->x(0);
  const size_t oldYlength = inputWorkspace->blocksize();
  const size_t oldVerticalAxislength = yAxis->length();

  // The input Y axis may be binned so the new X data should be too
  size_t newNhist(oldYlength), newXsize(oldVerticalAxislength),
      newYsize(oldNhist);
  MatrixWorkspace_sptr outputWorkspace = inputWorkspace->cloneEmpty();
  outputWorkspace->initialize(newNhist, newXsize, newYsize);
  outputWorkspace->setTitle(inputWorkspace->getTitle());
  outputWorkspace->setComment(inputWorkspace->getComment());
  outputWorkspace->copyExperimentInfoFrom(inputWorkspace.get());
  outputWorkspace->setYUnit(inputWorkspace->YUnit());
  outputWorkspace->setYUnitLabel(inputWorkspace->YUnitLabel());
  outputWorkspace->setDistribution(inputWorkspace->isDistribution());

  // Create a new numeric axis for Y the same length as the old X array
  // Values come from input X
  API::NumericAxis *newYAxis(nullptr);
  if (inputWorkspace->isHistogramData()) {
    newYAxis = new API::BinEdgeAxis(inX.rawData());
  } else {
    newYAxis = new API::NumericAxis(inX.rawData());
  }

  newYAxis->unit() = inputWorkspace->getAxis(0)->unit();
  outputWorkspace->getAxis(0)->unit() = inputWorkspace->getAxis(1)->unit();
  outputWorkspace->replaceAxis(1, newYAxis);
  setProperty("OutputWorkspace", outputWorkspace);
  return outputWorkspace;
}

/**
 * Creates the output workspace for this algorithm
 * @param inputWorkspace A parent workspace to initialize from.
 * @return A pointer to the output workspace.
 */
API::MatrixWorkspace_sptr Transpose::createOutputWorkspace(
    API::MatrixWorkspace_const_sptr inputWorkspace) {
  Mantid::API::Axis *yAxis = getVerticalAxis(inputWorkspace);
  const size_t oldNhist = inputWorkspace->getNumberHistograms();
  const MantidVec &inX = inputWorkspace->readX(0);
  const size_t oldYlength = inputWorkspace->blocksize();
  const size_t oldVerticalAxislength = yAxis->length();

  // The input Y axis may be binned so the new X data should be too
  size_t newNhist(oldYlength), newXsize(oldVerticalAxislength),
      newYsize(oldNhist);
  MatrixWorkspace_sptr outputWorkspace = WorkspaceFactory::Instance().create(
      inputWorkspace, newNhist, newXsize, newYsize);

  // Create a new numeric axis for Y the same length as the old X array
  // Values come from input X
  API::NumericAxis *newYAxis(nullptr);
  if (inputWorkspace->isHistogramData()) {
    newYAxis = new API::BinEdgeAxis(inX);
  } else {
    newYAxis = new API::NumericAxis(inX);
  }

  newYAxis->unit() = inputWorkspace->getAxis(0)->unit();
  outputWorkspace->getAxis(0)->unit() = inputWorkspace->getAxis(1)->unit();
  outputWorkspace->replaceAxis(1, newYAxis);
  setProperty("OutputWorkspace", outputWorkspace);
  return outputWorkspace;
}

/** Performs a simple check to see if the sizes of two workspaces are compatible
 *for a binary operation
 *  In order to be size compatible then the larger workspace
 *  must divide be the size of the smaller workspace leaving no remainder
 *
 *  @param lhs :: the first workspace to compare
 *  @param rhs :: the second workspace to compare
 *  @retval true The two workspaces are size compatible
 *  @retval false The two workspaces are NOT size compatible
 */
std::string Multiply::checkSizeCompatibility(
    const API::MatrixWorkspace_const_sptr lhs,
    const API::MatrixWorkspace_const_sptr rhs) const {
  if (!m_keepEventWorkspace && !m_AllowDifferentNumberSpectra) {
    // Fallback on the default checks
    return CommutativeBinaryOperation::checkSizeCompatibility(lhs, rhs);
  } else {

    // A SingleValueWorkspace on the right, or matches anything
    if (rhs->size() == 1)
      return "";

    // A SingleValueWorkspace on the left only matches if rhs was single value
    // too. Why are you using mantid to do simple math?!?
    if (lhs->size() == 1)
      return "The left side cannot contain a single value if the right side "
             "isn't also a single value.";

    // RHS only has one value (1D vertical), so the number of histograms needs
    // to match.
    // Each lhs spectrum will be divided by that scalar
    if (rhs->blocksize() == 1 &&
        lhs->getNumberHistograms() == rhs->getNumberHistograms())
      return "";

    if (m_matchXSize) {
      // Past this point, for a 2D WS operation, we require the X arrays to
      // match. Note this only checks the first spectrum
      if (!WorkspaceHelpers::matchingBins(*lhs, *rhs, true)) {
        return "X arrays must match when multiplying 2D workspaces.";
      }
    }

    // We don't need to check for matching bins for events. Yay events!
    const size_t rhsSpec = rhs->getNumberHistograms();

    // If the rhs has a single spectrum, then we can divide. The block size does
    // NOT need to match,
    if (rhsSpec == 1)
      return "";

    // Are we allowing the division by different # of spectra, using detector
    // IDs to match up?
    if (m_AllowDifferentNumberSpectra) {
      return "";
    }

    // Otherwise, the number of histograms needs to match, but the block size of
    // each does NOT need to match.

    if (lhs->getNumberHistograms() == rhs->getNumberHistograms()) {
      return "";
    } else {
      return "Number of histograms not identical.";
    }
  }
}

/** Performs a simple check to see if the sizes of two workspaces are compatible
 *for a binary operation
 *  In order to be size compatible then the larger workspace
 *  must divide be the size of the smaller workspace leaving no remainder
 *
 *  @param lhs :: the first workspace to compare
 *  @param rhs :: the second workspace to compare
 *  @retval "" The two workspaces are size compatible
 *  @retval "<reason why not compatible>" The two workspaces are NOT size
 *compatible
 */
std::string Divide::checkSizeCompatibility(
    const API::MatrixWorkspace_const_sptr lhs,
    const API::MatrixWorkspace_const_sptr rhs) const {
  // --- Check for event workspaces - different than workspaces 2D! ---

  // A SingleValueWorkspace on the right matches anything
  if (rhs->size() == 1)
    return "";

  // A SingleValueWorkspace on the left only matches if rhs was single value
  // too. Why are you using mantid to do simple math?!?
  if (lhs->size() == 1)
    return "The left side cannot contain a single value if the right side "
           "isn't also a single value.";

  // If RHS only has one value (1D vertical), the number of histograms needs to
  // match.
  // Each lhs spectrum will be divided by that scalar
  // std::cout << "rhs->blocksize() " << rhs->blocksize() << std::endl;
  // Are we allowing the division by different # of spectra, using detector IDs
  // to match up?
  if (m_AllowDifferentNumberSpectra ||
      (rhs->blocksize() == 1 &&
       lhs->getNumberHistograms() == rhs->getNumberHistograms())) {
    return "";
  }

  if (m_matchXSize) {
    // Past this point, for a 2D WS operation, we require the X arrays to match.
    // Note this only checks the first spectrum
    if (!WorkspaceHelpers::matchingBins(lhs, rhs, true)) {
      return "X arrays must match when dividing 2D workspaces.";
    }
  }

  // We don't need to check for matching bins for events. Yay events!
  const size_t rhsSpec = rhs->getNumberHistograms();

  // If the rhs has a single spectrum, then we can divide. The block size does
  // NOT need to match,
  if (rhsSpec == 1)
    return "";

  // Otherwise, the number of histograms needs to match, but the block size of
  // each does NOT need to match.

  if (lhs->getNumberHistograms() == rhs->getNumberHistograms()) {
    return "";
  } else {
    return "Number of histograms not identical.";
  }
}

    /** Performs a simple check to see if the sizes of two workspaces are compatible for a binary operation
     *  In order to be size compatible then the larger workspace
     *  must divide be the size of the smaller workspace leaving no remainder
     *  @param lhs :: the first workspace to compare
     *  @param rhs :: the second workspace to compare
     *  @retval true The two workspaces are size compatible
     *  @retval false The two workspaces are NOT size compatible
     */
    bool BinaryOperation::checkSizeCompatibility(const API::MatrixWorkspace_const_sptr lhs,const API::MatrixWorkspace_const_sptr rhs) const
    {
      const size_t lhsSize = lhs->size();
      const size_t rhsSize = rhs->size();
      // A SingleValueWorkspace on the right matches anything
      if ( rhsSize == 1 ) return true;
      // The rhs must not be smaller than the lhs
      if ( lhsSize < rhsSize ) return false;

      //Did checkRequirements() tell us that the X histogram size did not matter?
      if (!m_matchXSize)
        //If so, only the vertical # needs to match
        return (lhs->getNumberHistograms() == rhs->getNumberHistograms());

      // Otherwise they must match both ways, or horizontally or vertically with the other rhs dimension=1
      if ( rhs->blocksize() == 1 && lhs->getNumberHistograms() == rhs->getNumberHistograms() ) return true;
      // Past this point, we require the X arrays to match. Note this only checks the first spectrum
      if ( !WorkspaceHelpers::matchingBins(lhs,rhs,true) ) return false;
      
      const size_t rhsSpec = rhs->getNumberHistograms();

      return ( lhs->blocksize() == rhs->blocksize() && ( rhsSpec==1 || lhs->getNumberHistograms() == rhsSpec ) );
    }

/** Checks if workspaces input to Q1D or Qxy are reasonable
  @param dataWS data workspace
  @param binWS (WavelengthAdj) workpace that will be checked to see if it has one spectrum and the same number of bins as dataWS
  @param detectWS (PixelAdj) passing NULL for this wont raise an error, if set it will be checked this workspace has as many histograms as dataWS each with one bin
  @throw invalid_argument if the workspaces are not mututially compatible
*/
void Qhelper::examineInput(API::MatrixWorkspace_const_sptr dataWS, 
     API::MatrixWorkspace_const_sptr binAdj, API::MatrixWorkspace_const_sptr detectAdj)
{
  if ( dataWS->getNumberHistograms() < 1 )
  {
    throw std::invalid_argument("Empty data workspace passed, can not continue");
  }

  //it is not an error for these workspaces not to exist
  if (binAdj)
  {
    if ( binAdj->getNumberHistograms() != 1 )
    {
      throw std::invalid_argument("The WavelengthAdj workspace must have one spectrum");
    }
    if ( binAdj->readY(0).size() != dataWS->readY(0).size() )
    {
      throw std::invalid_argument("The WavelengthAdj workspace's bins must match those of the detector bank workspace");
    }
    MantidVec::const_iterator reqX = dataWS->readX(0).begin();
    MantidVec::const_iterator testX = binAdj->readX(0).begin();
    for ( ; reqX != dataWS->readX(0).end(); ++reqX, ++testX)
    {
      if ( *reqX != *testX )
      {
        throw std::invalid_argument("The WavelengthAdj workspace must have matching bins with the detector bank workspace");
      }
    }
    if ( binAdj->isDistribution() != dataWS->isDistribution() )
    {
      throw std::invalid_argument("The distrbution/raw counts status of the wavelengthAdj and DetBankWorkspace must be the same, use ConvertToDistribution");
    }
  }
  else if( ! dataWS->isDistribution() )
  {
    //throw std::invalid_argument("The data workspace must be a distrbution if there is no Wavelength dependent adjustment");
  }
  
  if (detectAdj)
  {
    if ( detectAdj->blocksize() != 1 )
    {
      throw std::invalid_argument("The PixelAdj workspace must point to a workspace with single bin spectra, as only the first bin is used");
    }
    if ( detectAdj->getNumberHistograms() != dataWS->getNumberHistograms() )
    {
      throw std::invalid_argument("The PixelAdj workspace must have one spectrum for each spectrum in the detector bank workspace");
    }
  }
}

 void Divide::setOutputUnits(const API::MatrixWorkspace_const_sptr lhs,const API::MatrixWorkspace_const_sptr rhs,API::MatrixWorkspace_sptr out)
 {
   if ( rhs->YUnit().empty() || !WorkspaceHelpers::matchingBins(lhs,rhs,true) )
   {
     // Do nothing
   }
   // If the Y units match, then the output will be a distribution and will be dimensionless
   else if ( lhs->YUnit() == rhs->YUnit() && rhs->blocksize() > 1 )
   {
     out->setYUnit("");
     out->isDistribution(true);
   }
   // Else we need to set the unit that results from the division
   else
   {
     if ( ! lhs->YUnit().empty() ) out->setYUnit(lhs->YUnit() + "/" + rhs->YUnit());
     else out->setYUnit("1/" + rhs->YUnit());
   }
 }

/** Checks if workspaces input to Q1D or Qxy are reasonable
  @param dataWS data workspace
  @param binAdj (WavelengthAdj) workpace that will be checked to see if it has
  one spectrum and the same number of bins as dataWS
  @param detectAdj (PixelAdj) passing NULL for this wont raise an error, if set
  it will be checked this workspace has as many histograms as dataWS each with
  one bin
  @throw invalid_argument if the workspaces are not mututially compatible
*/
void Qhelper::examineInput(API::MatrixWorkspace_const_sptr dataWS,
                           API::MatrixWorkspace_const_sptr binAdj,
                           API::MatrixWorkspace_const_sptr detectAdj) {
  if (dataWS->getNumberHistograms() < 1) {
    throw std::invalid_argument(
        "Empty data workspace passed, can not continue");
  }

  // it is not an error for these workspaces not to exist
  if (binAdj) {
    if (binAdj->getNumberHistograms() != 1) {
      throw std::invalid_argument(
          "The WavelengthAdj workspace must have one spectrum");
    }
    if (binAdj->readY(0).size() != dataWS->readY(0).size()) {
      throw std::invalid_argument("The WavelengthAdj workspace's bins must "
                                  "match those of the detector bank workspace");
    }
    MantidVec::const_iterator reqX = dataWS->readX(0).begin();
    MantidVec::const_iterator testX = binAdj->readX(0).begin();
    for (; reqX != dataWS->readX(0).end(); ++reqX, ++testX) {
      if (*reqX != *testX) {
        throw std::invalid_argument("The WavelengthAdj workspace must have "
                                    "matching bins with the detector bank "
                                    "workspace");
      }
    }
    if (binAdj->isDistribution() != dataWS->isDistribution()) {
      throw std::invalid_argument("The distrbution/raw counts status of the "
                                  "wavelengthAdj and DetBankWorkspace must be "
                                  "the same, use ConvertToDistribution");
    }
  } else if (!dataWS->isDistribution()) {
    // throw std::invalid_argument("The data workspace must be a distrbution if
    // there is no Wavelength dependent adjustment");
  }

  // Perform tests on detectAdj

  if (detectAdj) {
    if (detectAdj->blocksize() != 1) {
      throw std::invalid_argument("The PixelAdj workspace must point to a "
                                  "workspace with single bin spectra, as only "
                                  "the first bin is used");
    }
    if (detectAdj->getNumberHistograms() != dataWS->getNumberHistograms()) {
      throw std::invalid_argument("The PixelAdj workspace must have one "
                                  "spectrum for each spectrum in the detector "
                                  "bank workspace");
    }

    // test that when detector adjustment value less than or equal to zero that
    // the corresponding detector
    // in the workspace is masked

    size_t num_histograms = dataWS->getNumberHistograms();
    for (size_t i = 0; i < num_histograms; i++) {
      double adj = (double)detectAdj->readY(i)[0];
      if (adj <= 0.0) {
        bool det_is_masked;

        try {
          det_is_masked = dataWS->getDetector(i)->isMasked();
        } catch (...) {
          // just ignore. There are times, when the detector is not masked
          // because it does not exist at all.
          det_is_masked = true;
        }
        if (!det_is_masked) {
          throw std::invalid_argument(
              "Every detector with non-positive PixelAdj value must be masked");
        }
      }
    }
  }
}

    /**
     * A map detector ID and Q ranges
     * This method looks unnecessary as it could be calculated on the fly but
     * the parallelization means that lazy instantation slows it down due to the 
     * necessary CRITICAL sections required to update the cache. The Q range 
     * values are required very frequently so the total time is more than
     * offset by this precaching step
     */
    void SofQW2::initQCache(API::MatrixWorkspace_const_sptr workspace)
    {
      Mantid::Kernel::Timer clock;
      const size_t nhist(workspace->getNumberHistograms());
      const size_t nxpoints = workspace->blocksize();
      const MantidVec & X = workspace->readX(0);
      m_qcached.clear();

      PARALLEL_FOR1(workspace)
      for(int64_t i = 0 ; i < (int64_t)nhist; ++i)
      {
        PARALLEL_START_INTERUPT_REGION

        IDetector_const_sptr det;
        try
        {
           det = workspace->getDetector(i);
           if( det->isMonitor() ) det.reset();
        }
        catch(Kernel::Exception::NotFoundError&)
        {
          // Catch if no detector. Next line tests whether this happened - test placed
          // outside here because Mac Intel compiler doesn't like 'continue' in a catch
          // in an openmp block.
        }
        // If no detector found, skip onto the next spectrum
        if ( !det ) continue;

        std::vector<QValues> qvalues(nxpoints);
        DetectorGroup_const_sptr detGroup = boost::dynamic_pointer_cast<const DetectorGroup>(det);
        if( detGroup )
        {
          std::vector<IDetector_const_sptr> dets = detGroup->getDetectors();
          const size_t ndets(dets.size());
          for( size_t j = 0; j < ndets; ++j )
          {
            IDetector_const_sptr det_j = dets[j];
            QRangeCache qrange(static_cast<size_t>(i), 1.0/(double)ndets);
            for( size_t k = 0; k < nxpoints; ++k)
            {
              qvalues[k] = calculateQValues(det_j, X[k], X[k+1]);
            }
            qrange.qValues = qvalues;
            PARALLEL_CRITICAL(qcache_a)
            {
              m_qcached.insert(m_qcached.end(), qrange);
            }
          }
        }
        else
        {
          QRangeCache qrange(static_cast<size_t>(i), 1.0);
          for( size_t k = 0; k < nxpoints; ++k)
          {
            qvalues[k] = calculateQValues(det, X[k], X[k+1]);
          }
          qrange.qValues = qvalues;
          PARALLEL_CRITICAL(qcache_b)
          {
            m_qcached.insert(m_qcached.end(), qrange);
          }
        }

        PARALLEL_END_INTERUPT_REGION
      }