1 /* 2 misc.c 3 4 *%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 5 * 6 * Part of: SExtractor 7 * 8 * Author: E.BERTIN (IAP) 9 * 10 * Contents: miscellaneous functions. 11 * 12 * Last modify: 13/12/2002 13 * 14 *%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 15 */ 16 17 #ifdef HAVE_CONFIG_H 18 #include "config.h" 19 #endif 20 21 #include "define.h" 22 #include "globals.h" 23 24 25 /******************************** hmedian ***********************************/ 26 /* 27 Median using Heapsort algorithm (for float arrays) (based on Num.Rec algo.). 28 Warning: changes the order of data! 29 */ hmedian(float * ra,int n)30float hmedian(float *ra, int n) 31 32 { 33 int l, j, ir, i; 34 float rra; 35 36 37 if (n<2) 38 return *ra; 39 ra--; 40 for (l = ((ir=n)>>1)+1;;) 41 { 42 if (l>1) 43 rra = ra[--l]; 44 else 45 { 46 rra = ra[ir]; 47 ra[ir] = ra[1]; 48 if (--ir == 1) 49 { 50 ra[1] = rra; 51 return n&1? ra[n/2+1] : (float)((ra[n/2]+ra[n/2+1])/2.0); 52 } 53 } 54 for (j = (i=l)<<1; j <= ir;) 55 { 56 if (j < ir && ra[j] < ra[j+1]) 57 ++j; 58 if (rra < ra[j]) 59 { 60 ra[i] = ra[j]; 61 j += (i=j); 62 } 63 else 64 j = ir + 1; 65 } 66 ra[i] = rra; 67 } 68 69 /* (the 'return' is inside the loop!!) */ 70 } 71 72 73