Like in everyday photography, it is crucial in (electronic) microscopy to get images that are not blurred. The focus, is the process of adapting the optic conditions to get clear images. On an electron microscope, the focus mainly depends on the position of the lens on the z axis, the working distance (distance from electron gun to specimen) and the magnification. Out of those 3 dimensions, we will typically adjust the working distance to focus the image.

Some microscope interfaces like Zeiss’ already provide an autofocus, but it is typically slow (more than 15 seconds), yields bad result with images not focused or even worse, and they have to be triggered manually. To automate specimen analysis, we cannot afford to only focus at the beginning of the analysis, because samples can be very irregular in height, so the that the distance between the lens and the specimen cannot be constant, leading to blurred samples. The analysis algorithm should detect when the image is out of focus and automatically trigger the focus process.

A good autofocus should be fast, robust to noise and other factors, and produce sharp images. As most traditional methods to autofocus an image are iterative, the speed of an autofocus is mostly dictated by the number of tries it needs, rather than on complexity or other factors. New machine learning based methods are able to make much better initial guesses on the ideal working distance but here we will focus on a basic optimization based method.

Prerequisites/setup

As we are concerned with electronic microscopes, we need to have access to the microscope API to implement an autofocus. Depending on the model, this may be a TCP API or a SDK (Zeiss provides them in C/C#, Thermofischer provides it in Python). The API should provide a function to set (and get) the working distance, and a function to take images. (We will assume that the magnification and position of the stage are fixed).

There are ultimately 2 techniques to autofocus. Techniques that rely on the statistical variance of pixel values or spatial frequency content of the image. We will focus on the first one.

Detection

In order to trigger our autofocus, we first need to be able to tell when an image is out of focus. Intuitively, the more an image is out focus or blurred, the less it will have sharp edges.

A known tool in computer vision to detect edges is a Laplacian filter, a variant of convolutions, one of the most important operators in deep learning. To understand this, we can reason that a sharp image will have a high variance in pixel intensities, whereas blurring an image will by definition bring the pixel values to the mean.

To detect edges, we should thus aggregate the Laplacian filter in such a way to reflect the number of edges in the picture. One method to achieve this is to use the variance of the Laplacian. Mathematically it is computed as:

$$ \text{{Var}} = \frac{1}{{\text{{ksize}}^2}} \sum_{i=0}^{{\text{{ksize}}-1}} \sum_{j=0}^{{\text{{ksize}}-1}} ( \text{{Lap}}(\text{{src}})(i,j) - \text{{mean}} )^2 $$

Where the Laplacian definition is:

$$ \text{{Lap}}(f(x, y)) = \frac{{\partial^2 f(x, y)}}{{\partial x^2}} + \frac{{\partial^2 f(x, y)}}{{\partial y^2}} $$

This is the theoretical case in a continuous domain. For image processing, the Laplacian is approximated with the Laplacian filter mask, which is a convolution operator of the following shape:

$$ \begin{bmatrix} 0 & 1 & 0 \ 1 & -4 & 1 \ 0 & 1 & 0 \end{bmatrix} $$

This kernel is then simply convolved with the image to obtain the second derivative. The values in the kernel reflect the weights that the neighboring pixels (surrounding the center pixel) get in this operation.

The variance of the Laplacian is pretty much constant across images from the same distribution, so we can find the variance for images we know are not blurred, and use that value to detect out of focus pictures in the future. Practically, this can either be done with manual intervention, or by systematically autofocusing the first image to get a baseline value.

detect() {
	if laplacian(var(x)) > threshold:
		focus()
}

Optimization

To find the optimal working distance, we can just treat the microscope as a discrete (noisy) function. As computing the variance of the Laplacian would be pretty hard, we will simply use Brent’s method, which does not require computing derivative, but just evaluating the function at various points. We will use the minimize function from the optimize scipy module, so we will need to minimize $1 / Var$.

def autofocus(self, max_iter: int = 6) -> float:
	# assume self points to a class of a Microscope API
    def f(z: float) -> np.ndarray:      
				# get image from microscope and evaluate Laplacian variance
        img = self.scan_image(z)     
        return 1/lap_var(img)
    wd = optimize.minimize_scalar(
        f, 
        method='bounded', 
        bounds = [10, 100], # typically it is possible to give good bounds for the working distance
        options= {'maxiter':max_iter} # this will determine the performance vs. quality tradeoff of the autofocus
    )
    self.set_wd(wd.x)    
    return wd

Practically we may need to pass the image through a high pass filter as seen in the literature, but here we will assume this is not needed.

Simulation

Let’s test our “algorithm” now by simulating an electron microscope. Let’s assume a fixed position and image. We will need to model the amount of blur based on the working distance. This is done by mapping a working distance to the radius of disk shaped filter.

import cv2
import numpy as np
from scipy import optimize

def create_disk_kernel(radius):
    """Create a disk-shaped kernel of a given radius."""
    kernel_size = 2*radius + 1
    kernel = np.zeros((kernel_size, kernel_size))

    y, x = np.ogrid[-radius:radius+1, -radius:radius+1]
    mask = x**2 + y**2 <= radius**2
    kernel[mask] = 1

    kernel /= np.sum(kernel)
    
    return kernel

def lap_var(image):
    return cv2.Laplacian(image, cv2.CV_64F).var()

class SimMicroscope:
    def __init__(self, url="399471.png"):
        self.image = np.array(Image.open(url))
        # Optimal working distance
        self.wd_optimal = 50
        self.wd=self.wd_optimal
        # Proportionality constant that relates PSF radius to working distance difference
        self.blur_scale = 0.1 

    def set_wd(self, wd):
        self.wd=wd

    def scan_image(self, wd=None) -> np.ndarray:
        if wd is None: 
            wd=self.wd
        else:
            self.wd=wd
        # Calculate the difference in working distance from the optimal value
        d_diff = np.abs(wd - self.wd_optimal)

        # Calculate the radius of the disk-shaped PSF
        radius = int(round(d_diff * self.blur_scale))
        psf = create_disk_kernel(radius)
        
        # Convolve the image with the PSF
        blurred = cv2.filter2D(self.image, -1, psf)
        
        return blurred

Let’s see if this is a realistic model. We would expect the variance of the Laplacian to be a bell shaped curve around wd_optimal . By plotting it with matplotlib, we get the following distribution, which is pretty good for now:

X = np.linspace(10, 100, 10)
y = [cv2.Laplacian(microscope.scan_image(wd=x), cv2.CV_64F).var() for x in X]
plt.plot(X, y)

Then, we can simply add our autofocus function by extending this class with the method.

class AutoSimMicroscope(SimMicroscope):
    def autofocus(self, max_iter: int = 6) -> float:
            def lap_v_at_wd(z: float) -> np.ndarray:      
                img = self.scan_image(z)     
                return 1/lap_var(img) # we want to maximise lap_var so we take 1/x
            wd = optimize.minimize_scalar(
                lap_v_at_wd, 
                method='bounded', 
                bounds = [10, 100],
                options= {'maxiter':max_iter}
            )
            self.set_wd(wd.x)    
            return wd

m2 = AutoSimMicroscope()
m2.set_wd(10) # put the microscope out of focus on purpose
Image.fromarray(m2.scan_image()) # => the resulting image is blurred
m2.autofocus(max_iter=5) # => wd is now ~50 as expected
Image.fromarray(m2.scan_image()) # => the resulting image is in focus

There you have it, we have an autofocus that already works better than most algorithms that come with the microscope (such as Zeiss’s). Depending on the image settings and manufacturer, we can run this in ~5 seconds vs. 15 seconds with a much better quality assurance. If you’re interested in learning more, either about the context or about new machine learning based methods, a good paper to read is Deepfocus.