%% load the image and convert it to double in [0,1]
f = im2double(imread('cameraman.tif'));
imshow(f);

%% take the fourier transform
F = fft2(f);

%% move the origin to the centre
F = fftshift(F);

%% show the spectrum (amplitude)
% the '[]' argument does automatic scaling of the dynamic range for display
figure,imshow(abs(F),[]);

%% show the log of the spectrum
figure,imshow(log(abs(F)),[]);

%% create a lowpass filter in frequency domain
H = fspecial('gaussian',256,20);
H = H/max(H(:));    % normalize so the maximum (DC component) is 1
figure,imshow(H);

%% perform filtering
Ft = F.*H;
figure,imshow(log(abs(Ft)),[]);

%% transform back into spatial domain
ft = real(ifft2(ifftshift(Ft)));
figure,imshow(ft);