scipy.ndimage.filters.convolve 和乘以傅里叶变换给出不同的结果

Question

这是我的代码：

from scipy.ndimage import filters
import numpy

a = numpy.array([[2,43,42,123,461],[453,12,111,123,55] ,[123,112,233,12,255]])
b = numpy.array([[0,2,2,3,0],[0,15,12,100,0],[0,45,32,22,0]])

ab = filters.convolve(a,b, mode='constant', cval=0)

af = numpy.fft.fftn(a)
bf = numpy.fft.fftn(b)

abf = af*bf

abif = numpy.fft.ifftn(abf)

print numpy.around(ab)
print numpy.around(abif)

结果是：

[[ 1599  2951  7153 13280 18311]
 [ 8085 51478 13028 40239 30964]
 [18192 32484 23527 36122  8726]]

[[ 37416.+0.j  32251.+0.j  46375.+0.j  32660.+0.j  23986.+0.j]
 [ 30265.+0.j  33206.+0.j  62450.+0.j  19726.+0.j  17613.+0.j]
 [ 40239.+0.j  38095.+0.j  24492.+0.j  51478.+0.j  13028.+0.j]]

如何修正我使用 FFT 进行卷积的方式，以保证它给出与 scipy.ndimage.filters.convolve 相同的结果？

谢谢。

Answer 1

事实证明这是一个有趣的问题。似乎使用离散傅里叶变换的卷积（由 numpy.fft.fftn 实现）相当于 circular convolution。所以我们需要做的就是使用'wrap'卷积模式，并适当设置原点：

>>> filters.convolve(a, b, mode='wrap', origin=(-1, -2))
array([[37416, 32251, 46375, 32660, 23986],
       [30265, 33206, 62450, 19726, 17613],
       [40239, 38095, 24492, 51478, 13028]])
>>> numpy.fft.ifftn(numpy.fft.fftn(a) * numpy.fft.fftn(b))
array([[ 37416.+0.j,  32251.+0.j,  46375.+0.j,  32660.+0.j,  23986.+0.j],
       [ 30265.+0.j,  33206.+0.j,  62450.+0.j,  19726.+0.j,  17613.+0.j],
       [ 40239.+0.j,  38095.+0.j,  24492.+0.j,  51478.+0.j,  13028.+0.j]])
>>> (filters.convolve(a, b, mode='wrap', origin=(-1, -2)) ==
...  numpy.around(numpy.fft.ifftn(numpy.fft.fftn(a) * numpy.fft.fftn(b))))
array([[ True,  True,  True,  True,  True],
       [ True,  True,  True,  True,  True],
       [ True,  True,  True,  True,  True]], dtype=bool)

差异仅与 filters.convolve 处理边缘的方式有关。一种使用 fftn 在其他模式下执行卷积的方法并没有立即打动我；有关解决该问题的聪明（事后诸葛亮）方法，请参阅的出色答案。

Answer 2

正如@senderle 指出的那样，当您使用 FFT 实现卷积时，您会得到 circular 卷积。 @senderle 的回答显示了如何调整 filters.convolve 的参数以进行循环卷积。要修改 FFT 计算以生成与您最初使用 filters.convolve 相同的结果，您可以用 0 填充参数，然后提取结果的适当部分：

from scipy.ndimage import filters
import numpy

a = numpy.array([[2.0,43,42,123,461], [453,12,111,123,55], [123,112,233,12,255]])
b = numpy.array([[0.0,2,2,3,0], [0,15,12,100,0], [0,45,32,22,0]])

ab = filters.convolve(a,b, mode='constant', cval=0)

print numpy.around(ab)
print

nrows, ncols = a.shape
# Assume b has the same shape as a.
# Pad the bottom and right side of a and b with zeros.
pa = numpy.pad(a, ((0, nrows-1), (0, ncols-1)), mode='constant')
pb = numpy.pad(b, ((0, nrows-1), (0, ncols-1)), mode='constant')
paf = numpy.fft.fftn(pa)
pbf = numpy.fft.fftn(pb)
pabf = paf*pbf
p0 = nrows // 2
p1 = ncols // 2
pabif = numpy.fft.ifftn(pabf).real[p0:p0+nrows, p1:p1+ncols]

print pabif

输出：

[[  1599.   2951.   7153.  13280.  18311.]
 [  8085.  51478.  13028.  40239.  30964.]
 [ 18192.  32484.  23527.  36122.   8726.]]

[[  1599.   2951.   7153.  13280.  18311.]
 [  8085.  51478.  13028.  40239.  30964.]
 [ 18192.  32484.  23527.  36122.   8726.]]

scipy.ndimage.filters.convolve 和乘以傅里叶变换给出不同的结果

scipy.ndimage.filters.convolve and multiplying Fourier Transforms give different results

python

numpy

convolution

scipy