Monday, July 13, 2009

ACTIVITY 6 - Properties of the 2D Fourier Transform

Activity 6.A Familiarization with FT of Different 2D Patterns

This activity is similar to the previous one. We created different patterns, took their respective Fourier Transform and then display their FT modulus. This patterns were done with the aid of MS Paint and Scilab made it possible to compute their Fourier Transform.

The first pattern is a square aperture and beside it is its FT. The FT of a square aperture is similar to the FT of a circle (from Activity 5.A) but instead of an "airy disk" pattern, it has a central square maximum and maxima and minima intensities along its four sides. This pattern is know as the sinc function.

The next pattern is an annular aperture. Again, this is similar to the pattern obtained from a circular aperture except that at some fringes are not present due to the circular blocking at the center.
Then, a square annular aperture. This was analogous to the circular and annular aperture. The pattern obtained was almost the same with the square aperture and as expected, there are fringes missing. This is similar to a diffraction from a slit only that it is along x and y.

Next is a two-slit pattern. Along its side is its FT which yields from the cancellation of the upper and lower sections from the aperture. Thus, we see a horizontal pattern which is due to the interference coming from two slits.

The last pattern is a two-dot pattern or a double circular aperture. As has been said, the result looks like that of a circular aperture with maxima and minima intensities.

Activity 6.B Anamorphic Property of the Fourier Transform

Using Scilab, a 2D sinusoid in the x direction was created, with frequency that can be varied.

nx = 100; ny = 100;
x = linspace(-1,1,nx);
y = linspace(-1,1,ny);
[X,Y] = ndgrid(x,y);
f = 4 //frequency
z = sin(2*%pi*f*X);
imshow(z,[]);


And then, just like in previous part its Fourier Transform was taken and its FT modulus was displayed.

The figure above shows the source signal and its FT in 2D and 3D, respectively. The first figure shows a sinusoidal signal with frequency equals to 5 and below it is a signal with frequency 10. They corresponding FT shows that the peaks moves farther apart as the frequency increases. This is expected because the FT of a sinusoidal signal will have its peak on its frequency. Thus if you increase it, the peaks will farther along its axis.

-------------------------------------------------------------------------------------------------

Then, a real image was simulated by adding a constant bias to the sinusoid created earlier. Its FT was also taken.

z = sin(2*%pi*20*X)+0.95
The figure above shows the sinusoidal signal with constant bias and its FT in 2D and 3D. On the contrary with the previous result, the FT of a sinusoidal signal with constant bias has 3 peaks. The two peaks with equal amplitude are from the frequency of the sinusoid and the peak at the center comes from the constant bias. We recall that the FT of a constant is a dirac delta. Also note that the sinusoid signal used has a frequency 20.
--------------------------------------------------------------------------------------------------

If we suppose to have an image of an interferogram in a Young's Double Slit experiment, we can use this routine to extract its actual frequencies. We just have to find where the FT peaked. If a non-constant bias was added like those very low frequency sinusoids, we just neglect the peaks at low frequencies.

--------------------------------------------------------------------------------------------------

The sinusoid was then rotated.

theta = 30;
z = sin(2*%pi*f*(Y*sin(theta) + X*cos(theta)));


As expected, the peaks of its FT are also rotated.

-------------------------------------------------------------------------------------------------

Then a combination of sinusoids along x and y are investigated.

z = sin(2*%pi*4*X).*sin(2*%pi*4*Y);


Again, we expect this results, a four-peak plot. The first two peaks are due to the sinusoid along x and the last two are due to the sinusoid along y.

-------------------------------------------------------------------------------------------------

Aside from a combination of sinusoids, we added rotated sinusoid to the signal.

z = sin(2*%pi*f*X).*cos(2*%pi*f*Y)+sin(2*%pi*f*(Y*sin(30) + X*cos(30)));

Before calculating and displaying its FT, I had a prediction that the rotated FT of the additional signal will add up to the four-peak FR of the two-sinusoid combination.

As I predicted, six peaks appear on the FT.

-------------------------------------------------------------------------------------------------

I will give myself a grade of 8/10 for finishing the blog late. Although I was able to finish the activity during class hours, I'm having hard time doing blog because I have no access to the internet at home. Thanks again to Raffy and Gilbert for always answering my questions and for always having a good discussions.

Monday, July 6, 2009

ACTIVITY 5 - Fourier Transform Model of Image Formation

Activity 5.A - Familiarization with Discrete FFT

The function fft2() was demonstrated by utilizing a 128x128 images. This FFT was used for an image of a 'circle' which was created with the given code:

x = [-1:0.01:1];
[X,Y] = meshgrid(x);
r = sqrt(X.^2 + Y.^2);
circle = zeros(size(X,1), size(X,2));
circle(find (r <=0.5)) = 1.0; imshow(circle,[]);


Because of being complex, the abs() function was used for the output of the FFT to be able to view the intensity values. The first image shows the output of fft2() which has interchanged quadrants along the diagonal.

image = imread('imagefilename');
gray = im2gray(image);
FT = fft2(gray);
a = abs(FT);


To compensate, the function fftshift() was used to bring the interchanged quadrants back which is show in the second figure. With a closer look, this image is consistent to the analytical fourier transform of a circle which is an airy disk pattern.

b = fftshift(abs(FT));

Then, fft2() was again applied and the intensity of the image was displayed.

c = abs(fft2(FT)
This procedure was also done for an image of letter 'A'.

The inversion of the image after FFT being applied was clearly seen from this sample than the image of a circle. Thus, applying FFT twice will make the image inverted.



Activity 6.B - Simulation of an Imaging Device

A 128X128 image of the letters "VIP" which represents an object and an image of a white circle against black background which represents the "aperture" of a circular lens were created using MS Paint. The images were opened as grayscale in Scilab, taking the 2D FFT of the "VIP" image using fft2() function and using fftshift() function for the aperture image (the aperture is already in Fourier space).

object = gray_imread('vip.jpg');
aperture = gray_imread('circle.jpg');
offt = fft2(object);
afftshift = fftshift(aperture);


The product their FFT was inversed to get the convolved image.

FOA = offt.*(afftshift);
IOA = fft2(FOA);
The images above show the results for different radii of white circles. From left to right, the radius of white circle decreases. Thus, for the largest circle resulted for the sharpest image and for the smallest circle resulted for the blurriest image. This shows that the system acts just like a microscope.



Activity 6.C Template Matching Using Correlation

A 128X128 image with the text "THE RAIN IN SPAIN STAYS MAINLY IN THE PLAIN" on it was created in MS Paint. Then another 128X128 image with letter "A" at the center was also created with the same type of font as the first image. These images were opened as grayscale in Scilab and then got their 2D FFT.

text = gray_imread('text.bmp');
A = gray_imread('atext.bmp');
TFFT = fft2(text);
TA = fft2(A);


Then the Fourier transform of the image with letter "A" was multiplied element per element to the conjugate of the Fourier transform of the text image. Lastly, the inverse FFT of the product was computed.

P = TA.*(conj(TFFT));
R = fftshift((fft2(P)));

The picture above shows the output image besides the original image. It noticeable that the output image shows hazy counterpart of each letter on the original image. There is a significant peaks (white dots) on the letter "A" counterpart. This template matching technique recognizes the pattern (letter "A") from the template (text image). However, this can be problematic when the pattern looks like any part of each letter on the text. Take for example, when we use letter "I" as our pattern the output image will shows peaks on parts of each letters with straight vertical line just like the pattern.



Activity 6.D Edge Detection Using the Convolution Integral

Using Scilab, a 3X3 matrix pattern of an edge which has a total sum of zero was created. Three variants of this pattern were created.

horizontal:
template = [-1 -1 -1; 2 2 2; -1 -1 -1];

vertical:
template = [-1 2 -1; -1 2 -1; -1 2 -1];

spot:
template = [-1 -1 -1; -1 8 -1; -1 -1 -1];

This pattern was convolve with a "VIP" image using the function imcorrcoef().

C = imcorrcoef(vip, template);

The image above shows the results. We can see there are sharp edges on the output image depending on the pattern used. The pattern was detected along the edge of the template.


It felt good that I finished this activity during class hours. However, I was not able to finish the blog right away because of lack of time. I will give myself 10/10 for this activity because I was able to produce expected output although I have published the blog late. Thanks again for the help of Gilbert and Raffy in discussing the activity.