Commit c3d1dd66 authored by Matthew Bowles's avatar Matthew Bowles
Browse files

Merge branch '19341_crystal_field_multisite' of...

Merge branch '19341_crystal_field_multisite' of https://github.com/mantidproject/mantid into 20358_CrystalFieldFunction_wrapper
re #20398
parents dfb6324a 08289732
......@@ -150,6 +150,7 @@ set ( SRC_FILES
src/GetEi.cpp
src/GetEi2.cpp
src/GetEiMonDet2.cpp
src/GetQsInQENSData.cpp
src/GetTimeSeriesLogInformation.cpp
src/GravitySANSHelper.cpp
src/GroupWorkspaces.cpp
......@@ -469,6 +470,7 @@ set ( INC_FILES
inc/MantidAlgorithms/GetEi.h
inc/MantidAlgorithms/GetEi2.h
inc/MantidAlgorithms/GetEiMonDet2.h
inc/MantidAlgorithms/GetQsInQENSData.h
inc/MantidAlgorithms/GetTimeSeriesLogInformation.h
inc/MantidAlgorithms/GravitySANSHelper.h
inc/MantidAlgorithms/GroupWorkspaces.h
......@@ -799,6 +801,7 @@ set ( TEST_FILES
GetEiMonDet2Test.h
GetEiTest.h
GetEiV1Test.h
GetQsInQENSDataTest.h
GetTimeSeriesLogInformationTest.h
GroupWorkspacesTest.h
HRPDSlabCanAbsorptionTest.h
......
......@@ -80,11 +80,20 @@ private:
/// Creates an output workspace.
API::MatrixWorkspace_sptr
createOutputWorkspace(API::Progress &progress, const std::string &targetUnit,
API::MatrixWorkspace_sptr &inputWS, size_t nHist);
API::MatrixWorkspace_sptr &inputWS);
// Map to which the conversion to the unit is stored.
/// Map to which the conversion to the unit is stored.
std::multimap<double, size_t> m_indexMap;
/// Vector of axis in case ordering is not asked.
std::vector<double> m_axis;
/// Flag whether ordering is needed.
bool m_toOrder;
/// Emplaces to value and the index pair into the map.
void emplaceIndexMap(double value, size_t wsIndex);
/// Getting Efixed
double getEfixed(const size_t detectorIndex,
const Mantid::Geometry::DetectorInfo &detectorInfo,
......
#ifndef MANTID_ALGORITHMS_GETQSINQENSDATA_H
#define MANTID_ALGORITHMS_GETQSINQENSDATA_H
#include "MantidAPI/Algorithm.h"
#include "MantidAPI/MatrixWorkspace.h"
namespace Mantid {
namespace Algorithms {
/** Extracts Q-values from a workspace containing QENS data
Required Properties:
<UL>
<LI> InputWorkspace - The name of the workspace to extract the Q-values from
</LI>
<LI> RaiseMode - If set to true, exceptions will be raised, instead of an
empty list being returned </LI>
<LI> Qvalues - The Q-values extracted from the input workspace </LI>
</UL>
Copyright &copy; 2008 ISIS Rutherford Appleton Laboratory, NScD Oak Ridge
National Laboratory & European Spallation Source
This file is part of Mantid.
Mantid is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 3 of the License, or
(at your option) any later version.
Mantid is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with this program. If not, see <http://www.gnu.org/licenses/>.
File change history is stored at: <https://github.com/mantidproject/mantid>.
Code Documentation is available at: <http://doxygen.mantidproject.org>
*/
class DLLExport GetQsInQENSData : public API::Algorithm {
public:
/// Algorithm's category for identification. @see Algorithm::category
const std::string category() const override { return "Inelastic\\Indirect"; }
/// Algorithm's version for identification. @see Algorithm::version
int version() const override { return 1; }
/// Algorithm's summary for use in the GUI and help. @see Algorithm::summary
const std::string summary() const override {
return "Get Q-values in the vertical axis of a MatrixWorkspace containing "
"QENS S(Q,E) of S(theta,E) data.";
}
/// Algorithms name for identification. @see Algorithm::name
const std::string name() const override { return "GetQsInQENSData"; }
/// Cross-check properties with each other @see IAlgorithm::validateInputs
std::map<std::string, std::string> validateInputs() override;
private:
/// Initialization code
void init() override;
/// Execution code
void exec() override;
/// Extracts Q-values from the specified matrix workspace
MantidVec extractQValues(const Mantid::API::MatrixWorkspace_sptr workspace);
};
}
}
#endif /* MANTID_ALGORITHMS_GETQSINQENSDATA_H */
\ No newline at end of file
......@@ -49,6 +49,7 @@ public:
std::unique_ptr<Kernel::Property> property() const;
std::string propertyDoc() const;
std::string validateInputSize(const size_t size) const;
void applyInplace(HistogramData::Histogram &inOut, size_t stepSize) const;
void applyInPlace(const HistogramData::Histogram &in,
......
......@@ -62,6 +62,7 @@ public:
private:
void init() override;
void exec() override;
std::map<std::string, std::string> validateInputs() override;
API::MatrixWorkspace_uptr
doSimulation(const API::MatrixWorkspace &inputWS, const size_t nevents,
......
......@@ -55,6 +55,7 @@ private:
API::MatrixWorkspace_sptr load(const std::string filename);
void loadOldCalibration();
void createNewCalTable();
void createInformationWorkspaces();
std::function<double(double)> getDSpacingToTof(const detid_t detid);
std::vector<double> dSpacingWindows(const std::vector<double> &centres,
const double widthMax);
......@@ -62,13 +63,14 @@ private:
const double tzero);
void setCalibrationValues(const detid_t detid, const double difc,
const double difa, const double tzero);
void fitDIFCtZeroDIFA(const std::vector<double> &d,
const std::vector<double> &tof, double &difc,
double &t0, double &difa);
void fitDIFCtZeroDIFA_LM(const std::vector<double> &d,
const std::vector<double> &tof,
const std::vector<double> &height2, double &difc,
double &t0, double &difa);
API::MatrixWorkspace_sptr m_uncalibratedWS;
API::ITableWorkspace_sptr m_calibrationTable;
API::ITableWorkspace_sptr m_peakPositionTable;
std::vector<double> m_peaksInDspacing;
std::string calParams;
std::map<detid_t, size_t> m_detidToRow;
double m_tofMin{0.};
double m_tofMax{0.};
......@@ -77,6 +79,7 @@ private:
double m_difaMin{0.};
double m_difaMax{0.};
bool m_hasDasIds{false};
size_t m_numberMaxParams{0};
};
} // namespace Algorithms
......
#include "MantidAlgorithms/ConvertSpectrumAxis2.h"
#include "MantidTypes/SpectrumDefinition.h"
#include "MantidGeometry/Instrument/DetectorInfo.h"
#include "MantidAPI/InstrumentValidator.h"
#include "MantidAPI/NumericAxis.h"
#include "MantidAPI/Run.h"
......@@ -8,12 +6,14 @@
#include "MantidAPI/SpectrumInfo.h"
#include "MantidAPI/WorkspaceFactory.h"
#include "MantidGeometry/Instrument.h"
#include "MantidKernel/CompositeValidator.h"
#include "MantidGeometry/Instrument/DetectorInfo.h"
#include "MantidKernel/BoundedValidator.h"
#include "MantidKernel/CompositeValidator.h"
#include "MantidKernel/ListValidator.h"
#include "MantidKernel/Unit.h"
#include "MantidKernel/UnitConversion.h"
#include "MantidKernel/UnitFactory.h"
#include "MantidTypes/SpectrumDefinition.h"
#include <cfloat>
......@@ -65,18 +65,33 @@ void ConvertSpectrumAxis2::init() {
declareProperty("EFixed", EMPTY_DBL(), mustBePositive,
"Value of fixed energy in meV : EI (EMode=Direct) or EF "
"(EMode=Indirect))");
declareProperty("OrderAxis", true, "Whether or not to sort the resulting"
" spectrum axis.");
}
void ConvertSpectrumAxis2::exec() {
// Get the input workspace.
API::MatrixWorkspace_sptr inputWS = getProperty("InputWorkspace");
// Assign value to the member variable storing the number of histograms.
size_t nHist = inputWS->getNumberHistograms();
const size_t nHist = inputWS->getNumberHistograms();
// The unit to convert to.
const std::string unitTarget = getProperty("Target");
Progress progress(this, 0.0, 1.0, inputWS->getNumberHistograms());
// Whether needs to be ordered
m_toOrder = getProperty("OrderAxis");
size_t nProgress = nHist;
if (m_toOrder) {
// we will need to loop twice, once to build the indexMap,
// once to copy over the spectra and set the output
nProgress *= 2;
} else {
m_axis.reserve(nHist);
}
Progress progress(this, 0.0, 1.0, nProgress);
// Call the functions to convert to the different forms of theta or Q.
if (unitTarget == "theta" || unitTarget == "Theta" ||
......@@ -88,7 +103,7 @@ void ConvertSpectrumAxis2::exec() {
// Create an output workspace and set the property for it.
MatrixWorkspace_sptr outputWS =
createOutputWorkspace(progress, unitTarget, inputWS, nHist);
createOutputWorkspace(progress, unitTarget, inputWS);
setProperty("OutputWorkspace", outputWS);
}
......@@ -121,11 +136,11 @@ void ConvertSpectrumAxis2::createThetaMap(API::Progress &progress,
}
if (!spectrumInfo.isMonitor(i)) {
if (signedTheta)
m_indexMap.emplace(spectrumInfo.signedTwoTheta(i) * rad2deg, i);
emplaceIndexMap(spectrumInfo.signedTwoTheta(i) * rad2deg, i);
else
m_indexMap.emplace(spectrumInfo.twoTheta(i) * rad2deg, i);
emplaceIndexMap(spectrumInfo.twoTheta(i) * rad2deg, i);
} else {
m_indexMap.emplace(0.0, i);
emplaceIndexMap(0.0, i);
}
progress.report("Converting to theta...");
......@@ -172,16 +187,17 @@ void ConvertSpectrumAxis2::createElasticQMap(
}
// Convert to MomentumTransfer
double elasticQInAngstroms = Kernel::UnitConversion::run(theta, efixed);
double elasticQInAngstroms =
Kernel::UnitConversion::convertToElasticQ(theta, efixed);
if (targetUnit == "ElasticQ") {
m_indexMap.emplace(elasticQInAngstroms, i);
emplaceIndexMap(elasticQInAngstroms, i);
} else if (targetUnit == "ElasticQSquared") {
// The QSquared value.
double elasticQSquaredInAngstroms =
elasticQInAngstroms * elasticQInAngstroms;
m_indexMap.emplace(elasticQSquaredInAngstroms, i);
emplaceIndexMap(elasticQSquaredInAngstroms, i);
}
progress.report("Converting to Elastic Q...");
......@@ -194,22 +210,29 @@ void ConvertSpectrumAxis2::createElasticQMap(
* @param progress :: Progress indicator
* @param targetUnit :: Target conversion unit
* @param inputWS :: Input workspace
* @param nHist :: Stores the number of histograms
*/
MatrixWorkspace_sptr ConvertSpectrumAxis2::createOutputWorkspace(
API::Progress &progress, const std::string &targetUnit,
API::MatrixWorkspace_sptr &inputWS, size_t nHist) {
// Create the output workspace. Can not re-use the input one because the
// spectra are re-ordered.
MatrixWorkspace_sptr outputWorkspace = WorkspaceFactory::Instance().create(
inputWS, m_indexMap.size(), inputWS->x(0).size(), inputWS->y(0).size());
// Now set up a new numeric axis holding the theta values corresponding to
// each spectrum.
auto const newAxis = new NumericAxis(m_indexMap.size());
outputWorkspace->replaceAxis(1, newAxis);
API::MatrixWorkspace_sptr &inputWS) {
progress.setNumSteps(nHist + m_indexMap.size());
MatrixWorkspace_sptr outputWorkspace = nullptr;
NumericAxis *newAxis = nullptr;
if (m_toOrder) {
// Can not re-use the input one because the spectra are re-ordered.
outputWorkspace = WorkspaceFactory::Instance().create(
inputWS, m_indexMap.size(), inputWS->x(0).size(), inputWS->y(0).size());
std::vector<double> axis;
axis.reserve(m_indexMap.size());
for (const auto &it : m_indexMap) {
axis.emplace_back(it.first);
}
newAxis = new NumericAxis(std::move(axis));
} else {
// If there is no reordering we can simply clone.
outputWorkspace = inputWS->clone();
newAxis = new NumericAxis(m_axis);
}
outputWorkspace->replaceAxis(1, newAxis);
// Set the units of the axis.
if (targetUnit == "theta" || targetUnit == "Theta" ||
......@@ -221,20 +244,23 @@ MatrixWorkspace_sptr ConvertSpectrumAxis2::createOutputWorkspace(
newAxis->unit() = UnitFactory::Instance().create("QSquared");
}
std::multimap<double, size_t>::const_iterator it;
size_t currentIndex = 0;
for (it = m_indexMap.begin(); it != m_indexMap.end(); ++it) {
// Set the axis value.
newAxis->setValue(currentIndex, it->first);
// Copy over the data.
outputWorkspace->setHistogram(currentIndex, inputWS->histogram(it->second));
// We can keep the spectrum numbers etc.
outputWorkspace->getSpectrum(currentIndex)
.copyInfoFrom(inputWS->getSpectrum(it->second));
++currentIndex;
progress.report("Creating output workspace...");
// Note that this is needed only for ordered case
if (m_toOrder) {
size_t currentIndex = 0;
std::multimap<double, size_t>::const_iterator it;
for (it = m_indexMap.begin(); it != m_indexMap.end(); ++it) {
// Copy over the data.
outputWorkspace->setHistogram(currentIndex,
inputWS->histogram(it->second));
// We can keep the spectrum numbers etc.
outputWorkspace->getSpectrum(currentIndex)
.copyInfoFrom(inputWS->getSpectrum(it->second));
++currentIndex;
progress.report("Setting output spectrum #" +
std::to_string(currentIndex));
}
}
return outputWorkspace;
}
......@@ -278,5 +304,17 @@ double ConvertSpectrumAxis2::getEfixed(
return efixed;
}
/** Emplaces inside the ordered or unordered index registry
* @param value :: value to insert
* @param wsIndex :: workspace index
*/
void ConvertSpectrumAxis2::emplaceIndexMap(double value, size_t wsIndex) {
if (m_toOrder) {
m_indexMap.emplace(value, wsIndex);
} else {
m_axis.emplace_back(value);
}
}
} // namespace Algorithms
} // namespace Mantid
......@@ -309,6 +309,8 @@ MatrixWorkspace_sptr CreateSampleWorkspace::createHistogramWorkspace(
bool isRandom) {
BinEdges x(numBins + 1, LinearGenerator(x0, binDelta));
// there is a oddity here that y is evaluated from x=0, and x is from XMin
// changing it requires changing unit tests that use this algorithm
std::vector<double> xValues(cbegin(x), cend(x) - 1);
Counts y(evalFunction(functionString, xValues, isRandom ? 1 : 0));
......@@ -432,6 +434,7 @@ CreateSampleWorkspace::evalFunction(const std::string &functionString,
std::string replaceStr = boost::lexical_cast<std::string>(replace_val);
replaceAll(parsedFuncString, token, replaceStr);
}
g_log.information(parsedFuncString);
IFunction_sptr func_sptr =
FunctionFactory::Instance().createInitialized(parsedFuncString);
......
......@@ -118,14 +118,15 @@ void EQSANSTofStructure::execEvent(
const size_t numHists = inputWS->getNumberHistograms();
Progress progress(this, 0.0, 1.0, numHists);
// Get the nominal sample-to-detector distance (in mm)
// This now points to the correct distance and makes the naming clearer
// Get the nominal sample flange-to-detector distance (in mm)
Mantid::Kernel::Property *prop =
inputWS->run().getProperty("sample_detector_distance");
inputWS->run().getProperty("sampleflange_detector_distance");
auto dp = dynamic_cast<Mantid::Kernel::PropertyWithValue<double> *>(prop);
if (!dp) {
throw std::runtime_error("sample_detector_distance log not found.");
throw std::runtime_error("sampleflange_detector_distance log not found.");
}
const double SDD = *dp / 1000.0;
const double SFDD = *dp / 1000.0;
const auto &spectrumInfo = inputWS->spectrumInfo();
const auto l1 = spectrumInfo.l1();
......@@ -141,7 +142,7 @@ void EQSANSTofStructure::execEvent(
continue;
}
const auto l2 = spectrumInfo.l2(ispec);
double tof_factor = (l1 + l2) / (l1 + SDD);
double tof_factor = (l1 + l2) / (l1 + SFDD);
// Get the pointer to the output event list
std::vector<TofEvent> &events = inputWS->getSpectrum(ispec).getEvents();
......@@ -431,6 +432,8 @@ double EQSANSTofStructure::getTofOffset(EventWorkspace_const_sptr inputWS,
else
det_name = temp[0];
// Checked 8/11/2017 here detector_z is sfdd which has been updated
// in eqsansload.cpp
double source_z = inputWS->getInstrument()->getSource()->getPos().Z();
double detector_z =
inputWS->getInstrument()->getComponentByName(det_name)->getPos().Z();
......@@ -452,6 +455,7 @@ double EQSANSTofStructure::getTofOffset(EventWorkspace_const_sptr inputWS,
g_log.information() << i << " " << chopper_actual_phase[i] << " "
<< chopper_wl_1[i] << " " << chopper_wl_2[i] << '\n';
// Checked 8/10/2017
double low_wl_discard = 3.9560346 * low_tof_cut / source_to_detector;
double high_wl_discard = 3.9560346 * high_tof_cut / source_to_detector;
......
#include "MantidAlgorithms/GetQsInQENSData.h"
#include "MantidAPI/NumericAxis.h"
#include "MantidGeometry/IDetector.h"
#include "MantidGeometry/Instrument/DetectorInfo.h"
#include "MantidKernel/ArrayProperty.h"
#include "MantidKernel/UnitConversion.h"
#include "MantidKernel/Unit.h"
#include <stdexcept>
namespace {
Mantid::Kernel::Logger g_log("ConvolutionFitSequential");
}
namespace Mantid {
namespace Algorithms {
using namespace API;
using namespace Kernel;
using namespace Geometry;
// Register the Algorithm into the AlgorithmFactory
DECLARE_ALGORITHM(GetQsInQENSData)
// Initializes the Algorithm
void GetQsInQENSData::init() {
declareProperty(make_unique<WorkspaceProperty<MatrixWorkspace>>(
"InputWorkspace", "", Direction::Input),
"Input QENS data as MatrixWorkspace");
declareProperty("RaiseMode", false, "Set to True if an Exception, instead of "
"any empty list of Q values, is "
"desired.");
declareProperty(
make_unique<ArrayProperty<double>>("Qvalues", Direction::Output));
}
/*
* Validates the input properties
*/
std::map<std::string, std::string> GetQsInQENSData::validateInputs() {
std::map<std::string, std::string> issues;
MatrixWorkspace_sptr inputWs = getProperty("InputWorkspace");
// Check whether the input workspace could be found
if (!inputWs) {
issues["InputWorkspace"] = "InputWorkspace is not a MatrixWorkspace";
}
return issues;
}
/*
* Executes the Algorithm
*/
void GetQsInQENSData::exec() {
MatrixWorkspace_sptr inputWs = getProperty("InputWorkspace");
try {
setProperty("Qvalues", extractQValues(inputWs));
} catch (std::exception &e) {
g_log.error(e.what());
// If the RaiseMode property has been set to true, raise any
// exception which is thrown.
bool inRaiseMode = getProperty("RaiseMode");
if (inRaiseMode) {
throw;
}
setProperty("Qvalues", std::vector<double>());
}
}
/*
* Extracts Q-values from the specified workspace.
*
* @param workspace The workspace from which to extract Q-values.
* @return The extracted Q-values as a vector.
*/
MantidVec GetQsInQENSData::extractQValues(
const Mantid::API::MatrixWorkspace_sptr workspace) {
size_t numSpectra = workspace->getNumberHistograms();
Axis *qAxis;
try {
qAxis = workspace->getAxis(1);
} catch (std::exception &) {
throw std::runtime_error("Vertical axis is empty");
}
MantidVec qValues(qAxis->length());
// Check if the specified workspace is already in Q-space.
if (qAxis->unit()->unitID() == "MomentumTransfer") {
// Add axis values to vector of Q-values
for (size_t i = 0; i < qAxis->length(); i++) {
qValues[i] = qAxis->getValue(i);
}
// Check if the Q-values are stored as histogram data.
if (qValues.size() == numSpectra + 1) {
// Convert Q-values to point values.
qValues.pop_back();
qValues.erase(qValues.begin());
std::transform(qValues.begin(), qValues.end(), qValues.begin(),
std::bind2nd(std::divides<double>(), 2.0));
}
} else {
// Iterate over all spectrum in the specified workspace.
try {
for (size_t i = 0; i < numSpectra; i++) {
IDetector_const_sptr detector = workspace->getDetector(i);
double efixed = workspace->getEFixed(detector);
double theta = 0.5 * workspace->detectorTwoTheta(*detector);
qValues[i] = UnitConversion::convertToElasticQ(theta, efixed);
}
} catch (std::exception &) {
throw std::runtime_error(
"Detectors are missing from the input workspace");
}
}
return qValues;
}
}
}
......@@ -18,6 +18,8 @@ std::vector<std::string> OPTIONS{LINEAR_OPT, CSPLINE_OPT};
namespace Mantid {
using HistogramData::interpolateLinearInplace;
using HistogramData::interpolateCSplineInplace;
using HistogramData::minSizeForCSplineInterpolation;
using HistogramData::minSizeForLinearInterpolation;
using Kernel::Property;
namespace Algorithms {
......@@ -66,6 +68,32 @@ std::string InterpolationOption::propertyDoc() const {
return "Method of interpolation used to compute unsimulated values.";
}
/**
* Validate the size of input histogram
* @param size number of points in the input histogram
* @return an error message if there was one, otherwise an empty string
*/
std::string InterpolationOption::validateInputSize(const size_t size) const {
size_t nMin;
switch (m_value) {
case Value::Linear:
nMin = minSizeForLinearInterpolation();
if (size < nMin) {
return "Linear interpolation requires at least " + std::to_string(nMin) +
" points.";
}
break;
case Value::CSpline: