声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 5524|回复: 11

[FFT] FFT变换数据补零的问题

[复制链接]
发表于 2015-10-21 14:55 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
看资料说FFT分析的时候数据尽量去2的幂次方,这是为什么呢?
如果不够2的幂次方可以采用补零的办法补到2的幂次方,请问这样处理是否会给结果带来误差?
这个误差的大小怎么分析?
回复
分享到:

使用道具 举报

发表于 2015-10-22 09:02 | 显示全部楼层
补零这个问题网上讨论很多,先给你转一个总结性的帖子,你先看看

1.转自http://www.ilovematlab.cn/viewthread.php?tid=36349&sid=lqshA1

补零后的FFT和补零前的FFT,两者会有较大的不同,但两频谱的包络还是一致的。设补零前数据长N,补零后数据长M(补了M-N个零值),则补零前的FFT有N条谱线,分别代表的频率点是(0,1,...,N-1)*df1,df1=fs/N;补零前的FFT有M条谱线,分别代表的频率点是(0,1,...,M-1)*df2,df2=fs/M。由于补零前后数据长度不一样,它们的分辨率(分别为df1,df2)不一样,在频域中谱线所代表的频率也不一样,所以这两个频谱所描述的对象也不相同。这里举一个例子,fs=1000HZ,补零前后数据长度N=500和M=800,对应的df1=2,df2=1.25。补零前的频谱是对应于0,2,4,...,500HZ的频谱,而补零后的频谱是对应于0,1.25,2.5,...,500HZ的频谱,所以两频谱中对应频率不两同,描述当然不同。
但是当特殊情况:M=(2^n)N时,补零后的频谱相当于在补零前的频谱中插入(2^n)-1条谱线。与补零前的频谱中相重合的谱线,它们的幅值和相位完全一致。

2.下面这部分仅作为参考,暂时还没有弄懂是什么意思。转自:http://savasxch.spaces.live.com/ ... CC09016DB!449.entry

fft补零的提高分辨率吗? 不知你注意到没有?

评论一篇文章:一篇关于fft补零提高频率分辨率的讨论

这是一篇值得讨论的问题,作者认为补零fft可以提高频率分辨率,并给出了试验结果,可以看出确实提高了对频率细节的观察能力,本人可以肯定这个试验是真实的试验。但是所有的数字信号教课书上都认为补零fft并不能提高频率分辨率,是不是有矛盾?

1 从分析角度, 设fs为采样频率,fft长度为N, 那么频率分辨率为fs/N, 如果N增加那么频率分辨率增加。这是下面一篇文章的用的论据。

2 从另一角度,设fs为采样频率,fft长度为N, 则频率分辨率为fs/N, 我们引进另一个概念:时间长度DT(duration of time),  可以看出DT = 1/频率分辨率.     则频率分辨率=1/DT 。 从这一角度看只要DT不变,频率分辨率就不会变。因此尽管补零或插值,都不会提高分辨率。这是所有目前信号处理教科书的观点,但这些教科书都没有给出原因,不知道为什么,我发现这个问题是曾经找过不少教科书,没有一本给出原因,问老师也答的含糊不清。后来我反复考虑,感觉应该如此解释,如若有意见,欢迎讨论。

为什么两个角度看竟然矛盾?同样一个问题为什么有不同的解释

从1 我们看出,增加的值全为零,不是原信号内容,这就造成了特殊性,我们的信号变了不是原来信号了!而是新的补零信号的周期延拓。但是可以证明两个信号在对应点上的频谱值相同(直接利用定义即可推出)。至于补零后其它多点处的频谱是否是原信号的内容,这是问题的关键。

事实上,用于实用的方法不是下文里提到的方法,而是利用采样数据抽取,降低采样频率的方法来实现。因为数据长度一般在使用时都是最大长度,尤其是这种应用,肯定已经采用最大数据处理长度,不用问的。

我有一个试验,是多年前和一位同学讨论此类问题的试验,有兴趣的不妨试试.

%% 用于检验补零FFT是否提高分辨率

%%  结论:  1。 补零fft提高分辨率是指信号加窗后的合成信号的分辨率。

%     这种情况下fft可以帮助分辨真实的峰值,但分辨率你可以计算一下应该改是不变的。   

%           2。 如果信号=加窗后的合成信号   提高分辨率,如果加入点为真实的数据,当然提高分辨率。  

%           4。 提高分辨细节的本质是由于窗的展宽,窗分辨细节率的提高引起的。这是一个用大窗观 察含有小窗信号的小窗的过程。补零而看到的频率细节不是信号本身的细节,而是窗的细节。

  1.   clear;

  2. Nfft=16;

  3. span=[0:Nfft-1];

  4. omga=[1:3]*pi/8;

  5. x=exp(j*omga(1)*span)+exp(j*omga(2)*span)+exp(j*omga(3)*span);

  6. stem(pi/Nfft*[0:Nfft-1],abs(fft(x,Nfft)),'r');  figure(2);

  7. stem(pi/Nfft/2*[0:Nfft*2-1],abs(fft(x,Nfft*2)),'b');  figure(3);

  8. stem(pi/Nfft/16*[0:Nfft*16-1],abs(fft(x,Nfft*16)),'g');
复制代码


原文连接:http://www.ed-china.com/ART_8800009929_400014_500012_TS.HTM

原文简介:

采用C语言程序处理FFT算法无能为力的扩大频谱的问题

如果你只是对计算小频率范围上的频谱有兴趣,那HRFT可能比FFT更有效。频率范围越窄,HRFT的吸引力就越大。当必须利用有限的存储资源在嵌入式系统中执行运算时,这一点尤为明显。

许多科学和工程应用都要求信号准确的频谱或傅立叶变换。式1以无量纲形式(dimensionless form)给出,其中Ω = ΩT,T是采样间隔,q[n] = q(nT),代表信号q(t)的第n次采样。该式也假设信号为有限长序列,故总共仅有N 个连续采用点。Q(Ω)是连续变量Ω的周期函数,周期为2e。

常规求解方法是在区间Ω = [0, 2e]内的N个均匀间隔点上计算Q(Ω)的值。这个过程叫做离散傅立叶变换(DFT),通常利用快速傅利叶变换(FFT)的算法来完成。DFT给出点Ω = 2ek/N (k = 0...N -1)的傅立叶变换,从而按照实际频率给出了Ω=2e/(NT) 或f=1/(NT)的分辨率。对许多应用来说,这种分辨率可能已经足够。

对于某些应用,必须非常准确地确定频谱峰值(spectral peak)的位置。在这种情况下,DFT提供的分辨率可能不够。提高分辨率的方法之一就是简单地用额外的零来增加采样点。例如,为加倍频率分辨率,你可以在q[n]序列末尾增添N个零,使总的序列长度为2N。DFT将在2N个点上计算Q(Ω),从而使分辨率提高一倍。

为提高分辨率而增加的零的个数是有限制的。通常,只有在很小的频率范围内才需要提高分辨率。这里提供了一种简洁明了的办法,那就是直接在你感兴趣的频率范围上计算式1。下面给出有效算法。

采用欧拉恒等式,我们可以把式1分为实数和虚数两部分:

通过契比雪夫多项式(Chebyshev polynomials),我们从sinΩ和cosΩ开始,有效地进行sin(nΩ)和cos(nΩ)的递归运算。
在式4和式5中,已知T0 =1、T1=cosΩ、U0=1以及U1=2cosΩ,我们可以利用下列公式计算每一个Tn和Un现在我们只需要cosΩ和sinΩ,就能够以任意的高分辨率(远大于标准FFT)在一个非常小的频率范围上计算子频谱了。这一程序hrft用C语言编写,执行如下:

hrft f1 f2 Nf N sps infile outfile

其中,f1=开始频率(Hz),f2=终止频率(Hz),Nf=计算FT的频率点数, N=从输入文件读取的采样点数,sps=每秒采样点数,infile=输入文件,outfile=输出文件。

图中显示了由FFT算法产生的从2260Hz到2268 Hz的部分频谱,以及由高分辨率傅立叶变换(HRFT)在该频率范围上产生的频谱。

两种方法都采用了相同的2048点数据集。FFT具有4 Hz/bin的分辨率,HRFT计算了41个点,分辨率为0.2 Hz/bin。显然,41点的HRFT3比3点FFT的图更平滑。数据集的主要频率为2265 Hz,据此我们可以看出HRFT方法准确定位了频率峰值,而用FFT定位的峰值则偏离了1Hz。为了利用FFT获得同样好的分辨率,必须添加零,以达到65,536个点的序列长度。
如果你只是对计算小频率范围上的频谱有兴趣,那HRFT可能比FFT更有效。频率范围越窄,HRFT的吸引力就越大。当必须利用有限的存储资源在嵌入式系统中执行运算时,这一点尤为明显。

3.下面这篇文章的讲述容易理解一些

转自:http://xialulee.spaces.live.com/ ... ACFA82DB!1041.entry

时域补零对于DFT谱的影响

      节前总是疯狂加班,今天晚上还在实验室干到十点。其实这篇文章早就写好了,只是太忙,没时间上网。

      前几天上数字信号处理(本以为第二次上这个课只是简单地重复过去学习过的内容,但是这次有了很多新的发现),书上说对时域信号补零之后再作DFT并不能提高频谱的频率分辨率,提高采样频率也不能提高DFT谱的频率分辨率。这个很新鲜,以前上课时没有考虑过这个问题,以前的课本好像也没有开辟专门的章节论述这个问题。

      提高采样频率不能提高频率分辨率的原因其实很简单,因为提高了采样频率,虽然在相同的观察时长那的点数增多了,但与此同时采样频率也变大了,点数增加几倍采样频率增加几倍,所以不改变观察时长而仅仅提高采样频率并不能提高DFT谱的频率分辨率。

      但是时域补零呢?采样频率没有变化,而点数增加无疑会减小DFT谱的相邻谱线间隔,相邻谱线间隔的缩小为什么不能提高频率的分辨率呢?书上是这样写的:“错把‘计算分辨率’当成了‘物理分辨率’……补零没有对原信号增加任何新的信息,因此不可能提高分辨率。但补零……补零还可以对原X(k)做插值。”(《数字信号处理——理论、算法与实现(第二版)》清华大学出版社,胡广书)

      为了更好地理解这个问题,我又一次借用MATLAB的强大力量,写了一个简单的程序如下:

  1. %点数
  2. n=0:127;
  3. %频率
  4. f=0.1;
  5. %信号,正弦叠加矩形
  6. y1=sin(2*pi*f*n);
  7. y1(1:16)=y1(1:16)+1;
  8. %绘制y1的fft谱幅度
  9. %谱线较多,直接画的包络
  10. figure;
  11. plot(abs(fft(y1)));
  12. %对信号进行截短
  13. y2s=y1(1:32);
  14. %绘制y1截断后没有补零的fft谱幅度
  15. figure;
  16. fy2s=abs(fft(y2s));
  17. stem(fy2s);
  18. %然后补零使y1和y2一样长
  19. y2=[y2s zeros(1,128-32)];
  20. %打开一个绘图窗口
  21. figure;
  22. %绘制y1的fft谱幅度
  23. %谱线较多,直接画包络
  24. plot(abs(fft(y1)));
  25. %在同一个figure中继续绘图
  26. hold on;
  27. %绘制y2的fft谱幅度(红色)
  28. %谱线较多,直接画包络
  29. plot(abs(fft(y2)),'r');
  30. %绘制y2s的fft幅度谱
  31. stem(1:4:128,fy2s,'k');
  32. hold off;
复制代码


      程序生成的图像如下:
abc1.jpg
图1(原信号的谱,因为点数较多,绘制的是包络)
abc2.jpg
图2(截短信号的谱)
abc3.jpg
图3(比较原信号【蓝】,截短信号【黑】,补零信号【红】三者谱的关系,补零信号的谱由于点数较多,绘制的是包络)

      程序的大意是:首先生成了一个长度为128点的信号,绘制了它的DFT谱,然后将该信号截短,求其DFT谱,然后对截短的信号补零,使其长度为128点再求DFT谱,并将原信号的谱与补零信号的谱进行比较,结果一目了然。

      从图三中我们可以发现,红色的补零信号的谱,仅仅是对黑色的截短信号的谱的插值,也就是说补零信号的谱,是通过截短信号的谱进行了推测(插值算法)得来的,它并不能反映原信号的谱(因为原信号在截短的过程中部分信息丢失了,而补零并没有将这些丢失的信息找回来),所以虽然补零信号的谱线间隔变小了,但是除了从截短信号的谱中取出来的32根谱线以外,其余的96根谱线都是无效的。去掉这些无效的谱线,采样频率不变,有效的谱线数不变,所以其物理频率分辨率自然没有改变。

      在时域补零起到了对频谱的插值作用,考虑到傅立叶变换的对称性,考虑到变换与反变换的表达式除了系数和符号的差别外并没有根本的不同,推断出这样一个结果:在频域补零也会造成时域的插值(当然在频域进行操作一定要注意谱的对称性,否则反变换回来得到的将是一个复信号,至于为什么会这样参看《xialulee来说明什么是负频率》 ),于是我又写了如下的程序:  

  1. f=0.1;

  2. n=0:31;

  3. %产生正弦信号

  4. y1=sin(2*pi*f*n);

  5. %y1的fft谱

  6. f1=fft(y1);

  7. %给y1的谱补零

  8. f2=[f1(1:end/2+1) zeros(1,31) f1(end/2+1:end)];

  9. stem(abs(f2));

  10. %对补零后的谱作反变换

  11. y2=ifft(f2);

  12. %打开绘图窗口

  13. figure;

  14. %绘制多幅图形

  15. hold on;

  16. %绘制y2

  17. stem(y2);

  18. %绘制原信号(除以2以比较插值效果)

  19. stem(1:2:64,y1/2,'r');

  20. hold off;
复制代码


      程序的输出图形如下:
abc4.jpg
图4      仿佛证明了这个推测。这是不是就是所谓的傅立叶插值?

转自:http://zhichengma2007.blog.163.c ... 292009102933630710/

评分

1

查看全部评分

发表于 2015-10-22 14:18 | 显示全部楼层
FFT方法的前提是数据数量是2的整数冪,才能得到快速的效果。加零以后当然和原来的不一样,于是很多人做了不少的研究来弥补。
但是这已经是过去的事了,现在电脑的速度飞速发展,如果不是2的整数冪,Matlab自动用DFT慢速做,也不是多大的问题,所以加零问题现在已经不那么必须的了,我建议不必再在这个问题花脑筋了!
 楼主| 发表于 2015-10-23 07:45 | 显示全部楼层
simon21 发表于 2015-10-22 09:02
补零这个问题网上讨论很多,先给你转一个总结性的帖子,你先看看

1.转自http://www.ilovematlab.cn/view ...

谢谢,我好好看一下
 楼主| 发表于 2015-10-23 07:46 | 显示全部楼层
hcharlie 发表于 2015-10-22 14:18
FFT方法的前提是数据数量是2的整数冪,才能得到快速的效果。加零以后当然和原来的不一样,于是很多人做了不 ...

因为要自己编程实现相关的功能,程序要嵌入到我们的硬件系统中,不能用matlab
发表于 2015-10-26 21:40 | 显示全部楼层
对于整周期采样的一段数据,补零虽然能加快计算但会产生泄漏,只能通过加窗弥补或采用不补零的DFT。
对于非整周期采样的一段数据,反正都有泄漏,不如补零计算快点。
发表于 2015-10-26 21:43 | 显示全部楼层
“程序的大意是:首先生成了一个长度为128点的信号,绘制了它的DFT谱,然后将该信号截短,求其DFT谱,然后对截短的信号补零,使其长度为128点再求DFT谱,并将原信号的谱与补零信号的谱进行比较,结果一目了然。”----------截短原始信号来比较太妥当吧,应该是原始信号的DFT(或FFT)与原始信号补零后的FFT比较吧。
 楼主| 发表于 2015-11-7 18:32 | 显示全部楼层
TestGuru 发表于 2015-10-26 21:43
“程序的大意是:首先生成了一个长度为128点的信号,绘制了它的DFT谱,然后将该信号截短,求其DFT谱,然后 ...

谢谢指教
发表于 2015-11-8 10:48 | 显示全部楼层
本帖最后由 TestGuru 于 2015-11-8 10:51 编辑

补零虽然不能提高实际频率分辨率,但能让可能让躲在栅栏后的谱线露出来。

下图为时域1024点频域1024点hanning窗的频谱图。
NoZeroPadding.png

下图为时域1024点频域65536点(即:补零65536-1024个点)hanning窗的频谱图。
ZeroPadding.png

显然通过补零,hanning窗的真实频谱才显露出来了。
 楼主| 发表于 2015-12-7 13:58 | 显示全部楼层
TestGuru 发表于 2015-11-8 10:48
补零虽然不能提高实际频率分辨率,但能让可能让躲在栅栏后的谱线露出来。

下图为时域1024点频域1024点ha ...

很好的例子,谢谢
发表于 2015-12-10 23:53 | 显示全部楼层
TestGuru 发表于 2015-11-8 10:48
补零虽然不能提高实际频率分辨率,但能让可能让躲在栅栏后的谱线露出来。

下图为时域1024点频域1024点ha ...

请问这是在什么软件里面弄得?

点评

直接下载地址: www.multi-tech.cn/MIsetup.exe . 安装后,启动软件,然后按左上方的绿色按钮停止示波器运行,去[文件]>[打开]然后选C:\VIRTINS Multi-Instrument 3.5\wav\window\...里面有几十种窗函数的时域波形可  详情 回复 发表于 2015-12-11 17:55
发表于 2015-12-11 17:55 | 显示全部楼层
hustshiyi 发表于 2015-12-10 23:53
请问这是在什么软件里面弄得?

直接下载地址: www.multi-tech.cn/MIsetup.exe .  安装后,启动软件,然后按左上方的绿色按钮停止示波器运行,去[文件]>[打开]然后选C:\VIRTINS Multi-Instrument 3.5\wav\window\...里面有几十种窗函数的时域波形可选。
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-11-24 21:51 , Processed in 0.111596 second(s), 27 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表