Home > NoiseTools > nt_sgram.m

# nt_sgram

## PURPOSE

[s,f,t]=nt_sgram(x,window,noverlap,nfft,sr,flags) - spectrogram

## SYNOPSIS

function [s,f,t]=nt_sgram(x,window,noverlap,nfft,sr,flags)

## DESCRIPTION

    [s,f,t]=nt_sgram(x,window,noverlap,nfft,sr,flags) - spectrogram

     Without output arguments: plot cubic root of power spectrogram.
     With output arguments: return power spectrogram.
      s: spectrogram matrix (frequency X time)
      f: array of frequencies
      t: array of times

      x: data (if multidimensional, power is averaged over higher dimensions)
      window: size of hamming window in samples, if vector, shape of window
      noverlap: frame period in samples
      nfft: number of points in FFT (default 2^nextpow2(window size))
      sr: sampling rate.
      flags: 'nomean': remove mean of each window

      NoiseTools

## CROSS-REFERENCE INFORMATION

This function calls:

- nt_demean [y,mn]=nt_demean(x,w) - remove weighted mean over cols
- nt_greetings nt_greetings - display message the first time the toolbox is used
- nt_imagescc nt_imagescc - plot image with symmetric scaling
- nt_vecmult y=nt_vecmult(x,v) - multiply all rows or columns of matrix by vector

This function is called by:

- nt_bias_cluster [c0,c1,A,todss,pwr0,pwr1]=nt_bias_cluster(x,dsr,flags) - cluster covariance
- nt_cluster_jd [IDX,todss,SCORE,COVS]=nt_cluster_jd(x,dsr,smooth,flags,init,verbose) - cluster with joint diagonalization

## SOURCE CODE

    0001 function [s,f,t]=nt_sgram(x,window,noverlap,nfft,sr,flags)
    0002 %[s,f,t]=nt_sgram(x,window,noverlap,nfft,sr,flags) - spectrogram
    0003 %
    0004 % Without output arguments: plot cubic root of power spectrogram.
    0005 % With output arguments: return power spectrogram.
    0006 %  s: spectrogram matrix (frequency X time)
    0007 %  f: array of frequencies
    0008 %  t: array of times
    0009 %
    0010 %  x: data (if multidimensional, power is averaged over higher dimensions)
    0011 %  window: size of hamming window in samples, if vector, shape of window
    0012 %  noverlap: frame period in samples
    0013 %  nfft: number of points in FFT (default 2^nextpow2(window size))
    0014 %  sr: sampling rate.
    0015 %  flags: 'nomean': remove mean of each window
    0016 %
    0017 %  NoiseTools
    0018 nt_greetings;
    0019
    0020 if nargin<1; error('!'); end
    0021 if nargin<2; window=[]; end
    0022 if nargin<3; noverlap=[]; end
    0023 if nargin<4; nfft=[]; end
    0024 if nargin<5; sr=[]; end
    0025 if nargin<6; flags=[]; end
    0026
    0027 if isempty(window)
    0028     window=2^nextpow2(round(size(x,1)/8));
    0029 end
    0030 window=window(:);
    0031 if numel(window)==1;
    0032     window=hanning(window);
    0033 end
    0034 if isempty(noverlap)
    0035     noverlap=floor(numel(window)/2);
    0036 end
    0037 if isempty(nfft)
    0038     nfft=2^nextpow2(numel(window));
    0039 end
    0040 if isempty(sr)
    0041     sr=1;
    0042 end
    0043
    0044
    0045 x=reshape(x,[size(x,1), prod(size(x))/size(x,1)]); % fold all higher dimensions
    0046
    0047 nframes=floor( (size(x,1)-numel(window))/noverlap )  +  1;
    0048 if nframes<1; error('data too short for window size'); end
    0049
    0050 s=zeros(nfft/2+1,nframes);
    0051
    0052 for k=1:nframes
    0053     start=(k-1)*noverlap;
    0054     xx=nt_vecmult(x(start+1:start+size(window),:),window);
    0055     if ~isempty(flags) & strcmp(flags,'demean');
    0056         xx=nt_demean(xx);
    0057     end
    0058     XX=abs(fft(xx,nfft)).^2;
    0059     s(:,k)=sum(XX(1:(nfft/2+1),:),2);
    0060 end
    0061
    0062 f=(0:nfft/2)*sr/nfft;
    0063 t= (numel(window)/2 + (0:nframes-1)*noverlap) / sr;
    0064
    0065
    0066
    0067 if nargout==0;
    0068     ss=s.^(1/5);
    0069     %ss=ss-repmat(mean(ss,2),1,size(ss,2));
    0070
    0071     nt_imagescc(ss);
    0072     ytick=niceticks(sr/2);
    0073     set(gca,'ytick',1+ytick*nfft/sr, 'yticklabel',num2str(ytick'));
    0074     ylabel('Hz');
    0075     xtick=niceticks(size(x,1)/sr);
    0076     set(gca,'xtick',1+xtick/noverlap*sr - nfft/2/noverlap, 'xticklabel',num2str(xtick'));
    0077     xlabel('s');
    0078     set(gca,'tickdir','out')
    0079     s=[];t=[];f=[];
    0080     axis xy
    0081 end
    0082
    0083