/** Parse table workspace
  */
void RemovePeaks::parsePeakTableWorkspace(TableWorkspace_sptr peaktablews,
                                          vector<double> &vec_peakcentre,
                                          vector<double> &vec_peakfwhm) {
  // Get peak table workspace information
  vector<string> colnames = peaktablews->getColumnNames();
  int index_centre = -1;
  int index_fwhm = -1;
  for (int i = 0; i < static_cast<int>(colnames.size()); ++i) {
    string colname = colnames[i];
    if (colname.compare("TOF_h") == 0)
      index_centre = i;
    else if (colname.compare("FWHM") == 0)
      index_fwhm = i;
  }

  if (index_centre < 0 || index_fwhm < 0) {
    throw runtime_error(
        "Input Bragg peak table workspace does not have TOF_h and/or FWHM");
  }

  // Get values
  size_t numrows = peaktablews->rowCount();
  vec_peakcentre.resize(numrows, 0.);
  vec_peakfwhm.resize(numrows, 0.);

  for (size_t i = 0; i < numrows; ++i) {
    double centre = peaktablews->cell<double>(i, index_centre);
    double fwhm = peaktablews->cell<double>(i, index_fwhm);
    vec_peakcentre[i] = centre;
    vec_peakfwhm[i] = fwhm;
  }

  return;
}

  /** Select background points via a given background function
    */
  void ProcessBackground::execSelectBkgdPoints2()
  {
    // Process properties
    BackgroundFunction_sptr bkgdfunc = createBackgroundFunction(m_bkgdType);
    TableWorkspace_sptr bkgdtablews = getProperty("BackgroundTableWorkspace");

    // Set up background function from table
    size_t numrows = bkgdtablews->rowCount();
    map<string, double> parmap;
    for (size_t i = 0; i < numrows; ++i)
    {
      TableRow row = bkgdtablews->getRow(i);
      string parname;
      double parvalue;
      row >> parname >> parvalue;
      if (parname[0] == 'A')
        parmap.insert(make_pair(parname, parvalue));
    }

    int bkgdorder = static_cast<int>(parmap.size()-1); // A0 - A(n) total n+1 parameters
    bkgdfunc->setAttributeValue("n", bkgdorder);
    for (map<string, double>::iterator mit = parmap.begin(); mit != parmap.end(); ++mit)
    {
      string parname = mit->first;
      double parvalue = mit->second;
      bkgdfunc->setParameter(parname, parvalue);
    }

    // Filter out
    m_outputWS = filterForBackground(bkgdfunc);

    return;
  }

/** Select background points via a given background function
  */
void ProcessBackground::selectFromGivenFunction() {
  // Process properties
  BackgroundFunction_sptr bkgdfunc = createBackgroundFunction(m_bkgdType);
  TableWorkspace_sptr bkgdtablews = getProperty("BackgroundTableWorkspace");

  // Set up background function from table
  size_t numrows = bkgdtablews->rowCount();
  map<string, double> parmap;
  for (size_t i = 0; i < numrows; ++i) {
    TableRow row = bkgdtablews->getRow(i);
    string parname;
    double parvalue;
    row >> parname >> parvalue;
    if (parname[0] == 'A')
      parmap.emplace(parname, parvalue);
  }

  int bkgdorder =
      static_cast<int>(parmap.size() - 1); // A0 - A(n) total n+1 parameters
  bkgdfunc->setAttributeValue("n", bkgdorder);
  for (auto &mit : parmap) {
    string parname = mit.first;
    double parvalue = mit.second;
    bkgdfunc->setParameter(parname, parvalue);
  }

  // Filter out
  m_outputWS = filterForBackground(bkgdfunc);

  return;
}

  /** Executes the algorithm
     *
     *  @throw Exception::RuntimeError If ... ...
     */
  void GetDetOffsetsMultiPeaks::importFitWindowTableWorkspace(TableWorkspace_sptr windowtablews)
  {
    // Check number of columns matches number of peaks
    size_t numcols = windowtablews->columnCount();
    size_t numpeaks = m_peakPositions.size();
    if (numcols != 2*numpeaks+1)
      throw std::runtime_error("Number of columns is not 2 times of number of referenced peaks. ");

    // Check number of spectra should be same to input workspace
    size_t numrows = windowtablews->rowCount();
    if (numrows != m_inputWS->getNumberHistograms())
      throw std::runtime_error("Number of spectra in fit window workspace does not match input workspace. ");

    // Create workspace
    m_vecFitWindow.clear();
    m_vecFitWindow.resize(numrows);

    for (size_t i = 0; i < numrows; ++i)
    {
      // spectrum number
      int spec = windowtablews->cell<int>(i, 0);
      if (spec < 0 || spec >= static_cast<int>(numrows))
      {
        std::stringstream ess;
        ess << "Peak fit windows at row " << i << " has spectrum " << spec
            << ", which is out of allowed range! ";
        throw std::runtime_error(ess.str());
      }
      else if (m_vecFitWindow[spec].size() != 0)
      {
        std::stringstream ess;
        ess << "Peak fit windows at row " << i << " has spectrum " << spec
            << ", which appears before in fit window table workspace. ";
        throw std::runtime_error(ess.str());
      }

      // fit windows
      std::vector<double> fitwindows(numcols-1);
      for (size_t j = 1; j < numcols; ++j)
      {
        double dtmp = windowtablews->cell<double>(i, j);
        fitwindows[j-1] = dtmp;
      }

      // add to vector of fit windows
      m_vecFitWindow[spec] = fitwindows;
    }

    return;
  }

bool PoldiPeakCollection::checkColumns(
    const TableWorkspace_sptr &tableWorkspace) {
  if (tableWorkspace->columnCount() != 9) {
    return false;
  }

  std::vector<std::string> shouldNames{"HKL", "d", "delta d", "Q", "delta Q",
                                       "Intensity", "delta Intensity",
                                       "FWHM (rel.)", "delta FWHM (rel.)"};

  std::vector<std::string> columnNames = tableWorkspace->getColumnNames();

  return columnNames == shouldNames;
}

TableWorkspace_sptr PoldiPeakSummary::getSummaryTable(
    const PoldiPeakCollection_sptr &peakCollection) const {
  if (!peakCollection) {
    throw std::invalid_argument(
        "Cannot create summary of a null PoldiPeakCollection.");
  }

  TableWorkspace_sptr peakResultWorkspace = getInitializedResultWorkspace();

  for (size_t i = 0; i < peakCollection->peakCount(); ++i) {
    storePeakSummary(peakResultWorkspace->appendRow(), peakCollection->peak(i));
  }

  return peakResultWorkspace;
}

/**
 * Creates Dead Time Table using all the data between begin and end.
 * @param specToLoad :: vector containing the spectrum numbers to load
 * @param deadTimes :: vector containing the corresponding dead times
 * @return Dead Time Table create using the data
 */
TableWorkspace_sptr
LoadMuonNexus1::createDeadTimeTable(std::vector<int> specToLoad,
                                    std::vector<double> deadTimes) {
  TableWorkspace_sptr deadTimeTable =
      boost::dynamic_pointer_cast<TableWorkspace>(
          WorkspaceFactory::Instance().createTable("TableWorkspace"));

  deadTimeTable->addColumn("int", "spectrum");
  deadTimeTable->addColumn("double", "dead-time");

  for (size_t i = 0; i<specToLoad.size(); i++) {
    TableRow row = deadTimeTable->appendRow();
    row << specToLoad[i] << deadTimes[i];
  }

  return deadTimeTable;
}

void PoldiPeakCollection::recoverDataFromLog(
    const TableWorkspace_sptr &tableWorkspace) {
  LogManager_sptr tableLog = tableWorkspace->logs();

  m_intensityType = intensityTypeFromString(getIntensityTypeFromLog(tableLog));
  m_profileFunctionName = getProfileFunctionNameFromLog(tableLog);
  m_pointGroup = pointGroupFromString(getPointGroupStringFromLog(tableLog));
  m_unitCell = unitCellFromString(getUnitCellStringFromLog(tableLog));
}

void PoldiPeakCollection::dataToTableLog(const TableWorkspace_sptr &table) {
  LogManager_sptr tableLog = table->logs();
  tableLog->addProperty<std::string>("IntensityType",
                                     intensityTypeToString(m_intensityType));
  tableLog->addProperty<std::string>("ProfileFunctionName",
                                     m_profileFunctionName);
  tableLog->addProperty<std::string>("PointGroup",
                                     pointGroupToString(m_pointGroup));
  tableLog->addProperty<std::string>("UnitCell",
                                     Geometry::unitCellToStr(m_unitCell));
}

void PoldiPeakCollection::peaksToTable(const TableWorkspace_sptr &table) {
  for (std::vector<PoldiPeak_sptr>::const_iterator peak = m_peaks.begin();
       peak != m_peaks.end(); ++peak) {
    TableRow newRow = table->appendRow();

    newRow << MillerIndicesIO::toString((*peak)->hkl()) << (*peak)->d().value()
           << (*peak)->d().error() << (*peak)->q().value()
           << (*peak)->q().error() << (*peak)->intensity().value()
           << (*peak)->intensity().error()
           << (*peak)->fwhm(PoldiPeak::Relative).value()
           << (*peak)->fwhm(PoldiPeak::Relative).error();
  }
}

void PoldiPeakCollection::constructFromTableWorkspace(
    const TableWorkspace_sptr &tableWorkspace) {
  if (checkColumns(tableWorkspace)) {
    size_t newPeakCount = tableWorkspace->rowCount();
    m_peaks.resize(newPeakCount);

    recoverDataFromLog(tableWorkspace);

    for (size_t i = 0; i < newPeakCount; ++i) {
      TableRow nextRow = tableWorkspace->getRow(i);
      std::string hklString;
      double d, deltaD, q, deltaQ, intensity, deltaIntensity, fwhm, deltaFwhm;
      nextRow >> hklString >> d >> deltaD >> q >> deltaQ >> intensity >>
          deltaIntensity >> fwhm >> deltaFwhm;

      PoldiPeak_sptr peak = PoldiPeak::create(
          MillerIndicesIO::fromString(hklString), UncertainValue(d, deltaD),
          UncertainValue(intensity, deltaIntensity),
          UncertainValue(fwhm, deltaFwhm));
      m_peaks[i] = peak;
    }
  }
}

/// Create a TableWorkspace for the statistics with appropriate columns or get
/// one from the ADS.
ITableWorkspace_sptr
SortHKL::getStatisticsTable(const std::string &name) const {
  TableWorkspace_sptr tablews;

  // Init or append to a table workspace
  bool append = getProperty("Append");
  if (append && AnalysisDataService::Instance().doesExist(name)) {
    tablews = AnalysisDataService::Instance().retrieveWS<TableWorkspace>(name);
  } else {
    tablews = boost::make_shared<TableWorkspace>();
    tablews->addColumn("str", "Resolution Shell");
    tablews->addColumn("int", "No. of Unique Reflections");
    tablews->addColumn("double", "Resolution Min");
    tablews->addColumn("double", "Resolution Max");
    tablews->addColumn("double", "Multiplicity");
    tablews->addColumn("double", "Mean ((I)/sd(I))");
    tablews->addColumn("double", "Rmerge");
    tablews->addColumn("double", "Rpim");
    tablews->addColumn("double", "Data Completeness");
  }

  return tablews;
}

/** Copy table workspace content from one workspace to another
  * @param sourceWS :: table workspace from which the content is copied;
  * @param targetWS :: table workspace to which the content is copied;
  */
void ExtractMaskToTable::copyTableWorkspaceContent(
    TableWorkspace_sptr sourceWS, TableWorkspace_sptr targetWS) {
  // Compare the column names.  They must be exactly the same
  vector<string> sourcecolnames = sourceWS->getColumnNames();
  vector<string> targetcolnames = targetWS->getColumnNames();
  if (sourcecolnames.size() != targetcolnames.size()) {
    stringstream errmsg;
    errmsg << "Soruce table workspace " << sourceWS->name()
           << " has different number of columns (" << sourcecolnames.size()
           << ") than target table workspace's (" << targetcolnames.size()
           << ")";
    throw runtime_error(errmsg.str());
  }
  for (size_t i = 0; i < sourcecolnames.size(); ++i) {
    if (sourcecolnames[i].compare(targetcolnames[i])) {
      stringstream errss;
      errss << "Source and target have incompatible column name at column " << i
            << ". "
            << "Column name of source is " << sourcecolnames[i] << "; "
            << "Column name of target is " << targetcolnames[i];
      throw runtime_error(errss.str());
    }
  }

  // Copy over the content
  size_t numrows = sourceWS->rowCount();
  for (size_t i = 0; i < numrows; ++i) {
    double xmin, xmax;
    string speclist;
    TableRow tmprow = sourceWS->getRow(i);
    tmprow >> xmin >> xmax >> speclist;

    TableRow newrow = targetWS->appendRow();
    newrow << xmin << xmax << speclist;
  }

  return;
}

int PeakIntegration::fitneighbours(int ipeak, std::string det_name, int x0,
                                   int y0, int idet, double qspan,
                                   PeaksWorkspace_sptr &Peaks,
                                   const detid2index_map &pixel_to_wi) {
  UNUSED_ARG(ipeak);
  UNUSED_ARG(det_name);
  UNUSED_ARG(x0);
  UNUSED_ARG(y0);
  Geometry::IPeak &peak = Peaks->getPeak(ipeak);
  // Number of slices
  int TOFmax = 0;

  IAlgorithm_sptr slice_alg = createChildAlgorithm("IntegratePeakTimeSlices");
  slice_alg->setProperty<MatrixWorkspace_sptr>("InputWorkspace", inputW);
  std::ostringstream tab_str;
  tab_str << "LogTable" << ipeak;

  slice_alg->setPropertyValue("OutputWorkspace", tab_str.str());
  slice_alg->setProperty<PeaksWorkspace_sptr>("Peaks", Peaks);
  slice_alg->setProperty("PeakIndex", ipeak);
  slice_alg->setProperty("PeakQspan", qspan);

  int nPixels = std::max<int>(0, getProperty("NBadEdgePixels"));

  slice_alg->setProperty("NBadEdgePixels", nPixels);
  slice_alg->executeAsChildAlg();
  Mantid::API::MemoryManager::Instance().releaseFreeMemory();

  MantidVec &Xout = outputW->dataX(idet);
  MantidVec &Yout = outputW->dataY(idet);
  MantidVec &Eout = outputW->dataE(idet);
  TableWorkspace_sptr logtable = slice_alg->getProperty("OutputWorkspace");

  peak.setIntensity(slice_alg->getProperty("Intensity"));
  peak.setSigmaIntensity(slice_alg->getProperty("SigmaIntensity"));

  TOFmax = static_cast<int>(logtable->rowCount());
  for (int iTOF = 0; iTOF < TOFmax; iTOF++) {
    Xout[iTOF] = logtable->getRef<double>(std::string("Time"), iTOF);
    if (m_IC) // Ikeda-Carpenter fit
    {
      Yout[iTOF] = logtable->getRef<double>(std::string("TotIntensity"), iTOF);
      Eout[iTOF] =
          logtable->getRef<double>(std::string("TotIntensityError"), iTOF);
    } else {
      Yout[iTOF] = logtable->getRef<double>(std::string("ISAWIntensity"), iTOF);
      Eout[iTOF] =
          logtable->getRef<double>(std::string("ISAWIntensityError"), iTOF);
    }
  }

  outputW->getSpectrum(idet)->clearDetectorIDs();
  // Find the pixel ID at that XY position on the rectangular detector
  int pixelID = peak.getDetectorID(); // det->getAtXY(x0,y0)->getID();

  // Find the corresponding workspace index, if any
  auto wiEntry = pixel_to_wi.find(pixelID);
  if (wiEntry != pixel_to_wi.end()) {
    size_t wi = wiEntry->second;
    // Set detectorIDs
    outputW->getSpectrum(idet)
        ->addDetectorIDs(inputW->getSpectrum(wi)->getDetectorIDs());
  }

  return TOFmax - 1;
}

  /** Parse table workspace to a map of Parameters
    */
  void RefinePowderInstrumentParameters2::parseTableWorkspace(TableWorkspace_sptr tablews,
                                                              map<string, Parameter>& parammap)
  {
    // 1. Process Table column names
    std::vector<std::string> colnames = tablews->getColumnNames();
    map<string, size_t> colnamedict;
    convertToDict(colnames, colnamedict);

    int iname = getStringIndex(colnamedict, "Name");
    int ivalue = getStringIndex(colnamedict, "Value");
    int ifit = getStringIndex(colnamedict, "FitOrTie");
    int imin = getStringIndex(colnamedict, "Min");
    int imax = getStringIndex(colnamedict, "Max");
    int istep = getStringIndex(colnamedict, "StepSize");

    if (iname < 0 || ivalue < 0 || ifit < 0)
      throw runtime_error("TableWorkspace does not have column Name, Value and/or Fit.");

    // 3. Parse
    size_t numrows = tablews->rowCount();
    for (size_t irow = 0; irow < numrows; ++irow)
    {
      string parname = tablews->cell<string>(irow, iname);
      double parvalue = tablews->cell<double>(irow, ivalue);
      string fitq = tablews->cell<string>(irow, ifit);

      double minvalue;
      if (imin >= 0)
        minvalue = tablews->cell<double>(irow, imin);
      else
        minvalue = -DBL_MAX;

      double maxvalue;
      if (imax >= 0)
        maxvalue = tablews->cell<double>(irow, imax);
      else
        maxvalue = DBL_MAX;

      double stepsize;
      if (istep >= 0)
        stepsize = tablews->cell<double>(irow, istep);
      else
        stepsize = 1.0;

      Parameter newpar;
      newpar.name = parname;
      newpar.value = parvalue;
      newpar.minvalue = minvalue;
      newpar.maxvalue = maxvalue;
      newpar.stepsize = stepsize;

      // If empty string, fit is default to be false
      bool fit = false;
      if (fitq.size() > 0)
      {
        if (fitq[0] == 'F' || fitq[0] == 'f')
          fit = true;
      }
      newpar.fit = fit;

      parammap.insert(make_pair(parname, newpar));
    }

    return;
  }

/** Fits each spectrum in the workspace to f(x) = A * sin( w * x + p)
 * @param ws :: [input] The workspace to fit
 * @param freq :: [input] Hint for the frequency (w)
 * @param groupName :: [input] The name of the output workspace group
 * @param resTab :: [output] Table workspace storing the asymmetries and phases
 * @param resGroup :: [output] Workspace group storing the fitting results
 */
void CalMuonDetectorPhases::fitWorkspace(const API::MatrixWorkspace_sptr &ws,
                                         double freq, std::string groupName,
                                         API::ITableWorkspace_sptr resTab,
                                         API::WorkspaceGroup_sptr &resGroup) {

  int nhist = static_cast<int>(ws->getNumberHistograms());

  // Create the fitting function f(x) = A * sin ( w * x + p )
  // The same function and initial parameters are used for each fit
  std::string funcStr = createFittingFunction(freq, true);

  // Set up results table
  resTab->addColumn("int", "Spectrum number");
  resTab->addColumn("double", "Asymmetry");
  resTab->addColumn("double", "Phase");

  const auto &indexInfo = ws->indexInfo();

  // Loop through fitting all spectra individually
  const static std::string success = "success";
  for (int wsIndex = 0; wsIndex < nhist; wsIndex++) {
    reportProgress(wsIndex, nhist);
    const auto &yValues = ws->y(wsIndex);
    auto emptySpectrum = std::all_of(yValues.begin(), yValues.end(),
                                     [](double value) { return value == 0.; });
    if (emptySpectrum) {
      g_log.warning("Spectrum " + std::to_string(wsIndex) + " is empty");
      TableWorkspace_sptr tab = boost::make_shared<TableWorkspace>();
      tab->addColumn("str", "Name");
      tab->addColumn("double", "Value");
      tab->addColumn("double", "Error");
      for (int j = 0; j < 4; j++) {
        API::TableRow row = tab->appendRow();
        if (j == PHASE_ROW) {
          row << "dummy" << 0.0 << 0.0;
        } else {
          row << "dummy" << ASYMM_ERROR << 0.0;
        }
      }

      extractDetectorInfo(*tab, *resTab, indexInfo.spectrumNumber(wsIndex));

    } else {
      auto fit = createChildAlgorithm("Fit");
      fit->initialize();
      fit->setPropertyValue("Function", funcStr);
      fit->setProperty("InputWorkspace", ws);
      fit->setProperty("WorkspaceIndex", wsIndex);
      fit->setProperty("CreateOutput", true);
      fit->setPropertyValue("Output", groupName);
      fit->execute();

      std::string status = fit->getProperty("OutputStatus");
      if (!fit->isExecuted()) {
        std::ostringstream error;
        error << "Fit failed for spectrum at workspace index " << wsIndex;
        error << ": " << status;
        throw std::runtime_error(error.str());
      } else if (status != success) {
        g_log.warning("Fit failed for spectrum at workspace index " +
                      std::to_string(wsIndex) + ": " + status);
      }

      API::MatrixWorkspace_sptr fitOut = fit->getProperty("OutputWorkspace");
      resGroup->addWorkspace(fitOut);
      API::ITableWorkspace_sptr tab = fit->getProperty("OutputParameters");
      // Now we have our fitting results stored in tab
      // but we need to extract the relevant information, i.e.
      // the detector phases (parameter 'p') and asymmetries ('A')
      extractDetectorInfo(*tab, *resTab, indexInfo.spectrumNumber(wsIndex));
    }
  }
}

/** 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';
//.........这里部分代码省略.........

/** Add a list of spectra (detector IDs) to the output table workspace.
  * If a detector is masked in input MaskTableWorkspace, then it will not be
 * added to a new row
  * @param outws :: table workspace to write
  * @param maskeddetids :: vector of detector IDs of which detectors masked
  * @param xmin :: minumim x
  * @param xmax :: maximum x
  * @param prevmaskedids :: vector of previous masked detector IDs
  */
void ExtractMaskToTable::addToTableWorkspace(TableWorkspace_sptr outws,
                                             vector<detid_t> maskeddetids,
                                             double xmin, double xmax,
                                             vector<detid_t> prevmaskedids) {
  // Sort vector of detectors ID
  size_t numdetids = maskeddetids.size();
  if (numdetids == 0) {
    stringstream warnss;
    warnss << "Attempting to add an empty vector of masked detectors IDs to "
              "output workspace.  Operation failed.";
    g_log.warning(warnss.str());
    return;
  } else {
    sort(maskeddetids.begin(), maskeddetids.end());
  }

  // Exclude previously masked detectors IDs from masked detectors IDs
  if (prevmaskedids.size() > 0) {
    sort(prevmaskedids.begin(), prevmaskedids.end());
    maskeddetids = subtractVector(maskeddetids, prevmaskedids);
    numdetids = maskeddetids.size();
  } else {
    g_log.debug() << "[DB] There is no previously masked detectors."
                  << ".\n";
  }

  if (numdetids == 0) {
    // I don't know what should be done here
    throw std::runtime_error("Empty detector ID list");
  }

  // Convert vector to string
  stringstream spectralist;
  detid_t previd = maskeddetids[0];
  detid_t headid = maskeddetids[0];
  for (size_t i = 1; i < numdetids; ++i) {
    detid_t tmpid = maskeddetids[i];
    if (tmpid == previd + 1) {
      // Continuous ID
      previd = tmpid;
    } else if (tmpid > previd + 1) {
      // Skipped ID: make a pair
      if (previd == headid) {
        // Single item
        spectralist << " " << headid << ", ";
      } else {
        // Multiple items
        spectralist << " " << headid << "-" << previd << ", ";
      }

      // New head
      headid = tmpid;
      previd = tmpid;
    } else {
      g_log.error() << "Current ID = " << tmpid << ", Previous ID = " << previd
                    << ", Head ID = " << headid << ".\n";
      throw runtime_error("Impossible!  Programming logic error!");
    }
  } // ENDFOR (i)

  // Last one
  if (previd == headid)
    spectralist << " " << headid;
  else
    spectralist << " " << headid << "-" << previd;

  // Add to table workspace
  string specliststr = spectralist.str();
  TableRow newrow = outws->appendRow();
  newrow << xmin << xmax << specliststr;

  return;
}

/** Convert the workspace units using TOF as an intermediate step in the
 * conversion
 * @param fromUnit :: The unit of the input workspace
 * @param outputWS :: The output workspace
 */
void ConvertUnitsUsingDetectorTable::convertViaTOF(
    Kernel::Unit_const_sptr fromUnit, API::MatrixWorkspace_sptr outputWS) {
  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 read them into some vectors.
  auto l1Column = paramWS->getColVector<double>("l1");
  auto l2Column = paramWS->getColVector<double>("l2");
  auto twoThetaColumn = paramWS->getColVector<double>("twotheta");
  auto efixedColumn = paramWS->getColVector<double>("efixed");
  auto emodeColumn = paramWS->getColVector<int>("emode");
  auto spectraColumn = paramWS->getColVector<int>("spectra");

  EventWorkspace_sptr eventWS =
      boost::dynamic_pointer_cast<EventWorkspace>(outputWS);
  assert(static_cast<bool>(eventWS) == m_inputEvents); // Sanity check

  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 = outputWS->getAxis(0)->unit();

  std::vector<double> emptyVec;
  int failedDetectorCount = 0;

  // ConstColumnVector<int> spectraNumber = paramWS->getVector("spectra");

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

    // Lets find what row this spectrum ID 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 specid = det->getID();

      // int spectraNumber = static_cast<int>(spectraColumn->toDouble(i));
      // wsid = outputWS->getIndexFromSpectrumNumber(spectraNumber);
      g_log.debug() << "###### Spectra #" << specid
                    << " ==> Workspace ID:" << wsid << std::endl;

      // Now we need to find the row that contains this spectrum
      std::vector<int>::iterator specIter;

      specIter = std::find(spectraColumn.begin(), spectraColumn.end(), specid);
      if (specIter != spectraColumn.end()) {
        size_t detectorRow = std::distance(spectraColumn.begin(), specIter);
        double l1 = l1Column[detectorRow];
        double l2 = l2Column[detectorRow];
        double twoTheta = twoThetaColumn[detectorRow] * deg2rad;
        double efixed = efixedColumn[detectorRow];
        int emode = emodeColumn[detectorRow];

        g_log.debug() << "specId from detector table = "
                      << spectraColumn[detectorRow] << std::endl;

        // l1 = l1Column->toDouble(detectorRow);
        // l2 = l2Column->toDouble(detectorRow);
        // twoTheta = deg2rad * twoThetaColumn->toDouble(detectorRow);
        // efixed = efixedColumn->toDouble(detectorRow);
        // emode = static_cast<int>(emodeColumn->toDouble(detectorRow));

        g_log.debug() << "###### Spectra #" << specid
                      << " ==> Det Table Row:" << detectorRow << std::endl;
//.........这里部分代码省略.........

/** Executes the algorithm
   *
   *  @throw Exception::RuntimeError If ... ...
   */
void GetDetOffsetsMultiPeaks::importFitWindowTableWorkspace(
    TableWorkspace_sptr windowtablews) {
  // Check number of columns matches number of peaks
  size_t numcols = windowtablews->columnCount();
  size_t numpeaks = m_peakPositions.size();

  if (numcols != 2 * numpeaks + 1)
    throw std::runtime_error(
        "Number of columns is not 2 times of number of referenced peaks. ");

  // Check number of spectra should be same to input workspace
  size_t numrows = windowtablews->rowCount();
  bool needuniversal = false;
  if (numrows < m_inputWS->getNumberHistograms())
    needuniversal = true;
  else if (numrows > m_inputWS->getNumberHistograms())
    throw std::runtime_error(
        "Number of rows in table workspace is larger than number of spectra.");

  // Clear and re-size of the vector for fit windows
  m_vecFitWindow.clear();
  m_vecFitWindow.resize(m_inputWS->getNumberHistograms());

  std::vector<double> vec_univFitWindow;
  bool founduniversal = false;

  // Parse the table workspace
  for (size_t i = 0; i < numrows; ++i) {
    // spectrum number
    int spec = windowtablews->cell<int>(i, 0);
    if (spec >= static_cast<int>(numrows)) {
      std::stringstream ess;
      ess << "Peak fit windows at row " << i << " has spectrum " << spec
          << ", which is out of allowed range! ";
      throw std::runtime_error(ess.str());
    }
    if (spec < 0 && founduniversal) {
      th