Power spectrum of an image

Question

I have begun (a small project) to calculate the power spectrum of an image in the frequency domain.

So, what I have till now is the following:

%// close all; clear all; %// not generally appreciated
img   = imread('ajw_pic.jpg','jpg'); % it is a color image
img = rgb2gray(img); %// change to gray
psd = 10*log10(abs(fftshift(fft2(img))).^2 );
figure(2); clf
mesh(psd)

So far it looks good; I get the mesh plot which resembles the spectra I see in various academic papers.

However, what I am looking for is a graph plot of this power spectra versus the frequency, and I am not entirely sure how to get this frequency vector. I could do for instance:

N=400;        %// the image is 400 x 400
f=-N/2:N/2-1; %// possible frequencies?

but I am not convinced this is entirely correct as this gives rise to negative frequencies.

I'd really appreciate if someone could point me in the right direction to plot log frequency versus the power spectrum.

Answer 1

fft "splits" a signal into frequency "bins", where the maximum frequency you can observe is the Nyquist frequency or half your sampling frequency. This means that for:

Y = fft(X,N); %  (1D case)

the frequencies corresponding to fft-values in Y(1:N/2+1) will be:

f = [fs/2*linspace(0,1,N/2+1)]; % where fs is your sampling frequency

The other half of Y is just mirrored and comes from the intrinsics of the Fourier transform. If you don't want to understand it fully, I would say that there is no need to bother with it other than what you can find on wikipedia. But for the sake of curiosity you can check out a nice intuitive illustration of the origin of positive and negative frequencies: https://dsp.stackexchange.com/questions/431/what-is-the-physical-significance-of-negative-frequencies/449#449.

Some key differences for your 2D-image case:

1. Using fftshift you've moved the 0-frequencies to the center of the matrix, rather than having them at the edges as in the above 1D example. So you would actually get f = fs/2 * linspace(-1,1,N) (once again, don't mind the negative frequencies)

2. The next issue is to obtain the sampling frequency. Spatial frequencies are typically measured in [mm^-1], so in order to obtain it you actually need to know the physical distance between pixel centres (pixel pitch). But you could of course consider the spatial frequencies in [pixels^-1] in which case you are ready to go.

Answer 2

To plot the power spectra versus frequency of the image, one can use a process called 'radial averaging'. This calculates the average value of pixels that are a certain radial distance from the center of the image. Applying this to a power spectral density matrix results in a line plot of power versus frequency.

For more info and a a sample MATLAB code: http://www.mathworks.com/matlabcentral/fileexchange/46468-radialavg-zip/content/radialavg.m