C
C Subroutine to sort an integer array. The subroutine is a modification
C of subroutine "INDEXX" from Numerical Recipes.
C
C Version: 6.40
C Date: 2010, June 16
C
C Coded by Petr Bulant
C http://sw3d.cz/staff/bulant.htm
C
C-----------------------------------------------------------------------
C
C Usage:
C     CALL INDEXI(N,IARR,INDX)
C
C Input:
C     N...    Dimension of arrays IARR and INDX.
C     IARR... Array of integer values to be sorted.
C Output:
C     INDX... Array of the indices of sorted values in array IARR.
C             I.e., the lowest value in array IARR is IARR(INDX(1)),
C             the second lowest value is IARR(INDX(2)), and so on.
C Note that array IARR remains unchanged on output.
C
C-----------------------------------------------------------------------
C
      SUBROUTINE indexi(n,iarr,indx)
      INTEGER n,indx(n),M,NSTACK
      INTEGER iarr(n)
      PARAMETER (M=7,NSTACK=50)
      INTEGER i,indxt,ir,itemp,j,jstack,k,l,istack(NSTACK)
      INTEGER ia
      do 11 j=1,n
        indx(j)=j
11    continue
      jstack=0
      l=1
      ir=n
1     if(ir-l.lt.M)then
        do 13 j=l+1,ir
          indxt=indx(j)
          ia=iarr(indxt)
          do 12 i=j-1,1,-1
            if(iarr(indx(i)).le.ia)goto 2
            indx(i+1)=indx(i)
12        continue
          i=0
2         indx(i+1)=indxt
13      continue
        if(jstack.eq.0)return
        ir=istack(jstack)
        l=istack(jstack-1)
        jstack=jstack-2
      else
        k=(l+ir)/2
        itemp=indx(k)
        indx(k)=indx(l+1)
        indx(l+1)=itemp
        if(iarr(indx(l+1)).gt.iarr(indx(ir)))then
          itemp=indx(l+1)
          indx(l+1)=indx(ir)
          indx(ir)=itemp
        endif
        if(iarr(indx(l)).gt.iarr(indx(ir)))then
          itemp=indx(l)
          indx(l)=indx(ir)
          indx(ir)=itemp
        endif
        if(iarr(indx(l+1)).gt.iarr(indx(l)))then
          itemp=indx(l+1)
          indx(l+1)=indx(l)
          indx(l)=itemp
        endif
        i=l+1
        j=ir
        indxt=indx(l)
        ia=iarr(indxt)
3       continue
          i=i+1
        if(iarr(indx(i)).lt.ia)goto 3
4       continue
          j=j-1
        if(iarr(indx(j)).gt.ia)goto 4
        if(j.lt.i)goto 5
        itemp=indx(i)
        indx(i)=indx(j)
        indx(j)=itemp
        goto 3
5       indx(l)=indx(j)
        indx(j)=indxt
        jstack=jstack+2
        if(jstack.gt.NSTACK)then
C         INDEXI-01
          CALL ERROR('INDEXI-01: NSTACK too small in indexi')
        endif
        if(ir-i+1.ge.j-l)then
          istack(jstack)=ir
          istack(jstack-1)=i
          ir=j-1
        else
          istack(jstack)=j-1
          istack(jstack-1)=l
          l=i
        endif
      endif
      goto 1
      END
C
C=======================================================================
C