Code covered by the BSD License

# Single Pole Cepstrum

### Speech Processing (view profile)

29 Jan 2014 (Updated )

This exercise compares and contrasts three methods for computing the complex cepstrum of a signall.

[nseq,nplot,x,xhat1,xhat2,xhat3]=compute_single_pole_cepstrum(a)
function [nseq,nplot,x,xhat1,xhat2,xhat3]=compute_single_pole_cepstrum(a)
%
% compute complex cepstrum of minimum phase sequence, x[n]=a.^n u[n],
% in 3 ways
%   method 1--xhat1 is the analytical solution
%   method 2--xhat2 is the recursion solution
%   method 3--xhat3 is the data sequence solution via fft and ifft
%
% use sequence x[n]=a.^n u[n]
% Input:
%   a: radius of single pole
%
% Outputs:
%   nseq: length of input sequence in samples
%   nplot: number of cepstral coefficients to plot
%   x: input sequence
%   xhat1: analytical cepstrum
%   xhat2: recursion cepstrum
%   xhat3: cepstrum from fft, log and ifft

% define temporal and cepstral sequence parameters
% nseq: length of input sequence (samples)
% nplot: number of cepstral coefficients (quefrencies) to plot
nseq=1024;
nplot=100;
x(1:nseq)=a.^(0:nseq-1);

% method 1 -- calculate xhat1 using analytical formula
xhat1(1)=0;
xhat1(2:nseq)=a.^(1:nseq-1)./(1:nseq-1);

% method 2 -- calculate xhat2 using recursion relation
xhat2(1)=log(x(1));
for n=1:nseq-1
xhat2(n+1)=x(n+1)/x(1)-sum(xhat2(2:n).*x(n:-1:2).*(1:n-1))/((n)*x(1));
end

% method 3 -- compute xhat3 from data sequence using fft, log and ifft
xc=x(1:nseq)';
xcf=fft(xc,nseq);
xcfl=log(xcf);
xhat3=ifft(xcfl);

end