# DSP Trick: Magnitude Estimator

```From: Grant Griffin
Subject: DSP Trick: Magnitude Estimator
Date: 02 Oct 1999 00:00:00 GMT
THIS WORK IS PLACED IN THE PUBLIC DOMAIN
```

Name: Magnitude Estimator

Category: Algorithmic

Application: Signal strength measurement, AM demodulation, etc.

Advantages: This estimation algorithm is very quick compared to calculating magnitude to full precision using a square-root operation.

Introduction:

Given a complex number represented as "I + jQ", the magnitude is: sqrt(I**2 + Q**2). Unfortunately, the square-root operation in this formula is inherently slow and complicated to calculate, either in software or in hardware.

For applications which do not require a full-precision magnitude value, the use of a magnitude _estimation_ can save calculations. This well-known and widely used algorithm is a true gem of DSP because it provides considerable savings in calculation, at the cost of only a minimal loss of accuracy.

The Trick:

Given a complex number whose real part is "I" and whose imaginary part is "Q", the algorithm estimates the magnitude as:

Mag ~=Alpha * max(|I|, |Q|) + Beta * min(|I|, |Q|)

"Alpha" and "Beta" are two constants whose values can be chosen to trade among RMS error, peak error, and implementation complexity.

How It Works:

The absolute value operations folds the complex number into the range of 0-90 degrees, and the min, max operations further fold the complex number into the range of 0-45 degrees. Within this limited range, a linear combination of I and Q are a good approximation of magnitude.

Values for Alpha and Beta:

A program to demonstrate and test this algorithm is given at bottom. It includes a table of the most useful values for Alpha and Beta. The program prints out the following:

```=====================================================================
Alpha * Max + Beta * Min Magnitude Estimator

Name                  Alpha           Beta       Avg Err   RMS   Peak
(linear)  (dB)  (dB)
---------------------------------------------------------------------
Min RMS Err      0.947543636291 0.392485425092   0.000547 -32.6 -25.6
Min Peak Err     0.960433870103 0.397824734759  -0.013049 -31.4 -28.1
Min RMS w/ Avg=0 0.948059448969 0.392699081699   0.000003 -32.6 -25.7
1, Min RMS Err   1.000000000000 0.323260990000  -0.020865 -28.7 -23.8
1, Min Peak Err  1.000000000000 0.335982538000  -0.025609 -28.3 -25.1
1, 1/2           1.000000000000 0.500000000000  -0.086775 -20.7 -18.6
1, 1/4           1.000000000000 0.250000000000   0.006456 -27.6 -18.7
Frerking         1.000000000000 0.400000000000  -0.049482 -25.1 -22.3
1, 11/32         1.000000000000 0.343750000000  -0.028505 -28.0 -24.8
1, 3/8           1.000000000000 0.375000000000  -0.040159 -26.4 -23.4
15/16, 15/32     0.937500000000 0.468750000000  -0.018851 -29.2 -24.1
15/16, 1/2       0.937500000000 0.500000000000  -0.030505 -26.9 -24.1
31/32, 11/32     0.968750000000 0.343750000000  -0.000371 -31.6 -22.9
31/32, 3/8       0.968750000000 0.375000000000  -0.012024 -31.4 -26.1
61/64, 3/8       0.953125000000 0.375000000000   0.002043 -32.5 -24.3
61/64, 13/32     0.953125000000 0.406250000000  -0.009611 -31.8 -26.6
=====================================================================
```

Variations:

If you need a more accurate estimation than this algorithm provides, you can use some variation of it. For example, varying values of Alpha and Beta can be taken from a small lookup table, driven by the relative size of min and max values. Another possibility is to use this estimate as the "seed" of an iterative magnitude estimator.

Alternatives:

As an alternative, consider using the CORDIC algorithm, especially in hardware applications.

References:

This algorithm is described in Understanding Digital Signal Processing by Richard G. Lyons [Lyo97], in Digital Signal Processing in Communication Systems by Marvin E. Frerking [Fre94], and elsewhere.

Thanks to Clay S. Turner for providing some of the coefficients.

The Code:

```#include <math.h>
#include <stdio.h>
/*********************************************************************
*
*
* Name: mag_est.c
*
* Synopsis:
*
*   Demonstrates and tests the "Alpha * Min + Beta * Max" magnitude
*   estimation algorithm.
*
* Description:
*
*   This program demonstrates the "Alpha, Beta" algorithm for
*   estimating the magnitude of a complex number.  Compared to
*   calculating the magnitude directly using sqrt(I^2 + Q^2), this
*   estimation is very quick.
*
*   Various values of Alpha and Beta can be used to trade among RMS
*   error, peak error, and coefficient complexity.  This program
*   includes a table of the most useful values, and it prints out the
*   resulting RMS and peak errors.
*
* Copyright 1999  Grant R. Griffin
*
*                    The Wide Open License (WOL)
*
* Permission to use, copy, modify, distribute and sell this software
* and its documentation for any purpose is hereby granted without
* appear in all source copies. THIS SOFTWARE IS PROVIDED "AS IS"
* WITHOUT EXPRESS OR IMPLIED WARRANTY OF ANY KIND. See
*
*********************************************************************/
/********************************************************************/
double alpha_beta_mag(double alpha, double beta, double inphase,
{
/* magnitude ~= alpha * max(|I|, |Q|) + beta * min(|I|, |Q|) */
double abs_inphase = fabs(inphase);
return alpha * abs_inphase + beta * abs_quadrature;
} else {
return alpha * abs_quadrature + beta * abs_inphase;
}
}
/*********************************************************************/
double decibels(double linear)
{
#define SMALL 1e-20
if (linear <= SMALL) {
linear = SMALL;
}
return 20.0 * log10(linear);
}
/*********************************************************************/
void test_alpha_beta(char *name, double alpha, double beta,
int num_points)
{
#define PI 3.141592653589793
int ii;
double phase, real, imag, err, avg_err, rms_err;
double peak_err = 0.0;
double sum_err = 0.0;
double sum_err_sqrd = 0.0;
double delta_phase = (2.0 * PI) / num_points;
for (ii = 0; ii < num_points; ii++) {
phase = delta_phase * ii;
real = cos(phase);
imag = sin(phase);
err = sqrt(real * real + imag * imag)
- alpha_beta_mag(alpha, beta, real, imag);
sum_err += err;
sum_err_sqrd += err * err;
err = fabs(err);
if (err > peak_err) {
peak_err = err;
}
}
avg_err = sum_err / num_points;
rms_err = sqrt(sum_err_sqrd / num_points);
printf("%-16s %14.12lf %14.12lf  %9.6lf %4.1lf %4.1lf\n",
name, alpha, beta, avg_err, decibels(rms_err),
decibels(peak_err));
}
/*********************************************************************/
void main(void)
{
#define NUM_CHECK_POINTS 100000
typedef struct tagALPHA_BETA {
char *name;
double alpha;
double beta;
} ALPHA_BETA;
#define NUM_ALPHA_BETA 16
const ALPHA_BETA coeff[NUM_ALPHA_BETA] = {
{ "Min RMS Err",      0.947543636291, 0.3924854250920 },
{ "Min Peak Err",     0.960433870103, 0.3978247347593 },
{ "Min RMS w/ Avg=0", 0.948059448969, 0.3926990816987 },
{ "1, Min RMS Err",              1.0,     0.323260990 },
{ "1, Min Peak Err",             1.0,     0.335982538 },
{ "1, 1/2",                      1.0,      1.0 / 2.0  },
{ "1, 1/4",                      1.0,      1.0 / 4.0  },
{ "Frerking",                    1.0,            0.4  },
{ "1, 11/32",                    1.0,     11.0 / 32.0 },
{ "1, 3/8",                      1.0,      3.0 / 8.0  },
{ "15/16, 15/32",        15.0 / 16.0,     15.0 / 32.0 },
{ "15/16, 1/2",          15.0 / 16.0,      1.0 / 2.0  },
{ "31/32, 11/32",        31.0 / 32.0,     11.0 / 32.0 },
{ "31/32, 3/8",          31.0 / 32.0,      3.0 / 8.0  },
{ "61/64, 3/8",          61.0 / 64.0,      3.0 / 8.0  },
{ "61/64, 13/32",        61.0 / 64.0,     13.0 / 32.0 }
};
int ii;
printf("\n             Alpha * Max + Beta * Min Magnitude
Estimator\n\n");
printf("Name                  Alpha           Beta       Avg Err
RMS   Peak\n");
printf("                                                 (linear)
(dB)  (dB)\n");

printf("---------------------------------------------------------------------\n");
for (ii = 0; ii < NUM_ALPHA_BETA; ii++) {
test_alpha_beta(coeff[ii].name, coeff[ii].alpha, coeff[ii].beta,
1024);
}
}

```