python - 使用 fftconvolve 计算包裹的二维相关性

Question

我有一组二维数组，我必须计算它们的二维相关性。我一直在尝试许多不同的事情(甚至用 Fortran 进行编程)，但我认为最快的方法是使用 FFT 进行计算。

基于我的测试和 this answer我可以使用scipy.signal.fftconvolve如果我尝试重现 scipy.signal.correlate2d 的输出，它工作得很好。与 boundary='fill' 。所以基本上是这样的

scipy.signal.fftconvolve(a, a[::-1, ::-1], mode='same')

等于这个(除了轻微的偏移)

scipy.signal.correlate2d(a, a, boundary='fill', mode='same')

问题是数组应该以包装模式计算，因为它们是二维周期性数组(即 boundary='wrap' )。因此，如果我尝试重现输出

scipy.signal.correlate2d(a, a, boundary='wrap', mode='same')

我不能，或者至少我不知道该怎么做。 (我想使用 FFT 方法，因为它更快。)

显然 Scipy 曾经有 something like that这可能已经成功了，但显然它被遗忘了，我找不到它，所以我认为 Scipy 可能已经放弃了对它的支持。

无论如何，有没有办法使用scipy或numpy的 FFT 例程来计算周期数组的相关性？

Answer 1

包裹相关性可以使用FFT来实现。下面是一些演示如何操作的代码:

In [276]: import numpy as np

In [277]: from scipy.signal import correlate2d

创建一个随机数组a来使用:

In [278]: a = np.random.randn(200, 200)

使用scipy.signal.correlate2d计算二维相关性:

In [279]: c = correlate2d(a, a, boundary='wrap', mode='same')

现在使用 numpy.fft 中的 2D FFT 函数计算相同的结果。 (此代码假设 a 是正方形。)

In [280]: from numpy.fft import fft2, ifft2

In [281]: fc = np.roll(ifft2(fft2(a).conj()*fft2(a)).real, (a.shape[0] - 1)//2, axis=(0,1))

验证两种方法是否给出相同的结果:

In [282]: np.allclose(c, fc)
Out[282]: True

正如您所指出的，使用 FFT 速度快得多。对于此示例，速度大约快 1000 倍:

In [283]: %timeit c = correlate2d(a, a, boundary='wrap', mode='same')
1 loop, best of 3: 3.2 s per loop

In [284]: %timeit fc = np.roll(ifft2(fft2(a).conj()*fft2(a)).real, (a.shape[0] - 1)//2, axis=(0,1))
100 loops, best of 3: 3.19 ms per loop

其中包括 fft2(a) 的重复计算。当然，fft2(a) 应该只计算一次:

In [285]: fta = fft2(a)

In [286]: fc = np.roll(ifft2(fta.conj()*fta).real, (a.shape[0] - 1)//2, axis=(0,1))