### Matlab listing: oboeanal.m

The following Matlab script illustrates use of the findpeaks function above to determine the pitch of an oboe tone (given the general location of the correct spectral peak) and configure a spectrum analysis using the rectangular, Hamming, and Blackman windows. This script was used to create Figures 3.16 through 3.18.

```diary './dia/oboeanal.dia';
diary on;

zpf = 5; % desired min zero-padding factor in all FFTs
fn = '../wav/oboe.ff.C4B4.wav';
nx = length(x);
dur = nx/fs;
disp(sprintf(['Read %s, fs = %f, nbits = %d, ',...
'length = %d samples = %0.1f sec'],...
fn,fs,nbits,nx,dur));
soundsc(x,fs);
M = fs/10; % 100 ms window for pitch estimation
N = 2^nextpow2(M*zpf) % zero pad
Xzp = fft(x,N);

figure(K-1);
fmax = 500; % maximum freq in Hz
kmax = fmax*N/fs; % maximum freq in bins
prange = 1:kmax+1;
Xdb = dbn(Xzp(prange)); % dB normalized to 0 max
f = fs*[0:N-1]/N/1000; % frequency axis in kHz
plot(f(prange),Xdb);
ylabel('Magnitude (dB)');
xlabel('Normalized Frequency (cycles/sample)');
[amps,bins] = findpeaks(Xdb,1);
k0 = bins(1)-1;
f0 = k0*fs/N;
nP = fs/f0; % samples per period (used below)
disp(sprintf('Estimated pitch = %f Hz',f0));

% Test pitch estimate:
pause(ceil(dur)); % let previous sound finish
T = 1/fs;
t = 0:T:dur;
y = sawtooth(2*pi*f0*t);
sound(y,fs)

winnames={'none','boxcar','none','hamming','none','blackman'};

% Configure window analysis

for K=[2 4 6]
wintype = winnames{K};
M = ceil(K*nP); % min Blackman length
cmd = sprintf('w = %s(M);',wintype);
eval(cmd);
Nw = 2^nextpow2(M*zpf) % zero pad
xw = x(1:M) .* w;
Xwzp = fft(xw,Nw);
%fmax = 5*f0; % maximum freq in Hz
fmax = fs/4;
kmax = fmax*Nw/fs; % maximum freq in bins
prange = 1:kmax+1;
Xwdb = dbn(Xwzp(prange)); % dB normalized to 0 max

figure(K);
subplot(2,1,1);
plot(xw,'-k');
axis tight;
title(sprintf('Oboe data - %s window, K=%d',wintype,K));
ylabel('Amplitude');
xlabel('Time (samples)');
subplot(2,1,2);
f = fs*[0:Nw-1]/Nw/1000; % frequency axis in kHz
plot(f(prange),Xwdb,'-k'); grid on;
axis([f(prange(1)) f(prange(end)) -90 0]);
ylabel('Magnitude (dB)');
xlabel('Frequency (kHz)');

plotfile = sprintf('../eps/oboe%s.eps',wintype);
saveplot(plotfile);

end

diary off;
```

Next Section:
Matlab listing: myspectrogram.m
Previous Section:
Matlab listing: zpfmin.m