0

Good morning!

I'm coding up a project on a microcontroller to read in some analog audio (specifically, the sound of someone whistling: a near perfect sine wave) and determine which piano note tones are present in the signal.

Although the FFT is the obvious first choice, the frequencies I wish to detect in the signal don't line up with the FFT frequencies. Furthermore, I only want to detect around 10-15 frequencies, not the whole spectrum.

So my thoughts quickly turn to the DTFT, which can be implemented with a matrix multiplication (here is my prototype MATLAB code. Details provided on request). However, surprisingly, it still takes more computation than an entire FFT!

This brings me to my question. I've made all the obvious google searches ("fast DTFT", "FFT for specific frequencies") without much luck. Does anyone here know of an efficient DTFT method?

Mahkoe
  • 131
  • 4
  • search for "pitch detection". Loads and loads of questions of this kind have already been asked! – Marcus Müller Mar 14 '17 at 13:23
  • 3
    regarding your question: yes, for single tones you can do Goertzel instead of FFT, but as soon as you do a significant amount of Goertzels, your complexity quickly exceeds that of a complete FFT. Anyway, neither is the way to go for pitch detection without post-processing. – Marcus Müller Mar 14 '17 at 13:24

1 Answers1

4

According to Oppenheim and Schafer's "Discrete Time Signal Processing", the Goertzel algorithm will be more efficient than the FFT in computing an N point DFT if less than $2 Log_2 N $ DFT coefficients are needed. So if you need to compute 15 frequency points, the Goertzel will be more efficient for total number of samples $N > 65536$. Anything below that and the FFT will be more efficient than Goertzel as Marcus has implied.

As you probably discovered, one point in the DFT requires N phase rotations and N complex additions (assuming a complex input). (As the DFT output is computed by simply multiplying your input signal through successive rotations and accumulating the result-- in this manner you are correlating to each output frequency of interest).

Each phase rotation requires 4 multiplications and 2 additions (if you have time to iterate with little multipliers, see the CORDIC algorithm as an alternative phase rotator that requires no multipliers). The phase rotations that are +/-j are simplified, since a rotation by j is done by simply changing the sign of Q and swapping I and Q, but for a larger FFT this will be an insignificant number of the total computation.

enter image description here

So therefore for each point using the DFT, there are approximately 4N total multiplications and 4N total additions. So for an N pt DFT that would be $4N^2$ of each - wow!

In comparison the FFT requires $2Nlog_2N$ multiplications and $2Nlog_2N$ additions. When considering a complete DFT, the FFT offers a dramatic reduction in required computations. (Pause here for some respectful bowing to Cooley-Tukey).

I would be interested in seeing how the Goertzel algorithm lands on the graph below but haven't muddled through those details yet.

enter image description here

Regarding your other challenge mentioned that your frequencies do not line up; for that I recommend windowing first to reduce the effect of distant frequencies on each tone you are measuring and then interpolate between the bins to more accurately estimate your frequency result for each tone. (A simple linear interpolation may be sufficient). The windowing will reduce the frequency resolution of each bin but will also reduce the sidelobes from the other bins, so you would need to confirm that your frequency separation for your tones of interest are sufficiently larger than the equivalent bandwidth of your window+FFT. For a good summary of equivalent bandwidths vs window refer to fred harris' paper on windowing: https://www.utdallas.edu/~cpb021000/EE%204361/Great%20DSP%20Papers/Harris%20on%20Windows.pdf

Dan Boschen
  • 50,942
  • 2
  • 57
  • 135
  • Off Topic: Nice, and I just wonder why do you want to replace the long time established nomenclature of the complex multiplication, with a phase rotation, and is there really a benefit from this change, for understanding of the DFT computation? And if one day someone comes and says "I'm bored of calling a multiplication with $W_N^{kn}$ as a complex multiplication, what can I else call it?" I would suggest, okay call it a correlation with a complex exponential or better call it an inner product with a base function that span some space,,, but what's the phasor rotation there? – Fat32 Mar 14 '17 at 22:38
  • I explain many DSP processes with "spinning phasors" as I found to a larger audience that isn't as disciplined in signal processing with complex signals it offers a more intuitive understanding (it is very visual). Mileage may vary ;). (I also did not make it up and have seen phase rotation used by others) – Dan Boschen Mar 14 '17 at 22:47
  • yes it's very visual, I agree :-) – Fat32 Mar 14 '17 at 22:51
  • Off Topic-2: I really need to ask you a question (albeit this isn't the place). I really cannot explain how a single wire antenna can pick-up electric field of an EMW and create a flowing current I of electrons through it into some receiving device, dispite the wire having no return path or any capacitance to temporarily pile up some electrons? You should know the answer I assume (but I do not want to go to electronics se to ask this, my laziness!) Do you know the answer? (without spinning phasors please! :-)) – Fat32 Mar 14 '17 at 23:01
  • @Fat32 Not sure how to explain it in a few sentences. Basically the antenna is resonant with the incoming wave, so the electrons are continuously going back and forth based on the length of the antenna. The current and voltage in the wire is different along the length at the ratio is the impedance at that position, the impedance is matched to the receiver at the feedpoint for maximum power transfer. So the radiated signal (electromagnetic wave) is picked up by the resonant antenna and then optimally transferred to the receiver. Lenz's law and Faraday's law may offer insight. – Dan Boschen Mar 14 '17 at 23:47
  • @Fat32, now back to spinning phasors- what I think really shows well what I am getting at is the DFT itself-- consider each row in the DFT matrix, the first row is all ones--so we multiply our signal by DC (a phasor magnitude 1 angle 0 that does not rotate with time) to get our first DFT bin. The next one is a phasor rotating clockwise at the minimum rate; once over the signal duration, the next one is twice the rate, etc. Now look at the correlation formula- The DFT is correlating our signal with each complex frequency, that I visualize as a spinning phasor. – Dan Boschen Mar 14 '17 at 23:51
  • Let us continue this discussion in chat. – Dan Boschen Mar 14 '17 at 23:52