0
$\begingroup$

I am trying to do some filtering with a gray scale image in the frequency domain. I am just getting into matlab after some time away from any signal processing or coding and can't quite seem to get my code to give me the proper output.

Questions:

  1. What is the correct data type to use for the centered image? All the odd elements inside the zero-padded array will be negative. Leaving it as a double plots this totally white for me. Making it an integer works better due to negative values becoming 0.

  2. When I use FFT2 on my centered image i do not see a periodic frequency spectrum like I was expecting. Why am I not seeing a larger repeating pattern? IF I zoom into the center of my DFT there is a tiny pattern that is period but I was expecting a larger pattern. Not sure if this has to do with my data type selection from question 1. But I am unsure I am am doing something in incorrect steps or using the wrong data type?

I need to complete the following steps:

  1. Given an input f(x,y) of size M x N obtain the padding sizes P and Q. P = 2M, Q = 2N.
  2. Form a padded image fp(x,y) of size P x Q using zero padding.
  3. Center the image using fp(x,y) * (-1)^(x+y)
  4. Compute the DFT F(u,v) of the image from step 3.
  5. Construct a real symmetric filter transfer function H(u,v) of size P x Q with center at (P/2, Q/2).
  6. Form the product G(u,v)=G(u,v)*H(u,v) using element wise multiplication.
  7. Obtain the filtered image by computing the IDFT of G(u,v), and decenter the image. gp(x,y) = (real[IDFT{G(u,v)}])*(-1)^(x+y)
  8. Obtain the final filtered result g(x,y) of the same size as the input image by extracting the M x N region from the top, left quadrant of gp(x,y).

I have the following code (edited with correct version):

%Imported an image of NxM size
%N = number of rows
%M = number of columns
%imagePart3 = imread('Fig2.19(a).jpg');
%imagePart3 = imread('Fig2.21(a).jpg');
imagePart3 = imread('Fig0431(d)(blown_ic_crop).tif');
f_xy_Image = double(imagePart3);

%zero image of size PxQ which is P=2*M, Q=2*M
[m,n]=size(f_xy_Image);
fp_xy_PaddedImage=zeros(2*m,2*n);

%creating looping variables to index the original image into the new zero
%padded image
[p,q]=size(fp_xy_PaddedImage);

%Looping through the zeropadded array at every row and column position to
%copy the values of the original image into the zero padded image.
%The loop compares each row level j from 1to q against the total number of
%rows in the originalimage. As long as the zero padded array is within the
%bounds of the original image the values of that image are copied into this
%new array. Otherwise once the bounds have been exceed zeros are placed in
%the remaining spots
for i = 1:p
    for j = 1:q
        if i <= m && j <= n
            fp_xy_PaddedImage(i,j) = f_xy_Image(i,j);
        else
            fp_xy_PaddedImage(i,j) = 0;
        end
    end
end

%Centering the image by creating a new variable to store the resulting
%center. Then multiply the padded image by (-1)^(x+y)
%In order to multiply every element of the array we need to loop through
%all (i,j) coordinates of the padded image. Therefore the loops cover from
%1 to the column size p, and 1 to the row size q.
centalizedImage = zeros(p,q);
for i = 1:p
    for j = 1:q
        centalizedImage(i,j) = fp_xy_PaddedImage(i,j).*((-1).^(i+j));
    end
end

centalizedImage1= uint8(centalizedImage);
%Computing the 2D DFT (FFT), F(u,v) on the image
F_uv_dftImage = fft2(centalizedImage);

%Implementing a centeredGaussian lowpass filter transferfunction using
%freqspace and meshgrid

[a,b] = freqspace(p,'meshgrid');
lpf_Mask = zeros(p,q);
for i = 1:p
    for j = 1:q
        lpf_Mask(i,j) = sqrt(a(i,j).^2 + b(i,j).^2);
    end
end
%Selecting the cut off frequency for the filter
%This will cut out any array elements below the selected frequency
%attenuating our image
H_uv = size(p,q);
for i = 1:p
    for j = 1:q
        if abs(lpf_Mask(i,j)) <= 0.1
            H_uv(i,j) = 1;
        else
            H_uv(i,j) = 0;
        end
    end
end

%Multiplying the the computed DFT F(u,v) with the low pass filter transfer
%function H(u,v)

G_uv = F_uv_dftImage.*H_uv;

%Computing the IDFT of product G

G_uv_realImage = real(ifft2(G_uv));

%Multiplying the IDFT by (-1)^(x+y) to complete the decentralize and obtain
%gp(x,y)

gp_xy_decentralizedImage = size(p,q);
for i = 1:p
    for j = 1:q
        gp_xy_decentralizedImage(i,j) = G_uv_realImage(i,j).*((-1).^(i+j));
    end
end

%Obtaining original image from the idftby copying out the eact emelent
%locations (i,j) from the top left of the idft image

g_ux = zeros(m,n);

for i = 1:m
    for j = 1:n
        g_ux(i,j) = gp_xy_decentralizedImage(i,j);
    end
end
g_ux1 = uint8(g_ux);

figure,subplot(241),imshow(imagePart3,[]);title('Original Image - f(x,y)');axis on;
subplot(242),imshow(fp_xy_PaddedImage,[]);title('Padded Image - fp(x,y)');axis on;
subplot(243), imshow(centalizedImage,[]);title('Centered Image - fp(x,y) * (-1)^(x+y)');axis on;
subplot(244), imshow(F_uv_dftImage,[]);title('DFT Image using Matlab FFT');axis on;
subplot(245), imshow(H_uv,[]);title('Transfer Function H(u,v) for a low pass filter');axis on;
subplot(246), imshow(G_uv,[]);title('Product G(u,v) = H(u,v)*F(u,v)');axis on;
subplot(247), imshow(gp_xy_decentralizedImage,[]);title('Filtered Image - Inverse DFT (IDFT)');axis on;
subplot(248), imshow(g_ux,[]);title('Final Filtered Result - g(x,y)');axis on;


``````````````````````````````````````````````````````````````````````````````
$\endgroup$
  • 1
    $\begingroup$ you seem to have forgotten to ask a precise signal processing question! $\endgroup$ – Marcus Müller Sep 27 '19 at 20:44
  • $\begingroup$ Sorry about that I'm not really sure how to narrow it down to a precise question. The problem is when this code runs i generate 8 pictures. The original image picture and the image with zero padding are fine (picutres1 & 2). I guess my first question is once I center the picture it will make every odd element in the array negative. If I keep this as a double matlab plots everything as white. If I make it an unsigned integer then all the negatives become 0 which I believe is the correct action. I wanted help to unerdstand which method would produce a better result. $\endgroup$ – Olek Sep 27 '19 at 21:15
  • $\begingroup$ I found out that for the DFT imageto display properly you have to apply a log transformation to it in the form of: s = c*log(1+r). Where r is your normalized DFT image (DFT_Image/255), and c is a constant. Then you plot s. This works great! $\endgroup$ – Olek Sep 28 '19 at 1:15
2
$\begingroup$

I figured it ouy. In order to display DFT image properly they need to be log scaled. I used the log transformation s = c*log(1+r) where r is your normalized input image, and c is a constant. Plot s and you will be able to see the DFT correctly.

%Log transform of DFT Image
normalized_dftImage_F_uv = F_uv_dftImage/255;
c1 = 1;
s_LogTransform_F_uv_dftImage = c1*log(1 + (normalized_dftImage_F_uv));
| improve this answer | |
$\endgroup$

Not the answer you're looking for? Browse other questions tagged or ask your own question.