      subroutine int_qsort(data,order,length)
C 
C Description:
C
C This routine sorts a list of integer numbers using the quick sort algorithm. 
C It returns a list of integers which are the indices of the sorted numbers 
C in ascending sorted order.
C 
C Name    Type     Structure  In/out  Purpose
C
C data    integer  (length)   in      The list of numbers to be sorted
C order   integer  (length)   out     The list of indices to the sorted 
C                                     numbers.
C length  integer  scalar     in      The length of the arrays 'data' and 
C                                     'order'
C
C Errors:
C
C If the length of the array to be sorted is one or less, the routine
C returns without attempting to sort. The size of the input array cannot
C exceed 2**32.
C
      implicit none
      integer data,key
      integer order,length,i,stack_pointer,stack,first,last,
     $     low,high,key_index
      dimension data(length),order(length),stack(0:31,2)
      if (length .gt. 1) then
         do i=1,length
            order(i)=i
         end do
         stack_pointer=0
         stack(stack_pointer,1)=1
         stack(stack_pointer,2)=length
         stack_pointer=1
         do while (stack_pointer .gt. 0)
            stack_pointer=stack_pointer-1
            first=stack(stack_pointer,1)
            last=stack(stack_pointer,2)
            do while (first. lt. last)
               low=first
               high=last
               key=data(order(first))
               key_index=order(first)
               do while (high .gt. low)
                  do while ((data(order(high)) .gt. key) .and.
     $                 (high .gt. low))
                     high=high-1
                  end do
                  if (high .gt. low) then
                     order(low)=order(high)
                     low=low+1
                  end if
                  do while ((data(order(low)) .lt. key) .and.
     $                 (high .gt. low))
                     low=low+1
                  end do
                  if (high .gt. low) then
                     order(high)=order(low)
                     high=high-1
                  end if
               end do
               order(high)=key_index
               if (high-first .ge. last-high) then
                  stack(stack_pointer,1)=first
                  stack(stack_pointer,2)=high-1
                  stack_pointer=stack_pointer+1
                  first=high+1
               else
                  stack(stack_pointer,1)=high+1
                  stack(stack_pointer,2)=last
                  stack_pointer=stack_pointer+1
                  last=high-1             
               end if
            end do
         end do
      end if
      return
      end
