$$\renewcommand\AA{\unicode{x212B}}$$

# Exercise 6 Solutions¶

## Simple 1D¶

The aim of this exercise is to implement a Lorentz function, that has no derivative.

from mantid.api import *
import math
import numpy as np

INVERSE_PI = 1.0/math.pi

class Lorentz(IFunction1D):

def init(self):
# Tell Mantid about the 3 parameters that are involved in the function
self.declareParameter("Amplitude", 0.0)
self.declareParameter("PeakCentre", 0.0)
self.declareParameter("Gamma", 0.0)

def function1D(self, xvals):
# xvals is a 1D numpy array that contains the X values for the defined fitting range.

# Get the current values of the 3 parameters
amp = self.getParameterValue("Amplitude") # Access current value during the fit
half_gamma = 0.5*self.getParameterValue("Gamma")
c = self.getParameterValue("PeakCentre")

denom = (xvals - c)**2 + half_gamma*half_gamma
return amp*INVERSE_PI*half_gamma/denom

# Register with Mantid
FunctionFactory.subscribe(Lorentz)


## Analytical Derivative¶

The aim of this exercise is to implement a Lorentz function, that contains a derivative.

from mantid.api import *
import math
import numpy as np

INVERSE_PI = 1.0/math.pi

class LorentzWithDeriv(IFunction1D):

def init(self):
# Tell Mantid about the 3 parameters that are involved in the function
self.declareParameter("Amplitude", 1.0)
self.declareParameter("PeakCentre", 0.0)
self.declareParameter("Gamma", 0.1)

def function1D(self, xvals):
# xvals is a 1D numpy array that contains the X values for the defined fitting range.

# Get the current values of the 3 parameters
amp = self.getParameterValue("Amplitude") # Access current value during the fit
half_gamma = 0.5*self.getParameterValue("Gamma")
peak_centre = self.getParameterValue("PeakCentre")

denom = (xvals - peak_centre)**2 + half_gamma*half_gamma
return amp*INVERSE_PI*half_gamma/denom

def functionDeriv1D(self, xvals, out):
# xvals is a 1D numpy array that contains the X values for the defined fitting range.
# out is a Jacobian matrix object. Mantid expects the partial derivatives
# w.r.t the parameters and x values to be stored here
# Get the current values of the 3 parameters

amplitude = self.getParameterValue("Amplitude")
peakCentre = self.getParameterValue("PeakCentre")
gamma = self.getParameterValue("Gamma")
halfGamma = 0.5*gamma

for i, xval in enumerate(xvals):
diff = xval-peakCentre
invDen1 = 1.0/(gamma*gamma + 4.0*diff*diff)
dfda = 2.0*INVERSE_PI*gamma*invDen1
out.set(i, 0, dfda)

invDen2 =  1/(diff*diff + halfGamma*halfGamma)
dfdxo = amplitude*INVERSE_PI*gamma*diff*invDen2*invDen2
out.set(i, 1, dfdxo);

dfdg = -2.0*amplitude*INVERSE_PI*(gamma*gamma - 4.0*diff*diff)*invDen1*invDen1
out.set(i, 2, dfdg)

# Register with Mantid
FunctionFactory.subscribe(LorentzWithDeriv)


## Peak Function¶

The aim of this exercise is to implement a peak fit function function to fit.

from mantid.api import *
import numpy as np
import math

INVERSE_PI = 1.0/math.pi

class LorentzPeak(IPeakFunction):

def init(self):
# Tell Mantid about the 3 parameters that are involved in the function
self.declareParameter("Amplitude", 0.0)
self.declareParameter("PeakCentre", 0.0)
self.declareParameter("Gamma", 0.0)

def functionLocal(self, xvals):
# xvals is a 1D numpy array that contains the X values for the defined fitting range.
half_gamma = 0.5*self.getParameterValue("Gamma")
denom = (xvals - self.getParameterValue("PeakCentre"))**2 + half_gamma*half_gamma
return self.getParameterValue("Amplitude")*INVERSE_PI*half_gamma/denom

def functionDerivLocal(self, xvals, out):
# xvals is a 1D numpy array that contains the X values for the defined fitting range.
# out is a Jacobian matrix object. Mantid expects the partial derivatives
# w.r.t the parameters and x values to be stored here

# Get the current parameter values
amplitude = self.getParameterValue("Amplitude")
peakCentre = self.getParameterValue("PeakCentre")
gamma = self.getParameterValue("Gamma")
halfGamma = 0.5*gamma

for i, xval in enumerate(xvals):
diff = xval-peakCentre
invDen1 = 1.0/(gamma*gamma + 4.0*diff*diff)
dfda = 2.0*INVERSE_PI*gamma*invDen1
out.set(i, 0, dfda)

invDen2 =  1/(diff*diff + halfGamma*halfGamma)
dfdxo = amplitude*INVERSE_PI*gamma*diff*invDen2*invDen2
out.set(i, 1, dfdxo);

dfdg = -2.0*amplitude*INVERSE_PI*(gamma*gamma - 4.0*diff*diff)*invDen1*invDen1
out.set(i, 2, dfdg)

def centre(self):
# Return a guess at the centre
return self.getParameterValue("PeakCentre")

def height(self):
# Return a guess at the height
return self.getParameterValue("Amplitude")

def fwhm(self):
# Return a guess at the FWHM
return 2*self.getParameterValue("Gamma")

def setCentre(self, new_centre):
# Update centre guess when a new value is chosen from GUI
self.setParameter("PeakCentre", new_centre)

def setHeight(self, new_height):
# Update Amplitude guess when a new height is chosen from GUI
self.setParameter("Amplitude", new_height)

def setFwhm(self, new_fwhm):
# Update Gamma guess when a new width is chosen from GUI
self.setParameter("Gamma", new_fwhm/2.0)

# Register function with Mantid
FunctionFactory.subscribe(LorentzPeak)