Spectral peak location estimation using complex DFT

Question

In this paper a simple method to estimate a spectral peak is proposed, by using quadratic interpolation between three samples of the DFT of the signal. Namely, the position of the peak relative to the index $k$ (they call this difference $\delta$) is calculated with this formula in the paper:

$$\delta = \frac{|X_{k+1}|-|X_{k-1}|}{4|X_{k}|-2|X_{k-1}|-2|X_{k+1}|}$$

This picture will make the formula easier to understand:

Then in the paper it is stated that:

The expression is simple, but it is statistically biased and performs poorly in the presence of noise. Some simple changes [...] improve its accuracy dramatically, for instance by using the complex DFT values rather than the magnitudes as follows: $$\delta = -\mathrm{Re}\left[\frac{X_{k+1}-X_{k-1}}{2X_{k}-X_{k-1}-X_{k+1}}\right]$$

I understand where the first equation comes from. One can get to that equation by simply finding the coefficients of a parabola such that $y(m) = |X_m|$ for $m\in \{ k-1, k, k+1\}$, and then one can find the maximum for that parabola and get to the first expression for $\delta$. However, I don't understand why when using the complex DFT, the result the paper shows is that one.

Why does that negative sign appear? Why does the denominator seem to be divided by 2? Why do they take the real part instead of the magnitude, completely ignoring the imaginary part?

Answer 1

I admit that I failed my first attempt.

The text suggests that we would take $g(m) = am^2 + b m + c$, interpolating $g(k) = X(k+l)$, for $k \in \{-1, 0, 1\}$, interpolating this would give us polynomial parameterized in $X[k], X[k-1], X[k+1]$, derive, solve and it gives the optimum value at $m = k+\delta$. Not that easy, and it gave me three roots, none of them seemed to be the desired value :(

The comment on the paper is not helpful but they give us a reference [4], from the one of the authors, and there it gives more details. There you will find the equation (3.7) and they give a concise derivation on the Appendix A.

Derivation

Here I will fill some details of the proof.

You express the interpolator polynomial centered at $k + \delta$, so delta will maximize it by construction (more clever the attempt described above).

$$g(m) = a (m - k - \delta)^2 + b$$

Then $g(m)$, determine delta so that interpolates the three samples. Since the quadratic polynomial interpolating three points is unique, we know that $\delta$ here is unique as well.

$$\begin{eqnarray} X[k-1] &=& g(k-1) &=& b+((k-1) - k - \delta)^2 &=& b + a\cdot(\delta^2+2\delta + 1) \\ X[k ] &=& g(k ) &=& b+((k ) - k - \delta)^2 &=& b + a\cdot(\delta^2) \\ X[k+1] &=& g(k+1) &=& b+((k+1) - k - \delta)^2 &=& b + a\cdot(\delta^2-2\delta + 1) \end{eqnarray}$$

Start looking at the coefficient of the rightmost terms, and note that $$2 X[k] - X[k-1] - X[k+1] = 2 a$$

Also note that

$$ X[k-1] - X[k+1] = 4 a \delta$$

$$\delta = \frac{X[k-1] - X[k+1]}{4a} = \frac{X[k-1] - X[k+1]}{2(2X[k] - X[k-1] - X[k+1])}$$

With the sentence

A simple adaptation provides an estimator suitable for estimating tone frequencies from complex DFT outputs. They introduce the formula we see in the first paper you mentioned.

The adaptation they mention is because $\delta$ is complex, so he takes the real value. Maybe a justification for that would be that $\min_{f} |f - \delta|$ is $Re[\delta]$.