关于python的fftconvolve，我有一个问题。在我目前的研究中，我被要求计算两个函数之间的卷积。为此，我使用fourier变换（我使用numpy.fft并对其进行规范化）计算它。问题是，如果我想使用fftconvolve包对其进行比较，它无法给出正确的结果。这是我的代码：
#!/usr/bin/python
import numpy as np
from scipy.signal import fftconvolve , convolve
def FFT(array , sign):
if sign==1:
return np.fft.fftshift(np.fft.fft(np.fft.fftshift(array))) * dw / (2.0 * np.pi)
elif sign==-1:
return np.fft.fftshift(np.fft.ifft(np.fft.fftshift(array))) * dt * len(array)
def convolve_arrays(array1,array2,sign):
sign = int(sign)
temp1 = FFT(array1 , sign,)
temp2 = FFT(array2 , sign,)
temp3 = np.multiply(temp1 , temp2)
return FFT(temp3 , -1 * sign) / (2. * np.pi)
""" EXAMPLE """
dt = .1
N = 2**17
t_max = N * dt / 2
time = dt * np.arange(-N / 2 , N / 2 , 1)
dw = 2. * np.pi / (N * dt)
w_max = N * dw / 2.
w = dw * np.arange(-N / 2 , N / 2 , 1)
eta_fourier = 1e-10
Gamma = 1.
epsilon = .5
omega = .5
G = zeros(N , complex)
G[:] = 1. / (w[:] - epsilon + 1j * eta_fourier)
D = zeros(N , complex)
D[:] = 1. / (w[:] - omega + 1j * eta_fourier) - 1. / (w[:] + omega + 1j * eta_fourier)
H = convolve_arrays(D , G , 1)
J = fftconvolve(D , G , mode = 'same') * np.pi / (2. * N)
如果你画出
H，J的实部/虚部，你会看到w轴的移动，而且我必须乘以J的结果，以便以某种方式接近（但仍然不是）正确的结果。有什么建议吗？
谢谢！
最佳答案
计算卷积时，边界条件很重要。
当您卷积两个信号时，结果的边缘取决于您在输入边缘之外假设的值。fftconvolve假设零填充边界计算卷积。
看看source code of fftconvolve。注意他们为了达到零填充边界条件而进行的恶作剧，特别是这些线：
size = s1 + s2 - 1
...
fsize = 2 ** np.ceil(np.log2(size)).astype(int) #For speed; often suboptimal!
fslice = tuple([slice(0, int(sz)) for sz in size])
...
ret = ifftn(fftn(in1, fsize) * fftn(in2, fsize))[fslice].copy()
...
return _centered(ret, s1) #strips off padding
这是好东西！如果你想理解基于傅里叶的卷积，那么仔细阅读
fftconvolve的代码可能是值得的，而且是一个很好的教育。小品
前向FFT零填充每个信号以防止周期性边界条件：
a = np.array([3, 4, 5])
b = np.fft.ifftn(np.fft.fftn(a, (5,))).real
print b #[ 3. 4. 5. 0. 0.]
前向FFTs的乘积的逆FFT给出了一个填充结果：
a = np.array([3, 4, 5])
b = np.array([0., 0.9, 0.1])
b = np.fft.ifftn(np.fft.fftn(a, (5,))*
np.fft.fftn(b, (5,))
).real
print b #[ 0. 2.7 3.9 4.9 0.5]
并且
_centered函数在末尾去掉额外的填充像素（假设使用mode='same'选项）。关于python - scipy.signal.fftconvolve没有给出所需的结果，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/19241021/