diff --git a/Framework/CurveFitting/inc/MantidCurveFitting/Functions/CrystalFieldControl.h b/Framework/CurveFitting/inc/MantidCurveFitting/Functions/CrystalFieldControl.h
index c13ca3b7e4e49806613f992378217c508bb0faf9..2f75b05900291d8235f0efb86188464fa576bec8 100644
--- a/Framework/CurveFitting/inc/MantidCurveFitting/Functions/CrystalFieldControl.h
+++ b/Framework/CurveFitting/inc/MantidCurveFitting/Functions/CrystalFieldControl.h
@@ -88,10 +88,6 @@ private:
std::vector<double> m_temperatures;
/// Cache the default peak FWHMs
std::vector<double> m_FWHMs;
- /// Cache number of fitted peaks
- // mutable std::vector<size_t> m_nPeaks;
- /// Cache the list of "spectra" corresponding to physical properties
- std::vector<int> m_physprops;
/// Caches of the width functions
std::vector<std::vector<double>> m_fwhmX;
std::vector<std::vector<double>> m_fwhmY;
diff --git a/Framework/CurveFitting/inc/MantidCurveFitting/Functions/CrystalFieldFunction.h b/Framework/CurveFitting/inc/MantidCurveFitting/Functions/CrystalFieldFunction.h
index 29080f77220e60f7118d96e4bf3d808d93ece112..4168512286dc9117c90f8c45d4693e0682b36575 100644
--- a/Framework/CurveFitting/inc/MantidCurveFitting/Functions/CrystalFieldFunction.h
+++ b/Framework/CurveFitting/inc/MantidCurveFitting/Functions/CrystalFieldFunction.h
@@ -241,6 +241,7 @@ private:
m_mapPrefixes2PhysProps;
/// Temporary cache for parameter values during source function resetting.
mutable std::vector<double> m_parameterResetCache;
+ mutable std::vector<bool> m_fixResetCache;
};
} // namespace Functions
diff --git a/Framework/CurveFitting/src/Functions/CrystalFieldControl.cpp b/Framework/CurveFitting/src/Functions/CrystalFieldControl.cpp
index b017bee64d023f445d6d516c859df8417a33d253..ff604962cb3b2be6c93d07fbc6a2c4517561df5a 100644
--- a/Framework/CurveFitting/src/Functions/CrystalFieldControl.cpp
+++ b/Framework/CurveFitting/src/Functions/CrystalFieldControl.cpp
@@ -55,10 +55,31 @@ void CrystalFieldControl::setAttribute(const std::string &name,
m_temperatures = attr.asVector();
buildControls();
} else if (name == "FWHMs") {
+ const size_t nSpec = m_temperatures.size();
m_FWHMs = attr.asVector();
- if (m_FWHMs.size() == 1 && m_FWHMs.size() != m_temperatures.size()) {
- m_FWHMs.assign(m_temperatures.size(), m_FWHMs.front());
+ if (m_FWHMs.size() == 1 && m_FWHMs.size() != nSpec) {
+ m_FWHMs.assign(nSpec, m_FWHMs.front());
+ }
+ if (!m_FWHMs.empty()) {
+ m_fwhmX.resize(nSpec);
+ m_fwhmY.resize(nSpec);
+ for (size_t i = 0; i < nSpec; ++i) {
+ m_fwhmX[i].clear();
+ m_fwhmY[i].clear();
+ if (nSpec > 1) {
+ auto &control = *getFunction(i);
+ control.setAttribute("FWHMX", Attribute(std::vector<double>()));
+ control.setAttribute("FWHMY", Attribute(std::vector<double>()));
+ } else {
+ setAttribute("FWHMX", Attribute(std::vector<double>()));
+ setAttribute("FWHMY", Attribute(std::vector<double>()));
+ }
+ }
}
+ } else if ((name.compare(0, 5, "FWHMX") == 0 ||
+ name.compare(0, 5, "FWHMY") == 0) &&
+ !attr.asVector().empty()) {
+ m_FWHMs.clear();
}
API::IFunction::setAttribute(name, attr);
}
diff --git a/Framework/CurveFitting/src/Functions/CrystalFieldFunction.cpp b/Framework/CurveFitting/src/Functions/CrystalFieldFunction.cpp
index df242ad42379f308cfd3986aa84597715baa8991..99c8063ad8e508a98fbbbb2e88b7ef10da61b6d7 100644
--- a/Framework/CurveFitting/src/Functions/CrystalFieldFunction.cpp
+++ b/Framework/CurveFitting/src/Functions/CrystalFieldFunction.cpp
@@ -492,6 +492,10 @@ void CrystalFieldFunction::setAttribute(const std::string &attName,
m_source.reset();
}
attRef.first->setAttribute(attRef.second, attr);
+ if (attName.find("FWHM") != std::string::npos ||
+ attName.find("Background") != std::string::npos) {
+ m_dirtyTarget = true;
+ }
}
/// Check if attribute attName exists
@@ -632,9 +636,12 @@ void CrystalFieldFunction::buildSourceFunction() const {
m_parameterResetCache.size() == m_source->nParams()) {
for (size_t i = 0; i < m_parameterResetCache.size(); ++i) {
m_source->setParameter(i, m_parameterResetCache[i]);
+ if (m_fixResetCache[i])
+ m_source->fix(i);
}
}
m_parameterResetCache.clear();
+ m_fixResetCache.clear();
}
/// Update spectrum function if necessary.
@@ -741,9 +748,9 @@ void CrystalFieldFunction::buildSingleSiteMultiSpectrum() const {
for (size_t i = 0; i < nSpec; ++i) {
auto intensityScaling =
m_control.getFunction(i)->getParameter("IntensityScaling");
- fun->addFunction(buildSpectrum(nre, energies, waveFunctions,
- temperatures[i], FWHMs[i], i, addBackground,
- intensityScaling));
+ fun->addFunction(buildSpectrum(
+ nre, energies, waveFunctions, temperatures[i],
+ FWHMs.size() > i ? FWHMs[i] : 0., i, addBackground, intensityScaling));
fun->setDomainIndex(i, i);
}
auto &physProps = m_control.physProps();
@@ -843,9 +850,10 @@ void CrystalFieldFunction::buildMultiSiteMultiSpectrum() const {
for (size_t i = 0; i < nSpec; ++i) {
auto spectrumIntensityScaling =
m_control.getFunction(i)->getParameter("IntensityScaling");
- spectra[i]->addFunction(buildSpectrum(
- nre, energies, waveFunctions, temperatures[i], FWHMs[i], i,
- addBackground, ionIntensityScaling * spectrumIntensityScaling));
+ spectra[i]->addFunction(
+ buildSpectrum(nre, energies, waveFunctions, temperatures[i],
+ FWHMs.size() > i ? FWHMs[i] : 0., i, addBackground,
+ ionIntensityScaling * spectrumIntensityScaling));
}
size_t i = 0;
@@ -934,7 +942,7 @@ API::IFunction_sptr CrystalFieldFunction::buildSpectrum(
}
auto xVec = m_control.getFunction(iSpec)->getAttribute("FWHMX").asVector();
- auto yVec = m_control.getFunction(iSpec)->getAttribute("FWHMX").asVector();
+ auto yVec = m_control.getFunction(iSpec)->getAttribute("FWHMY").asVector();
CrystalFieldUtils::buildSpectrumFunction(*spectrum, peakShape, values, xVec,
yVec, fwhmVariation, fwhm,
nRequiredPeaks, fixAllPeaks);
@@ -1052,7 +1060,7 @@ void CrystalFieldFunction::updateSingleSiteSingleSpectrum() const {
auto peakShape = m_control.getAttribute("PeakShape").asString();
bool fixAllPeaks = m_control.getAttribute("FixAllPeaks").asBool();
auto xVec = m_control.getAttribute("FWHMX").asVector();
- auto yVec = m_control.getAttribute("FWHMX").asVector();
+ auto yVec = m_control.getAttribute("FWHMY").asVector();
auto &FWHMs = m_control.FWHMs();
auto defaultFWHM = FWHMs.empty() ? 0.0 : FWHMs[0];
size_t indexShift = hasBackground() ? 1 : 0;
@@ -1085,7 +1093,8 @@ void CrystalFieldFunction::updateSingleSiteMultiSpectrum() const {
auto &FWHMs = m_control.FWHMs();
for (size_t iSpec = 0; iSpec < temperatures.size(); ++iSpec) {
updateSpectrum(*fun.getFunction(iSpec), nre, energies, waveFunctions,
- temperatures[iSpec], FWHMs[iSpec], iSpec, iFirst);
+ temperatures[iSpec],
+ FWHMs.size() > iSpec ? FWHMs[iSpec] : 0., iSpec, iFirst);
}
for (auto &prop : m_mapPrefixes2PhysProps) {
@@ -1099,7 +1108,7 @@ void CrystalFieldFunction::updateMultiSiteSingleSpectrum() const {
auto peakShape = getAttribute("PeakShape").asString();
bool fixAllPeaks = getAttribute("FixAllPeaks").asBool();
auto xVec = m_control.getAttribute("FWHMX").asVector();
- auto yVec = m_control.getAttribute("FWHMX").asVector();
+ auto yVec = m_control.getAttribute("FWHMY").asVector();
auto &FWHMs = m_control.FWHMs();
auto defaultFWHM = FWHMs.empty() ? 0.0 : FWHMs[0];
@@ -1142,7 +1151,8 @@ void CrystalFieldFunction::updateMultiSiteMultiSpectrum() const {
auto &ionSpectrum =
dynamic_cast<CompositeFunction &>(*spectrum.getFunction(ionIndex));
updateSpectrum(ionSpectrum, nre, energies, waveFunctions,
- temperatures[iSpec], FWHMs[iSpec], iSpec, iFirst);
+ temperatures[iSpec],
+ FWHMs.size() > iSpec ? FWHMs[iSpec] : 0., iSpec, iFirst);
}
std::string prefix("ion");
@@ -1174,7 +1184,7 @@ void CrystalFieldFunction::updateSpectrum(
const auto peakShape = getAttribute("PeakShape").asString();
const bool fixAllPeaks = getAttribute("FixAllPeaks").asBool();
auto xVec = m_control.getFunction(iSpec)->getAttribute("FWHMX").asVector();
- auto yVec = m_control.getFunction(iSpec)->getAttribute("FWHMX").asVector();
+ auto yVec = m_control.getFunction(iSpec)->getAttribute("FWHMY").asVector();
auto intensityScaling =
m_control.getFunction(iSpec)->getParameter("IntensityScaling");
FunctionValues values;
@@ -1390,8 +1400,10 @@ void CrystalFieldFunction::cacheSourceParameters() const {
}
auto np = m_source->nParams();
m_parameterResetCache.resize(np);
+ m_fixResetCache.resize(np);
for (size_t i = 0; i < np; ++i) {
m_parameterResetCache[i] = m_source->getParameter(i);
+ m_fixResetCache[i] = m_source->isFixed(i);
}
}
diff --git a/Framework/CurveFitting/src/Functions/CrystalFieldMultiSpectrum.cpp b/Framework/CurveFitting/src/Functions/CrystalFieldMultiSpectrum.cpp
index 246fec09c7f9e1398a8f8b857c20006f2e891e37..f1af30eea7407b205972ccd4e73adc63fcabba0a 100644
--- a/Framework/CurveFitting/src/Functions/CrystalFieldMultiSpectrum.cpp
+++ b/Framework/CurveFitting/src/Functions/CrystalFieldMultiSpectrum.cpp
@@ -15,6 +15,7 @@
#include "MantidAPI/ParameterTie.h"
#include "MantidKernel/Exception.h"
+#include <boost/regex.hpp>
#include <iostream>
namespace Mantid {
@@ -29,6 +30,10 @@ DECLARE_FUNCTION(CrystalFieldMultiSpectrum)
namespace {
+// Regex for the FWHMX# type strings (single-site mode)
+const boost::regex FWHMX_ATTR_REGEX("FWHMX([0-9]+)");
+const boost::regex FWHMY_ATTR_REGEX("FWHMY([0-9]+)");
+
/// Define the source function for CrystalFieldMultiSpectrum.
/// Its function() method is not needed.
class Peaks : public CrystalFieldPeaksBase, public API::IFunctionGeneral {
@@ -51,11 +56,15 @@ public:
std::vector<size_t> m_IntensityScalingIdx;
std::vector<size_t> m_PPLambdaIdxChild;
std::vector<size_t> m_PPLambdaIdxSelf;
+ std::vector<size_t> m_PPChi0IdxChild;
+ std::vector<size_t> m_PPChi0IdxSelf;
/// Declare the intensity scaling parameters: one per spectrum.
void declareIntensityScaling(size_t nSpec) {
m_IntensityScalingIdx.clear();
m_PPLambdaIdxChild.resize(nSpec, -1);
m_PPLambdaIdxSelf.resize(nSpec, -1);
+ m_PPChi0IdxChild.resize(nSpec, -1);
+ m_PPChi0IdxSelf.resize(nSpec, -1);
for (size_t i = 0; i < nSpec; ++i) {
auto si = std::to_string(i);
try { // If parameter has already been declared, don't declare it.
@@ -71,6 +80,8 @@ public:
if (m_PPLambdaIdxSelf.size() <= iSpec) {
m_PPLambdaIdxSelf.resize(iSpec + 1, -1);
m_PPLambdaIdxChild.resize(iSpec + 1, -1);
+ m_PPChi0IdxSelf.resize(iSpec + 1, -1);
+ m_PPChi0IdxChild.resize(iSpec + 1, -1);
}
auto si = std::to_string(iSpec);
try { // If parameter has already been declared, don't declare it.
@@ -78,7 +89,13 @@ public:
"Effective exchange coupling of dataset " + si);
} catch (std::invalid_argument &) {
}
+ try { // If parameter has already been declared, don't declare it.
+ declareParameter("Chi0" + si, 0.0,
+ "Effective exchange coupling of dataset " + si);
+ } catch (std::invalid_argument &) {
+ }
m_PPLambdaIdxSelf[iSpec] = parameterIndex("Lambda" + si);
+ m_PPChi0IdxSelf[iSpec] = parameterIndex("Chi0" + si);
}
};
}
@@ -123,6 +140,7 @@ CrystalFieldMultiSpectrum::createEquivalentFunctions() const {
/// Perform custom actions on setting certain attributes.
void CrystalFieldMultiSpectrum::setAttribute(const std::string &name,
const Attribute &attr) {
+ boost::smatch match;
if (name == "Temperatures") {
// Define (declare) the parameters for intensity scaling.
const auto nSpec = attr.asVector().size();
@@ -130,13 +148,23 @@ void CrystalFieldMultiSpectrum::setAttribute(const std::string &name,
m_nOwnParams = m_source->nParams();
m_fwhmX.resize(nSpec);
m_fwhmY.resize(nSpec);
+ std::vector<double> new_fwhm = getAttribute("FWHMs").asVector();
+ const auto nWidths = new_fwhm.size();
+ if (nWidths != nSpec) {
+ new_fwhm.resize(nSpec);
+ if (nWidths > nSpec) {
+ for (size_t iSpec = nWidths; iSpec < nSpec; ++iSpec) {
+ new_fwhm[iSpec] = new_fwhm[0];
+ }
+ }
+ }
+ FunctionGenerator::setAttribute("FWHMs", Attribute(new_fwhm));
for (size_t iSpec = 0; iSpec < nSpec; ++iSpec) {
const auto suffix = std::to_string(iSpec);
declareAttribute("FWHMX" + suffix, Attribute(m_fwhmX[iSpec]));
declareAttribute("FWHMY" + suffix, Attribute(m_fwhmY[iSpec]));
}
- }
- if (name == "PhysicalProperties") {
+ } else if (name == "PhysicalProperties") {
const auto physpropId = attr.asVector();
const auto nSpec = physpropId.size();
auto &source = dynamic_cast<Peaks &>(*m_source);
@@ -162,6 +190,12 @@ void CrystalFieldMultiSpectrum::setAttribute(const std::string &name,
m_nOwnParams = m_source->nParams();
}
}
+ } else if (boost::regex_match(name, match, FWHMX_ATTR_REGEX)) {
+ auto iSpec = std::stoul(match[1]);
+ m_fwhmX[iSpec].clear();
+ } else if (boost::regex_match(name, match, FWHMY_ATTR_REGEX)) {
+ auto iSpec = std::stoul(match[1]);
+ m_fwhmY[iSpec].clear();
}
FunctionGenerator::setAttribute(name, attr);
}
@@ -332,6 +366,8 @@ API::IFunction_sptr CrystalFieldMultiSpectrum::buildPhysprop(
spectrum.setAttribute("powder", getAttribute("powder" + suffix));
dynamic_cast<Peaks &>(*m_source).m_PPLambdaIdxChild[iSpec] =
spectrum.parameterIndex("Lambda");
+ dynamic_cast<Peaks &>(*m_source).m_PPChi0IdxChild[iSpec] =
+ spectrum.parameterIndex("Chi0");
return retval;
}
case Magnetisation: {
@@ -408,6 +444,8 @@ void CrystalFieldMultiSpectrum::updateSpectrum(
auto &source = dynamic_cast<Peaks &>(*m_source);
suscept.setParameter(source.m_PPLambdaIdxChild[iSpec],
getParameter(source.m_PPLambdaIdxSelf[iSpec]));
+ suscept.setParameter(source.m_PPChi0IdxChild[iSpec],
+ getParameter(source.m_PPChi0IdxSelf[iSpec]));
break;
}
case Magnetisation: {
diff --git a/Framework/CurveFitting/src/Functions/CrystalFieldSusceptibility.cpp b/Framework/CurveFitting/src/Functions/CrystalFieldSusceptibility.cpp
index 123c6ef59a6ce8622eee2bfe981d38957fa62ac5..26d83385c0c99c0183cbd8861b513543c087d8a8 100644
--- a/Framework/CurveFitting/src/Functions/CrystalFieldSusceptibility.cpp
+++ b/Framework/CurveFitting/src/Functions/CrystalFieldSusceptibility.cpp
@@ -139,6 +139,19 @@ void CrystalFieldSusceptibilityBase::function1D(double *out,
} else {
calculate(out, xValues, nData, m_en, m_wf, m_nre, H, convfact);
}
+ const double lambda = getParameter("Lambda");
+ const double EPS = 1.e-6;
+ if (fabs(lambda) > EPS) {
+ for (size_t i = 0; i < nData; i++) {
+ out[i] /= (1. - lambda * out[i]); // chi = chi_cf/(1 - lambda.chi_cf)
+ }
+ }
+ const double chi0 = getParameter("Chi0");
+ if (fabs(lambda) > 1.e-6) {
+ for (size_t i = 0; i < nData; i++) {
+ out[i] += chi0;
+ }
+ }
if (getAttribute("inverse").asBool()) {
for (size_t i = 0; i < nData; i++) {
out[i] = 1. / out[i];
@@ -150,12 +163,6 @@ void CrystalFieldSusceptibilityBase::function1D(double *out,
out[i] *= fact;
}
}
- const double lambda = getParameter("Lambda");
- if (fabs(lambda) > 1.e-6) {
- for (size_t i = 0; i < nData; i++) {
- out[i] /= (1. - lambda * out[i]); // chi = chi0/(1 - lambda.chi0)
- }
- }
}
DECLARE_FUNCTION(CrystalFieldSusceptibility)
@@ -164,6 +171,7 @@ CrystalFieldSusceptibility::CrystalFieldSusceptibility()
: CrystalFieldPeaksBase(), CrystalFieldSusceptibilityBase(),
m_setDirect(false) {
declareParameter("Lambda", 0.0, "Effective exchange interaction");
+ declareParameter("Chi0", 0.0, "Background or remnant susceptibility");
}
// Sets the eigenvectors / values directly
@@ -186,6 +194,7 @@ void CrystalFieldSusceptibility::function1D(double *out, const double *xValues,
CrystalFieldSusceptibilityCalculation::CrystalFieldSusceptibilityCalculation() {
declareParameter("Lambda", 0.0, "Effective exchange interaction");
+ declareParameter("Chi0", 0.0, "Background or remnant susceptibility");
}
// Sets the eigenvectors / values directly
diff --git a/Framework/PythonInterface/mantid/api/src/Exports/IFunction.cpp b/Framework/PythonInterface/mantid/api/src/Exports/IFunction.cpp
index 7e132662437cd9b759b352623aa24684a488dc16..0406931e6735830975e7cd13522739a56909c6a0 100644
--- a/Framework/PythonInterface/mantid/api/src/Exports/IFunction.cpp
+++ b/Framework/PythonInterface/mantid/api/src/Exports/IFunction.cpp
@@ -209,6 +209,9 @@ void export_IFunction() {
(void (IFunction::*)(const std::string &)) & IFunction::removeTie,
(arg("self"), arg("name")), "Remove the tie of the named parameter")
+ .def("getTies", &IFunction::writeTies, arg("self"),
+ "Returns the list of current ties as a string")
+
.def("addConstraints", &IFunction::addConstraints,
addConstraints_Overloads(
(arg("self"), arg("constraints"), arg("isDefault")),
@@ -218,6 +221,9 @@ void export_IFunction() {
(arg("self"), arg("name")),
"Remove the constraint on the named parameter")
+ .def("getConstraints", &IFunction::writeConstraints, arg("self"),
+ "Returns the list of current constraints as a string")
+
.def("getNumberDomains", &IFunction::getNumberDomains, (arg("self")),
"Get number of domains of a multi-domain function")
diff --git a/Testing/Data/SystemTest/Sm2O3.cif.md5 b/Testing/Data/SystemTest/Sm2O3.cif.md5
new file mode 100644
index 0000000000000000000000000000000000000000..4ca349a92ca0f0a33fefdf7e4e9a73c8310bbcf6
--- /dev/null
+++ b/Testing/Data/SystemTest/Sm2O3.cif.md5
@@ -0,0 +1 @@
+1d127c4d48f67c2260e06e1844e87725
diff --git a/Testing/SystemTests/tests/analysis/CrystalFieldPythonInterface.py b/Testing/SystemTests/tests/analysis/CrystalFieldPythonInterface.py
new file mode 100644
index 0000000000000000000000000000000000000000..7dbf920e6b907049e3f99fbc7bcd66e371abd162
--- /dev/null
+++ b/Testing/SystemTests/tests/analysis/CrystalFieldPythonInterface.py
@@ -0,0 +1,391 @@
+#pylint: disable=no-init
+from stresstesting import MantidStressTest
+from mantid.simpleapi import CreateWorkspace, LinearBackground, FlatBackground, Gaussian, CalculateChiSquared, mtd
+from CrystalField.fitting import makeWorkspace
+from PyChop import PyChop2
+import numpy as np
+
+
+class CrystalFieldPythonInterface(MantidStressTest):
+ """ Runs all the commands in the crystal field python interface documentation file
+ This test only checks that the syntax of all the commands are still valid.
+ Unit tests for the individual functions / fitting procedure check for correctness.
+ """
+
+ def my_create_ws(self, outwsname, x, y):
+ jitter = (np.random.rand(np.shape(y)[0])-0.5)*np.max(y)/100
+ CreateWorkspace(x, y + jitter, y*0+np.max(y)/100, Distribution=True, OutputWorkspace=outwsname)
+ return mtd[outwsname]
+
+ def my_func(self, en):
+ return (25-en)**(1.5) / 200 + 0.1
+
+ def runTest(self):
+ from CrystalField import CrystalField, CrystalFieldFit, CrystalFieldMultiSite, Background, Function, ResolutionModel
+
+ cf = CrystalField('Ce', 'C2v')
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770)
+ cf['B40'] = -0.031
+ b = cf['B40']
+
+ # Calculate and return the Hamiltonian matrix as a 2D numpy array.
+ h = cf.getHamiltonian()
+ print(h)
+ # Calculate and return the eigenvalues of the Hamiltonian as a 1D numpy array.
+ e = cf.getEigenvalues()
+ print(e)
+ # Calculate and return the eigenvectors of the Hamiltonian as a 2D numpy array.
+ w = cf.getEigenvectors()
+ print(w)
+ # Using the keyword argument
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, Temperature=44)
+ # Using the property
+ cf.Temperature = 44
+
+ print(cf.getPeakList())
+ #[[ 0.00000000e+00 2.44006198e+01 4.24977124e+01 1.80970926e+01 -2.44006198e+01]
+ # [ 2.16711565e+02 8.83098530e+01 5.04430056e+00 1.71153708e-01 1.41609425e-01]]
+ cf.ToleranceIntensity = 1
+ print(cf.getPeakList())
+ #[[ 0. 24.40061976 42.49771237]
+ # [ 216.71156467 88.30985303 5.04430056]]
+ cf.PeakShape = 'Gaussian'
+ cf.FWHM = 0.9
+ sp = cf.getSpectrum()
+ print(cf.function)
+ CrystalField_Ce = CreateWorkspace(*sp)
+ print(CrystalField_Ce)
+
+ # If the peak shape is Gaussian
+ cf.peaks.param[1]['Sigma'] = 2.0
+ cf.peaks.param[2]['Sigma'] = 0.01
+
+ # If the peak shape is Lorentzian
+ cf.PeakShape = 'Lorentzian'
+ cf.peaks.param[1]['FWHM'] = 2.0
+ cf.peaks.param[2]['FWHM'] = 0.01
+
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.12544,
+ Temperature=44.0, FWHM=1.1)
+ cf.background = Background(peak=Function('Gaussian', Height=10, Sigma=1),
+ background=Function('LinearBackground', A0=1.0, A1=0.01))
+ h = cf.background.peak.param['Height']
+ a1 = cf.background.background.param['A1']
+ print(h)
+ print(a1)
+
+ cf.ties(B20=1.0, B40='B20/2')
+ cf.constraints('1 < B22 <= 2', 'B22 < 4')
+
+ print(cf.function[1].getTies())
+ print(cf.function[1].getConstraints())
+
+ cf.background.peak.ties(Height=10.1)
+ cf.background.peak.constraints('Sigma > 0')
+ cf.background.background.ties(A0=0.1)
+ cf.background.background.constraints('A1 > 0')
+
+ print(cf.function[0][0].getConstraints())
+ print(cf.function[0][1].getConstraints())
+ print(cf.function.getTies())
+ print(cf.function.getConstraints())
+
+ cf.peaks.ties({'f2.FWHM': '2*f1.FWHM', 'f3.FWHM': '2*f2.FWHM'})
+ cf.peaks.constraints('f0.FWHM < 2.2', 'f1.FWHM >= 0.1')
+
+ cf.PeakShape = 'Gaussian'
+ cf.peaks.tieAll('Sigma=0.1', 3)
+ cf.peaks.constrainAll('0 < Sigma < 0.1', 4)
+ cf.peaks.tieAll('Sigma=f0.Sigma', 1, 3)
+ cf.peaks.ties({'f1.Sigma': 'f0.Sigma', 'f2.Sigma': 'f0.Sigma', 'f3.Sigma': 'f0.Sigma'})
+
+ rm = ResolutionModel(([1, 2, 3, 100], [0.1, 0.3, 0.35, 2.1]))
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, Temperature=44.0, ResolutionModel=rm)
+
+ rm = ResolutionModel(self.my_func, xstart=0.0, xend=24.0, accuracy=0.01)
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, Temperature=44.0, ResolutionModel=rm)
+
+ marires = PyChop2('MARI')
+ marires.setChopper('S')
+ marires.setFrequency(250)
+ marires.setEi(30)
+ rm = ResolutionModel(marires.getResolution, xstart=0.0, xend=29.0, accuracy=0.01)
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, Temperature=44.0, ResolutionModel=rm)
+
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, Temperature=44.0, ResolutionModel=rm, FWHMVariation=0.1)
+
+ # ---------------------------
+
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.12544,
+ Temperature=[44.0, 50], FWHM=[1.1, 0.9])
+ cf.PeakShape = 'Lorentzian'
+ cf.peaks[0].param[0]['FWHM'] = 1.11
+ cf.peaks[1].param[1]['FWHM'] = 1.12
+ cf.background = Background(peak=Function('Gaussian', Height=10, Sigma=0.3),
+ background=Function('FlatBackground', A0=1.0))
+ cf.background[1].peak.param['Sigma'] = 0.8
+ cf.background[1].background.param['A0'] = 1.1
+
+ # The B parameters are common for all spectra - syntax doesn't change
+ cf.ties(B20=1.0, B40='B20/2')
+ cf.constraints('1 < B22 <= 2', 'B22 < 4')
+
+ # Backgrounds and peaks are different for different spectra - must be indexed
+ cf.background[0].peak.ties(Height=10.1)
+ cf.background[0].peak.constraints('Sigma > 0.1')
+ cf.background[1].peak.ties(Height=20.2)
+ cf.background[1].peak.constraints('Sigma > 0.2')
+ cf.peaks[1].tieAll('FWHM=2*f1.FWHM', 2, 5)
+ cf.peaks[0].constrainAll('FWHM < 2.2', 1, 4)
+
+ rm = ResolutionModel([self.my_func, marires.getResolution], 0, 100, accuracy = 0.01)
+ cf.ResolutionModel = rm
+
+ # Calculate second spectrum, use the generated x-values
+ sp = cf.getSpectrum(1)
+ # Calculate second spectrum, use the first spectrum of a workspace
+ sp = cf.getSpectrum(1, 'CrystalField_Ce')
+ # Calculate first spectrum, use the i-th spectrum of a workspace
+ i=0
+ sp = cf.getSpectrum(0, 'CrystalField_Ce', i)
+
+ print(cf.function)
+
+ cf.Temperature = [5, 50, 150]
+
+ print()
+ print(cf.function)
+
+ ws = 'CrystalField_Ce'
+ ws1 = 'CrystalField_Ce'
+ ws2 = 'CrystalField_Ce'
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.12544, Temperature=5)
+
+ # In case of a single spectrum (ws is a workspace)
+ fit = CrystalFieldFit(Model=cf, InputWorkspace=ws)
+
+ # Or for multiple spectra
+ fit = CrystalFieldFit(Model=cf, InputWorkspace=[ws1, ws2])
+ cf.Temperature = [5, 50]
+ fit.fit()
+
+ params = {'B20': 0.377, 'B22': 3.9, 'B40': -0.03, 'B42': -0.116, 'B44': -0.125,
+ 'Temperature': [44.0, 50], 'FWHM': [1.1, 0.9]}
+ cf1 = CrystalField('Ce', 'C2v', **params)
+ cf2 = CrystalField('Pr', 'C2v', **params)
+ cfms = cf1 + cf2
+ cf = 2*cf1 + cf2
+
+ cfms = CrystalFieldMultiSite(Ions=['Ce', 'Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44.0], FWHMs=[1.1])
+ cfms['ion0.B40'] = -0.031
+ cfms['ion1.B20'] = 0.37737
+ b = cfms['ion0.B22']
+
+ print(b)
+ print(cfms.function)
+
+ cfms = CrystalFieldMultiSite(Ions=['Ce', 'Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44.0], FWHMs=[1.1],
+ parameters={'ion0.B20': 0.37737, 'ion0.B22': 3.9770, 'ion1.B40':-0.031787,
+ 'ion1.B42':-0.11611, 'ion1.B44':-0.12544})
+ cfms = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHMs=[1.0],
+ Background='name=Gaussian,Height=0,PeakCentre=1,Sigma=0;name=LinearBackground,A0=0,A1=0')
+ cfms = CrystalFieldMultiSite(Ions=['Ce'], Symmetries=['C2v'], Temperatures=[50], FWHMs=[0.9],
+ Background=LinearBackground(A0=1.0), BackgroundPeak=Gaussian(Height=10, Sigma=0.3))
+ cfms = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHMs=[1.0],
+ Background= Gaussian(PeakCentre=1) + LinearBackground())
+ cfms = CrystalFieldMultiSite(Ions=['Ce','Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44, 50], FWHMs=[1.1, 0.9],
+ Background=FlatBackground(), BackgroundPeak=Gaussian(Height=10, Sigma=0.3),
+ parameters={'ion0.B20': 0.37737, 'ion0.B22': 3.9770, 'ion1.B40':-0.031787,
+ 'ion1.B42':-0.11611, 'ion1.B44':-0.12544})
+ cfms.ties({'sp0.bg.f0.Height': 10.1})
+ cfms.constraints('sp0.bg.f0.Sigma > 0.1')
+ cfms.constraints('ion0.sp0.pk1.FWHM < 2.2')
+ cfms.ties({'ion0.sp1.pk2.FWHM': '2*ion0.sp1.pk1.FWHM', 'ion1.sp1.pk3.FWHM': '2*ion1.sp1.pk2.FWHM'})
+
+ # --------------------------
+
+ params = {'ion0.B20': 0.37737, 'ion0.B22': 3.9770, 'ion1.B40':-0.031787, 'ion1.B42':-0.11611, 'ion1.B44':-0.12544}
+ cf = CrystalFieldMultiSite(Ions=['Ce', 'Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44.0, 50.0],
+ FWHMs=[1.0, 1.0], ToleranceIntensity=6.0, ToleranceEnergy=1.0, FixAllPeaks=True,
+ parameters=params)
+
+ cf.fix('ion0.BmolX', 'ion0.BmolY', 'ion0.BmolZ', 'ion0.BextX', 'ion0.BextY', 'ion0.BextZ', 'ion0.B40',
+ 'ion0.B42', 'ion0.B44', 'ion0.B60', 'ion0.B62', 'ion0.B64', 'ion0.B66', 'ion0.IntensityScaling',
+ 'ion1.BmolX', 'ion1.BmolY', 'ion1.BmolZ', 'ion1.BextX', 'ion1.BextY', 'ion1.BextZ', 'ion1.B40',
+ 'ion1.B42', 'ion1.B44', 'ion1.B60', 'ion1.B62', 'ion1.B64', 'ion1.B66', 'ion1.IntensityScaling',
+ 'sp0.IntensityScaling', 'sp1.IntensityScaling')
+
+ chi2 = CalculateChiSquared(str(cf.function), InputWorkspace=ws1, InputWorkspace_1=ws2)[1]
+
+ fit = CrystalFieldFit(Model=cf, InputWorkspace=[ws1, ws2], MaxIterations=10)
+ fit.fit()
+
+ print(chi2)
+
+ cfms = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2', Temperatures=[25], FWHMs=[1.0], PeakShape='Gaussian',
+ BmolX=1.0, B40=-0.02)
+ print(str(cfms.function).split(',')[0])
+
+ # --------------------------
+
+ # Create some crystal field data
+ origin = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.12544,
+ Temperature=44.0, FWHM=1.1)
+ x, y = origin.getSpectrum()
+ ws = makeWorkspace(x, y)
+
+ # Define a CrystalField object with parameters slightly shifted.
+ cf = CrystalField('Ce', 'C2v', B20=0, B22=0, B40=0, B42=0, B44=0,
+ Temperature=44.0, FWHM=1.0, ResolutionModel=([0, 100], [1, 1]), FWHMVariation=0)
+
+ # Set any ties on the field parameters.
+ cf.ties(B20=0.37737)
+ # Create a fit object
+ fit = CrystalFieldFit(cf, InputWorkspace=ws)
+ # Find initial values for the field parameters.
+ # You need to define the energy splitting and names of parameters to estimate.
+ # Optionally additional constraints can be set on tied parameters (eg, peak centres).
+ fit.estimate_parameters(EnergySplitting=50,
+ Parameters=['B22', 'B40', 'B42', 'B44'],
+ Constraints='20<f1.PeakCentre<45,20<f2.PeakCentre<45',
+ NSamples=1000)
+ print('Returned', fit.get_number_estimates(), 'sets of parameters.')
+ # The first set (the smallest chi squared) is selected by default.
+ # Select a different parameter set if required
+ fit.select_estimated_parameters(3)
+ print(cf['B22'], cf['B40'], cf['B42'], cf['B44'])
+ # Run fit
+ fit.fit()
+
+ # --------------------------
+
+ from CrystalField import PointCharge
+ axial_pc_model = PointCharge([[-2, 0, 0, -4], [-2, 0, 0, 4]], 'Nd')
+ axial_blm = axial_pc_model.calculate()
+ print(axial_blm)
+
+ from CrystalField import PointCharge
+ from mantid.geometry import CrystalStructure
+ perovskite_structure = CrystalStructure('4 4 4 90 90 90', 'P m -3 m', 'Ce 0 0 0 1 0; Al 0.5 0.5 0.5 1 0; O 0.5 0.5 0 1 0')
+ cubic_pc_model = PointCharge(perovskite_structure, 'Ce', Charges={'Ce':3, 'Al':3, 'O':-2}, MaxDistance=7.5)
+
+ cubic_pc_model = PointCharge(perovskite_structure, 'Ce', Charges={'Ce':3, 'Al':3, 'O':-2}, Neighbour=2)
+ print(cubic_pc_model)
+
+ cif_pc_model = PointCharge('Sm2O3.cif')
+ print(cif_pc_model.getIons())
+
+ cif_pc_model.Charges = {'O1':-2, 'O2':-2, 'Sm1':3, 'Sm2':3, 'Sm3':3}
+ cif_pc_model.IonLabel = 'Sm2'
+ cif_pc_model.Neighbour = 1
+ cif_blm = cif_pc_model.calculate()
+ print(cif_blm)
+ bad_pc_model = PointCharge('Sm2O3.cif', MaxDistance=7.5, Neighbour=2)
+ print(bad_pc_model.Neighbour)
+ print(bad_pc_model.MaxDistance)
+
+ cif_pc_model.Charges = {'O':-2, 'Sm':3}
+ cif_blm = cif_pc_model.calculate()
+ print(cif_blm)
+
+ cf = CrystalField('Sm', 'C2', Temperature=5, FWHM=10, **cif_pc_model.calculate())
+ fit = CrystalFieldFit(cf, InputWorkspace=ws)
+
+ fit = CrystalFieldFit(cf, InputWorkspace=ws)
+ fit.fit()
+
+ # --------------------------
+
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, Temperature=44.0)
+ Cv = cf.getHeatCapacity() # Calculates Cv(T) for 1<T<300K in 1K steps (default)
+
+ T = np.arange(1,900,5)
+ Cv = cf.getHeatCapacity(T) # Calculates Cv(T) for specified values of T (1 to 900K in 5K steps here)
+
+ # Temperatures from a single spectrum workspace
+ ws = CreateWorkspace(T, T, T)
+ Cv = cf.getHeatCapacity(ws) # Use the x-values of a workspace as the temperatures
+ ws_cp = CreateWorkspace(*Cv)
+
+ # Temperatures from a multi-spectrum workspace
+ ws = CreateWorkspace(T, T, T, NSpec=2)
+ Cv = cf.getHeatCapacity(ws, 1) # Uses the second spectrum's x-values for T (e.g. 450<T<900)
+
+ chi_v = cf.getSusceptibility(T, Hdir=[1, 1, 1])
+ chi_v_powder = cf.getSusceptibility(T, Hdir='powder')
+ chi_v_cgs = cf.getSusceptibility(T, Hdir=[1, 1, 0], Unit='SI')
+ chi_v_bohr = cf.getSusceptibility(T, Unit='bohr')
+ print(type([chi_v, chi_v_powder, chi_v_cgs, chi_v_bohr]))
+ moment_t = cf.getMagneticMoment(Temperature=T, Hdir=[1, 1, 1], Hmag=0.1) # Calcs M(T) with at 0.1T field||[111]
+ H = np.linspace(0, 30, 121)
+ moment_h = cf.getMagneticMoment(Hmag=H, Hdir='powder', Temperature=10) # Calcs M(H) at 10K for powder sample
+ moment_SI = cf.getMagneticMoment(H, [1, 1, 1], Unit='SI') # M(H) in Am^2/mol at 1K for H||[111]
+ moment_cgs = cf.getMagneticMoment(100, Temperature=T, Unit='cgs') # M(T) in emu/mol in a field of 100G || [001]
+ print(type([moment_t, moment_h, moment_SI, moment_cgs]))
+
+ # --------------------------
+
+ from CrystalField import CrystalField, CrystalFieldFit, PhysicalProperties
+ # Fits a heat capacity dataset - you must have subtracted the phonon contribution by some method already
+ # and the data must be in J/mol/K.
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.12544,
+ PhysicalProperty=PhysicalProperties('Cv'))
+ fitcv = CrystalFieldFit(Model=cf, InputWorkspace=ws)
+ fitcv.fit()
+
+ params = {'B20':0.37737, 'B22':3.9770, 'B40':-0.031787, 'B42':-0.11611, 'B44':-0.12544}
+ cf = CrystalField('Ce', 'C2v', **params)
+ cf.PhysicalProperty = PhysicalProperties('Cv')
+ fitcv = CrystalFieldFit(Model=cf, InputWorkspace=ws)
+ fitcv.fit()
+
+ # Fits a susceptibility dataset. Data is the volume susceptibility in SI units
+ cf = CrystalField('Ce', 'C2v', **params)
+ cf.PhysicalProperty = PhysicalProperties('susc', Hdir='powder', Unit='SI')
+ fit_chi = CrystalFieldFit(Model=cf, InputWorkspace=ws)
+ fit_chi.fit()
+
+ # Fits a magnetisation dataset. Data is in emu/mol, and was measured at 5K with the field || [111].
+ cf = CrystalField('Ce', 'C2v', **params)
+ cf.PhysicalProperty = PhysicalProperties('M(H)', Temperature=5, Hdir=[1, 1, 1], Unit='cgs')
+ fit_mag = CrystalFieldFit(Model=cf, InputWorkspace=ws)
+ fit_mag.fit()
+
+ # Fits a magnetisation vs temperature dataset. Data is in Am^2/mol, measured with a 0.1T field || [110]
+ cf = CrystalField('Ce', 'C2v', **params)
+ cf.PhysicalProperty = PhysicalProperties('M(T)', Hmag=0.1, Hdir=[1, 1, 0], Unit='SI')
+ fit_moment = CrystalFieldFit(Model=cf, InputWorkspace=ws)
+ fit_moment.fit()
+
+ # --------------------------
+
+ # Pregenerate the required workspaces
+ for tt in [10, 44, 50]:
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.12544, Temperature=tt, FWHM=0.5)
+ x, y = cf.getSpectrum()
+ self.my_create_ws('ws_ins_'+str(tt)+'K', x, y)
+ ws_ins_10K = mtd['ws_ins_10K']
+ ws_ins_44K = mtd['ws_ins_44K']
+ ws_ins_50K = mtd['ws_ins_50K']
+ ws_cp = self.my_create_ws('ws_cp', *cf.getHeatCapacity())
+ ws_chi = self.my_create_ws('ws_chi', *cf.getSusceptibility(np.linspace(1,300,100), Hdir='powder', Unit='cgs'))
+ ws_mag = self.my_create_ws('ws_mag', *cf.getMagneticMoment(Hmag=np.linspace(0, 30, 100), Hdir=[1,1,1], Unit='bohr', Temperature=5))
+
+ # --------------------------
+
+ # Fits an INS spectrum (at 10K) and the heat capacity simultaneously
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.12544)
+ cf.Temperature = 10
+ cf.FWHM = 1.5
+ cf.PhysicalProperty = PhysicalProperties('Cv')
+ fit = CrystalFieldFit(Model=cf, InputWorkspace=[ws_ins_10K, ws_cp])
+ fit.fit()
+
+ # Fits two INS spectra (at 44K and 50K) and the heat capacity, susceptibility and magnetisation simultaneously.
+ PPCv = PhysicalProperties('Cv')
+ PPchi = PhysicalProperties('susc', 'powder', Unit='cgs')
+ PPMag = PhysicalProperties('M(H)', [1, 1, 1], 5, 'bohr')
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.12544,
+ Temperature=[44.0, 50], FWHM=[1.1, 0.9], PhysicalProperty=[PPCv, PPchi, PPMag] )
+ fit = CrystalFieldFit(Model=cf, InputWorkspace=[ws_ins_44K, ws_ins_50K, ws_cp, ws_chi, ws_mag])
+ fit.fit()
diff --git a/docs/source/fitfunctions/CrystalFieldSusceptibility.rst b/docs/source/fitfunctions/CrystalFieldSusceptibility.rst
index 9820e65ea8f47e8696a077b2073d188cf923a77c..f06644e26df6d1765bfdc3857a0e222771559408 100644
--- a/docs/source/fitfunctions/CrystalFieldSusceptibility.rst
+++ b/docs/source/fitfunctions/CrystalFieldSusceptibility.rst
@@ -31,14 +31,15 @@ be along the quantisation axis of the crystal field (which is usually defined to
:math:`B_x`, :math:`B_y`, and :math:`B_z` are the components of the unit vector pointing in the direction of the applied
magnetic field in this coordinate system.
-Finally, in order to account for the effect of any exchange interactions in the system which will shift the susceptibility curve
-up or down (analogous to the Curie-Weiss temperature), the actual magnetic susceptibility calculated by this function is:
+Finally, in order to account for the effect of any exchange interactions in the system which will shift the susceptiblity curve
+up or down (analogous to the Curie-Weiss temperature), and any residual (background) susceptibility in the sample (perhaps from
+an impurity), the actual magnetic susceptibility calculated by this function is:
-.. math:: \chi^{\mathrm{eff}} = \frac{\chi(T)}{1 - \lambda \chi(T)}
+.. math:: \chi^{\mathrm{eff}} = \frac{\chi(T)}{1 - \lambda \chi(T)} + \chi_0
-where :math:`\lambda` parameterises an effective exchange interaction and :math:`\chi` is the bare (paramagnetic Crystal Field)
-susceptibility. A negative :math:`\lambda` indicates overall antiferromagnetic interactions, whilst a positive :math:`\lambda`
-corresponds to overall ferromagnetic interactions.
+where :math:`\lambda` parameterises an effective exchange interaction with :math:`\chi` the bare (paramagnetic Crystal Field)
+susceptibility, and :math:`\chi_0` the residual susceptibility. A negative :math:`\lambda` indicates overall antiferromagnetic
+interactions, whilst a positive :math:`\lambda` corresponds to overall ferromagnetic interactions.
Example
-------
diff --git a/docs/source/interfaces/Crystal Field Python Interface.rst b/docs/source/interfaces/Crystal Field Python Interface.rst
index 461865af1b369b6d0a2ae94c0b93d8f7f2b09ec8..2dc1d1bf1a38283327bcbb37c8182965098d2b27 100644
--- a/docs/source/interfaces/Crystal Field Python Interface.rst
+++ b/docs/source/interfaces/Crystal Field Python Interface.rst
@@ -7,18 +7,18 @@ Crystal Field Python Interface
:local:
The python facilities for Crystal Field calculations are available in Mantid from module `CrystalField`.
-There are two main classes the module provides: `CrystalFiled` that defines various properties of a crystal
-field and `CrystalFieldFit` that manages the fitting process.
+The module provides two main classes: `CrystalField` defines various properties of a crystal field and
+`CrystalFieldFit` manages the fitting process.
Setting up crystal field parameters
-----------------------------------
A crystal field computation starts with creating an instance of the `CrystalField` class. The constructor
-has two mandatory arguments: `Ion` - a symbolic name of the ion, and `Symmetry` - a name of the symmetry group
-of the field. The rest of the parameters are optional.
+has two mandatory arguments: `Ion` - the symbolic name of the ion, and `Symmetry` - the name of the point symmetry
+group of the field. The rest of the parameters are optional.
-Possible values for the `Ion` argument::
+Possible values for the `Ion` argument are::
Ce, Pr, Nd, Pm, Sm, Eu, Gd, Tb, Dy, Ho, Er, Tm, Yb
@@ -30,7 +30,7 @@ or half-integer value, e.g. `Ion=S2` or `Ion=S1.5`. In these cases, the g-factor
The prefix letter can also be `J` instead of `S`, and lower case letters are also supported. (e.g. `Ion=j1`,
`Ion=s2.5` and `Ion=J0.5` are all valid).
-Allowed values for `Symmetry`::
+Allowed values for `Symmetry` are::
C1, Ci, C2, Cs, C2h, C2v, D2, D2h, C4, S4, C4h, D4, C4v, D2d, D4h, C3,
S6, D3, C3v, D3d, C6, C3h, C6h, D6, C6v, D3h, D6h, T, Td, Th, O, Oh
@@ -41,7 +41,7 @@ The minimum code to create a crystal field object is::
cf = CrystalField('Ce', 'C2v')
Names of the crystal field parameters have the form `Bnn` and `IBnn` where `nn` are two digits between 0 and 6.
-The `Bnn` is a real and the `IBnn` is an imaginary part of a complex parameter. If a parameter isn't set explicitly
+`Bnn` is the real and `IBnn` is the imaginary part of a complex parameter. If a parameter isn't set explicitly
its default value is 0. To set a parameter pass it to the `CrystalField` constructor as a keyword argument, e.g.::
cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770)
@@ -58,7 +58,7 @@ Which can also be used to query the value of a parameter::
Calculating the Eigensystem
---------------------------
-`CrystalField` class has methods to calculate the Hamiltonian and its eigensystem::
+The `CrystalField` class has methods to calculate the Hamiltonian and its eigensystem::
# Calculate and return the Hamiltonian matrix as a 2D numpy array.
h = cf.getHamiltonian()
@@ -91,7 +91,7 @@ Knowing the temperature allows us to calculate a peak list: a list of transition
print cf.getPeakList()
-The output is::
+Which produces the output::
[[ 0.00000000e+00 2.44006198e+01 4.24977124e+01 1.80970926e+01 -2.44006198e+01]
[ 2.16711565e+02 8.83098530e+01 5.04430056e+00 1.71153708e-01 1.41609425e-01]]
@@ -112,9 +112,9 @@ The new output::
[[ 0. 24.40061976 42.49771237]
[ 216.71156467 88.30985303 5.04430056]]
-To calculate a spectrum we need to define a shape of each peak (peak profile function) and its default width (`FWHM`).
+To calculate a spectrum we need to define the shape of each peak (peak profile function) and its default width (`FWHM`).
The width can be set either via a keyword argument or a property with name `FWHM`. If the peak shape isn't set the default
-of Lorentzian is assumed. To set a different shape use the `PeakShape` property::
+of `Lorentzian` is assumed. To set a different shape use the `PeakShape` property::
cf.PeakShape = 'Gaussian'
cf.FWHM = 0.9
@@ -128,8 +128,9 @@ After the peak shape is defined a spectrum can be calculated::
The output is a tuple of two 1d numpy arrays (x, y) that can be used with `matplotlib` to plot::
- pyplot.plot(*sp)
- pyplot.show()
+ import matplotlib.pyplot as plt
+ plt.plot(*sp)
+ plt.show()
.. image:: /images/CrystalFieldSpectrum1.png
:height: 300
@@ -170,7 +171,13 @@ Plotting in MantidPlot
----------------------
To plot a spectrum using MantidPlot's graphing facilities `CrystalField` has method `plot`. It has the same arguments as `getSpectrum`
-and opens a window with a plot.
+and opens a window with a plot, e.g.::
+
+ cf.plot()
+
+In addition to plotting, the `plot` method creates a workspace named `CrystalField_<Ion>` with the plot data. Subsequent calls to `plot`
+for the same `CrystalField` object will use the same plot window as created by the first call unless this window has been closed in the
+mean time.
Adding a Background
@@ -194,7 +201,7 @@ Setting Ties and Constraints
----------------------------
Setting ties and constraints are done by calling the `ties` and `constraints` methods of the `CrystalField` class or its components.
-To `Bnn` parameters are tied by the `CrystalField` class directly specifying the tied parameter as a keyword argument::
+The `Bnn` parameters are tied by the `CrystalField` class directly specifying the tied parameter as a keyword argument::
cf.ties(B20=1.0, B40='B20/2')
@@ -238,14 +245,16 @@ which is equivalent to::
Setting Resolution Model
------------------------
-Resolution model here is a way to constrain widths of the peaks to realistic numbers which agree with a measured or
+A resolution model is a way to constrain the widths of the peaks to realistic numbers which agree with a measured or
calculated instrument resolution function. A model is a function that returns a FWHM for a peak centre. The Crystal
-Field python interface defines helper class `ResolutionModel` to help define and set resolution models.
+Field python interface defines the helper class `ResolutionModel` to help define and set resolution models.
To construct an instance of `ResolutionModel` one needs to provide up to four input parameters. The first parameter, `model`, is
-mandatory and can be either of the two
+mandatory and can be either of:
-1. A tuple containing two arrays (lists) of real numbers which will be interpreted as tabulated values of the model function. The first element of the tuple is a list of increasing values for peak centres, and the second element is a list of corresponding widths. Values between the tabulated peak positions will be linearly interpolated.
+1. A tuple containing two arrays (lists) of real numbers which will be interpreted as tabulated values of the model function.
+ The first element of the tuple is a list of increasing values for peak centres, and the second element is a list of corresponding
+ widths. Values between the tabulated peak positions will be linearly interpolated.
2. A python function that takes a :class:`numpy.ndarray` of peak positions and returns a numpy array of widths.
@@ -254,21 +263,38 @@ function. `xstart` and `xend` define the interval of interpolation which must in
:math:`10^{-4}` and defines an approximate desired accuracy of the approximation. The interval will be split until the largest error of the interpolation
is smaller than `accuracy`. Note that subdivision cannot go on to infinity as the number of points is limited by the class member `ResolutionModel.max_model_size`.
-Example of setting a resolution model::
+Example of setting a resolution model using a tuple of two arrays::
+ from CrystalField import CrystalField, ResolutionModel
rm = ResolutionModel(([1, 2, 3, ...., 100], [0.1, 0.3, 0.35, ..., 2.1]))
cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, ..., Temperature=44.0, ResolutionModel=rm)
- ...
+Or using an arbitrary function `my_func`::
+
+ def my_func(en):
+ return (25-en)**(1.5) / 200 + 0.1
+
+ rm = ResolutionModel(my_func, xstart=0.0, xend=24.0, accuracy=0.01)
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, ..., Temperature=44.0, ResolutionModel=rm)
+
+Finally, the :ref:`PyChop` interface may be used to generate the resolution function for a particular spectrometer::
- rm = ResolutionModel(my_func, xstart=0.0, xend=120.0, accuracy=0.01)
+ from PyChop import PyChop2
+ marires = PyChop2('MARI')
+ marires.setChopper('S')
+ marires.setFrequency(250)
+ marires.setEi(30)
+ rm = ResolutionModel(marires.getResolution, xstart=0.0, xend=29.0, accuracy=0.01)
cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, ..., Temperature=44.0, ResolutionModel=rm)
-When a resolution model is set the peak width will be constrained to have a value close to the model. The degree of deviation is controlled by the
-`FWHMVariation` parameter. It has the default of 0.1 and is an absolute maximum difference a width can have. If set to 0 the widths will be fixed
-to their calculated values (depending on the instant values of their peak centres). For example::
+When a resolution model is set, the peak width will be constrained to have a value close to the model. The degree of deviation is controlled by the
+`FWHMVariation` parameter. It has the default of 0.1 and is the maximum difference from the value given by the resolution model a width can have.
+If set to 0 the widths will be fixed to their calculated values (depending on the instant values of their peak centres). For example::
- cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, ..., Temperature=44.0, ResolutionModel=rm, FWHMVariation=0.001)
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, ..., Temperature=44.0, ResolutionModel=rm, FWHMVariation=0.1)
+
+will allow the peak widths to vary between :math:`\Delta(E)-0.1` and :math:`\Delta(E)+0.1` where :math:`\Delta(E)` is the value of the
+resolution model at the peak position :math:`E`.
@@ -311,22 +337,22 @@ The resolution model also needs to be initialised from a list::
# or
rm = ResolutionModel([func0, func1], 0, 100, accuracy = 0.01)
-
+ cf.ResolutionModel = rm
To calculate a spectrum call the same method `getSpectrum` but pass the spectrum index as its first parameter::
- # Calculate second spectrum, use the generated x-values
- sp = cf.getSpectrum(1)
+ # Calculate second spectrum, use the generated x-values
+ sp = cf.getSpectrum(1)
- # Calculate third spectrum, use a list for x-values
- x = [0, 1, 2, 3, ...]
- sp = cf.getSpectrum(2, x)
-
- # Calculate second spectrum, use the first spectrum of a workspace
- sp = cf.getSpectrum(1, ws)
-
- # Calculate first spectrum, use the i-th spectrum of a workspace
- sp = cf.getSpectrum(0, ws, i)
+ # Calculate third spectrum, use a list for x-values
+ x = [0, 1, 2, 3, ...]
+ sp = cf.getSpectrum(2, x)
+
+ # Calculate second spectrum, use the first spectrum of a workspace
+ sp = cf.getSpectrum(1, ws)
+
+ # Calculate first spectrum, use the i-th spectrum of a workspace
+ sp = cf.getSpectrum(0, ws, i)
Note that the attributes `Temperature`, `FWHM`, `peaks` and `background` may be set separately from the constructor, e.g.::
@@ -380,9 +406,11 @@ means that the intensity of `cf1` should be twice that of `cf2`.
Alternatively, you can create a `CrystalFieldMultiSite` object directly. This takes Ions, Symmetries, Temperatures and peak widths as lists::
+ from CrystalField import CrystalFieldMultiSite
cfms = CrystalFieldMultiSite(Ions=['Ce', 'Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44.0], FWHMs=[1.1])
-To access parameters of a CrystalFieldMultiSite object, prefix with the ion index::
+Note that `Temperature` and `FWHM` (without plural) can also be used in place of the equivalent plural parameters.
+To access parameters of a CrystalFieldMultiSite object, prefix them with the ion index::
cfms['ion0.B40'] = -0.031
cfms['ion1.B20'] = 0.37737
@@ -393,7 +421,7 @@ Parameters can be set when creating the object by passing in a dictionary using
cfms = CrystalFieldMultiSite(Ions=['Ce', 'Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44.0], FWHMs=[1.1],
parameters={'ion0.B20': 0.37737, 'ion0.B22': 3.9770, 'ion1.B40':-0.031787,
- 'ion1.B42':-0.11611, 'ion1.B44':-0.12544})
+ 'ion1.B42':-0.11611, 'ion1.B44':-0.12544})
A background can also be set this way, or using `cfms.background.` It can be passed as a string, a Function object(s), or a
CompositeFunction object::
@@ -407,7 +435,7 @@ CompositeFunction object::
cfms = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHMs=[1.0],
Background= Gaussian(PeakCentre=1) + LinearBackground())
-Ties and constraints are set similiarly to `CrystalField` objects. `f` prefixes have been changed to be more descriptive::
+Ties and constraints are set similarly to `CrystalField` objects. `f` prefixes have been changed to be more descriptive::
cfms = CrystalFieldMultiSite(Ions=['Ce','Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44, 50], FWHMs=[1.1, 0.9],
Background=FlatBackground(), BackgroundPeak=Gaussian(Height=10, Sigma=0.3),
@@ -418,8 +446,11 @@ Ties and constraints are set similiarly to `CrystalField` objects. `f` prefixes
cfms.constraints('ion0.sp0.pk1.FWHM < 2.2')
cfms.ties({'ion0.sp1.pk2.FWHM': '2*ion0.sp1.pk1.FWHM', 'ion1.sp1.pk3.FWHM': '2*ion1.sp1.pk2.FWHM'})
-When fitting, all parameters are assumed to be free. Parameters must be fixed explicitly. Fitting example::
+Parameters which are not allowed by the specified symmetry will be fixed to be zero, but unlike for the single-site case,
+all other parameters are assumed to be free (in the single-site case, parameters which are unset are assumed to be fixed
+to be zero). For the multi-site case, parameters must be fixed explicitly. For example::
+ params = {'ion0.B20': 0.37737, 'ion0.B22': 3.9770, 'ion1.B40':-0.031787, 'ion1.B42':-0.11611, 'ion1.B44':-0.12544}
cf = CrystalFieldMultiSite(Ions=['Ce', 'Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44.0, 50.0],
FWHMs=[1.0, 1.0], ToleranceIntensity=6.0, ToleranceEnergy=1.0, FixAllPeaks=True,
parameters=params)
@@ -567,14 +598,14 @@ gives::
{'O1': [0.125, 0.125, 0.375],
'O2': [0.125, 0.375, 0.375],
- 'Yb1': [0.25, 0.25, 0.25],
- 'Yb2': [0.021, 0.0, 0.25],
- 'Yb3': [0.542, 0.0, 0.25]}
+ 'Sm1': [0.25, 0.25, 0.25],
+ 'Sm2': [0.021, 0.0, 0.25],
+ 'Sm3': [0.542, 0.0, 0.25]}
You can then define the charges for each site, the magnetic ion and the maximum distance, and calculate::
- cif_pc_model.Charges = {'O1':-2, 'O2':-2, 'Yb1':3, 'Yb2':3, 'Yb3':3}
- cif_pc_model.IonLabel = 'Yb2'
+ cif_pc_model.Charges = {'O1':-2, 'O2':-2, 'Sm1':3, 'Sm2':3, 'Sm3':3}
+ cif_pc_model.IonLabel = 'Sm2'
cif_pc_model.Neighbour = 1
cif_blm = cif_pc_model.calculate()
print(cif_blm)
@@ -591,14 +622,14 @@ then the value for ``MaxDistance`` will be used regardless of where it appears i
For ``Charges``, instead of listing the charges of each site, you can just give the charge for each element, e.g.::
- cif_pc_model.Charges = {'O':-2, 'Yb':3}
+ cif_pc_model.Charges = {'O':-2, 'Sm':3}
cif_blm = cif_pc_model.calculate()
The result of the ``calculate()`` method can be put directly into a ``CrystalField`` object and used either
to calculate a spectrum or as the starting parameters in a fit::
- cf = CrystalField('Yb', 'C2', Temperature=5, FWHM=10, **cif_pc_model.calculate())
- plot(**cf.getSpectrum())
+ cf = CrystalField('Sm', 'C2', Temperature=5, FWHM=10, **cif_pc_model.calculate())
+ plot(*cf.getSpectrum())
fit = CrystalFieldFit(cf, InputWorkspace=ws)
fit.fit()
@@ -620,12 +651,14 @@ susceptibility, and magnetisation. The calculated values can be invoked using th
To calculate the heat capacity use::
+ import matplotlib.pyplot as plt
+ cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, Temperature=44.0)
Cv = cf.getHeatCapacity() # Calculates Cv(T) for 1<T<300K in 1K steps (default)
- pyplot.plot(*Cv) # Returns a tuple of (x, y) values
+ plt.plot(*Cv) # Returns a tuple of (x, y) values
T = np.arange(1,900,5)
Cv = cf.getHeatCapacity(T) # Calculates Cv(T) for specified values of T (1 to 900K in 5K steps here)
- pyplot.plot(T, Cv[1])
+ plt.plot(T, Cv[1])
# Temperatures from a single spectrum workspace
ws = CreateWorkspace(T, T, T)
@@ -710,22 +743,26 @@ class as the `PhysicalProperty` attribute of `CrystalField`, either as a keyword
or separately after construction::
- cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.12544)
+ params = {'B20':0.37737, 'B22':3.9770, 'B40':-0.031787, 'B42':-0.11611, 'B44':-0.12544}
+ cf = CrystalField('Ce', 'C2v', **params)
cf.PhysicalProperty = PhysicalProperties('Cv')
fitcv = CrystalFieldFit(Model=cf, InputWorkspace=ws)
fitcv.fit()
# Fits a susceptibility dataset. Data is the volume susceptibility in SI units
+ cf = CrystalField('Ce', 'C2v', **params)
cf.PhysicalProperty = PhysicalProperties('susc', Hdir='powder', Unit='SI')
fit_chi = CrystalFieldFit(Model=cf, InputWorkspace=ws)
fit_chi.fit()
# Fits a magnetisation dataset. Data is in emu/mol, and was measured at 5K with the field || [111].
+ cf = CrystalField('Ce', 'C2v', **params)
cf.PhysicalProperty = PhysicalProperties('M(H)', Temperature=5, Hdir=[1, 1, 1], Unit='cgs')
fit_mag = CrystalFieldFit(Model=cf, InputWorkspace=ws)
fit_mag.fit()
# Fits a magnetisation vs temperature dataset. Data is in Am^2/mol, measured with a 0.1T field || [110]
+ cf = CrystalField('Ce', 'C2v', **params)
cf.PhysicalProperty = PhysicalProperties('M(T)', Hmag=0.1, Hdir=[1, 1, 0], Unit='SI')
fit_moment = CrystalFieldFit(Model=cf, InputWorkspace=ws)
fit_moment.fit()
@@ -754,9 +791,9 @@ properties should be a list. E.g.::
fit.fit()
# Fits two INS spectra (at 44K and 50K) and the heat capacity, susceptibility and magnetisation simultaneously.
- PPCv = PhysicalProperty('Cv')
- PPchi = PhysicalProperty('susc', 'powder', Unit='cgs')
- PPMag = PhysicalProperty('M(H)', 5, [1, 1, 1], 'bohr')
+ PPCv = PhysicalProperties('Cv')
+ PPchi = PhysicalProperties('susc', 'powder', Unit='cgs')
+ PPMag = PhysicalProperties('M(H)', [1, 1, 1], 5, 'bohr')
cf = CrystalField('Ce', 'C2v', B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.12544,
Temperature=[44.0, 50], FWHM=[1.1, 0.9], PhysicalProperty=[PPCv, PPchi, PPMag] )
fit = CrystalFieldFit(Model=cf, InputWorkspace=[ws_ins_44K, ws_ins_50K, ws_cp, ws_chi, ws_mag])
@@ -764,7 +801,7 @@ properties should be a list. E.g.::
Note that `PhysicalProperty` requires the type of physical property (either `'Cv'` or `'Cp'` or `'heatcap'` for
heat capacity; `'susc'` or `'chi'` for susceptibility; `'mag'` or `'M(H)'` for magnetic moment vs applied field;
-or `'mom'` or `'M(T)'` for moment vs temperature) as the first argument. subsequent arguments are optional, and
+or `'mom'` or `'M(T)'` for moment vs temperature) as the first argument. Subsequent arguments are optional, and
are in the following order::
PhysicalProperties('Cp') # No further parameters required for heat capacity
diff --git a/docs/source/release/v3.12.0/direct_inelastic.rst b/docs/source/release/v3.12.0/direct_inelastic.rst
index 1949d7a42fb508f15d9e3e929330f700345edfdf..3db959ac18f8b2ae6732015484109151c990c6c8 100644
--- a/docs/source/release/v3.12.0/direct_inelastic.rst
+++ b/docs/source/release/v3.12.0/direct_inelastic.rst
@@ -27,6 +27,12 @@ Algorithms
Crystal Field
#############
+Multi-site calculations and fitting are now supported by the crystal field (Python commandline) interface.
+
+Calculation of dipole transition matrix elements has been added, together with the addition of a :math:`\chi_0` term in the :ref:`CrystalFieldSusceptibility <func-CrystalFieldSusceptibility>` function.
+
+Several bugs in the Python and C++ code has been fixed - see the `github page <https://github.com/mantidproject/mantid/pull/21604>`_ for details.
+
Features Removed
----------------
diff --git a/docs/source/release/v3.12.0/spectroscopy.rst b/docs/source/release/v3.12.0/spectroscopy.rst
index c327a21767bbd88294b14be30836de85e2202f16..e19120bf466500cda86a7a1eb4fb099ed38830c6 100644
--- a/docs/source/release/v3.12.0/spectroscopy.rst
+++ b/docs/source/release/v3.12.0/spectroscopy.rst
@@ -11,6 +11,7 @@ Spectroscopy Changes
- The algorithms :ref:`algm-SofQWCentre`, :ref:`algm-SofQWPolygon` and :ref:`algm-SofQWNormalisedPolygon`, which rebin an inelastic workspace (has a `DeltaE` axis) from spectrum numbers (angle) to `MomentumTransfer` may now rebin the energy (`DeltaE`) axis as well as the :math:`|Q|` (`MomentumTransfer`) axes.
- :ref:`algm-SofQWNormalisedPolygon` now has uses a faster method for calculating the polygon intersections.
+- The crystal field computation and fitting engine is now feature complete. It can now handle multi-site computation and simultaneous fitting of inelastic spectra and physical properties dataset. See the :ref:`Crystal Field Python Interface` help page for details, and `<http://www.mantidproject.org/Crystal_Field_Examples>`_ for examples of use.
Direct Geometry
---------------
diff --git a/scripts/Inelastic/CrystalField/CrystalFieldMultiSite.py b/scripts/Inelastic/CrystalField/CrystalFieldMultiSite.py
index 5c693ea580b2229c067402a52b5db65e63ba5541..2eb89e593e877772a820f123ec393f6244c3a1d2 100644
--- a/scripts/Inelastic/CrystalField/CrystalFieldMultiSite.py
+++ b/scripts/Inelastic/CrystalField/CrystalFieldMultiSite.py
@@ -1,6 +1,7 @@
import numpy as np
from CrystalField import CrystalField, Function
+from .fitting import islistlike
def makeWorkspace(xArray, yArray, child=True, ws_name='dummy'):
@@ -66,9 +67,13 @@ class CrystalFieldMultiSite(object):
self.Symmetries = Symmetries
self._plot_window = {}
self.chi2 = None
+ self._resolutionModel = None
parameter_dict = kwargs.pop('parameters', None)
attribute_dict = kwargs.pop('attributes', None)
+ ties_dict = kwargs.pop('ties', None)
+ constraints_list = kwargs.pop('constraints', None)
+ fix_list = kwargs.pop('fixedParameters', None)
kwargs = self._setMandatoryArguments(kwargs)
@@ -78,22 +83,32 @@ class CrystalFieldMultiSite(object):
self._setRemainingArguments(kwargs)
+ self.default_spectrum_size = 200
+
if attribute_dict is not None:
for name, value in attribute_dict.items():
self.function.setAttributeValue(name, value)
if parameter_dict is not None:
for name, value in parameter_dict.items():
self.function.setParameter(name, value)
+ if ties_dict:
+ for name, value in parameter_dict.items():
+ self.function.tie(name, value)
+ if constraints_list:
+ self.function.addConstraints(','.join(constraints_list))
+ if fix_list:
+ for param in fix_list:
+ self.function.fixParameter(param)
def _setMandatoryArguments(self, kwargs):
- if 'Temperatures' in kwargs:
- self.Temperatures = kwargs.pop('Temperatures')
- if 'FWHMs' in kwargs:
- self.FWHMs = kwargs.pop('FWHMs')
+ if 'Temperatures' in kwargs or 'Temperature' in kwargs:
+ self.Temperatures = kwargs.pop('Temperatures') if 'Temperatures' in kwargs else kwargs.pop('Temperature')
+ if 'FWHM' in kwargs or 'FWHMs' in kwargs:
+ self.FWHM = kwargs.pop('FWHMs') if 'FWHMs' in kwargs else kwargs.pop('FWHM')
elif 'ResolutionModel' in kwargs:
self.ResolutionModel = kwargs.pop('ResolutionModel')
else:
- raise RuntimeError("If temperatures are set, must also set FWHMs or ResolutionModel")
+ raise RuntimeError("If temperatures are set, must also set FWHM or ResolutionModel")
return kwargs
def _setRemainingArguments(self, kwargs):
@@ -137,6 +152,16 @@ class CrystalFieldMultiSite(object):
return self._calcSpectrum(i, ws, ws_index)
+ def calc_xmin_xmax(self, i=0):
+ peaks = np.array([])
+ for idx in range(len(self.Ions)):
+ blm = {}
+ for bparam in CrystalField.field_parameter_names:
+ blm[bparam] = self.function.getParameterValue('ion{}.'.format(idx) + bparam)
+ _cft = CrystalField(self.Ions[idx], 'C1', Temperature=self.Temperatures[i], **blm)
+ peaks = np.append(peaks, _cft.getPeakList()[0])
+ return np.min(peaks), np.max(peaks)
+
def getSpectrum(self, *args):
"""
Get a specified spectrum calculated with the current field and peak parameters.
@@ -144,14 +169,16 @@ class CrystalFieldMultiSite(object):
Alternatively can be called getSpectrum(workspace, ws_index). Spectrum index is assumed zero.
Examples:
+ cf.getSpectrum() # Calculate the first spectrum using automatically generated x-values
+ cf.getSpectrum(1) # Calculate the second spectrum using automatically generated x-values
cf.getSpectrum(1, ws, 5) # Calculate the second spectrum using the x-values from the 6th spectrum
# in workspace ws.
- cf.getSpectrum(ws) # Calculate the first spectrum using the x-values from the 1st spectrum
- # in workspace ws.
- cf.getSpectrum(ws, 3) # Calculate the first spectrum using the x-values from the 4th spectrum
- # in workspace ws.
- cf.getSpectrum(3, ws) # Calculate the third spectrum using the x-values from the 1st spectrum
- # in workspace ws.
+ cf.getSpectrum(ws) # Calculate the first spectrum using the x-values from the 1st spectrum
+ # in workspace ws.
+ cf.getSpectrum(ws, 3) # Calculate the first spectrum using the x-values from the 4th spectrum
+ # in workspace ws.
+ cf.getSpectrum(2, ws) # Calculate the third spectrum using the x-values from the 1st spectrum
+ # in workspace ws.
@return: A tuple of (x, y) arrays
"""
@@ -160,11 +187,18 @@ class CrystalFieldMultiSite(object):
raise RuntimeError('You must first define a temperature for the spectrum')
return self._calcSpectrum(args[0], args[1], args[2])
elif len(args) == 1:
- return self._calcSpectrum(0, args[0], 0)
+ if isinstance(args[0], int):
+ x_min, x_max = self.calc_xmin_xmax(args[0])
+ xArray = np.linspace(x_min, x_max, self.default_spectrum_size)
+ return self._calcSpectrum(args[0], xArray, 0)
+ else:
+ return self._calcSpectrum(0, args[0], 0)
elif len(args) == 2:
return self._getSpectrumTwoArgs(*args)
else:
- raise RuntimeError('getSpectrum expected 1-3 arguments, got {}s'.format(len(args)))
+ x_min, x_max = self.calc_xmin_xmax()
+ xArray = np.linspace(x_min, x_max, self.default_spectrum_size)
+ return self._calcSpectrum(0, xArray, 0)
def _convertToWS(self, wksp_list):
"""
@@ -350,7 +384,7 @@ class CrystalFieldMultiSite(object):
params = get_parameters_for_add_from_multisite(self, 0)
params.update(get_parameters_for_add_from_multisite(other, len(self.Ions)))
new_cf = CrystalFieldMultiSite(Ions=ions, Symmetries=symmetries, Temperatures=self.Temperatures,
- FWHMs=self.FWHMs, parameters=params, abundances=abundances)
+ FWHM=self.FWHM, parameters=params, abundances=abundances)
return new_cf
def __getitem__(self, item):
@@ -378,7 +412,7 @@ class CrystalFieldMultiSite(object):
params = get_parameters_for_add_from_multisite(self, 0)
params.update(get_parameters_for_add(other, len(self.Ions)))
new_cf = CrystalFieldMultiSite(Ions=ions, Symmetries=symmetries, Temperatures=self.Temperatures,
- FWHMs=self.FWHMs, parameters=params, abundances=abundances)
+ FWHM=self.FWHM, parameters=params, abundances=abundances)
return new_cf
def __radd__(self, other):
@@ -395,13 +429,31 @@ class CrystalFieldMultiSite(object):
params = get_parameters_for_add(other, 0)
params.update(get_parameters_for_add_from_multisite(self, 1))
new_cf = CrystalFieldMultiSite(Ions=ions, Symmetries=symmetries, Temperatures=self.Temperatures,
- FWHMs=self.FWHMs, parameters=params, abundances=abundances)
+ FWHM=self.FWHM, parameters=params, abundances=abundances)
return new_cf
@property
def background(self):
return self._background
+ @background.setter
+ def background(self, value):
+ if hasattr(value, 'peak') and hasattr(value, 'background'):
+ # Input is a CrystalField.Background object
+ if value.peak and value.background:
+ self._setBackground(peak=str(value.peak.function), background=str(value.background.function))
+ elif value.peak:
+ self._setBackground(peak=str(value.peak.function))
+ else:
+ self._setBackground(background=str(value.background.function))
+ elif hasattr(value, 'function'):
+ self._setBackground(background=str(value.function))
+ else:
+ self._setBackground(background=value)
+ # Need this for a weird python bug: "IndexError: Function index (2) out of range (2)"
+ # if user calls print(self.function) after setting background
+ _ = self.function.getTies() # noqa: F841
+
@property
def Ions(self):
string_ions = self.function.getAttributeValue('Ions')
@@ -456,6 +508,14 @@ class CrystalFieldMultiSite(object):
def Temperatures(self, value):
self.function.setAttributeValue('Temperatures', value)
+ @property
+ def Temperature(self):
+ return list(self.function.getAttributeValue("Temperatures"))
+
+ @Temperature.setter
+ def Temperatures(self, value):
+ self.function.setAttributeValue('Temperatures', value)
+
@property
def FWHMs(self):
fwhm = self.function.getAttributeValue('FWHMs')
@@ -466,10 +526,46 @@ class CrystalFieldMultiSite(object):
@FWHMs.setter
def FWHMs(self, value):
- if len(value) == 1:
- value = value[0]
+ if islistlike(value):
+ if len(value) != len(self.Temperatures):
+ value = [value[0]] * len(self.Temperatures)
+ else:
value = [value] * len(self.Temperatures)
self.function.setAttributeValue('FWHMs', value)
+ self._resolutionModel = None
+
+ @property
+ def FWHM(self):
+ return self.FWHMs
+
+ @FWHM.setter
+ def FWHM(self, value):
+ self.FWHMs = value
+
+ @property
+ def ResolutionModel(self):
+ return self._resolutionModel
+
+ @ResolutionModel.setter
+ def ResolutionModel(self, value):
+ from .function import ResolutionModel
+ if hasattr(value, 'model'):
+ self._resolutionModel = value
+ else:
+ self._resolutionModel = ResolutionModel(value)
+ nSpec = len(self.Temperatures)
+ if nSpec > 1:
+ if not self._resolutionModel.multi or self._resolutionModel.NumberOfSpectra != nSpec:
+ raise RuntimeError('Resolution model is expected to have %s functions, found %s' %
+ (nSpec, self._resolutionModel.NumberOfSpectra))
+ for i in range(nSpec):
+ model = self._resolutionModel.model[i]
+ self.function.setAttributeValue('sp%i.FWHMX' % i, model[0])
+ self.function.setAttributeValue('sp%i.FWHMY' % i, model[1])
+ else:
+ model = self._resolutionModel.model
+ self.function.setAttributeValue('FWHMX', model[0])
+ self.function.setAttributeValue('FWHMY', model[1])
@property
def FWHMVariation(self):
diff --git a/scripts/Inelastic/CrystalField/__init__.py b/scripts/Inelastic/CrystalField/__init__.py
index 17cf614835b81aeb3b6b1c7d4b8bafe7bc7a5449..20fd13f37bb8abddefd4308b4d230dd7d7264fc4 100644
--- a/scripts/Inelastic/CrystalField/__init__.py
+++ b/scripts/Inelastic/CrystalField/__init__.py
@@ -1,6 +1,7 @@
from __future__ import (absolute_import, division, print_function)
-from .fitting import CrystalField, CrystalFieldFit, CrystalFieldMulti
+from .fitting import CrystalField, CrystalFieldFit
from .function import PeaksFunction, Background, Function, ResolutionModel, PhysicalProperties
from .pointcharge import PointCharge
-__all__ = ['CrystalField', 'CrystalFieldFit', 'CrystalFieldMulti', 'PeaksFunction',
+from .CrystalFieldMultiSite import CrystalFieldMultiSite
+__all__ = ['CrystalField', 'CrystalFieldFit', 'CrystalFieldMultiSite', 'PeaksFunction',
'Background', 'Function', 'ResolutionModel', 'PhysicalProperties', 'PointCharge']
diff --git a/scripts/Inelastic/CrystalField/fitting.py b/scripts/Inelastic/CrystalField/fitting.py
index aab99ebca585ddecb1347b76fc8a463a861b26af..f07526d4cfd5b9ca3a1a6fed0a54e474269757c4 100644
--- a/scripts/Inelastic/CrystalField/fitting.py
+++ b/scripts/Inelastic/CrystalField/fitting.py
@@ -2,7 +2,7 @@ from __future__ import (absolute_import, division, print_function)
import numpy as np
import re
import warnings
-from six import string_types
+from six import string_types, iteritems
# RegEx pattern matching a composite function parameter name, eg f2.Sigma.
@@ -26,20 +26,68 @@ def makeWorkspace(xArray, yArray):
def islistlike(arg):
- return (not hasattr(arg, "strip")) and (hasattr(arg, "__getitem__") or hasattr(arg, "__iter__"))
+ return (not hasattr(arg, "strip")) and (hasattr(arg, "__getitem__") or hasattr(arg, "__iter__")) and hasattr(arg, "__len__")
+
+
+def ionname2Nre(ionname):
+ ion_nre_map = {'Ce': 1, 'Pr': 2, 'Nd': 3, 'Pm': 4, 'Sm': 5, 'Eu': 6, 'Gd': 7,
+ 'Tb': 8, 'Dy': 9, 'Ho': 10, 'Er': 11, 'Tm': 12, 'Yb': 13}
+ if ionname not in ion_nre_map.keys():
+ msg = 'Value %s is not allowed for attribute Ion.\nList of allowed values: %s' % \
+ (ionname, ', '.join(list(ion_nre_map.keys())))
+ arbitraryJ = re.match('[SJsj]([0-9\.]+)', ionname)
+ if arbitraryJ and (float(arbitraryJ.group(1)) % 0.5) == 0:
+ nre = int(-float(arbitraryJ.group(1)) * 2.)
+ if nre < -99:
+ raise RuntimeError('J value ' + str(-nre / 2) + ' is too large.')
+ else:
+ raise RuntimeError(msg+', S<n>, J<n>')
+ else:
+ nre = ion_nre_map[ionname]
+ return nre
+
+
+def cfpstrmaker(x, pref='B'):
+ return [pref+str(k)+str(x) for k in [2, 4, 6] if x <= k]
+
+
+def getSymmAllowedParam(sym_str):
+ if 'T' in sym_str or 'O' in sym_str:
+ return ['B40', 'B44', 'B60', 'B64']
+ if any([sym_str == val for val in ['C1', 'Ci']]):
+ return sum([cfpstrmaker(i) for i in range(7)] +
+ [cfpstrmaker(i, 'IB') for i in range(1, 7)],[])
+ retval = cfpstrmaker(0)
+ if '6' in sym_str or '3' in sym_str:
+ retval += cfpstrmaker(6)
+ if any([sym_str == val for val in ['C6', 'C3h', 'C6h']]):
+ retval += cfpstrmaker(6, 'IB')
+ if ('3' in sym_str and '3h' not in sym_str) or 'S6' in sym_str:
+ retval += cfpstrmaker(3)
+ if any([sym_str == val for val in ['C3', 'S6']]):
+ retval += cfpstrmaker(3, 'IB') + cfpstrmaker(6, 'IB')
+ if '4' in sym_str or '2' in sym_str:
+ retval += cfpstrmaker(4)
+ if any([sym_str == val for val in ['C4', 'S4', 'C4h']]):
+ retval += cfpstrmaker(4, 'IB')
+ if ('2' in sym_str and '2d' not in sym_str) or 'Cs' in sym_str:
+ retval += cfpstrmaker(2)
+ if any([sym_str == val for val in ['C2', 'Cs', 'C2h']]):
+ retval += cfpstrmaker(2, 'IB') + cfpstrmaker(4, 'IB')
+ return retval
#pylint: disable=too-many-instance-attributes,too-many-public-methods
class CrystalField(object):
"""Calculates the crystal fields for one ion"""
- ion_nre_map = {'Ce': 1, 'Pr': 2, 'Nd': 3, 'Pm': 4, 'Sm': 5, 'Eu': 6, 'Gd': 7,
- 'Tb': 8, 'Dy': 9, 'Ho': 10, 'Er': 11, 'Tm': 12, 'Yb': 13}
-
allowed_symmetries = ['C1', 'Ci', 'C2', 'Cs', 'C2h', 'C2v', 'D2', 'D2h', 'C4', 'S4', 'C4h',
'D4', 'C4v', 'D2d', 'D4h', 'C3', 'S6', 'D3', 'C3v', 'D3d', 'C6', 'C3h',
'C6h', 'D6', 'C6v', 'D3h', 'D6h', 'T', 'Td', 'Th', 'O', 'Oh']
+ lande_g = [6.0 / 7., 4.0 / 5., 8.0 / 11., 3.0 / 5., 2.0 / 7., 0.0, 2.0,
+ 3.0 / 2., 4.0 / 3., 5.0 / 4., 6.0 / 5., 7.0 / 6., 8.0 / 7.]
+
default_peakShape = 'Gaussian'
default_background = 'FlatBackground'
default_spectrum_size = 200
@@ -132,6 +180,12 @@ class CrystalField(object):
free_parameters = {key: kwargs[key] for key in kwargs if key in CrystalField.field_parameter_names}
+ if 'ResolutionModel' in kwargs and 'FWHM' in kwargs:
+ msg = 'Both ''ResolutionModel'' and ''FWHM'' specified but can only accept one width option.'
+ msg += ' Prefering to use ResolutionModel, and ignoring FWHM.'
+ kwargs.pop('FWHM')
+ warnings.warn(msg, SyntaxWarning)
+
for key in kwargs:
if key == 'ToleranceEnergy':
self.ToleranceEnergy = kwargs[key]
@@ -157,7 +211,7 @@ class CrystalField(object):
for param in CrystalField.field_parameter_names:
if param in free_parameters:
self.function.setParameter(param, free_parameters[param])
- else:
+ elif param not in getSymmAllowedParam(self.Symmetry):
self.function.fixParameter(param)
self._setPeaks()
@@ -173,8 +227,6 @@ class CrystalField(object):
self._peakList = None
# Spectra
- self._dirty_spectra = True
- self._spectra = {}
self._plot_window = {}
# self._setDefaultTies()
@@ -248,6 +300,7 @@ class CrystalField(object):
"""
if hasattr(i, 'toString'):
out = i.toString()
+ ppobj = i
else:
if self._physprop is None:
raise RuntimeError('Physical properties environment not defined.')
@@ -261,6 +314,11 @@ class CrystalField(object):
if len(fieldParams) > 0:
out += ',%s' % ','.join(['%s=%s' % item for item in fieldParams.items()])
ties = self._getFieldTies()
+ from .function import PhysicalProperties
+ if ppobj.TypeID == PhysicalProperties.SUSCEPTIBILITY:
+ ties += ',' if ties else ''
+ ties += 'Lambda=0' if ppobj.Lambda == 0. else ''
+ ties += ',Chi0=0' if ppobj.Chi0 == 0. else ''
if len(ties) > 0:
out += ',ties=(%s)' % ties
constraints = self._getFieldConstraints()
@@ -269,7 +327,6 @@ class CrystalField(object):
return out
def makeMultiSpectrumFunction(self):
- import re
return re.sub(r'FWHM[X|Y]\d+=\(\),', '', str(self.function))
@property
@@ -290,14 +347,10 @@ class CrystalField(object):
...
cf.Ion = 'Pr'
"""
- if value not in self.ion_nre_map.keys():
- msg = 'Value %s is not allowed for attribute Ion.\nList of allowed values: %s' % \
- (value, ', '.join(list(self.ion_nre_map.keys())))
- raise RuntimeError(msg)
+ self._nre = ionname2Nre(value)
self.crystalFieldFunction.setAttributeValue('Ion', value)
self._dirty_eigensystem = True
self._dirty_peaks = True
- self._dirty_spectra = True
@property
def Symmetry(self):
@@ -324,7 +377,6 @@ class CrystalField(object):
self.crystalFieldFunction.setAttributeValue('Symmetry', value)
self._dirty_eigensystem = True
self._dirty_peaks = True
- self._dirty_spectra = True
@property
def ToleranceEnergy(self):
@@ -336,7 +388,6 @@ class CrystalField(object):
"""Set energy tolerance"""
self.crystalFieldFunction.setAttributeValue('ToleranceEnergy', float(value))
self._dirty_peaks = True
- self._dirty_spectra = True
@property
def ToleranceIntensity(self):
@@ -348,7 +399,6 @@ class CrystalField(object):
"""Set intensity tolerance"""
self.crystalFieldFunction.setAttributeValue('ToleranceIntensity', float(value))
self._dirty_peaks = True
- self._dirty_spectra = True
@property
def IntensityScaling(self):
@@ -378,7 +428,6 @@ class CrystalField(object):
self.crystalFieldFunction.setParameter(paramName, value[i])
self._dirty_peaks = True
- self._dirty_spectra = True
@property
def Temperature(self):
@@ -401,7 +450,6 @@ class CrystalField(object):
return
self.crystalFieldFunction.setAttributeValue('Temperature', float(value))
self._dirty_peaks = True
- self._dirty_spectra = True
@property
def FWHM(self):
@@ -420,12 +468,30 @@ class CrystalField(object):
if self._isMultiSpectrum:
if not islistlike(value):
value = [value] * self.NumberOfSpectra
+ if len(value) != len(self.Temperature):
+ if self.PhysicalProperty is not None and len(value) == len(self.Temperature) - len(self.PhysicalProperty):
+ value = value + [0] * len(self.PhysicalProperty)
+ else:
+ raise RuntimeError('Vector of FWHMs must either have same size as '
+ 'Temperatures (%i) or have size 1.' % (len(self.Temperature)))
self.crystalFieldFunction.setAttributeValue('FWHMs', value)
else:
if islistlike(value):
raise ValueError('FWHM is expected to be a single floating point value')
self.crystalFieldFunction.setAttributeValue('FWHM', float(value))
- self._dirty_spectra = True
+ # If both FWHM and ResolutionModel is set, may cause runtime errors
+ self._resolutionModel = None
+ if self._isMultiSpectrum:
+ for i in range(self.NumberOfSpectra):
+ if self.crystalFieldFunction.getAttributeValue('FWHMX%s' % i):
+ self.crystalFieldFunction.setAttributeValue('FWHMX%s' % i, [])
+ if self.crystalFieldFunction.getAttributeValue('FWHMY%s' % i):
+ self.crystalFieldFunction.setAttributeValue('FWHMY%s' % i, [])
+ else:
+ if self.crystalFieldFunction.getAttributeValue('FWHMX'):
+ self.crystalFieldFunction.setAttributeValue('FWHMX', [])
+ if self.crystalFieldFunction.getAttributeValue('FWHMY'):
+ self.crystalFieldFunction.setAttributeValue('FWHMY', [])
@property
def FWHMVariation(self):
@@ -434,13 +500,11 @@ class CrystalField(object):
@FWHMVariation.setter
def FWHMVariation(self, value):
self.crystalFieldFunction.setAttributeValue('FWHMVariation', float(value))
- self._dirty_spectra = True
def __getitem__(self, item):
return self.crystalFieldFunction.getParameterValue(item)
def __setitem__(self, key, value):
- self._dirty_spectra = True
self.crystalFieldFunction.setParameter(key, value)
@property
@@ -450,7 +514,7 @@ class CrystalField(object):
@ResolutionModel.setter
def ResolutionModel(self, value):
from .function import ResolutionModel
- if isinstance(value, ResolutionModel):
+ if hasattr(value, 'model'):
self._resolutionModel = value
else:
self._resolutionModel = ResolutionModel(value)
@@ -466,6 +530,11 @@ class CrystalField(object):
model = self._resolutionModel.model
self.crystalFieldFunction.setAttributeValue('FWHMX', model[0])
self.crystalFieldFunction.setAttributeValue('FWHMY', model[1])
+ # If FWHM is set, it overrides resolution model, so unset it
+ if self._isMultiSpectrum and any(self.crystalFieldFunction.getAttributeValue('FWHMs')):
+ self.crystalFieldFunction.setAttributeValue('FWHMs', [0.] * self.NumberOfSpectra)
+ elif not self._isMultiSpectrum and self.crystalFieldFunction.getAttributeValue('FWHM'):
+ self.crystalFieldFunction.setAttributeValue('FWHM', 0.)
@property
def FixAllPeaks(self):
@@ -511,11 +580,10 @@ class CrystalField(object):
value: an instance of function.Background class or a list of instances
in a multi-spectrum case
"""
- from .function import Background
from mantid.simpleapi import FunctionFactory
if self._background is not None:
raise ValueError('Background has been set already')
- if not isinstance(value, Background):
+ if not hasattr(value, 'toString'):
raise TypeError('Expected a Background object, found %s' % str(value))
if not self._isMultiSpectrum:
fun_str = value.toString() + ';' + str(self.function)
@@ -548,9 +616,8 @@ class CrystalField(object):
@PhysicalProperty.setter
def PhysicalProperty(self, value):
- from .function import PhysicalProperties
vlist = value if islistlike(value) else [value]
- if all([isinstance(pp, PhysicalProperties) for pp in vlist]):
+ if all([hasattr(pp, 'TypeID') for pp in vlist]):
nOldPP = len(self._physprop) if islistlike(self._physprop) else (0 if self._physprop is None else 1)
self._physprop = value
else:
@@ -576,7 +643,14 @@ class CrystalField(object):
self.function.setAttributeValue('PhysicalProperties', [0]*len(tt)+ppids)
for attribs in [pp.getAttributes(i+len(tt)) for i, pp in enumerate(vlist)]:
for item in attribs.items():
- self.function.setAttributeValue(item[0], item[1])
+ if 'Lambda' in item[0] or 'Chi0' in item[0]:
+ self.function.setParameter(item[0], item[1])
+ if item[1] == 0.:
+ self.function.tie(item[0], '0.')
+ else:
+ self.function.removeTie(item[0])
+ else:
+ self.function.setAttributeValue(item[0], item[1])
@property
def isPhysicalPropertyOnly(self):
@@ -642,10 +716,6 @@ class CrystalField(object):
@param ws_index: An index of a spectrum from workspace to use.
@return: A tuple of (x, y) arrays
"""
- if self._dirty_spectra:
- self._spectra = {}
- self._dirty_spectra = False
-
wksp = workspace
# Allow to call getSpectrum with a workspace as the first argument.
if not isinstance(i, int):
@@ -668,16 +738,12 @@ class CrystalField(object):
return self._calcSpectrum(i, wksp, ws_index)
if xArray is None:
- if i in self._spectra:
- return self._spectra[i]
- else:
- x_min, x_max = self.calc_xmin_xmax(i)
- xArray = np.linspace(x_min, x_max, self.default_spectrum_size)
+ x_min, x_max = self.calc_xmin_xmax(i)
+ xArray = np.linspace(x_min, x_max, self.default_spectrum_size)
yArray = np.zeros_like(xArray)
wksp = makeWorkspace(xArray, yArray)
- self._spectra[i] = self._calcSpectrum(i, wksp, 0)
- return self._spectra[i]
+ return self._calcSpectrum(i, wksp, 0)
def getHeatCapacity(self, workspace=None, ws_index=0):
"""
@@ -822,6 +888,22 @@ class CrystalField(object):
return self._getPhysProp(PhysicalProperties(pptype, *args, **kwargs), workspace, ws_index)
+ def getDipoleMatrix(self):
+ """Returns the dipole transition matrix as a numpy array"""
+ from scipy.constants import physical_constants
+ from CrystalField.energies import energies
+ self._calcEigensystem()
+ _, _, hx = energies(self._nre, BextX=1.0)
+ _, _, hy = energies(self._nre, BextY=1.0)
+ _, _, hz = energies(self._nre, BextZ=1.0)
+ ix = np.dot(np.conj(np.transpose(self._eigenvectors)), np.dot(hx, self._eigenvectors))
+ iy = np.dot(np.conj(np.transpose(self._eigenvectors)), np.dot(hy, self._eigenvectors))
+ iz = np.dot(np.conj(np.transpose(self._eigenvectors)), np.dot(hz, self._eigenvectors))
+ gj = 2. if (self._nre < 1) else self.lande_g[self._nre - 1]
+ gJuB = gj * physical_constants['Bohr magneton in eV/T'][0] * 1000.
+ trans = np.multiply(ix, np.conj(ix)) + np.multiply(iy, np.conj(iy)) + np.multiply(iz, np.conj(iz))
+ return trans / (gJuB ** 2)
+
def plot(self, i=0, workspace=None, ws_index=0, name=None):
"""Plot a spectrum. Parameters are the same as in getSpectrum(...)"""
from mantidplot import plotSpectrum
@@ -898,12 +980,32 @@ class CrystalField(object):
symmetries = [self.Symmetry, other.Symmetry]
params = {}
temperatures = [self._getTemperature(x) for x in range(self.NumberOfSpectra)]
- fwhms = [self._getFWHM(x) for x in range(self.NumberOfSpectra)]
for bparam in CrystalField.field_parameter_names:
params['ion0.' + bparam] = self[bparam]
params['ion1.' + bparam] = other[bparam]
- return CrystalFieldMultiSite(Ions=ions, Symmetries=symmetries, Temperatures=temperatures,
- FWHMs = fwhms, parameters=params, abundances=[1.0, 1.0])
+ ties = {}
+ fixes = []
+ for prefix, obj in iteritems({'ion0.':self, 'ion1.':other}):
+ tiestr = obj.function.getTies()
+ if tiestr:
+ for tiepair in [tie.split('=') for tie in tiestr.split(',')]:
+ ties[prefix + tiepair[0]] = tiepair[1]
+ for par_id in [id for id in range(obj.function.nParams()) if obj.function.isFixed(id)]:
+ parName = obj.function.getParamName(par_id)
+ if obj.background is not None:
+ parName = parName.split('.')[-1]
+ if parName not in self.field_parameter_names:
+ continue
+ fixes.append(prefix + parName)
+ if self._resolutionModel is None:
+ fwhms = [self._getFWHM(x) for x in range(self.NumberOfSpectra)]
+ return CrystalFieldMultiSite(Ions=ions, Symmetries=symmetries, Temperatures=temperatures,
+ FWHM = fwhms, parameters=params, abundances=[1.0, 1.0],
+ ties = ties, fixedParameters = fixes)
+ else:
+ return CrystalFieldMultiSite(Ions=ions, Symmetries=symmetries, Temperatures=temperatures,
+ ResolutionModel=self._resolutionModel, parameters=params,
+ abundances=[1.0, 1.0], ties = ties, fixedParameters = fixes)
def __mul__(self, factor):
ffactor = float(factor)
@@ -970,13 +1072,11 @@ class CrystalField(object):
return params
def _getFieldTies(self):
- import re
- ties = re.match('ties=\((.*?)\)', str(self.crystalFieldFunction))
+ ties = re.search('ties=\((.*?)\)', str(self.crystalFieldFunction))
return ties.group(1) if ties else ''
def _getFieldConstraints(self):
- import re
- constraints = re.match('constraints=\((.*?)\)', str(self.crystalFieldFunction))
+ constraints = re.search('constraints=\((.*?)\)', str(self.crystalFieldFunction))
return constraints.group(1) if constraints else ''
def _getPhysProp(self, ppobj, workspace, ws_index):
@@ -1012,9 +1112,8 @@ class CrystalField(object):
"""
if self._dirty_eigensystem:
import CrystalField.energies as energies
- nre = self.ion_nre_map[self.Ion]
self._eigenvalues, self._eigenvectors, self._hamiltonian = \
- energies.energies(nre, **self._getFieldParameters())
+ energies.energies(self._nre, **self._getFieldParameters())
self._dirty_eigensystem = False
def _calcPeaksList(self, i):
@@ -1066,6 +1165,7 @@ class CrystalFieldSite(object):
def __add__(self, other):
from CrystalField.CrystalFieldMultiSite import CrystalFieldMultiSite
+ from CrystalField import CrystalField
if isinstance(other, CrystalField):
abundances = [self.abundance, 1.0]
elif isinstance(other, CrystalFieldSite):
@@ -1079,233 +1179,18 @@ class CrystalFieldSite(object):
ions = [self.crystalField.Ion, other.Ion]
symmetries = [self.crystalField.Symmetry, other.Symmetry]
temperatures = [self.crystalField._getTemperature(x) for x in range(self.crystalField.NumberOfSpectra)]
- FWHMs = [self.crystalField._getFWHM(x) for x in range(self.crystalField.NumberOfSpectra)]
params = {}
for bparam in CrystalField.field_parameter_names:
params['ion0.' + bparam] = self.crystalField[bparam]
params['ion1.' + bparam] = other[bparam]
- return CrystalFieldMultiSite(Ions=ions, Symmetries=symmetries, Temperatures=temperatures, FWHMs = FWHMs,
- abundances=abundances, parameters=params)
-
-
-class CrystalFieldMulti(object):
- """CrystalFieldMulti represents crystal field of multiple ions."""
-
- def __init__(self, *args):
- self.sites = []
- self.abundances = []
- for arg in args:
- if isinstance(arg, CrystalField):
- self.sites.append(arg)
- self.abundances.append(1.0)
- elif isinstance(arg, CrystalFieldSite):
- self.sites.append(arg.crystalField)
- self.abundances.append(arg.abundance)
- else:
- raise RuntimeError('Cannot include an object of type %s into a CrystalFieldMulti' % type(arg))
- self._ties = {}
-
- def makeSpectrumFunction(self):
- fun = ';'.join([a.makeSpectrumFunction() for a in self.sites])
- fun += self._makeIntensityScalingTies()
- ties = self.getTies()
- if len(ties) > 0:
- fun += ';ties=(%s)' % ties
- return 'composite=CompositeFunction,NumDeriv=1;' + fun
-
- def makePhysicalPropertiesFunction(self):
- # Handles relative intensities. Scaling factors here a fixed attributes not
- # variable parameters and we require the sum to be unity.
- factors = np.array(self.abundances)
- sum_factors = np.sum(factors)
- factstr = [',ScaleFactor=%s' % (str(factors[i] / sum_factors)) for i in range(len(self.sites))]
- fun = ';'.join([a.makePhysicalPropertiesFunction()+factstr[i] for a,i in enumerate(self.sites)])
- ties = self.getTies()
- if len(ties) > 0:
- fun += ';ties=(%s)' % ties
- return fun
-
- def makeMultiSpectrumFunction(self):
- fun = ';'.join([a.makeMultiSpectrumFunction() for a in self.sites])
- fun += self._makeIntensityScalingTiesMulti()
- ties = self.getTies()
- if len(ties) > 0:
- fun += ';ties=(%s)' % ties
- return 'composite=CompositeFunction,NumDeriv=1;' + fun
-
- def ties(self, **kwargs):
- """Set ties on the parameters."""
- for tie in kwargs:
- self._ties[tie] = kwargs[tie]
-
- def getTies(self):
- ties = ['%s=%s' % item for item in self._ties.items()]
- return ','.join(ties)
-
- def getSpectrum(self, i=0, workspace=None, ws_index=0):
- tt = []
- for site in self.sites:
- tt = tt + (list(site.Temperature) if islistlike(site.Temperature) else [site.Temperature])
- if any([val < 0 for val in tt]):
- raise RuntimeError('You must first define a temperature for the spectrum')
- largest_abundance= max(self.abundances)
- if workspace is not None:
- xArray, yArray = self.sites[0].getSpectrum(i, workspace, ws_index)
- yArray *= self.abundances[0] / largest_abundance
- ia = 1
- for arg in self.sites[1:]:
- _, yyArray = arg.getSpectrum(i, workspace, ws_index)
- yArray += yyArray * self.abundances[ia] / largest_abundance
- ia += 1
- return xArray, yArray
- x_min = 0.0
- x_max = 0.0
- for arg in self.sites:
- xmin, xmax = arg.calc_xmin_xmax(i)
- if xmin < x_min:
- x_min = xmin
- if xmax > x_max:
- x_max = xmax
- xArray = np.linspace(x_min, x_max, CrystalField.default_spectrum_size)
- _, yArray = self.sites[0].getSpectrum(i, xArray, ws_index)
- yArray *= self.abundances[0] / largest_abundance
- ia = 1
- for arg in self.sites[1:]:
- _, yyArray = arg.getSpectrum(i, xArray, ws_index)
- yArray += yyArray * self.abundances[ia] / largest_abundance
- ia += 1
- return xArray, yArray
-
- def update(self, func):
- nFunc = func.nFunctions()
- assert nFunc == len(self.sites)
- for i in range(nFunc):
- self.sites[i].update(func[i])
-
- def update_multi(self, func):
- nFunc = func.nFunctions()
- assert nFunc == len(self.sites)
- for i in range(nFunc):
- self.sites[i].update_multi(func[i])
-
- def _makeIntensityScalingTies(self):
- """
- Make a tie string that ties IntensityScaling's of the sites according to their abundances.
- """
- n_sites = len(self.sites)
- if n_sites < 2:
- return ''
- factors = np.array(self.abundances)
- i_largest = np.argmax(factors)
- largest_factor = factors[i_largest]
- tie_template = 'f%s.IntensityScaling=%s*' + 'f%s.IntensityScaling' % i_largest
- ties = []
- for i in range(n_sites):
- if i != i_largest:
- ties.append(tie_template % (i, factors[i] / largest_factor))
- s = ';ties=(%s)' % ','.join(ties)
- return s
-
- def _makeIntensityScalingTiesMulti(self):
- """
- Make a tie string that ties IntensityScaling's of the sites according to their abundances.
- """
- n_sites = len(self.sites)
- if n_sites < 2:
- return ''
- factors = np.array(self.abundances)
- i_largest = np.argmax(factors)
- largest_factor = factors[i_largest]
- tie_template = 'f{1}.IntensityScaling{0}={2}*f%s.IntensityScaling{0}' % i_largest
- ties = []
- n_spectra = self.sites[0].NumberOfSpectra
- for spec in range(n_spectra):
- for i in range(n_sites):
- if i != i_largest:
- ties.append(tie_template.format(spec, i, factors[i] / largest_factor))
- s = ';ties=(%s)' % ','.join(ties)
- return s
-
- @property
- def isPhysicalPropertyOnly(self):
- return all([a.isPhysicalPropertyOnly for a in self.sites])
-
- @property
- def PhysicalProperty(self):
- return [a.PhysicalProperty for a in self.sites]
-
- @PhysicalProperty.setter
- def PhysicalProperty(self, value):
- for a in self.sites:
- a.PhysicalProperty = value
-
- @property
- def NumberOfSpectra(self):
- """ Returns the number of expected workspaces """
- num_spec = []
- for site in self.sites:
- num_spec.append(site.NumberOfSpectra)
- if len(set(num_spec)) > 1:
- raise ValueError('Number of spectra for each site not consistent with each other')
- return num_spec[0]
-
- @property
- def Temperature(self):
- tt = []
- for site in self.sites:
- tt.append([val for val in (site.Temperature if islistlike(site.Temperature) else [site.Temperature])])
- if len(set([tuple(val) for val in tt])) > 1:
- raise ValueError('Temperatures of spectra for each site not consistent with each other')
- return tt[0]
-
- @Temperature.setter
- def Temperature(self, value):
- for site in self.sites:
- site.Temperature = value
-
- def __add__(self, other):
- if isinstance(other, CrystalFieldMulti):
- cfm = CrystalFieldMulti()
- cfm.sites += self.sites + other.sites
- cfm.abundances += self.abundances + other.abundances
- return cfm
- elif isinstance(other, CrystalField):
- cfm = CrystalFieldMulti()
- cfm.sites += self.sites + [other]
- cfm.abundances += self.abundances + [1]
- return cfm
- elif isinstance(other, CrystalFieldSite):
- cfm = CrystalFieldMulti()
- cfm.sites += self.sites + [other.crystalField]
- cfm.abundances += self.abundances + [other.abundance]
- return cfm
+ if self.crystalField.ResolutionModel is None:
+ FWHM = [self.crystalField._getFWHM(x) for x in range(self.crystalField.NumberOfSpectra)]
+ return CrystalFieldMultiSite(Ions=ions, Symmetries=symmetries, Temperatures=temperatures, FWHM = FWHM,
+ abundances=abundances, parameters=params)
else:
- raise TypeError('Cannot add %s to CrystalFieldMulti' % other.__class__.__name__)
-
- def __radd__(self, other):
- if isinstance(other, CrystalFieldMulti):
- cfm = CrystalFieldMulti()
- cfm.sites += other.sites + self.sites
- cfm.abundances += other.abundances + self.abundances
- return cfm
- elif isinstance(other, CrystalField):
- cfm = CrystalFieldMulti()
- cfm.sites += [other] + self.sites
- cfm.abundances += [1] + self.abundances
- return cfm
- elif isinstance(other, CrystalFieldSite):
- cfm = CrystalFieldMulti()
- cfm.sites += [other.crystalField] + self.sites
- cfm.abundances += [other.abundance] + self.abundances
- return cfm
- else:
- raise TypeError('Cannot add %s to CrystalFieldMulti' % other.__class__.__name__)
-
- def __len__(self):
- return len(self.sites)
-
- def __getitem__(self, item):
- return self.sites[item]
+ return CrystalFieldMultiSite(Ions=ions, Symmetries=symmetries, Temperatures=temperatures,
+ ResolutionModel=self.crystalField.ResolutionModel,
+ abundances=abundances, parameters=params)
#pylint: disable=too-few-public-methods
@@ -1314,12 +1199,15 @@ class CrystalFieldFit(object):
Object that controls fitting.
"""
- def __init__(self, Model=None, Temperature=None, FWHM=None, InputWorkspace=None, **kwargs):
+ def __init__(self, Model=None, Temperature=None, FWHM=None, InputWorkspace=None,
+ ResolutionModel=None, **kwargs):
self.model = Model
if Temperature is not None:
self.model.Temperature = Temperature
if FWHM is not None:
self.model.FWHM = FWHM
+ if ResolutionModel is not None:
+ self.model.ResolutionModel = ResolutionModel
self._input_workspace = InputWorkspace
self._output_workspace_base_name = 'fit'
self._fit_properties = kwargs
@@ -1347,11 +1235,22 @@ class CrystalFieldFit(object):
def estimate_parameters(self, EnergySplitting, Parameters, **kwargs):
from CrystalField.normalisation import split2range
+ from CrystalField.CrystalFieldMultiSite import CrystalFieldMultiSite
from mantid.api import mtd
self.check_consistency()
- ranges = split2range(Ion=self.model.Ion, EnergySplitting=EnergySplitting,
- Parameters=Parameters)
- constraints = [('%s<%s<%s' % (-bound, parName, bound)) for parName, bound in ranges.items()]
+ if isinstance(self.model, CrystalFieldMultiSite):
+ constraints = []
+ for ni in range(len(self.model.Ions)):
+ pars = Parameters[ni] if islistlike(Parameters[ni]) else Parameters
+ ion = self.model.Ions[ni]
+ ranges = split2range(Ion=ion, EnergySplitting=EnergySplitting, Parameters=pars)
+ constraints += [('%s<ion%d.%s<%s' % (-bound, ni, parName, bound))
+ for parName, bound in ranges.items()]
+ else:
+ ranges = split2range(Ion=self.model.Ion, EnergySplitting=EnergySplitting,
+ Parameters=Parameters)
+ constraints = [('%s<%s<%s' % (-bound, parName, bound))
+ for parName, bound in ranges.items()]
self.model.constraints(*constraints)
if 'Type' not in kwargs or kwargs['Type'] == 'Monte Carlo':
if 'OutputWorkspace' in kwargs and kwargs['OutputWorkspace'].strip() != '':
@@ -1378,12 +1277,14 @@ class CrystalFieldFit(object):
ne = self.get_number_estimates()
if ne == 0:
raise RuntimeError('There are no estimated parameters.')
- if index >= ne:
+ if index > ne:
raise RuntimeError('There are only %s sets of estimated parameters, requested set #%s' % (ne, index))
+ from CrystalField.CrystalFieldMultiSite import CrystalFieldMultiSite
for row in range(self._estimated_parameters.rowCount()):
name = self._estimated_parameters.cell(row, 0)
value = self._estimated_parameters.cell(row, index)
- self.model[name] = value
+ model_pname = name if isinstance(self.model, CrystalFieldMultiSite) else name.split('.')[-1]
+ self.model[model_pname] = value
if self._function is not None:
self._function.setParameter(name, value)
diff --git a/scripts/Inelastic/CrystalField/function.py b/scripts/Inelastic/CrystalField/function.py
index 5d6f24e20778cd609b540287efacb76ca709709d..c92547e806ad2d4cbd3b5a842cb8fe38368c2121 100644
--- a/scripts/Inelastic/CrystalField/function.py
+++ b/scripts/Inelastic/CrystalField/function.py
@@ -433,13 +433,18 @@ class ResolutionModel:
x = self._mergeArrays(x, xx)
y = self._mergeArrays(y, yy)
n = len(x)
- return x, y
+ return list(x), list(y)
class PhysicalProperties(object):
"""
Contains information about measurement conditions of physical properties
"""
+ HEATCAPACITY = 1
+ SUSCEPTIBILITY = 2
+ MAGNETISATION = 3
+ MAGNETICMOMENT = 4
+
def _str2id(self, typeid):
mappings = [['cp', 'cv', 'heatcap'], ['chi', 'susc'], ['mag', 'm(h)'], ['mom', 'm(t)']]
for id in range(4):
@@ -461,19 +466,21 @@ class PhysicalProperties(object):
:param temperature: the temperature in Kelvin of measurements of M(H)
:param inverse: a boolean indicating whether susceptibility is chi or 1/chi or M(T) or 1/M(T)
:param unit: the unit the data was measured in. Either: 'bohr', 'SI' or 'cgs'.
+ :param lambda: (susceptibility only) the value of the exchange constant in inverse susc units
+ :param chi0: (susceptibility only) the value of the residual (background) susceptibility
typeid is required in all cases, and all other parameters may be specified as keyword arguments.
otherwise the syntax is:
PhysicalProperties('Cp') # No further parameters required for heat capacity
- PhysicalProperties('chi', hdir, inverse, unit)
+ PhysicalProperties('chi', hdir, inverse, unit, lambda, chi0)
PhysicalProperties('chi', unit)
PhysicalProperties('mag', hdir, temp, unit)
PhysicalProperties('mag', unit)
PhysicalProperties('M(T)', hmag, hdir, inverse, unit)
PhysicalProperties('M(T)', unit)
- Defaults are: hdir=[0, 0, 1]; hmag=1; temp=1; inverse=False; unit='cgs'.
+ Defaults are: hdir=[0, 0, 1]; hmag=1; temp=1; inverse=False; unit='cgs'; lambda=chi0=0.
"""
self._physpropUnit = 'cgs'
self._suscInverseFlag = False
@@ -481,6 +488,7 @@ class PhysicalProperties(object):
self._hmag = 1.
self._physpropTemperature = 1.
self._lambda = 0. # Exchange parameter (for susceptibility only)
+ self._chi0 = 0. # Residual/background susceptibility (for susceptibility only)
self._typeid = self._str2id(typeid) if isinstance(typeid, string_types) else int(typeid)
try:
initialiser = getattr(self, 'init' + str(self._typeid))
@@ -532,11 +540,11 @@ class PhysicalProperties(object):
@property
def Inverse(self):
- return self._suscInverseFlag if (self._typeid == 2 or self._typeid == 4) else None
+ return self._suscInverseFlag if (self._typeid == self.SUSCEPTIBILITY or self._typeid == self.MAGNETICMOMENT) else None
@Inverse.setter
def Inverse(self, value):
- if (self._typeid == 2 or self._typeid == 4):
+ if (self._typeid == self.SUSCEPTIBILITY or self._typeid == self.MAGNETICMOMENT):
if isinstance(value, string_types):
self._suscInverseFlag = value.lower() in ['true', 't', '1', 'yes', 'y']
else:
@@ -546,40 +554,49 @@ class PhysicalProperties(object):
@property
def Hdir(self):
- return self._hdir if (self._typeid > 1) else None
+ return self._hdir if (self._typeid != self.HEATCAPACITY) else None
@Hdir.setter
def Hdir(self, value):
- if (self._typeid > 1):
+ if (self._typeid != self.HEATCAPACITY):
self._hdir = self._checkhdir(value)
@property
def Hmag(self):
- return self._hmag if (self._typeid == 4) else None
+ return self._hmag if (self._typeid == self.MAGNETICMOMENT) else None
@Hmag.setter
def Hmag(self, value):
- if (self._typeid == 4):
+ if (self._typeid == self.MAGNETICMOMENT):
self._hmag = float(value)
@property
def Temperature(self):
- return self._physpropTemperature if (self._typeid == 3) else None
+ return self._physpropTemperature if (self._typeid == self.MAGNETISATION) else None
@Temperature.setter
def Temperature(self, value):
- if (self._typeid == 3):
+ if (self._typeid == self.MAGNETISATION):
self._physpropTemperature = float(value)
@property
def Lambda(self):
- return self._lambda if (self._typeid == 2) else None
+ return self._lambda if (self._typeid == self.SUSCEPTIBILITY) else None
@Lambda.setter
def Lambda(self, value):
- if (self._typeid == 2):
+ if (self._typeid == self.SUSCEPTIBILITY):
self._lambda = float(value)
+ @property
+ def Chi0(self):
+ return self._chi0 if (self._typeid == self.SUSCEPTIBILITY) else None
+
+ @Chi0.setter
+ def Chi0(self, value):
+ if (self._typeid == self.SUSCEPTIBILITY):
+ self._chi0 = float(value)
+
def init1(self, *args, **kwargs):
""" Initialises environment for heat capacity data """
if len(args) > 0:
@@ -602,7 +619,7 @@ class PhysicalProperties(object):
def init2(self, *args, **kwargs):
""" Initialises environment for susceptibility data """
- mapping = ['Hdir', 'Inverse', 'Unit', 'Lambda']
+ mapping = ['Hdir', 'Inverse', 'Unit', 'Lambda', 'Chi0']
self._parseargs(mapping, *args, **kwargs)
def init3(self, *args, **kwargs):
@@ -620,35 +637,38 @@ class PhysicalProperties(object):
types = ['CrystalFieldHeatCapacity', 'CrystalFieldSusceptibility',
'CrystalFieldMagnetisation', 'CrystalFieldMoment']
out = 'name=%s' % (types[self._typeid - 1])
- if self._typeid > 1:
+ if self._typeid != self.HEATCAPACITY:
out += ',Unit=%s' % (self._physpropUnit)
if 'powder' in self._hdir:
out += ',powder=1'
else:
out += ',Hdir=(%s)' % (','.join([str(hh) for hh in self._hdir]))
- if self._typeid == 3: # magnetisation M(H)
+ if self._typeid == self.MAGNETISATION:
out += ',Temperature=%s' % (self._physpropTemperature)
else: # either susceptibility or M(T)
out += ',inverse=%s' % (1 if self._suscInverseFlag else 0)
- out += (',Hmag=%s' % (self._hmag)) if self._typeid==3 else ''
- if self._typeid == 2 and self._lambda != 0:
+ out += (',Hmag=%s' % (self._hmag)) if self._typeid == self.MAGNETISATION else ''
+ if self._typeid == self.SUSCEPTIBILITY and self._lambda != 0:
out += ',Lambda=%s' % (self._lambda)
+ if self._typeid == self.SUSCEPTIBILITY and self._chi0 != 0:
+ out += ',Chi0=%s' % (self._chi0)
return out
def getAttributes(self, dataset=None):
"""Returns a dictionary of PhysicalProperties attributes for use with IFunction"""
dataset = '' if dataset is None else str(dataset)
out = {}
- if self._typeid > 1:
+ if self._typeid != self.HEATCAPACITY:
out['Unit%s' % (dataset)] = self._physpropUnit
if 'powder' in self._hdir:
out['powder%s' % (dataset)] = 1
else:
out['Hdir%s' % (dataset)] = [float(hh) for hh in self._hdir] # needs to be list
- if self._typeid != 3: # either susceptibility or M(T)
+ if self._typeid != self.MAGNETISATION: # either susceptibility or M(T)
out['inverse%s' % (dataset)] = 1 if self._suscInverseFlag else 0
- if self._typeid==3:
+ if self._typeid == self.MAGNETICMOMENT:
out['Hmag%s' % (dataset)] = self._hmag
- if self._typeid == 2 and self._lambda != 0:
- out['Lambda%s=' % (dataset)] = self._lambda
+ if self._typeid == self.SUSCEPTIBILITY:
+ out['Lambda%s' % (dataset)] = self._lambda
+ out['Chi0%s' % (dataset)] = self._chi0
return out
diff --git a/scripts/Inelastic/CrystalField/normalisation.py b/scripts/Inelastic/CrystalField/normalisation.py
index cd2f9178697fbb906f32440c10419028d70c6ef1..f1df0597b84e9211e7e6ebd6792d23513d495548 100644
--- a/scripts/Inelastic/CrystalField/normalisation.py
+++ b/scripts/Inelastic/CrystalField/normalisation.py
@@ -1,6 +1,7 @@
from __future__ import (absolute_import, division, print_function)
import numpy as np
from CrystalField.energies import energies as CFEnergy
+from CrystalField.fitting import ionname2Nre
from six import iteritems
from six import string_types
@@ -9,13 +10,14 @@ def _get_normalisation(nre, bnames):
""" Helper function to calculate the normalisation factor.
Defined as: ||Blm|| = sum_{Jz,Jz'} |<Jz|Blm|Jz'>|^2 / (2J+1)
"""
- J = [0, 5./2, 4, 9./2, 4, 5./2, 0, 7./2, 6, 15./2, 8, 15./2, 6, 7./2]
+ Jvals = [0, 5./2, 4, 9./2, 4, 5./2, 0, 7./2, 6, 15./2, 8, 15./2, 6, 7./2]
+ J = (-nre / 2.) if (nre < 0) else Jvals[nre]
retval = {}
for bname in bnames:
bdict = {bname: 1}
ee, vv, ham = CFEnergy(nre, **bdict)
Omat = np.mat(ham)
- norm = np.trace(np.real(Omat * np.conj(Omat))) / (2*J[nre]+1)
+ norm = np.trace(np.real(Omat * np.conj(Omat))) / (2*J+1)
retval[bname] = np.sqrt(np.abs(norm)) * np.sign(norm)
return retval
@@ -25,12 +27,11 @@ def _parse_args(**kwargs):
"""
# Some definitions
Blms = ['B20', 'B21', 'B22', 'B40', 'B41', 'B42', 'B43', 'B44', 'B60', 'B61', 'B62', 'B63', 'B64', 'B65', 'B66']
- Ions = ['Ce', 'Pr', 'Nd', 'Pm', 'Sm', 'Eu', 'Gd', 'Tb', 'Dy', 'Ho', 'Er', 'Tm', 'Yb']
# Some Error checking
if 'Ion' not in kwargs.keys() and 'IonNum' not in kwargs.keys():
raise NameError('You must specify the ion using either the ''Ion'', ''IonNum'' keywords')
if 'Ion' in kwargs.keys():
- nre = [id for id, val in enumerate(Ions) if val == kwargs['Ion']][0] + 1
+ nre = ionname2Nre(kwargs['Ion'])
else:
nre = kwargs['IonNum']
# Now parses the list of input crystal field parameters
@@ -184,9 +185,8 @@ def split2range(*args, **kwargs):
# Error checking
if 'Ion' not in argin.keys() and 'IonNum' not in argin.keys():
raise NameError('You must specify the ion using either the ''Ion'', ''IonNum'' keywords')
- Ions = ['Ce', 'Pr', 'Nd', 'Pm', 'Sm', 'Eu', 'Gd', 'Tb', 'Dy', 'Ho', 'Er', 'Tm', 'Yb']
if 'Ion' in argin.keys():
- nre = [ id for id,val in enumerate(Ions) if val==kwargs['Ion'] ][0] + 1
+ nre = ionname2Nre(kwargs['Ion'])
else:
nre = argin['IonNum']
if 'EnergySplitting' not in argin.keys():
diff --git a/scripts/test/CrystalFieldMultiSiteTest.py b/scripts/test/CrystalFieldMultiSiteTest.py
index 96b7af7b73e404343ff8a3d1e8337ca482e174c7..dc8c8404c535a7f156ede63623c8f2113d0765c3 100644
--- a/scripts/test/CrystalFieldMultiSiteTest.py
+++ b/scripts/test/CrystalFieldMultiSiteTest.py
@@ -9,19 +9,19 @@ c_mbsr = 79.5774715459 # Conversion from barn to mb/sr
class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_init_single_ion(self):
- cfms = CrystalFieldMultiSite(Ions='Pm', Symmetries='D2', Temperatures=[20, 52], FWHMs=[1.0, 1.0])
+ cfms = CrystalFieldMultiSite(Ions='Pm', Symmetries='D2', Temperatures=[20, 52], FWHM=[1.0, 1.0])
self.assertEqual(cfms.Ions, ['Pm'])
self.assertEqual(cfms.Symmetries, ['D2'])
self.assertEqual(cfms.Temperatures, [20, 52])
- self.assertEqual(cfms.FWHMs, [1.0, 1.0])
+ self.assertEqual(cfms.FWHM, [1.0, 1.0])
def test_init_multi_ions(self):
cfms = CrystalFieldMultiSite(Ions=('Pm', 'Eu'), Symmetries=('D2', 'C3v'), Temperatures=[20, 52],
- FWHMs=[1.0, 1.0])
+ FWHM=[1.0, 1.0])
self.assertEqual(cfms.Ions, ['Pm', 'Eu'])
self.assertEqual(cfms.Symmetries, ['D2', 'C3v'])
self.assertEqual(cfms.Temperatures, [20, 52])
- self.assertEqual(cfms.FWHMs, [1.0, 1.0])
+ self.assertEqual(cfms.FWHM, [1.0, 1.0])
def test_toleranceEnergy_property(self):
cfms = CrystalFieldMultiSite(('Pm', 'Eu'), ('D2', 'C3v'), ToleranceEnergy=1.0)
@@ -48,28 +48,28 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
self.assertEqual(cfms.FWHMVariation, 0.3)
def test_init_parameters(self):
- cfms = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2', Temperatures=[15], FWHMs=[1.0],
+ cfms = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2', Temperatures=[15], FWHM=[1.0],
parameters={'BmolX': 1.0, 'B40': -0.02})
self.assertEqual(cfms.function.getParameterValue('BmolX'), 1.0)
self.assertEqual(cfms.function.getParameterValue('B40'), -0.02)
def test_init_parameters_multi_ions(self):
cfms = CrystalFieldMultiSite(Ions=('Pm', 'Eu'), Symmetries=('D2', 'C3v'), Temperatures=[20, 52],
- FWHMs=[1.0, 1.0], parameters={'ion0.B40': -0.02, 'ion0.B42': -0.11,
- 'ion1.B42': -0.12})
+ FWHM=[1.0, 1.0], parameters={'ion0.B40': -0.02, 'ion0.B42': -0.11,
+ 'ion1.B42': -0.12})
self.assertEqual(cfms.function.getParameterValue('ion0.B42'), -0.11)
self.assertEqual(cfms.function.getParameterValue('ion0.B40'), -0.02)
self.assertEqual(cfms.function.getParameterValue('ion1.B42'), -0.12)
def test_peak_values(self):
- cfms = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2', Temperatures=[25], FWHMs=[1.5],
+ cfms = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2', Temperatures=[25], FWHM=[1.5],
BmolX=1.0, B40=-0.02)
self.assertEqual(int(cfms.function.getParameterValue('pk0.Amplitude')), 48)
self.assertEqual(cfms.function.getParameterValue('pk0.FWHM'), 1.5)
self.assertEqual(cfms.function.getParameterValue('pk0.PeakCentre'), 0)
def test_peak_values_gaussian(self):
- cfms = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2', Temperatures=[25], FWHMs=[1.0], PeakShape='Gaussian',
+ cfms = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2', Temperatures=[25], FWHM=[1.0], PeakShape='Gaussian',
BmolX=1.0, B40=-0.02)
self.assertEqual(int(cfms.function.getParameterValue('pk0.Height')), 45)
self.assertAlmostEqual(cfms.function.getParameterValue('pk0.Sigma'), 0.42, 2)
@@ -77,7 +77,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_peak_values_multi_ions(self):
cfms = CrystalFieldMultiSite(Ions=('Pr', 'Nd'), Symmetries=('D2', 'C3v'), Temperatures=[20],
- FWHMs=[1.5], parameters={'ion0.B60': -0.02, 'ion1.B62': -0.12})
+ FWHM=[1.5], parameters={'ion0.B60': -0.02, 'ion1.B62': -0.12})
self.assertEqual(int(cfms.function.getParameterValue('ion0.pk0.Amplitude')), 278)
self.assertEqual(cfms.function.getParameterValue('ion0.pk0.FWHM'), 1.5)
self.assertEqual(cfms.function.getParameterValue('ion0.pk1.PeakCentre'), 982.8)
@@ -87,7 +87,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_peak_values_multi_ions_and_spectra(self):
cfms = CrystalFieldMultiSite(Ions=('Pm', 'Ce'), Symmetries=('D2', 'C3v'), Temperatures=[20, 52],
- FWHMs=[1.0, 1.0], parameters={'ion0.B40': -0.02, 'ion1.B42': -0.12})
+ FWHM=[1.0, 1.0], parameters={'ion0.B40': -0.02, 'ion1.B42': -0.12})
self.assertEqual(int(cfms.function.getParameterValue('ion0.sp0.pk0.Amplitude')), 308)
self.assertEqual(cfms.function.getParameterValue('ion0.sp1.pk0.FWHM'), 1.0)
self.assertEqual(cfms.function.getParameterValue('ion0.sp0.pk1.PeakCentre'), 0)
@@ -97,7 +97,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_peak_values_multi_gaussian(self):
cfms = CrystalFieldMultiSite(Ions=('Pm', 'Dy'), Symmetries=('D2', 'C3v'), Temperatures=[20, 52],
- FWHMs=[1.0, 1.5], parameters={'ion0.B40': -0.02, 'ion1.B42': -0.12},
+ FWHM=[1.0, 1.5], parameters={'ion0.B40': -0.02, 'ion1.B42': -0.12},
PeakShape='Gaussian')
self.assertEqual(int(cfms.function.getParameterValue('ion0.sp0.pk0.Height')), 289)
self.assertAlmostEqual(cfms.function.getParameterValue('ion0.sp1.pk0.Sigma'), 0.64, 2)
@@ -107,7 +107,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
self.assertAlmostEqual(cfms.function.getParameterValue('ion1.sp0.pk1.PeakCentre'), 40.957515, 6)
def test_get_spectrum_from_list(self):
- cfms = CrystalFieldMultiSite(Ions=['Ce'], Symmetries=['C2v'], Temperatures=[4.0], FWHMs=[0.1],
+ cfms = CrystalFieldMultiSite(Ions=['Ce'], Symmetries=['C2v'], Temperatures=[4.0], FWHM=[0.1],
B20=0.035, B40=-0.012, B43=-0.027, B60=-0.00012, B63=0.0025, B66=0.0068,
ToleranceIntensity=0.001*c_mbsr)
r = [0.0, 1.45, 2.4, 3.0, 3.85]
@@ -120,7 +120,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_get_spectrum_ws(self):
from mantid.simpleapi import CreateWorkspace
cfms = CrystalFieldMultiSite(['Ce'], ['C2v'], B20=0.035, B40=-0.012, B43=-0.027, B60=-0.00012, B63=0.0025,
- B66=0.0068, Temperatures=[4.0], FWHMs=[0.1], ToleranceIntensity=0.001*c_mbsr)
+ B66=0.0068, Temperatures=[4.0], FWHM=[0.1], ToleranceIntensity=0.001*c_mbsr)
x = np.linspace(0.0, 2.0, 30)
y = np.zeros_like(x)
@@ -146,7 +146,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
self.assertAlmostEqual(y[1], 0.054428, 6)
def test_get_spectrum_from_list_multi_spectra(self):
- cfms = CrystalFieldMultiSite(Ions=['Ce'], Symmetries=['C2v'], Temperatures=[4.0, 50.0], FWHMs=[0.1, 0.2],
+ cfms = CrystalFieldMultiSite(Ions=['Ce'], Symmetries=['C2v'], Temperatures=[4.0, 50.0], FWHM=[0.1, 0.2],
B20=0.035, B40=-0.012, B43=-0.027, B60=-0.00012, B63=0.0025, B66=0.0068)
r = [0.0, 1.45, 2.4, 3.0, 3.85]
x, y = cfms.getSpectrum(0, r)
@@ -164,7 +164,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_get_spectrum_ws_multi_spectra(self):
from mantid.simpleapi import CreateWorkspace
cfms = CrystalFieldMultiSite(['Ce'], ['C2v'], B20=0.035, B40=-0.012, B43=-0.027, B60=-0.00012, B63=0.0025, B66=0.0068,
- Temperatures=[4.0, 50.0], FWHMs=[0.1, 0.2])
+ Temperatures=[4.0, 50.0], FWHM=[0.1, 0.2])
x = np.linspace(0.0, 2.0, 30)
y = np.zeros_like(x)
@@ -201,7 +201,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
'ion0.B44': -0.12544, 'ion1.B20': 0.37737, 'ion1.B22': 3.9770, 'ion1.B40': -0.031787,
'ion1.B42': -0.11611, 'ion1.B44': -0.12544}
cfms = CrystalFieldMultiSite(Ions=['Ce', 'Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44.0, 50.0],
- FWHMs=[1.1, 1.2], parameters=params)
+ FWHM=[1.1, 1.2], parameters=params)
r = [0.0, 1.45, 2.4, 3.0, 3.85]
x, y = cfms.getSpectrum(0, r)
y = y / c_mbsr
@@ -233,7 +233,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
params = {'ion0.B20': 0.37737, 'ion0.B22': 3.9770, 'ion0.B40': -0.031787, 'ion0.B42': -0.11611,
'ion0.B44': -0.12544, 'ion1.B20': 0.37737, 'ion1.B22': 3.9770, 'ion1.B40': -0.031787,
'ion1.B42': -0.11611, 'ion1.B44': -0.12544}
- cf = CrystalFieldMultiSite(Ions=['Ce', 'Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44.0], FWHMs=[1.1],
+ cf = CrystalFieldMultiSite(Ions=['Ce', 'Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44.0], FWHM=[1.1],
ToleranceIntensity=6.0, ToleranceEnergy=1.0, FixAllPeaks=True, parameters=params)
cf.fix('ion0.BmolX', 'ion0.BmolY', 'ion0.BmolZ', 'ion0.BextX', 'ion0.BextY', 'ion0.BextZ', 'ion0.B40',
@@ -267,7 +267,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
'ion0.B44': -0.12544, 'ion1.B20': 0.37737, 'ion1.B22': 3.9770, 'ion1.B40': -0.031787,
'ion1.B42': -0.11611, 'ion1.B44': -0.12544}
cf = CrystalFieldMultiSite(Ions=['Ce', 'Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44.0, 50.0],
- FWHMs=[1.0, 1.0], ToleranceIntensity=6.0, ToleranceEnergy=1.0, FixAllPeaks=True,
+ FWHM=[1.0, 1.0], ToleranceIntensity=6.0, ToleranceEnergy=1.0, FixAllPeaks=True,
parameters=params)
cf.fix('ion0.BmolX', 'ion0.BmolY', 'ion0.BmolZ', 'ion0.BextX', 'ion0.BextY', 'ion0.BextZ', 'ion0.B40',
@@ -285,7 +285,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
self.assertTrue(cf.chi2 < chi2)
def test_set_background(self):
- cf = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHMs=[1.0], Background='name=LinearBackground,A0=1')
+ cf = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHM=[1.0], Background='name=LinearBackground,A0=1')
self.assertEquals('"name=LinearBackground,A0=1"', cf['Background'])
self.assertEquals(cf.background.param['A0'], 1)
self.assertEquals(cf.background.background.param['A0'], 1)
@@ -294,7 +294,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_set_background_as_function(self):
from mantid.fitfunctions import LinearBackground
- cf = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHMs=[1.0],
+ cf = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHM=[1.0],
Background=LinearBackground(A0=1))
self.assertEquals('"name=LinearBackground,A0=1,A1=0"', cf['Background'])
self.assertEquals(cf.background.param['A0'], 1)
@@ -304,7 +304,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_set_background_with_peak(self):
from mantid.fitfunctions import Gaussian
- cf = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHMs=[1.0], Background='name=LinearBackground', BackgroundPeak=Gaussian(Height=1))
+ cf = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHM=[1.0], Background='name=LinearBackground', BackgroundPeak=Gaussian(Height=1))
self.assertEquals('"name=Gaussian,Height=1,PeakCentre=0,Sigma=0;name=LinearBackground"', cf['Background'])
self.assertEquals(cf.background.peak.param['Height'], 1)
self.assertEquals(cf.background.param['f0.Height'], 1)
@@ -316,7 +316,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_set_background_peak_only(self):
from mantid.fitfunctions import Gaussian
- cf = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHMs=[1.0], BackgroundPeak=Gaussian(Sigma=1))
+ cf = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHM=[1.0], BackgroundPeak=Gaussian(Sigma=1))
self.assertEquals('"name=Gaussian,Height=0,PeakCentre=0,Sigma=1"', cf['Background'])
self.assertEquals(cf.background.peak.param['Sigma'], 1)
self.assertEquals(cf.background.param['Sigma'], 1)
@@ -325,7 +325,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_set_background_composite(self):
from mantid.fitfunctions import Gaussian, LinearBackground
- cf = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHMs=[1.0],
+ cf = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHM=[1.0],
Background= Gaussian(PeakCentre=1) + LinearBackground())
self.assertEquals('"name=Gaussian,Height=0,PeakCentre=1,Sigma=0;name=LinearBackground,A0=0,A1=0"', cf['Background'])
cf.background.param['f1.A0'] = 1
@@ -334,7 +334,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
self.assertEquals(cf.background.param['f0.PeakCentre'], 0.5)
def test_set_background_composite_as_string(self):
- cf = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHMs=[1.0],
+ cf = CrystalFieldMultiSite(Ions='Ce', Symmetries='C2v', Temperatures=[20], FWHM=[1.0],
Background='name=Gaussian,Height=0,PeakCentre=1,Sigma=0;name=LinearBackground,A0=0,A1=0')
self.assertEquals('"name=Gaussian,Height=0,PeakCentre=1,Sigma=0;name=LinearBackground,A0=0,A1=0"', cf['Background'])
self.assertEquals(cf.background.param['f0.Sigma'], 0)
@@ -345,7 +345,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
from mantid.fitfunctions import Gaussian, LinearBackground
from mantid.simpleapi import FunctionFactory
- cf = CrystalFieldMultiSite(Ions=['Ce'], Symmetries=['C2v'], Temperatures=[50], FWHMs=[0.9],
+ cf = CrystalFieldMultiSite(Ions=['Ce'], Symmetries=['C2v'], Temperatures=[50], FWHM=[0.9],
B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.12544,
Background=LinearBackground(A0=1.0), BackgroundPeak=Gaussian(Height=10, Sigma=0.3))
@@ -377,7 +377,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_constraints_multi_spectrum(self):
from mantid.fitfunctions import FlatBackground, Gaussian
- cf = CrystalFieldMultiSite(Ions=['Ce'], Symmetries=['C2v'], Temperatures=[44, 50], FWHMs=[1.1, 0.9],
+ cf = CrystalFieldMultiSite(Ions=['Ce'], Symmetries=['C2v'], Temperatures=[44, 50], FWHM=[1.1, 0.9],
Background=FlatBackground(), BackgroundPeak=Gaussian(Height=10, Sigma=0.3),
B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.12544)
@@ -408,7 +408,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_constraints_multi_spectrum_and_ion(self):
from mantid.fitfunctions import FlatBackground, Gaussian
- cf = CrystalFieldMultiSite(Ions=['Ce','Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44, 50], FWHMs=[1.1, 0.9],
+ cf = CrystalFieldMultiSite(Ions=['Ce','Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44, 50], FWHM=[1.1, 0.9],
Background=FlatBackground(), BackgroundPeak=Gaussian(Height=10, Sigma=0.3),
parameters={'ion0.B20': 0.37737, 'ion0.B22': 3.9770, 'ion1.B40':-0.031787,
'ion1.B42':-0.11611, 'ion1.B44':-0.12544})
@@ -440,7 +440,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_add_CrystalField(self):
from CrystalField import CrystalField
- cfms = CrystalFieldMultiSite(Ions=['Pr'], Symmetries=['C2v'], Temperatures=[44, 50], FWHMs=[1.1, 0.9],
+ cfms = CrystalFieldMultiSite(Ions=['Pr'], Symmetries=['C2v'], Temperatures=[44, 50], FWHM=[1.1, 0.9],
B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.15)
params = {'B20': 0.37737, 'B22': 3.9770, 'B40': -0.031787, 'B42': -0.11611, 'B44': -0.12544}
@@ -463,7 +463,7 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
def test_add_CrystalFieldSite(self):
from CrystalField import CrystalField
- cfms = CrystalFieldMultiSite(Ions=['Pr'], Symmetries=['C2v'], Temperatures=[44, 50], FWHMs=[1.1, 0.9],
+ cfms = CrystalFieldMultiSite(Ions=['Pr'], Symmetries=['C2v'], Temperatures=[44, 50], FWHM=[1.1, 0.9],
B20=0.37737, B22=3.9770, B40=-0.031787, B42=-0.11611, B44=-0.15)
params = {'B20': 0.37737, 'B22': 3.9770, 'B40': -0.031787, 'B42': -0.11611, 'B44': -0.12544}
@@ -482,11 +482,11 @@ class CrystalFieldMultiSiteTests(unittest.TestCase):
self.assertTrue('ion0.IntensityScaling=0.125*ion1.IntensityScaling' in s)
def test_add_CrystalFieldMultiSite(self):
- cfms1 = CrystalFieldMultiSite(Ions=['Pm'], Symmetries=['D2'], Temperatures=[44, 50], FWHMs=[1.1, 0.9],
+ cfms1 = CrystalFieldMultiSite(Ions=['Pm'], Symmetries=['D2'], Temperatures=[44, 50], FWHM=[1.1, 0.9],
B20=0.37737, B40=-0.031787, B42=-0.11611, B44=-0.15)
cfms2 = CrystalFieldMultiSite(Ions=['Ce', 'Pr'], Symmetries=['C2v', 'C2v'], Temperatures=[44, 50],
- FWHMs=[1.1, 0.9], parameters={'ion0.B20': 0.34, 'ion0.B22': 3.9770,
- 'ion1.B40': -0.03})
+ FWHM=[1.1, 0.9], parameters={'ion0.B20': 0.34, 'ion0.B22': 3.9770,
+ 'ion1.B40': -0.03})
cfms3 = cfms1 + cfms2
self.assertEqual(len(cfms3.Ions), 3)
self.assertEqual(len(cfms3.Symmetries), 3)
diff --git a/scripts/test/CrystalFieldTest.py b/scripts/test/CrystalFieldTest.py
index 5fa0706c89f93135f46d73dbf903b3b4ec54c0bc..b704f56075dc1feff282ccc94a2b3bfb17abcf06 100644
--- a/scripts/test/CrystalFieldTest.py
+++ b/scripts/test/CrystalFieldTest.py
@@ -1105,14 +1105,14 @@ class CrystalFieldFitTest(unittest.TestCase):
def test_ResolutionModel_func_multi(self):
from CrystalField import ResolutionModel
def func0(x):
- return np.sin(x)
+ return np.sin(np.array(x))
class CalcWidth:
def __call__(self, x):
- return np.cos(x / 2)
+ return np.cos(np.array(x) / 2)
def model(self, x):
- return np.tan(x / 2)
+ return np.tan(np.array(x) / 2)
func1 = CalcWidth()
func2 = func1.model
rm = ResolutionModel([func0, func1, func2], -np.pi/2, np.pi/2, accuracy = 0.01)