DCT and IDCT Without Using Function (Matlab 2013a)

Program:

clc
clear all
x = input('Enter the sequence = ');
l = length(x);
for k = 1:l
    y(k) = 0;
    if (k==1)
        w(k) = 1/sqrt(l);
    else
        w(k) = sqrt(2/l);
    end       
    for n = 1:l
        y(k) = y(k)+x(n)*cos(pi*(2*n-1)*(k-1)/(2*l));
    end
    y(k) = y(k)*w(k);
end
t = 0:l-1;
subplot(2,2,1);
stem(t,x);
ylabel('Amplitude ---->');
xlabel('n ---->');
title('Input Sequence');
grid on;
magnitude = abs(y); % Find the magnitudes of individual DCT points
disp('DCT Sequence = ');
disp(magnitude);
%code block to plot the DCT sequence
t = 0:l-1;
subplot(2,2,2);
stem(t,magnitude);
ylabel('Amplitude ---->');
xlabel('K ---->');
title('DCT Sequence');
grid on;
phase = angle(y); % Find the phases of individual DCT points
disp('Phase = ');
disp(phase);
%code block to plot the phase response
t = 0:l-1;
subplot(2,2,3);
stem(t,phase);
ylabel('Phase ---->');
xlabel('K ---->');
title('Phase Response');
grid on;
for n = 1:l
    X(n) = 0;
    if (n==1)
        w(n) = 1/sqrt(l);
    else
        w(n) = sqrt(2/l);
    end        
    for k = 1:l
        X(n) = X(n)+w(k)*y(k)*cos(pi*(2*n-1)*(k-1)/(2*l));
    end
end
t = 0:l-1;
subplot(2,2,4);
stem(t,X);
disp('IDCT Sequence = ');
disp(X);
ylabel('Amplitude ---->');
xlabel('n ---->');
title('IDCT sequence');
grid on;

Output:

Enter the sequence = [2 4 3 1 2]
DCT Sequence =
    5.3666    1.1152    0.8279    1.8045    0.1208

Phase =
         0         0    3.1416    3.1416         0

IDCT Sequence =
    2.0000    4.0000    3.0000    1.0000    2.0000

Waveforms: