FFT Bin Estimation - Quadratic Interpolation Equation Differences

75 Views Asked by At

Regarding FFT peak Bin estimation, excerpts from Julius O. Smith's SPECTRAL AUDIO SIGNAL PROCESSING and PARSHL detail a derivation for curve fitting a magnitude-maximum FFT bin with its two neighbors, totaling 3 points to form a parabola. Equations follow.

$$y(x) = a(x - p)^2 + b \tag{1. ref. parabola}$$

\begin{align} \alpha &= y(-1) \\\\ \beta &= y(0) \\\\ \gamma &= y(1) \end{align}

$$p = \frac 12 \frac{\alpha - \gamma}{\alpha - 2\beta + \gamma} \in [-0.5, 0.5] \tag{2. Quadratic Interpolator}$$

Shortly after stumbling upon this, I came across an adaption of this equation for complex numbers on Eric Jacobsen's page.

$$p_c = \Re \bigg[ -\frac{\alpha - \gamma}{\alpha - 2\beta + \gamma} \bigg] \tag{3. Complex Numbers Quadterp}$$

(If you have noticed the equation is polarity flipped from what's listed in Eric's PDF, their document is in error. You may verify with a quick glance at Eric's MATLab implementation.)

(3) can be interpreted as (2) multipled by $-2$. I would like to figure out why this simple multiplier is so important. Is it compensation for complex number math behavior?

I have observed in GNU Octave that using the 1/2 multiplier from (2) on (3) does indeed throw the output off, as illustrated in the graph below.

Eric stresses the point that their interpolator is not the same as fitting a curve to the DFT output as (2) does; that theirs "uses the quadratic formula to derive an expression that is then applied to the three complex outputs around the DFT magnitude maximizer and the real part of this complex calculation is used as the estimate." Unfortunately that derivation is not further explained on the website or any document linked. My attempts to "use the quadratic equation" only brought (2) back out.

Please explain how (3) was derived, and specifically why it is essentially EQ. (2) multiplied by $-2$. Moreover why is that modification useful?

P.S. As a bonus, here is a graph I've coded based on Eric's MATLab files. 1000 peak bin location estimates are averaged over their actual bin locations using Eric's quadratic estimator. The reference bin locations are black X's, Eric's original estimator is the blue line. Interpret bin "0" as an arbitrary bin with adjacent bins surrounding it.

Multiplying Eric's estimator by 0.5 to make it more akin to EQ. (2) is graphed in red circles, and then a polarity flip to match it exactly to EQ. (2) is plotted as green diamonds.

If anything, this clearly demonstrates that the performance degrades as the equation is converted to EQ. (2).

Graph of bin estimates vs Actual