view functions/dceps.m @ 0:b4e26b53072f tip

Initial commit.

function c = dceps(X, w, numCoeffs, lambda)
% computes the discrete cepstral coefficients of a sequence of partials X.
%
% c = dceps(X, w, numCoeffs)
%
% w contains the normalized frequencies of the partials. Normalized
% frequencies should be in the range [0 .. pi]
% X and w must be column vectors.
% numCoeffs specifies the number of cepstral coefficients.
%
% Further details can be found in:
% [1] Diemo Schwarz. Spectral envelopes in sound analysis and synthesis.
% Master's thesis, Universitaet Stuttgart, 1998, pp. 36-38.
% [2] T. Galas and X. Rodet. An improved cepstral method for deconvolution
% of source filter systems with discrete spectra: Application to musical
% sound signals. In International Computer Music Conference, 1990.


p = numCoeffs-1;
n = length(X); % number of partials
S = ones(n,1); % amplitude of partials of source spectrum
%lambda = 0.0005; % regularization term

if sum(X == 0) %if one or more element is zero
    error('All elements of X must be nonzero');
end

% compute all r_i (i from 0 to 2p)
r = .5 * sum(cos(w * (0:2*p)));

% compute matrix A
A = zeros(p+1);
B = zeros(p+1);
for i = 0:p
    for j = 0:p
        A(i+1,j+1) = r(i+j+1) + r(abs(i-j)+1); % +1 because 1st idx in r is 0!
    end
    B(i+1,i+1) = 8 * (pi*(i-1))^2;
end

A = A + lambda * B;

% compute b
b = sum( repmat(log(X./S),1,p+1) .* cos(w *(0:p)) )';

% compute c
c = A \ b; % equals inv(A) * b