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

 /**
  * Checks if the spectra at the given index of either input workspace is masked. If so then the output spectra has zeroed data
  * and is also masked. 
  * @param lhs :: A pointer to the left-hand operand
  * @param rhs :: A pointer to the right-hand operand
  * @param index :: The workspace index to check
  * @param out :: A pointer to the output workspace
  * @returns True if further processing is not required on the spectra, false if the binary operation should be performed.
  */
 bool BinaryOperation::propagateSpectraMask(const API::MatrixWorkspace_const_sptr lhs, const API::MatrixWorkspace_const_sptr rhs, 
     const int64_t index, API::MatrixWorkspace_sptr out)
 {
   bool continueOp(true);
   IDetector_const_sptr det_lhs, det_rhs;
   try
   {
     det_lhs = lhs->getDetector(index);
     det_rhs = rhs->getDetector(index);
   }
   catch(std::runtime_error &)
   {
   }
   catch(std::domain_error &)
   {
     // try statement will throw a domain_error when the axis is not a spectra axis.
     return continueOp;
   }
   if( (det_lhs && det_lhs->isMasked()) || ( det_rhs && det_rhs->isMasked()) )
   {
     continueOp = false;
     out->maskWorkspaceIndex(index);
   }
   return continueOp;
 }

    /**
     * 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::initThetaCache(API::MatrixWorkspace_const_sptr workspace)
    {
      const size_t nhist = workspace->getNumberHistograms();
      m_thetaPts = std::vector<double>(nhist);
      size_t ndets(0);
      double minTheta(DBL_MAX), maxTheta(-DBL_MAX);

      for(int64_t i = 0 ; i < (int64_t)nhist; ++i) //signed for OpenMP
      {

        m_progress->report("Calculating detector angles");
        IDetector_const_sptr det;
        try
        {
          det = workspace->getDetector(i);
          // Check to see if there is an EFixed, if not skip it
          try
          {
            m_EmodeProperties.getEFixed(det);
          }
          catch(std::runtime_error&)
          {
            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 || det->isMonitor() )
        {
          m_thetaPts[i] = -1.0; // Indicates a detector to skip
        }
        else
        {
          ++ndets;
          const double theta = workspace->detectorTwoTheta(det);
          m_thetaPts[i] = theta;
          if( theta < minTheta )
          {
            minTheta = theta;
          }
          else if( theta > maxTheta )
          {
            maxTheta = theta;
          }
        }
      }

      m_thetaWidth = (maxTheta - minTheta)/static_cast<double>(ndets);
      g_log.information() << "Calculated detector width in theta=" << (m_thetaWidth*180.0/M_PI) << " degrees.\n";
    }

double ConvertEmptyToTof::getL2(API::MatrixWorkspace_const_sptr workspace,
                                int detId) {
  // Get a pointer to the instrument contained in the workspace
  Geometry::Instrument_const_sptr instrument = workspace->getInstrument();
  // Get the distance between the source and the sample (assume in metres)
  Geometry::IComponent_const_sptr sample = instrument->getSample();
  // Get the sample-detector distance for this detector (in metres)
  double l2 =
      workspace->getDetector(detId)->getPos().distance(sample->getPos());
  return l2;
}

/** Method updates the column, which describes if current detector/spectra is
   masked
    It is used if one tries to process multiple workspaces obtained from a
   series of experiments  where the masked detectors can change */
void PreprocessDetectorsToMD::updateMasksState(
    const API::MatrixWorkspace_const_sptr &inputWS,
    DataObjects::TableWorkspace_sptr &targWS) {
  int *pMasksArray = targWS->getColDataArray<int>("detMask");
  if (!pMasksArray)
    throw std::invalid_argument(
        "target workspace " + targWS->getName() +
        " does not have defined masks column to update");

  size_t nHist = targWS->rowCount();
  const size_t nRows = inputWS->getNumberHistograms();
  if (nHist != nRows)
    throw std::invalid_argument(
        " source workspace " + inputWS->getName() + " and target workspace " +
        targWS->getName() +
        " are inconsistent as have different numner of detectors");

  uint32_t liveDetectorsCount(0);
  for (size_t i = 0; i < nHist; i++) {
    // get detector or detector group which corresponds to the spectra i
    Geometry::IDetector_const_sptr spDet;
    try {
      spDet = inputWS->getDetector(i);
    } catch (Kernel::Exception::NotFoundError &) {
      continue;
    }

    // Check that we aren't dealing with monitor...
    if (spDet->isMonitor())
      continue;

    // if masked detectors state is not used, masked detectors just ignored;
    bool maskDetector = spDet->isMasked();
    *(pMasksArray + liveDetectorsCount) = maskDetector ? 1 : 0;

    liveDetectorsCount++;
  }
}

/** Finds the first index number of the first wavelength bin that should
* included based on the
*  the calculation: W = Wcut (Rcut-R)/Rcut
*  @param dataWS data workspace
*  @param RCut the radius cut off, should be value of the property RadiusCut
* (unit is mm)
*  @param WCut this wavelength cut off, should be equal to the value WaveCut
*  @param specInd spectrum that is being analysed
*  @return index number of the first bin to include in the calculation
*/
size_t Qhelper::waveLengthCutOff(API::MatrixWorkspace_const_sptr dataWS,
                                 const double RCut, const double WCut,
                                 const size_t specInd) const {
  double l_WCutOver = 0.0;
  double l_RCut = 0.0; // locally we store RCut in units of meters
  if (RCut > 0 && WCut > 0) {
    l_RCut = RCut / 1000.0; // convert to meters
    // l_RCut = RCut;  // convert to meters
    l_WCutOver = WCut / l_RCut;
  }

  if (!(l_RCut > 0)) {
    return 0;
  }
  // get the distance of between this detector and the origin, which should be
  // the along the beam center
  const V3D posOnBank = dataWS->getDetector(specInd)->getPos();
  double R = (posOnBank.X() * posOnBank.X()) + (posOnBank.Y() * posOnBank.Y());
  R = std::sqrt(R);

  const double WMin = l_WCutOver * (l_RCut - R);
  const MantidVec &Xs = dataWS->readX(specInd);
  return std::lower_bound(Xs.begin(), Xs.end(), WMin) - Xs.begin();
}

/**
*  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::formGroups( API::MatrixWorkspace_const_sptr inputWS, API::MatrixWorkspace_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_GroupSpecInds.size()), 1, 1);

  g_log.debug() << name() << ": Preparing to group spectra into " << m_GroupSpecInds.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_GroupSpecInds.begin(); it != m_GroupSpecInds.end() ; ++it )
  {
    // This is the grouped spectrum
    ISpectrum * outSpec = outputWS->getSpectrum(outIndex);

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

    // Copy over X data from first spectrum, the bin boundaries for all spectra are assumed to be the same here
    outSpec->dataX() = inputWS->readX(0);

    // 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( std::vector<size_t>::const_iterator wsIter = it->second.begin(); wsIter != it->second.end(); ++wsIter)
    {
      const size_t originalWI = *wsIter;

      // detectors to add to firstSpecNum
      const ISpectrum * fromSpectrum = inputWS->getSpectrum(originalWI);

      // Add up all the Y spectra and store the result in the first one
      // Need to keep the next 3 lines inside loop for now until ManagedWorkspace mru-list works properly
      MantidVec &firstY = outSpec->dataY();
      MantidVec::iterator fYit;
      MantidVec::iterator fEit = outSpec->dataE().begin();
      MantidVec::const_iterator Yit = fromSpectrum->dataY().begin();
      MantidVec::const_iterator Eit = fromSpectrum->dataE().begin();
      for (fYit = firstY.begin(); fYit != firstY.end(); ++fYit, ++fEit, ++Yit, ++Eit)
      {
        *fYit += *Yit;
        // Assume 'normal' (i.e. Gaussian) combination of errors
        *fEit = std::sqrt( (*fEit)*(*fEit) + (*Eit)*(*Eit) );
      }

      // detectors to add to the output spectrum
      outSpec->addDetectorIDs(fromSpectrum->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();
    }
    outIndex ++;
  }
  
  // Refresh the spectraDetectorMap
  outputWS->generateSpectraMap();

  if ( bhv == 1 && requireDivide )
  {
    g_log.debug() << "Running Divide algorithm to perform averaging.\n";
    Mantid::API::IAlgorithm_sptr divide = createChildAlgorithm("Divide");
    divide->initialize();
//.........这里部分代码省略.........

/** 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
      }

  /**
   * Function that retrieves the two-theta and azimuthal angles from a given
   * detector. It then looks up the nearest neighbours. Using those detectors,
   * it calculates the two-theta and azimuthal angle widths.
   * @param workspace : the workspace containing the needed detector information
   */
  void SofQW3::getValuesAndWidths(API::MatrixWorkspace_const_sptr workspace)
  {
    // Trigger a build of the nearst neighbors outside the OpenMP loop
    const int numNeighbours = 4;
    const size_t nHistos = workspace->getNumberHistograms();
    g_log.debug() << "Number of Histograms: " << nHistos << std::endl;

    this->m_theta = std::vector<double>(nHistos);
    this->m_thetaWidths = std::vector<double>(nHistos);
    this->m_phi = std::vector<double>(nHistos);
    this->m_phiWidths = std::vector<double>(nHistos);

    for (size_t i = 0; i < nHistos; ++i)
    {
      m_progress->report("Calculating detector angular widths");
      DetConstPtr detector = workspace->getDetector(i);
      g_log.debug() << "Current histogram: " << i << std::endl;
      specid_t inSpec = workspace->getSpectrum(i)->getSpectrumNo();
      SpectraDistanceMap neighbours = workspace->getNeighboursExact(inSpec,
                                                                    numNeighbours,
                                                                    true);

      g_log.debug() << "Current ID: " << inSpec << std::endl;
      // Convert from spectrum numbers to workspace indices
      double thetaWidth = -DBL_MAX;
      double phiWidth = -DBL_MAX;

      // Find theta and phi widths
      double theta = workspace->detectorTwoTheta(detector);
      double phi = detector->getPhi();

      specid_t deltaPlus1 = inSpec + 1;
      specid_t deltaMinus1 = inSpec - 1;
      specid_t deltaPlusT = inSpec + this->m_detNeighbourOffset;
      specid_t deltaMinusT = inSpec - this->m_detNeighbourOffset;

      for (SpectraDistanceMap::iterator it = neighbours.begin();
           it != neighbours.end(); ++it)
      {
        specid_t spec = it->first;
        g_log.debug() << "Neighbor ID: " << spec << std::endl;
        if (spec == deltaPlus1 || spec == deltaMinus1 ||
            spec == deltaPlusT || spec == deltaMinusT)
        {
          DetConstPtr detector_n = workspace->getDetector(spec - 1);
          double theta_n = workspace->detectorTwoTheta(detector_n);
          double phi_n = detector_n->getPhi();

          double dTheta = std::fabs(theta - theta_n);
          double dPhi = std::fabs(phi - phi_n);
          if (dTheta > thetaWidth)
          {
            thetaWidth = dTheta;
            //g_log.debug() << "Current ThetaWidth: " << thetaWidth << std::endl;
          }
          if (dPhi > phiWidth)
          {
            phiWidth = dPhi;
            //g_log.debug() << "Current PhiWidth: " << phiWidth << std::endl;
          }
        }
      }
      this->m_theta[i] = theta;
      this->m_phi[i] = phi;
      this->m_thetaWidths[i] = thetaWidth;
      this->m_phiWidths[i] = phiWidth;
    }
  }

/** method does preliminary calculations of the detectors positions to convert
results into k-dE space ;
and places the results into static cash to be used in subsequent calls to this
algorithm */
void PreprocessDetectorsToMD::processDetectorsPositions(
    const API::MatrixWorkspace_const_sptr &inputWS,
    DataObjects::TableWorkspace_sptr &targWS) {
  g_log.information()
      << "Preprocessing detector locations in a target reciprocal space\n";
  //
  Geometry::Instrument_const_sptr instrument = inputWS->getInstrument();
  // this->pBaseInstr                = instrument->baseInstrument();
  //
  Geometry::IComponent_const_sptr source = instrument->getSource();
  Geometry::IComponent_const_sptr sample = instrument->getSample();
  if ((!source) || (!sample)) {
    g_log.error() << " Instrument is not fully defined. Can not identify "
                     "source or sample\n";
    throw Kernel::Exception::InstrumentDefinitionError(
        "Instrument not sufficiently defined: failed to get source and/or "
        "sample");
  }

  // L1
  try {
    double L1 = source->getDistance(*sample);
    targWS->logs()->addProperty<double>("L1", L1, true);
    g_log.debug() << "Source-sample distance: " << L1 << std::endl;
  } catch (Kernel::Exception::NotFoundError &) {
    throw Kernel::Exception::InstrumentDefinitionError(
        "Unable to calculate source-sample distance for workspace",
        inputWS->getTitle());
  }
  // Instrument name
  std::string InstrName = instrument->getName();
  targWS->logs()->addProperty<std::string>(
      "InstrumentName", InstrName,
      true); // "The name which should unique identify current instrument");
  targWS->logs()->addProperty<bool>("FakeDetectors", false, true);

  // get access to the workspace memory
  auto &sp2detMap = targWS->getColVector<size_t>("spec2detMap");
  auto &detId = targWS->getColVector<int32_t>("DetectorID");
  auto &detIDMap = targWS->getColVector<size_t>("detIDMap");
  auto &L2 = targWS->getColVector<double>("L2");
  auto &TwoTheta = targWS->getColVector<double>("TwoTheta");
  auto &Azimuthal = targWS->getColVector<double>("Azimuthal");
  auto &detDir = targWS->getColVector<Kernel::V3D>("DetDirections");

  // Efixed; do we need one and does one exist?
  double Efi = targWS->getLogs()->getPropertyValueAsType<double>("Ei");
  float *pEfixedArray(nullptr);
  const Geometry::ParameterMap &pmap = inputWS->constInstrumentParameters();
  if (m_getEFixed)
    pEfixedArray = targWS->getColDataArray<float>("eFixed");

  // check if one needs to generate masked detectors column.
  int *pMasksArray(nullptr);
  if (m_getIsMasked)
    pMasksArray = targWS->getColDataArray<int>("detMask");

  //// progress message appearance
  size_t div = 100;
  size_t nHist = targWS->rowCount();
  Mantid::API::Progress theProgress(this, 0, 1, nHist);
  //// Loop over the spectra
  uint32_t liveDetectorsCount(0);
  for (size_t i = 0; i < nHist; i++) {
    sp2detMap[i] = std::numeric_limits<uint64_t>::quiet_NaN();
    detId[i] = std::numeric_limits<int32_t>::quiet_NaN();
    detIDMap[i] = std::numeric_limits<uint64_t>::quiet_NaN();
    L2[i] = std::numeric_limits<double>::quiet_NaN();
    TwoTheta[i] = std::numeric_limits<double>::quiet_NaN();
    Azimuthal[i] = std::numeric_limits<double>::quiet_NaN();
    //     detMask[i]  = true;

    // get detector or detector group which corresponds to the spectra i
    Geometry::IDetector_const_sptr spDet;
    try {
      spDet = inputWS->getDetector(i);
    } catch (Kernel::Exception::NotFoundError &) {
      continue;
    }

    // Check that we aren't dealing with monitor...
    if (spDet->isMonitor())
      continue;

    // if masked detectors state is not used, masked detectors just ignored;
    bool maskDetector = spDet->isMasked();
    if (m_getIsMasked)
      *(pMasksArray + liveDetectorsCount) = maskDetector ? 1 : 0;
    else if (maskDetector)
      continue;

    // calculate the requested values;
    sp2detMap[i] = liveDetectorsCount;
    detId[liveDetectorsCount] = int32_t(spDet->getID());
    detIDMap[liveDetectorsCount] = i;
    L2[liveDetectorsCount] = spDet->getDistance(*sample);
//.........这里部分代码省略.........

/** Convert the workspace units using TOF as an intermediate step in the
* conversion
* @param fromUnit :: The unit of the input workspace
* @param inputWS :: The input workspace
* @returns A shared pointer to the output workspace
*/
MatrixWorkspace_sptr ConvertUnitsUsingDetectorTable::convertViaTOF(
    Kernel::Unit_const_sptr fromUnit, API::MatrixWorkspace_const_sptr inputWS) {
  using namespace Geometry;

  // Let's see if we are using a TableWorkspace to override parameters
  TableWorkspace_sptr paramWS = getProperty("DetectorParameters");

  // See if we have supplied a DetectorParameters Workspace
  // TODO: Check if paramWS is NULL and if so throw an exception

  //      const std::string l1ColumnLabel("l1");

  // Let's check all the columns exist and are readable
  try {
    auto spectraColumnTmp = paramWS->getColumn("spectra");
    auto l1ColumnTmp = paramWS->getColumn("l1");
    auto l2ColumnTmp = paramWS->getColumn("l2");
    auto twoThetaColumnTmp = paramWS->getColumn("twotheta");
    auto efixedColumnTmp = paramWS->getColumn("efixed");
    auto emodeColumnTmp = paramWS->getColumn("emode");
  } catch (...) {
    throw Exception::InstrumentDefinitionError(
        "DetectorParameter TableWorkspace is not defined correctly.");
  }

  // Now let's take a reference to the vectors.
  const auto &l1Column = paramWS->getColVector<double>("l1");
  const auto &l2Column = paramWS->getColVector<double>("l2");
  const auto &twoThetaColumn = paramWS->getColVector<double>("twotheta");
  const auto &efixedColumn = paramWS->getColVector<double>("efixed");
  const auto &emodeColumn = paramWS->getColVector<int>("emode");
  const auto &spectraColumn = paramWS->getColVector<int>("spectra");

  Progress prog(this, 0.2, 1.0, m_numberOfSpectra);
  int64_t numberOfSpectra_i =
      static_cast<int64_t>(m_numberOfSpectra); // cast to make openmp happy

  // Get the unit object for each workspace
  Kernel::Unit_const_sptr outputUnit = m_outputUnit;
  std::vector<double> emptyVec;
  int failedDetectorCount = 0;

  // Perform Sanity Validation before creating workspace
  size_t checkIndex = 0;
  int checkSpecNo = inputWS->getDetector(checkIndex)->getID();
  auto checkSpecIter =
      std::find(spectraColumn.begin(), spectraColumn.end(), checkSpecNo);
  if (checkSpecIter != spectraColumn.end()) {
    size_t detectorRow = std::distance(spectraColumn.begin(), checkSpecIter);
    // copy the X values for the check
    auto checkXValues = inputWS->readX(checkIndex);
    // Convert the input unit to time-of-flight
    auto checkFromUnit = std::unique_ptr<Unit>(fromUnit->clone());
    auto checkOutputUnit = std::unique_ptr<Unit>(outputUnit->clone());
    double checkdelta = 0;
    checkFromUnit->toTOF(checkXValues, emptyVec, l1Column[detectorRow],
                         l2Column[detectorRow], twoThetaColumn[detectorRow],
                         emodeColumn[detectorRow], efixedColumn[detectorRow],
                         checkdelta);
    // Convert from time-of-flight to the desired unit
    checkOutputUnit->fromTOF(checkXValues, emptyVec, l1Column[detectorRow],
                             l2Column[detectorRow], twoThetaColumn[detectorRow],
                             emodeColumn[detectorRow],
                             efixedColumn[detectorRow], checkdelta);
  }

  // create the output workspace
  MatrixWorkspace_sptr outputWS = this->setupOutputWorkspace(inputWS);
  EventWorkspace_sptr eventWS =
      boost::dynamic_pointer_cast<EventWorkspace>(outputWS);
  assert(static_cast<bool>(eventWS) == m_inputEvents); // Sanity check

  // TODO: Check why this parallel stuff breaks
  // Loop over the histograms (detector spectra)
  // PARALLEL_FOR_IF(Kernel::threadSafe(*outputWS))
  for (int64_t i = 0; i < numberOfSpectra_i; ++i) {

    // Lets find what row this spectrum Number appears in our detector table.

    // PARALLEL_START_INTERUPT_REGION

    std::size_t wsid = i;

    try {

      double deg2rad = M_PI / 180.;

      auto det = outputWS->getDetector(i);
      int specNo = det->getID();

      // int spectraNumber = static_cast<int>(spectraColumn->toDouble(i));
      // wsid = outputWS->getIndexFromSpectrumNumber(spectraNumber);
      g_log.debug() << "###### Spectra #" << specNo
                    << " ==> Workspace ID:" << wsid << '\n';
//.........这里部分代码省略.........