Dave Mitchell's DigitalMicrograph™ Scripting Website

// Function to fit a 2nd order polynomial to the three points of a peak (pre maximum, maximum and post maximum)

// and return the position of the maximum (or minimum). Values entered are the x1, x2, x3 coordinates and their

// corresponding values f1, f2, f3. The position of the maximum (or minimum) and its value are passed into the

// function as references to xmax and fmax respectively.

// D. R. G. Mitchell, adminnospam@dmscripting.com (remove the nospam to make this email address work)

// version:20070528

// v1.0, May 28th 2007.

// Thanks to Dr Bill Bertram of ANSTO for providing the algorithm

void interpolatemaximum(number x1, number f1, number x2, number f2, number x3, number f3, number &xmax, number &fmax)

{

number a, b, c // The three variables in the polynomial y=a+bx+cx^2

// Error trap - if the three intensity values are the same, there is no inflection point in the curve, and an error would

// occur. If this happens, the inflection point/intensity value are set to the middle point

if(f1==f2 && f2==f3)

{

xmax=x2

fmax=f2

return

}

// Calculate a

number temp1=(f3*x1*(x1-x2)*x2)

number temp2=(f1*x2*(x2-x3))+(f2*x1*(-x1+x3))

number temp3=((x1-x2)*(x1-x3)*(x2-x3))

a=(temp1+(x3*temp2))/temp3

// calculate b

b=-((f2*x1**2)-(f3*x1**2)-(f1*x2**2)+(f3*x2**2)+(f1*x3**2)-(f2*x3**2))/((x1-x2)*(x1-x3)*(-x2+x3))

//calculate c

c=-((f2*x1)-(f3*x1)-(f1*x2)+(f3*x2)+(f1*x3)-(f2*x3))/((x2-x3)*(x1**2-(x1*x2)-(x1*x3)+(x2*x3)))

// Calculate the position of the curve inflection and value

xmax=-b/(2*c)

fmax=a+(b*xmax)+(c*xmax**2)

}

// Main script starts here - some test data

number x1, x2, x3, f1, f2, f3, xmax, fmax

x1=10

x2=20

x3=30

f1=2

f2=3

f3=3.1

interpolatemaximum(x1, f1, x2, f2, x3, f3, xmax, fmax)

result("\n xmax="+xmax+" ymax="+fmax)