2009-08-24

FastMath :: atan2 lookup

Math.atan2() is very convenient, yet very slow. A lookup table will be roughly 8x faster. You lose accuracy, but heck, that's obvious.

   public static final float atan2Deg(float y, float x)
   {
      return FastMath.atan2(y, x) * DEG;
   }

   public static final float atan2DegStrict(float y, float x)
   {
      return (float) Math.atan2(y, x) * DEG;
   }

   public static final float atan2(float y, float x)
   {
      float add, mul;

      if (x < 0.0f)
      {
         if (y < 0.0f)
         {
            x = -x;
            y = -y;

            mul = 1.0f;
         }
         else
         {
            x = -x;
            mul = -1.0f;
         }

         add = -3.141592653f;
      }
      else
      {
         if (y < 0.0f)
         {
            y = -y;
            mul = -1.0f;
         }
         else
         {
            mul = 1.0f;
         }

         add = 0.0f;
      }

      float invDiv = ATAN2_DIM_MINUS_1 / ((x < y) ? y : x);

      int xi = (int) (x * invDiv);
      int yi = (int) (y * invDiv);

      return (atan2[yi * ATAN2_DIM + xi] + add) * mul;
   }


private static final int     ATAN2_BITS        = 7;

   private static final int     ATAN2_BITS2       = ATAN2_BITS << 1;
   private static final int     ATAN2_MASK        = ~(-1 << ATAN2_BITS2);
   private static final int     ATAN2_COUNT       = ATAN2_MASK + 1;
   private static final int     ATAN2_DIM         = (int) Math.sqrt(ATAN2_COUNT);

   private static final float   ATAN2_DIM_MINUS_1 = (ATAN2_DIM - 1);

   private static final float[] atan2             = new float[ATAN2_COUNT];

   static
   {
      for (int i = 0; i < ATAN2_DIM; i++)
      {
         for (int j = 0; j < ATAN2_DIM; j++)
         {
            float x0 = (float) i / ATAN2_DIM;
            float y0 = (float) j / ATAN2_DIM;

            atan2[j * ATAN2_DIM + i] = (float) Math.atan2(y0, x0);
         }
      }
   }

And some test code:
int dim = 633 * 2; // random number times 2

      float maxDiff = 0.0f;
      float sumDiff = 0.0f;

      for (int i = 0; i < dim * dim; i++)
      {
         float x = (float) ((i % dim) - (dim / 2)) / (dim / 2);
         float y = (float) ((i / dim) - (dim / 2)) / (dim / 2);
         float slow = (float) Math.atan2(y, x);
         float fast = FastMath.atan2(y, x);
         float diff = Math.abs(slow - fast);
         if (diff > maxDiff)
            maxDiff = diff;
         sumDiff += diff;
      }

      float avgDiff = sumDiff / (dim * dim);

      System.out.println("maxDiff=" + maxDiff); // 0.007858515
      System.out.println("avgDiff=" + avgDiff); // 0.002910751

12 comments:

  1. Hello there. I added this code to Eclipse, and it immediately complained about INV_ATAN2_DIM_MINUS_1 needing to be cast to int. However, with the cast added, INV_ATAN2_DIM_MINUS_1 is always zero, which I doubt is the intention. What should go there?
    Cheers
    Ryanm

    ReplyDelete
  2. You are correct. It should have been a float field. I really should stop modifying code in a webpage...

    ReplyDelete
  3. Also, should the line
    float invDiv = INV_ATAN2_DIM_MINUS_1 / ((x < y) ? y : x);
    be
    float invDiv = ((x < y) ? y : x) /INV_ATAN2_DIM_MINUS_1;
    ?

    ReplyDelete
  4. A lovely bit of code this, I'm measuring it at about 5 times faster that Math.atan2, even with the extra code to normalise down to the -1->1 range. Max error seems to be about 0.4 degrees, mean of 0.17 degrees.
    Cheers!

    ReplyDelete
  5. THAT WILL TEACH ME POSTING UNTESTED CODE =D

    I took the actual code from Eclipse, and optimized it a bit more. I provided a test-case so you can run (test...) it for yourself.

    ReplyDelete
  6. "even with the extra code to normalise down to the -1->1 range"

    No need to normalize the values.

    For example:
    Fastmath.atan2(137, -53) ~~ Math.atan2(137, -53)

    ReplyDelete
  7. "Max error seems to be about 0.4 degrees, mean of 0.17 degrees"

    Just change the "ATAN2_BITS=7" field, if you want to modify the accuracy.

    ReplyDelete
  8. I must have broken it somehow, I was getting array out-of-bounds exceptions when out of -1>1 range. Working lovely now though, at about 6x faster.

    ReplyDelete
  9. Keep in mind that due to:
    ATAN2_DIM_MINUS_1 / ((x < y) ? y : x)
    your will get out-of-bounds exceptions when both x and y are tiny, near Math.ulp().

    ReplyDelete
  10. <3 this code, I'm getting an average of 11 times faster than the default math version :D

    ReplyDelete
  11. Very useful and fast :) ! Under which license do you publish this code? I notice you use cc-by for posts but Apache would be nice for code :D ? Would like to integrate it into an existing apache project.

    ReplyDelete
  12. Regarding the index exceptions. Can't we solve this ala:

    int xi = Math.min((int) (x * invDiv), ATAN2_DIM);
    int yi = Math.min((int) (y * invDiv), ATAN2_DIM);

    Still identical fast for my usecase. (even slightly faster??)

    ReplyDelete