DSPRelated.com
Code

2D Linear Convolution for 2D Image Processing

Senthilkumar February 27, 20112 comments Coded in Scilab

Performing 2D linear convolution involves taking 2D FFT of given 2D sequences and 2D IFFT of the resulting sequences.

function [y,X,H] = conv2d2(x,h)

[x1,x2] = size(x);
[h1,h2] = size(h);
X = zeros(x1+h1-1,x2+h2-1); 
H = zeros(x1+h1-1,x2+h2-1); 
Y = zeros(x1+h1-1,x2+h2-1); 
for i = 1:x1
    for j = 1:x2
        X(i,j)=x(i,j);
    end
end
for i =1:h1
    for j = 1:h2
        H(i,j)=h(i,j);
    end
end
disp(X,'x=') 
disp(H,'h =')
Y = fft2d(X).*fft2d(H);
disp(Y)
y = ifft2d(Y);
endfunction
/////////////////////////////////////////////
function [a2] = fft2d(a) 
//a = any real or complex 2D matrix
//a2 = 2D-DFT of 2D matrix  'a'
m=size(a,1)
n=size(a,2)
// fourier transform along the rows
for i=1:n
a1(:,i)=exp(-2*%i*%pi*(0:m-1)'.*.(0:m-1)/m)*a(:,i) 
end
// fourier transform along the columns
for j=1:m
a2temp=exp(-2*%i*%pi*(0:n-1)'.*.(0:n-1)/n)*(a1(j,:)).' 
a2(j,:)=a2temp.'
end
for i = 1:m
    for j = 1:n
        if((abs(real(a2(i,j)))<0.0001)&(abs(imag(a2(i,j)))<0.0001))
            a2(i,j)=0;
        elseif(abs(real(a2(i,j)))<0.0001)
            a2(i,j)= 0+%i*imag(a2(i,j));
        elseif(abs(imag(a2(i,j)))<0.0001)
            a2(i,j)= real(a2(i,j))+0;
        end
    end
end
/////////////////////////////////////////////////
function [a] =ifft2d(a2) 
//a2 = 2D-DFT of any real or complex 2D matrix
//a = 2D-IDFT of a2  
m=size(a2,1)
n=size(a2,2)
//Inverse Fourier transform along the rows
for i=1:n
a1(:,i)=exp(2*%i*%pi*(0:m-1)'.*.(0:m-1)/m)*a2(:,i) 
end
//Inverse fourier transform along the columns
for j=1:m
atemp=exp(2*%i*%pi*(0:n-1)'.*.(0:n-1)/n)*(a1(j,:)).' 
a(j,:)=atemp.'
end
a = a/(m*n)
a = real(a)
endfunction