00001 /*****************************************************************************
00002    Major portions of this software are copyrighted by the Medical College
00003    of Wisconsin, 1994-2000, and are released under the Gnu General Public
00004    License, Version 2.  See the file README.Copyright for details.
00005 ******************************************************************************/
00006 
00007 #include "cs.h"
00008 
00009 /********************************************************************************/
00010 /* insertion_sort : sort an array of float + void*                              */
00011 
00012 static void isort_floatstuff( int n , float * ar , void ** iar )
00013 {
00014    register int  j , p ;  /* array indices */
00015    register float temp ;  /* a[j] holding place */
00016    register void * itemp ;
00017    register float  * a = ar ;
00018    register void ** ia = iar ;
00019 
00020    if( n < 2 ) return ;
00021 
00022    for( j=1 ; j < n ; j++ ){
00023 
00024      if( a[j] < a[j-1] ){   /* out of order */
00025         p    = j ;
00026         temp = a[j] ; itemp = ia[j] ;
00027 
00028        do{
00029            a[p] =  a[p-1] ; /* at this point, a[p-1] > temp, so move it up */
00030           ia[p] = ia[p-1] ;
00031           p-- ;
00032         } while( p > 0 && temp < a[p-1] ) ;
00033 
00034         a[p] = temp ;       /* finally, put temp in its place */
00035        ia[p] = itemp ;
00036      }
00037    }
00038 }
00039 
00040 /********************************************************************************/
00041 /* qsrec : recursive part of quicksort (stack implementation)                   */
00042 
00043 #define QS_STACK  1024  /* stack size */
00044 #define QS_SWAPF(x,y) ( temp=(x),(x)=(y),(y)= temp)
00045 #define QS_SWAPI(i,j) (itemp=(i),(i)=(j),(j)=itemp)
00046 #define QS_SWAPV(i,j) (vtemp=(i),(i)=(j),(j)=vtemp)
00047 
00048 static void qsrec_floatstuff( int n , float * ar , void ** iar , int cutoff )
00049 {
00050    register int i , j ;           /* scanning indices */
00051    register float temp , pivot ;  /* holding places */
00052    register void * ipivot ;
00053    register float * a = ar ;
00054    register void ** ia = iar ;
00055    int itemp ;
00056    void * vtemp ;
00057 
00058    int left , right , mst , stack[QS_STACK] , nnew ;
00059 
00060    /* return if too short (insertion sort will clean up) */
00061 
00062    if( cutoff < 3 ) cutoff = 3 ;
00063    if( n < cutoff ) return ;
00064 
00065    /* initialize stack to start with whole array */
00066 
00067    stack[0] = 0   ;
00068    stack[1] = n-1 ;
00069    mst      = 2   ;
00070 
00071    /* loop while the stack is nonempty */
00072 
00073    while( mst > 0 ){
00074       right = stack[--mst] ;  /* work on subarray from left -> right */
00075       left  = stack[--mst] ;
00076 
00077       i = ( left + right ) / 2 ;           /* middle of subarray */
00078 
00079       /* sort the left, middle, and right a[]'s */
00080 
00081       if( a[left] > a[i]     ){ QS_SWAPF(a[left] ,a[i]    ); QS_SWAPV(ia[left] ,ia[i]    ); }
00082       if( a[left] > a[right] ){ QS_SWAPF(a[left] ,a[right]); QS_SWAPV(ia[left] ,ia[right]); }
00083       if( a[i] > a[right]    ){ QS_SWAPF(a[right],a[i]    ); QS_SWAPV(ia[right],ia[i]    ); }
00084 
00085       pivot  = a[i] ;                      /* a[i] is the median-of-3 pivot! */
00086       a[i]   = a[right] ;
00087       ipivot = ia[i] ;
00088       ia[i]  = ia[right] ;
00089 
00090       i = left ;                           /* initialize scanning */
00091       j = right ;
00092 
00093       /*----- partition:  move elements bigger than pivot up and elements
00094                           smaller than pivot down, scanning in from ends -----*/
00095 
00096       do{
00097         for( ; a[++i] < pivot ; ) ;  /* scan i up,   until a[i] >= pivot */
00098         for( ; a[--j] > pivot ; ) ;  /* scan j down, until a[j] <= pivot */
00099 
00100         if( j <= i ) break ;         /* if j meets i, quit */
00101 
00102         QS_SWAPF( a[i] , a[j] ) ; QS_SWAPV( ia[i] , ia[j] ) ;
00103       } while( 1 ) ;
00104 
00105       /*----- at this point, the array is partitioned -----*/
00106 
00107       a[right]  = a[i] ;           /*restore the pivot*/
00108       a[i]      = pivot ;
00109       ia[right] = ia[i] ;
00110       ia[i]     = ipivot ;
00111 
00112       /*----- push subarrays [left..i-1] and [i+1..right] onto stack, if big -----*/
00113 
00114       nnew = 0 ;
00115       if( (i-left)  > cutoff ){ stack[mst++] = left ; stack[mst++] = i-1   ; nnew++ ; }
00116       if( (right-i) > cutoff ){ stack[mst++] = i+1  ; stack[mst++] = right ; nnew++ ; }
00117 
00118       /* if just added two subarrays to stack, make sure shorter one comes first */
00119 
00120       if( nnew == 2 && stack[mst-3] - stack[mst-4] > stack[mst-1] - stack[mst-2] ){
00121          QS_SWAPI( stack[mst-4] , stack[mst-2] ) ;
00122          QS_SWAPI( stack[mst-3] , stack[mst-1] ) ;
00123       }
00124 
00125    }  /* end of while stack is non-empty */
00126 
00127 }
00128 
00129 /********************************************************************************/
00130 /* quick_sort :  sort an array partially recursively, and partially insertion   */
00131 
00132 #ifndef QS_CUTOFF
00133 #define QS_CUTOFF 10
00134 #endif
00135 
00136 void qsort_floatstuff( int n , float * a , void ** ia )
00137 {
00138    qsrec_floatstuff( n , a , ia , QS_CUTOFF ) ;
00139    isort_floatstuff( n , a , ia ) ;
00140    return ;
00141 }