PROGRAM TRY
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
      PARAMETER (MAXM=20)
      DIMENSION X(MAXM)
      integer IINDEX(MAXM), IRANK(MAXM)
      
      ISEED = 54356
      DO I = 1, 54356
      A = RAN1(ISEED)
      ENDDO
      
      DO I = 1, MAXM
      X(I) = (ran1(iseed)-0.5)*5
      print *, X(I), round(X(I),1), round(X(I),2), x(i)-nint(X(I))
      ENDDO

      
      STOP
      END

      double precision function round(tt,mm)
      implicit double precision (a-h,o-z)
      t = tt 
      n = t 
      if( t .lt. 0.d0 ) n = n - 1
      round = t - n
      if(mm.eq.2.and. abs(round-1.d0) .lt. 1.d-6) round = 0.d0
      return
      end
      
      
      DOUBLE PRECISION FUNCTION RAN1(IDUM)
C     ---------------------------------------------------------
      INTEGER IDUM,IA,IM,IQ,IR,NTAB,NDIV
      DOUBLE PRECISION AM,EPS,RNMX
      PARAMETER (IA=16807,IM=2147483647,AM=1.D0/IM,IQ=127773,IR=2836,
     A           NTAB=32,NDIV=1+(IM-1)/NTAB,EPS=1.2D-7,RNMX=1.D0-EPS)
      INTEGER J,K,IV(NTAB),IY
      SAVE IV,IY
      DATA IV /NTAB*0/, IY /0/
 54   IF (IDUM.LE.0.OR.IY.EQ.0) THEN
      IDUM=MAX(-IDUM,1)
      DO J=NTAB+8,1,-1
      K=IDUM/IQ
      IDUM=IA*(IDUM-K*IQ)-IR*K
      IF (IDUM.LT.0) IDUM=IDUM+IM
      IF (J.LE.NTAB) IV(J)=IDUM
      ENDDO
      IY=IV(1)
      ENDIF
      K=IDUM/IQ
      IDUM=IA*(IDUM-K*IQ)-IR*K
      IF (IDUM.LT.0) IDUM=IDUM+IM
      J=1+IY/NDIV
      IY=IV(J)
      IV(J)=IDUM
      RAN1=MIN(AM*IY,RNMX)
C     MAKE IT (0,1)
      IF (RAN1.LE.0.D0) GOTO 54
      IF (RAN1.GE.1.D0) GOTO 54      
      RETURN
      END
C     ---------------------------------------------------------

       
      SUBROUTINE SORT (ARRAY, IINDEX, N, IRANK)
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C     **************************************************
C     Sort an array with the HEAPSORT algorithm 
C     (cf. W. H. Press et al., "Numerical Recipes", 
C     Cambridge(1986), p. 229 - 235)
C     **************************************************
      DIMENSION ARRAY(N), IINDEX(N), IRANK(N)
C     Initialize the index array
      DO J = 1,N
      IINDEX(J) = J
      ENDDO
C     Index an array ARRAY of length N, i. e. outputs the 
C     array IINDEX such that ARRAY(IINDEX(J)) is in ascending 
C     order for J=1,2,...,N.  ARRAY itself is not changed.
      IF (N.EQ.1) GOTO 30
      L = N/2+1
      IR = N
 10   IF (L.GT.1) THEN
      L = L-1
      INDEXT = IINDEX(L)
      Q = ARRAY(INDEXT)
      ELSE
      INDEXT = IINDEX(IR)
      Q = ARRAY(INDEXT)
      IINDEX(IR) = IINDEX(1)
      IR = IR-1
      IF (IR.EQ.1) THEN
      IINDEX(1) = INDEXT
      GOTO 30
      ENDIF
      ENDIF
      I = L
      J = L+L
 20   IF (J.LE.IR) THEN
      IF (J.LT.IR) THEN
      IF (ARRAY(IINDEX(J)).LT.ARRAY(IINDEX(J+1))) J=J+1
      ENDIF
      IF (Q.LT.ARRAY(IINDEX(J))) THEN
      IINDEX(I) = IINDEX(J)
      I = J
      J = J+J
      ELSE
      J = IR+1
      ENDIF
      GOTO 20
      ENDIF
      IINDEX(I) = INDEXT
      GOTO 10
 30   DO J=1,N
      IRANK(IINDEX(J)) = J
      ENDDO
      RETURN
      END


C
      subroutine tri(array,iindex,n,irank)
      implicit double precision (a-h,o-z)

c*******************************************************************************
c     Sorting of an array with the HEAPSORT algorithm (cf. W. H. Press et al., 
c     "Numerical Recipes", Cambridge(1986), p. 229 - 235.
c*******************************************************************************

      dimension array(5),iindex(5),irank(5)
      
      irank(1) = 1
      print *, 'irank = ', irank(1)

c *** "indexing"
      call sindex(n,array,iindex)
      
      do i=1,n
      print *, i, iindex(i)
      enddo
      print *, ' '

c *** "ranking"
      call srank(n,iindex,irank)

      return
      end

*###############################################################################

      subroutine sindex(n,array,iindex)

      implicit double precision (a-h,o-z)
c*******************************************************************************
c     indexes an array ARRAY of length N, i. e. outputs the array IINDEX such 
c     that ARRAY(IINDEX(J)) is in ascending order for J=1,2,...,N. The input 
c     quantities N and ARRAY are not changed.
c*******************************************************************************

      dimension array(5),iindex(5)

c *** initialize the index array with consecutive integers
      do 1, j=1,n
        iindex(j) = j
    1 continue
      
ctc
      if(n.eq.1) return

c *** heapsort algorithm
      l = n/2+1
      ir = n
   10 continue
        if (l.gt.1) then
          l = l-1
          indext = iindex(l)
          q = array(indext)
        else
          indext = iindex(ir)
          q = array(indext)
          iindex(ir) = iindex(1)
          ir = ir-1
          if (ir.eq.1) then
            iindex(1) = indext
            return
          end if
        end if
        i = l
        j = l+l
   20   if (j.le.ir) then
            if (j.lt.ir) then
              if (array(iindex(j)).lt.array(iindex(j+1))) j = j+1
            end if
            if (q.lt.array(iindex(j))) then
            iindex(i) = iindex(j)
            i = j
            j = j+j
          else
            j = ir+1
          end if
          go to 20
        end if
        iindex(i) = indext
      go to 10

      end

*###############################################################################

      subroutine srank (n,iindex,irank)
      implicit double precision (a-h,o-z)

c*******************************************************************************
c     Given IINDEX of length N as output from the routine INDEX, this routine
c     returns an array IRANK, the corresponding table of ranks
c*******************************************************************************

      dimension iindex(5)
      dimension irank(5)

      print *, 'coming in'
      
      do i=1,n
      print *, i, iindex(i)
      irank(i) = i
      enddo
      print *, ' '
      
      do i=1,n
      print *, i, irank(i)
      enddo

      print *, 'tulu '
            
      do j=1,n
      print *, 'j = ', j, iindex(j), irank(iindex(j)), iindex(3)
      irank(iindex(j)) = j
      print *, j, iindex(j), irank(iindex(j)) , iindex(3)
      enddo

      do i=1,n
      print *, i, irank(i)
      enddo
      print *, ' '

      end