% spectrum.m - compute fft of X(4848:5359) (512 samples around 5104)
 
clear Pyy
W = hamming(512);        % 32 ms window at 16kHz SR
 
Y = fft(W .* X(4848:5359),512);
Pyy = (Y .* conj(Y)/512)';
 
S = log(Pyy(1:256));
f = 16000 * (0:255)/512;
plot(f,S(1:256))
 
xlabel('f (Hz)')