% KAIST CS891 Special Lecture on Music Information Retrieval (MIR)
% Lab2.3: Time-domain Onset Detection
% Kyogu Lee (August 2008)

clear all;
close all;

% read audio file
[x,fs] = wavread('../wav/T39-piano-8k.wav'); % vowel 'ooh'
if size(x,2) > size(x,1)
    x = x'; % make a column vector
end
x = mean(x,2); % make mono
x = x(1:2*fs);

% smoothing window size
M = 256;
win = hamming(M);

% compute E(n)
for i=1:length(x)-M+1
    e(i) = sum(win.*x(i:i+M-1).^2)/M;
end

% compute dE(n)
diffe = diff(e); % first-order difference
for i=1:length(diffe)-M+1
  de(i) = sum(win'.*diffe(i:i+M-1))/M;
end

% compute log(dE)
diffloge = diff(log(e+0.00001)); % +0.00001 - to prevent -Inf
for i=1:length(diffloge)-M+1
  dloge(i) = sum(win'.*diffloge(i:i+M-1))/M;
end

% display
subplot(411);
plot([0:length(x)-1]/fs,x); title('signal'); ylabel('x(t)');
subplot(412);
plot(([0:length(e)-1]+M/2)/fs,e); title('energy'); ylabel('E(t)');
subplot(413);
plot(([0:length(de)-1]+M/2)/fs,de); title('energy difference'); ylabel('dE(t)');
subplot(414);
plot(([0:length(dloge)-1]+M/2)/fs,dloge); title('log energy difference'); ylabel('dlogE(t)');