%% Calculate IMDCT
% Read in an audio file, convert it to mono, and then plot it.

audioIn = audioread('FunkyDrums-44p1-stereo-25secs.mp3');
audioIn = mean(audioIn,2);

figure(1)
plot(audioIn,'bo')
ylabel('Amplitude')
xlabel('Sample Number')

%
% Calculate the MDCT using a 4096-point sine window. Plot the power of the MDCT coefficients over time.

N = 4096;
wdw = sin(pi*((1:N)-0.5)/N);

C = mdct(audioIn,wdw);

figure(2)
surf(20*log10(C.*conj(C)),'EdgeColor','none');
view([0 90])
xlabel('Frame')
ylabel('Frequency')
axis([0 size(C,2) 0 size(C,1)])
colorbar

%
% Transform the representation back to the time domain. Verify the perfect reconstruction property by computing the mean squared error. Plot the reconstructed signal over the original signal.

audioReconstructed = imdct(C,wdw);
err = mean((audioIn-audioReconstructed(1:size(audioIn,1),:)).^2)

figure(1)
hold on
plot(audioReconstructed,'r.')
ylabel('Amplitude')
xlabel('Sample Number')

%
% Copyright 2018 The MathWorks, Inc.