      subroutine heapst (n,rb,is,iclust,idim,rrb)                       HEAP001
c ----------------------------------------------------------------------HEAP002
c     This program module is part of the author's X-ray Library         HEAP003
c     and may not be sold or reproduced for commercial purposes.        HEAP004
c     Although it has been tested carefully, there is no warranty       HEAP005
c     regarding accuracy or freedom from side-effects.                  HEAP006
c                   (C) 1986-1994 Bernhard Rupp                         HEAP007
c ----------------------------------------------------------------------HEAP008
c     HEAPST sorts did after distances did(i,4) is denotes the element  HEAP009
c     after which is sorted (is=4=distance)                             HEAP010
c ----------------------------------------------------------------------HEAP011
c      parameter (m=8)                                                  HEAP012
      dimension rb(iclust,idim),rrb(idim)                               HEAP013
                                                                        HEAP014
      l=n/2+1                                                           HEAP015
      ir=n                                                              HEAP016
10    continue                                                          HEAP017
        if(l.gt.1)then                                                  HEAP018
          l=l-1                                                         HEAP019
          rra=rb(l,is)                                                  HEAP020
          do 1001 ia=1,idim                                                HEAP021
 1001        rrb(ia)=rb(l,ia)                                           HEAP022
        else                                                            HEAP023
          rra=rb(ir,is)                                                 HEAP024
          do 1002 ia=1,idim                                                HEAP025
 1002        rrb(ia)=rb(ir,ia)                                          HEAP026
          rb(ir,is)=rb(1,is)                                            HEAP027
          do 1003 ia=1,idim                                                HEAP028
 1003        rb(ir,ia)=rb(1,ia)                                         HEAP029
          ir=ir-1                                                       HEAP030
          if(ir.eq.1)then                                               HEAP031
            rb(1,is)=rra                                                HEAP032
            do 1004 ia=1,idim                                              HEAP033
 1004          rb(1,ia)=rrb(ia)                                         HEAP034
            return                                                      HEAP035
          endif                                                         HEAP036
        endif                                                           HEAP037
        i=l                                                             HEAP038
        j=l+l                                                           HEAP039
20      if(j.le.ir)then                                                 HEAP040
          if(j.lt.ir)then                                               HEAP041
            if(rb(j,is).lt.rb(j+1,is))j=j+1                             HEAP042
          endif                                                         HEAP043
          if(rra.lt.rb(j,is))then                                       HEAP044
            rb(i,is)=rb(j,is)                                           HEAP045
            do 1005 ia=1,idim                                              HEAP046
 1005          rb(i,ia)=rb(j,ia)                                        HEAP047
            i=j                                                         HEAP048
            j=j+j                                                       HEAP049
          else                                                          HEAP050
            j=ir+1                                                      HEAP051
          endif                                                         HEAP052
        go to 20                                                        HEAP053
        endif                                                           HEAP054
        rb(i,is)=rra                                                    HEAP055
        do 1006 ia=1,idim                                                  HEAP056
 1006      rb(i,ia)=rrb(ia)                                             HEAP057
      go to 10                                                          HEAP058
      return                                                            HEAP059
      end                                                               HEAP060
