Skip to content

Commit

Permalink
Re-implemented Voigt as a peak function. Refs #5923
Browse files Browse the repository at this point in the history
  • Loading branch information
martyngigg committed Oct 30, 2012
1 parent 02ca561 commit cf7de44
Show file tree
Hide file tree
Showing 3 changed files with 132 additions and 17 deletions.
26 changes: 19 additions & 7 deletions Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/Voigt.h
Original file line number Diff line number Diff line change
Expand Up @@ -2,8 +2,8 @@
#define MANTID_CURVEFITTING_VOIGT_H_

#include "MantidCurveFitting/DllConfig.h"
#include "MantidAPI/IFunction1D.h"
#include "MantidAPI/ParamFunction.h"
#include "MantidAPI/IPeakFunction.h"


namespace Mantid
{
Expand Down Expand Up @@ -34,23 +34,35 @@ namespace Mantid
File change history is stored at: <https://svn.mantidproject.org/mantid/trunk/Code/Mantid>
Code Documentation is available at: <http://doxygen.mantidproject.org>
*/
class MANTID_CURVEFITTING_DLL Voigt :
public virtual API::ParamFunction, public virtual API::IFunction1D
class MANTID_CURVEFITTING_DLL Voigt : public API::IPeakFunction
{
private:
/// Return a string identifier for the function
std::string name() const { return "Voigt"; }
/// Declare parameters
void declareParameters();

/// Fill out with function values at given x points
void function1D(double *out, const double *xValues, const size_t nData) const;
void functionLocal(double *out, const double *xValues, const size_t nData) const;
/// Derivatives of function with respect to active parameters
void functionDeriv1D(API::Jacobian* out, const double* xValues, const size_t nData);
void functionDerivLocal(API::Jacobian* out, const double* xValues, const size_t nData);

/// Calculate both function & derivative together
void calculateFunctionAndDerivative(const double *xValues, const size_t nData,
double *functionValues, API::Jacobian * derivatives) const;

/// Return value of centre of peak
double centre()const;
/// Return value of height of peak
double height()const;
/// Return value of FWHM of peak
double fwhm()const;
/// Set the centre of the peak
void setCentre(const double value);
/// Set the height of the peak
void setHeight(const double value);
/// Set the FWHM of the peak
void setFwhm(const double value);
};


Expand Down
80 changes: 70 additions & 10 deletions Code/Mantid/Framework/CurveFitting/src/Voigt.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -24,6 +24,11 @@ namespace Mantid
const double COEFFC[NLORENTZIANS] = {-0.3085, 0.5906, -0.3085, 0.5906};
const double COEFFD[NLORENTZIANS] = {0.0210, -1.1858, -0.0210, 1.1858};

const char * LORENTZ_AMP = "LorentzAmp";
const char * LORENTZ_POS = "LorentzPos";
const char * LORENTZ_FWHM = "LorentzFWHM";
const char * GAUSSIAN_FWHM = "GaussianFWHM";

const double SQRTLN2 = std::sqrt(std::log(2.0));
const double SQRTPI = std::sqrt(M_PI);
///@endcond
Expand All @@ -34,10 +39,10 @@ namespace Mantid
*/
void Voigt::declareParameters()
{
declareParameter("LorentzAmp", 0.0, "Value of the Lorentzian amplitude");
declareParameter("LorentzPos", 0.0, "Position of the Lorentzian peak");
declareParameter("LorentzFWHM", 0.0, "Value of the full-width half-maximum for the Lorentzian");
declareParameter("GaussianFWHM", 0.0, "Value of the full-width half-maximum for the Gaussian");
declareParameter(LORENTZ_AMP, 0.0, "Value of the Lorentzian amplitude");
declareParameter(LORENTZ_POS, 0.0, "Position of the Lorentzian peak");
declareParameter(LORENTZ_FWHM, 0.0, "Value of the full-width half-maximum for the Lorentzian");
declareParameter(GAUSSIAN_FWHM, 0.0, "Value of the full-width half-maximum for the Gaussian");
}

/**
Expand All @@ -46,7 +51,7 @@ namespace Mantid
* @param xValues :: The X values
* @param nData :: The number of X values to evaluate
*/
void Voigt::function1D(double *out, const double *xValues, const size_t nData) const
void Voigt::functionLocal(double *out, const double *xValues, const size_t nData) const
{
calculateFunctionAndDerivative(xValues,nData,out,NULL);
}
Expand All @@ -57,7 +62,7 @@ namespace Mantid
* @param xValues :: The X values
* @param nData :: The number of X values to evaluate
*/
void Voigt::functionDeriv1D(API::Jacobian* out, const double* xValues, const size_t nData)
void Voigt::functionDerivLocal(API::Jacobian* out, const double* xValues, const size_t nData)
{
calculateFunctionAndDerivative(xValues, nData, NULL, out);
}
Expand All @@ -72,10 +77,10 @@ namespace Mantid
void Voigt::calculateFunctionAndDerivative(const double *xValues, const size_t nData,
double *functionValues, API::Jacobian * derivatives) const
{
const double a_L = getParameter("LorentzAmp");
const double lorentzPos = getParameter("LorentzPos");
const double gamma_L = getParameter("LorentzFWHM");
const double gamma_G = getParameter("GaussianFWHM");
const double a_L = getParameter(LORENTZ_AMP);
const double lorentzPos = getParameter(LORENTZ_POS);
const double gamma_L = getParameter(LORENTZ_FWHM);
const double gamma_G = getParameter(GAUSSIAN_FWHM);

const double rtln2oGammaG = SQRTLN2/gamma_G;
const double prefactor = (a_L*SQRTPI*gamma_L*SQRTLN2/gamma_G);
Expand Down Expand Up @@ -113,6 +118,61 @@ namespace Mantid
}
}

/**
* Returns the value of the "LorentzPos" parameter.
* @return value of centre of peak
*/
double Voigt::centre()const
{
return getParameter(LORENTZ_POS);
}

/**
* Return the value of the "LorentzAmp" parameter
* @return value of height of peak
*/
double Voigt::height()const
{
return 2.0*getParameter(LORENTZ_AMP);
}

/**
* Gives the FWHM of the peak. This is estimated as
* 0.5*(LorentzFWHM + GaussianFWHM)
* @return value of FWHM of peak
*/
double Voigt::fwhm()const
{
return (getParameter(LORENTZ_FWHM) + getParameter(GAUSSIAN_FWHM));
}

/**
* Set the centre of the peak, the LorentzPos parameter
* @param value :: The new value for the centre of the peak
*/
void Voigt::setCentre(const double value)
{
this->setParameter(LORENTZ_POS, value);
}

/**
* Set the height of the peak, the LorentzAmp parameter
* @param value :: The new value for the centre of the peak
*/
void Voigt::setHeight(const double value)
{
this->setParameter(LORENTZ_AMP, 0.5*value);
}

/**
* Set the FWHM of the peak
* @param value :: The new value for the FWHM of the peak
*/
void Voigt::setFwhm(const double value)
{
this->setParameter(LORENTZ_FWHM, 0.5*value);
this->setParameter(GAUSSIAN_FWHM, 0.5*value);
}

} // namespace CurveFitting
} // namespace Mantid
43 changes: 43 additions & 0 deletions Code/Mantid/Framework/CurveFitting/test/VoigtTest.h
Original file line number Diff line number Diff line change
Expand Up @@ -51,6 +51,7 @@ class VoigtTest : public CxxTest::TestSuite
{
const double a_L(5), pos(-1), gamma_L(0.9), gamma_G(0.1);
auto voigtFn = createFunction(a_L,pos,gamma_L,gamma_G);


Mantid::API::FunctionValues outputs(*m_domain);
voigtFn->function(*m_domain, outputs);
Expand Down Expand Up @@ -92,7 +93,46 @@ class VoigtTest : public CxxTest::TestSuite
TS_ASSERT_DELTA(dxDa[i][j], jacobian.get(i,j), 1e-10);
}
}
}

void test_function_is_a_peak_function()
{
const double a_L(5), pos(-1), gamma_L(0.9), gamma_G(0.1);
auto voigtFn = createFunction(a_L,pos,gamma_L,gamma_G);

auto peakFn = boost::dynamic_pointer_cast<Mantid::API::IPeakFunction>(voigtFn);
TSM_ASSERT("Voigt function should be a PeakFunction", peakFn);
}

void test_peak_functions_return_expected_results()
{
const double a_L(5), pos(-1), gamma_L(0.9), gamma_G(0.1);
auto voigtFn = createFunction(a_L,pos,gamma_L,gamma_G);
auto peakFn = boost::dynamic_pointer_cast<Mantid::API::IPeakFunction>(voigtFn);

TS_ASSERT_DELTA(peakFn->centre(), pos, 1e-12);
TS_ASSERT_DELTA(peakFn->height(), 2.0*a_L, 1e-12);
TS_ASSERT_DELTA(peakFn->fwhm(), (gamma_L+gamma_G), 1e-12);
}

void test_setting_peak_functions_set_expected_parameters()
{
double a_L(5), pos(-1), gamma_L(0.9), gamma_G(0.1);
auto voigtFn = createFunction(a_L,pos,gamma_L,gamma_G);
auto peakFn = boost::dynamic_pointer_cast<Mantid::API::IPeakFunction>(voigtFn);

pos = 1.2;
peakFn->setCentre(pos);
TS_ASSERT_DELTA(peakFn->centre(), pos, 1e-12);

a_L = 3.5;
peakFn->setHeight(a_L);
TS_ASSERT_DELTA(peakFn->height(), a_L, 1e-12);

gamma_L = 1.2;
gamma_G = 0.4;
peakFn->setFwhm(gamma_L+gamma_G);
TS_ASSERT_DELTA(peakFn->fwhm(), (gamma_L+gamma_G), 1e-12);
}

private:
Expand All @@ -101,6 +141,9 @@ class VoigtTest : public CxxTest::TestSuite
const double gamma_L, const double gamma_G)
{
boost::shared_ptr<IFunction> voigtFn = boost::make_shared<Voigt>();
auto peakFn = boost::dynamic_pointer_cast<Mantid::API::IPeakFunction>(voigtFn);
// Set a fairly wide radius for simple tests
peakFn->setPeakRadius(10);
voigtFn->initialize();

voigtFn->setParameter("LorentzAmp", a_L);
Expand Down

0 comments on commit cf7de44

Please sign in to comment.