# spectral factorization using cepstral deconvolution

Started by December 23, 2006
```Hi there,

I want to find the minimum phase spectral factor of a real autocorrelation
sequence using cepstral deconvolution. However I run into problems when the
spectrum has a null. Is there a way to get around this problem, or is the
cepstral deconvolution doomed to fail in case of a null? Could anyone
suggest me another efficient method that will take about the same time for
a sequence of length ~2^13-2^14. (Could I use Levinson algorithm?)

In case anyone wonders if I implement this correctly in Matlab, I run the
following code:

R = fft(r);
Rhat = 1/2*log(R);
rcep = ifft(Rhat);
rcc(1) = rcep(1)/2;
rcc(2:N/2) =  rcep(2:N/2);
rcc(N/2+1:N) = 0;
Rh = exp(fft(rcc));
x = real(ifft(Rh));

P.S. r is the autocorrelation sequence, i.e. r = [1; 0.5; zeros(61,1);
0.5], for which the above problem occurs. N = 64 in this case.

Thanks 1/eps :)

Emre
```
```"emre" <eguven@ece.neu.edu> wrote in message
news:SNudnbv_UrDlLRHYnZ2dnUVZ_hy3nZ2d@giganews.com...
> Hi there,
>
> I want to find the minimum phase spectral factor of a real autocorrelation
> sequence using cepstral deconvolution. However I run into problems when
the
> spectrum has a null. Is there a way to get around this problem, or is the
> cepstral deconvolution doomed to fail in case of a null? Could anyone
> suggest me another efficient method that will take about the same time for
> a sequence of length ~2^13-2^14. (Could I use Levinson algorithm?)
>
> In case anyone wonders if I implement this correctly in Matlab, I run the
> following code:
>
> R = fft(r);
> Rhat = 1/2*log(R);
> rcep = ifft(Rhat);
> rcc(1) = rcep(1)/2;
> rcc(2:N/2) =  rcep(2:N/2);
> rcc(N/2+1:N) = 0;
> Rh = exp(fft(rcc));
> x = real(ifft(Rh));
>
> P.S. r is the autocorrelation sequence, i.e. r = [1; 0.5; zeros(61,1);
> 0.5], for which the above problem occurs. N = 64 in this case.
>
> Thanks 1/eps :)
>
> Emre

I doubt it should ever really be zero in the  real world since there is also
noise of some description.
Besides, spectral factorization doesn't mean anything unless you have noise
in the first place!! You should always consider there to be some  residual
white noise besides any coloured additive noise.

Tam

--
Posted via a free Usenet account from http://www.teranews.com

```
```"emre" <eguven@ece.neu.edu> wrote in message
news:SNudnbv_UrDlLRHYnZ2dnUVZ_hy3nZ2d@giganews.com...
> Hi there,
>
> I want to find the minimum phase spectral factor of a real autocorrelation
> sequence using cepstral deconvolution. However I run into problems when
the
> spectrum has a null. Is there a way to get around this problem, or is the
> cepstral deconvolution doomed to fail in case of a null? Could anyone
> suggest me another efficient method that will take about the same time for
> a sequence of length ~2^13-2^14. (Could I use Levinson algorithm?)
>
> In case anyone wonders if I implement this correctly in Matlab, I run the
> following code:
>
> R = fft(r);
> Rhat = 1/2*log(R);
> rcep = ifft(Rhat);
> rcc(1) = rcep(1)/2;
> rcc(2:N/2) =  rcep(2:N/2);
> rcc(N/2+1:N) = 0;
> Rh = exp(fft(rcc));
> x = real(ifft(Rh));
>
> P.S. r is the autocorrelation sequence, i.e. r = [1; 0.5; zeros(61,1);
> 0.5], for which the above problem occurs. N = 64 in this case.
>
> Thanks 1/eps :)
>
> Emre

I miss-read your original post - yes - you don't need the autocorrelation -

Tam

--
Posted via a free Usenet account from http://www.teranews.com

```
```Well, this really is a filter design problem, not spectrum estimation, etc.
I already assume perfect knowledge of discrete power spectrum at the
beginning. My starting point is actually the autocorrelation of a
deterministic (but unknown) sequence, from which I want to find the
(minimum-phase) signal that results in the given autocorrelation.

I hope this clarifies the problem.

Thanks again,
Emre

>I miss-read your original post - yes - you don't need the autocorrelation
-
>
>
>Tam
```
```emre skrev:
> Hi there,
>
> I want to find the minimum phase spectral factor of a real autocorrelation
> sequence using cepstral deconvolution. However I run into problems when the
> spectrum has a null. Is there a way to get around this problem, or is the
> cepstral deconvolution doomed to fail in case of a null?

Cepstrum processing is notoriously unstable, for two reasons.
One reason is nulls in the spectra, as you have discovered, beacuse
computing the cepstrum involves computing the logarithm of the
spectrum magnitude. The second reason is that one some times
need to unwrap the phase of the spectrum, which is not trivial.

Rune

```
```>Cepstrum processing is notoriously unstable, for two reasons.
>One reason is nulls in the spectra, as you have discovered, beacuse
>computing the cepstrum involves computing the logarithm of the
>spectrum magnitude. The second reason is that one some times
>need to unwrap the phase of the spectrum, which is not trivial.
>
>Rune

Thanks a lot Rune. In my specific case I have a real spectrum, so I don't
need phase unwrapping, but the former issue is apparently troublesome
enough.. I wish there was a computationally efficient alternative.

Any suggestions?

Emre
```
```emre skrev:
> >Cepstrum processing is notoriously unstable, for two reasons.
> >One reason is nulls in the spectra, as you have discovered, beacuse
> >computing the cepstrum involves computing the logarithm of the
> >spectrum magnitude. The second reason is that one some times
> >need to unwrap the phase of the spectrum, which is not trivial.
> >
> >Rune
>
> Thanks a lot Rune. In my specific case I have a real spectrum, so I don't
> need phase unwrapping, but the former issue is apparently troublesome
> enough.. I wish there was a computationally efficient alternative.
>
> Any suggestions?

If the minimum hase representation is what you want, exploiting
an AR(p) model springs to mind...

Rune

```
```emre wrote:
> I want to find the minimum phase spectral factor of a real autocorrelation
> sequence using cepstral deconvolution. However I run into problems when the
> spectrum has a null.

You can't take the log of a null.  What you can try is to add some
decreasing noise floors, take the cepstrums, and interpolate the
limit as the noise floor decreases, if a limit seems to exist or to
just converge somewhere within your precision requirements.

IMHO. YMMV.
--
rhn A.T nicholson d.0.t C-o-M

```
```Hi Emre,

IIRC from >20 years ago, there was an approach to cepstral
deconvolution of f(n), n>=0, where there were nulls in the spectrum,
that multiplied f(n) with a decaying exponential g(n)=e^(-an), for
n>=0, a>0.  Again, IIRC, the multiplication shifted the zeroes away
from the unit circle so the the nulls were no longer present AND
modified the cepstrum in a very predictable way.  IIRC the cepstrum
could then be modified, inverted, and the scaling by g(n) removed.

I think you can work through the math (I did at the time), but I
believe I saw this in Ron Schafer's MS (or is it S.M.?) disertation.

If anyone remembers more detail, feel free to contibute.

Dirk

Dirk Bell
DSP Consultant

emre wrote:
> Hi there,
>
> I want to find the minimum phase spectral factor of a real autocorrelation
> sequence using cepstral deconvolution. However I run into problems when the
> spectrum has a null. Is there a way to get around this problem, or is the
> cepstral deconvolution doomed to fail in case of a null? Could anyone
> suggest me another efficient method that will take about the same time for
> a sequence of length ~2^13-2^14. (Could I use Levinson algorithm?)
>
> In case anyone wonders if I implement this correctly in Matlab, I run the
> following code:
>
> R = fft(r);
> Rhat = 1/2*log(R);
> rcep = ifft(Rhat);
> rcc(1) = rcep(1)/2;
> rcc(2:N/2) =  rcep(2:N/2);
> rcc(N/2+1:N) = 0;
> Rh = exp(fft(rcc));
> x = real(ifft(Rh));
>
> P.S. r is the autocorrelation sequence, i.e. r = [1; 0.5; zeros(61,1);
> 0.5], for which the above problem occurs. N = 64 in this case.
>
> Thanks 1/eps :)
>
> Emre

```
```dbell skrev:
> Hi Emre,
>
> IIRC from >20 years ago, there was an approach to cepstral
> deconvolution of f(n), n>=0, where there were nulls in the spectrum,
> that multiplied f(n) with a decaying exponential g(n)=e^(-an), for
> n>=0, a>0.  Again, IIRC, the multiplication shifted the zeroes away
> from the unit circle so the the nulls were no longer present AND
> modified the cepstrum in a very predictable way.  IIRC the cepstrum
> could then be modified, inverted, and the scaling by g(n) removed.
>
> I think you can work through the math (I did at the time), but I
> believe I saw this in Ron Schafer's MS (or is it S.M.?) disertation.

I wouldn't be surprised if this was investigated by Schafer; I know
Oppenheim has written a thesis over cepstrum processing.

However, nulls is only one of the problems with cepstra. Another
big problem is phase unwrapping.

In order to extract the interesting pulse shapes, one needs to
both unwrap the phase of the spectrum unambiguously, and also
keep the FT pair of the log magnitude and unwrapped phase
in sync with the requirements to causal signals.

Not very easy.

Rune

```