短时傅里叶变换matlab程序

合集下载
  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

function [Spec,Freq]=STFT(Sig,nLevel,WinLen,SampFreq)

%计算离散信号的短时傅里叶变换;

% Sig 待分析信号;

% nLevel 频率轴长度划分(默认值512);

% WinLen 汉宁窗长度(默认值64);

% SampFreq 信号的采样频率(默认值1);

if (nargin <1),

error('At least one parameter required!');

end;

Sig=real(Sig);

SigLen=length(Sig);

if (nargin <4),

SampFreq=1;

end

if (nargin <3),

WinLen=64;

end

if (nargin <2),

nLevel=513;

end

nLevel=ceil(nLevel/2)*2+1;

WinLen=ceil(WinLen/2)*2+1;

WinFun=exp(-6*linspace(-1,1,WinLen).^2);

WinFun=WinFun/norm(WinFun);

Lh=(WinLen-1)/2;

Ln=(nLevel-1)/2;

Spec=zeros(nLevel,SigLen);

wait=waitbar(0,'Under calculation,please wait...');

for iLoop=1:SigLen,

waitbar(iLoop/SigLen,wait);

iLeft=min([iLoop-1,Lh,Ln]);

iRight=min([SigLen-iLoop,Lh,Ln]);

iIndex=-iLeft:iRight;

iIndex1=iIndex+iLoop;

iIndex2=iIndex+Lh+1;

Index=iIndex+Ln+1;

Spec(Index,iLoop)=Sig(iIndex1).*conj(WinFun(iIndex2)); end;

close(wait);

Spec=fft(Spec);

Spec=abs(Spec(1:(end-1)/2,:));

Freq=linspace(0,0.5,(nLevel-1)/2)*SampFreq; t=(0:(SigLen-1))/SampFreq;

clf

set(gcf,'Position',[20 100 500 430]);

set(gcf,'Color','w');

axes('Position',[0.1 0.45 0.53 0.5]);

mesh(t,Freq,Spec);

axis([min(t) max(t) 0 max(Freq)]);

colorbar

xlabel('t/s');

ylabel('f/Hz');

title('STFT时频谱图');

axes('Position',[0.1 0.1 0.55 0.25]);

plot(t,Sig);

axis tight

ylabel('x(t)');

title('时域波形');

axes('Position',[0.73 0.45 0.24 0.5]);

PSP=abs(fft(Sig));

Freq=linspace(0,1,SigLen)*SampFreq;

plot(PSP(1:end/2),Freq(1:end/2));

title('频谱');

相关文档
最新文档