Activity 8: Enhancement in the Frequency Domain

The Fourier Transform is so powerful that it allows enhancing images in the frequency domain. The convolution theorem can be used to filter and mask unwanted frequencies in the frequency domain. The following activities show the prowess of the Fourier Transform:

First I will showcase the two dots symmetrically along the x- axis and I take their Fourier Transform:

Figure 1

To generate the distance, I used the following code as a function of distance:

dist = 50
MAT = zeros(256,256);
MAT([128],[128 - dist, 128 +dist]) = 255,255;
Im = mat2gray(abs(fft2(MAT)));

As one can see, the frequency of the vertical lines increase as the distance between the 2 pixels increase. If I wanted horizontal lines, I would simply edit my code as follows:

dist = 50
MAT = zeros(256,256);
MAT([128 - dist, 128 + dist],[128,128]) = 255,255;
Im = mat2gray(abs(fft2(MAT)));

and get the following results:

Figure 2

I could make a checkerboard by changing MAT to:
MAT([128 - dist, 128 + dist],[128 + dist,128 + dist])

The next part of the activity shows an original (hand-drawn) image of two circles with a certain radius:

Figure 3

I take its FFT and get this:

Figure 4

The patterns look something from the Victorian era huh?

Now instead of circles, let's deal with squares, and see how the FFT vary in the frequency domain:

3 x 3 square (Left) and 5 x 5 square (right)

I manually input the squares locations at the center of the image. 

//Three by three Squares

N = 256

MAT2 = zeros(N,N);

MAT2([N/2 + 1,N/2 ,N/2 + 1],[N/2 - 11,N/2 - 10,N/2 - 9]) = [255,255,255;255,255,255;255,255,255];

MAT2([[N/2 + 1,N/2 ,N/2 + 1],[N/2 + 11,N/2 + 10,N/2 + 9]) = [255,255,255;255,255,255;255,255,255];

IM2 = mat2gray(abs(fft2(MAT2)));

imwrite(IM2, 'Squaredots.bmp');

imshow(IM2);

//Five by Five Squares

MAT3 = zeros(256,256);

MAT3([126,127,128,129,130],[116,117,118,119,120]) = 255*ones(5,5);

MAT3([126,127,128,129,130],[136,137,138,139,140]) = 255*ones(5,5);

IM3 = mat2gray(abs(fft2(MAT3)));

The code snippet above shows exactly how I manually input the squares. Note as one increases the size of the square, the greater the frequency of its FFT. On the next process, I will perform the FFT on ten randomly placed pixels on the image and convolve it with a square mask. The following shows the result:

Figure 6

As expected the results show the square mask displayed in the randomly generated area.

In the next part, 

//Convolution

A = zeros(200,200);

d = zeros(200,200);

Randomat = grand(10,2,"uin",0,200); //Generates the random locations of the numbers

for i = 1:1:10

     A(Randomat(i,1), Randomat(i,2) ) = 1;

end

five = rand(5,5);

d([98,99,100,101,102],[98,99,100,101,102])=five;

FTA = fft2(A);

Ftd = fft2(d);

FTAd = FTA.*Ftd;

IM4 = mat2gray(abs(fft2(FTAd)));

imwrite(IM4, 'RandConvol.bmp');

imshow(IM4);

Now we can do much more complicated filtering as in Figure 7. Here, we have a Lunar Orbiter Image with horizontal lines (left image):

Figure 7: Left image(Lunar Finding without Filter) Right image (Lunar Finding with Filter)

Figure 7: Left image(FFT of the image without the filter) Right image (FFT of the image with the filter)

The code snippet below shows how "elegantly," I generated the filter mask. Line 6 defines how much of the FFT image I remove with values ranging from 0 to 0.5. If I put 0.5, it removes everything, if I put 0, nothing happens. The best image is 0.49. The reason is that I apply the filter mask before the log, so the threshold is very small, as 0.48 and 0.49 are very distinct. The strength of this filter mask is that it can be applied on any image.

[1] //Lunar Landing Scanned Pictures: Line Removal

[2] Lunar = double(imread('C:\Users\Phil\Desktop\Academic Folder\Academic Folder 13-14 First Sem\AP 186\Activity 8\lun.jpg'));

[3] R = size(Lunar,1);

[4] C = size(Lunar,2);

[5] fftLunar = fftshift(fft2(Lunar));

[6] remov = 0.49;

[7] fftLunar(R/2 +1, 1:remov*C) = zeros(1,remov*C);

[8] fftLunar(R/2 +1,(1-remov)*C:C)= zeros(1,remov*C +1);

[9] fftLunar(1:remov*R,C/2 +1) = rand(remov*R,1);

[10] fftLunar((1-remov)*R:R,1+ C/2) = rand(remov*R+1,1);

[11] Im6 = mat2gray(abs(fft2(fftLunar)));

[12] imshow(Im6);

Note: Due to SPP and exams, I was unable to update the blog. The blog will continually be updated. Thank you!

Activity 7: Fourier Transform Model of Image Formation

Unlike the last activity in which I had such a hard time, this was more lenient, because Fourier Transforms are easy to deal with. Given a certain image, the Fourier Transform of the image in matrix form can be used to obtain its image in spatial frequency domain. The equation of the Fourier Transform can be divided into real and imaginary parts given by:

Figure 1

\

So  here Figure 1 shows an image of a circle with radius of 50% of total length of image or 64-bit radius. We take the Fourier transform and shift the image. Note we have to take the FFTshift or else the image will be just a dark figure. Figure 2 is the resulting image.

Figure 2

 Next we take the Fourier Transform of the letter "A." (Figure 3) The resulting Fourier Transform of Figure 3 is Figure 4. Take the Fourier Transform of Figure 4 again which is essentially taking the inverse Fourier Transform of figure 3 and we come back to the original image.

Figure 3

Figure 4

Figure 5

Next we perform convolution of two images. Convolution of two images is essentially the product of element per element of the two images in matrix form in inverse spatial frequency domain. Here the Fourier Transform is performed on the image VIP. It is then multiplied to a circle with an aperture of different radius. Since we are dealing with 128 bits, the circle is forced to have 128 x 128 elements. The elements are then multiplied to the Fourier transformed image. Figure 6 shows the resulting images.

Figure 6

The resulting image becomes unidentifiable as the aperture becomes smaller and smaller. Next we perform correlation. In correlation before multiplying two images element per element, the conjugate of one of the images in matrix form is taken in its frequency domain. Here Figure 7' s matrix is applied the function conj().

before element per element multiplication of the letter A. 

Figure 7

Figure 8

The result shows a sharp bright dot in which the letter A can be found. This is a great application in processing large image files and facial recognition. Its only downside is that the original image can no longer be used.

Edge Detection:

Another important application of Fourier Transform, is edge detection. By performing a convolution of 3 x 3 unit matrix and  an image, the edges of the image can be seen. At points where there is no change, the matrix gives a value of zero. When the matrix does change, the value is no longer zero. The following matrices are convolved against the same VIP image. The more variation there are in the matrix, the better the edge detection.

Figure 9

Here, I try to recreate the Sobel edge detection. Since most edge detection split their matrix operator into a horizontal and vertical component, I do exactly the same method.  The Gradient magnitude is then computed by [3]:

where Gx and Gy are the matrix operators shown in Figure 10. As you can see the resulting image gradient: the Sobel Gradient shows a very sharp line on the edge of the letters VIP!

Figure 10

I would give myself a 12 for this activity!!!

Code Snippet:

//Set-up

x = [-1: 0.0157480:1];

[X,Y] = meshgrid(x);

r = sqrt(X.^2 + Y.^2);

circle = zeros(size(X,1), size(X,2));

circle(find (r <=0.05)) = 1.0;

//imshow(circle);

A = double(rgb2gray(imread('C:\Users\Phil\Desktop\A.png')));

//imshow(A);

VIP = double(rgb2gray(imread('C:\Users\Phil\Desktop\VIP.png')));

//imshow(VIP);

Spain = double(rgb2gray(imread('C:\Users\Phil\Desktop\CorrelationC.png')));

//imshow(Spain)

//FFT Circle

r = fftshift(fft2(circle));

//imwrite(mat2gray(abs(r)), 'FFTCircle.bmp');

//FFT twice

s = fft2(fft2(circle));

//imwrite(mat2gray(abs(s)), 'FFTtwice.bmp');

//FFT Letter A

z = fftshift(fft2(A));

//imwrite(mat2gray(log(abs(z))),'FFTA.bmp');

t = fft2(fft2(A));

//imwrite(mat2gray(log(abs(t))),'InverseFFT.bmp');

//Convolution

FFTVIP = fft2(VIP); 

FFTCircle = fftshift(circle);

convolved = fft2(FFTVIP.*FFTCircle);

G = mat2gray(log(abs(convolved)));

//imwrite(G,'FFTConvolved0.05.bmp');

//Correlation

FFTA = fft2(A);

FFTSpain = fft2(Spain);

F = mat2gray(log(abs(fftshift(fft2(FFTA.*conj(FFTSpain))))));

//imwrite(F,"FFTCorrelation.bmp");

//Edge Detection

MAT = zeros(128,128);

MAT2 = zeros(128,128);

randomat = [-1,0,3;5,-2,-4;2,-1,-2];

submat = [-1,-1,-1;2,2,2;-1,-1,-1];

submat2 = [-2,1,-2;1,4,1;-2,1,-2];

submat3 = [0,1,0;1,-4,1;0,1,0];

submat4 = [3,-1,3;-1,-8,-1;3,-1,3];

submat5 = [-4,7,2;-8,-1,3;2,-2,1];

sobelx = [-1,0,1;-2,0,2;-1,0,1];

submat7 = [-0.25,0,0.25;0,0,0;0.25,0,-0.25];

sobely = [1,2,1;0,0,0;-1,-2,-1]

MAT([63 64 65], [63 64 65]) = sobely;

MAT2([63 64 65], [63 64 65]) = sobelx;

FFTVIP = fft2(VIP); 

FFTMAT = fft2(MAT);

FFTMAT2 = fft2(MAT2);

convolvedmatrx = fftshift(fft2(FFTVIP.*FFTMAT));

convolvedmatrx2 =fftshift(fft2(FFTVIP.*FFTMAT2));

H = mat2gray(abs(convolvedmatrx));

I = mat2gray(abs(convolvedmatrx2));

J = sqrt((H.*H) + (I.*I))

imwrite(J,'sobel.bmp');

imwrite(edge(VIP,'sobel'),'sobel2.bmp')

References:

1. Soriano, Jing Activity 7 Activity 7: Fourier Transform Model of Image Formation

2. http://www.csce.uark.edu/~jgauch/5683/notes/ch04a.pdf

http://en.wikipedia.org/wiki/Sobel_operator

Activity 6 Enhancement by Histogram Manipulation

First and foremost, I really, really had a hard time with this activity. I was messing up with the program left and right. Why? There were two parts that I had trouble with:

(Ignore my ranting here and skip to the introduction if you don't want to read this part)

1. The back-projecting of the pixels between the desired CDF and original CDF
2. Faster way of switching pixels in the image

Don't worry. This blog contains all the cheat codes so you won't spend the midnight oil like I did, so if you want some challenge DO NOT READ my blog... I'm just joking.

Introduction:

Images believe it or not have a distribution function. The pixels, specifically, make up the probability distribution function in which a histogram of the image can be created. Pixel colors range from 0 to 255. 0 being totally black and 255 being sharpest white. Between these two values, 254 shades of gray (no pun intended) can be technically observed.

Now let us talk about Histogram Manipulation!
Histogram manipulation is one way of improving the quality of an image, enhancing certain image features, or mimicking the responds of different systems such as the human eye.[1]

To this extent, one can take the probability density function (PDF) of an image by simply counting the number of times each pixel number ranging from 0-255 appeared in the image and plotting No. of pixels vs pixel (0-255). Then one can take the cumulative density function (CDF), which is essentially similar to the probability density function except for two things:

1. It's a cumulative, which means the pixel numbers are increasing.
2. The pixel number is normalized by dividing the y-axis value by the total number of pixel. This means the y-axis always ranges from 0 to 1.

This process can easily be done in scilab, since scilab allows easy conversion of images to a matrix of pixel numbers.
Another process which will be used is the back-projecting of a desired CDF to the orignal CDF of the image. It's similar to reshuffling the image pixels of your choosing by simply deciding what CDF to use.
This back-projection can be understood easily thanks to Mum Jing's figure shown below:

Figure 1 [1]

Now let's get started!

Since I recently watched Monsters University in theaters, I'm going to take the image of Michael "Mike" Wazowski, best friend of Sully. 

 Figure 2 [2]

 Figure 2

 

Figure 3: PDF (left) and CDF (right)

Note that the x-axis is Pixel from 0-255 and y-axis for the PDF is the occurrence of the specific pixel in the image and the y-axis for the CDF ranges from 0 to 1.

Now I'm going to make a CDF that is a straight line like this:

Figure 4

To do this we make use of the code that I tirelessly worked on:

(1)Image = imread('C:\Users\Phil\Desktop\MU.jpg');

(2)im = RGB2Gray(Image);

(3)imwrite(im,'MU.bmp');

(4)pixel = 255

(5) //CODE SNIPPET from: Mum Jing's blogspot of obtaining pixel color (0-255) from an image

(6)val=[];

(7)num=[];

(8)counter=1;

(9)for i=0:1:pixel 

(10)    [x,y]=find(im==i); //finds where im==i 

(11)   val(counter)=i; num(counter)=length(x); //find how many pixels of im have value of i 

(12)    counter=counter+1;

(13)end

(14)//End of Code Snippet

(15)[n,m] = size(im);

(16)nor = num/(n*m);

(17)mat1 = val';

(18)mat2 = nor';

(19)//plot(val, nor);

(20)mat3 = cumsum(mat2);//creates matrix containing both 0-255 and corresponding no.of pixels

(21)oldp = linspace(0,1,pixel+1); // creates a straight line y-axis value (desired CDF)

(22)//input1 = linspace(-10,10,pixel+1)

(23)//oldp = (1+exp(-1*input1))^-1

(24)newp = interp1(oldp, mat1, mat3,'nearest'); //

(25)newm = [mat1;mat3;mat1;oldp;newp];//0-255;0-1(currentCDF);0-255;0-1 (new CDF);0-255

(26)editim = im;

(27)// This part of the code backprojects/replaces the old pixels of the image to the new one.

(28)for i = 1:pixel  

(29)    editim(find(im==newm(1,i)))= newm(5,i); //

(30)end

(31)imwrite(editim,'EditedMU.bmp');

(32)[edity, editx] = imhist(editim);

(33)clf()

(34)plot(editx,edity/(n*m));

(35)sumedit = cumsum(edity/(n*m));

(36)clf();

(37)plot(editx,sumedit);

Tada!! Here is the resulting PDF:

Figure 5

And Figure 6 shows the resulting image. As you can see, there is a distinct difference between this image and the original. The contrast between white and black is greater. This can be attributed to the fact that linearity creates pixels counts which drop to zero (As seen from the PDF).

Therefore there is a loss of information...

Now let's try to create a sigmoid function, defined by:

Here is my poor attempt at trying to replicate the PDF curve of the sigmoid function using GIMP:

(Will still update. Discuss Code and Gimp.....)

References:

1. Soriano, Jing, Activity 6 Enhancement by Histogram Manipulation

2. http://blogcritics.org/movie-review-monsters-university-3/