Mogens Kjaer wrote:
...
Padé approximation:
http://list.dprg.org/archive/2007-January/030202.html
You could try this piece of code:
static double atan2x(double y, double x)
{
double z;
if(x==0)
{
if(y<0) return -M_PI_2;
else if(y>0) return M_PI_2;
else return 0;
}
else
{
z = y/x;
if(z>=-1 && z<=1)
{
z=z*(15.0+4.0*z*z)/(15.0+9.0*z*z);
if(x>0) return z;
else
{
if(y>=0) return M_PI+z;
else return z-M_PI;
}
}
else
{
z=x/y;
z=M_PI_2-z*(15.0+4.0*z*z)/(15.0+9.0*z*z);
if(y<=0) return z-M_PI;
else return z;
}
}
}
Testing it with:
maxe=0;
for(i=0;i<10001;i++)
for(j=0;j<10001;j++)
{
x=0.0001*i-0.5;
y=0.0001*j-0.5;
a1=atan2(y,x);
a2=atan2x(y,x);
e=fabs(a1-a2);
if(e>maxe) maxe=e;
}
printf("maxe: %g\n", maxe);
gives maxe: 0.0062685
Timing (Xeon 5160 3 GHz) for the loop with only one function call:
atan2: 7.028u
atan2x: 1.586u
I'm not an expert on this, so use it at your own risk :-)
Mogens
--
Mogens Kjaer, Carlsberg A/S, Computer Department
Gamle Carlsberg Vej 10, DK-2500 Valby, Denmark
Phone: +45 33 27 53 25, Fax: +45 33 27 47 08
Email: mk@xxxxxx Homepage: http://www.crc.dk
