
% KAIST CS891 Special Lecture on Music Information Retrieval (MIR)
% Lab1.3: Spectral leaking, windowing
% 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
nfft = 1024; % fft size
xs = x(1:nfft); % snippet of x
X = fft(xs); % complex numbers
Xdb = 20*log10(abs(X)/max(abs(X))); % normalized magnitude in dB scale

% windowing
win = hann(nfft)'; % Hann window; be careful about transpose operator '
xw = win.*xs; % windowed signal; be careful about point-wise operator .*
Xw = fft(xw);
Xwdb = 20*log10(abs(Xw)/max(abs(Xw)));


% display
% time-domain wavform
subplot(221);
ti = 0:T:dur-T; % in seconds
plot(xs); xlim([0 nfft]);
hold on; plot(win,'r'); xlim([0 nfft]);
title('Original signal');
xlabel('time (samples)');
ylabel('amplitude');

% windowed signal
subplot(222);
plot(xw); xlim([0 nfft]);
hold on; plot(win,'r'); xlim([0 nfft]);
title('Windowed signal');
xlabel('time (samples)');

% fft of original signal
subplot(223);
fi = 0:fs/nfft:fs-fs/nfft;
plot(fi,Xdb); xlim([0 fs/2]); ylim([-80 0]); 
title('Magnitude response of original signal');
xlabel('frequency (Hz)');
ylabel('normalized magnitude (dB)');

% fft of windowed signal
subplot(224);
plot(fi,Xwdb); xlim([0 fs/2]); ylim([-80 0]);
title('Magnitude response of windowed signal');
xlabel('frequency (Hz)');
