scipy.signal.fftconvolve没有给出所需的结果

我有一个关于python的fftconvolve的问题。在我目前的研究中，我被要求计算两个函数之间的一些卷积。为此，我使用傅立叶变换（我使用numpy.fft并对其进行归一化）进行计算。问题是，如果我想用fftconvolve包进行比较，它不能给出正确的结果。这里是我的代码：scipy.signal.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的结果关闭（但仍然没有）到正确的结果。

有什么建议吗？

谢谢！

来源

2013-10-08 user2857354

请转到['scipy.fftconvolve']（https://github.com/scipy/scipy/blob/master/scipy/signal/signaltools.py#L153）并观察算法没有你奇怪的变化或缩放你想用FFT功能实现什么？ –

当计算卷积时，边界条件很重要。

当您对两个信号进行卷积时，结果的边缘取决于您假定输入的边缘以外的值。 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.]

正向FFT的乘积的逆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'选项）。

来源

2014-01-04 23:33:20 Andrew

它看起来像更新版本的fftconvolve [做得更好，选择填补多少以获得快速的FFT计算]（https://github.com/scipy/scipy/blob/v0.15.1/scipy/signal/signaltools.py ＃L349） – Andrew

scipy.signal.fftconvolve没有给出所需的结果

回答

相关问题