function y = dct1( u )
% The discrete cosine transform of type 1.
%
% This code also works for complex inputs. 
%
% To obtain the eigenvectors of 
% 
% [ 2  -2             ]
% [ -1  2  -1         ]
% [                   ]
% [           -1  2 -1]
% [              -2  2]
% 
% try C = dct1( eye( N ) );
% 
% Alex Townsend, September 2014.

% Fold out to a FFT:
n = size(u, 1);
% cos(kz) = (e^(ikz)+e^(-ikz))/2 
% The term e^(-ikz)  =>  roots of unity in reverse order, flip!
y = fft( [ u ; u(end-1:-1:2 ,: ) ]/2 );
% Take one of the halves: 
y = y( 1:n, : );

end