2011-02-16 128 views
7

我想了解FFT算法,到目前为止我认为我理解它背后的主要概念。不过,我对'framesize'和'window'之间的区别感到困惑。与FFT算法混淆

根据我的理解,它们似乎彼此是多余的?例如,我将输入的一个样本块的帧大小设置为1024.因此,我将字节[1024]显示为输入。

那么窗口函数的目的是什么?自从开始以来,我认为窗口函数的目的是从原始数据中选择样本块。

谢谢!

回答

13

那么窗口函数的目的是什么?

这是要处理所谓的“频谱泄漏”:FFT假定无限序列,重复给定的样本帧一遍又一遍。如果你的样本帧中有一个正弦波是整数个周期,那么所有的都是好的,FFT在适当的频率给你一个很好的窄峰。但是如果你有一个不是整数个周期的正弦波,最后一个和第一个样本之间会有一个不连续点,而FFT会给你提供假的谐波。

窗口函数降低了采样帧开始和结束时的幅度,以减少由这种不连续性引起的谐波。

National Instruments webpage on windowing一些图:

周期积分#:循环

enter image description here

非整数#:

enter image description here

的其他信息:

http://www.tmworld.com/article/322450-Windowing_Functions_Improve_FFT_Results_Part_I.php

http://zone.ni.com/reference/en-XX/help/371361B-01/lvanlsconcepts/char_smoothing_windows/

http://www.physik.uni-wuerzburg.de/~praktiku/Anleitung/Fremde/ANO14.pdf

+0

稀释我明白了。所以基本上,窗口函数只是为了改进FFT算法的结果而修改或改进样本中的数据。我对吗?谢谢! – user488792 2011-02-16 04:00:04

+0

@ user488792:是的,这几乎就是它 - 窗口函数*平滑了否则在采样窗口边缘会出现瞬态。 – 2011-02-16 09:11:14

+0

是的,我明白了。非常感谢您的澄清和提供一个很好的答案! – user488792 2011-02-16 13:54:52

2

如果不修改的样本值,并选择数据作为FFT长度的相同的长度,这等同于使用矩形窗口,在其中情况下框架和窗口是相同的。然而,将输入数据乘以时域中的矩形窗口与将输入信号的频谱与频域中的Sinc函数进行卷积相同,这将在频域中的FFT孔径中不完全周期性地频散任何频谱峰值整个频谱。

通常使用非矩形窗口,因此得到的FFT频谱与Sinc函数相比有更多“聚焦”的卷积。

您还可以使用与FFT长度或光圈不同的矩形窗口。在数据窗口较短的情况下,FFT帧可以被填充为零,这可以导致更平滑的内插FFT结果频谱。您甚至可以使用矩形窗口,该窗口的长度大于FFT的长度,方法是通过以循环合计的方式将数据包裹在FFT窗口周围,以获得频率分辨率的一些有趣效果。

ADDED由于请求:

通过在时域中的窗口相乘产生相同的结果作为与该窗口的在频域中变换进行卷积。

一般来说,一个更窄的时域窗口产生更宽的频域变换。这是零填充产生更平滑频率图的原因。较窄的时域窗口产生更宽的Sinc,与帧宽度相比,相对于帧宽度更丰富和更平滑的曲线,从而使内插频率结果看起来比其相同的非零填充FFT更平滑帧长度。

反过来在一定程度上也是如此。更宽的矩形窗口将产生更窄的Sinc,空值更接近峰值。因此,您可能可以使用精心选择的较宽窗口来产生更窄的外观Sinc,以便将频率更接近感兴趣的频段的频率设置为低于1频率范围。你如何使用更宽的窗户?将数据包裹起来并进行求和,这与使用未被截断为1FFT帧长的FT基向量相同。但是,由于这样做时FFT结果向量比数据更短,这是一个有损耗的过程,它将引入伪像,并引入一些新的混叠。但它会给你在每个bin上更尖锐的频率选择峰值,并且可以放置在小于1 bin的位置上的陷波滤波器,比如在bin之间的中间等。

6

长度为M的矩形窗口的频率响应为sin(ω*M/2)/sin(ω/2),对于k≠0,当ω = 2*π*k/M为零。对于长度为N的DFT,其中ω = 2*π*n/N,在n = k * N/M处存在空值。比率N/M不一定是整数。例如,如果N = 40,M = 32,那么在1.25的倍数处有零点,但在DFT中将仅出现整数倍数,在这种情况下是DFT 5,10,15和20。

这里的一个32点的矩形窗口的1024点DFT的曲线图:

M = 32 
N = 1024 
w = ones(M) 
W = rfft(w, N) 
K = N/M 
nulls = abs(W[K::K]) 
plot(abs(W)) 
plot(r_[K:N/2+1:K], nulls, 'ro') 
xticks(r_[:512:64]) 
grid(); axis('tight') 

DFT of a 32-point rectangular window

注意在每一个N/M = 32的二进制位的空值。如果N = M(即窗口长度等于DFT长度),则除了n = 0之外的所有窗口中都有空值。

当您将窗口乘以信号时,频域中的对应操作是窗口频谱与信号频谱的循环卷积。例如,正弦曲线的DTFT是位于正弦曲线的正频率和负频率的加权德尔塔函数(即具有无限高度,无穷小延伸和有限区域的脉冲)。使用delta函数来描述频谱只是将其转换到三角洲的位置并按三角洲的权重进行缩放。因此,当您在样本域中用正弦曲线乘以一个窗口时,窗口的频率响应会缩放并转换为正弦曲线的频率。

有几种情况可以检查矩形窗口的长度。首先让我们看看窗口长度是正弦曲线周期的整数倍的情况,例如,余弦的32个样本的矩形窗口周期为32/8 = 4个样本:

x1 = cos(2*pi*8*r_[:32]/32) # ω0 = 8π/16, bin 8/32 * 1024 = 256 
X1 = rfft(x1 * w, 1024) 
plot(abs(X1)) 
xticks(r_[:513:64]) 
grid(); axis('tight') 

DFT of cosine at w0 = 8π/16 with a 32-point rectangular window

如前,有在N/M = 32,但该窗口的倍数空频谱已经转移到正弦曲线的二进制256,并按正弦频率和负频率(我只绘制正频率)之间的0.5倍进行缩放。如果DFT长度为32,那么空值将在每个垃圾箱中排队,这提示外观没有泄漏。但是这种误导性的外观只是DFT长度的一个函数。如果你用零填充窗口信号(如上),你将会看到空值之间频率的类似于sinc的响应。

现在让我们看看窗口长度不是正弦曲线周期的整数倍的情况,例如,用的7.5π/ 16的角频率的余弦(周期为64个样品):

x2 = cos(2*pi*15*r_[:32]/64) # ω0 = 7.5π/16, bin 15/64 * 1024 = 240 
X2 = rfft(x2 * w, 1024) 
plot(abs(X2)) 
xticks(r_[-16:513:64]) 
grid(); axis('tight') 

DFT of cosine at w0 = 7.5π/16 with a 32-point rectangular window

中心仓位置不再在32的整数倍,但由一个半错开下降到仓240.所以让我们看看相应的32点DFT是什么样子的(推断32点矩形窗口)。我会计算并绘制X2 [N]的32点DFT,并叠加1024点DFT的32倍抽取副本:

X2_32 = rfft(x2, 32) 
X2_sample = X2[::32] 
stem(r_[:17],abs(X2_32)) 
plot(abs(X2_sample), 'rs') # red squares 
grid(); axis([0,16,0,11]) 

32-point DFT of cosine at w0 = 7.5π/16

正如你在前面的情节看,零点不再以32的倍数对齐,所以32点DFT的大小在每个bin处都是非零的。在32点DFT中,窗口的零点仍然以N/M = 32/32 = 1 bin为间隔,但由于ω0=7.5π/ 16,中心位于'bin'7.5,零点为0.5,1.5等等,所以他们不在32点DFT中。

的一般信息是,一个窗信号的频谱泄漏是总是存在,但可以在DFT被掩蔽,如果信号specrtum,窗口长度,并且DFT长度在只是以正确的方式来一起排队空。除此之外,您应该忽略这些DFT伪像,并专注于信号的DTFT(即填充零点以更高分辨率对DTFT进行采样,以便您可以清楚地检查泄漏)。

由窗口光谱的卷积引起的光谱泄漏将始终存在,这就是为什么制作特殊形状窗户的艺术如此重要。每个窗口类型的频谱都针对特定任务进行了定制,例如动态范围或灵敏度。

这里的一个矩形窗的输出进行比较Vs的汉明窗的示例:

from pylab import * 
import wave 

fs = 44100 
M = 4096 
N = 16384 
# load a sample of guitar playing an open string 6 
# with a fundamental frequency of 82.4 Hz 
g = fromstring(wave.open('dist_gtr_6.wav').readframes(-1), 
       dtype='int16') 
L = len(g)/4 
g_t = g[L:L+M] 
g_t = g_t/float64(max(abs(g_t))) 
# compute the response with rectangular vs Hamming window 
g_rect = rfft(g_t, N) 
g_hamm = rfft(g_t * hamming(M), N) 

def make_plot(): 
    fmax = int(82.4 * 4.5/fs * N) # 4 harmonics 
    subplot(211); title('Rectangular Window') 
    plot(abs(g_rect[:fmax])); grid(); axis('tight') 
    subplot(212); title('Hamming Window') 
    plot(abs(g_hamm[:fmax])); grid(); axis('tight') 

if __name__ == "__main__": 
    make_plot() 

rectangular vs hamming window