Below is some code for filtering out certain spatial frequencies from a grayscale image using python. The process is rather straightforward in theory. First, perform the two-dimensional discrete Fourier transform (using numpy's fft2 function) on the data from the original image. This produces a new representation of the image in which each pixel represents a spatial frequency and orientation, rather than an xy coordinate. Second, design a filter. The filter should be an array with the same shape and size as the image, with values between 0 and 1. The filter describes by how much to attenuate each frequency/orientation component; since a pixel's radial distance determines its frequency, a filter with a filled circle of 1s at its centre and 0s everywhere else will produce an image containing all frequencies below a certain value, and no frequencies above that value, irrespective of orientation. Third, multiply the Fourier image by the filter, and perform the inverse Fourier transform on the result. This produces an array of complex numbers whose absolute values correspond to the filtered image.
There are a couple of things to watch out for, however. When Fourier-transformed images are plotted graphically, the low frequencies are found at the centre; this is not what fft2 actually produces, so we need to also apply numpy's fftshift—and later, ifftshift—functions. Another caveat when plotting Fourier images is that, quite often, the zero-frequency component (i.e., the DC-offset of the original image) is enormous compared to all other components. Therefore, Fourier images (such as those found in other online tutorials) are sometimes log-transformed prior to plotting.
A third issue is ringing, wherein sharp transitions within the filter causes distortions in the final image. The simple on/off (or 'ideal') filter described in first paragraph causes a lot of ringing, so it is usually better to use a filter design with smooth transitions, such as the Butterworth filter. A final issue is that spatial filtering may significantly change the image's overall luminance pattern, so some form of equalisation is often necessary.
The image below shows lowpass, highpass, and bandpass Butterworth filters applied to a grayscale image. The code requires numpy, scipy, and matplotlib. It also uses the scikits image-processing toolbox for histogram equalisation. I define my frequencies in degrees of visual angle, assuming I'm sitting 57 cm from my monitor.
Banner image is The Little Owl by Albrecht Durer.
Below are some links on image processing I found useful:
- This site contains a quite thorough guide to image processing (including filtering).
- A complete description (along with some C code) of the Discrete Fourier transform can be found here.
- An intuitive explanation of Fourier analysis can be found here.
- I also found ImageMagick's documentation quite informative.