
% KAIST CS891 Special Lecture on Music Information Retrieval (MIR)
% Lab1.2: Frequency resolution
% Kyogu Lee (August 2008)

clear all;
close all;

% Generate a sinusoid
fs = 44100; % sampling rate
amp = 1; % amplitude
freq = 440; % frequency in Hz
T = 1/fs; % sampling period
dur = 2; % duration in seconds
phi = 0; % phase in radian

% time-domain signal
x = amp*sin(2*pi*freq*[0:T:dur-T]+phi);

% Fourier transform
nfft1 = 1024; % fft size
xs = x(1:nfft1); % snippet of x
X = fft(xs); % complex numbers
Xdb = 20*log10(abs(X)); % magnitude in dB scale

% zero-padding
nfft2 = 8192;
xz = [xs,zeros(1,nfft2-length(xs))];
Xz = fft(xz);
Xzdb = 20*log10(abs(Xz));

% display
% time-domain wavform
subplot(511);
ti = 0:T:dur-T; % in seconds
plot(ti(1:nfft1),xs);

% zero-padded signal
subplot(512);
plot(ti(1:nfft2),xz);

% fft magnitude
subplot(513);
fi = 0:fs/nfft1:fs-fs/nfft1; % freq. resolution = fs/nfft1 = 43.07 Hz
plot(fi,Xdb,'bo-'); xlim([400 500]);

subplot(514);
fi = 0:fs/nfft2:fs-fs/nfft2; % freq. resolution = fs/nfft2 = 5.38
plot(fi,Xzdb,'bo-'); xlim([400 500]);

% fft(x,nfft) does zeropadding automatically when nfft > length(x)
Xz2 = fft(xs,nfft2);
Xz2db = 20*log10(abs(Xz2));
subplot(515);
plot(fi,Xz2db,'bo-'); xlim([400 500]);