// This macro shows how to estimate the resolution of an image of a sharp edge
// from a line profile across that edge and fitting it with the error function
// (This assumes that the point spread function is a Gaussian).
// It is based on the fact that the convolution of a Gaussian with a step
// function is the integral over the Gaussian, i.e., the error function erf(x).
// Note that the error function erf(x) is defined for a Gaussian exp(-x^2),
// while a Gaussian with a standard deviation of unity is given by exp(-x^2/2).
// Thus, one has to divide the result of the erf fit by a factor of sqrt2 to get
// the standard deviation sigma of the Gaussian.
//
// For a sum of two Gaussians (of equal amplitude), to have two separate maxima
// (and not blend into one maximum), their distance must be >2 sigma.
// Thus, 2 sigma may be taken as an approximation of the resolution.
// (This value is also close to the Rayleigh limit for Airy discs of similar
// size as the Gaussians. The Rayleigh criterion is not defined for Gaussians
// since Gaussians, in contrast to Airy discs, do not have a minimum.)
//
// For the human eye+brain it often seems that the resolution is actually
// somewhat better: The brain automatically does some deconvolution, so we
// can also recognize an image with two Gaussians with a distance somewhat
// below 2 sigma as resulting from two sources (although there are not two
// separate maxima any more).
// The macro assumes that the objects of the 'Blobs' image have sharp borders
// that are blurred only by the resolution of the optics, and that the border
// of the object it the line selection runs perpendicular to the line.
// The line width of 5 is chosen to reduce the noise.
requires("1.52f")
run("Blobs (25K)");
setLineWidth(5);
makeLine(211, 227, 232, 227); // line perpendicular to a sharp edge
y=getProfile();
x=Array.getSequence(y.length);
Fit.doFit("Error Function", x,y); // y = a+b*erf((x-c)/d)
Fit.plot;
if (Fit.rSquared < 0.95)
print("Warning: poor fit, uncertain result. R^2="+Fit.rSquared);
dParam = Fit.p(3); //fit parameter 'd'
sigma = dParam/sqrt(2);
resolution=2*sigma;
print("Point spread function: sigma = "+d2s(sigma,2)+" pixels");
print("Resolution limit = "+d2s(resolution,2)+" pixels");