Understanding the output of fftfreq function and the fft plot for a single row in an image

Question

I am trying to understand the function fftfreq and the resulting plot generated by adding real and imaginary components for one row in the image. Here is what I did:

import numpy as np
import cv2
import matplotlib.pyplot as plt

image = cv2.imread("images/construction_150_200_background.png", 0)
image_fft = np.fft.fft(image)
real = image_fft.real
imag = image_fft.imag

real_row_bw = image_fft[np.ceil(image.shape[0]/2).astype(np.int),0:image.shape[1]]
imag_row_bw = image_fft[np.ceil(image.shape[0]/2).astype(np.int),0:image.shape[1]]

sum = real_row_bw + imag_row_bw

plt.plot(np.fft.fftfreq(image.shape[1]), sum)
plt.show()

Here is image of the plot generated :

I read the image from the disk, calculate the Fourier transform and extract the real and imaginary parts. Then I sum the sine and cosine components and plot using the pyplot library.

Could someone please help me understand the fftfreq function? Also what does the peak represent in the plot for the following image:

I understand that Fourier transform maps the image from spatial domain to the frequency domain but I cannot make much sense from the graph.

Note: I am unable to upload the images directly here, as at the moment of asking the question, I am getting an upload error.

Answer 1

I don't think that you really need fftfreq to look for frequency-domain information in images, but I'll try to explain it anyway.

fftfreq is used to calculate the frequencies that correspond to each bin in an FFT that you calculate. You are using fftfreq to define the x coordinates on your graph.

fftfreq has two arguments: one mandatory, one optional. The mandatory first argument is an integer, the window length you used to calculate an FFT. You will have the same number of frequency bins in the FFT as you had samples in the window. The optional second argument is the time period per window. If you don't specify it, the default is a period of 1. I don't know whether a sample rate is a meaningful quantity for an image, so I can understand you not specifying it. Maybe you want to give the period in pixels? It's up to you.

Your FFT's frequency bins start at the negative Nyquist frequency, which is half the sample rate (default = -0.5), or a little higher; and it ends at the positive Nyquist frequency (+0.5), or a little lower.

The fftfreq function returns the frequencies in a funny order though. The zero frequency is always the zeroth element. The frequencies count up to the maximum positive frequency, and then flip to the maximum negative frequency and count upwards towards zero. The reason for this strange ordering is that if you're doing FFT's with real-valued data (you are, image pixels do not have complex values), the negative frequency data is exactly equal to the corresponding positive frequency data and is redundant. This ordering makes it easy to throw the negative frequencies away: just take the first half of the array. Since you aren't doing that, you're plotting the negative frequencies too. If you should choose to ignore the second half of the array, the negative frequencies will be removed.

As for the strong spike that you see at the zero frequency in your image, this is probably because your image data is RGB values which range from 0 to 255. There's a huge "DC offset" in your data. It looks like you're using Matplotlib. If you are plotting in an interactive window, you can use the zoom rectangle to look at that horizontal line. If you push the DC offset off scale, setting the Y axis scale to perhaps ±500, I bet you will start to see that the horizontal line isn't exactly horizontal after all.

Once you know which bin contains your DC offset, if you don't want to see it, you can just assign the value of the fft in that bin to zero. Then the graph will scale automatically.

By the way, these two lines of code perform identical calculations, so you aren't actually taking the sine and cosine components like your text says:

real_row_bw = image_fft[np.ceil(image.shape[0]/2).astype(np.int),0:image.shape[1]]
imag_row_bw = image_fft[np.ceil(image.shape[0]/2).astype(np.int),0:image.shape[1]]

And one last thing: to sum the sine and cosine components properly (once you have them), since they're at right angles, you need to use a vector sum rather than a scalar sum. Look at the function numpy.linalg.norm.