## 4 Proposed Method: XQIFFT
### Interpolation Bias  
We know that the MQIFFT and the LQIFFT present systematic fractional bin index and amplitude errors that are strongly dependent [30] on the distance of the sinusoid's true fractional bin index \(\mathcal{K}\) from a bin center \(\lfloor\mathcal{K}\rceil\). \(\lfloor\cdot\rceil\) indicates the _rounding_ function, which returns its argument rounded to the nearest integer. This distance, denoted by \(\Delta_{\mathcal{K}}\), is defined as  
\[\Delta_{\mathcal{K}}=\mathcal{K}-\lfloor\mathcal{K}\rceil\,. \tag{31}\]  
By construction, bin estimates have odd symmetry around \(\Delta_{\mathcal{K}}=0\) and magnitude estimates have even symmetry around \(\Delta_{\mathcal{K}}=0\) (symmetries exploited in the previous error mean and worst-case definitions). The errors in estimates associated with the XIQFFT have all these properties as well. These can be seen in the _bias curves_ in Figure 2, which show \(e_{\mathcal{K}}\) and \(e_{\mathcal{X}}\) as a function of bin offset \(\Delta_{\mathcal{K}}\) for a length-4096 Hann window.  
Figure 2: Mean and worst-case of absolute value of errors for each method as a function of bin offset. (**a**) nearest bin method; (**b**) MQIFFT; (**c**) LQIFFT; (**d**) XQIFFT, \(p=0.1\); (**e**) nearest bin method; (**f**) MQIFFT; (**g**) LQIFFT; (**h**) XQIFFT, \(p=0.1\).  
Note that while it may sometimes (for certain windows) be possible to obtain closed form solutions for those two statistics, these solutions are generally very complex, even for simple smoothing windows. In the context of this article, we propose instead the following heuristic to evaluate those statistics numerically with arbitrary precision  
Finding the worst-case errors \(e_{\mathcal{K}}^{\text{w.c.}}\) and \(e_{\mathcal{X}}^{\text{w.c.}}\) requires finding the global maximum on the appropriate bias curve in Figure 2. Doing so numerically requires first to identify local maxima as the bias curves do not necessarily exhibit a single local maximum. We identify the regions containing a local maximum by computing the empirical error at a coarse sampling of the bin index space and finding the different local peaks that we associate with a single neighboring local maximum. We then run a Fibonacci search [37; 38; 39] (pp. 275-278, [40]) in the neighborhood of each found peak to find a local maximum. Using the Fibonacci search algorithm allows us to compute the location of that maximum up to a preset precision. We finally pick the global maximum as the largest local maximum.  
Finding the mean errors \(\text{E}[|e_{\mathcal{K}}(\mathcal{K})|]\) and \(\text{E}[|e_{\mathcal{X}}(\mathcal{K})|]\) requires to evaluate the integral of the error function across the range of all possible bin indices (i.e., between \(-0.5\) and \(+0.5\)). To do so numerically, we run an adaptive integration algorithm, the Simpson's quadrature algorithm [41], on the empirical error function. This algorithm allows us to control the desired numerical precision of the final result.  
For both methods, the procedure described above is aimed at obtaining estimates with arbitrary numerical precision (allowed by the Fibonacci search and the adaptive integration algorithms), while minimizing the number of empirical error function evaluations required during the extremum search (allowed by the Fibonacci search algorithm) and not requiring the derivation of a closed form error function.