1 /*****************************************************************************
2 Major portions of this software are copyrighted by the Medical College
3 of Wisconsin, 1994-2000, and are released under the Gnu General Public
4 License, Version 2. See the file README.Copyright for details.
5 ******************************************************************************/
6
7 #include "cs.h"
8
9 /********************************************************************************/
10 /* insertion_sort : sort an array of float + int64_t */
11
isort_floatint64_t(int64_t n,float * ar,int64_t * iar)12 static void isort_floatint64_t( int64_t n , float * ar , int64_t * iar )
13 {
14 register int64_t j , p ; /* array indices */
15 register float temp ; /* a[j] holding place */
16 register int64_t itemp ;
17 register float * a = ar ;
18 register int64_t * ia = iar ;
19
20 if( n < 2 ) return ;
21
22 for( j=1 ; j < n ; j++ ){
23
24 if( a[j] < a[j-1] ){ /* out of order */
25 p = j ;
26 temp = a[j] ; itemp = ia[j] ;
27
28 do{
29 a[p] = a[p-1] ; /* at this point, a[p-1] > temp, so move it up */
30 ia[p] = ia[p-1] ;
31 p-- ;
32 } while( p > 0 && temp < a[p-1] ) ;
33
34 a[p] = temp ; /* finally, put temp in its place */
35 ia[p] = itemp ;
36 }
37 }
38 }
39
40 /********************************************************************************/
41 /* qsrec : recursive part of quicksort (stack implementation) */
42
43 #undef QS_SWAPF
44 #undef QS_SWAPI
45 #define QS_SWAPF(x,y) ( temp=(x),(x)=(y),(y)= temp)
46 #define QS_SWAPI(i,j) (itemp=(i),(i)=(j),(j)=itemp)
47 #ifndef QS_STACK
48 # define QS_STACK 9999
49 #endif
50
qsrec_floatint64_t(int64_t n,float * ar,int64_t * iar,int64_t cutoff)51 static void qsrec_floatint64_t( int64_t n , float * ar , int64_t * iar , int64_t cutoff )
52 {
53 register int64_t i , j ; /* scanning indices */
54 register float temp , pivot ; /* holding places */
55 register int64_t itemp , ipivot ;
56 register float * a = ar ;
57 register int64_t * ia = iar ;
58
59 int64_t left , right , mst , stack[QS_STACK] , nnew ;
60
61 /* return if too short (insertion sort will clean up) */
62
63 if( cutoff < 3 ) cutoff = 3 ;
64 if( n < cutoff ) return ;
65
66 /* initialize stack to start with whole array */
67
68 stack[0] = 0 ;
69 stack[1] = n-1 ;
70 mst = 2 ;
71
72 /* loop while the stack is nonempty */
73
74 while( mst > 0 ){
75 right = stack[--mst] ; /* work on subarray from left -> right */
76 left = stack[--mst] ;
77
78 i = ( left + right ) / 2 ; /* middle of subarray */
79
80 /* sort the left, middle, and right a[]'s */
81
82 if( a[left] > a[i] ){ QS_SWAPF(a[left] ,a[i] ); QS_SWAPI(ia[left] ,ia[i] ); }
83 if( a[left] > a[right] ){ QS_SWAPF(a[left] ,a[right]); QS_SWAPI(ia[left] ,ia[right]); }
84 if( a[i] > a[right] ){ QS_SWAPF(a[right],a[i] ); QS_SWAPI(ia[right],ia[i] ); }
85
86 pivot = a[i] ; /* a[i] is the median-of-3 pivot! */
87 a[i] = a[right] ;
88 ipivot = ia[i] ;
89 ia[i] = ia[right] ;
90
91 i = left ; /* initialize scanning */
92 j = right ;
93
94 /*----- partition: move elements bigger than pivot up and elements
95 smaller than pivot down, scanning in from ends -----*/
96
97 do{
98 for( ; a[++i] < pivot ; ) ; /* scan i up, until a[i] >= pivot */
99 for( ; a[--j] > pivot ; ) ; /* scan j down, until a[j] <= pivot */
100
101 if( j <= i ) break ; /* if j meets i, quit */
102
103 QS_SWAPF( a[i] , a[j] ) ; QS_SWAPI( ia[i] , ia[j] ) ;
104 } while( 1 ) ;
105
106 /*----- at this point, the array is partitioned -----*/
107
108 a[right] = a[i] ; /*restore the pivot*/
109 a[i] = pivot ;
110 ia[right] = ia[i] ;
111 ia[i] = ipivot ;
112
113 /*----- push subarrays [left..i-1] and [i+1..right] onto stack, if big -----*/
114
115 nnew = 0 ;
116 if( (i-left) > cutoff ){ stack[mst++] = left ; stack[mst++] = i-1 ; nnew++ ; }
117 if( (right-i) > cutoff ){ stack[mst++] = i+1 ; stack[mst++] = right ; nnew++ ; }
118
119 /* if just added two subarrays to stack, make sure shorter one comes first */
120
121 if( nnew == 2 && stack[mst-3] - stack[mst-4] > stack[mst-1] - stack[mst-2] ){
122 QS_SWAPI( stack[mst-4] , stack[mst-2] ) ;
123 QS_SWAPI( stack[mst-3] , stack[mst-1] ) ;
124 }
125
126 } /* end of while stack is non-empty */
127
128 }
129
130 /********************************************************************************/
131 /* quick_sort : sort an array partially recursively, and partially insertion */
132
133 #ifndef QS_CUTOFF
134 #define QS_CUTOFF 10
135 #endif
136
qsort_floatint64_t(int64_t n,float * a,int64_t * ia)137 void qsort_floatint64_t( int64_t n , float * a , int64_t * ia )
138 {
139 qsrec_floatint64_t( n , a , ia , QS_CUTOFF ) ;
140 isort_floatint64_t( n , a , ia ) ;
141 return ;
142 }
143