xref: /netbsd-src/lib/libm/src/s_ilogb.c (revision d9158b13b5dfe46201430699a3f7a235ecf28df3)
1 /* @(#)s_ilogb.c 5.1 93/09/24 */
2 /*
3  * ====================================================
4  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
5  *
6  * Developed at SunPro, a Sun Microsystems, Inc. business.
7  * Permission to use, copy, modify, and distribute this
8  * software is freely granted, provided that this notice
9  * is preserved.
10  * ====================================================
11  */
12 
13 #ifndef lint
14 static char rcsid[] = "$Id: s_ilogb.c,v 1.4 1994/03/03 17:04:38 jtc Exp $";
15 #endif
16 
17 /* ilogb(double x)
18  * return the binary exponent of non-zero x
19  * ilogb(0) = 0x80000001
20  * ilogb(inf/NaN) = 0x7fffffff (no signal is raised)
21  */
22 
23 #include <math.h>
24 #include <machine/endian.h>
25 
26 #if BYTE_ORDER == LITTLE_ENDIAN
27 #define n0	1
28 #else
29 #define n0	0
30 #endif
31 
32 #ifdef __STDC__
33 	int ilogb(double x)
34 #else
35 	int ilogb(x)
36 	double x;
37 #endif
38 {
39 	int hx,lx,ix;
40 
41 	hx  = (*(n0+(unsigned*)&x))&0x7fffffff;	/* high word of x */
42 	if(hx<0x00100000) {
43 	    lx = *(1-n0+(int*)&x);
44 	    if((hx|lx)==0)
45 		return 0x80000001;	/* ilogb(0) = 0x80000001 */
46 	    else			/* subnormal x */
47 		if(hx==0) {
48 		    for (ix = -1043; lx>0; lx<<=1) ix -=1;
49 		} else {
50 		    for (ix = -1022,hx<<=11; hx>0; hx<<=1) ix -=1;
51 		}
52 	    return ix;
53 	}
54 	else if (hx<0x7ff00000) return (hx>>20)-1023;
55 	else return 0x7fffffff;
56 }
57