Goertzel Algorithm
9th September 2013I’ve always wondered how complicated it would be to write the code that listens to the numbers you dial on your phone, and works out what numbers you pressed.
I often enter my credit card details over the phone, and think: “If someone had a good enough microphone to record the sounds, they could probably get my information!”
The fancy name for this is: Dual-tone multi-frequency signalling. When you press the buttons, two sine waves are generated with different frequencies. The receiver then has to measure these two frequencies and can work out what number you’ve pressed. I do love the similarity of this with formants (the way we create the different vowels sounds).
Whilst an FFT would easily do the trick, after a bit of research, it’s obvious that the most common algorithm to use is the Goertzel algorithm.
My time pocket principals say I should have just done an FFT and been done with it. The point here wasn’t to get a working decoder, but to understand how it is done, and to actually write some code that I can have.
I decided to work in scilab. I was tempted by C/C++, but I fancied learning a bit more about scilab (and not fiddling with audio libraries).
I created a couple of simple functions to generate some pure “tones” to test the algorithm, as well as adding noise, and gaps.
There’s plenty of examples of the Goertzel code out there. I started off doing different overlapping partitions to sample the audio, but in the end took a far simpler approach.
I’m not exactly satisfied with the speed, or with the way I’ve implemented the divisions. Ironically the first time I wrote it was probably a better method, but the second is definitely faster. So I certainly made some mistakes. Creating a matrix for all the frequencies at all the sample points probably was a bit silly… but I ran out of time.
In addition to the Goertzel I looked at dealing with noise, silence and series of numbers. The algorithm will sample the frequencies (by default) every 0.05 seconds, and will adjust to the sample rate of the audio data. If there isn’t a strong signal, then it will assume there’s no number being pressed. If the value is the same as the previous check, no new number is written out.
I tested it by recording my pressing buttons on a phone, and it actually worked out great! I was quite surprised.
All in all a good way to spend a few hours!
SciLab code here:
(I would imagine it would be fairly easy to port to Matlab).
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 | // This creates an array of the goertzel's response. // Use: data x, at frequency f, sampled at points defined by at. function [res] = goertzelArr(x, f, sampleRate, at) a=0; b=0; c=0; val = 2 * cos(2 * %pi * f / sampleRate); res = zeros(1,length(at)); nextAt = 1; nextTrig = at(1); for i = 1:length(x), a = val * b - c + x(i); c = b; b = a; if i == nextTrig then res(nextAt) = (b*b + c*c - b*c*val); nextAt = nextAt + 1; b = 0; c = 0; if (nextAt > length(at) ) then return; end nextTrig = at(nextAt); end end endfunction // Main decoding function. // Call with: evalTone(AUDIO_VECTOR, SAMPLE_RATE_OF_AUDIO [, time between samples]) function [t] = evalToneAt(x,sampleRate,timebetween) // Default time between checks is 0.05 [out, inp] = argn(0); if inp < 3 then, timebetween = 0.05, end // Work out how many samples per check, and build array spaceSamples = timebetween * sampleRate; at = [spaceSamples:spaceSamples:length(x)]; // Relevant frequencies for number checkComb = [1209, 1336, 1477, 1633,697, 770, 852, 941]; chars = ['1','2','3','A','4','5','6','B','7','8','9','C','*','0','#','D']; // Create results array r = zeros(8,length(at)); // Get results from goertzels for i=1:8 r(i,:) = goertzelArr(x,checkComb(i),sampleRate,at); end lastVal = 0; t = ''; // Loop through measurements for i=1:length(at) // Find best tone pair col = 1; row = 1; colval = r(1,i); rowval = r(5,i); for j=2:4 if (r(j,i) > colval) then colval = r(j,i); col = j; end if (r(j+4,i) > rowval) then rowval = r(j+4,i); row = j; end end // This is the best candidate tone pair v = col + (row - 1) * 4; // If the magnitude doesn't satisfy contraints ignore. if (colval < spaceSamples | rowval < spaceSamples) then lastVal = 0; continue; end // If this is same tone as before. Ignore if (v ~= lastVal) then t = strcat([t,chars(v)]); lastVal = v; end end endfunction |
These functions create the sample tones:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 | // Assistant function to create the tones function [x] = soundChar(c, len, sampleRate) [out, inp] = argn(0); if inp < 3 then, sampleRate = 22050, end if inp < 2 then, len = 0.5, end checkCol = [1209, 1336, 1477, 1633]; checkRow = [697, 770, 852, 941]; chars = ['1','2','3','A';'4','5','6','B';'7','8','9','C';'*','0','#','D']; sig1 = 0; sig2 = 0; for i = 1:size(chars,1) for j = 1:size(chars,2) if (chars(i,j) == c) then sig2 = checkRow(i); sig1 = checkCol(j); end end end v1 = len * 2 * %pi; v2 = len * sampleRate; s1 = sin(linspace(0, v1*sig1, v2)); s2 = sin(linspace(0, v1*sig2, v2)); x = s1 + s2; endfunction // Create a tone sequence eg: genTone('0123456789#*ABCD',0.2,44100,0.05); // Will return an audio vector with each tone playing for 0.2 seconds, // with gaps in between of 0.05 seconds at a sample rate of 44100 function [x] = genTone(s, len, sampleRate, gap) [out, inp] = argn(0); if inp < 4 then, gap = 0.2, end if inp < 3 then, sampleRate = 22050, end if inp < 2 then, len = 0.5, end x = []; for i = 1:length(s) x = [x,soundChar(part(s,i), len, sampleRate),zeros(1,round(gap*sampleRate))]; end endfunction |