1*29957Smckusick /* 2*29957Smckusick * @(#)zabs.c 5.1 11/03/86 3*29957Smckusick */ 4*29957Smckusick 5*29957Smckusick /* 6*29957Smckusick * "@(#)zabs.c 1.1" 7*29957Smckusick */ 8*29957Smckusick 9*29957Smckusick /* THIS IS BASED ON TAHOE FP REPRESENTATION */ 10*29957Smckusick #include <tahoemath/FP.h> 11*29957Smckusick 12*29957Smckusick double zabs(real, imag) 13*29957Smckusick double real, imag; 14*29957Smckusick { 15*29957Smckusick double temp, sqrt(); 16*29957Smckusick 17*29957Smckusick if(real < 0) 18*29957Smckusick *(long int *)&real ^= SIGN_BIT; 19*29957Smckusick if(imag < 0) 20*29957Smckusick *(long int *)&imag ^= SIGN_BIT; 21*29957Smckusick if(imag > real){ 22*29957Smckusick temp = real; 23*29957Smckusick real = imag; 24*29957Smckusick imag = temp; 25*29957Smckusick } 26*29957Smckusick if(imag == 0.) /* if((real+imag) == real) */ 27*29957Smckusick return(real); 28*29957Smckusick 29*29957Smckusick temp = imag/real; 30*29957Smckusick temp = real*sqrt(1.0 + temp*temp); /*overflow!!*/ 31*29957Smckusick return(temp); 32*29957Smckusick } 33