Matlab中的2D卷积 - 代码优化

Question

这是我们的图像处理作业练习。我的代码可以正常工作。我希望能够得到一些关于代码优化的帮助。

function C = convolve_slow(A,B)
(file name is accordingly convolve_slow.m ) 
This routine performs convolution between an image A and a mask B.
Input:      A - a grayscale image (values in [0,255]) 
            B - a grayscale image (values in [0,255]) serves as a mask in the convolution.
Output:     C - a grayscale image (values in [0,255]) - the output of the convolution. 
                      C is the same size as A.

Method:  Convolve A with mask B using zero padding. Assume the origin of B is at 
     floor(size(B)/2)+1.
Do NOT use matlab convolution routines (conv,conv2,filter2 etc). 
Make the routine as efficient as possible: Restrict usage of for loops which are expensive (use matrix multiplications and matlab routines such as dot etc).
To simplify and reduce ifs, you should pad the image with zeros before starting your convolution loop.
Do not assume the size of A nor B (B might actually be larger than A sometimes).

这是我们的解决方案

function [ C ] = convolve_slow( A,B )
%This routine performs convolution between an image A and a mask B.
% Input:      A - a grayscale image (values in [0,255])
%             B - a grayscale image (values in [0,255]) serves as a mask in the convolution.
% Output:     C - a grayscale image (values in [0,255]) - the output of the convolution. 
%             C is the same size as A.
% 
% Method:  Convolve A with mask B using zero padding. Assume the origin of B is at floor(size(B)/2)+1.
% init C to size A with zeros
C = zeros(size(A));
% make b xy-reflection and vector
vectB = reshape(flipdim(flipdim(B,1),2)' ,[] , 1);
% padding A with zeros
paddedA = padarray(A, [floor(size(B,1)/2) floor(size(B,2)/2)]);
% Loop over A matrix:
for i = 1:size(A,1)
    for j = 1:size(A,2)
        startAi = i;
        finishAi = i + size(B,1) - 1;
        startAj = j;
        finishAj = j + size(B,2) - 1;
        vectPaddedA = reshape(paddedA(startAi :finishAi,startAj:finishAj)',1,[]);
        C(i,j) = vectPaddedA* vectB;
    end
end
end  

因为我是新手，对于图像处理和Matlab还不熟悉。你能帮助我优化代码吗？特别是关于矩阵运算方面。是否有可能不使用循环语句？

- Gilad

---

你可以以任何方式使用fft2和ifft2函数吗？如果可以的话，这将是一个相当快速的解决方案。 - MarkV

不好意思，没有快速傅里叶变换（FFT），主要的概念是我们实现卷积以了解它的工作原理。好吧，它能工作！现在我想看看是否有更好的方法来做到这一点。 - Gilad

Answer 1

不需要明确编写代码，我可以想到一种方法将其简化为一个主要的for循环。基本上，将矩阵A和B视为列向量，通过将A和B的每列展开成向量（这是在MATLAB中内部存储的方式）来实现。然后，可以使用函数sub2ind将A的每个坐标映射到线性索引k。接下来，在A的主体（忽略填充）中的每个线性索引处计算与该线性索引周围子矩阵对应的线性索引列表（这可能是最难的部分）。然后计算A( theseIndices )和B(:)的点积。使用此方法，您只需循环遍历A的每个线性索引。

- MarkV

Answer 2

我不确定这个方法是否更快，但至少它没有使用 for 循环（也就是说，在最近的 MATLAB版本中不一定意味着它一定更快）。

function A = tmpConv(A,B)

    filterSize = size(B,1);
    filterSize2 = floor(filterSize/2);
    inputSize = size(A);

    A = padarray(A,[filterSize2 filterSize2]);

    f = repmat(B(:),[1 inputSize(1)*inputSize(2)]);
    A = im2col(A,[filterSize filterSize]);
    A = reshape(sum(A.*f),inputSize);

- user1958650

---

这个计算的是相关性，而不是卷积。如果要进行卷积，您需要将滤波器旋转180度。此外，构建“f”矩阵并不是必要的，因为使用“A'”与（反转的）滤波器向量进行矩阵乘法可以得到相同的结果。 - Joost