wolffd@0: % pic = CorrelogramFrame(data, picWidth, start, winLen)
wolffd@0: % Compute one frame of a correlogram.  The input data is a
wolffd@0: % two-dimensional array of cochlear data, each row representing
wolffd@0: % firing probabilities from one cochlear channel.  The output
wolffd@0: % picture is a two dimensional array of width "picWidth".
wolffd@0: %
wolffd@0: % The correlogram is computed with autocorrelation using 
wolffd@0: % data from the input array.  For each channel, the data from
wolffd@0: % is extracted starting at column "start" and extending for 
wolffd@0: % "winLength" time steps.  
wolffd@0: 
wolffd@0: % (c) 1998 Interval Research Corporation
wolffd@0: 
wolffd@0: function pic = CorrelogramFrame(data, picWidth, start, winLen)
wolffd@0: 
wolffd@0: if nargin < 2
wolffd@0: 	disp('Syntax: pic=CorrelogramFrame(data, picWidth[, start, len])');
wolffd@0: 	return
wolffd@0: end
wolffd@0: 
wolffd@0: if nargin < 3
wolffd@0: 	start = 1;
wolffd@0: end
wolffd@0: 
wolffd@0: if nargin < 4
wolffd@0: 	[channels, winLen] = size(data);
wolffd@0: end
wolffd@0: 	
wolffd@0: [channels, dataLen] = size(data);
wolffd@0: start = max(1, start);
wolffd@0: last = min(dataLen, start+winLen-1);
wolffd@0: 
wolffd@0: pic = zeros(channels, picWidth);
wolffd@0: fftSize = 2^(nextpow2(max(picWidth, winLen))+1);
wolffd@0: % disp(['CorrelogramFrame fftSize is ' int2str(fftSize)]);
wolffd@0: 
wolffd@0: a = .54;
wolffd@0: b = -.46;
wolffd@0: wr = sqrt(64/256);
wolffd@0: phi = pi/winLen;
wolffd@0: ws = 2*wr/sqrt(4*a*a+2*b*b)*(a + b*cos(2*pi*(0:winLen-1)/winLen + phi));
wolffd@0: 
wolffd@0: for i=1:channels
wolffd@0: 	f = zeros(1, fftSize);
wolffd@0: %	d = zeros(1, winLen);
wolffd@0: %	d(1:(last-start+1)) = data(i,start:last) .* ws(1:(last-start+1));
wolffd@0: %	f(1:(winLen/2)) = d(winLen/2+1:winLen);
wolffd@0: %	f(fftSize-(winLen/2)+1:fftSize) = d(1:min(length(d),winLen/2));
wolffd@0: 	f(1:(last-start+1)) = data(i,start:last) .* ws(1:(last-start+1));
wolffd@0: 	f = fft(f);
wolffd@0: 	f = ifft(f.*conj(f));
wolffd@0: 	pic(i,:) = real(f(1:picWidth));
wolffd@0: 	if pic(i,1) > pic(i,2) & pic(i,1) > pic(i,3)
wolffd@0: 		pic(i,:) = pic(i,:)/sqrt(pic(i,1));
wolffd@0: 	else
wolffd@0: 		pic(i,:) = zeros(1,picWidth);
wolffd@0: 	end
wolffd@0: end