应用窗函数后FFT发生意外的频移 [英] FFT unexpected frequency shift after window function application
问题描述
我得到了这个python代码，用于声音信号的FFT计算:
I got this python code for FFT calculation of a sound signal:
from math import *
from cmath import exp, pi
def fft(x):
N = len(x)
if N <= 1: return x
even = fft(x[0::2])
odd = fft(x[1::2])
return ([even[k] + exp(-2j * pi * k / N) * odd[k] for k in xrange(N / 2)] +
[even[k] - exp(-2j * pi * k / N) * odd[k] for k in xrange(N / 2)])
N = 64
res = [sin(k) for k in xrange(N)]
# Window function
a = 2*pi/(N-1)
for k in xrange(N):
z = a * res[k]
res[k] = 0.42659 - 0.49656*cos(z) + 0.076849*cos(2*z)
res = fft(res)
for k in xrange(N/2):
# get the amplitude...
sqr = sqrt(res[k].real * res[k].real + res[k].imag * res[k].imag)
if sqr > 0:
print 20 * log10(sqr) # ...in decibels
else:
print "-INF"
我得到以下结果:
没有窗口功能(注释掉):
WITHOUT window function (commented-out):
-20.3017238269
-16.9192604108
-12.5089302395
-8.97999530657
-5.96033201086
-3.12975820108
-0.242090896634
2.97021879504
6.95134203457
12.8752188937
29.5096108632 <-- PEAK
17.1668404562
10.6485650284
7.24321329787
4.98448122464
3.3242079033
2.03154022635
0.987966110459
0.124898554197
-0.600705355004
-1.21748302238
-1.74534177237
-2.1985940834
-2.5878009699
-2.92091067118
-3.20399051424
-3.44171254421
-3.63768393032
-3.79467588076
-3.91478386211
-3.99953964778
-4.04998822971
具有窗口功能:
-6.55943077129
-65.8567720414
-65.7987645827
-65.7012678903
-65.5625673034
-65.380788761
-65.1529344157
-64.8750852394
-64.5420675211
-64.1470597764
-63.6810080181
-63.131731575
-62.4825087571
-61.7097419947
-60.7788888801
-59.6368610687
-58.1964601495
-56.3001921054
-53.6185951634
-49.2554491173
-38.3322646561 <-- PEAK
-43.3318138698
-52.0838904305
-56.7277347745
-60.2038755771
-62.9772322874
-65.442363488
-67.7550361967
-70.0212827894
-72.3056579688
-74.5822818952
-76.5522909937
由于某种原因，峰出现了偏移.这是2倍频移！
为了检查结果，我尝试了以下Java小程序:
http://www.random-science-tools.com/maths/FFT.htm
似乎没有任何窗口函数的结果都是正确的(峰值位于频谱的三分之一处).相反，如果我在python脚本中应用window函数，则峰值显示在频谱的2/3处.
The peak appears shifted for some reason. It is a 2x frequency shift!
To check the results, I tried this Java applet:
http://www.random-science-tools.com/maths/FFT.htm
And it appears that the results WITHOUT any window function are the correct ones (peak at 1 third of the spectrum). Instead, if I apply the window function in my python script the peak shows at 2/3 of the spectrum.
这应该发生吗?我该如何解决?
Is this supposed to happen? How do I fix it?
推荐答案
好的，与此同时，我意识到出了什么问题.我在问题中写的window函数完全没有意义.
Ok, In the meanwhile I realized what was wrong. The window function as I wrote it in the question was totally meaningless.
这是正确的:
a = 2*pi/(N-1)
for k in xrange(N):
z = a * k
res[k] *= 0.42659 - 0.49656*cos(z) + 0.076849*cos(2*z) # Blackman
结果:
-63.8888312044
-62.1859660802
-59.4560808775
-57.5235455007
-57.0010514385
-59.4284419437
-66.6535724743
-46.1441434426
-2.31562840406
16.0873761957
22.4136439765 <-- PEAK
19.5784749467
6.43274013629
-28.3842042716
-55.5273291654
-68.8982705127
-53.3843989911
-49.731974213
-48.3131204305
-47.6953570892
-47.4386151256
-47.361972079
-47.3787962267
-47.4434419084
-47.530228024
-47.6240076874
-47.7155325706
-47.799012933
-47.870764286
-47.9284264139
-47.9705003855
-47.9960714351
现在峰正好位于其应有的位置.
The peak is now exactly where it is supposed to be.
您可能想尝试其他一些窗口:
Some other windows you may want to try:
res[k] *= 0.355768 - 0.487396*cos(z) + 0.144232*cos(2*z) - 0.012604*cos(3*z)
res[k] *= 1 - 1.93*cos(z) + 1.29*cos(2*z) - 0.388*cos(3*z) + 0.028*cos(4*z)
res[k] *= 1 - 1.985844164102*cos(z) + 1.791176438506*cos(2*z) - 1.282075284005*cos(3*z) + 0.667777530266*cos(4*z) - 0.240160796576*cos(5*z) + 0.056656381764*cos(6*z) - 0.008134974479*cos(7*z) + 0.000624544650*cos(8*z) - 0.000019808998*cos(9*z) + 0.000000132974*cos(10*z)
顺序:Nuttall，FTSRS，HFT248D.
In order: Nuttall, FTSRS, HFT248D.
https://holometer.fnal.gov/GH_FFT.pdf
这篇关于应用窗函数后FFT发生意外的频移的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
