代码之家 › 专栏 › 技术社区 › ggkmath

为什么卷积结果在时域和频域中有不同的长度?

fftw fft signal-processing

ggkmath · 技术社区 · 15 年前

例如,如果滤波器和波形都是N个点长,则第一种方法(即卷积)产生N+N-1个点长的结果,其中该响应的前半部分是滤波器填充,第二部分是滤波器清空。为了获得稳态响应,滤波器需要比待滤波波形具有更少的点。

用第二种方法继续这个例子,假设离散时域波形数据都是实的(不复杂),滤波器的FFT和波形都产生N个点长的FFT。如果将两个谱相乘得到的结果相乘,得到的时域结果也是N点长。在这里,滤波器填充和清空的响应在时域中相互重叠,没有稳态响应。这是循环卷积的效果。为了避免这种情况,通常滤波器尺寸将小于波形尺寸,并且两者都将被零填充,以便在两个频谱的乘积IFFT之后,允许频率卷积在时间上扩展的空间。

我的问题是,我经常在文献中看到来自知名专家/公司的工作,他们有一个离散(实时)时域波形(N个点),他们对其进行FFT处理,然后乘以一些滤波器(也就是N个点),然后将结果进行后续处理。我天真的想法是,这个结果不应该包含稳态响应,因此应该包含来自过滤器填充/清空的工件,这些工件会导致解释结果数据时出错,但是我一定遗漏了一些东西。在什么情况下,这是一种有效的方法?

4 回复 | 直到 10 年前

tom10 15 年前

基本问题不是关于零填充与假定周期性的关系,而是傅里叶分析将信号分解为正弦波,在最基本的水平上,假设正弦波的范围是无限的。 两种方法都是正确的 两种方法都不正确

回到你的第一段:通常,在DSP中,我遇到的最大问题是FFT是非因果的,因此我通常喜欢停留在时域中,例如使用FIR和IIR滤波器。

更新:

在问题陈述中,OP正确地指出了在使用fft来过滤信号时可能出现的一些问题,例如,边缘效应,当进行长度(在时域中)与采样波形相比较的卷积时,这些问题可能特别有问题。需要注意的是,并不是所有的过滤都是用fft完成的,而且 在OP引用的论文中,他们没有使用FFT滤波器,并且FFT滤波器实现中可能出现的问题也不会使用他们的方法出现

例如,考虑一个过滤器,它使用两个不同的实现实现128个采样点的简单平均值。

快速傅里叶变换 最后的 127个(即未来的)wfm样本(跳过wfm的中间部分),或者127个零(如果您将填充零)。两者都不正确。

:另一种方法是用FIR滤波器实现平均值。例如,这里可以使用128寄存器FIFO队列中值的平均值。也就是说,当每个采样点进入时,1)将其放入队列中,2)将最旧的项目出列,3)平均队列中剩余的128个项目;这是您对此采样点的结果。这种方法连续运行,每次处理一个点,并在每个采样后返回滤波结果,并且在应用于有限的样本块时,没有FFT所出现的问题。每个结果只是当前样本和之前127个样本的平均值。

this free book

Eryk Sun 15 年前

这里的例子是没有填充的线性卷积。这是长度M=32序列与长度L=128序列的卷积(使用Numpy/Matplotlib):

f = rand(32); g = rand(128)
h1 = convolve(f, g)
h2 = real(ifft(fft(f, 128)*fft(g)))
plot(h1); plot(h2,'r')
grid()

alt text 第一个M-1点是不同的,它比M-1点短,因为它不是零填充的。如果要进行块卷积,这些差异是一个问题,但是使用诸如overlap和save或overlap和add之类的技术来克服这个问题。否则,如果只计算一次性过滤操作,有效结果将从索引M-1开始,在索引L-1结束,长度为L-M+1。

编辑:关于共轭传递函数,锁相环有一个实值脉冲响应,每个实值信号都有一个共轭对称谱。在代码中你可以看到,他们只是使用Hfinal[N-i]来获得负频率,而不需要共轭。我绘制了-50 MHz到50 MHz的传输函数:

N = 150000                    # number of samples. Need >50k to get a good spectrum. 
res = 100e6/N                 # resolution of single freq point  
f = res * arange(-N/2, N/2)   # set the frequency sweep [-50MHz,50MHz), N points
s = 2j*pi*f                   # set the xfer function to complex radians 

f1 = 22e6       # define 3dB corner frequency for H1 
zeta1 = 0.54    # define peaking for H1 
f2 = 7e6        # define 3dB corner frequency for H2 
zeta2 = 0.54    # define peaking for H2    
f3 = 1.0e6      # define 3dB corner frequency for H3 

# w1 = natural frequency   
w1 = 2*pi*f1/((1 + 2*zeta1**2 + ((1 + 2*zeta1**2)**2 + 1)**0.5)**0.5)  
# H1 transfer function 
H1 = ((2*zeta1*w1*s + w1**2)/(s**2 + 2*zeta1*w1*s + w1**2))            

# w2 = natural frequency 
w2 = 2*pi*f2/((1 + 2*zeta2**2 + ((1 + 2*zeta2**2)**2 + 1)**0.5)**0.5)  
# H2 transfer function  
H2 = ((2*zeta2*w2*s + w2**2)/(s**2 + 2*zeta2*w2*s + w2**2))            

w3 = 2*pi*f3        # w3 = 3dB point for a single pole high pass function. 
H3 = s/(s+w3)       # the H3 xfer function is a high pass

Ht = 2*(H1-H2)*H3   # Final transfer based on the difference functions

subplot(311); plot(f, abs(Ht)); ylabel("abs")
subplot(312); plot(f, real(Ht)); ylabel("real")
subplot(313); plot(f, imag(Ht)); ylabel("imag")

alt text

hotpaw2 15 年前

虽然假设数据的矩形窗口在FFT孔径宽度处是周期性的,这是圆卷积在没有足够的零填充的情况下所起作用的一种解释,但是这种差异可能大到足以淹没所讨论的数据分析。

Peter Tillemans 15 年前

我可以解释一下为什么在应用FFT之前应用“加窗”的原因。

如前所述,FFT假设我们有一个无穷大的信号。当我们在有限时间T上取一个样本,这在数学上相当于用一个矩形函数乘以信号。

在时域中的乘法在频域中变成了卷积。矩形的频率响应是同步函数,即sin(x)/x。分子中的x是kicker,因为它衰减为O(1/N)。

如果频率分量正好是1/T的倍数,这并不重要,因为同步函数在所有点上都是零,除了频率是1。

然而,如果你有一个正弦在2点之间,你会看到同步功能采样在频率点。它像是一个放大版的同步功能和'鬼'信号造成的卷积下降到1/N或6dB/倍频程。如果你有一个信号在噪音底限以上60分贝,你将看不到主信号左右1000个频率的噪音,它将被同步功能的“裙子”淹没。

如果你使用不同的时间窗口,你会得到不同的频率响应,比如余弦会以1/x^2消失,对于不同的测量有专门的窗口。汉宁窗常用作通用窗。

关键在于,不应用任何“窗口函数”时使用的矩形窗口会产生比选择良好的窗口更糟糕的人工制品。i、通过“扭曲”时间样本,我们可以在更接近“真实”的频率域中获得更好的图像,或者更接近于我们期望和希望看到的“真实”。