Go to notes and algorithms index. ~ Go to home page.

auto_arc_tan() version 1


Integer Version. This version uses an arc tangent emulation formula to calculate the arc tangent values. For a similar function that uses a look-up table instead see version 5 of this function.


long auto_arc_tan( long v_dx, long v_dy )
   {
   // V1b
   // (C) Jon P, 1998.
   //
   // ** INTEGER VERSION **
   //
   // Returns an angle from 0 to 4095 (which equates to: 0 to nearly 360 degrees)
   // Angles are measured clockwise from vertical.
   //
   // This function will find the angle of orientation of a line specified
   // by parameters dx and dy.  Imagine this line as the hypotenuse of a
   // right angled triangle, dx and dy would be the lengths of the other
   // two sides.  Although dx and dy aren't strictly just lengths -- they
   // need a sign as well.... decide which end of your line is the origin
   // and calculate dx and dy as follows:-
   //   o Measure dx from the origin towards the right.
   //     (If the line extends left of the origin then dx should be negative.)
   //   o Measure dy from the origin upwards.
   //     (If the line extends below the origin then dy should be negative.)
   // In other words this function assumes that you're using a co-ordinate
   // system where Y values increase in the upward direction and X values
   // increase towards the right.
   //
   // "The 1 o'clock Example".
   // For the most straight-forward example imagine a triangle where
   // the hypotenuse was the hour hand of a clock pointing to 1 o'clock and
   // the second side was a vertical line from the centre of the clock face
   // upwards, and the third side was a horizontal line from the tip
   // of the hour hand projecting to the left intersecting with the
   // vertical side at right angles.  The origin is the centre of the clock
   // face.  This function would give the answer of the equivelent of 30
   // degrees (approx. 341 since the results are 0 to 4095 and not in degrees)
   //
   // This function works by emulating an 'arc tangent' (the opposite
   // of a tangent) where a ratio of dx against dy is created and from that
   // an angle is returned.  This ratio is called the gradient (for which
   // the variable name in this code is 'grad').
   //
   // Arc tangent is accurate and easy to use to calculate angles from 0 to
   // 45 degrees (gradients 0.0 to 1.0) but arc tangent becomes increasingly
   // problematic to use when the gradient becomes larger...and for 90 degrees
   // the ratio is infinity which computers cannot simply handle.
   // And to calculate angles approaching 90 degrees dy has to be measured
   // with very high accuracy to get good results.
   //
   // ... so that is one reason why this function only emulates arc tangent
   // for the ratios 0.0 to 1.0.  The other reason is that this region of the
   // arc tangent curve is much easier to emulate. Actually, for extremly
   // rough requirements this part of the arc tangent curve could even be
   // replaced with a straight line... but this function approximates the
   // curve.
   //
   // "The 2 o'clock Example".
   // However this function will calculate any angle from 0.0 to 360.0 degrees
   // (0 to 4095 inclusive).  So how can it do this when it only emulates
   // arc tangents that result in 0 to 45 degrees?  This is done by dividing
   // up the full 360 arc into 8 regions of 45 degrees each.  Up until #5# all
   // gradients and angles are for 0 to 45 degrees.   To demonstrate this
   // imagine the hour hand of the clock was pointing to 2 o'clock (instead
   // of 1 o'clock)... this would result in 60 degrees -- an angle greater
   // than 45 degrees and hence a gradient greater than 1.0.  However if we
   // were to do dx/dy (instead of dy/dx) the calculation would result in the
   // same gradient & angle as our first example: 30 degrees. The trick is then
   // to simply take 30 degrees away from 90 degrees... the result is the
   // required answer of 60 degrees (#6# is the line of code that would
   // perform this subtraction for this example).  All we've done is measured
   // the angle from 90 degrees anti-clockwise (instead of the usual way by
   // measuring from zero degrees clockwise).  When an angle is beyond 90
   // degrees the same process applies yet an additional amount 90, 180 or
   // 270 degrees has to be added to get the end result (this is done in the
   // code after #5#).
   //
   // The first process in dividing the calculateion into 8 regions is in the
   // 'if' statement marked #1#.  This makes sure the following code deals
   // only with a gradient between 0.0 and 1.0 (which is represented in
   // this code as 0 to 4096). It does this by deciding wheather to do
   // dy/dx (as would be for the 1 o'clock example) or dx/dy (as it would be
   // for the 2 o'clock example).  By the time execution has reached #2# or
   // #3# 'grad' holds a number between 0 and 4096 (representing a gradient
   // between 0.0 and 1.0).
   //
   // The line just after #2# and #3# sets a flag called 'half'.  If half
   // is set to 1 this will indicate to the code following #5# that your line
   // is in the first half of a quadrant.  In the 1 o'clock example
   // the #2# branch is executed, the gradient is calculated by dy/dx,
   // and 'half' is set to 1.
   //
   // The 4 lines of code just after the #4# comment is the formula that
   // emulates the arc tangent curve between 0.0 and 1.0.
   // There is little that can be explained for this since I found the formula
   // via trial and error by trying various curves and adjusting their
   // parameters until I found one that would approximate the arc tangent
   // curve.  A look-up table could instead be used here and would be slightly
   // quicker yet use more memory.  This section of code calculates 'ang' which
   // is an angle between 0 and 45 degrees (an actual value of 0 to 512).
   //
   // The code after #5# reconfigures the 0 to 45 degree angle stored in 'ang'
   // to be the full angle 0 to 360 required. as described above.
   //
   // A similar function to this one is supplied as a standard part of C
   // and is called atan2() (or something similar in other computer languages).
   // However that function is usually a floating-point function and I needed
   // an integer function because the platform I was working on didn't
   // have a floating-point processor... so that's why I wrote this function.
      
   long grad, ang, abs_dx, abs_dy, half;
   
   abs_dx = v_dx;
   if ( abs_dx < 0 ) abs_dx = 0 - abs_dx;
   abs_dy = v_dy;
   if ( abs_dy < 0 ) abs_dy = 0 - abs_dy;
   // ... make some unsigned versions of dx and dy.
   
   // find the gradient between 0 and one (or 0 and 4096 in this impelmentation)....
   if ( abs_dy <= abs_dx )  // #1#
      {
      if ( abs_dx == 0 )
         { //exact angle can be returned... also this prevents divide by 0 error...
         if ( v_dy >= 0 )
            return 0;
         else
            return 2048;
         }
      grad = abs_dy;
      grad <<= 12;     // This line multiplies by 4096 making 'grad' a fraction in 4096ths.
      grad /= abs_dx;  // 'grad' is 0 to 4096 (for gradients 0.0 to 1.0)
      // #2#
      half = 1;
      }
   else
      {
      if ( abs_dy == 0 )
         { //exact angle can be returned... also this prevents divide by 0 error...
         if ( v_dx >= 0 )
            return 1024;
         else
            return 3072;
         }
      grad = abs_dx;
      grad <<= 12;     // This line multiplies by 4096 making 'grad' a fraction in 4096ths.
      grad /= abs_dy;  // 'grad' is 0 to 4096 (for gradients 0.0 to 1.0)
      // #3#
      half = 0;
      }

   // #4#  This is the formula that emulates arc tangent for gradients
   //      0.0 to 1.0  ('grad' = 0 to 4096) ...
   ang = ( grad * 2409 ) >> 12;
   ang *= 4096 - ang;
   ang += grad * 3097;
   ang >>= 15;


   // #5#   Find which quadrant the line is in by examining the
   // signs of dx and dy.   Note that this function assumes Y ordinates
   // increase in the upwards direction, and X ordinates increase
   // in the right-hand direction.
   //
   if ( v_dx >= 0 )
      { // handle 0 to 180 degrees...
         
      if ( v_dy >= 0 )
         { // handle 0 to 90 degrees...
     
         if ( half ) ang = 1024 - ang;  // #6#
         }
      else
         { // handle over 90 degrees to 180 degrees...
            
         if ( half )
            ang += 1024;
         else
            ang = 2048 - ang;
         }
      }
   else // handle over 180 degrees to nearly 360 degrees...
      {
      if ( v_dy >= 0 )
         {  // handle 270 degrees to nearly 360 degrees...
            
         if ( half )
            ang += 3072;
         else
            ang = 4096 - ang;
         }
      else
         {  // handle over 180 degrees to nearly 270 degrees...
            
         if ( half )
            ang = 3072 - ang;
         else
            ang += 2048;
         }
      }
   
   return ang;
   }

I developed a new more acurate (but much slower) arc tangent emulation formula based on the curve "n * ( 1.0 - n )". The code for it is listed below and can be substituted for the section marked #4# above.

   ang = grad * grad;
   ang = ( grad * -9378 + ang ) / -5282 *
         ( 4096 - ( ( grad * -921 + ang ) / 3175 ) );
   ang /= 25885;
   ang += grad;
   ang >>= 3;