      subroutine heap_iper(n,ix,iper) ! bernd berg, feb 12 2003.  
c version of heap_per, which sorts an integer array ix.
      include '../../libs/fortran/implicit.sta'
      dimension ix(n),iper(n)
c
      do i=1,n
        iper(i)=i ! initialize iper.
      end do
c
c building the heap:
c ==================
      nhalf=n/2
      do m=nhalf,1,-1
        itemp=ix(m) 
        itmp=iper(m)
        i1=m
        i2=m+m
1       continue
c         index the larger of ix(i2) and ix(i2+1) as i2:
          if(i2.lt.n .and. ix(i2+1).gt.ix(i2)) i2=i2+1
          if(ix(i2).gt.itemp) then 
            ix(i1)=ix(i2)
            iper(i1)=iper(i2)
            i1=i2
            i2=i1+i1
          else  
            i2=n+1
          endif
        if(i2.le.n) go to 1
        ix(i1)=itemp ! final position for itemp=ix(m).
        iper(i1)=itmp
      end do
c
c sorting the heap:
c =================
      do i=n,3,-1 ! i is the number of still competing elements.
        itemp=ix(i)
        itmp=iper(i)
        ix(i)=ix(1) ! store top of the heap.
        iper(i)=iper(1)
        i1=1
        i2=2
2       continue
c         index the larger of ix(i2) and ix(i2+1) as i2:
          if((i2+1).lt.i .and. ix(i2+1).gt.ix(i2)) i2=i2+1
          if(ix(i2).gt.itemp) then  
            ix(i1)=ix(i2)
            iper(i1)=iper(i2)
            i1=i2
            i2=i1+i1
          else 
            i2=i
          endif 
        if(i2.lt.i) go to 2
        ix(i1)=itemp ! final position for itemp=ix(i).
        iper(i1)=itmp
      end do 
      itemp=ix(2)
      itmp=iper(2)
      ix(2)=ix(1)
      iper(2)=iper(1)
      ix(1)=itemp
      iper(1)=itmp
      return
      end
