my_DCT(x)

function X = my_DCT(x)
% Discrete Cosine Transform Computation.
% This function computes the DCT of a real or complex input signal x[n] 
% stored in vector x by use of it's relationship with the DFT.

% First convert input vector x to a row vector for subsequent processing:
x = x(:).';
N = length(x);

% Create the following symmetric sequence y[n] = x[n] + x[2N-n-1].
x1 = [x fliplr(x)];

Y = fft(x1);  % Take the 2N-point DFT of the new sequence.

k = 0:N-1;
phasor = exp(-1i*pi*k/(2*N));

% Multiply the first N elements of Y with the complex phasor.
if isreal(x) 
       X = real(phasor.*Y(1:N));
else
       X = phasor.*Y(1:N);
end

% Scale as necessary to make it consistent with the MatLab dct.m function.
X = [1/(2*sqrt(N))*X(1) 1/sqrt(2*N)*X(2:end)];