声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 10866|回复: 59

[FFT] 这是我写的一个关于频谱细化的程序!希望和大家探讨!

[复制链接]
发表于 2006-8-11 21:56 | 显示全部楼层 |阅读模式

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

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

x
这是我写的一个关于频谱细化的程序!希望和大家探讨!
大家帮我看看有什么问题,谢谢!
回复
分享到:

使用道具 举报

 楼主| 发表于 2006-8-11 21:58 | 显示全部楼层

接上

clf;
N=512;
fs=5120;
n=0:(N-1);
t=n/fs;
x=cos(2*pi*110*t)-2*cos(2*pi*104*t)+3*cos(2*pi*108*t);
f0=100;%移动的频率
h=exp(-j*2*pi*f0*t);
y=x.*h;%数字变频,把w=100点移动到原点
y1=resample(y,1,10);%重新采样,采样频率为fs/N
k1=0:1:511;
h1=exp((j*2*pi*f0*k1)/512);
y1=y1.*h1;%反移频,使得前后频率一致
y1=abs(fft(y1));%快速傅里叶变换
subplot(2,2,1);plot(k1,y1);
axis([100,120,0,1000]);
title('细化10倍后的频率特性');
xlabel('频率');
ylabel('幅度值');

[ 本帖最后由 VibInfo 于 2006-8-13 06:57 编辑 ]
发表于 2006-8-13 07:01 | 显示全部楼层
我这里没有装matlab看不了,不过有个程序你可以参考一下

  1. % 三种频率成分:20Hz,20.5Hz,40Hz
  2. % 采样频率Fs=100Hz (满足NiQust定理)

  3. % 要想分辨出信号的所有频谱成分,则分辨率为 0.5Hz
  4. % 若采样率为Fs=100;最少点数N应满足: Fs/N=0.5 N 最少为200,才能满足0.5Hz的分辨率要求

  5. % N<200 则无法分辨 20与20.5

  6. N=256; %分别用128,210,256,512 ,比较频谱形状,比较分辨率,比较方差(即平滑程度)
  7. % 结论 : 增大N可以提高分辨率
  8. Fs=100;
  9. n=0:1:N-1;
  10. f1=20;
  11. f2=20.5;
  12. f3=40;
  13. %xn=sin(2*pi*f1*n/Fs)+sin(2*pi*f2*n/Fs)+sin(2*pi*f3*n/Fs);多个信号
  14. xn=sin(2*pi*f1*n/Fs);%单个信号
  15. subplot(4,1,1);
  16. plot(n,xn);

  17. Xk=fft(xn);
  18. AXk=abs(Xk(1:N/2));
  19. k=(0:N/2-1)*Fs/N;
  20. subplot(4,1,2);
  21. stem(k,AXk); %可分别用plot和stem

  22. % 补零 先分别另上面N=256和N=200,N=128, 补零到下面M=512,M=1024观察 谱线变化(细化),分辨率(不提高)

  23. M=512;
  24. xn2=[xn zeros(1,M-N)];
  25. Xk2=fft(xn2);
  26. AXk2=abs(Xk2(1:M/2));
  27. m=0:1:M-1;
  28. k2=(0:M/2-1)*Fs/M;
  29. subplot(4,1,3);
  30. plot(m,xn2);
  31. subplot(4,1,4);
  32. stem(k2,AXk2);%可分别用plot和stem
复制代码

评分

1

查看全部评分

发表于 2006-8-13 08:33 | 显示全部楼层
楼主在编ZFFT时,最大的错误是x用了512个数。x的512个数经频率移动和下采样(10:1)后只留下52个点,而你在程序中继续作512点处理(y1=y1.*h1)便产生错误 。我把它改为如下:
clf;
N=512;
fs=5120;
n=0:(N*10-1);
t=n/fs;
x=cos(2*pi*110*t)-2*cos(2*pi*104*t)+3*cos(2*pi*108*t);
f0=100;%移动的频率
h=exp(-j*2*pi*f0*t);
y=x.*h;%数字变频,把w=100点移动到原点
y1=resample(y,1,10);%重新采样,采样频率为fs/N
k1=0:1:511;
h1=exp((j*2*pi*f0*k1)/512);
y1=y1.*h1;%反移频,使得前后频率一致
y1=abs(fft(y1));%快速傅里叶变换
plot(k1,y1);
axis([100,120,0,1000]);
set(gca, 'XTickMode', 'manual', 'XTick', [100, 104, 108, 110, 115, 120]); grid;
title('细化10倍后的频率特性');
xlabel('频率');
ylabel('幅度值');
得图有(达到细化目的)
luo2a.jpg

评分

1

查看全部评分

 楼主| 发表于 2006-8-13 08:56 | 显示全部楼层
非常感谢,我试试看能不能运行
发表于 2006-8-13 08:57 | 显示全部楼层
一篇高怀钢“一种分析频谱局部特性的快速算法” ,其中详细介绍了ZFFT的方法和实现步骤。
发表于 2006-9-3 12:14 | 显示全部楼层
h1=exp((j*2*pi*f0*k1)/512);
y1=y1.*h1;%反移频,使得前后频率一致
我不明白这个地方是什么意思,能不能给我解释一下 啊
发表于 2006-9-3 13:28 | 显示全部楼层
频率细化的方法,提高了频率分辨率吗?好象没有啊,开始采样点数5120,采样率5120,分辨率为1。细化后采样率512,采样点数512,还是1啊?
发表于 2006-9-4 13:53 | 显示全部楼层
原帖由 happytaotao 于 2006-9-3 13:28 发表
频率细化的方法,提高了频率分辨率吗?好象没有啊,开始采样点数5120,采样率5120,分辨率为1。细化后采样率512,采样点数512,还是1啊?

分辨率的提高,是对初始时按512点作FFT来讲的,而不是对5120点。这样在同样点数的FFT下,分辨率提高了10倍。
发表于 2006-9-4 13:59 | 显示全部楼层
原帖由 happytaotao 于 2006-9-3 12:14 发表
h1=exp((j*2*pi*f0*k1)/512);
y1=y1.*h1;%反移频,使得前后频率一致
我不明白这个地方是什么意思,能不能给我解释一下 啊

在程序的前部,对信号作了移频:
h=exp(-j*2*pi*f0*t);
y=x.*h;%数字变频,把w=100点移动到原点
把100Hz频率移到原点。而在下采样以后,楼主又经上反移频,把100Hz频率移回原处。
发表于 2006-9-5 14:18 | 显示全部楼层

回复 #10 songzy41 的帖子

我想我还是不明白,说把100Hz频率移回原处是在细化完后作的工作,还是这也是细化的一部分呢?可是我看到的运行结果,100Hz的频谱成分还是在原点啊
发表于 2006-9-5 15:46 | 显示全部楼层
原帖由 happytaotao 于 2006-9-5 14:18 发表
我想我还是不明白,说把100Hz频率移回原处是在细化完后作的工作,还是这也是细化的一部分呢?可是我看到的运行结果,100Hz的频谱成分还是在原点啊

你为什么会认为“100Hz的频谱成分还是在原点”,是否因为作的图,看上去100Hz在原点?这是因为用了
axis([100,120,0,1000]);
只显示100~120Hz之间的频谱,使看得更清楚。而实际原点还是0Hz。
这并不是一定要这样做的。你不是看过高怀钢的“一种分析频谱局部特性的快速算法”一文,他给出是ZFFT的通用算法:频移(调制)-滤波-下采样-再FFT,这样一个过程。

评分

1

查看全部评分

发表于 2006-9-5 15:58 | 显示全部楼层
不错!
发表于 2006-9-6 09:16 | 显示全部楼层

谢谢

你好:songzy41
为什么你写的程序我在 matlab中运行不了呢?
错误的原因为:??? Error using ==> fft
Not enough input arguments.
是什么意思?
发表于 2006-9-6 14:14 | 显示全部楼层
这错误表示调用fft函数中的变量有问题。这程序我在Matlab v6.5中通过没有问题,不知对不同的Matlab 版本是否有影响。
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

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

GMT+8, 2024-5-19 18:56 , Processed in 0.073548 second(s), 22 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

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