Evaluate Formants Using LPC in Python

I am new to signal processing (and numpy, scipy and matlab for that matter). I'm trying to evaluate vowel formants with LPC in Python by adapting this matlab code:

http://www.mathworks.com/help/signal/ug/formant-estimation-with-lpc-coefficients.html

Here is my code:

#!/usr/bin/env python
import sys
import numpy
import wave
import math
from scipy.signal import lfilter, hamming
from scikits.talkbox import lpc

"""
Estimate formants using LPC.
"""

def get_formants(file_path):

    # Read from file.
    spf = wave.open(file_path, 'r') # http://www.linguistics.ucla.edu/people/hayes/103/Charts/VChart/ae.wav

    # Get file as numpy array.
    x = spf.readframes(-1)
    x = numpy.fromstring(x, 'Int16')

    # Get Hamming window.
    N = len(x)
    w = numpy.hamming(N)

    # Apply window and high pass filter.
    x1 = x * w
    x1 = lfilter([1., -0.63], 1, x1)

    # Get LPC.
    A, e, k = lpc(x1, 8)

    # Get roots.
    rts = numpy.roots(A)
    rts = [r for r in rts if numpy.imag(r) >= 0]

    # Get angles.
    angz = numpy.arctan2(numpy.imag(rts), numpy.real(rts))

    # Get frequencies.
    Fs = spf.getframerate()
    frqs = sorted(angz * (Fs / (2 * math.pi)))

    return frqs

print get_formants(sys.argv[1])

      

Using this file as input, my script returns this list:

[682.18960189917243, 1886.3054773107765, 3518.8326108511073, 6524.8112723782951]

      

I didn't even get to the last steps where they are filtering frequencies by bandwidth because the frequencies in the list are wrong. According to Praat, I should get something like this (this is a formant list for the middle of a vowel):

Time_s     F1_Hz        F2_Hz         F3_Hz         F4_Hz
0.164969   731.914588   1737.980346   2115.510104   3191.775838 

      

What am I doing wrong?

Many thanks

UPDATE:

I changed this

x1 = lfilter([1., -0.63], 1, x1)

to

x1 = lfilter([1], [1., 0.63], x1)

as suggested by Warren Walkers and now I get

[631.44354635609318, 1815.8629524985781, 3421.8288991389031, 6667.5030877036006]

I feel like I'm missing something because F3 is not working very much.

UPDATE 2:

I realized that I was being order

passed on scikits.talkbox.lpc

disabled due to the difference in sample rate. Changed:

Fs = spf.getframerate()
ncoeff = 2 + Fs / 1000
A, e, k = lpc(x1, ncoeff)

      

Now I am getting:

[257.86573127888488, 774.59006835496086, 1769.4624576002402, 2386.7093679399809, 3282.387975973973, 4413.0428174593926, 6060.8150432549655, 6503.3090645887842, 7266.5069407315023]

Much closer to evaluating Praat!

+3


source to share


3 answers


The problem is with the order passed to the lpc function. 2 + fs / 1000

where fs

is the sampling rate is the rule of thumb according to:



http://www.phon.ucl.ac.uk/courses/spsci/matlab/lect10.html

+2


source


I was not able to get the expected results, but I notice two things that may cause some differences:

  • Your code uses [1, -0.63]

    where the MATLAB code from the link you specified has [1 0.63]

    .
  • Your processing is applied to the entire vector at x

    once instead of smaller segments of it (see where the MATLAB code does this :) x = mtlb(I0:Iend);

    .



Hope it helps.

+1


source


There are at least two problems:

  • According to the link, the pre-focus filter is a full pole filter (AR (1)). Coefficient values given here are correct: [1, 0.63]

    . If you use [1, -0.63]

    , you get a low pass filter.

  • You have the first two arguments for scipy.signal.lfilter

    in reverse order.

So, try changing this:

x1 = lfilter([1., -0.63], 1, x1)

      

:

x1 = lfilter([1.], [1., 0.63], x1)

      

I haven't tried using your code yet, so I don't know if these are the only problems.

+1


source







All Articles