13e12c5d1SDavid du Colombier /* 23e12c5d1SDavid du Colombier floating-point arctangent 33e12c5d1SDavid du Colombier 43e12c5d1SDavid du Colombier atan returns the value of the arctangent of its 53e12c5d1SDavid du Colombier argument in the range [-pi/2,pi/2]. 63e12c5d1SDavid du Colombier 73e12c5d1SDavid du Colombier atan2 returns the arctangent of arg1/arg2 83e12c5d1SDavid du Colombier in the range [-pi,pi]. 93e12c5d1SDavid du Colombier 103e12c5d1SDavid du Colombier there are no error returns. 113e12c5d1SDavid du Colombier 123e12c5d1SDavid du Colombier coefficients are #5077 from Hart & Cheney. (19.56D) 133e12c5d1SDavid du Colombier */ 143e12c5d1SDavid du Colombier 153e12c5d1SDavid du Colombier #include <u.h> 163e12c5d1SDavid du Colombier #include <libc.h> 173e12c5d1SDavid du Colombier 183e12c5d1SDavid du Colombier #define sq2p1 2.414213562373095048802e0 193e12c5d1SDavid du Colombier #define sq2m1 .414213562373095048802e0 203e12c5d1SDavid du Colombier #define p4 .161536412982230228262e2 213e12c5d1SDavid du Colombier #define p3 .26842548195503973794141e3 223e12c5d1SDavid du Colombier #define p2 .11530293515404850115428136e4 233e12c5d1SDavid du Colombier #define p1 .178040631643319697105464587e4 243e12c5d1SDavid du Colombier #define p0 .89678597403663861959987488e3 253e12c5d1SDavid du Colombier #define q4 .5895697050844462222791e2 263e12c5d1SDavid du Colombier #define q3 .536265374031215315104235e3 273e12c5d1SDavid du Colombier #define q2 .16667838148816337184521798e4 283e12c5d1SDavid du Colombier #define q1 .207933497444540981287275926e4 293e12c5d1SDavid du Colombier #define q0 .89678597403663861962481162e3 303e12c5d1SDavid du Colombier 313e12c5d1SDavid du Colombier 323e12c5d1SDavid du Colombier /* 333e12c5d1SDavid du Colombier xatan evaluates a series valid in the 34*219b2ee8SDavid du Colombier range [-0.414...,+0.414...]. (tan(pi/8)) 353e12c5d1SDavid du Colombier */ 363e12c5d1SDavid du Colombier 373e12c5d1SDavid du Colombier static 383e12c5d1SDavid du Colombier double xatan(double arg)393e12c5d1SDavid du Colombierxatan(double arg) 403e12c5d1SDavid du Colombier { 413e12c5d1SDavid du Colombier double argsq, value; 423e12c5d1SDavid du Colombier 433e12c5d1SDavid du Colombier argsq = arg*arg; 443e12c5d1SDavid du Colombier value = ((((p4*argsq + p3)*argsq + p2)*argsq + p1)*argsq + p0); 453e12c5d1SDavid du Colombier value = value/(((((argsq + q4)*argsq + q3)*argsq + q2)*argsq + q1)*argsq + q0); 463e12c5d1SDavid du Colombier return value*arg; 473e12c5d1SDavid du Colombier } 483e12c5d1SDavid du Colombier 493e12c5d1SDavid du Colombier /* 503e12c5d1SDavid du Colombier satan reduces its argument (known to be positive) 513e12c5d1SDavid du Colombier to the range [0,0.414...] and calls xatan. 523e12c5d1SDavid du Colombier */ 533e12c5d1SDavid du Colombier 543e12c5d1SDavid du Colombier static 553e12c5d1SDavid du Colombier double satan(double arg)563e12c5d1SDavid du Colombiersatan(double arg) 573e12c5d1SDavid du Colombier { 583e12c5d1SDavid du Colombier 593e12c5d1SDavid du Colombier if(arg < sq2m1) 603e12c5d1SDavid du Colombier return xatan(arg); 613e12c5d1SDavid du Colombier if(arg > sq2p1) 623e12c5d1SDavid du Colombier return PIO2 - xatan(1/arg); 633e12c5d1SDavid du Colombier return PIO2/2 + xatan((arg-1)/(arg+1)); 643e12c5d1SDavid du Colombier } 653e12c5d1SDavid du Colombier 663e12c5d1SDavid du Colombier /* 673e12c5d1SDavid du Colombier atan makes its argument positive and 683e12c5d1SDavid du Colombier calls the inner routine satan. 693e12c5d1SDavid du Colombier */ 703e12c5d1SDavid du Colombier 713e12c5d1SDavid du Colombier double atan(double arg)723e12c5d1SDavid du Colombieratan(double arg) 733e12c5d1SDavid du Colombier { 743e12c5d1SDavid du Colombier 753e12c5d1SDavid du Colombier if(arg > 0) 763e12c5d1SDavid du Colombier return satan(arg); 773e12c5d1SDavid du Colombier return -satan(-arg); 783e12c5d1SDavid du Colombier } 79