% Function to construct n X n DCT matrix
%     for two dimensional analysis

function C = dctMat(n)
    % tic
    i = 0:n-1;
    j = 0:n-1;
    ai = (1/sqrt(2))*(i==0) + (i~=0);
    C = (sqrt(2/n))*diag(ai)*cos((i'*(2*j+1)*pi)/(2*n));
    % toc
        % one line representation of code below
        %   (dctMat(5000) takes ~1.7 seconds)
        % diag(ai)*cos() allows multiplication of corresponding rows
        %   by transforming n x 1 matrix ai to diagonal matrix with
        %       diagonal entries equal to entries of ai
        %   i.e. (1st row of ai)*(first row of cos())
        % i'*(2*j+1) ensures a 5x5 result matrix
    
    % Easier to read/understand (dctMat(5000) takes 42-44 seconds)
    % tic
    % for i=1:n
    %   for j=1:n
    %     C(i,j)=cos((i-1)*(2*j-1)*pi/(2*n));
    %   end
    % end
    % C=sqrt(2/n)*C;
    % C(1,:)=C(1,:)/sqrt(2);
    % toc