LoadFullprofResolution.cpp

Go to the documentation of this file.

// Mantid Repository : https://github.com/mantidproject/mantid
//
// Copyright © 2018 ISIS Rutherford Appleton Laboratory UKRI,
//   NScD Oak Ridge National Laboratory, European Spallation Source,
//   Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
// SPDX - License - Identifier: GPL - 3.0 +
#include "MantidDataHandling/LoadFullprofResolution.h"
#include "MantidAPI/FileProperty.h"
#include "MantidAPI/InstrumentDataService.h"
#include "MantidAPI/MatrixWorkspace.h"
#include "MantidAPI/TableRow.h"
#include "MantidAPI/WorkspaceGroup.h"
#include "MantidAPI/WorkspaceProperty.h"
#include "MantidGeometry/Instrument.h"
#include "MantidGeometry/Instrument/InstrumentDefinitionParser.h"
#include "MantidKernel/ArrayProperty.h"
 
#include <boost/algorithm/string.hpp>
#include <boost/algorithm/string/finder.hpp>
#include <boost/algorithm/string/iter_find.hpp>
#include <boost/algorithm/string/trim.hpp>
 
#include <Poco/DOM/AutoPtr.h>
#include <Poco/DOM/DOMWriter.h>
#include <Poco/DOM/Element.h>
#include <Poco/XML/XMLWriter.h>
 
#include <fstream>
 
using namespace Mantid;
using namespace Mantid::API;
using namespace Mantid::DataObjects;
using namespace Mantid::Kernel;
using namespace std;
using namespace Poco::XML;
 
using Geometry::Instrument_const_sptr;
using Geometry::Instrument_sptr;
 
namespace Mantid::DataHandling {
 
DECLARE_ALGORITHM(LoadFullprofResolution)
 
std::map<std::string, size_t> LoadFullprofResolution::m_rowNumbers;
 
//----------------------------------------------------------------------------------------------
void LoadFullprofResolution::init() {
  // Input file name
  declareProperty(std::make_unique<FileProperty>("Filename", "", FileProperty::Load, ".irf"),
                  "Path to an Fullprof .irf file to load.");
 
  // Output table workspace
  auto wsprop = std::make_unique<WorkspaceProperty<API::ITableWorkspace>>("OutputTableWorkspace", "", Direction::Output,
                                                                          PropertyMode::Optional);
  declareProperty(std::move(wsprop), "Name of the output TableWorkspace containing "
                                     "profile parameters or bank information. ");
 
  // Use bank numbers as given in file
  declareProperty(std::make_unique<PropertyWithValue<bool>>("UseBankIDsInFile", true, Direction::Input),
                  "Use bank IDs as given in file rather than ordinal number of bank."
                  "If the bank IDs in the file are not unique, it is advised to set this "
                  "to false.");
 
  // Bank to import
  declareProperty(std::make_unique<ArrayProperty<int>>("Banks"), "ID(s) of specified bank(s) to load, "
                                                                 "The IDs are as specified by UseBankIDsInFile."
                                                                 "Default is all banks contained in input .irf file.");
 
  // Workspace to put parameters into. It must be a workspace group with one
  // workpace per bank from the IRF file
  declareProperty(
      std::make_unique<WorkspaceProperty<WorkspaceGroup>>("Workspace", "", Direction::InOut, PropertyMode::Optional),
      "A workspace group with the instrument to which we add the "
      "parameters from the Fullprof .irf file with one workspace "
      "for each bank of the .irf file");
 
  // Workspaces for each bank
  declareProperty(std::make_unique<ArrayProperty<int>>("WorkspacesForBanks"),
                  "For each Fullprof bank,"
                  " the ID of the corresponding workspace in same order as the "
                  "Fullprof banks are specified. "
                  "ID=1 refers to the first workspace in the workspace group, "
                  "ID=2 refers to the second workspace and so on. "
                  "Default is all workspaces in numerical order."
                  "If default banks are specified, they too are taken to be in "
                  "numerical order");
}
 
//----------------------------------------------------------------------------------------------
void LoadFullprofResolution::exec() {
  // Get input
  string datafile = getProperty("Filename");
  const bool useBankIDsInFile = getProperty("UseBankIDsInFile");
  vector<int> outputbankids = getProperty("Banks");
  WorkspaceGroup_sptr wsg = getProperty("Workspace");
  vector<int> outputwsids = getProperty("WorkspacesForBanks");
 
  // Import data
  vector<string> lines;
  loadFile(datafile, lines);
 
  // Get Prof number
  int nProf = getProfNumber(lines);
 
  // Examine bank information
  vector<int> vec_bankinirf;
  map<int, int> bankstartindexmap, bankendindexmap;
  scanBanks(lines, useBankIDsInFile, vec_bankinirf, bankstartindexmap, bankendindexmap);
  if (useBankIDsInFile)
    sort(vec_bankinirf.begin(), vec_bankinirf.end());
 
  for (auto bank : vec_bankinirf)
    g_log.debug() << "Irf containing bank " << bank << ".\n";
 
  // Bank-workspace correspondence
  map<int, size_t> workspaceOfBank;
 
  vector<int> vec_bankids; // bank IDs to output to table workspace
 
  if (vec_bankinirf.empty()) {
    throw runtime_error("No Bank is found in input file.");
  } else if (outputbankids.empty()) {
    vec_bankids = vec_bankinirf;
    // If workspaces, set up Bank-Workpace correspondence
    if (wsg)
      createBankToWorkspaceMap(vec_bankids, outputwsids, workspaceOfBank);
  } else {
    // If workspaces, set up Bank-Workpace correspondece, before banks are
    // sorted.
    if (wsg)
      createBankToWorkspaceMap(outputbankids, outputwsids, workspaceOfBank);
 
    // Deal with banks
    sort(outputbankids.begin(), outputbankids.end());
    for (auto outputbankid : outputbankids) {
      if (outputbankid < 0) {
        g_log.warning() << "Input bank ID (" << outputbankid << ") is negative.  It is not allowed and is  ignored. "
                        << ".\n";
      } else {
        auto fiter = lower_bound(vec_bankinirf.begin(), vec_bankinirf.end(), outputbankid);
        if (fiter == vec_bankinirf.end() || *fiter != outputbankid) {
          // Specified bank ID does not exist.
          stringstream errmsg;
          errmsg << "Specified output bank (ID = " << outputbankid << ") cannot be found in input " << datafile
                 << ", which includes bank : ";
          for (size_t i = 0; i < vec_bankinirf.size(); ++i) {
            errmsg << vec_bankinirf[i];
            if (i != vec_bankinirf.size() - 1)
              errmsg << ",";
          }
          g_log.error(errmsg.str());
          throw runtime_error(errmsg.str());
        } else {
          vec_bankids.emplace_back(outputbankid);
        }
      }
    }
    if (vec_bankids.empty()) {
      g_log.error("There is no valid specified bank IDs for output.");
      throw runtime_error("There is no valid specified bank IDs for output.");
    }
  }
 
  // Parse .irf and export profile parameters
  map<int, map<string, double>> bankparammap;
  for (int bankid : vec_bankids) {
    g_log.debug() << "Parse bank " << bankid << " of total " << vec_bankids.size() << ".\n";
    map<string, double> parammap;
    parseResolutionStrings(parammap, lines, useBankIDsInFile, bankid, bankstartindexmap[bankid],
                           bankendindexmap[bankid], nProf);
    bankparammap.emplace(bankid, parammap);
  }
 
  // Generate output table workspace
  API::ITableWorkspace_sptr outTabWs = genTableWorkspace(bankparammap);
 
  if (!getPropertyValue("OutputTableWorkspace").empty()) {
    // Output the output table workspace
    setProperty("OutputTableWorkspace", outTabWs);
  }
 
  // If workspace, put parameters there
  if (wsg) {
    // First check that number of workspaces in group matches number of banks,
    // if no WorkspacesForBanks is specified.
    if ((outputwsids.empty()) && (wsg->size() != vec_bankids.size())) {
      std::ostringstream mess;
      mess << "Number of banks (" << vec_bankids.size() << ") does not match number of workspaces (" << wsg->size()
           << ") in group. Parameters not put into workspaces.";
      g_log.error(mess.str());
    } else // Numbers match, so put parameters into workspaces.
    {
      getTableRowNumbers(outTabWs, LoadFullprofResolution::m_rowNumbers);
      for (size_t i = 0; i < vec_bankids.size(); ++i) {
        int bankId = vec_bankids[i];
        size_t wsId = workspaceOfBank[bankId];
        Workspace_sptr wsi = wsg->getItem(wsId - 1);
        auto workspace = std::dynamic_pointer_cast<MatrixWorkspace>(wsi);
        // Get column from table workspace
        API::Column_const_sptr OutTabColumn = outTabWs->getColumn(i + 1);
        std::string parameterXMLString;
        putParametersIntoWorkspace(OutTabColumn, workspace, nProf, parameterXMLString);
 
        // Load the string into the workspace
        Algorithm_sptr loadParamAlg = createChildAlgorithm("LoadParameterFile");
        loadParamAlg->setProperty("ParameterXML", parameterXMLString);
        loadParamAlg->setProperty("Workspace", workspace);
        loadParamAlg->execute();
      }
    }
  } else if (getPropertyValue("OutputTableWorkspace").empty()) {
    // We don't know where to output
    throw std::runtime_error("Either the OutputTableWorkspace or Workspace property must be set.");
  }
}
 
//----------------------------------------------------------------------------------------------
void LoadFullprofResolution::loadFile(const string &filename, vector<string> &lines) {
  ifstream myfile(filename.c_str());
 
  if (myfile.is_open()) {
    while (!myfile.eof()) {
      string line;
      getline(myfile, line);
 
      boost::algorithm::trim(line);
      if (!line.empty())
        lines.emplace_back(line);
    }
 
    myfile.close();
  } else {
    stringstream errmsg;
    errmsg << "Input .irf file " << filename << " cannot be open. ";
    g_log.error(errmsg.str());
    throw runtime_error(errmsg.str());
  }
}
 
//----------------------------------------------------------------------------------------------
int LoadFullprofResolution::getProfNumber(const vector<string> &lines) {
  // Assume the NPROF number is on the second line
  if (lines[1].find("NPROF") != string::npos) {
    // Split line to get the NPROF number
    size_t nStart = lines[1].find("NPROF");
    size_t nEq = lines[1].find('=', nStart);
    size_t nEnd = lines[1].find(' ', nStart); // Assume the NRPOF number is followed by space
    if (nEq == string::npos || nEnd == string::npos)
      return (-1);
    size_t nNumber = nEq + 1;
    return (boost::lexical_cast<int>(lines[1].substr(nNumber, nEnd - nNumber)));
  }
 
  return (0);
}
 
//----------------------------------------------------------------------------------------------
void LoadFullprofResolution::scanBanks(const vector<string> &lines, const bool useFileBankIDs, vector<int> &banks,
                                       map<int, int> &bankstartindexmap, map<int, int> &bankendindexmap) {
  int startindex = -1;
  int endindex;
  int bankid = 0;
  for (size_t i = 0; i < lines.size(); ++i) {
    string line = lines[i];
    if (line.find("Bank") != string::npos) {
      // A new line found
      if (startindex >= 0) {
        // Previous line is in a bank range.  Then finish the previous bank
        // range
        endindex = static_cast<int>(i) - 1;
        bankstartindexmap.emplace(banks.back(), startindex);
        bankendindexmap.emplace(banks.back(), endindex);
      }
 
      // Start the new pair
      startindex = static_cast<int>(i);
 
      // Get bank ID
      if (useFileBankIDs) { // Get bank ID from line
        vector<string> level1s;
        boost::split(level1s, line, boost::is_any_of("Bank"));
        vector<string> level2s;
        string bankterm = level1s.back();
        boost::algorithm::trim(bankterm);
        boost::split(level2s, bankterm, boost::is_any_of(" "));
        bankid = std::stoi(level2s[0]);
      } else { // Get bank ID as ordinal number of bank
        bankid++;
      }
      banks.emplace_back(bankid);
    }
  }
  if (startindex >= 0) {
    endindex = static_cast<int>(lines.size()) - 1;
    bankstartindexmap.emplace(banks.back(), startindex);
    bankendindexmap.emplace(banks.back(), endindex);
  }
 
  g_log.debug() << "[DB1112] Number of bank IDs = " << banks.size() << ", "
                << "Number of ranges = " << bankstartindexmap.size() << '\n';
  for (const auto &bank : banks) {
    g_log.debug() << "Bank " << bank << " From line " << bankstartindexmap[bank] << " to " << bankendindexmap[bank]
                  << '\n';
  }
}
 
//----------------------------------------------------------------------------------------------
void LoadFullprofResolution::parseResolutionStrings(map<string, double> &parammap, const vector<string> &lines,
                                                    const bool useFileBankIDs, int bankid, int startlineindex,
                                                    int endlineindex, int profNumber) {
  string bankline = lines[startlineindex];
  double cwl;
  int tmpbankid;
  parseBankLine(bankline, cwl, tmpbankid);
  if (tmpbankid != -1) {
    g_log.debug() << "Found CWL = " << cwl << ", Bank ID = " << tmpbankid << "\n";
  } else {
    g_log.warning() << "No CWL found for bank " << bankid;
    tmpbankid = bankid;
  }
  if (useFileBankIDs && (bankid != tmpbankid)) {
    stringstream errss;
    errss << "Input bank ID (" << bankid << ") is not same as the bank ID (" << tmpbankid
          << ") found in the specified region from input. ";
    throw runtime_error(errss.str());
  }
 
  parammap["NPROF"] = profNumber;
  parammap["CWL"] = cwl;
 
  for (int i = startlineindex + 1; i <= endlineindex; ++i) {
    string line = lines[i];
 
    // Skip information line
    if (line[0] == '!')
      continue;
 
    // skip NPROF line, which is processed by getProfNumber
    if (line.find("NPROF") != string::npos)
      continue;
 
    // Parse
    g_log.debug() << "Parse Line " << i << "\t\t" << line << "\n";
 
    if (line.starts_with("TOFRG")) {
      // TOFRG tof-min step tof-max
      vector<string> terms;
      boost::split(terms, line, boost::is_any_of(" "), boost::token_compress_on);
      if (terms.size() != 4) {
        stringstream errmsg;
        errmsg << "Line TOFRG has " << terms.size() << " terms.  Different from 4 terms in definition.";
        g_log.error(errmsg.str());
        throw runtime_error(errmsg.str());
      } else {
        parammap["tof-min"] = parseDoubleValue(terms[1], "tof-min");
        parammap["step"] = parseDoubleValue(terms[2], "step");
        parammap["tof-max"] = parseDoubleValue(terms[3], "tof-max");
      }
    } else if (line.starts_with("D2TOF")) {
      // D2TOF Dtt1 Dtt2 Zero
      vector<string> terms;
      boost::split(terms, line, boost::is_any_of(" "), boost::token_compress_on);
      if (terms.size() != 2 && terms.size() != 4) {
        stringstream errmsg;
        errmsg << "Line D2TOF has " << terms.size() << " terms.  Different from 2/4 terms in definition.";
        g_log.error(errmsg.str());
        throw runtime_error(errmsg.str());
      } else {
        parammap["Dtt1"] = parseDoubleValue(terms[1], "Dtt1");
        if (terms.size() == 4) {
          parammap["Dtt2"] = parseDoubleValue(terms[2], "Dtt2");
          parammap["Zero"] = parseDoubleValue(terms[3], "Zero");
        } else {
          parammap["Dtt2"] = 0.0;
          parammap["Zero"] = 0.0;
        }
      }
    } // "D2TOF"
    else if (line.starts_with("ZD2TOF")) {
      vector<string> terms;
      boost::split(terms, line, boost::is_any_of(" "), boost::token_compress_on);
      if (terms.size() != 3) {
        stringstream errmsg;
        errmsg << "Line ZD2TOF has " << terms.size() << " terms.  Different from 4 terms in definition.";
        g_log.error(errmsg.str());
        throw runtime_error(errmsg.str());
      } else {
        parammap["Zero"] = parseDoubleValue(terms[1], "Zero");
        parammap["Dtt1"] = parseDoubleValue(terms[2], "Dtt1");
        parammap["Dtt2"] = 0.;
      }
    } // "ZD2TOF"
    else if (line.starts_with("D2TOT")) {
      vector<string> terms;
      boost::split(terms, line, boost::is_any_of(" "), boost::token_compress_on);
      if (terms.size() != 6) {
        stringstream errmsg;
        errmsg << "Line TOFRG has " << terms.size() << " terms.  Different from 4 terms in definition.";
        g_log.error(errmsg.str());
        throw runtime_error(errmsg.str());
      } else {
        parammap["Dtt1t"] = parseDoubleValue(terms[1], "Dtt1t");
        parammap["Dtt2t"] = parseDoubleValue(terms[2], "Dtt2t");
        parammap["Tcross"] = parseDoubleValue(terms[3], "Tcross");
        parammap["Width"] = parseDoubleValue(terms[4], "Width");
        parammap["Zerot"] = parseDoubleValue(terms[5], "Zerot");
      }
    } // "D2TOT"
    else if (line.starts_with("ZD2TOT")) {
      vector<string> terms;
      boost::split(terms, line, boost::is_any_of(" "), boost::token_compress_on);
      if (terms.size() != 6) {
        stringstream errmsg;
        errmsg << "Line TOFRG has " << terms.size() << " terms.  Different from 4 terms in definition.";
        g_log.error(errmsg.str());
        throw runtime_error(errmsg.str());
      } else {
        parammap["Zerot"] = parseDoubleValue(terms[1], "Zerot");
        parammap["Dtt1t"] = parseDoubleValue(terms[2], "Dtt1t");
        parammap["Dtt2t"] = parseDoubleValue(terms[3], "Dtt2t");
        parammap["Tcross"] = parseDoubleValue(terms[4], "Tcross");
        parammap["Width"] = parseDoubleValue(terms[5], "Width");
      }
    } // "ZD2TOT"
    else if (line.starts_with("TWOTH")) {
      vector<string> terms;
      boost::split(terms, line, boost::is_any_of(" "), boost::token_compress_on);
      if (terms.size() != 2) {
        stringstream errmsg;
        errmsg << "Line TOFRG has " << terms.size() << " terms.  Different from 4 terms in definition.";
        g_log.error(errmsg.str());
        throw runtime_error(errmsg.str());
      } else {
        parammap["twotheta"] = parseDoubleValue(terms[1], "twotheta");
      }
    } // "TWOTH"
    else if (line.starts_with("SIGMA")) {
      vector<string> terms;
      boost::split(terms, line, boost::is_any_of(" "), boost::token_compress_on);
      if (terms.size() != 4) {
        stringstream errmsg;
        errmsg << "Line TOFRG has " << terms.size() << " terms.  Different from 4 terms in definition.";
        g_log.error(errmsg.str());
        throw runtime_error(errmsg.str());
      } else {
        parammap["Sig2"] = sqrt(parseDoubleValue(terms[1], "Sig2"));
        parammap["Sig1"] = sqrt(parseDoubleValue(terms[2], "Sig1"));
        parammap["Sig0"] = sqrt(parseDoubleValue(terms[3], "Sig0"));
      }
    } // "SIGMA"
    else if (line.starts_with("GAMMA")) {
      vector<string> terms;
      boost::split(terms, line, boost::is_any_of(" "), boost::token_compress_on);
      if (terms.size() != 4) {
        stringstream errmsg;
        errmsg << "Line TOFRG has " << terms.size() << " terms.  Different from 4 terms in definition.";
        g_log.error(errmsg.str());
        throw runtime_error(errmsg.str());
      } else {
        parammap["Gam2"] = parseDoubleValue(terms[1], "Gam2");
        parammap["Gam1"] = parseDoubleValue(terms[2], "Gam1");
        parammap["Gam0"] = parseDoubleValue(terms[3], "Gam0");
      }
    } // "GAMMA"
    else if (line.starts_with("ALFBE")) {
      // ALFBE alph0 beta0 alph1 beta1
      vector<string> terms;
      boost::split(terms, line, boost::is_any_of(" "), boost::token_compress_on);
      if (terms.size() != 5) {
        stringstream errmsg;
        errmsg << "Line ALFBE has " << terms.size() << " terms.  Different from 4 terms in definition.";
        g_log.error(errmsg.str());
        throw runtime_error(errmsg.str());
      } else {
        parammap["Alph0"] = parseDoubleValue(terms[1], "Alph0");
        parammap["Beta0"] = parseDoubleValue(terms[2], "Beta0");
        parammap["Alph1"] = parseDoubleValue(terms[3], "Alph1");
        parammap["Beta1"] = parseDoubleValue(terms[4], "Beta1");
      }
    } // "ALFBE"
    else if (line.starts_with("ALFBT")) {
      vector<string> terms;
      boost::split(terms, line, boost::is_any_of(" "), boost::token_compress_on);
      if (terms.size() != 5) {
        stringstream errmsg;
        errmsg << "Line ALFBT has " << terms.size() << " terms.  Different from 4 terms in definition.";
        g_log.error(errmsg.str());
        throw runtime_error(errmsg.str());
      } else {
        parammap["Alph0t"] = parseDoubleValue(terms[1], "Alph0t");
        parammap["Beta0t"] = parseDoubleValue(terms[2], "Beta0t");
        parammap["Alph1t"] = parseDoubleValue(terms[3], "Alph1t");
        parammap["Beta1t"] = parseDoubleValue(terms[4], "Beta1t");
      }
    } // "ALFBT"
    else if (line.starts_with("END")) {
      // Ignore END line
      g_log.debug("END line of bank.");
    } else {
      g_log.warning() << "Line '" << line << "' is not processed.\n";
    } // END -IF StartWith
 
  } // For-all-line
}
 
//----------------------------------------------------------------------------------------------
void LoadFullprofResolution::parseBankLine(string line, double &cwl, int &bankid) {
  // 1. Split along 'Bank'
  std::vector<std::string> v;
  iter_split(v, line, boost::algorithm::first_finder("Bank"));
 
  // 2. Split the rest around CWL if there is
  string infostr = v[1];
  boost::algorithm::trim(infostr);
 
  cwl = -1.0;
  bankid = -1;
 
  if (infostr.find("CWL") != string::npos) {
    // There is CWL
    v.clear();
    iter_split(v, infostr, boost::algorithm::first_finder("CWL"));
 
    // Bank ID
    bankid = std::stoi(v[0]);
 
    // CWL
    infostr = v[1];
    v.clear();
    boost::split(v, infostr, boost::is_any_of("=A"));
    for (size_t i = 0; i < v.size(); ++i) {
      g_log.debug() << "Last CWL splitted.  Term " << i << ": \t\t"
                    << "'" << v[i] << "'\n";
      string candidate = v[i];
      boost::algorithm::trim(candidate);
      if (!candidate.empty()) {
        cwl = parseDoubleValue(candidate);
        break;
      }
    }
  }
}
 
double LoadFullprofResolution::parseDoubleValue(const std::string &value, const std::string &label) {
  if (!value.empty()) {
    try {
      return std::stod(value);
    } catch (...) {
      std::stringstream msg;
      msg << "Failed to convert \"" << value << "\" to a double";
      if (!label.empty()) {
        msg << " for \"" << label << "\"";
      }
      msg << " using 0. instead";
      g_log.warning(msg.str());
    }
  }
  return 0.;
}
 
//----------------------------------------------------------------------------------------------
TableWorkspace_sptr LoadFullprofResolution::genTableWorkspace(map<int, map<string, double>> bankparammap) {
  g_log.notice() << "Start to generate table workspace ...."
                 << ".\n";
 
  // Retrieve some information
  size_t numbanks = bankparammap.size();
  if (numbanks == 0)
    throw runtime_error("Unable to generate a table from an empty map!");
 
  auto bankmapiter = bankparammap.begin();
  size_t numparams = bankmapiter->second.size();
 
  // vector of all parameter name
  vector<string> vec_parname;
  vector<int> vec_bankids;
 
  map<string, double>::iterator parmapiter;
  for (parmapiter = bankmapiter->second.begin(); parmapiter != bankmapiter->second.end(); ++parmapiter) {
    string parname = parmapiter->first;
    vec_parname.emplace_back(parname);
  }
 
  for (bankmapiter = bankparammap.begin(); bankmapiter != bankparammap.end(); ++bankmapiter) {
    int bankid = bankmapiter->first;
    vec_bankids.emplace_back(bankid);
  }
 
  g_log.debug() << "[DBx240] Number of imported parameters is " << numparams
                << ", Number of banks = " << vec_bankids.size() << "."
                << "\n";
 
  // Create TableWorkspace
  auto tablews = std::make_shared<TableWorkspace>();
 
  // set columns :
  // Any 2 columns cannot have the same name.
  tablews->addColumn("str", "Name");
  for (size_t i = 0; i < numbanks; ++i) {
    stringstream colnamess;
    int bankid = vec_bankids[i];
    colnamess << "Value_" << bankid;
    tablews->addColumn("double", colnamess.str());
  }
 
  g_log.debug() << "Number of column = " << tablews->columnCount() << ".\n";
 
  // add BANK ID row
  TableRow newrow = tablews->appendRow();
  newrow << "BANK";
  for (size_t i = 0; i < numbanks; ++i)
    newrow << static_cast<double>(vec_bankids[i]);
 
  g_log.debug() << "Number of row now = " << tablews->rowCount() << ".\n";
 
  // add profile parameter rows
  for (size_t i = 0; i < numparams; ++i) {
    newrow = tablews->appendRow();
 
    string parname = vec_parname[i];
    newrow << parname;
 
    for (size_t j = 0; j < numbanks; ++j) {
      int bankid = vec_bankids[j];
 
      // Locate map of bank 'bankid'
      map<int, map<string, double>>::iterator bpmapiter;
      bpmapiter = bankparammap.find(bankid);
      if (bpmapiter == bankparammap.end()) {
        throw runtime_error("Bank cannot be found in map.");
      }
 
      // Locate parameter
      parmapiter = bpmapiter->second.find(parname);
      if (parmapiter == bpmapiter->second.end()) {
        throw runtime_error("Parameter cannot be found in a bank's map.");
      } else {
        double pvalue = parmapiter->second;
        newrow << pvalue;
      }
 
    } // END(j)
  } // END(i)
 
  return tablews;
}
 
//----------------------------------------------------------------------------------------------
void LoadFullprofResolution::createBankToWorkspaceMap(const std::vector<int> &banks, const std::vector<int> &workspaces,
                                                      std::map<int, size_t> &workspaceOfBank) {
  if (workspaces.empty()) {
    for (size_t i = 0; i < banks.size(); i++) {
      workspaceOfBank.emplace(banks[i], i + 1);
    }
  } else {
    for (size_t i = 0; i < banks.size(); i++) {
      workspaceOfBank.emplace(banks[i], workspaces[i]);
    }
  }
}
 
//----------------------------------------------------------------------------------------------
void LoadFullprofResolution::putParametersIntoWorkspace(const API::Column_const_sptr &column,
                                                        const API::MatrixWorkspace_sptr &ws, int nProf,
                                                        std::string &parameterXMLString) {
 
  // Get instrument name from matrix workspace
  std::string instrumentName = ws->getInstrument()->getName();
 
  // Convert table workspace column into DOM XML document
  //   Set up writer to Paremeter file
  DOMWriter writer;
  writer.setNewLine("\n");
  writer.setOptions(XMLWriter::PRETTY_PRINT);
 
  //   Get current time
  Types::Core::DateAndTime date = Types::Core::DateAndTime::getCurrentTime();
  std::string ISOdate = date.toISO8601String();
  std::string ISOdateShort = ISOdate.substr(0, 19); // Remove fraction of seconds
 
  //   Create document
  AutoPtr<Document> mDoc = new Document();
  AutoPtr<Element> rootElem = mDoc->createElement("parameter-file");
  rootElem->setAttribute("date", ISOdateShort);
  mDoc->appendChild(rootElem);
 
  //   Add instrument
  AutoPtr<Element> instrumentElem = mDoc->createElement("component-link");
  instrumentElem->setAttribute("name", instrumentName);
  rootElem->appendChild(instrumentElem);
  if (nProf == 9) // put parameters into BackToBackExponential function
  {
    addBBX_S_Parameters(column, mDoc, instrumentElem);
    addBBX_A_Parameters(column, mDoc, instrumentElem);
    addBBX_B_Parameters(column, mDoc, instrumentElem);
  } else // Assume IkedaCarpenter PV
  {
    addALFBEParameter(column, mDoc, instrumentElem, "Alph0");
    addALFBEParameter(column, mDoc, instrumentElem, "Beta0");
    addALFBEParameter(column, mDoc, instrumentElem, "Alph1");
    addALFBEParameter(column, mDoc, instrumentElem, "Beta1");
    addSigmaParameters(column, mDoc, instrumentElem);
    addGammaParameters(column, mDoc, instrumentElem);
  }
 
  // Convert DOM XML document into string
  std::ostringstream outFile;
  writer.writeNode(outFile, mDoc);
  parameterXMLString = outFile.str();
 
  // Useful code for testing upgrades commented out for production use
  // std::ofstream outfileDebug("C:/Temp/test4_fullprof.xml");
  // outfileDebug << parameterXMLString;
  // outfileDebug.close();
}
 
/* Add an Ikeda Carpenter PV ALFBE parameter to the XML document according to
 *the table workspace
 *
 *  paramName is the name of the parameter as it appears in the table workspace
 */
void LoadFullprofResolution::addALFBEParameter(const API::Column_const_sptr &column, Poco::XML::Document *mDoc,
                                               Element *parent, const std::string &paramName) {
  AutoPtr<Element> parameterElem = mDoc->createElement("parameter");
  parameterElem->setAttribute("name", getXMLParameterName(paramName));
  parameterElem->setAttribute("type", "fitting");
 
  AutoPtr<Element> formulaElem = mDoc->createElement("formula");
  formulaElem->setAttribute("eq", getXMLEqValue(column, paramName));
  if (paramName != "Beta1")
    formulaElem->setAttribute("result-unit", "TOF");
  parameterElem->appendChild(formulaElem);
 
  AutoPtr<Element> fixedElem = mDoc->createElement("fixed");
  parameterElem->appendChild(fixedElem);
 
  parent->appendChild(parameterElem);
}
 
/* Add a set of Ikeda Carpenter PV SIGMA parameters to the XML document
 * according to the table workspace
 * for the bank at the given column of the table workspace
 */
void LoadFullprofResolution::addSigmaParameters(const API::Column_const_sptr &column, Poco::XML::Document *mDoc,
                                                Poco::XML::Element *parent) {
  AutoPtr<Element> parameterElem = mDoc->createElement("parameter");
  parameterElem->setAttribute("name", "IkedaCarpenterPV:SigmaSquared");
  parameterElem->setAttribute("type", "fitting");
 
  AutoPtr<Element> formulaElem = mDoc->createElement("formula");
  std::string eqValue = getXMLSquaredEqValue(column, "Sig1") + "*centre^2+" + getXMLSquaredEqValue(column, "Sig0");
  formulaElem->setAttribute("eq", eqValue);
  formulaElem->setAttribute("unit", "dSpacing");
  formulaElem->setAttribute("result-unit", "TOF^2");
  parameterElem->appendChild(formulaElem);
 
  parent->appendChild(parameterElem);
}
 
/* Add a set of Ikeda Carpenter PV GAMMA parameters to the XML document
 * according to the table workspace
 * for the bank at the given column of the table workspace
 */
void LoadFullprofResolution::addGammaParameters(const API::Column_const_sptr &column, Poco::XML::Document *mDoc,
                                                Poco::XML::Element *parent) {
  AutoPtr<Element> parameterElem = mDoc->createElement("parameter");
  parameterElem->setAttribute("name", "IkedaCarpenterPV:Gamma");
  parameterElem->setAttribute("type", "fitting");
 
  AutoPtr<Element> formulaElem = mDoc->createElement("formula");
  std::string eqValue = getXMLEqValue(column, "Gam1") + "*centre";
  formulaElem->setAttribute("eq", eqValue);
  formulaElem->setAttribute("unit", "dSpacing");
  formulaElem->setAttribute("result-unit", "TOF");
  parameterElem->appendChild(formulaElem);
 
  parent->appendChild(parameterElem);
}
 
/* Add a set of BackToBackExponential S parameters to the XML document according
 * to the table workspace
 * for the bank at the given column of the table workspace
 */
void LoadFullprofResolution::addBBX_S_Parameters(const API::Column_const_sptr &column, Poco::XML::Document *mDoc,
                                                 Poco::XML::Element *parent) {
  AutoPtr<Element> parameterElem = mDoc->createElement("parameter");
  parameterElem->setAttribute("name", "BackToBackExponential:S");
  parameterElem->setAttribute("type", "fitting");
 
  AutoPtr<Element> formulaElem = mDoc->createElement("formula");
  std::string eqValue = "sqrt(" + getXMLSquaredEqValue(column, "Sig2") + "*centre^4 + " +
                        getXMLSquaredEqValue(column, "Sig1") + "*centre^2 + " + getXMLSquaredEqValue(column, "Sig0") +
                        ")";
  formulaElem->setAttribute("eq", eqValue);
  formulaElem->setAttribute("unit", "dSpacing");
  formulaElem->setAttribute("result-unit", "TOF");
  parameterElem->appendChild(formulaElem);
 
  parent->appendChild(parameterElem);
}
 
/* Add a set of BackToBackExponential A parameters to the XML document according
 * to the table workspace
 * for the bank at the given column of the table workspace
 */
void LoadFullprofResolution::addBBX_A_Parameters(const API::Column_const_sptr &column, Poco::XML::Document *mDoc,
                                                 Poco::XML::Element *parent) {
  AutoPtr<Element> parameterElem = mDoc->createElement("parameter");
  parameterElem->setAttribute("name", "BackToBackExponential:A");
  parameterElem->setAttribute("type", "fitting");
 
  AutoPtr<Element> formulaElem = mDoc->createElement("formula");
  std::string eqValue = "(" + getXMLEqValue(column, "Alph1") + "/centre) + " + getXMLEqValue(column, "Alph0");
  formulaElem->setAttribute("eq", eqValue);
  formulaElem->setAttribute("unit", "dSpacing");
  formulaElem->setAttribute("result-unit", "TOF");
  parameterElem->appendChild(formulaElem);
 
  AutoPtr<Element> fixedElem = mDoc->createElement("fixed");
  parameterElem->appendChild(fixedElem);
 
  parent->appendChild(parameterElem);
}
 
/* Add a set of BackToBackExponential B parameters to the XML document according
 * to the table workspace
 * for the bank at the given column of the table workspace
 */
void LoadFullprofResolution::addBBX_B_Parameters(const API::Column_const_sptr &column, Poco::XML::Document *mDoc,
                                                 Poco::XML::Element *parent) {
  AutoPtr<Element> parameterElem = mDoc->createElement("parameter");
  parameterElem->setAttribute("name", "BackToBackExponential:B");
  parameterElem->setAttribute("type", "fitting");
 
  AutoPtr<Element> formulaElem = mDoc->createElement("formula");
  std::string eqValue = "(" + getXMLEqValue(column, "Beta1") + "/centre^4) + " + getXMLEqValue(column, "Beta0");
  formulaElem->setAttribute("eq", eqValue);
  formulaElem->setAttribute("unit", "dSpacing");
  formulaElem->setAttribute("result-unit", "TOF");
  parameterElem->appendChild(formulaElem);
 
  AutoPtr<Element> fixedElem = mDoc->createElement("fixed");
  parameterElem->appendChild(fixedElem);
 
  parent->appendChild(parameterElem);
}
 
/*
 *  Get the XML name of a parameter given its Table Workspace name
 */
std::string LoadFullprofResolution::getXMLParameterName(const std::string &name) {
  // Only used for ALFBE parameters
  std::string prefix = "IkedaCarpenterPV:";
  if (name == "Alph0")
    return prefix + "Alpha0";
  if (name == "Beta0")
    return prefix + "Beta0";
  if (name == "Alph1")
    return prefix + "Alpha1";
  if (name == "Beta1")
    return prefix + "Kappa";
  return "?" + name;
}
 
/*
 * Get the value string to put in the XML eq attribute of the formula element of
 * the paramenter element
 * given the name of the parameter in the table workspace.
 */
std::string LoadFullprofResolution::getXMLEqValue(const API::Column_const_sptr &column, const std::string &name) {
  size_t paramNumber = LoadFullprofResolution::m_rowNumbers[name];
  // API::Column_const_sptr column = tablews->getColumn( columnIndex );
  double eqValue = column->cell<double>(paramNumber);
  return boost::lexical_cast<std::string>(eqValue);
}
 
/*
 * Get the value string to put in the XML eq attribute of the formula element of
 * the SQUARED paramenter element
 * given the name of the parameter in the table workspace.
 */
std::string LoadFullprofResolution::getXMLSquaredEqValue(const API::Column_const_sptr &column,
                                                         const std::string &name) {
  size_t paramNumber = LoadFullprofResolution::m_rowNumbers[name];
  // API::Column_const_sptr column = tablews->getColumn( columnIndex );
  double eqValue = column->cell<double>(paramNumber);
  return boost::lexical_cast<std::string>(eqValue * eqValue);
}
 
/* This function fills in a list of the row numbers starting 0 of the parameters
   in the table workspace, so one can find the position in a column of
   the value of the given parameter.
*/
void LoadFullprofResolution::getTableRowNumbers(const API::ITableWorkspace_sptr &tablews,
                                                std::map<std::string, size_t> &parammap) {
  parammap.clear();
 
  size_t numrows = tablews->rowCount();
  for (size_t i = 0; i < numrows; ++i) {
    TableRow row = tablews->getRow(i);
    std::string paramName;
    row >> paramName;
    parammap.emplace(paramName, i);
  }
}
 
} // namespace Mantid::DataHandling