Float vs fix math

[iMo]

As Jeoren3 wrote above for a real world app you have to consider the "environment" the mult lives in.
You have to feed the the the data in, make the multiply, write the result somewhere, increment some pointers or loop counters and/or to make other operations with the result for each mult.
So when you say "I need to make 1million mults within 5 secs" it says little about the reality.
Let us consider for example there is an overhead of 100 cycles per one mult (4 cycles).
So 104 cycles at 240MHz is 433ns and 1million of such events is 0.43secs.
Also look at the FPU spec of the ESP, it could be there is a MAC operation as well (mult and accumulate in one shot, like X = X + A*B) which is usually much faster than mults and additions made separately (when your math needs it).

[radiolistener]

The Goertzel algorithm is essentially a single-frequency FFT calculation.  It uses cos(), division, sqrt(), and lots of multiplication.  However, the cos() and division can be avoided by pre-computing the eight DTMF frequency coefficients and hard-coding them.  And sqrt() is only required to determine the magnitude of the signal to compare to the pre-set threshold, and magnitude squared would work just as well, compared to a squared threshold.

I remember, some time ago (more than 25 years ago  ), I implemented a DTMF detector using a simplified correlation approach with only 8-bit integer operations — basically XOR, ADD, and small lookup tables for precomputed sin/cos values with 1-bit samples (8 samples per byte). It worked surprisingly well even over noisy phone lines.

It was first implemented and debugged on a Z80, later I ported it to a PIC16F84. The entire code for detecting DTMF sequences and displaying them on an I2C LCD, including the display and keypad drivers — fit within the PIC16F84 1K program memory.

So, there is actually no need for sin(), cos(), sqrt(), or even multiplication for good DTMF detection.

[Peabody]

If you could make that available, or a link to the process you used, that would be very helpful.

[langwadt]

If you could make that available, or a link to the process you used, that would be very helpful.

sounds like some variation of the Walsh–Hadamard transform, basically a Fourier transform with square waves instead of sin/cos

[radiolistener]

The code probably still exists somewhere, but it’s stored on a CD-R, and I haven’t had a CD drive for years, so there’s no way to read it now. As far as I remember, the algorithm worked roughly as follows (I might have forgotten some details): the input signal from a comparator was fed into a GPIO pin, which sampled it at a fixed rate. Lookup tables with 1-bit samples for the detectable frequencies were generated using the same sample rate. The input bits were XORed with the corresponding sin/cos table values and accumulated. The number of accumulators therefore equaled twice the number of detectable frequencies. All processing was done in real time. After an accumulation period (about 50 ms, if I recall correctly), the algorithm evaluated the accumulated values — summing the squares of the sin and cos components (I don’t remember the exact evaluation operation to get magnitude, but everything was implemented entirely using ADD instruction), checking against a detection threshold, and selecting the two frequencies with the highest levels. The detected frequency pair was then mapped to the corresponding DTMF symbol.

[tggzzz]

The code probably still exists somewhere, but it’s stored on a CD-R, and I haven’t had a CD drive for years, so there’s no way to read it now. As far as I remember, the algorithm worked roughly as follows (I might have forgotten some details): the input signal from a comparator was fed into a GPIO pin, which sampled it at a fixed rate. Lookup tables with 1-bit samples for the detectable frequencies were generated using the same sample rate. The input bits were XORed with the corresponding sin/cos table values and accumulated. The number of accumulators therefore equaled twice the number of detectable frequencies. All processing was done in real time. After an accumulation period (about 50 ms, if I recall correctly), the algorithm evaluated the accumulated values — summing the squares of the sin and cos components (I don’t remember the exact evaluation operation to get magnitude, but everything was implemented entirely using ADD instruction), checking against a detection threshold, and selecting the two frequencies with the highest levels. The detected frequency pair was then mapped to the corresponding DTMF symbol.

Sounds like you had an IQ mixer, mixed the incoming waveform with local oscillator of an interesting frequency, calculated the magnitude of the output, and integrated.

[nfmax]

Not really related to the float/fixed topic, but here is my experience with a related problem.

I wanted to obtain the amplitude of a narrow-band signal: basically looking for echo peaks from an ultrasonic transducer. The 'proper' way to do this is to use either the Hilbert transform to generate the analytic signal, and thence the envelope, or to use a quadrature bandpass filter pair to create what is effectively the analytic signal over a limited range of frequencies. The first technique works best in the FFT domain, using relatively short records, the second is more suited to continuous operation with long-term signals.

There was a lot of low-frequency interference (basically amplifier/switch settling) that had to be filtered out to isolate the actual echo.

My fastest 'quick and dirty' technique was based on the observation that the sampling rate gave me about four samples per cycle of the transducer signal. What I did was to sum the absolute value of the differences between four consecutive pairs of samples (five samples in all). Dividing by two gives an estimate of the signal amplitude at frequencies around one-fourth of the sample rate, at the instant of the middle sample. You then move along to the next set of 5 samples, dropping the oldest and picking up the next.

This technique gives bandpass filtering and envelope calculation at the expense of eight fixed-point add/subtracts, four ANDs (to clear the sign bit) and a shift - and even the shift wasn't necessary for my needs. Address calculations probably require as many instructions as the actual computation (as they so often do).

To be fair, this is only a rough approximation to the true signal envelope, but for finding the echo peak, it was perfectly adequate.

[Simon]

The most intense part here is the DSP filter that will run a calculation.  It could be around 1M multiplies to do every 5 sec or so.  I am also realizing the heart of these equations for the DSP rely on floating math, so I may have little choice, or at least I could prob do some conversions. 

I appreciate and agree this is hard to really 'answer' without a specific use case.  Was just trying to get ideas in my head as I work if I should be thinking of ways to use integers for this and that. 

find out if you have a 64 bit result register for multiplication if you are using 32 bit values. I got a bit burnt on my M0+ MCU that claimed single cycle multiplication. I followed advice that using 32 bit values in my code would be faster as the processor would not be going to extra work to "ignore" upper bits. The result was that at 48 MHz it was taking 5+µs to do a multiplication (250 cycles). If you are using 16 bit values then you will be fine.

I don't know anything about DSP stuff but once I started working with an M4F based MCU I noticed a lot of talk about saturating math for DSP use. I don't think an ARM M4F is meant for actual DSP but they have included specific instructions that you don't normally expect to help. So might be worth looking at what the processor can actually do.

[Peabody]

@radiolistener, Perplexity knows of no publicly available AVR code using the Walsh–Hadamard transform to detect DTMF (if that's what you used).  But in general, if you used a comparator to convert the incoming audio to square waves, how did you deal with the audio level?  How did you prevent normal speech from triggering DTMF detections?  Well, maybe we're moving a bit off-topic with this.  I did get my fixed point Goertzel to work, including two >>14 operations.  Now I need to see how the float version compares.

[SteveThackery]

Presumably not possible in the OP's scenario, but has anyone used the DTMF encoder and decoder chips that are available?  Any comments on their performance, ease of use, etc?

I ask because I'm considering DTMF for a low bit rate project I have in mind.

[radiolistener]

Perplexity knows of no publicly available AVR code using the Walsh–Hadamard transform to detect DTMF (if that's what you used).  But in general, if you used a comparator to convert the incoming audio to square waves, how did you deal with the audio level?  How did you prevent normal speech from triggering DTMF detections?

The algorithm operated on a binary, discretized representation of the signal — essentially it responded to the timing of zero crossings rather than to amplitude. So the signal level wasn’t really relevant. In practice it was used for remote control of a telephone answering machine, and I never noticed it being triggered by normal speech, although I suppose certain sounds could have caused false detections if you tried hard enough.

As for the algorithm itself, I originally found it by disassembling the firmware of that answering machine. I’m not sure what the underlying math was — ChatGPT suggested it might be a “1-bit correlation detector with a Bayesian decision criterion under unknown a priori frequency probabilities”.

I remember that the manual for the answering machine mentioned that DTMF recognition was based on a “new decision-making algorithm under conditions of a priori uncertainty”. There was no other information about the algorithm provided in the documentation.

[tggzzz]

As for the algorithm itself, I originally found it by disassembling the firmware of that answering machine. I’m not sure what the underlying math was — ChatGPT suggested it might be a “1-bit correlation detector with a Bayesian decision criterion under unknown a priori frequency probabilities”

That answer response is content-free word salad.

Not surprising, given the source.

[radiolistener]

That answer response is content-free word salad.

Not surprising, given the source.

I tried feeding your criticism back to ChatGPT, and it replied:
“The phrasing does sound overcomplicated. In simpler terms, it was just a 1-bit correlation detector: the signal from the comparator was compared against precomputed sine/cosine patterns for each frequency, and the two with the highest correlation were selected as the detected DTMF tones.”

I also tried to mention the similarity to an XOR-based phase detector, and it responded:
“In more general terms, this approach can be regarded as a form of nonlinear correlation or binary phase detection.”

That still doesn’t fully clarify the underlying math, though. Perhaps someone familiar with the theory behind such methods could shed some light on the mathematical basis. I used this algorithm quite a bit back in the day for frequency detection on low-end processors, and it performed surprisingly well. Intuitively, I see a connection to the Fourier transform and correlation, but I’m not sure how exactly this approach would be formally classified.

[nctnico]

Presumably not possible in the OP's scenario, but has anyone used the DTMF encoder and decoder chips that are available?  Any comments on their performance, ease of use, etc?

Back in the old days I have used these chips. Typically from Mitel but those are probably long obsolete.

[tggzzz]

That answer response is content-free word salad.

Not surprising, given the source.

I tried feeding your criticism back to ChatGPT, and it replied:
“The phrasing does sound overcomplicated. In simpler terms, it was just a 1-bit correlation detector: the signal from the comparator was compared against precomputed sine/cosine patterns for each frequency, and the two with the highest correlation were selected as the detected DTMF tones.”

I also tried to mention the similarity to an XOR-based phase detector, and it responded:
“In more general terms, this approach can be regarded as a form of nonlinear correlation or binary phase detection.”

That still doesn’t fully clarify the underlying math, though. Perhaps someone familiar with the theory behind such methods could shed some light on the mathematical basis. I used this algorithm quite a bit back in the day for frequency detection on low-end processors, and it performed surprisingly well. Intuitively, I see a connection to the Fourier transform and correlation, but I’m not sure how exactly this approach would be formally classified.

Keep in feeding criticism back to it, and it will generate more answers. Rinse and repeat.

I wonder if it would eventually blunder upon a respectable answer, and whether anyone would be able to determine that that specific response was accurate and useful.

[jheissjr]

The Cortex-M85 has the new Helium floating point and fixed point extension. An ARM engineer wrote a great reference guide about how it works and how to use it. https://armkeil.blob.core.windows.net/developer/Files/pdf/ebook/arm-helium-technology-mve.pdf

[nctnico]

The code probably still exists somewhere, but it’s stored on a CD-R, and I haven’t had a CD drive for years, so there’s no way to read it now. As far as I remember, the algorithm worked roughly as follows (I might have forgotten some details): the input signal from a comparator was fed into a GPIO pin, which sampled it at a fixed rate. Lookup tables with 1-bit samples for the detectable frequencies were generated using the same sample rate. The input bits were XORed with the corresponding sin/cos table values and accumulated. The number of accumulators therefore equaled twice the number of detectable frequencies. All processing was done in real time. After an accumulation period (about 50 ms, if I recall correctly), the algorithm evaluated the accumulated values — summing the squares of the sin and cos components (I don’t remember the exact evaluation operation to get magnitude, but everything was implemented entirely using ADD instruction), checking against a detection threshold, and selecting the two frequencies with the highest levels. The detected frequency pair was then mapped to the corresponding DTMF symbol.

Sounds like you had an IQ mixer, mixed the incoming waveform with local oscillator of an interesting frequency, calculated the magnitude of the output, and integrated.

Yep. XOR function works as a mixer and having two mixers fed with two LO frequencies 90 degrees out of phase gets you an IQ demodulator.

[tggzzz]

The code probably still exists somewhere, but it’s stored on a CD-R, and I haven’t had a CD drive for years, so there’s no way to read it now. As far as I remember, the algorithm worked roughly as follows (I might have forgotten some details): the input signal from a comparator was fed into a GPIO pin, which sampled it at a fixed rate. Lookup tables with 1-bit samples for the detectable frequencies were generated using the same sample rate. The input bits were XORed with the corresponding sin/cos table values and accumulated. The number of accumulators therefore equaled twice the number of detectable frequencies. All processing was done in real time. After an accumulation period (about 50 ms, if I recall correctly), the algorithm evaluated the accumulated values — summing the squares of the sin and cos components (I don’t remember the exact evaluation operation to get magnitude, but everything was implemented entirely using ADD instruction), checking against a detection threshold, and selecting the two frequencies with the highest levels. The detected frequency pair was then mapped to the corresponding DTMF symbol.

Sounds like you had an IQ mixer, mixed the incoming waveform with local oscillator of an interesting frequency, calculated the magnitude of the output, and integrated.

Yep. XOR function works as a mixer and having two mixers fed with two LO frequencies 90 degrees out of phase gets you an IQ demodulator.

I wonder how much time and prompting it would take to get an LLM to emit that 

[langwadt]

The code probably still exists somewhere, but it’s stored on a CD-R, and I haven’t had a CD drive for years, so there’s no way to read it now. As far as I remember, the algorithm worked roughly as follows (I might have forgotten some details): the input signal from a comparator was fed into a GPIO pin, which sampled it at a fixed rate. Lookup tables with 1-bit samples for the detectable frequencies were generated using the same sample rate. The input bits were XORed with the corresponding sin/cos table values and accumulated. The number of accumulators therefore equaled twice the number of detectable frequencies. All processing was done in real time. After an accumulation period (about 50 ms, if I recall correctly), the algorithm evaluated the accumulated values — summing the squares of the sin and cos components (I don’t remember the exact evaluation operation to get magnitude, but everything was implemented entirely using ADD instruction), checking against a detection threshold, and selecting the two frequencies with the highest levels. The detected frequency pair was then mapped to the corresponding DTMF symbol.

Sounds like you had an IQ mixer, mixed the incoming waveform with local oscillator of an interesting frequency, calculated the magnitude of the output, and integrated.

Yep. XOR function works as a mixer and having two mixers fed with two LO frequencies 90 degrees out of phase gets you an IQ demodulator.

sorta related, a Tayloe mixer

[SiliconWizard]

As others have said, Goertzel has been used a lot for DTMF detection, and it's adequate when used properly. So, I agree that claiming it's "useless" is, at best, odd.
Now, it's not a magical tool and will work better with a "clean" input signal, ideally with a reasonably "normalized" amplitude, so some amount of pre-filtering and a basic AGC may be required, unless those properties are guaranteed with the source signal. So, yes, depending on the source signal, some other approaches may turn out more effective.

Not too surprised "cargo cult" engineering may rear its head, as in "I've been bitten once, so I'm never gonna use this again".

[Peabody]

The question of AGC is an interesting one.  In playing with Goertzel on a Nano, I've found that with 5V ADC values, you pretty much need the signal to fill up most of that 5V range (centered on 2.5V) for Goertzel to work well.  And even after picking the highest level row and column frequencies, both of their levels must also be greater than some predetermined threshold.  If the threshold is too low, there would be no way to prevent voice audio from generating false positives.  So it seems that DTMF depends on the tones being relatively loud.  I'm not sure AGC works for that.

My project is for a home VOIP system such as Ooma and MagicJack.  I'll need to detect DTMF both from another extension phone on the same line and from someone at the other end of the conversation.  Not surprisingly, the levels produced by those sources are very different.  But it's still just two different levels.  Well, there is a third non-DTMF signal I need to detect, which is dial tone.  So three different levels.  But I can switch between three gain levels depending on what I'm listening for, so it doesn't have to be automatic gain control, it just has to be gain control.  At least that's my plan.

[nctnico]

As others have said, Goertzel has been used a lot for DTMF detection, and it's adequate when used properly. So, I agree that claiming it's "useless" is, at best, odd.

It is not. Look at the relatively wide frequency and detection signal level specs of DTMF and how you want to meet those specs with Goertzel. You simply can't because the frequency response of Goertzel is like a parabole. At frequencies which are 'off' but still in spec, you will need signals which are a lot stronger. In addition to that, Goertzel will be sensitive for signals outside the frequency range if they are strong enough. Actually, using Goertzel is the definition of cargo cult engineering. People read a paper written by a few students using Goertzel for DTMF decoding and then think using Goertzel is the defacto standard for DTMF decoding without considering which specs actually need to be met.

[radiolistener]

I attempted to reconstruct the frequency detection method I mentioned earlier. The following Octave code demonstrates it.

For simplicity, it uses the CCITT No.5 multi-frequency signaling scheme, which operates similarly to DTMF but with a different set of frequencies.

Code: [Select]

% Binary frequency detector simulation with multiple tones

Fs = 48000;                  % sample rate
pause_len = 0.1;             % 100 ms pause
tone_len  = 0.1;             % 100 ms per tone
freqs = [700 900 1100 1300 1500 1700];  % CCITT5 - Signaling System No. 5 frequencies

% --- Generate composite test signal ---
sig = [];
% signle tones with pause
for k = 1:length(freqs)
    % 100 ms pause before first tone only
    if k == 1
        sig = [sig, zeros(1, round(pause_len*Fs))];
    end

    % Tone segment
    t = (0:1/Fs:tone_len - 1/Fs);
    sig = [sig, sin(2*pi*freqs(k)*t)];

    % Pause after each tone
    sig = [sig, zeros(1, round(pause_len*Fs))];
end
% two tones
for k = 2:length(freqs)
    t = (0:1/Fs:tone_len - 1/Fs);
    tone = (sin(2*pi*freqs(k)*t) + sin(2*pi*freqs(k-1)*t)) / 2;
    sig = [sig, tone];
end


t_sig = (0:length(sig)-1)/Fs;

% --- Comparator simulation ---
bin_sig = uint8(sig > 0) & 1;


% --- Precompute binary sin/cos references ---
win = round(0.01*Fs);        % 10 ms window
for k = 1:length(freqs)
    f = freqs(k);
    ref_sin{k} = uint8(sin(2*pi*f*(0:win-1)/Fs) > 0) & 1;  % square wave sin/cos ref table
    ref_cos{k} = uint8(cos(2*pi*f*(0:win-1)/Fs) > 0) & 1;
end


% --- Sliding detection ---
step = round(win/4);
nsteps = floor((length(bin_sig)-win)/step);
det = zeros(length(freqs), nsteps);

for i = 1:nsteps
    pos = (i-1)*step + 1;
    seg = bin_sig(pos:pos+win-1);
    for k = 1:length(freqs)
        sxor = bitxor(seg, ref_sin{k});
        cxor = bitxor(seg, ref_cos{k});
        sxor = 1 - 2*sxor;    % 0/1 => -1/+1
        cxor = 1 - 2*cxor;
        xs = sum(sxor);  % count bit matches
        xc = sum(cxor);
        det(k,i) = min((xs^2 + xc^2)/256, 256);
    end
end

t_det = (0:(nsteps-1))*step/Fs;

% --- Plot result ---
figure(1); clf;
subplot(2,1,1);
plot(t_sig, sig);
xlabel('Time (s)');
ylabel('Signal');
title('Test signal');
grid on; grid minor;

subplot(2,1,2);
plot(t_det, det');
xlabel('Time (s)');
ylabel('Detection amplitude');
legend(arrayfun(@(f) sprintf('%d Hz', f), freqs, 'UniformOutput', false));
title('Frequency detector output');
grid on; grid minor;

Note: min((xs^2 + xc^2)/256, 256) can be replaced with 256 byte lookup table and ADD, so there is no need for multiplication operation, just XOR and signed ADD

Sample rate can be decreased to 12 kHz, so the reference sin/cos table for 10 ms window will need 120*2 bits = 30 bytes per frequency (180 bytes for CCITT5 or 240 bytes for DTMF).

For optimization, you can accumulate eight bits of the input signal and then perform XOR and ADD operations for all eight samples in a single step.

[Peabody]

Thanks for posting the code.  I'll look it over.

[radiolistener]

Also you can use Chebyshev distance as a simple magnitude approximation without MUL and SQRT operations, so there is no need for 256 byte LUT:

Code: [Select]

for i = 1:nsteps
    pos = (i-1)*step + 1;
    seg = bin_sig(pos:pos+win-1);
    for k = 1:length(freqs)
        sxor = bitxor(seg, ref_sin{k});
        cxor = bitxor(seg, ref_cos{k});
        sxor = 1 - 2*sxor;    % 0/1 => -1/+1
        cxor = 1 - 2*cxor;
        xs = sum(sxor);  % count bit matches
        xc = sum(cxor);

        %mag = sqrt(xs^2 + xc^2);       % Classic Magnitude
        %mag = xs^2 + xc^2;             % Squared Magnitude
        %mag = abs(xs) + abs(xc);       % Manhattan norm (L1)
        mag = max(abs(xs), abs(xc));   % Chebyshev norm (L∞)
       
        det(k,i) = min(mag, 256);
    end
end

This allows the DTMF detector to operate entirely with integer arithmetic, eliminating the need of ADC, amplitude scaling, floating-point calculations, and costly MUL, SQRT, or trigonometric operations, relying solely on integer XOR and ADD.

On 32-bit MCU you can process 32 samples at once with just 3-4 processor instructions.

Result: