Furrier Cats

Last week I gave an introductory lecture on imaging with radio interferometers. I started off with a Python demonstration of Fourier filtering using cat pictures, which seemed pretty popular and lead to one of my fellow lecturers coining the term “Furrier Cats” (this is genius and I wish I’d thought of it).

Since it was so well received, I thought I’d repeat it here. If you’re already au fait with Fourier transforms then this is perhaps not for you, but for all those new starters out there who (1) would like a visual introduction to Fourier analysis and (2) like cat pictures, read on.

Inspiration

This description was inspired by the ever insightful xkcd comic strip, and in particular this cartoon:

Which lead to an enjoyable afternoon for me, looking at cat pictures on google and selecting one that I thought would illustrate the points I wanted to make in an easily accessible way. In the end I actually used the first picture I found:

And this was what I read into my Python notebook:

%matplotlib inline

import numpy as np   # for array manipulation and the fft
import pylab as pl   # for plotting
import cv2           # for image file handling

cat = cv2.imread('./FIGURES/cat1.jpg',0)

pl.imshow(cat,cmap='gray')
pl.show()

Fourier (“Furrier”) Cats

The first thing I did was to take Fourier transform of the cat. Fourier transforms of functions are typically indicated by putting a tilde over the symbol denoting the function, i.e. becomes .

So obviously, we should label the Fourier transform of cat as cat_squiggle.

cat_squiggle = np.fft.fft2(cat)
cat_squiggle_shifted = np.fft.fftshift(cat_squiggle)
cat_spectrum = 20*np.log(np.abs(cat_squiggle_shifted))

There are a couple of things to notice here:

1. After transforming the cat, I have applied a shift to the output data. This shift puts the zero-point of Fourier space at the centre of our data array;
2. I have taken the logarithm of the magnitude of the shifted data using: 20*np.log(np.abs(data)). This converts the complex-valued output of the Fourier transform into decibels because I just want to visualise the amount of power at different frequencies.

Let’s see how it looks:

pl.subplot(121),pl.imshow(cat, cmap = 'gray')
pl.title('Cat'), pl.xticks([]), pl.yticks([])
pl.subplot(122),pl.imshow(cat_spectrum, cmap = 'gray')
pl.title('Fourier Cat'), pl.xticks([]), pl.yticks([])
pl.show()

The Fourier image is showing us how much power there is in the image on different spatial scales (Fourier frequencies). The key to understanding the pattern in Fourier space is that large spatial scales (i.e. big things in the image) correspond to small Fourier frequencies and that small spatial scales (i.e. small details in the image) correspond to large Fourier frequencies. i.e.

Big is Small & Small is Big

If we want to pick out information on particular scales, it can help to filter out the information on other scales that might confuse things.

Filtering Fourier Cats

Now we can set things up for filtering our Fourier Cat. We need to know (1) the dimensions of the image and (2) where the centre is.

rows, cols = cat.shape
crow,ccol = rows/2 , cols/2

To start with, let’s make a simple filter function that separates the inner most 40 x 40 pixels from everything else.

filter_fnc = np.zeros(cat_squiggle_shifted.shape)
filter_fnc[crow-20:crow+20, ccol-20:ccol+20] = 1.0

We can then use this filter to, firstly, mask out the inner most 40 x 40 pixels in Fourier space. This removes our small Fourier frequencies, i.e. removes the large scale information in our image.

cat_squiggle_hpf = np.copy(cat_squiggle_shifted)
cat_squiggle_hpf[np.where(filter_fnc==1.)] = 0.0+0*1j

…and, see how it looks.

pl.subplot(121),pl.imshow(20*np.log(np.abs(cat_squiggle_hpf)), cmap = 'gray')
pl.title('Filtered Fourier Cat'), pl.xticks([]), pl.yticks([])
pl.subplot(122),pl.imshow(cat_filtered_hpf)
pl.title('HPF Cat'), pl.xticks([]), pl.yticks([])
pl.show()

We can see that the large scale components of the image have gone and we’ve been left with small-scale features, which are often associated with the edges of a smooth area  (which are narrow changes in intensity and therefore small scale) and the details like whiskers and fur patterns.

Now let’s filter out the large Fourier frequencies:

cat_squiggle_lpf = np.copy(cat_squiggle_shifted)
cat_squiggle_lpf[np.where(filter_fnc==0.)] = 0.+0.*1j

and Fourier transform that back into image space:

cat_filtered = np.fft.ifftshift(cat_squiggle_lpf)
cat_filtered_lpf = np.fft.ifft2(cat_filtered)
cat_filtered_lpf = np.abs(cat_filtered_lpf)

pl.subplot(121),pl.imshow(20*np.log(np.abs(cat_squiggle_lpf)), cmap = 'gray')
pl.title('Filtered Fourier Cat'), pl.xticks([]), pl.yticks([])
pl.subplot(122),pl.imshow(cat_filtered_lpf)
pl.title('LPF Cat'), pl.xticks([]), pl.yticks([])
pl.show()

Here we can see that the large features in the image are retained, but we’ve lost the small scale detail such as the the whiskers and the pattern of the fur.

In both cases our output is a convolution. Mathematically this is because the filter we’ve applied in Fourier space is a logical (i.e. one or zero) multiplication and one of the most useful principles of Fourier analysis is that a multiplication in one space produces a convolution in the corresponding transform space. This is known as the Convolution Theorem.

Et voilà, Fourier – or Furrier – cats.

Now for the blog this.