Logo Search packages:      
Sourcecode: aubio version File versions  Download package

smpl_t vec_median ( fvec_t input  ) 

returns the median of the vector

This Quickselect routine is based on the algorithm described in "Numerical recipes in C", Second Edition, Cambridge University Press, 1992, Section 8.5, ISBN 0-521-43108-5

This code by Nicolas Devillard - 1998. Public domain, available at http://ndevilla.free.fr/median/median/

Definition at line 260 of file mathutils.c.

References _fvec_t::data, and _fvec_t::length.

                                  {
  uint_t n = input->length;
  smpl_t * arr = (smpl_t *) input->data[0];
  uint_t low, high ;
  uint_t median;
  uint_t middle, ll, hh;

  low = 0 ; high = n-1 ; median = (low + high) / 2;
  for (;;) {
    if (high <= low) /* One element only */
      return arr[median] ;

    if (high == low + 1) {  /* Two elements only */
      if (arr[low] > arr[high])
        ELEM_SWAP(arr[low], arr[high]) ;
      return arr[median] ;
    }

    /* Find median of low, middle and high items; swap into position low */
    middle = (low + high) / 2;
    if (arr[middle] > arr[high])    ELEM_SWAP(arr[middle], arr[high]);
    if (arr[low]    > arr[high])    ELEM_SWAP(arr[low],    arr[high]);
    if (arr[middle] > arr[low])     ELEM_SWAP(arr[middle], arr[low]) ;

    /* Swap low item (now in position middle) into position (low+1) */
    ELEM_SWAP(arr[middle], arr[low+1]) ;

    /* Nibble from each end towards middle, swapping items when stuck */
    ll = low + 1;
    hh = high;
    for (;;) {
      do ll++; while (arr[low] > arr[ll]) ;
      do hh--; while (arr[hh]  > arr[low]) ;

      if (hh < ll)
        break;

      ELEM_SWAP(arr[ll], arr[hh]) ;
    }

    /* Swap middle item (in position low) back into correct position */
    ELEM_SWAP(arr[low], arr[hh]) ;

    /* Re-set active partition */
    if (hh <= median)
      low = ll;
    if (hh >= median)
      high = hh - 1;
  }
}


Generated by  Doxygen 1.6.0   Back to index