DSP Trick: Fixed-Point Atan2 With Self Normalization

Subject: DSP Trick:  Fixed-pt. Atan2 with Self Normalization
From: Jim Shima
Date: 1999/04/23
Newsgroups: comp.dsp
THIS WORK IS PLACED IN THE PUBLIC DOMAIN

Name: Fast fixed-pt. atan2 calculation with self normalization

Category: Algorithmic

Application: Used when one wants to compute the 4-quadrant arctangent of a complex number (or any number with x-y coordinates) with a self-normalizing function.

Example Applications: digital FM demodulation, phase angle computations

Advantages: Fast

Introduction:

Computing a 4-quadrant arctangent on DSPs has been the subject of many discussions. Several techniques such as table lookup and polynomial expansion are well known.

In fixed-point DSPs, some normalization of the complex number may be necessary, effectively implementing a hard limiter or amplitude invariance function.

In fact, computing:

theta = atan(y/x)

includes the necessary normalization, but in a fixed-pt. DSP, the division can result in values outside the fixed-pt. range of [-1,1).

Also, for certain applications, such as digital FM demodulation, any amplitude fluctuations must be removed before computing the phase angle.

The Trick:

I computed a self-normalizing ratio depending on the quadrant that the complex number resides in. For a complex number z, let x=Re(z) and y=Im(z).

For a complex number in quadrant I, compute the ratio:

    x-y
r = ---     (1)
    x+y

To get the phase angle, compute:

theta1 = pi/4 - pi/4*r (2)

Likewise, if the complex number resides in quadrant II, compute the ratio:

    x+y
r = ---     (3)
    y-x

And to get the quadrant II phase angle, compute:

theta2=3*pi/4 - pi/4*r (4)

If it turns out that the complex number was really in quad IV instead of quad I, just negate the answer resulting from (2).

Likewise, do the same if the number was in quad III instead of quad II. By doing this, you have a 4-quadrant arctan function.

The max error using equations (2) or (4) is a little less than 0.07 rads (only at a few angles though). The accuracy of the estimator is actually quite good considering using a 1st-order polynomial to estimate the phase angle.

If you use a higher degree polynomial, it turns out that the even powers of the poly will disappear (due to the odd function), thus relaxing some of the computational load.

QUICK NOTE FOR BETTER ACCURACY:

To obtain better accuracy (a max error of .01 rads), one can replace equations (2) and (4) with:

theta1=0.1963 * r^3 - 0.9817 * r + pi/4 (2a)
theta2=0.1963 * r^3 - 0.9817 * r + 3*pi/4 (4a)

equations (2a) or (4a) can be computed using 2 MACs, which does not involve much more computation for a 7x increase in accuracy.

Here is some C pseudocode (not optimized in any fashion) using equations (1)-(4):

//-----------------------------------------------
// Fast arctan2
float arctan2(float y, float x)
{
   coeff_1 = pi/4;
   coeff_2 = 3*coeff_1;
   abs_y = fabs(y)+1e-10      // kludge to prevent 0/0 condition
   if (x>=0)
   {
      r = (x - abs_y) / (x + abs_y);
      angle = coeff_1 - coeff_1 * r;
   }
   else
   {
      r = (x + abs_y) / (abs_y - x);
      angle = coeff_2 - coeff_1 * r;
   }
   if (y < 0)
   return(-angle);     // negate if in quad III or IV
   else
   return(angle);
}