Mailing Lists: Apple Mailing Lists
Image of Mac OS face in stamp
 
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: Fast atan2


On Jan 20, 2005, at 11:46 AM, rick wrote:


http://www.dspguru.com/comp.dsp/tricks/alg/fxdatan2.htm

On 18 Jan 2005, at 23:56, Robert Purves wrote:

Shark indicates that an important bottleneck in my app is
angle = atan2( y, x );
Does anyone know of a fast approximate atan2 algorithm (like those delightfully sleazy hacks for sqrt and its reciprocal)?
I don't need anything like full double, or even single, precision. An absolute error of 1e-4 would be acceptable.

Robert P. 

I don't know that that method is really saving you much time:

ollmia:/tmp iano$ gcc -O3 main4.c -Wmost
ollmia:/tmp iano$ ./a.out
best libm time for 1000 calls: 3.06011e-06 seconds (792544)
best cheesy time for 1000 calls: 2.97011e-06 seconds (1.56224e+06)
ollmia:/tmp iano$ cat main4.c
#include <math.h>
#include <stdint.h>
#include <mach/mach_time.h>
#include <stdio.h>
#include <string.h>

static double MySubtractTime( uint64_t endTime, uint64_t startTime )
{
    uint64_t    difference = endTime - startTime;
    static double conversion = 0.0;

    if( 0.0 == conversion )
    {
        mach_timebase_info_data_t       timebase;
        kern_return_t                   err;

        memset( &timebase, 0, sizeof( timebase ));

        err = mach_timebase_info( &timebase );
        if( 0 != err )
            return err;

conversion = 1e-9 * (double) timebase.numer / (double) timebase.denom;
}

    return (double) difference * conversion;
}

float arctan2(float y, float x)
{
   float  angle;
   float coeff_1 = M_PI/4;
   float coeff_2 = 3*coeff_1;
   float abs_y = fabsf(y)+1e-10;      // kludge to prevent 0/0 condition
   if (x>=0)
   {
      double r = (x - abs_y) / (x + abs_y);
      angle = coeff_1 - coeff_1 * r;
   }
   else
   {
      float 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);
}

int main( void )
{
        uint64_t        startTime, endTime;
        double          currentTime, bestTime;
        int i, j;
        float sum = 0.0;

bestTime = 1e20;
for( i = 0; i < 1000; i++ )
{
startTime = mach_absolute_time();
for( j = 0; j < 1000; j++ )
sum += atan2f(1.0f, 1.0f );
endTime = mach_absolute_time();
currentTime = MySubtractTime( endTime, startTime );
if( currentTime < bestTime )
bestTime = currentTime;
}
printf("best libm time for 1000 calls: %g seconds (%g)\n", bestTime, sum );

bestTime = 1e20;
for( i = 0; i < 1000; i++ )
{
startTime = mach_absolute_time();
for( j = 0; j < 1000; j++ )
sum += arctan2(1.0f, 1.0f );
endTime = mach_absolute_time();
currentTime = MySubtractTime( endTime, startTime );
if( currentTime < bestTime )
bestTime = currentTime;
}
printf("best cheesy time for 1000 calls: %g seconds (%g)\n", bestTime, sum );


        return 0;
}

Note: I didn't spend any time to check to make sure that I translated the pseudo-code presented in that article correctly.

Ian

_______________________________________________
Do not post admin requests to the list. They will be ignored.
PerfOptimization-dev mailing list      (email@hidden)
Help/Unsubscribe/Update your Subscription:
http://lists.apple.com/mailman/options/perfoptimization-dev/email@hidden

This email sent to email@hidden
References: 
 >Fast atan2 (From: Robert Purves <email@hidden>)
 >Re: Fast atan2 (From: rick <email@hidden>)



Visit the Apple Store online or at retail locations.
1-800-MY-APPLE

Contact Apple | Terms of Use | Privacy Policy

Copyright © 2007 Apple Inc. All rights reserved.