tfrview.m

来自「时频分析工具箱」· M 代码 · 共 571 行 · 第 1/2 页

  else
   freqs=freqr;
   spec=spec(1:Nf2,:);
  end
  Maxsp=max(max(spec));
  if linlogspec==0,
   plot(freqs,spec);
   title('Linear scale');  
   set(axspec,'ytick',[0 round(Maxsp/2) round(Maxsp)]);
   set(gca,'Ylim',[0 Maxsp*1.2]);
  else
   plot(freqs,10*log10(spec/Maxsp)); 
   title('Log. scale [dB]');
   set(gca,'ytick',[fix(10^(log10(Maxsp)/2));fix(Maxsp)]);
   str1=['   ']; str2=['   '];
   st1=num2str(fix(10*log10(Maxsp)/2)); str1(1:length(st1))=st1;
   st2=num2str(fix(10*log10(Maxsp)));   str2(1:length(st2))=st2; 
   set(axspec,['Y' TickLabelStr],[str1;str2]);
   set(gca,'Ylim',[10*log10(min(min(spec))/Maxsp) 0]);
  end
  xlabel('Energy spectral density');
  Nsp=length(spec); 
  set(gca,'Xlim',[freqs(1) freqs(Nsp)]);
  set(gca,['X' TickLabelStr],[],'view',[-90 90]);

 end;

  if iscolorbar,                     % Is there a color bar ?
   PositionTfr=get(axtfr,'Position');
   set(axtfr,'Position',PositionTfr-[0 0 0.03 0]);
   set(axcb, 'Position',[PositionTfr(1)+PositionTfr(3)-0.01,...
                         PositionTfr(2) 0.01 PositionTfr(4)]);
   axes(axcb); 
   Ncolors=length(colormap); 
   [cmin,cmax]=caxis; colorvect=linspace(cmax,cmin,Ncolors);
   imagesc(colorvect'); axis('off'); 
   if issig,                          % there is a signal
    PositionSig=get(axsig,'Position');
    set(axsig,'Position',PositionSig-[0 0 0.03 0]);
   end 
  end

 axes(axtfr);                         % Display the tfr
 if (tcol==1),                                          % 时间长度为1
  plot(freqr,tfr(1:Nf2));
 else                                                   % 时间长度为tcol,来自t(trow,tcol)
  if strcmp(computer,'MAC'),
   tfr=flipud(tfr(1:Nf2,:));
  else
   tfr=tfr(1:Nf2,:);                                    % *** 关键之处 ***，若时长大于1，则按需选择正频率对应的tfr值。
  end;
 
  indmin=find(tfr<mini);                                % 找tfr中比最小值mini小的元素索引，将其置最小值。在自己编成时，需引用。
  tfr(indmin)=mini*ones(1,length(indmin));
 
  indmax=find(tfr>maxi);                                % 找tfr中比最大值mani大的元素索引，将其置最大值。
  tfr(indmax)=maxi*ones(1,length(indmax));
  if isaffine & (display==2),
   fprintf('imagesc does not support non-linear scales for axes. Replaced by pcolor.\n');
   display=3; 
  end;

        % 主要的contour，imagesc，pcolor，surf和mesh显示工作。
  if display==1,                  % contour
   contour(ts,freqr,tfr,levels);  % contour(tfr,levels,ts,freqr);       % ts：时间轴刻度（如：1:512），freqr：频率轴刻度（如：0:0.5，不是间隔1）
   DisplayStr=', contour';
  elseif display==2,              % imagesc
   if linlogtfr==0,
    imagesc(ts,freqr,tfr); axis('xy');
   else
    imagesc(ts,freqr,log10(tfr));axis('xy');
   end
    DisplayStr=', imagesc';

  elseif display==3,              % pcolor
   if linlogtfr==0,
    pcolor(ts,freqr,tfr); shading interp;
   else
    pcolor(ts,freqr,log10(tfr)); shading interp;
   end
   DisplayStr=', pcolor';

  elseif display==4,              % surf
   if linlogtfr==0,
    surf(ts,freqr,tfr); shading interp;
    zlabel('Amplitude');
    axis([ts(1) ts(tcol) freqr(1) freqr(Nf2) mini maxi]);
   else
    surf(ts,freqr,log10(tfr)); shading interp;
    zlabel('Positive values');
    axis([ts(1) ts(tcol) freqr(1) freqr(Nf2) log10(mini) log10(maxi)]);
   end
   DisplayStr=', surf';

  elseif display==5,              % mesh
   if linlogtfr==0,
    mesh(ts,freqr,tfr);  shading interp;
    zlabel('Amplitude');
    axis([ts(1) ts(tcol) freqr(1) freqr(Nf2) mini maxi]);
   else
    mesh(ts,freqr,log10(tfr));shading interp;
    zlabel('Positive values');
    axis([ts(1) ts(tcol) freqr(1) freqr(Nf2) log10(mini) log10(maxi)]);
   end
   DisplayStr=', mesh';
  end
 
 % Define the title and check the input arguments depending on 'method'
    % 拼接在图上方显示的分布名称，刻度类型，绘图方式，tfr最小值的门限
 method=method(4:length(method));
 if nargin==5, % if there is no additional parameters, do the best.
  title([method, LinLogStr,DisplayStr,...
         ', Threshold=',num2str(threshold),'%']);

 elseif strcmp(method,'WV'  ) | strcmp(method,'MH') | ...
        strcmp(method,'PAGE'), % no parameters
  title([method,', Nf=',num2str(Nf2), LinLogStr,DisplayStr,...
        ', Threshold=',num2str(threshold),'%']);

 elseif strcmp(method,'PWV'  )|strcmp(method,'PMH'  )| ...
        strcmp(method,'SP'   )|strcmp(method,'PPAGE')| ...
        strcmp(method,'RSP'  )|strcmp(method,'RPPAG')| ...
        strcmp(method,'RPWV' )|strcmp(method,'RPMH' ),
  h=p1;[hrow,hcol]=size(h); Lh=(hrow-1)/2; % one parameter
  if (hcol~=1)|(rem(hrow,2)==0),
   error('h must be a smoothing window with odd length'); end;
  title([method, ', Lh=',num2str(Lh), ', Nf=',num2str(Nf2),...
        LinLogStr, DisplayStr,', Threshold=',num2str(threshold),'%']);
 
 elseif strcmp(method,'STFT'), % short-time fourier transform case
  h=p1;[hrow,hcol]=size(h); Lh=(hrow-1)/2;
  if (hcol~=1)|(rem(hrow,2)==0),
   error('h must be a smoothing window with odd length'); end;
  title(['|',method, '|^2, Lh=',num2str(Lh),...
        ', Nf=',num2str(Nf2), LinLogStr, DisplayStr,', Thld=',...
        num2str(threshold),'%']);
 
 elseif strcmp(method,'SPWV' ) | strcmp(method,'MHS'  )| ...
        strcmp(method,'RSPWV') | strcmp(method,'RMHS' )| ...
        strcmp(method,'ZAM'  ) | strcmp(method,'RIDBN')|...
        strcmp(method,'BJ'   ) | strcmp(method,'RIDB' )| ...
        strcmp(method,'RIDH' ) | strcmp(method,'RIDT' ),
  g=p1; [grow,gcol]=size(g); Lg=(grow-1)/2;
  if (gcol~=1)|(rem(grow,2)==0),
   error('g must be a smoothing window with odd length'); end;
  h=p2; [hrow,hcol]=size(h); Lh=(hrow-1)/2; 
  if (hcol~=1)|(rem(hrow,2)==0),
   error('h must be a smoothing window with odd length'); end;
  title([method,', Lg=',num2str(Lg),', Lh=',num2str(Lh),...
         ', Nf=',num2str(Nf2),LinLogStr, DisplayStr,...
         ', Threshold=',num2str(threshold),'%']);
 
 elseif strcmp(method,'MMCE'),
  h=p1;[hrow,hcol]=size(h); Lh=(hrow-1)/2;
  if (rem(hrow,2)==0),
   error('h must be a smoothing window with odd length'); end;
  title([method, ', Lh=',num2str(Lh), ', Nf=',num2str(Nf2),...
        LinLogStr, DisplayStr,', Threshold=',num2str(threshold),'%']);
 
 elseif strcmp(method,'CW' ) | strcmp(method,'BUD'),
  g=p1; [grow,gcol]=size(g); Lg=(grow-1)/2;
  if (gcol~=1)|(rem(grow,2)==0),
   error('g must be a smoothing window with odd length'); end;
  h=p2; [hrow,hcol]=size(h); Lh=(hrow-1)/2; 
  if (hcol~=1)|(rem(hrow,2)==0),
   error('h must be a smoothing window with odd length'); end;
  sigma=p3;
  title([method,', Lg=',num2str(Lg),', Lh=',num2str(Lh),...
         ' sigma=',num2str(sigma),', Nf=',num2str(Nf2),...
         LinLogStr, DisplayStr, ', Threshold=',num2str(threshold),'%']);
 
 elseif strcmp(method,'GRD')
  g=p1; [grow,gcol]=size(g); Lg=(grow-1)/2;
  if (gcol~=1)|(rem(grow,2)==0),
   error('g must be a smoothing window with odd length'); end;
  h=p2; [hrow,hcol]=size(h); Lh=(hrow-1)/2; 
  if (hcol~=1)|(rem(hrow,2)==0),
   error('h must be a smoothing window with odd length'); end;
  title([method,', Lg=',num2str(Lg),', Lh=',num2str(Lh),...
         ', rs =',num2str(p3), ', M/N =',num2str(p4), ...
         ', Nf =',num2str(Nf2), ...
         LinLogStr, DisplayStr, ', Threshold=',num2str(threshold),'%']);
 
 elseif strcmp(method,'MSC' ) | strcmp(method,'RMSC' )
  f0T=p1; if (f0T<=0), error('f0T must be positive'); end;
  title([method, ', f0T=',num2str(f0T), ', Nf=',num2str(Nf2),...
        LinLogStr, DisplayStr, ', Threshold=',num2str(threshold),'%']);
 
 elseif strcmp(method,'RGAB' )
  Nh=p1; if (Nh<=0), error('Nh must be positive'); end;
  title([method, ', Nh=',num2str(Nh), ', Nf=',num2str(Nf2),...
        LinLogStr, DisplayStr, ', Threshold=',num2str(threshold),'%']);
 
 
 elseif strcmp(method,'DFLA' ) | strcmp(method,'UNTER' )| ...
        strcmp(method,'BERT' ),
  N=p1; 
  if (N<=0),                error('N must be positive'); end;
  title([method, ', N=',num2str(N), LinLogStr, DisplayStr, ', Threshold=',...
         num2str(threshold), '%']);  
 
 elseif strcmp(method,'SCALO'),
  Nh0=p1; N=p2; 
  if (Nh0<0),               error('Nh0 must be positive'); end;
  if (N<=0),                error('N must be positive'); end;
  if (Nh0>0), 
   title([method, ', Morlet wavelet, Nh0=', num2str(Nh0), ...
         ', N=',num2str(N), LinLogStr, DisplayStr, ', Thld=',...
         num2str(threshold), '%']);  
  else 
   title([method, ', Mexican hat, N=',num2str(N), LinLogStr, DisplayStr, ...
          ', Thld=', num2str(threshold), '%']);  
  end
 
 elseif strcmp(method,'ASPW'),
  Nh0=p1; Ng0=p2; N=p3; 
  if (Nh0<0),               error('Nh0 must be positive'); end;
  if (Ng0<0),               error('Ng0 must be positive'); end;
  if (N<=0),                error('N must be positive'); end;
  if (Nh0>0), 
   title([method, ', Morlet wlt, Nh0=', num2str(Nh0), ', Ng0=',...
         num2str(Ng0), ', N=',num2str(N), LinLogStr, DisplayStr, ', Thld=',...
         num2str(threshold), '%']);  
  else 
   title([method, ', Mexican hat, Ng0=',num2str(Ng0),...  
         ', N=',num2str(N), LinLogStr, DisplayStr, ', Thld=',... 
         num2str(threshold), '%']);  
  end
 
 elseif strcmp(method,'SPAW'),
  K=p1; Nh0=p2; Ng0=p3; N=p4; 
  if (Nh0<0),               error('Nh0 must be positive'); end;
  if (Ng0<0),               error('Ng0 must be positive'); end;
  if (N<=0),                error('N must be positive'); end;
  if (Nh0>0), 
   title([method, ', K=', num2str(K), ', Morlet wlt, Nh0=',...
         num2str(Nh0), ', Ng0=',...
         num2str(Ng0), ', N=',num2str(N), LinLogStr, DisplayStr, ', Thld=',...
         num2str(threshold), '%']);  
  else 
   title([method, ', K=', num2str(K), ', Mexican hat, Ng0=',...
         num2str(Ng0),', N=',num2str(N), LinLogStr, DisplayStr, ', Thld=',...
         num2str(threshold), '%']);  
  end
 
 
 elseif strcmp(method,'GABOR'),
  N=p1; Q=p2; h=p3; [hrow,hcol]=size(h); Lh=(hrow-1)/2;
  if (hcol~=1)|(rem(hrow,2)==0),
   error('h must be a smoothing window with odd length'); end;
  title([method, ', Lh=',num2str(Lh), ', Nf=',...
         num2str(Nf2),', N=',num2str(N),', Q=',num2str(Q),...
         LinLogStr, DisplayStr, ', Thld=',num2str(threshold),'%']);
 
 end;
end

% add the correct legend on the axes            时间、频率显示图例
if unitHz==1,
 xlabel('Time [s]'); ylabel('Frequency [Hz]');
elseif unitHz==2,
 xlabel('Time [ms]'); ylabel('Frequency [kHz]');
elseif unitHz==3,
 xlabel('Time [祍]'); ylabel('Frequency [MHz]');
end


if (isgridsig & issig),		% Updating of the grids
 axes(axsig); grid on	
elseif (~isgridsig & issig),
 axes(axsig); grid off
end 
if (isgridspec & isspec),
 axes(axspec); grid on
elseif (~isgridspec & isspec),
 axes(axspec); grid off
end 

if (isgridtfr),			% upating of the grid on the tfr
 axes(axtfr); grid on  
elseif (~isgridtfr),
 axes(axtfr); grid off
end