2D Convolution Matrix in Matlab

matrix_image_conv.m

% Specify image
im_size = [128,128];
X = double(imread('peppers.png'));
X = imresize(X,im_size);

% Grab some kind of filter/PSF
PSF_dim = [4,4];
PSF = ones(PSF_dim);
PSF = PSF ./ norm(PSF(:));

% Generate the matrix
H = convmtx2(PSF,im_size);

% Calculate convolution
tic
Y = H*X(:);
matrix_time = toc;

% Reshape the solution
new_dim = [(im_size(1)+PSF_dim(1)-1) (im_size(2)+PSF_dim(2)-1)];
Y = reshape(Y,new_dim);

% Comparison -- Using conv2
tic
Y2 = conv2(X,PSF,'full');
function_time = toc;
Y2 = reshape(Y2,new_dim);

% Display the result
figure(1);
  subplot(1,3,1);
    imagesc(X); 
    axis image;
    title('Original');
  subplot(1,3,2);
    imagesc(Y);
    axis image;
    title('Matrix Convolution');
    xlabel(sprintf('Compute Time = %0.2e',matrix_time));
  subplot(1,3,3);
    imagesc(Y2);
    axis image;
    title('Functional Convolution');
    xlabel(sprintf('Compute Time = %0.2e',function_time));
colormap(gray);