annotate toolboxes/MIRtoolbox1.3.2/AuditoryToolbox/spectrogram.m @ 0:e9a9cd732c1e tip

first hg version after svn
author wolffd
date Tue, 10 Feb 2015 15:05:51 +0000
parents
children
rev   line source
wolffd@0 1 function [array,raw] = spectrogram(wave,segsize,nlap,ntrans);
wolffd@0 2 %function array = spectrogram(wave,segsize,nlap,ntrans);
wolffd@0 3 % defaults spectrogram(wave,128,8,4)
wolffd@0 4 % nlap is number of hamming windows overlapping a point;
wolffd@0 5 % ntrans is factor by which transform is bigger than segment;
wolffd@0 6 % returns a spectrogram 'array' with fourth root of power,
wolffd@0 7 % filter smoothed and formatted for display.
wolffd@0 8
wolffd@0 9 % Added option to return raw spectrogram.... Malcolm 5/26/95
wolffd@0 10 % Added code so that input could be any direction ... Malcolm 5/26/95
wolffd@0 11 % (c) 1998 Interval Research Corporation
wolffd@0 12
wolffd@0 13 if nargin < 4; ntrans=4; end
wolffd@0 14 if nargin < 3; nlap=8; end
wolffd@0 15 if nargin < 2; segsize=128; end
wolffd@0 16
wolffd@0 17 [r c] = size(wave);
wolffd@0 18 if (r < c)
wolffd@0 19 wave = filter([1 -0.95],[1],wave');
wolffd@0 20 else
wolffd@0 21 wave = filter([1 -0.95],[1],wave);
wolffd@0 22 end
wolffd@0 23
wolffd@0 24 s = length(wave);
wolffd@0 25 nsegs = floor(s/(segsize/nlap))-nlap+1;
wolffd@0 26 array = zeros(ntrans/2*segsize,nsegs);
wolffd@0 27 window = 0.54-0.46*cos(2*pi/(segsize+1)*(1:segsize)');
wolffd@0 28 for i = 1:nsegs
wolffd@0 29 seg = zeros(ntrans*segsize,1); % leave half full of zeroes
wolffd@0 30 seg(1:segsize) = ...
wolffd@0 31 window.*wave(((i-1)*segsize/nlap+1):((i+nlap-1)*segsize/nlap));
wolffd@0 32 seg = abs(fft(seg));
wolffd@0 33 % reverse for image display
wolffd@0 34 array(:,i) = seg(((ntrans/2*segsize)+1):(ntrans*segsize));
wolffd@0 35 end
wolffd@0 36
wolffd@0 37 if nargout > 1
wolffd@0 38 raw = array;
wolffd@0 39 end
wolffd@0 40
wolffd@0 41 array = array .* array; % back into power domain for smoothing
wolffd@0 42
wolffd@0 43 for i=1:nsegs % smooth the spectral slices
wolffd@0 44 array(:,i) = filter([.2 1 .2],[1],array(:,i));
wolffd@0 45 end
wolffd@0 46
wolffd@0 47 for i=1:ntrans/2*segsize % smooth the channels
wolffd@0 48 array(i,:) = filter([.2 1 .2],[1],array(i,:));
wolffd@0 49 end
wolffd@0 50
wolffd@0 51 % compress with square root of amplitude (fourth root of power)
wolffd@0 52 off = 0.0001*max(max(array)); % low end stabilization offset,
wolffd@0 53 array = (off+array).^0.25-off^0.25; % better than a threshold hack!
wolffd@0 54 array = 255/max(max(array))*array;
wolffd@0 55
wolffd@0 56
wolffd@0 57