function [DSS, lambdas] = getDiscreteSlepianSequence(N, a, b)
% I is the interval, [a,b] \subset [-1,1]. 

    CorrSeq = (b-a)/2*sinc((b-a)/2*(0:N-1)).*exp(1i*(b+a)/2*pi*(0:N-1));
    CorrMat = toeplitz(CorrSeq.'); 

    [DSS, Lambda] = eig(CorrMat); 

    lambdas = diag(Lambda); 

    [~, order] = sort(lambdas, 'descend'); 
    lambdas = lambdas(order); 

    DSS = DSS(:, order); 

end
