Problems implementing my own 2d DCT function. Smearing effect

Hi, I'm trying to implement my own Discrete Cosine transform and inverse transform functions. I have all the code written however the output images from the fuction seem to have a blurring effect.

EX:

Here is my code :

 function dct = my_dct_2d(image)
%% Calculate the cosine coefficients

image = double(image);

N = length(image);

dct = zeros(N);

for u=0: N-1
    for v=0: N-1
        if u==0
           cu = 1/sqrt(N); 
        else
           cu = sqrt(2/N);
        end
        if v==0
           cv = 1/sqrt(N); 
        else
           cv = sqrt(2/N);
        end
        
        mysum = 0;
        for x=0:N-1
            for y=0:N-1
                temp = image(x+1,y+1)*cos(pi*u*(2*x+1)/(2*N))*cos(pi*v*(2*y+1)/(2*N));
                
                %temp = image(k,l)*cos(((2*(k-1)+1)*i*pi)/(2*n))*cos(((2*(l-1)+1)*j*pi)/(2*n));
                mysum = mysum+ temp;
            end
        end
        dct(u+1,v+1) = cu*cv*mysum;
    end
end

end
 
function idct = my_idct_2d(image)
%% Calculate the cosine coefficients

N = length(image);

idct = zeros(N);

for u=0: N-1
    for v=0: N-1
        if u==0
           cu = 1/sqrt(N); 
        else
           cu = sqrt(2/N);
        end
        if v==0
           cv = 1/sqrt(N); 
        else
           cv = sqrt(2/N);
        end
        
        mysum = 0;
        for x=0:N-1
            for y=0:N-1
                temp = cu*cv*image(x+1,y+1)*cos(pi*u*(2*x+1)/(2*N))*cos(pi*v*(2*y+1)/(2*N));
                
                %temp = image(k,l)*cos(((2*(k-1)+1)*i*pi)/(2*n))*cos(((2*(l-1)+1)*j*pi)/(2*n));
                mysum = mysum+ temp;
            end
        end
        idct(u+1,v+1) = mysum;
    end
end

idct = uint8(idct);

end

Please let me know if you see an error or I am implementing the functions incorectly.

Thanks

-Matt

Edit: Here is the code that generated the images posted as an example:

I = imread('cameraman.tif');
temp = my_dct_2d(I);
newImage = my_idct_2d(temp);
imshow([I newImage]);

[Matt J]

You need to show the code that generated the smeared images as well

[speck]

This is now included

Answer 1

[Matt J]

Your forward and inverse DCT code look identical. In other words, you appear to be performing the forward DCT where you intended to do the inverse.

[Matt J]

Your code is also quite slow, making testing difficult. There are opportunities here to exploit vectorization and separability. Below is how I would implement it. You could probably do even better by writing it in terms of FFTs.

I = imread('cameraman.tif');
temp = my_dct_2d(I);
newImage = uint8( my_idct_2d(temp));
imshow([I newImage]);


function A=dctmatrix_1d(N)

     c(1:N,1)=sqrt(2/N);  c(1)=c(1)/sqrt(2);
     
     e=(0:N-1);
     
     A=c.*cos(  e(:)  *  (pi*(2*e+1)/(2*N))  );


end


 function dct = my_dct_2d(image)

  image = double(image);
  N = length(image);

  A=dctmatrix_1d(N);
  
  dct=A*image*A.';
  
  
end
 
function idct = my_idct_2d(image)

  image = double(image);
  N = length(image);

  A=dctmatrix_1d(N);
  
  idct=A\image/A.';

end

[Matt J]

Another, more compact way of implementing the above is by using my KronProd class.

I = double(imread('cameraman.tif'));

K=dct_kronprod(length(I));  %This can be re-used later on images of the same size

tic
 temp = K*I;
 newImage = uint8( K\temp );
toc %Elapsed time is 0.005276 seconds.

imshow([I newImage]);


function K=dct_kronprod(N)

     c(1:N,1)=sqrt(2/N);  c(1)=c(1)/sqrt(2);
     
     e=(0:N-1);
     
     A=c.*cos(  e(:)  *  (pi*(2*e+1)/(2*N))  );
     
     K=KronProd({A},[1,1]);
     
end