有关短时傅立叶的一个程序,帮我看看问题出在哪

wyb1552

我想通过短时傅立叶算法得到y函数的有效值曲线,但画出来的图形还多出了两堆曲线,不知道是什么回事,请高手们帮我看看


%采样波形
clear all
fs=3200;     %采样频率
t=[0:1/fs:0.1-1/fs];
a=1.414*ones(size(t));
index=find(0.035<=t&t<=0.075);
a(index)=1.414*0.7;
y=a.*sin(100*pi*t);
subplot(211);
plot(t,y);
xlabel('时间');
ylabel('幅值');
title('电压暂降信号 时域波形');

%短时傅立叶变换
N=length(y);
Nw=64;                     %窗函数长
L=16;                       %窗函数每次移动的样点数
Tn=(N-Nw)/L+1;              %计算把数据x共分成多少段
nfft=64;                   %FFT的长度
TF=zeros(Tn,nfft);          %将存放三维谱图，先清零
for i=1:Tn
    n1=L*(i-1)+1;
    n2=L*(i-1)+Nw;
    sf=fft(y(n1:n2),nfft);     %FFT变换
    TF(i,:)=sf;                %把谱图存放在TF中
end
subplot(212);
fnew=((1:nfft)-nfft/2)*fs/nfft;
tnew=(1:Tn)*L/fs;
[F,T]=meshgrid(fnew,tnew);
plot(T,2*abs(TF)/(sqrt(2)*Nw),'b');
xlabel('时间');
ylabel('有效值');
title('短时傅立叶跟踪电压暂降信号的有效值波形');
grid on;

[ 本帖最后由 wyb1552 于 2008-5-30 17:52 编辑 ]

songzy41

原帖由 wyb1552 于 2008-5-30 17:40 发表

                               
登录/注册后可看大图

我想通过短时傅立叶算法得到y函数的有效值曲线,但画出来的图形还多出了两堆曲线,不知道是什么回事,请高手们帮我看看
...
[F,T]=meshgrid(fnew,tnew);
plot(T,2*abs(TF)/(sqrt(2)*Nw),'b');

TF是一个17*64的数组，不能用plot命令来作图，同时T也是一个17*64的数组。把plot命令改一下，改为：
plot(tnew,2*abs(TF(:,2))/(sqrt(2)*Nw),'b');
就能得到下图。

wyb1552

谢谢您帮我解决这个问题,我知道该怎么做了

wyb1552

我想在问下,短时傅立叶是不是只能对一个周波积分?如何提高它的精度?

songzy41

短时傅立叶变换的基础还是傅立叶变换，每一帧的谱值相当于是该时间窗（在楼主的程序中取了一周）内的平均，所以要想用短时傅立叶变换来提高跟踪信号变化的精度是不可能的，建议楼主采用其它方法试试。作为傅立叶变换，窗长内的信号至少要一个周波，最好是多几个周波。

wyb1552

哦 ,好的,谢谢你的帮忙,学到了不少

和丹

学习了，很好的练习

和丹

你的短时傅立叶为啥不加窗函数？矩形窗......