Commit a0b240a9 authored by Roman Tolchenov's avatar Roman Tolchenov
Browse files

Re #15702. Multi-ion multi-spectrum fitting.

parent c39685e6
......@@ -160,6 +160,11 @@ public:
void setUpForFit() override;
/// Remove a constraint
void removeConstraint(const std::string &parName) override;
/// Get number of domains required by this function
size_t getNumberDomains() const override;
/// Split this function (if needed) into a list of independent functions.
std::vector<boost::shared_ptr<IFunction>>
createEquivalentFunctions() const override;
/* CompositeFunction own methods */
......
......@@ -746,5 +746,52 @@ CompositeFunction::getContainingFunction(const ParameterReference &ref) const {
return IFunction_sptr();
}
/// Get number of domains required by this function
size_t CompositeFunction::getNumberDomains() const {
auto n = nFunctions();
if (n == 0) {
return 1;
}
size_t nd = getFunction(0)->getNumberDomains();
for (size_t iFun = 1; iFun < n; ++iFun) {
if (getFunction(0)->getNumberDomains() != nd) {
throw std::runtime_error("CompositeFunction has members with "
"inconsistent domain numbers.");
}
}
return nd;
}
/// Split this function (if needed) into a list of independent functions.
/// The number of functions must be the number of domains this function is
/// working on (== getNumberDomains()). The result of evaluation of the
/// created functions on their domains must be the same as if this function
/// was evaluated on the composition of those domains.
std::vector<IFunction_sptr> CompositeFunction::createEquivalentFunctions() const {
auto nd = getNumberDomains();
if (nd == 1) {
return std::vector<IFunction_sptr>(
1, FunctionFactory::Instance().createInitialized(asString()));
}
auto nf = nFunctions();
std::vector<std::vector<IFunction_sptr>> equiv;
equiv.reserve(nf);
for(size_t i = 0; i < nf; ++i) {
equiv.push_back(getFunction(i)->createEquivalentFunctions());
}
std::vector<IFunction_sptr> funs;
funs.reserve(nd);
for(size_t i = 0; i < nd; ++i) {
auto comp = new CompositeFunction;
funs.push_back(IFunction_sptr(comp));
for(size_t j = 0; j < nf; ++j) {
comp->addFunction(equiv[j][i]);
}
}
return funs;
}
} // namespace API
} // namespace Mantid
......@@ -47,16 +47,20 @@ private:
/// Build a function for a single spectrum.
API::IFunction_sptr buildSpectrum(int nre, const DoubleFortranVector &en,
const ComplexFortranMatrix &wf,
double temperature, double fwhm) const;
double temperature, double fwhm, size_t i) const;
/// Update a function for a single spectrum.
bool updateSpectrum(API::IFunction &spectrum, int nre,
void updateSpectrum(API::IFunction &spectrum, int nre,
const DoubleFortranVector &en,
const ComplexFortranMatrix &wf, double temperature) const;
const ComplexFortranMatrix &wf, double temperature, size_t i) const;
/// Calculate excitations at given temperature
void calcExcitations(int nre, const DoubleFortranVector &en,
const ComplexFortranMatrix &wf, double temperature,
DoubleFortranVector &eExcitations,
DoubleFortranVector &iExcitations) const;
/// Cache number of fitted peaks
mutable std::vector<size_t> m_nPeaks;
/// Cache number of all peaks
mutable std::vector<size_t> m_maxNPeaks;
};
} // namespace Functions
......
#include "MantidCurveFitting/Functions/CrystalElectricField.h"
#include "MantidCurveFitting/Functions/CrystalFieldMultiSpectrum.h"
#include "MantidCurveFitting/Functions/CrystalElectricField.h"
#include "MantidCurveFitting/Functions/CrystalFieldPeaks.h"
#include "MantidAPI/FunctionFactory.h"
......@@ -110,8 +110,11 @@ void CrystalFieldMultiSpectrum::buildTargetFunction() const {
}
}
// Create the single-spectrum functions.
for (size_t i = 0; i < temperatures.size(); ++i) {
fun->addFunction(buildSpectrum(nre, en, wf, temperatures[i], fwhms[i]));
auto nSpec = temperatures.size();
m_nPeaks.resize(nSpec);
m_maxNPeaks.resize(nSpec);
for (size_t i = 0; i < nSpec; ++i) {
fun->addFunction(buildSpectrum(nre, en, wf, temperatures[i], fwhms[i], i));
fun->setDomainIndex(i, i);
}
}
......@@ -135,23 +138,27 @@ void CrystalFieldMultiSpectrum::calcExcitations(
/// Build a function for a single spectrum.
API::IFunction_sptr CrystalFieldMultiSpectrum::buildSpectrum(
int nre, const DoubleFortranVector &en, const ComplexFortranMatrix &wf,
double temperature, double fwhm) const {
double temperature, double fwhm, size_t i) const {
DoubleFortranVector eExcitations;
DoubleFortranVector iExcitations;
calcExcitations(nre, en, wf, temperature, eExcitations, iExcitations);
size_t nPeaks = eExcitations.size();
size_t maxNPeaks = nPeaks + nPeaks / 2 + 1;
auto nPeaks = eExcitations.size();
auto maxNPeaks = nPeaks + nPeaks / 2 + 1;
m_nPeaks[i] = nPeaks;
m_maxNPeaks[i] = maxNPeaks;
bool fixAll = getAttribute("FixAllPeakParameters").asBool();
auto peakShape = IFunction::getAttribute("PeakShape").asString();
auto bkgdShape = IFunction::getAttribute("Background").asUnquotedString();
if (!bkgdShape.empty() && bkgdShape.find("name=") != 0 && bkgdShape.front() != '(') {
if (!bkgdShape.empty() && bkgdShape.find("name=") != 0 &&
bkgdShape.front() != '(') {
bkgdShape = "name=" + bkgdShape;
}
auto spectrum = new CompositeFunction;
auto background = API::FunctionFactory::Instance().createInitialized(bkgdShape);
auto background =
API::FunctionFactory::Instance().createInitialized(bkgdShape);
spectrum->addFunction(background);
for (size_t i = 0; i < maxNPeaks; ++i) {
auto fun = API::FunctionFactory::Instance().createFunction(peakShape);
......@@ -194,45 +201,40 @@ void CrystalFieldMultiSpectrum::updateTargetFunction() const {
auto &fun = dynamic_cast<MultiDomainFunction &>(*m_target);
auto temperatures = getAttribute("Temperatures").asVector();
for (size_t i = 0; i < temperatures.size(); ++i) {
auto res =
updateSpectrum(*fun.getFunction(i), nre, en, wf, temperatures[i]);
if (!res) {
buildTargetFunction();
return;
}
updateSpectrum(*fun.getFunction(i), nre, en, wf, temperatures[i], i);
}
}
/// Update a function for a single spectrum.
/// If returns false the target must be rebuilt with buldTargetFunction()
bool CrystalFieldMultiSpectrum::updateSpectrum(API::IFunction &spectrum,
int nre,
const DoubleFortranVector &en,
const ComplexFortranMatrix &wf,
double temperature) const {
void CrystalFieldMultiSpectrum::updateSpectrum(
API::IFunction &spectrum, int nre, const DoubleFortranVector &en,
const ComplexFortranMatrix &wf, double temperature, size_t i) const {
DoubleFortranVector eExcitations;
DoubleFortranVector iExcitations;
calcExcitations(nre, en, wf, temperature, eExcitations, iExcitations);
size_t nPeaks = eExcitations.size();
size_t maxNPeaks = nPeaks + nPeaks / 2 +
1; // m_source->getAttribute("MaxPeakCount").asInt();
size_t nGoodPeaks = eExcitations.size();
size_t nOriginalPeaks = m_nPeaks[i];
size_t maxNPeaks = m_maxNPeaks[i];
auto &composite = dynamic_cast<CompositeFunction &>(spectrum);
if (maxNPeaks + 1 != composite.nFunctions()) {
return false;
throw std::logic_error(
"CrystalFieldMultiSpectrum target function size mismatch.");
}
for (size_t i = 0; i < maxNPeaks; ++i) {
auto fun = composite.getFunction(i + 1);
auto &peak = dynamic_cast<API::IPeakFunction &>(*fun);
if (i < nPeaks) {
if (i < nGoodPeaks) {
peak.setCentre(eExcitations.get(i));
peak.setIntensity(iExcitations.get(i));
} else {
peak.setHeight(0.0);
peak.fixAll();
if (i > nOriginalPeaks) {
peak.fixAll();
}
}
}
return true;
}
} // namespace Functions
......
......@@ -37,8 +37,6 @@ void CrystalFieldSpectrum::buildTargetFunction() const {
FunctionValues values;
m_source->function(domain, values);
size_t maxNPeaks = m_source->getAttribute("MaxPeakCount").asInt();
if (values.size() == 0) {
return;
}
......@@ -51,6 +49,7 @@ void CrystalFieldSpectrum::buildTargetFunction() const {
auto peakShape = IFunction::getAttribute("PeakShape").asString();
auto fwhm = IFunction::getAttribute("FWHM").asDouble();
auto nPeaks = values.size() / 2;
size_t maxNPeaks = nPeaks + nPeaks / 2 + 1;
for (size_t i = 0; i < maxNPeaks; ++i) {
auto fun = API::FunctionFactory::Instance().createFunction(peakShape);
auto peak = boost::dynamic_pointer_cast<API::IPeakFunction>(fun);
......@@ -83,12 +82,18 @@ void CrystalFieldSpectrum::updateTargetFunction() const {
m_source->function(domain, values);
auto nPeaks = values.size() / 2;
size_t maxNPeaks = m_source->getAttribute("MaxPeakCount").asInt();
size_t maxNPeaks = nPeaks + nPeaks / 2 + 1;
auto spectrum = dynamic_cast<CompositeFunction *>(m_target.get());
if (!spectrum || spectrum->nFunctions() != maxNPeaks) {
if (!spectrum) {
buildTargetFunction();
return;
}
if (spectrum->nFunctions() != maxNPeaks) {
maxNPeaks = spectrum->nFunctions();
if (nPeaks > maxNPeaks) {
nPeaks = maxNPeaks;
}
}
for (size_t i = 0; i < maxNPeaks; ++i) {
auto fun = spectrum->getFunction(i);
......
......@@ -147,6 +147,48 @@ public:
TS_ASSERT_EQUALS(fun->getParameter("f1.f0.f0.Sigma"), 0.2);
}
void test_composite_multispectral() {
std::string fun1 =
"name=CrystalFieldMultiSpectrum,Ion=Ce,Temperatures=(44, "
"50),ToleranceIntensity=0.001,B20=0.37737,B22=3.9770,"
"B40=-0.031787,B42=-0.11611,B44=-0.12544,"
"f0.f1.FWHM=1.6,f0.f2.FWHM=2.0,f0.f3.FWHM=2.3,f1.f1.FWHM=1.6,"
"f1.f2.FWHM=2.0,f1.f3.FWHM=2.3";
std::string fun2 =
"name=CrystalFieldMultiSpectrum,Ion=Pr,Temperatures=(44, "
"50),ToleranceIntensity=0.001,B20=0.37737,B22=3.9770,"
"B40=-0.031787,B42=-0.11611,B44=-0.12544,"
"f0.f1.FWHM=1.6,f0.f2.FWHM=2.0,f0.f3.FWHM=2.3,f1.f1.FWHM=1.6,"
"f1.f2.FWHM=2.0,f1.f3.FWHM=2.3";
auto fun = fun1 + ";" + fun2;
auto ws = createWorkspace();
auto alg = AlgorithmFactory::Instance().create("EvaluateFunction", -1);
alg->initialize();
alg->setPropertyValue("Function", fun);
alg->setProperty("InputWorkspace", ws);
alg->setProperty("InputWorkspace_1", ws);
alg->setProperty("OutputWorkspace", "out");
alg->execute();
TS_ASSERT(alg->isExecuted());
auto out = AnalysisDataService::Instance().retrieveWS<MatrixWorkspace>(
"Workspace_0");
TS_ASSERT(out);
TS_ASSERT_EQUALS(out->getNumberHistograms(), 3);
TS_ASSERT_DELTA(out->readY(1)[0], 2.9202, 0.001);
TS_ASSERT_DELTA(out->readY(1)[1], 2.4691, 0.001);
TS_ASSERT_DELTA(out->readY(1)[2], 1.3817, 0.001);
out = AnalysisDataService::Instance().retrieveWS<MatrixWorkspace>(
"Workspace_1");
TS_ASSERT(out);
TS_ASSERT_EQUALS(out->getNumberHistograms(), 3);
TS_ASSERT_DELTA(out->readY(1)[0], 2.9192, 0.001);
TS_ASSERT_DELTA(out->readY(1)[1], 2.4647, 0.001);
TS_ASSERT_DELTA(out->readY(1)[2], 1.3791, 0.001);
AnalysisDataService::Instance().clear();
}
private:
Workspace_sptr createWorkspace() {
auto ws = WorkspaceFactory::Instance().create("Workspace2D", 1, 100, 100);
......
......@@ -37,6 +37,8 @@ class CrystalField(object):
default_peakShape = 'Gaussian'
default_background = 'FlatBackground'
default_spectrum_size = 200
field_parameter_names = ['BmolX','BmolY','BmolZ','BextX','BextY','BextZ',
......@@ -116,9 +118,6 @@ class CrystalField(object):
self._intensityScaling = 1.0
self._resolutionModel = None
self.peaks = PeaksFunction()
self.background = None
for key in kwargs:
if key == 'ToleranceEnergy':
self._toleranceEnergy = kwargs[key]
......@@ -136,6 +135,12 @@ class CrystalField(object):
# Crystal field parameters
self._fieldParameters[key] = kwargs[key]
if isinstance(self._temperature, list) or isinstance(self._temperature, np.ndarray):
self.peaks = [PeaksFunction() for t in self._temperature]
else:
self.peaks = PeaksFunction()
self.background = None
# Eigensystem
self._dirty_eigensystem = True
self._eigenvalues = None
......@@ -541,6 +546,7 @@ class CrystalField(object):
Update values of the fitting parameters in case of a multi-spectrum function.
@param func: A IFunction object containing new parameter values.
"""
from .function import Function
n = func.nParams()
for i in range(n):
par = func.parameterName(i)
......@@ -551,6 +557,8 @@ class CrystalField(object):
ipeak = int(m.group(2))
par = m.group(3)
if ipeak == 0:
if self.background is None:
self.setBackground(background=Function(self.default_background))
background = self.background[ispec]
mb = re.match(fn_pattern, par)
if mb:
......@@ -689,6 +697,13 @@ class CrystalFieldMulti(object):
fun += ',ties=(%s)' % ties
return fun
def makeMultiSpectrumFunction(self):
fun = ';'.join([a.makeMultiSpectrumFunction() for a in self.args])
ties = self.getTies()
if len(ties) > 0:
fun += ',ties=(%s)' % ties
return fun
def ties(self, **kwargs):
"""Set ties on the parameters."""
for tie in kwargs:
......@@ -726,6 +741,12 @@ class CrystalFieldMulti(object):
for i in range(n):
self.args[i].update(func[i])
def update_multi(self, func):
n = func.nFunctions()
assert n == len(self.args)
for i in range(n):
self.args[i].update_multi(func[i])
def __add__(self, other):
if isinstance(other, CrystalFieldMulti):
return CrystalFieldMulti(*(self.args + other.args))
......@@ -744,7 +765,7 @@ class CrystalFieldFit(object):
Object that controls fitting.
"""
def __init__(self, Model=None, Temperature=None, FWHM=None, InputWorkspace=None):
def __init__(self, Model=None, Temperature=None, FWHM=None, InputWorkspace=None, **kwargs):
self.model = Model
if Temperature is not None:
self.model.Temperature = Temperature
......@@ -752,6 +773,7 @@ class CrystalFieldFit(object):
self.model.FWHM = FWHM
self._input_workspace = InputWorkspace
self._output_workspace_base_name = 'fit'
self._fit_properties = kwargs
def fit(self):
"""
......@@ -773,6 +795,7 @@ class CrystalFieldFit(object):
alg.setProperty('Function', fun)
alg.setProperty('InputWorkspace', self._input_workspace)
alg.setProperty('Output', 'fit')
self._set_fit_properties(alg)
alg.execute()
f = alg.getProperty('Function').value
self.model.update(f)
......@@ -792,6 +815,11 @@ class CrystalFieldFit(object):
alg.setProperty('InputWorkspace_%s' % i, ws)
i += 1
alg.setProperty('Output', 'fit')
self._set_fit_properties(alg)
alg.execute()
f = alg.getProperty('Function').value
self.model.update_multi(f)
def _set_fit_properties(self, alg):
for prop in self._fit_properties.items():
alg.setProperty(*prop)
......@@ -178,11 +178,11 @@ class CrystalFieldTests(unittest.TestCase):
self.assertAlmostEqual(y[63], 7.63981716, 8)
self.assertAlmostEqual(y[64], 4.08015236, 8)
x, y = cf.getSpectrum(1)
self.assertAlmostEqual(y[45], 0.29653544, 8)
self.assertAlmostEqual(y[46], 0.45996719, 8)
self.assertAlmostEqual(y[47], 0.65873995, 8)
self.assertAlmostEqual(y[48], 0.6924739, 8)
self.assertAlmostEqual(y[49], 0.51119472, 8)
self.assertAlmostEqual(y[45], 0.29822612216224065, 8)
self.assertAlmostEqual(y[46], 0.46181038787922241, 8)
self.assertAlmostEqual(y[47], 0.66075719314988057, 8)
self.assertAlmostEqual(y[48], 0.69469096259927476, 8)
self.assertAlmostEqual(y[49], 0.51364268980567007, 8)
def test_api_CrystalField_spectrum_from_list(self):
from CrystalField import CrystalField
......@@ -210,11 +210,11 @@ class CrystalFieldTests(unittest.TestCase):
self.assertEqual(x[3], 3.0)
self.assertEqual(x[4], 3.85)
self.assertAlmostEqual(y[0], 6.3037011243626706, 8)
self.assertAlmostEqual(y[1], 0.33097708354212135, 8)
self.assertAlmostEqual(y[2], 1.2245308604241509, 8)
self.assertAlmostEqual(y[3], 0.078438334846419058, 8)
self.assertAlmostEqual(y[4], 2.637989456640248, 8)
self.assertAlmostEqual(y[0], 6.3046701386938624, 8)
self.assertAlmostEqual(y[1], 0.33121919026244667, 8)
self.assertAlmostEqual(y[2], 1.2246681560002572, 8)
self.assertAlmostEqual(y[3], 0.078541076629159004, 8)
self.assertAlmostEqual(y[4], 2.6380618652343704, 8)
def test_api_CrystalField_spectrum_0(self):
from CrystalField import CrystalField
......@@ -247,11 +247,11 @@ class CrystalFieldTests(unittest.TestCase):
self.assertAlmostEqual(y[15], 0.050129858433581413, 8)
self.assertAlmostEqual(y[16], 0.054427788297191478, 8)
x, y = cf.getSpectrum(1, workspace)
self.assertAlmostEqual(y[0], 6.3037011243626706, 8)
self.assertAlmostEqual(y[1], 4.2744246158088393, 8)
self.assertAlmostEqual(y[2], 2.1770150365787457, 8)
self.assertAlmostEqual(y[3], 1.2003766255981201, 8)
self.assertAlmostEqual(y[4], 0.73968415199300741, 8)
self.assertAlmostEqual(y[0], 6.3046701386938624, 8)
self.assertAlmostEqual(y[1], 4.2753076741531455, 8)
self.assertAlmostEqual(y[2], 2.1778230746690772, 8)
self.assertAlmostEqual(y[3], 1.2011188019120242, 8)
self.assertAlmostEqual(y[4], 0.74036819427919942, 8)
x, y = cf.getSpectrum(workspace)
self.assertAlmostEqual(y[0], 12.474954833565066, 8)
self.assertAlmostEqual(y[1], 4.3004160689570403, 8)
......@@ -401,9 +401,6 @@ class CrystalFieldFitTest(unittest.TestCase):
self.assertAlmostEqual(cf.peaks.param[4]['PeakCentre'], 0.0, 8)
self.assertAlmostEqual(cf.peaks.param[4]['FWHM'], 1.0, 8)
self.assertAlmostEqual(cf.peaks.param[4]['Amplitude'], 0.0, 8)
self.assertAlmostEqual(cf.peaks.param[5]['PeakCentre'], 0.0, 8)
self.assertAlmostEqual(cf.peaks.param[5]['FWHM'], 1.0, 8)
self.assertAlmostEqual(cf.peaks.param[5]['Amplitude'], 0.0, 8)
def test_CrystalFieldFit_multi_spectrum(self):
......@@ -449,39 +446,39 @@ class CrystalFieldFitTest(unittest.TestCase):
fit = CrystalFieldFit(cf, InputWorkspace=[ws0, ws1])
fit.fit()
self.assertAlmostEqual(cf.background[0].peak.param['PeakCentre'], -0.0026088250802930577, 8)
self.assertAlmostEqual(cf.background[0].peak.param['Sigma'], 0.30300200254659215, 8)
self.assertAlmostEqual(cf.background[0].peak.param['Height'], 10.555088363800422, 8)
self.assertAlmostEqual(cf.background[0].background.param['A0'], 0.9933917156353519, 8)
self.assertAlmostEqual(cf.background[0].peak.param['PeakCentre'], -0.0026425783022528374, 8)
self.assertAlmostEqual(cf.background[0].peak.param['Sigma'], 0.30305576117004407, 8)
self.assertAlmostEqual(cf.background[0].peak.param['Height'], 10.55973474260564, 8)
self.assertAlmostEqual(cf.background[0].background.param['A0'], 0.9933285059936604, 8)
self.assertAlmostEqual(cf.background[1].peak.param['PeakCentre'], -0.5645843235671519, 8)
self.assertAlmostEqual(cf.background[1].peak.param['Sigma'], 0.990615576443698, 8)
self.assertAlmostEqual(cf.background[1].peak.param['Height'], 17.506721782939895, 8)
self.assertAlmostEqual(cf.background[1].background.param['A0'], 1.1006327347548288, 8)
self.assertAlmostEqual(cf.background[1].peak.param['PeakCentre'], -0.5641183841947577, 8)
self.assertAlmostEqual(cf.background[1].peak.param['Sigma'], 0.9904967068950424, 8)
self.assertAlmostEqual(cf.background[1].peak.param['Height'], 17.498648691583526, 8)
self.assertAlmostEqual(cf.background[1].background.param['A0'], 1.1005498193552903, 8)
self.assertAlmostEqual(cf.peaks[0].param[1]['PeakCentre'], 44.52027265331026, 8)
self.assertAlmostEqual(cf.peaks[0].param[1]['FWHM'], 2.9042743044799426, 8)
self.assertAlmostEqual(cf.peaks[0].param[1]['Amplitude'], 0.7468047008285558, 8)
self.assertAlmostEqual(cf.peaks[0].param[1]['PeakCentre'], 44.47522063599522, 8)
self.assertAlmostEqual(cf.peaks[0].param[1]['FWHM'], 2.9482185926005435, 8)
self.assertAlmostEqual(cf.peaks[0].param[1]['Amplitude'], 0.7481078193200019, 8)
self.assertAlmostEqual(cf.peaks[0].param[2]['PeakCentre'], 29.32253802971275, 8)
self.assertAlmostEqual(cf.peaks[0].param[2]['FWHM'], 1.1315246964379955, 8)
self.assertAlmostEqual(cf.peaks[0].param[2]['Amplitude'], 0.7452080584249358, 8)
self.assertAlmostEqual(cf.peaks[0].param[2]['PeakCentre'], 29.322977342784924, 8)
self.assertAlmostEqual(cf.peaks[0].param[2]['FWHM'], 1.1355719580493524, 8)
self.assertAlmostEqual(cf.peaks[0].param[2]['Amplitude'], 0.7471598262528436, 8)
self.assertAlmostEqual(cf.peaks[0].param[3]['PeakCentre'], 0.0, 8)
self.assertAlmostEqual(cf.peaks[0].param[3]['FWHM'], 0.0, 8)
self.assertAlmostEqual(cf.peaks[0].param[3]['PeakCentre'], 23.99290823852288, 8)
self.assertAlmostEqual(cf.peaks[0].param[3]['FWHM'], 0.4060486307010578, 8)
self.assertAlmostEqual(cf.peaks[0].param[3]['Amplitude'], 0.0, 8)
self.assertAlmostEqual(cf.peaks[1].param[1]['PeakCentre'], 44.52027265331026, 8)
self.assertAlmostEqual(cf.peaks[1].param[1]['FWHM'], 7.803607774897451, 8)
self.assertAlmostEqual(cf.peaks[1].param[1]['Amplitude'], 0.7462867984192054, 8)
self.assertAlmostEqual(cf.peaks[1].param[1]['PeakCentre'], 44.47522063599522, 8)
self.assertAlmostEqual(cf.peaks[1].param[1]['FWHM'], 7.83698643366528, 8)
self.assertAlmostEqual(cf.peaks[1].param[1]['Amplitude'], 0.7475888752255804, 8)
self.assertAlmostEqual(cf.peaks[1].param[2]['PeakCentre'], 29.32253802971275, 8)
self.assertAlmostEqual(cf.peaks[1].param[2]['FWHM'], 1.130386130370342, 8)
self.assertAlmostEqual(cf.peaks[1].param[2]['Amplitude'], 0.7446912632728735, 8)
self.assertAlmostEqual(cf.peaks[1].param[2]['PeakCentre'], 29.322977342784924, 8)
self.assertAlmostEqual(cf.peaks[1].param[2]['FWHM'], 1.1343917778015058, 8)
self.assertAlmostEqual(cf.peaks[1].param[2]['Amplitude'], 0.7466415397580233, 8)
self.assertAlmostEqual(cf.peaks[1].param[3]['PeakCentre'], 15.197734623597508, 8)
self.assertAlmostEqual(cf.peaks[1].param[3]['FWHM'], 25.956468289206597, 8)
self.assertAlmostEqual(cf.peaks[1].param[3]['Amplitude'], 0.0015178943479789188, 8)
self.assertAlmostEqual(cf.peaks[1].param[3]['PeakCentre'], 15.152243293210297, 8)
self.assertAlmostEqual(cf.peaks[1].param[3]['FWHM'], 25.947899288913575, 8)
self.assertAlmostEqual(cf.peaks[1].param[3]['Amplitude'], 0.0015192409710355985, 8)
def test_CrystalFieldFit_multi_spectrum_simple_background(self):
......@@ -562,13 +559,13 @@ class CrystalFieldFitTest(unittest.TestCase):
fit = CrystalFieldFit(cf, InputWorkspace=[ws0, ws1])
fit.fit()
self.assertAlmostEqual(cf.background[0].peak.param['PeakCentre'], -0.00831074804991939, 8)
self.assertAlmostEqual(cf.background[0].peak.param['Sigma'], 0.3121248098926919, 8)
self.assertAlmostEqual(cf.background[0].peak.param['Height'], 10.967831390443026, 8)
self.assertAlmostEqual(cf.background[0].peak.param['PeakCentre'], -0.00018532664109938181, 8)
self.assertAlmostEqual(cf.background[0].peak.param['Sigma'], 0.30014532208212547, 8)
self.assertAlmostEqual(cf.background[0].peak.param['Height'], 10.030054969751856, 8)
self.assertAlmostEqual(cf.background[1].peak.param['PeakCentre'], -0.4606116667468967, 8)
self.assertAlmostEqual(cf.background[1].peak.param['Sigma'], 0.9582966759931785, 8)
self.assertAlmostEqual(cf.background[1].peak.param['Height'], 15.963129385365477, 8)
self.assertAlmostEqual(cf.background[1].peak.param['PeakCentre'], -0.6266190281804616, 8)
self.assertAlmostEqual(cf.background[1].peak.param['Sigma'], 1.009442002821441, 8)
self.assertAlmostEqual(cf.background[1].peak.param['Height'], 18.497799303486225, 8)
def test_multi_ion_single_spectrum(self):
from CrystalField.fitting import MakeWorkspace
......@@ -591,17 +588,50 @@ class CrystalFieldFitTest(unittest.TestCase):
fit = CrystalFieldFit(Model=cf, InputWorkspace=ws)
fit.fit()
self.assertAlmostEqual(cf[0].param['B20'], 0.48824387792806034, 8)
self.assertAlmostEqual(cf[0].param['B22'], 3.977892421652792, 8)
self.assertAlmostEqual(cf[0].param['B40'], -0.024294997727440197, 8)
self.assertAlmostEqual(cf[0].param['B42'], -0.14604639435112998, 8)
self.assertAlmostEqual(cf[0].param['B44'], -0.072762568202823, 8)
self.assertAlmostEqual(cf[1].param['B20'], 0.43777220488713825, 8)
self.assertAlmostEqual(cf[1].param['B22'], 3.9635529299402585, 8)
self.assertAlmostEqual(cf[1].param['B40'], -0.03319808647891935, 8)
self.assertAlmostEqual(cf[1].param['B42'], -0.12160876404778258, 8)
self.assertAlmostEqual(cf[1].param['B44'], -0.11549947832481772, 8)
self.assertAlmostEqual(cf[0].param['B20'], 0.36568328044388737, 8)
self.assertAlmostEqual(cf[0].param['B22'], 3.9154384449144914, 8)
self.assertAlmostEqual(cf[0].param['B40'], -0.04489701572952994, 8)
self.assertAlmostEqual(cf[0].param['B42'], -0.07106796874986179, 8)
self.assertAlmostEqual(cf[0].param['B44'], -0.1236384086842752, 8)