cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine getfile(unit,question,inout,form) 
c      Opens files...filenames not recognized by VAX will generate 
c       an error message and request a second try. 
c      Sample call: 
c           call getfile(21,' Give input file','in','formatted') 
 
      integer unit 
      character filename*50, status*8 
      character*(*) question, inout, form 
      character bel*1 
      parameter (bel=char(7)) 
 
      common/savfilename/ filename 
 
      print *, ' ' 
c 10   write(*,101) question(1:lentrue(question))//': ' 
 10   write(*,101) question(1:lentrue(question)),': ' 
c 101  format('$',a,' ') 
 101  format(a,' ',$) 
      read (*,'(a12)') filename 
 
c     For INPUT, open the file with appropriate attributes... 
      if(inout.eq.'in') then 
        status='old'
        open(unit,file=filename,form=form,status=status,err=900,
     &  mode='read',share='denywr')
c     &     READONLY,SHARED) 
c         (Note that READONLY and SHARED are VMS/VAX specific.)
 
c     For OUTPUT, open the file with appropriate attributes... 
      else if(inout.eq.'out') then 
        status='unknown' 
        if(form.eq.'formatted') then 
          open(unit,file=filename,form=form,status=status, 
     &    err=900) 
c     &    carriagecontrol='list',err=900) 
c        else if(form.eq.'fortran') then 
c          open(unit,file=filename,form='formatted',status=status, 
c     &    carriagecontrol='fortran',err=900) 
        else if(form.eq.'unformatted') then 
          open(unit,file=filename,form=form,status=status,err=900) 
        endif

      else 
        print *, ' ERROR...inout variable in GETFILE not recognized.' 
        print *, '  inout =', inout 
        stop 
 
      endif
      return 
 
c.....Error message for bad open... 
 900  print *, bel 
      print *,' ** Error in opening  ',filename(1:lentrue(filename)), 
     &         ' .... Try again.' 
      close(unit) 
      print *, ' ' 
      go to 10 
 
      end 
 
 
 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine getfile2 (unit,question,inout,form,suffix) 
c       Opens files...filenames not recognized by VAX will generate 
c     an error message and request a second try. 
c       Giving a  blank filename or CR will cause the default filename 
c     appearing in square brackets to be used, unless the default is 
c     itself blank. 
c       The default filename is constructed from suffix in the 
c     following way.  If suffix starts with a '.' a filename root 
c     is searched for in a file called 'savename.tmp' and the new suffix 
c     is appended.  If suffix does not start with a '.', it is assumed 
c     to be the complete filename. 
 
c      Variables: 
c        unit     = io unit 
c        question = query requesting file name 
c        inout    = 'in' or 'out' 
c        form     = 'unformatted' (binary) or 
c                   'formatted' (ascii with no carriage control) or 
c                   'fortran' (ascii with fortran carriage control). 
c        suffix   = '.dat' for example or a complete file name 
 
c      Sample call: 
c           call getfile(21,' Give input file','in','formatted','.dat') 
 
      integer unit 
      character*(*) question, inout, form, suffix 
      character filename*50, defaultfile*50, status*8, suff*80 
      character bel*1 
      parameter (bel=char(7)) 
      common/savfilename/ filename 
 
      suff = suffix 
      call ljust(suff) 
      call lc2uc(suff)             
      lensuff = lentrue (suff)       

c     Suffix does not begin with '.' ... must be a complete filename... 
      if(lensuff.ne.0 .and. suff(1:1).ne.'.') then 
        defaultfile=suff 
 
c     Find filename and add suffix if not blank... 
      else 
        call makenewfn(defaultfile,suff) 
      endif 
 
 10   call lc2uc(defaultfile) 
      leng=lentrue(defaultfile) 
 
c     Ask for file, giving default filename in [] if there is one. 
      if(leng.ge.1) then 
        write(*,101) question(1:lentrue(question)), 
     &  ':  [', defaultfile(1:leng), ']  ' 
 101    format(/,4a,$) 
 
      else 
        write(*,101) question(1:lentrue(question)), ':  ' 
      endif 
 
      read '(a50)', filename 
      call ljust(filename)
      call lc2uc(filename) 
      lenf = lentrue(filename) 
 
      if(filename.eq.' ') then 
        filename = defaultfile 
 
      else if(filename.ne.' ') then 
        idot = index(filename,'.') 
        idotd= index(defaultfile,'.') 
        if (idotd.eq.0) then 
           lendef = lentrue (defaultfile) 
        else 
           lendef = idotd - 1 
        endif 
        if(idot.eq.0) then 
c         Only a new root is given.  Fix default and query again... 
          if(idotd.ne.0) then 
             defaultfile=filename(1:lenf)//defaultfile(idotd:leng) 
          else 
             defaultfile=filename(1:lenf) 
          endif 
          filename = ' ' 
        else if(idot.eq.1) then 
c         Only a new suffix is given.  Fix default and query again... 
          defaultfile=defaultfile(1:lendef)//filename(idot:lenf) 
          filename = ' ' 
        else 
c         Decent new filename is probably given so continue ... 
        endif 
      endif 
 
c     Repeat query if filename is still blank... 
      if(filename.eq.' ') goto 10 
 
c     Save the new filename unless a full filename was given in suffix. 
      if(lensuff.ne.0 .and. suff(1:1).ne.'.') then 
          continue 
      else 
          call putfn(filename) 
      endif 
 
c     Filename is in hand, so open INPUT file... 
      if(inout.eq.'in') then 
        status='old' 
        open(unit,file=filename,form=form,status=status,err=900, 
     &     mode='read',share='denywr')
c     &     READONLY,SHARED) 
c         (Note that READONLY and SHARED are VMS/VAX specific.) 
 
c     Filename is in hand, so open OUTPUT file... 
      else if(inout.eq.'out') then 
        status='unknown' 
        if(form.eq.'formatted') then 
          open(unit,file=filename,form=form,status=status, 
     &    err=900) 
c     &    carriagecontrol='list',err=900) 
c        else if(form.eq.'fortran') then 
c          open(unit,file=filename,form='formatted',status=status, 
c     &    carriagecontrol='fortran',err=900) 
        else if(form.eq.'unformatted') then 
          open(unit,file=filename,form=form,status=status,err=900) 
        endif 
 
      else 
        print *, ' ERROR...inout variable in GETFILE not recognized.' 
        print *, '  inout =', inout 
        stop 
 
      endif 
 
      return 
 
c.....Error message for bad open... 
 900  print *, bel 
      print *,' ** Error in opening  ',filename(1:lentrue(filename)), 
     &         ' .... Try again.' 
      close(unit) 
      go to 10 
 
      end 
 
 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine putfn(filename) 
c     Stores a filename in a file name SAVENAME.TMP for future reference. 
 
      character filename*(*), savename*80 
      parameter (savename='savename.tmp') 
 
      open(91,file=savename,form='formatted',status='unknown',err=90) 
      write(91,'(a)') filename (1:lentrue(filename)) 
      close(91) 
 
 90   return 
      end 
 
 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine makenewfn(filename,suffix) 
c       Finds the root part of a VAX filename in a save file and appends 
c     the given suffix... 
c       Given suffix should have the form '.dat' for example... 
c     If suffix is blank, and filename has a dot in it, then whole 
c     filename is returned as read. 
 
      character filename*(*), tempname*80, savename*80, suffix*(*) 
      parameter (savename='savename.tmp') 
 
      open(91,file=savename,form='formatted',status='old',
     1share='denywr',err=90)
      read(91,'(a)',err=90) filename 
      close(91) 
 
      leng = lentrue(filename) 
      lens = lentrue(suffix) 
 
      if(leng.ne.0 .and. lens.eq.0) then 
c       Filename but no suffix... 
        return 
 
      else if (leng.ne.0 .and. lens.ne.0) then 
c       Filename and new suffix 
        tempname = filename 
        j=index(tempname,']') 
        if(j.ne.0) tempname(1:)=tempname(j+1:) 
        i=index(tempname,'.') 
        i=i-1 
        if(i.le.0) i=lentrue (tempname) 
        filename = tempname(1:i)//suffix(1:lens) 
 
      else if(leng.eq.0) then 
        goto 90 
      endif 
 
      return 
 
c.....Error in filename obtained, so return query answer. 
 90   filename = '?'//suffix(1:lentrue(suffix)) 
      return 
 
      end 
 
 
C***** lentrue.f ***** 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      integer function lentrue (string) 
c       Gives position of last non-blank, non-tab, non-null 
c     character in a string. 
c       Returns 0 if no such beast exists in the string. 
c     Programmed by R.Simpson - U.S.G.S. 
 
      character*(*) string 
      character*1 blank, tab, null 
      parameter (blank=' ', tab=char(9), null=char(0)) 
 
      lentrue=0 
 
      do 100 i=len(string),1,-1 
      if(string(i:i).ne.blank  .and. 
     &   string(i:i).ne.tab    .and. 
     &   string(i:i).ne.null         ) then 
            lentrue=i 
            return 
      endif 
 100  continue 
 
      return 
      end 
 
 
C***** lc2uc.f ***** 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine lc2uc(string) 
c       Converts lower case letters to upper case... 
c     Assumes that the codes are in sequence for a-z and A-Z 
c     Programmed by R.Simpson - U.S.G.S. 
 
      character string*(*), s*1 
 
      lstring=len(string) 
 
      do 20 i=1,lstring 
      s=string(i:i) 
      if(s.ge.'a' .and. s.le.'z') then 
c       string(i:i)=char(ichar(s)-'a'+'A') 
       string(i:i)=char(ichar(s)-ichar('a')+ichar('A')) 
      endif 
 20   continue 
 
      return 
      end 
 
 
C***** uc2lc.f ***** 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine uc2lc(string) 
c       Converts upper case letters to lower case... 
c     Assumes that the codes are in sequence for a-z and A-Z 
c     Programmed by R.Simpson - U.S.G.S. 
 
      character string*(*), s*1 
 
      lstring=len(string) 
 
      do 20 i=1,lstring 
      s=string(i:i) 
      if(s.ge.'A'. and .s.le.'Z') then 
c       string(i:i)=char(ichar(s)-'A'+'a') 
       string(i:i)=char(ichar(s)-ichar('A')+ichar('a')) 
      endif 
 20   continue 
 
      return 
      end 
 
 
C***** ljust.f ***** 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine ljust(string) 
c       Left justifies a string by eliminating blanks, tabs, nulls 
c     at left end of string. 
c     Programmed by R.Simpson - U.S.G.S. 
 
      character*(*) string 
 
      ifirst=leftend(string) 
      if(ifirst.le.1) then 
        return 
      else 
        l2=len(string)-(ifirst-1) 
        string(1:l2)=string(ifirst:ifirst+l2-1) 
        string(l2+1:len(string))=' ' 
      endif 
 
      return 
      end 
 
 
C***** rjust.f ***** 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine rjust(string) 
c     Right justifies a string in its alloted length by eliminating 
c      blanks, tabs, and nulls on the right end. 
c     Programmed by R.Simpson - U.S.G.S. 
 
      character*(*) string 
 
      leng=len(string) 
      ilast=lentrue(string) 
      if(ilast.eq.0.or.ilast.eq.leng) return 
 
      do 100 i=1,ilast 
 100  string(leng+1-i:leng+1-i)=string(ilast+1-i:ilast+1-i) 
 
      do 200 i=ilast+1,leng 
 200  string(leng+1-i:leng+1-i)=' ' 
 
      return 
      end 
 
 
C***** ldeblank.f ***** 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine ldeblank(string1,string2,len2) 
c     Eliminates blanks at left end of string1, 
c      returns as string2 with new true length len2. 
c     Programmed by R.Simpson - U.S.G.S. 
 
      character*(*) string1, string2 
 
      ifirst = leftend(string1) 
      if (ifirst.eq.0) then 
        string2 = ' ' 
        len2 = 0 
      else 
        string2 = string1(ifirst:) 
        len2 = lentrue(string2) 
      endif 
 
      return 
      end 
 
 
C***** rdeblank.f ***** 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine rdeblank(string1,string2,len2) 
c     "Eliminates" blanks at right end of string1, 
c       returns as string2 with true length len2. 
c     Programmed by R.Simpson - U.S.G.S. 
 
      character*(*) string1, string2 
 
      string2 = string1 
      len2 = lentrue(string2) 
 
      return 
      end 
 
 
C***** ideblank.f ***** 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
C+ 
C    IDEBLANK- removes all blanks and returns length 
C 
C    ireturn = IDEBLANK (string) 
C 
C       ireturn = length of the string without blanks 
C        string = character string to be de-blanked 
C 
C- 
c     Programmed by Rick Saltus - U.S.G.S. 
C  ********************************************* 
      Integer Function ideblank(string) 
      Character*256 string*(*),temp 
      j=1 
      temp=' ' 
      itemp=len(string) 
      Do 10 i=1,itemp 
      If (string(i:i).EQ.' ')Go To 10 
      temp(j:j)=string(i:i) 
      j=j+1 
   10 Continue 
      string=temp 
      ideblank=j-1 
      Return 
      End 
 
 
C***** leftend.f ***** 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      integer function leftend(string) 
c     Gives position of first non-blank, non-tab, non-null 
c      character in string. 
c      Returns 0 if no such beast in string. 
c     Programmed by R.Simpson - U.S.G.S. 
 
      character*(*) string 
      character*1 blank, tab, null 
      parameter (blank=' ', tab=char(9), null=char(0)) 
 
      leftend=0 
 
      do 100 i=1,len(string) 
      if(string(i:i).ne.blank.and. 
     &   string(i:i).ne.tab.and. 
     &   string(i:i).ne.null) then 
            leftend=i 
            return 
      endif 
 100  continue 
 
      return 
      end 
 
 
C***** remblanks.f ***** 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine remblanks(string1,string2,len2) 
c     Programmed by R.Simpson - U.S.G.S. 
 
c     Removes all blanks from string1, returns result as string2. 
      character*(*) string1,string2 
 
      len1 = lentrue(string1) 
      len2 = len(string2) 
 
      j=0 
      do 50 i=1,len1 
      if(string1(i:i).ne.' ') then 
       j=j+1 
       string2(j:j)=string1(i:i) 
      endif 
 50   continue 
 
      do 60 i = j+1, len2 
 60   string2(i:i) = ' ' 
 
      len2=j 
 
      return 
      end 
 




 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine putpost (iout, id, post) 
 
c       Writes a single post record. 
 
      character id*8 
      real post(8) 
 
      write (iout) id, post 
 
      return 
      end 
 
 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine closeblp (iout) 
 
      close (iout) 
 
      return 
      end 
 
 
 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine getpost (inunit, id, post, *) 
 
c       Gets a single post record. 
 
      character id*8 
      real post(8) 
 
      read (inunit, end=90) id, post 
      return 
 
 90   return1 
      end 
 
 

cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine getint(ques,intdefault,int) 
c     Converts a string answer to an integer provided the string 
c      is not blank... 
c         ques = question to be asked 
c         intdefault = default answer 
c         int = integer answer 
c     Programmed by R.Simpson - U.S.G.S. 
 
      integer int, intdefault 
      character*(*) ques 
      character*10 default, ans, fmt 
      character bel*1
      parameter (bel=char(7)) 
 
      write(default,'(i10)') intdefault 
      call ljust(default) 
      lengdef=lentrue(default) 
 
      lengq=lentrue(ques) 
c 10   write(*,101) ques(1:lengq)//':  ['//default(1:lengdef)//']  ' 
 10   write(*,101) ques(1:lengq),':  [',default(1:lengdef),']  ' 
 101    format(/4a,$) 
      read (*,'(a)') ans 
      leng=lentrue(ans) 
      if(leng.eq.0) then 
        int=intdefault 
      else 
        write(fmt,'(a,i3,a)') '(i', leng, ')' 
        read(ans(1:leng),fmt,err=90) int 
      endif 
 
      return 
 
c.....Error message 
 90   print '(/,1x,a)', 
     &  '*** ERROR, expecting an integer answer...try again...' 
      print *, bel 
      goto 10 
 
      end 
 


 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine getopt(ques,opts,default,answer) 
c     Queries user for option desired--returns lower case answer. 
c       ques = question. 
c       opts = character string with options separated by /. 
c       default = default option to be used if a CR is the answer. 
c       answer = option selected (or default). 
c 
c     Sample call: 
c 
c       call getopt('Want tic marks?','y/n','n',answer) 
c 
c     Programmed by R.Simpson - U.S.G.S. 
 
      character*(*) ques, opts, default, answer 
      character*132 options, oo 
      character ans*10 
      character bel*1 
      parameter (bel=char(7))
    
 
      lengq=lentrue(ques) 
      lengo=lentrue(opts) 
      lengd=lentrue(default) 
 
c 10   write(*,101) ques(1:lengq)//' ('//opts(1:lengo)// 
c     &                     '):   '//'['//default(1:lengd)//']  ' 
 10   write(*,101) ques(1:lengq),' (',opts(1:lengo), 
     &                     '):   ','[',default(1:lengd),']  ' 
 101    format(/1x,7a$) 
      read (*,'(a10)') ans 
      call ljust(ans) 
      lenga=lentrue(ans)
 
      call uc2lc(opts) 
      call uc2lc(default) 
      call uc2lc(ans)           
 
c     Add slashes to front and back of option string to aid searches. 
c      options = '/'//opts(1:lengo)//'/' 
       options(1:1) = '/'
       options (2:lengo+1) = opts(1:lengo)
       options (lengo+2:lengo+2)= '/'
 
      if(lenga.eq.0) then 
        answer = default 
 
      else
         oo='/'//ans(1:lenga)//'/'
         lengoo = lentrue(oo)
 
        if (index(options, oo(1:lenga)) .eq.0) then
        print '(/,1x,a,a)', ' ****** ERROR...try again....',bel 
        goto 10 
        endif
      
        answer = ans 
 
      endif 
 
      return 
      end 
 
 
 




      SUBROUTINE ROSORT(CHAR8,IPOINT,N)
c      implicit  integer *2 (i-n)
c      CHARACTER*12 CHAR8 
        real char8
c       changed char8 to real to sort array of real #'s MRM
C***  THIS IS A FAST ROUTINE FOR SORTING A SET OF POINTERS
C***  CONTAINED IN ARRAY IPOINT BASED UPON INCREASING VALUES
C***  IN CHARACTER ARRAY CHAR8.  IT IS SIMILAR TO FSORT.FOR.
C***  ROUTINE ORIGINALLY WRITTEN BY CALVIN SMITH, U. S. ATOMIC
C***  ENERGY COMMISSION, GRAND JUNCTION OFFICE, FOR A CDC 3100
C***  COMPUTER.  THIS SORT IS NOT STABLE.
C***  RECODED FOR IBM PERSONAL COMPUTER BY W. D. GRUNDY APRIL 1984
C***
C***
C***  DIMENSION REQUIREMENTS: 
C***  IF CHAR8 IS DIMENSIONED AS CHAR8(M) IN THE MAIN PROGRAM, SET
C***  DIMENSIONS IU(N) AND IL(N) SUCH THAT 2**(N+1)-1 IS 
C***  GREATER THAN OR EQUAL TO M.
C***
      DIMENSION CHAR8(1),IPOINT(1),IU(17),IL(17)
      I=1
      M=1
      J=N
   10 IF(I.GE.J) GO TO 110
   20 K=I
      L=J
      IJ=(I+J)/2
      IF(CHAR8(IPOINT(I)).LE.CHAR8(IPOINT(IJ))) GO TO 30
      ITEMP=IPOINT(I)
      IPOINT(I)=IPOINT(IJ)
      IPOINT(IJ)=ITEMP
   30 IF(CHAR8(IPOINT(J)).GE.CHAR8(IPOINT(IJ))) GO TO 60
      ITEMP=IPOINT(J)
      IPOINT(J)=IPOINT(IJ)
      IPOINT(IJ)=ITEMP
      IF(CHAR8(IPOINT(IJ)).GE.CHAR8(IPOINT(I))) GO TO 60
      ITEMP=IPOINT(I)
      IPOINT(I)=IPOINT(IJ)
      IPOINT(IJ)=ITEMP
      GO TO 60
   40 ITEMP=IPOINT(L)
      IPOINT(L)=IPOINT(K)
      IPOINT(K)=ITEMP
      IF(IJ.NE.L) GO TO 50
      IJ=K
      GO TO 60
   50 IF(IJ.NE.K) GO TO 60
      IJ=L
   60 L=L-1
      IF(CHAR8(IPOINT(IJ)).LT.CHAR8(IPOINT(L))) GO TO 60
   70 K=K+1
      IF(CHAR8(IPOINT(K)).LT.CHAR8(IPOINT(IJ))) GO TO 70
      IF(K-L) 40,80,90
   80 K=K+1
      L=L-1
   90 IF((L-I).LE.(J-K)) GO TO 100
      IL(M)=I
      IU(M)=L
      I=K
      M=M+1
      GO TO 120
  100 IL(M)=K
      IU(M)=J
      J=L
      M=M+1
      GO TO 120
  110 M=M-1
      IF(M.EQ.0) RETURN
      I=IL(M)
      J=IU(M)
  120 IF((J-I).GE.11) GO TO 20
      IF(I.EQ.1) GO TO 10
      GO TO 140
  130 I=I+1
  140 IF(I.EQ.J) GO TO 110
      IF(CHAR8(IPOINT(I)).LE.CHAR8(IPOINT(I+1))) GO TO 130
      K=I
  150 ITEMP=IPOINT(K)
      IPOINT(K)=IPOINT(K+1)
      IPOINT(K+1)=ITEMP
      K=K-1
      IF(CHAR8(IPOINT(K+1)).LT.CHAR8(IPOINT(K))) GO TO 150
      GO TO 130
      END






 
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 
      subroutine getrow(iounit,ncol,row) 
      real row(ncol) 
      read(iounit) dummy,row 
      return 
      end 
 
 
 
