如何实现快速傅里叶变换以相关两个二维数组？

Question

我希望能够使用更大的数组来实现这个功能，但为了方便理解，我将使用一个较小的示例。

假设我有以下数组：

    A = [[0, 0, 100],
         [0, 0, 0],
         [0, 0, 0]] 

如果我想计算这个数组与另一个数组的相关性，可以通过乘以相应的条目来实现。例如，A * A 等于 A（1 * 1 = 1，其他地方都是零）。我读到快速傅里叶变换可以用于加速大型数组的操作。从我所读的内容中，如果我想要像 A*B 这样乘两个数组 A 和 B，我可以使用以下方法更快地完成（在 Python 中使用 numpy）：

a = np.conj(np.fft.fftn(A))
b = np.fft.fftn(B)
c = np.fft.ifft(a*b)

因此，实际上需要对A进行FFT，对B进行FFT，将两个结果相乘，然后获取该结果的反演。然而，我尝试了使用给定的A的情况，将其与自身相乘。我希望反向乘法会给我。

[[0, 0, 10000],
 [0, 0, 0    ],
 [0, 0, 0    ]]

然而，我得到了一些不同的东西，更接近于……

[[10000, 0, 0],
 [10000, 0, 0],
 [10000, 0, 0]]

有人知道发生了什么吗？抱歉，我猜我对快速傅里叶变换有些误解。

- Franz

---

我得到的输出结果与我预期的不同，我认为可以肯定地说这是因为我在这里做了一些我不理解的事情，而不是计算机出错了。 - Franz

你所看到的是边缘效应。为了正确处理，您需要在所有方向上将两个数组都用其大小的~1/2填充零。此外，我认为您想在最后一步（c）中使用 np.fft.ifftn。 - Joe Kington

抱歉，我本意是写ifftn。我会查看边缘效应，谢谢。 - Franz

等一下，我在结尾没有加上那个“n”...现在尝试了一下，看起来它的工作效果更好了。现在左上角只有一个1000的条目。感谢你指出这个问题。 - Franz

@Franz - 很高兴它能正常工作！然而，通常情况下你仍需要担心边缘效应。请注意，fft/fft2/fftn等函数都接受一个形状参数，如果你传递想要填充的形状，它将处理零填充。 - Joe Kington

Answer 1

你应该使用 scipy.signal.fftconvolve 代替。它已经被实现并且经过了广泛的测试，特别是关于边界处理。在二维卷积运算中，从卷积运算到相关运算所需要的唯一额外步骤是将滤波器数组旋转180°(参见这个答案)。

- rth

Answer 2

如果你必须自己实现，你应该知道在频域中的乘法对应于时域中的循环卷积。为了获得所需的线性卷积，你需要用零填充两个数组，使它们的长度至少是原始矩阵大小的两倍¹：

s = [2*x for x in np.shape(A)]
a = np.conj(np.fft.fftn(A,s))
b = np.fft.fftn(B,s)
c = np.fft.ifftn(a*b)

严格来说，大小为2n-1（而不是2n）就可以了，但是FFT在处理由小质因数的倍数构成的大小时性能更好。

- SleuthEye