/**
* Perform a call to nxgetslab, via the NexusClasses wrapped methods for a given
* block-size
* @param data :: The NXDataSet object
* @param blocksize :: The block-size to use
* @param period :: The period number
* @param start :: The index within the file to start reading from (zero based)
* @param hist :: The workspace index to start reading into
* @param spec_num :: The spectrum number that matches the hist variable
* @param local_workspace :: The workspace to fill the data with
*/
void LoadISISNexus2::loadBlock(NXDataSetTyped<int> &data, int64_t blocksize,
                               int64_t period, int64_t start, int64_t &hist,
                               int64_t &spec_num,
                               DataObjects::Workspace2D_sptr &local_workspace) {
  data.load(static_cast<int>(blocksize), static_cast<int>(period),
            static_cast<int>(start)); // TODO this is just wrong
  int *data_start = data();
  int *data_end = data_start + m_loadBlockInfo.numberOfChannels;
  int64_t final(hist + blocksize);
  while (hist < final) {
    m_progress->report("Loading data");
    MantidVec &Y = local_workspace->dataY(hist);
    Y.assign(data_start, data_end);
    data_start += m_detBlockInfo.numberOfChannels;
    data_end += m_detBlockInfo.numberOfChannels;
    MantidVec &E = local_workspace->dataE(hist);
    std::transform(Y.begin(), Y.end(), E.begin(), dblSqrt);
    // Populate the workspace. Loop starts from 1, hence i-1
    local_workspace->setX(hist, m_tof_data);
    if (m_load_selected_spectra) {
      // local_workspace->getAxis(1)->setValue(hist,
      // static_cast<specid_t>(spec_num));
      auto spec = local_workspace->getSpectrum(hist);
      specid_t specID = m_specInd2specNum_map.at(hist);
      // set detectors corresponding to spectra Number
      spec->setDetectorIDs(m_spec2det_map.getDetectorIDsForSpectrumNo(specID));
      // set correct spectra Number
      spec->setSpectrumNo(specID);
    }

    ++hist;
    ++spec_num;
  }
}

    /**  Load logs from Nexus file. Logs are expected to be in
    *   /raw_data_1/runlog group of the file. Call to this method must be done
    *   within /raw_data_1 group.
    *   @param ws :: The workspace to load the logs to.
    *   @param period :: The period of this workspace
    */
    void LoadISISNexus2::loadLogs(DataObjects::Workspace2D_sptr ws, int period)
    {
      IAlgorithm_sptr alg = createSubAlgorithm("LoadNexusLogs", 0.0, 0.5);
      alg->setPropertyValue("Filename", this->getProperty("Filename"));
      alg->setProperty<MatrixWorkspace_sptr>("Workspace", ws);
      try
      {
        alg->executeAsSubAlg();
      }
      catch(std::runtime_error&)
      {
        g_log.warning() << "Unable to load run logs. There will be no log "
                        << "data associated with this workspace\n";
        return;
      }
      ws->populateInstrumentParameters();
      // If we loaded an icp_event log then create the necessary period logs 
      if( ws->run().hasProperty("icp_event") )
      {
        Kernel::Property *log = ws->run().getProperty("icp_event");
        LogParser parser(log);
        ws->mutableRun().addProperty(parser.createPeriodLog(period));
        ws->mutableRun().addProperty(parser.createAllPeriodsLog());
      }

    }

/** Finalizes the calculation of the correlation spectrum
  *
  * This method offers a variable way of using the correlation spectrum
  *calculated previously.
  * The base version converts to Q and creates an appropriate output workspace.
  *
  * @param correctedCorrelatedIntensities :: Intensities of correlation
  *spectrum.
  * @param dValues :: d-spacings at which the spectrum was calculated.
  * @return A workspace containing the correlation spectrum.
  */
DataObjects::Workspace2D_sptr PoldiAutoCorrelationCore::finalizeCalculation(
    const std::vector<double> &correctedCorrelatedIntensities,
    const std::vector<double> &dValues) const {
  /* Finally, the d-Values are converted to q-Values for plotting etc. and
   * inserted into the output workspace. */
  size_t dCount = dValues.size();
  std::vector<double> qValues(dCount);

  PARALLEL_FOR_NO_WSP_CHECK()
  for (int i = 0; i < static_cast<int>(dCount); ++i) {
    qValues[dCount - i - 1] = Conversions::dToQ(dValues[i]);
  }

  m_logger.information() << "  Setting result..." << std::endl;
  DataObjects::Workspace2D_sptr outputWorkspace =
      boost::dynamic_pointer_cast<Mantid::DataObjects::Workspace2D>(
          WorkspaceFactory::Instance().create("Workspace2D", 1, dValues.size(),
                                              dValues.size()));

  outputWorkspace->getAxis(0)->setUnit("MomentumTransfer");

  outputWorkspace->dataY(0) = correctedCorrelatedIntensities;

  outputWorkspace->setX(0, qValues);

  return outputWorkspace;
}

void
SaveNXTomo::writeIntensityValue(const DataObjects::Workspace2D_sptr workspace,
                                ::NeXus::File &nxFile, int thisFileInd) {
  // Add Intensity to control if present, use 1 if not
  try {
    nxFile.openPath("/entry1/tomo_entry/control");
  } catch (...) {
    throw std::runtime_error("Unable to create a valid NXTomo file");
  }

  std::vector<double> intensityValue;
  intensityValue.push_back(1);

  if (workspace->run().hasProperty("Intensity")) {
    std::string tmpVal = workspace->run().getLogData("Intensity")->value();
    try {
      intensityValue[0] = boost::lexical_cast<double>(tmpVal);
    } catch (...) {
    }
    // Invalid Cast is handled below
  }

  nxFile.openData("data");
  nxFile.putSlab(intensityValue, thisFileInd, 1);
  nxFile.closeData();
}

/** Load in a single spectrum taken from a NeXus file
*  @param hist ::     The workspace index
*  @param i ::        The spectrum index
*  @param specNo ::   The spectrum number
*  @param nxload ::   A reference to the MuonNeXusReader object
*  @param lengthIn :: The number of elements in a spectrum
*  @param localWorkspace :: A pointer to the workspace in which the data will be
* stored
*/
void LoadMuonNexus1::loadData(size_t hist, specid_t &i, specid_t specNo, MuonNexusReader &nxload,
                              const int64_t lengthIn,
                              DataObjects::Workspace2D_sptr localWorkspace) {
  // Read in a spectrum
  // Put it into a vector, discarding the 1st entry, which is rubbish
  // But note that the last (overflow) bin is kept
  // For Nexus, not sure if above is the case, hence give all data for now
  MantidVec &Y = localWorkspace->dataY(hist);
  Y.assign(nxload.counts + i * lengthIn,
           nxload.counts + i * lengthIn + lengthIn);

  // Create and fill another vector for the errors, containing sqrt(count)
  MantidVec &E = localWorkspace->dataE(hist);
  typedef double (*uf)(double);
  uf dblSqrt = std::sqrt;
  std::transform(Y.begin(), Y.end(), E.begin(), dblSqrt);
  // Populate the workspace. Loop starts from 1, hence i-1

  // Create and fill another vector for the X axis  
  float *timeChannels = new float[lengthIn+1]();
  nxload.getTimeChannels(timeChannels, static_cast<const int>(lengthIn+1));
  // Put the read in array into a vector (inside a shared pointer)
  boost::shared_ptr<MantidVec> timeChannelsVec(
    new MantidVec(timeChannels, timeChannels + lengthIn+1));

  localWorkspace->setX(hist, timeChannelsVec);
  localWorkspace->getSpectrum(hist)->setSpectrumNo(specNo);

  // Clean up
  delete[] timeChannels;

}

/**
 * load vectors onto a Workspace2D with 3 bins (the three components of the
 * vectors)
 * dataX for the origin of the vector (assumed (0,0,0) )
 * dataY for the tip of the vector
 * dataE is assumed (0,0,0), no errors
 * @param h5file file identifier
 * @param gws pointer to WorkspaceGroup being filled
 * @param sorting_indexes permutation of qvmod indexes to render it in
 * increasing order of momemtum transfer
 */
const MantidVec LoadSassena::loadQvectors(const hid_t &h5file,
                                          API::WorkspaceGroup_sptr gws,
                                          std::vector<int> &sorting_indexes) {

  const std::string gwsName = this->getPropertyValue("OutputWorkspace");
  const std::string setName("qvectors");

  hsize_t dims[3];
  if (dataSetInfo(h5file, setName, dims) < 0) {
    throw Kernel::Exception::FileError(
        "Unable to read " + setName + " dataset info:", m_filename);
  }
  int nq = static_cast<int>(dims[0]); // number of q-vectors
  double *buf = new double[nq * 3];
  this->dataSetDouble(h5file, "qvectors", buf);

  MantidVec qvmod; // store the modulus of the vector
  double *curr = buf;
  for (int iq = 0; iq < nq; iq++) {
    qvmod.push_back(
        sqrt(curr[0] * curr[0] + curr[1] * curr[1] + curr[2] * curr[2]));
    curr += 3;
  }

  if (getProperty("SortByQVectors")) {
    std::vector<mypair> qvmodpair;
    for (int iq = 0; iq < nq; iq++)
      qvmodpair.push_back(mypair(qvmod[iq], iq));
    std::sort(qvmodpair.begin(), qvmodpair.end(), compare);
    for (int iq = 0; iq < nq; iq++)
      sorting_indexes.push_back(qvmodpair[iq].second);
    std::sort(qvmod.begin(), qvmod.end());
  } else
    for (int iq = 0; iq < nq; iq++)
      sorting_indexes.push_back(iq);

  DataObjects::Workspace2D_sptr ws =
      boost::dynamic_pointer_cast<DataObjects::Workspace2D>(
          API::WorkspaceFactory::Instance().create("Workspace2D", nq, 3, 3));
  std::string wsName = gwsName + std::string("_") + setName;
  ws->setTitle(wsName);

  for (int iq = 0; iq < nq; iq++) {
    MantidVec &Y = ws->dataY(iq);
    const int index = sorting_indexes[iq];
    curr = buf + 3 * index;
    Y.assign(curr, curr + 3);
  }
  delete[] buf;

  ws->getAxis(0)->unit() = Kernel::UnitFactory::Instance().create(
      "MomentumTransfer"); // Set the Units

  this->registerWorkspace(
      gws, wsName, ws, "X-axis: origin of Q-vectors; Y-axis: tip of Q-vectors");
  return qvmod;
}

    /// Run the Child Algorithm LoadInstrument (or LoadInstrumentFromNexus)
    void LoadISISNexus2::runLoadInstrument(DataObjects::Workspace2D_sptr localWorkspace)
    {

      IAlgorithm_sptr loadInst = createChildAlgorithm("LoadInstrument");

      // Now execute the Child Algorithm. Catch and log any error, but don't stop.
      bool executionSuccessful(true);
      try
      {
        loadInst->setPropertyValue("InstrumentName", m_instrument_name);
        loadInst->setProperty<MatrixWorkspace_sptr> ("Workspace", localWorkspace);
        loadInst->setProperty("RewriteSpectraMap", false);
        loadInst->execute();
      }
      catch( std::invalid_argument&)
      {
        g_log.information("Invalid argument to LoadInstrument Child Algorithm");
        executionSuccessful = false;
      }
      catch (std::runtime_error&)
      {
        g_log.information("Unable to successfully run LoadInstrument Child Algorithm");
        executionSuccessful = false;
      }
      if( executionSuccessful )
      {
        // If requested update the instrument to positions in the data file
        const Geometry::ParameterMap & pmap = localWorkspace->instrumentParameters();
        if( pmap.contains(localWorkspace->getInstrument()->getComponentID(),"det-pos-source") )
        {
          boost::shared_ptr<Geometry::Parameter> updateDets = pmap.get(localWorkspace->getInstrument()->getComponentID(),"det-pos-source");
          std::string value = updateDets->value<std::string>();
          if(value.substr(0,8)  == "datafile" )
          {
            IAlgorithm_sptr updateInst = createChildAlgorithm("UpdateInstrumentFromFile");
            updateInst->setProperty<MatrixWorkspace_sptr>("Workspace", localWorkspace);
            updateInst->setPropertyValue("Filename", m_filename);
            if(value  == "datafile-ignore-phi" )
            {
              updateInst->setProperty("IgnorePhi", true);
              g_log.information("Detector positions in IDF updated with positions in the data file except for the phi values");
            }
            else 
            {
              g_log.information("Detector positions in IDF updated with positions in the data file");
            }
            // We want this to throw if it fails to warn the user that the information is not correct.
            updateInst->execute();
          }
        }
      }

    }

/**
 * Convenience function to store a detector value into a given spectrum.
 * Note that this type of data doesn't use TOD, so that we use a single dummy
 * bin in X. Each detector is defined as a spectrum of length 1.
 * @param ws: workspace
 * @param specID: ID of the spectrum to store the value in
 * @param value: value to store [count]
 * @param error: error on the value [count]
 * @param wavelength: wavelength value [Angstrom]
 * @param dwavelength: error on the wavelength [Angstrom]
 */
void store_value(DataObjects::Workspace2D_sptr ws, int specID, double value,
                 double error, double wavelength, double dwavelength) {
  auto &X = ws->mutableX(specID);
  auto &Y = ws->mutableY(specID);
  auto &E = ws->mutableE(specID);
  // The following is mostly to make Mantid happy by defining a histogram with
  // a single bin around the neutron wavelength
  X[0] = wavelength - dwavelength / 2.0;
  X[1] = wavelength + dwavelength / 2.0;
  Y[0] = value;
  E[0] = error;
  ws->getSpectrum(specID).setSpectrumNo(specID);
}

/** Executes the algorithm. Reading in the file and creating and populating
 *  the output workspace
 *
 *  @throw Exception::NotFoundError Error when saving the PoldiDeadWires Results data to Workspace
 *  @throw std::runtime_error Error when saving the PoldiDeadWires Results data to Workspace
 */
void PoldiRemoveDeadWires::exec()
{

	////////////////////////////////////////////////////////////////////////
	// About the workspace
	////////////////////////////////////////////////////////////////////////

	DataObjects::Workspace2D_sptr localWorkspace = this->getProperty("InputWorkspace");

	this->m_channelsPerSpectrum		= localWorkspace.get()->blocksize();
	this->m_numberOfSpectra 		= localWorkspace.get()->size() / m_channelsPerSpectrum;

	g_log.debug() << "_poldi : m_numberOfSpectra     = " << m_numberOfSpectra     << std::endl;
	g_log.debug() << "_poldi : m_channelsPerSpectrum = " << m_channelsPerSpectrum << std::endl;

	////////////////////////////////////////////////////////////////////////
	// Load the data into the workspace
	////////////////////////////////////////////////////////////////////////


    //create table workspace
     try
     {
    	 ITableWorkspace_sptr outputws = WorkspaceFactory::Instance().createTable();

    	// remove the dead-declared wires
    	bool doRemoveExcludedWires = getProperty("RemoveExcludedWires");
    	if(doRemoveExcludedWires){
    		runExcludWires3(localWorkspace, outputws);
    	}


    	// remove the auto-detected dead wires
    	bool doAutoRemoveBadWires = getProperty("AutoRemoveBadWires");
    	if(doAutoRemoveBadWires){
    		autoRemoveDeadWires(localWorkspace, outputws);
    	}

    	setProperty("PoldiDeadWires",outputws);

     }
     catch(Mantid::Kernel::Exception::NotFoundError& )
     {
       throw std::runtime_error("Error when saving the PoldiDeadWires Results data to Workspace");
     }
     catch(std::runtime_error &)
     {
       throw std::runtime_error("Error when saving the PoldiDeadWires Results data to Workspace");
     }

}

 /**
  * Convenience function to store a detector value into a given spectrum.
  * Note that this type of data doesn't use TOD, so that we use a single dummy
  * bin in X. Each detector is defined as a spectrum of length 1.
  * @param ws: workspace
  * @param specID: ID of the spectrum to store the value in
  * @param value: value to store [count]
  * @param error: error on the value [count]
  * @param wavelength: wavelength value [Angstrom]
  * @param dwavelength: error on the wavelength [Angstrom]
  */
 void store_value(DataObjects::Workspace2D_sptr ws, int specID,
     double value, double error, double wavelength, double dwavelength)
 {
   MantidVec& X = ws->dataX(specID);
   MantidVec& Y = ws->dataY(specID);
   MantidVec& E = ws->dataE(specID);
   // The following is mostly to make Mantid happy by defining a histogram with
   // a single bin around the neutron wavelength
   X[0] = wavelength-dwavelength/2.0;
   X[1] = wavelength+dwavelength/2.0;
   Y[0] = value;
   E[0] = error;
   ws->getAxis(1)->setValue(specID, specID);
 }

/**  Log the run details from the file
* @param localWorkspace :: The workspace details to use
*/
void
LoadMuonNexus1::loadRunDetails(DataObjects::Workspace2D_sptr localWorkspace) {
  API::Run &runDetails = localWorkspace->mutableRun();

  runDetails.addProperty("run_title", localWorkspace->getTitle(), true);

  int numSpectra = static_cast<int>(localWorkspace->getNumberHistograms());
  runDetails.addProperty("nspectra", numSpectra);

  NXRoot root(m_filename);
  try {
    std::string start_time = root.getString("run/start_time");
    runDetails.addProperty("run_start", start_time);
  } catch (std::runtime_error &) {
    g_log.warning("run/start_time is not available, run_start log not added.");
  }

  try {
    std::string stop_time = root.getString("run/stop_time");
    runDetails.addProperty("run_end", stop_time);
  } catch (std::runtime_error &) {
    g_log.warning("run/stop_time is not available, run_end log not added.");
  }

  try {
    std::string dur = root.getString("run/duration");
    runDetails.addProperty("dur", dur);
    runDetails.addProperty("durunits", 1); // 1 means second here
    runDetails.addProperty("dur_secs", dur);
  } catch (std::runtime_error &) {
    g_log.warning("run/duration is not available, dur log not added.");
  }

  // Get sample parameters
  NXEntry runSample = root.openEntry("run/sample");

  if (runSample.containsDataSet("temperature")) {
    float temperature = runSample.getFloat("temperature");
    runDetails.addProperty("sample_temp", static_cast<double>(temperature));
  }

  if (runSample.containsDataSet("magnetic_field")) {
    float magn_field = runSample.getFloat("magnetic_field");
    runDetails.addProperty("sample_magn_field",
                           static_cast<double>(magn_field));
  }



}

/// Run the LoadLog Child Algorithm
void LoadMuonNexus1::runLoadLog(DataObjects::Workspace2D_sptr localWorkspace) {
  IAlgorithm_sptr loadLog = createChildAlgorithm("LoadMuonLog");
  // Pass through the same input filename
  loadLog->setPropertyValue("Filename", m_filename);
  // Set the workspace property to be the same one filled above
  loadLog->setProperty<MatrixWorkspace_sptr>("Workspace", localWorkspace);

  // Now execute the Child Algorithm. Catch and log any error, but don't stop.
  try {
    loadLog->execute();
  } catch (std::runtime_error &) {
    g_log.error("Unable to successfully run LoadMuonLog Child Algorithm");
  } catch (std::logic_error &) {
    g_log.error("Unable to successfully run LoadMuonLog Child Algorithm");
  }

  if (!loadLog->isExecuted())
    g_log.error("Unable to successfully run LoadMuonLog Child Algorithm");

  NXRoot root(m_filename);

  // Get main field direction
  std::string mainFieldDirection = "Longitudinal"; // default
  try {
    NXChar orientation = root.openNXChar("run/instrument/detector/orientation");
    // some files have no data there
    orientation.load();

    if (orientation[0] == 't') {
      auto p =
          Kernel::make_unique<Kernel::TimeSeriesProperty<double>>("fromNexus");
      std::string start_time = root.getString("run/start_time");
      p->addValue(start_time, -90.0);
      localWorkspace->mutableRun().addLogData(std::move(p));
      mainFieldDirection = "Transverse";
    }
  } catch (...) {
    // no data - assume main field was longitudinal
  }

  // set output property and add to workspace logs
  auto &run = localWorkspace->mutableRun();
  setProperty("MainFieldDirection", mainFieldDirection);
  run.addProperty("main_field_direction", mainFieldDirection);

  ISISRunLogs runLogs(run);
  runLogs.addStatusLog(run);
}

/**
   Read from the instrument file the dead wires and store the information in a TableWorkspace.
   If asked, the dead wires are removed from the data set.

   @param localWorkspace :: input raw data workspace, containing the information about the instrument
   @param outputws :: input dead wire liste workspace
  */
void PoldiRemoveDeadWires::runExcludWires3
(
    DataObjects::Workspace2D_sptr &localWorkspace,
    API::ITableWorkspace_sptr &outputws
)
{
    outputws->addColumn("int","DeadWires");

    boost::shared_ptr<const Mantid::Geometry::IComponent> comp = localWorkspace->getInstrument()->getComponentByName("holder");
    boost::shared_ptr<const Mantid::Geometry::ICompAssembly> bank = boost::dynamic_pointer_cast<const Mantid::Geometry::ICompAssembly>(comp);
    if (bank)
    {
        // Get a vector of children (recursively)
        std::vector<boost::shared_ptr<const Mantid::Geometry::IComponent> > children;
        bank->getChildren(children, true);

        std::vector<double> defaultDeadWires;
        int ewLine = 0;

        for (unsigned int it = 0; it < children.size(); ++it)
        {
            string wireName = children.at(it)->getName();
            std::vector<boost::shared_ptr<const Mantid::Geometry::IComponent> > tyty =
                localWorkspace.get()->getInstrument().get()->getAllComponentsWithName(wireName);

            std::vector<double> tempWire = tyty[0]->getNumberParameter("excluded");
            if(tempWire.size()>0) {
                int val = (int)tempWire[0];
                g_log.debug() << "_poldi : dead wires :" << val <<  std::endl;
                defaultDeadWires.push_back(val);
                for(unsigned int j=0; j<m_channelsPerSpectrum; j++) {
                    localWorkspace->maskBin(val-1,j,1);
                }
                ewLine++;


                TableRow t = outputws->appendRow();
                t << val ;

            }
        }
        g_log.information() << "_poldi : dead wires set to 0 (nb:" << ewLine << ")" <<  std::endl;
        setProperty("nbExcludedWires",ewLine);
    } else {
        g_log.information() << "_poldi : no dead wire removed" <<  std::endl;
    }
}

    /**
     * Populate spectra mapping to detector IDs
     *
     * TODO: Get the detector size information from the workspace directly
     *
     * @param localWorkspace: Workspace2D object
     * @param nxbins: number of bins in X
     * @param nybins: number of bins in Y
     */
    void LoadSpice2D::runLoadMappingTable(DataObjects::Workspace2D_sptr localWorkspace, int nxbins, int nybins)
    {
      // Get the number of monitor channels
      boost::shared_ptr<const Geometry::Instrument> instrument = localWorkspace->getInstrument();
      std::vector<detid_t> monitors = instrument->getMonitors();
      const int nMonitors = static_cast<int>(monitors.size());

      // Number of monitors should be consistent with data file format
      if( nMonitors != LoadSpice2D::nMonitors ) {
        std::stringstream error;
        error << "Geometry error for " << instrument->getName() <<
            ": Spice data format defines " << LoadSpice2D::nMonitors << " monitors, " << nMonitors << " were/was found";
        throw std::runtime_error(error.str());
      }

      const size_t ndet = nxbins*nybins + nMonitors;
      boost::shared_array<detid_t> udet(new detid_t[ndet]);
      boost::shared_array<specid_t> spec(new specid_t[ndet]);

      // Generate mapping of detector/channel IDs to spectrum ID

      // Detector/channel counter
      int icount = 0;

      // Monitor: IDs start at 1 and increment by 1
      for(int i=0; i<nMonitors; i++)
      {
        spec[icount] = icount;
        udet[icount] = icount+1;
        icount++;
      }

      // Detector pixels
      for(int ix=0; ix<nxbins; ix++)
      {
        for(int iy=0; iy<nybins; iy++)
        {
          spec[icount] = icount;
          udet[icount] = 1000000 + iy*1000 + ix;
          icount++;
        }
      }

      // Populate the Spectra Map with parameters
      localWorkspace->replaceSpectraMap(new API::SpectraDetectorMap(spec.get(), udet.get(), ndet));
    }

/**
 * Create workspace to store the structure factor.
 * First spectrum is the real part, second spectrum is the imaginary part
 * X values are the modulus of the Q-vectors
 * @param h5file file identifier
 * @param gws pointer to WorkspaceGroup being filled
 * @param setName string name of dataset
 * @param qvmod vector of Q-vectors' moduli
 * @param sorting_indexes permutation of qvmod indexes to render it in increasing order of momemtum transfer
 */
void LoadSassena::loadFQ(const hid_t& h5file, API::WorkspaceGroup_sptr gws, const std::string setName, const MantidVec &qvmod, const std::vector<int> &sorting_indexes)
{
  const std::string gwsName = this->getPropertyValue("OutputWorkspace");
  int nq = static_cast<int>( qvmod.size() ); //number of q-vectors

  DataObjects::Workspace2D_sptr ws = boost::dynamic_pointer_cast<DataObjects::Workspace2D>(API::WorkspaceFactory::Instance().create("Workspace2D", 2, nq, nq));
  const std::string wsName = gwsName + std::string("_") + setName;
  ws->setTitle(wsName);

  double* buf = new double[nq*2];
  this->dataSetDouble(h5file,setName,buf);
  MantidVec& re = ws->dataY(0); // store the real part
  ws->dataX(0) = qvmod;  //X-axis values are the modulus of the q vector
  MantidVec& im = ws->dataY(1); // store the imaginary part
  ws->dataX(1) = qvmod;
  double *curr = buf;
  for(int iq=0; iq<nq; iq++){
    const int index=sorting_indexes[iq];
    re[index]=curr[0];
    im[index]=curr[1];
    curr+=2;
  }
  delete[] buf;

  // Set the Units
  ws->getAxis(0)->unit() = Kernel::UnitFactory::Instance().create("MomentumTransfer");

  this->registerWorkspace(gws,wsName,ws, "X-axis: Q-vector modulus; Y-axis: intermediate structure factor");
}

    /** Read in a single spectrum from the raw file
     *  @param tcbs ::     The vector containing the time bin boundaries
     *  @param hist ::     The workspace index
     *  @param i ::        The spectrum number
     *  @param iraw ::     A reference to the ISISRAW object
     *  @param lengthIn :: The number of elements in a spectrum
     *  @param spectrum :: Pointer to the array into which the spectrum will be read
     *  @param localWorkspace :: A pointer to the workspace in which the data will be stored
     */
    void LoadRaw::loadData(const MantidVecPtr::ptr_type& tcbs, int32_t hist, specid_t& i, ISISRAW& iraw, const int& lengthIn, int* spectrum, DataObjects::Workspace2D_sptr localWorkspace)
    {
      // Read in a spectrum
      memcpy(spectrum, iraw.dat1 + i * lengthIn, lengthIn * sizeof(int));
      // Put it into a vector, discarding the 1st entry, which is rubbish
      // But note that the last (overflow) bin is kept
      MantidVec& Y = localWorkspace->dataY(hist);
      Y.assign(spectrum + 1, spectrum + lengthIn);
      // Create and fill another vector for the errors, containing sqrt(count)
      MantidVec& E = localWorkspace->dataE(hist);
      std::transform(Y.begin(), Y.end(), E.begin(), dblSqrt);
      // Populate the workspace. Loop starts from 1, hence i-1
      localWorkspace->setX(hist, tcbs);

      localWorkspace->getAxis(1)->setValue(hist, i);
      // NOTE: Raw numbers go straight into the workspace
      //     - no account taken of bin widths/units etc.
    }

/// Calculate the relative change in residuals with respect to the supplied
/// total number of counts
double PoldiAnalyseResiduals::relativeCountChange(
    const DataObjects::Workspace2D_sptr &sum, double totalMeasuredCounts) {
  const MantidVec &corrCounts = sum->readY(0);
  double csum = 0.0;
  for (double corrCount : corrCounts) {
    csum += fabs(corrCount);
  }

  return csum / totalMeasuredCounts * 100.0;
}

    /**
    * Load data about the sample
    *   @param local_workspace :: The workspace to load the logs to.
    *   @param entry :: The Nexus entry
    */
    void LoadISISNexus2::loadSampleData(DataObjects::Workspace2D_sptr local_workspace, NXEntry & entry)
    {
      /// Sample geometry
      NXInt spb = entry.openNXInt("isis_vms_compat/SPB");
      // Just load the index we need, not the whole block. The flag is the third value in
      spb.load(1, 2);
      int geom_id = spb[0];
      local_workspace->mutableSample().setGeometryFlag(spb[0]);

      NXFloat rspb = entry.openNXFloat("isis_vms_compat/RSPB");
      // Just load the indices we need, not the whole block. The values start from the 4th onward
      rspb.load(3, 3);
      double thick(rspb[0]), height(rspb[1]), width(rspb[2]);
      local_workspace->mutableSample().setThickness(thick);
      local_workspace->mutableSample().setHeight(height);
      local_workspace->mutableSample().setWidth(width);

      g_log.debug() << "Sample geometry -  ID: " << geom_id << ", thickness: " << thick << ", height: " << height << ", width: " << width << "\n";
    }

/// Calculate the relative change in residuals with respect to the supplied
/// total number of counts
double PoldiAnalyseResiduals::relativeCountChange(
    const DataObjects::Workspace2D_sptr &sum, double totalMeasuredCounts) {
  const MantidVec &corrCounts = sum->readY(0);
  double csum = 0.0;
  for (auto it = corrCounts.begin(); it != corrCounts.end(); ++it) {
    csum += fabs(*it);
  }

  return csum / totalMeasuredCounts * 100.0;
}

/// Adds the specified value to all spectra specified by the given workspace
/// indices.
void PoldiAnalyseResiduals::addValue(
    DataObjects::Workspace2D_sptr &workspace, double value,
    const std::vector<int> &workspaceIndices) const {
  for (size_t i = 0; i < workspaceIndices.size(); ++i) {
    MantidVec &counts = workspace->dataY(workspaceIndices[i]);
    for (size_t j = 0; j < counts.size(); ++j) {
      counts[j] += value;
    }
  }
}