Custom Image Spatial Filtering Code using loops not giving the same result as using the built in imtransform() function
    6 views (last 30 days)
  
       Show older comments
    
    Joshua Kearney
 on 12 Jun 2018
  
    
    
    
    
    Edited: Anton Semechko
      
 on 12 Jun 2018
            Hi, I'm working on trying to create a custom code to apply spatial filtering without Matlab functions for school. So I created a custom convolution function to be applied to an image and a kernel but the resultant image looks different for both of these images and I'm hitting a wall with why. My custom code is more blurred and I think my convolution function is incorrect but to me, it looks like I'm applying the equation correctly.
Any help would be appreciated. Thanks.
 clear all; close all; clc;
 % Read images 1 and 2
 im = imread('img1.png');
 A = rgb2gray(im);
    % Read Kernels
Kernel_1 = (1/9)*ones(3);
Kernel_2 = (1/49)*ones(7);
Kernel = Kernel_2;          % Set current kernel
 img_out = convolution(A,Kernel);   %Perform convolution on image and selected kernel
 img_out_filter = imfilter(A,Kernel,'same','conv');
 %%Display output images
subplot(1,3,1); imshow(A); title('Original') 
subplot(1,3,2); imshow(img_out,[]); title('Custom Created Function')
subplot(1,3,3); imshow(img_out_filter); title('Matlab imfilter Function')
    %%Rotate input matrix
function rot_mat = rot(mat,theta)
    theta = -25 *2*pi/360;
    R = [cos(theta)  sin(theta) 0;
         -sin(theta) cos(theta) 0;
         0           0          1];
    rot_mat = mat*R;
end
    %%Perform convolution on image and kernel
function B = convolution(A, k)
    [ky, kx] = size(k);             % Read kernel size
    im_pad = padarray(A, [kx ky]);  % Pad original image
    [y, x] = size(im_pad);          % Read image size
    B = zeros(x,y);         % Create empty matrix to store output image
    kr = rot90(k);          % Rotates kernel 180 deg for convolution
    kr = rot90(kr);
    for i=(1+ky):(y-ky)         % index through each image row
        for j=(1+kx):(x-kx)     % index through each image pixel        
            neigh = im_pad(i-floor(ky/2):i+floor(ky/2), j-floor(kx/2):j+floor(kx/2));   % Create local neighborhood of image
            accumulator = 0; 
            for u=1:ky      % index through each kernel row 
                for v=1:kx  % index through each kernel element
                    if(i>ky && i<y-ky && j>kx && j<y-kx)
                    temp = neigh(u,v)*kr(u,v);
                    accumulator = accumulator + temp;
                    end
                end
            end
            B(i,j) = accumulator;   %Set value of pixel in new image with convolution operation resultant
        end
    end
    B=B(1+ky:y-ky,1+kx:x-kx);       % Remove image padding
end
0 Comments
Accepted Answer
  Anton Semechko
      
 on 12 Jun 2018
        
      Edited: Anton Semechko
      
 on 12 Jun 2018
  
      Actually, both your 'convolution' function and built-in 'imfilter' function produce very similar results. The discrepancy you were observing was due to the fact that you forgot to cast A (the image) into double format before filtering. Here is the corrected part of your code.
 % Read images 1 and 2
 im = imread('img1.png');
 A = rgb2gray(im);
 A = double(A);
    % Read Kernels
Kernel_1 = (1/9)*ones(3);
Kernel_2 = (1/49)*ones(7);
Kernel = Kernel_2;          % Set current kernel
 img_out = convolution(A,Kernel);   %Perform convolution on image and selected kernel
 img_out_filter = imfilter(A,Kernel,'same','conv');
 %%Display output images
figure('color','w')
subplot(1,3,1); imshow(A,[]); title('Original') 
subplot(1,3,2); imshow(img_out,[]); title('Custom Created Function')
subplot(1,3,3); imshow(img_out_filter,[]); title('Matlab imfilter Function')
3 Comments
  Anton Semechko
      
 on 12 Jun 2018
				Glad I could help. You did a good job with your implementation of 2D convolution.
  Anton Semechko
      
 on 12 Jun 2018
				
      Edited: Anton Semechko
      
 on 12 Jun 2018
  
			Also note that there is not need to reflect the filtering kernel if its radially symmetric, like it is in your case.
More Answers (0)
See Also
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!