Monday, September 25, 2017

Prony Method

Given N equispaced data-points \(F_i = F(t = i \Delta t)\), where \(i = 0, 1, ..., N-1\), the Prony method can be used to fit a sum of m decaying exponenitals: \[F(t) = \sum_{i=1}^{m} a_i e^{b_i t}. \] The 2m unknowns are \(a_i\) and \(b_i\).

In the Prony method, the number of modes in the exponential (m) is pre-specified. There are other methods, which are more general.

Here is a python subprogram which implements the Prony method.
If you have arrays t and F, it can be called as:

a_est, b_est = prony(t, F, m)

4 comments:

  1. Dear Clueless Fundatma,
    Can you please post more information about the code or give a reference?
    The following line seems to be missing in the code
    N = len(t)
    probably after the beginning of the function Prony()
    Kind regards,
    João Henriques

    ReplyDelete
  2. Thanks for pointing that error out, Joao. There is a link to the method in the writeup above. I just implemented it.

    ReplyDelete
  3. Dear Clueless Fundatma,
    The method produces an error when you are computing a_est and b_est for a signal like:
    y[k]=math.exp(-3.0*t[k])+math.exp(-2.0*t[k])*math.cos(50.0*2.0*math.pi*t[k]);

    the correction is to assume a complex matrix for the least square solver.

    I modified the line:

    # Set up LLS problem to find the "a"s in step 3
    Amat = np.zeros((N, m),dtype=complex)

    This solves the problem when you need to compute the oscillation frequency of the sine or cosine waves.

    ReplyDelete