声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 3167|回复: 7

[FFT] 带扩展系数的zfft

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

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

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

x
一般来说,由于滤波器的过渡带的问题,ZFFT结果的两端幅值都会减小.
经过扩展后的ZFFT可解决这一问题.

  1. function [xz fz]=zoomffta(x,fs,N,fe,D,a)
  2. % 带扩展系数a的zoomfft
  3. % x---输入时间序列,注意x的长度不能小于(N-1)*D+2*M;
  4. % fs--x的采样频率;
  5. % N---做谱点数;
  6. % D---细化倍数;
  7. % a---扩展系数。a取大于0小于等于1的数,越大M越小,计算量越小

  8. M=4*D/a;                        %滤波器半阶数
  9. k=1:M;                          
  10. w=0.5+0.5*cos(pi*k/M);          %Hanning窗

  11. % 求取带通滤波器上下界;
  12. fl=max(fe-fs/(4*D),-fs/2);
  13. fh=min(fe+fs/(4*D),fs/2);

  14. %求取扩展带通滤波器上下界;
  15. hfl=fl-(fh-fl)*a/2;
  16. hfh=fh+(fh-fl)*a/2;

  17. %构造扩展带通滤波器;
  18. wl=2*pi*hfl/fs;
  19. wh=2*pi*hfh/fs;
  20. hr(1)=(wh-wl)/pi;
  21. hr(2:M+1)=(sin(wh*k)-sin(wl*k))./(pi*k).*w;
  22. hi(1)=0;
  23. hi(2:M+1)=(cos(wl*k)-cos(wh*k))./(pi*k).*w;

  24. %选抽滤波
  25. for k=1:N
  26.     kk=(k-1)*D+M;
  27.     xrz(k)=x(kk+1)*hr(1)+sum(hr(2:M+1).*(x(kk+2:kk+M+1)+x(kk:-1:kk-M+1)));
  28.     xiz(k)=x(kk+1)*hi(1)+sum(hi(2:M+1).*(-x(kk+2:kk+M+1)+x(kk:-1:kk-M+1)));
  29. end

  30. %移频,把fl移到0频
  31. yf=D*fl/fs;                       %移频量
  32. xz=(xrz+j*xiz).*exp(-j*2*pi*(0:N-1)*yf);

  33. xz=fft(xz);
  34. xz=xz(1:N/2)/N;                  %细化复数谱
  35. fz=(0:N/2-1)*fs/N/D+fl;          %细化谱对应的频率
复制代码

[ 本帖最后由 MVH 于 2006-10-4 15:51 编辑 ]

评分

1

查看全部评分

回复
分享到:

使用道具 举报

发表于 2007-4-12 23:23 | 显示全部楼层
%构造扩展带通滤波器;
wl=2*pi*hfl/fs;
wh=2*pi*hfh/fs;
hr(1)=(wh-wl)/pi;
hr(2:M+1)=(sin(wh*k)-sin(wl*k))./(pi*k).*w;
hi(1)=0;
hi(2:M+1)=(cos(wl*k)-cos(wh*k))./(pi*k).*w;

%选抽滤波
for k=1:N
    kk=(k-1)*D+M;
    xrz(k)=x(kk+1)*hr(1)+sum(hr(2:M+1).*(x(kk+2:kk+M+1)+x(kk:-1:kk-M+1)));
    xiz(k)=x(kk+1)*hi(1)+sum(hi(2:M+1).*(-x(kk+2:kk+M+1)+x(kk:-1:kk-M+1)));
end


请问hr(1)和hi(1)为什么要单独求呢?
for循环里是做卷积吗?能否稍微解释具体点。
谢谢
 楼主| 发表于 2007-4-13 11:31 | 显示全部楼层
单独求是因为滤波器在0点是0/0的形式,应该用洛比塔法则来求.
for 里可以说是做卷积,滤波就是一个卷积过程.
发表于 2007-5-15 21:13 | 显示全部楼层
%选抽滤波
for k=1:N
    kk=(k-1)*D+M;
    xrz(k)=x(kk+1)*hr(1)+sum(hr(2:M+1).*(x(kk+2:kk+M+1)+x(kk:-1:kk-M+1)));
    xiz(k)=x(kk+1)*hi(1)+sum(hi(2:M+1).*(-x(kk+2:kk+M+1)+x(kk:-1:kk-M+1)));
end

这几句是怎样实现抽取和滤波,实在是看不懂了,能否解释的详细一些?
谢谢!
 楼主| 发表于 2007-5-16 08:50 | 显示全部楼层
本帖最后由 wdhd 于 2016-9-14 10:31 编辑

  原帖由 yixiaofan 于 2007-5-15 21:13 发表
  %选抽滤波
  for k=1:N
  kk=(k-1)*D+M; %从第M点开始隔D点进行滤波
  xrz(k)=x(kk+1)*hr(1)+sum(hr(2:M+1).*(x(kk+2:kk+M+1)+x(kk:-1:kk-M+1))); %滤波即做一个卷积,卷积对每一个点来说又是对滤波器进行反转再做内积运算,具体参考卷积的步骤;这里考滤到滤波器实部的偶对称,只用了一半来表示滤波器的实部;
  xiz(k)=x(kk+1)*hi(1)+sum(hi(2:M+1).*(-x(kk+2:kk+M+1)+x(kk:-1:kk-M+1))); %虚部同样,只不过滤波器虚部是奇对称的
  e ...
发表于 2008-4-1 22:58 | 显示全部楼层

回复 5楼 的帖子

%选抽滤波
for k=1:N
    kk=(k-1)*D+M;
    xrz(k)=x(kk+1)*hr(1)+sum(hr(2:M+1).*(x(kk+2:kk+M+1)+x(kk:-1:kk-M+1)));
    xiz(k)=x(kk+1)*hi(1)+sum(hi(2:M+1).*(-x(kk+2:kk+M+1)+x(kk:-1:kk-M+1)));
end
为什么报错啊?
发表于 2008-4-2 11:37 | 显示全部楼层

回复 6楼 的帖子

是我自己搞错了,导入数据出的错。
发表于 2009-3-6 13:11 | 显示全部楼层
把抽取滤波和变频的顺序倒过来不会有什么问题吗?我看的很多的论文上面都是先变频后滤波的!这样是不是对滤波来说,要简单一些?先变频把关心的频率终点移到零频,这样应该需要的是两个低通滤波器而不是带通滤波器!如果先抽取滤波需要的是带通滤波器~

[ 本帖最后由 ruben052013 于 2009-3-6 15:28 编辑 ]
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

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

GMT+8, 2024-11-11 17:37 , Processed in 0.078750 second(s), 22 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

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