Laboratory 4

The aim of this laboratory session is to introduce you to MATLAB's Fourier transform routines, learn how to interpret Fourier transform data, and then to see how convolution is implemented in the frequency domain using simple averaging filters as an example.

Building up an Image from its Fourier Components


Original Image
im = imresize(im,0.5)


Partial reconstruction of your image using freqcomp.m,
Frequency components used 50

Partial reconstruction of your image using freqcomp.m
Frequency components used 200
Partial reconstruction of your image using freqcomp.m
Frequency components used 100

Gaussian Filter


Original Image im


imfft = fft2(im);

gafilt = fspecial('gaussian', size(im), 6)
gafiltfft = fft2(gafilt);


real(ifft2(newimgafft))
newimgafft = gafiltfft.*imfft;

Average Filter


Original Image im
imfft = fft2(im);

avfilt = averagefilter(size(im),21,21);
avfiltfft = fft2(avfilt);

real(ifft2(newimavfft))
newimavfft = avfiltfft.*imfft;

Listing of the code I used to obtain the above images: 

im = imread('d:\comp vision\images\travis.jpg');
im = rgb2gray(im);

avfilt = averagefilter(480,640,21,21);       % average filter
gafilt = fspecial('gaussian', [480 640], 6); % gaussian filter

imfft = fft2(im);
gafiltfft = fft2(gafilt);
avfiltfft = fft2(avfilt);
newimavfft = avfiltfft.*imfft;
newimgafft = gafiltfft.*imfft;

imwritesc(im, 'im.png');
imwritesc(gafilt, 'gafilt.png');
imwritesc(avfilt, 'avfilt.png');

imwritesc(fftshift(log(abs(imfft)+eps)), 'imfft.png');

imwritesc(fftshift(log(abs(gafiltfft)+eps)), 'gafiltfft.png');
imwritesc(fftshift(log(abs(avfiltfft)+eps)), 'avfiltfft.png');

imwritesc(fftshift(log(abs(newimavfft)+eps)), 'newimavfft.png');
imwritesc(fftshift(log(abs(newimgafft)+eps)), 'newimgafft.png');

imwritesc(fftshift(real(ifft2(newimavfft))), 'newimav.png');
imwritesc(fftshift(real(ifft2(newimgafft))), 'newimga.png');

 

Phase-only Images and Phase and Amplitude Swapping


im1
im2

im1 magnitude and im2 phase
im1 phase and im2 magnitude

im1 phase only
im2 phase only

Listing of the code I used to obtain the above images: 

im1 = imread('d:\comp vision\images\travis.jpg');
im2 = imread('d:\comp vision\images\travis1.jpg');

im1 = rgb2gray(im1); % convert to greyscale
im2 = rgb2gray(im2); % convert to greyscale

im1fft = fft2(im1);
im2fft = fft2(im2);

im1phase = angle(im1fft);
im2phase = angle(im2fft);

im1mag = abs(im1fft);
im2mag = abs(im2fft);

newimfft1 = im1mag.*(cos(im2phase) + i*sin(im2phase));
newimfft2 = im2mag.*(cos(im1phase) + i*sin(im1phase));

imwritesc(real(ifft2(newimfft1)), 'im1mag_im2phase.png');
imwritesc(real(ifft2(newimfft2)), 'im1phase_im2mag.png');

im1phaseonly = im1fft./im1mag;
im2phaseonly = im2fft./im2mag;

imwritesc(real(ifft2(im1phaseonly)), 'im1phaseonly.png');
imwritesc(real(ifft2(im2phaseonly)), 'im2phaseonly.png');

imwritesc(im1, 'im1.png');
imwritesc(im2, 'im2.png');