FFT - Calculating exact frequency between frequency bins

Question

I am using a nice FFT library I found online to see if I can write a pitch-detection program. So far, I have been able to successfully let the library do FFT calculation on a test audio signal containing a few sine waves including one at 440Hz (I'm using 16384 samples as the size and the sample rate at 44100Hz).

The FFT output looks like:

433.356Hz - Real: 590.644 - Imag: -27.9856 - MAG: 16529.5
436.047Hz - Real: 683.921 - Imag: 51.2798 - MAG: 35071.4
438.739Hz - Real: 4615.24 - Imag: 1170.8 - MAG: 5.40352e+006
441.431Hz - Real: -3861.97 - Imag: 2111.13 - MAG: 8.15315e+006
444.122Hz - Real: -653.75 - Imag: 341.107 - MAG: 222999
446.814Hz - Real: -564.629 - Imag: 186.592 - MAG: 105355

As you can see, the 441.431Hz and 438.739Hz bins both show equally high magnitude outputs (the right-most numbers following "MAG:"), so it's obvious that the target frequency 440Hz falls somewhere between. Increasing the resolution might be one way to close in, but that would add to the calculation time.

How do I calculate the exact frequency that falls between two frequency bins?

UPDATE:

I tried out Barry Quinn's "Second Estimator" discussed on the DSPGuru website and got excellent results. The following shows the result for 440Hz square wave - now I'm only off by 0.003Hz!

Here is the code I used. I simply adapted this example I found, which was for Swift. Thank you everyone for your very valuable input, this has been a great learning journey :)

The methods you are looking for are called DESA (discrete energy separation algorithm). — DaBler
Thank you for your tip, I will look into this. I also found this post showing several solutions... — morezoom
I’d expect the most computationally-inexpensive and most reliable method is oversampling the FFT (i.e., zero-padding the input, especially to some power-of-two). But here are a few other techniques you can look at: gist.github.com/fasiha/957035272009eb1c9eb370936a6af2eb — Ahmed Fasih
I did wonder if I could take advantage of the harmonic content to help with the estimation. I'll definitely explore these Harmonic Product Spectrum techniques more. Thank you! — morezoom

Danijel Danijel · Accepted Answer · 2017-01-23T07:18:10

To calculate the "true" frequency, once I used parabola fit algorithm. It worked very well for my use case.

This is the way I followed in order to find the fundamental frequency:

Calculate DFT (WOLA).
Find peaks in your DFT bins.
Find Harmonic Product Spectrum. Not the most reliable nor precise, but this is a very easy way of finding your fundamental frequency candidates.
Based on peaks and HPS, use parabola fit algorithm to find fundamental pitch frequency (and amplitude if needed).

For example, HPS says the fundamental (strongest) pitch is concentrated in bin x of your DFT; if bin x belongs to the peak y, then parabola fit frequency is taken from the peak y and that is the pitch you were looking for.

If you are not looking for fundamental pitch, but exact frequency in any bin, just apply parabola fit for that bin.

Some code to get you started:

struct Peak
{
    float         freq     ; // Peak frequency calculated by parabola fit algorithm. 
    float         amplitude; // True amplitude.   
    float         strength ; // Peak strength when compared to neighbouring bins.         
    uint16_t      startPos ; // Peak starting position (DFT bin).
    uint16_t      maxPos   ; // Peak location (DFT bin).
    uint16_t      stopPos  ; // Peak stop position (DFT bin).
}; 

void calculateTrueFrequency( Peak & peak, float const bins, uint32_t const fs, DFT_Magnitudes mags )
{
    // Parabola fit:
    float a = mags[ peak.maxPos - 1 ];
    float b = mags[ peak.maxPos     ];
    float c = mags[ peak.maxPos + 1 ];

    float p   = 0.5f * ( a - c ) / ( a - 2.0f * b + c );
    float bin = convert<float>( peak.maxPos ) + p;

    peak.freq      = convert<float>( fs ) * bin / bins / 2;
    peak.amplitude = b - 0.25f + ( a - c ) * p;
}

FFT - Calculating exact frequency between frequency bins

2 Answers