mengfuse 发表于 2012-8-28 18:59

平滑伪Wigner-Ville分布计算程序

% Wigner-Ville分布计算程序
function=SWWignerVille(Sig,SampFreq,FreqBins,GLen,HLen)
%Sig      : 输入信号
%FreqBins : 频率轴划分区间数(默认值:512)
%SampFreq : 信号的采集频率
%alfa   : 窗函数控制参数
%GLen   : 窗函数g的长度(默认值:FreqBins/5)
%HLen   : 窗函数h的长度(默认值:FreqBins/4)
%WVD      :   SWWigner-Ville计算结果

if (nargin<1),
    error('At least 1 input required');
end;
% Sig=hilbert(real(Sig));    % 计算信号的解析信号
SigLen=length(Sig);      % 获取信号的长度
if (nargin<2),
    SampFreq=1;
elseif (nargin<3),
    FreqBins=512;
elseif (nargin<4),
    GLen=floor(FreqBins/5);   
    HLen=floor(FreqBins/4);
end;
FreqBins=2^nextpow2(FreqBins);
FreqBins=min(FreqBins,SigLen);
GLen=GLen+1-rem(GLen,2);
HLen=HLen+1-rem(HLen,2);
GWin=window(@gausswin,GLen);
HWin=window(@gausswin,HLen);
Lg=(GLen-1)/2;
Lh=(HLen-1)/2;
HWin=HWin/HWin(Lh+1);
WVD=zeros(FreqBins,SigLen);
wait=waitbar(0,'Under calculation,please wait...');

for kk=1:SigLen,
    waitbar(kk/SigLen,wait);
    MTau=min();
    k=-min():min();
    SubG=GWin(Lg+1+k);
    SubG=SubG/sum(SubG);
    WVD(1,kk)=sum(SubG.*Sig(kk-k,1).*conj(Sig(kk-k)));%
    for tau=1:MTau,
      k=-min():min();
      SubG=GWin(Lg+1+k);
      SubG=SubG/sum(SubG);
      R=sum(SubG.*Sig(kk+tau-k,1).*conj(Sig(kk-tau-k)));
      WVD(1+tau,kk)=HWin(Lh+tau+1)*R;
      R=sum(SubG.*Sig(kk-tau-k,1).*conj(Sig(kk+tau-k)));
      WVD(FreqBins+1-tau,kk)=HWin(Lh-tau+1)*R;
    end;   
end;
close(wait);
WVD=fft(WVD);
f=linspace(0,0.5,FreqBins)*SampFreq;
t=(0:(SigLen-1))/SampFreq;
set(gcf,'position',);
set(gcf,'color','w');
axes('position',);
mesh(t,f,abs(WVD));
axis();
colorbar;
xlabel('t/s');
ylabel('f/Hz');
title('平滑伪Wigner-Ville分布');

输入以下:
%SampFreq=100;
%FreqBins=512;
%t=0:1/SampFreq:4;
%Sig=sin(20*pi*t)+sin(80*pi*t);
%WignerVille(Sig,SampFreq,FreqBins,Freqbins/5,FreqBins/4);

出现以下错误:
??? Index exceeds matrix dimensions.

Error in ==> SWWignerVille at 43
    WVD(Temp,kk)=sum(SubG.*Sig(kk-k,1).*conj(Sig(kk-k)));

谁能给修改下。谢谢

mengfuse 发表于 2012-9-6 08:59

{:{28}:}{:{28}:}

犟牛 发表于 2012-9-6 10:31

WVD(1,kk)=sum(SubG.*Sig(kk-k,1).*conj(Sig(kk-k)));改为:WVD(1,kk)=sum(SubG'.*Sig(kk-k).*conj(Sig(kk-k)));后面还有类似错误,自己改吧

ChaChing 发表于 2012-9-6 22:15

Ref: 12F, 常见的程序出错问题整理 http://forum.vibunion.com/thread-46001-1-1.html
from http://forum.vibunion.com/home-space-uid-63979-do-blog-id-18250.html

tf19830910 发表于 2012-10-20 16:54

想请教楼主一下,这个平滑算法是根据什么文献编写的?时频工具箱中的算法似乎和你是一样的,但是我找不到这个算法的理论文献。
页: [1]
查看完整版本: 平滑伪Wigner-Ville分布计算程序