声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

楼主: luoluo

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

[复制链接]
发表于 2008-6-14 09:06 | 显示全部楼层
呵呵,谢谢你!以后还要麻烦你!
这本书在那里可以直接下载,在线阅读吗?
呵呵,我搜到它的程序了.在
http://forum.vibunion.com/thread-23569-1-48.html
可是我想搞清楚它的原理.
回复 支持 反对
分享到:

使用道具 举报

发表于 2008-6-14 18:43 | 显示全部楼层
我对程序进行了改进,如下,有一个小问题,低通滤波后频域图形显示有错,这是怎么回事?怎样改进啊?谢谢了。
clear all; clc; close all;
fs=input('采样频率fs=');%fs的取值需满足采样定律
f0=input('需求的中心频谱f0=');
B=input('理想低通滤波器的带宽B=');
N=fs+1;
t=linspace(0,1,N);dt=1/N-1;%取信号长度Tp=1s,给出采样周期Ts=dt,采样频率fs=1/Ts=N-1
N1=(fs/B)*(N-1)+1;
t1=(1:N1)/fs;
x=10*sin(2*pi*64*t1)+10*sin(2*pi*250*t1)+20*sin(2*pi*256*t1)+30*sin(2*pi*258*t1)+20*sin(2*pi*500*t1);
%信号最高频率fc=512HZ
figure(1)
stem(t1(1:1024),x(1:1024));grid%画出原始信号
title('原始信号x(t1)')
xlabel('t1(s)');ylabel('x(t1)');

f=linspace(0,N-1,N);
X=fft(x,N);
figure(2)
plot(f,abs(X)/max(abs(X))),grid
xlabel('f(HZ)');ylabel('|X(jf)|');

x1=x.*exp(-j*2*pi*t1*f0);%'.*'为元素群运算
figure(3)
stem(t1(1:1024),x1(1:1024));grid
title('频移后的信号x1(t1)')
xlabel('t1(s)');ylabel('x1(t1)');
  
X1=fft(x1,N);
f11=(0:N-1)*fs/N+f0;
figure(4)
plot(f11,abs(X1)/max(abs(X1))),grid
xlabel('f11(HZ)');ylabel('|X1(jf11)|');
set(gca, 'XTickMode', 'manual', 'XTick', [f0, 250, 400, 500, 524, 762, 1024]);


XX1=fft(x1);
ff=(0:N1-1)*fs/N1;%频域取样点数N1
n=find(ff<B/2);%求出底通滤波器带宽内的下标
f1=f(n);%取出中段正频率
X2=XX1(n);%取出中段正频谱
n1=find(f<B/2);%求出底通滤波器带宽内的下标
f111=f(n1);%取出中段正频率
X21=X1(n1);%取出中段正频谱
figure(6) figure6.fig (5.64 KB, 下载次数: 18)
plot(f111,abs(X21)/max(abs(X21))),grid


x2=ifft(X2,N1);
figure(7)
stem(t1(1:1024),x2(1:1024));grid

y=resample(x2,1,fs/B);% 对信号进行fs/B倍的抽取
t2=resample(t1,1,fs/B);
fs1=B;%采样频率下降fs/B倍
f3=(0:N-1)*fs1/N+f0;
Y=fft(y,N);% 对信号进行抽取后的 FFT 运算
figure(8)
plot(f3,abs(Y)/max(abs(Y))),grid
xlabel('f(HZ)');ylabel('|Y(jf)|');
set(gca, 'XTickMode', 'manual', 'XTick', [f0, 250, 256, 258, 260, 270, 280, 300]);
%输入参数分别为  1024  240  51.2
发表于 2008-6-14 20:04 | 显示全部楼层

所谓的错误是否x轴的刻度不对?只要改为:
plot(f111+240,abs(X21)/max(abs(X21))),grid
就可以了。

如果想了解ZFFT的原理,我找到一些文献名,可参考一下:
refrence.JPG
发表于 2008-6-15 07:56 | 显示全部楼层
我想楼主是没有明白我的意思,移频后的频谱如图, figure4.fig (17.49 KB, 下载次数: 12) 经过低通滤波后,频谱显示应该是该图250附近处的谱线,其它谱线被滤除,即只应显示[f0,f0+B/2]之间频谱。
我想的是,在图7中,只显示图4中f0+B/2之前频谱。
用plot(f111+f0,abs(X21)/max(abs(X21))),grid
不是我想要的效果。
呵呵,不过还是要谢谢你。
发表于 2009-2-17 19:23 | 显示全部楼层

回复 9楼 songzy41 的帖子

小弟刚学ZOOMFFT,不知道说的对不对,我觉得在你的程序运行到x=cos(2*pi*110*t)-2*cos(2*pi*104*t)+3*cos(2*pi*108*t);这一行后,直接对x进行FFT,频率分辨率也是1啊,你下面细化后频率分辨率也是1啊,这点不太明白
发表于 2009-3-11 16:08 | 显示全部楼层
为什么看不到,怎么回事????
发表于 2009-3-11 17:11 | 显示全部楼层
感谢,先试试,正好要用!
发表于 2009-12-15 19:52 | 显示全部楼层
学习了。。。。。。。。。。。。。。。
发表于 2010-9-25 10:28 | 显示全部楼层
我想问一下,我现在的数据时512点,我做的是512 的FFT,按照前面各位的讨论,我根本不可能细化到十倍吧?因为降采样后的数据只有51个,做FFT时候,肯定是要补零的,也就没有所谓的细化吧?
发表于 2010-9-25 20:40 | 显示全部楼层
楼上各位的程序,为何我在运行中调试程序是出错的,而且打开的xdate.txt里面的数据应该是楼主自己给定的吧?
发表于 2010-9-26 09:57 | 显示全部楼层
suqiao 发表于 2010-9-25 20:40
楼上各位的程序,为何我在运行中调试程序是出错的,而且打开的xdate.txt里面的数据应该是楼主自己给定的吧? ...

请说明错误信息,以便大家分析
xdate.txt应该是你所有处理的时间序列
发表于 2010-9-26 10:18 | 显示全部楼层
我想问一下,我采样后的数据是512点,我做的是512 的FFT,按照前面各位的讨论,我根本不可能细化到十倍吧?因为降采样后的数据只有51个,做FFT时候,肯定是要补零的,也就没有所谓的细化吧?
发表于 2012-7-14 09:35 | 显示全部楼层
用matlab做这些东西确实方便不少!
发表于 2013-5-10 13:38 | 显示全部楼层
hnyanhua 发表于 2008-6-9 11:22
我想楼主是看过高怀钢的一文,在文中介绍ZFFT实现时,笫1步是频移,笫2步是滤波。在重采样之前为了防止混迭 ...

您说的是resample 函数不带滤波,这点我并不是赞同
resample采用多相滤波器对时间序列进行重采样,得到的序列y的长度为原来的序列x的长度的p/q倍,p和q都为正整数。此时,默认地采用使用FIR方法设计的抗混叠的低通滤波器。
发表于 2014-8-23 14:28 | 显示全部楼层
交流学习助人为乐
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

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

GMT+8, 2024-5-20 06:32 , Processed in 0.059005 second(s), 20 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

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