diff nonExposed/utils/powspec.m @ 14:b1901e8d8f5f

initial commit

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/nonExposed/utils/powspec.m	Tue Mar 17 09:34:13 2015 +0100
@@ -0,0 +1,65 @@
+function [y,e] = powspec(x, sr, wintime, steptime, dither)
+%[y,e] = powspec(x, sr, wintime, steptime, sumlin, dither)
+%
+% compute the powerspectrum and frame energy of the input signal.
+% basically outputs a power spectrogram
+%
+% each column represents a power spectrum for a given frame
+% each row represents a frequency
+%
+% default values:
+% sr = 8000Hz
+% wintime = 25ms (200 samps)
+% steptime = 10ms (80 samps)
+% which means use 256 point fft
+% hamming window
+%
+% $Header: /Users/dpwe/matlab/rastamat/RCS/powspec.m,v 1.3 2012/09/03 14:02:01 dpwe Exp dpwe $
+
+% for sr = 8000
+%NFFT = 256;
+%NOVERLAP = 120;
+%SAMPRATE = 8000;
+%WINDOW = hamming(200);
+
+if nargin < 2
+  sr = 8000;
+end
+if nargin < 3
+  wintime = 0.025;
+end
+if nargin < 4
+  steptime = 0.010;
+end
+if nargin < 5
+  dither = 1;
+end
+
+winpts = round(wintime*sr);
+steppts = round(steptime*sr);
+
+NFFT = 2^(ceil(log(winpts)/log(2)));
+%WINDOW = hamming(winpts);
+%WINDOW = [0,hanning(winpts)'];
+WINDOW = [hanning(winpts)'];
+% hanning gives much less noisy sidelobes
+NOVERLAP = winpts - steppts;
+SAMPRATE = sr;
+
+% Values coming out of rasta treat samples as integers, 
+% not range -1..1, hence scale up here to match (approx)
+y = abs(specgram(x*32768,NFFT,SAMPRATE,WINDOW,NOVERLAP)).^2;
+
+% imagine we had random dither that had a variance of 1 sample 
+% step and a white spectrum.  That's like (in expectation, anyway)
+% adding a constant value to every bin (to avoid digital zero)
+if (dither)
+  y = y + winpts;
+end
+% ignoring the hamming window, total power would be = #pts
+% I think this doesn't quite make sense, but it's what rasta/powspec.c does
+
+% that's all she wrote
+
+% 2012-09-03 Calculate log energy - after windowing, by parseval
+e = log(sum(y));