MATLAB中信号的内核二维卷积答案

【问题标题】：Kernel 2D Convolution of signal in MATLABMATLAB中信号的内核二维卷积
【发布时间】：2014-03-08 08:23:11
【问题描述】：

我想在 MATLAB 中实现信号的二维核卷积。这本质上是一个信号（不是图像）的非线性二次滤波器，如下式所示：

其中 k2(t1,t2) 是一个二维卷积核。 x 是 1d [N,1] 信号向量 & y 是输出 [N,1] 信号。

到目前为止，我一直在以一种非常粗暴、不优雅的方式来做这件事。我想知道我是否可以使用 MATLAB filter2/conv2 函数更有效地做到这一点！我知道这些功能是用于图像处理的，我不太了解，所以我希望有人能提供帮助！

【问题讨论】：

M 和 N 一样吗？
不，M 是内核大小（系统内存）。 N是信号长度。棒球场人物，M=30 个箱子，N=10,000 个箱子。

标签： matlab image-processing signal-processing convolution

【解决方案1】：

可能如下

y=diag(conv2(k,x*x.'));

或

y=diag(conv2(x,x,k));

【讨论】：

真正的优化需要知道 x 和 k 的相对大小。

【解决方案2】：

使用diag(conv2(x,x,k)); 时的问题是您正在计算更大的东西（整个二维矩阵），然后只保留对角线。根据信号的大小，它可能会很昂贵。你可以试试

n = 500; m = 50;

x = rand(n,1);
k2 = rand(m,m);

tic; res1 = diag(conv2(x,x,k2)); toc;

tic;
res2 = zeros(n+m-1,1);
for k = 1:n+m-1
    imin = k+1-min(k,m); imax = min(k,n);
    jmin = max(1,k-n+1); jmax = min(k,m);
    res2(k) = x(imax:-1:imin)'*k2(jmin:jmax,jmin:jmax)*x(imax:-1:imin);
end
toc;
norm(res1-res2)

在我尝试过的许多情况下，它比其他选项运行得更快。例如，一个输出可以是

>> script
Elapsed time is 0.012753 seconds.
Elapsed time is 0.006541 seconds.

ans =

   1.5059e-12

不知道你的信号或者内核有多大，你可以试试。

【讨论】：

棒球场人物，M=30 箱，N=10,000 箱。
@DankMasterDan 不要尝试使用这种尺寸的 conv2。它可能永远不会在普通计算机中结束。
@DankMasterDan 稍作修改就可以提高效率。只需添加括号即可使 matlab 矩阵向量与未转置的向量相乘。 x(imax:-1:imin)'*(k2(jmin:jmax,jmin:jmax)*x(imax:-1:imin)); 没有括号总是更快。
感谢塞巴斯！它不像我想象的那样优雅，但它把我的计算时间减少了一半！
不客气。下次你也可以用“算法”标记你的问题，你可能会引起比我更好的程序员的注意。

【解决方案3】：

对于您提到的大小，使用 FFT 可能是值得的，尽管这相当于循环卷积，因此您必须确保 N 包含足够的零填充。我假设在下面。

由于可分离性，您只需要 x 上的一维 fft

  X=fft(x);
  K=fft2(k,N,N);

  XX=X*X.';

  ytmp=ifft2(XX.*K);
  y=diag(ytmp);

正如 sebas 指出的那样，除了 ytmp 的对角线之外，您将丢弃所有内容，因此会浪费一些计算。如果您要对长度为 N 的许多信号 x 重复此操作，您可以通过预先计算一个 IDFT 矩阵来抑制这种浪费，

  Afft=ifft2(eye(N)).';

然后，不用计算ytmp，你可以直接计算ifft的对角线

  y=sum(Afft2.*ifft( XX.*(K.')) );

即使您不重复计算多个 x，FFT 的 NlogN 效率加上避免循环仍可能弥补浪费。

【讨论】：

反 fft 后需要得到 real(ytmp) 吗？
我不知道是否可以预期实值结果。这取决于 k 的对称性，而 OP 没有告诉我们这一点。但是，如果预期输出是实值，那么是的，您将不得不丢弃虚部中的任何浮点残差，就像使用 ifft() 等一样。