program make_vsop87
!
!     Read the VSOP87 ASCII data files, and produce separate head and
!     tail files of a Fortran BLOCK DATA subprogram containing the
!     cosine-only terms of the VSOP87D planetary coefficients for all
!     eight planets, for use by hcpvd.f
!
!     B. Knapp, 1998-05-06, 2000-09-01, 2006-08-15
!
C
C     RCS DATA
C     
C     $Header: /mizar/proj/sorce/razor_db/RAZOR_UNIVERSE/DOMAIN_01/sorce/Archive/RZ_VCS/src/ancillary/astronomy/make_vsop87.f,v 1.2 2006-08-15 20:57:54+00 knapp Exp $
C     
C     $Log: make_vsop87.f,v $
C     Revision 1.2  2006-08-15 20:57:54+00  knapp
C     Port to g77 compiler (fix some non-standard Fortran usage)
C
C     Revision 1.1  2001-05-22 21:24:57+00  knapp
C     Initial revision
C
C
      implicit none
!
      character*132 line1, line2
      character*16 infile
      character*8 pname(8)
     &   /' mercury', '   venus', '   earth', '    mars',
     &    ' jupiter', '  saturn', '  uranus', ' neptune' /
      character*18 c1
      character*14 c2
      character*20 c3
      character*1  c4
      integer*4 offset(0:5,3,8),nterms(0:5,3,8)
      integer*4 ios,i,j,k,len,m,mcheck,n,page
      character*1 ff
      ff=char(12)
!
!     Write the file containing the coefficients (tail)
      open(2,file='vsop87.tail',status='unknown',form='formatted')
      n = 0
      page = 0
      do k=1,8

!        Initialize (not all array elements are actually used)
         do j=1,3
            do i=0,5
               offset(i,j,k) = 0
               nterms(i,j,k) = 0
            enddo
         enddo
         m = 0
         mcheck = 0
!
         write(infile,'(a,a)') pname(k),'d.dat'
         call str_trim(3,infile,infile,len)
         write(*,*) infile
         open(1,file=infile,status='old',form='formatted')
         read(1,'(a)',iostat=ios) line2
         do while( ios .eq. 0 )
            line1 = line2
            read(1,'(a)',iostat=ios) line2
            if (index(line1,'VSOP87') .gt. 0) then
               if (mcheck .ne. m) stop 'Missing terms!'
               read(line1,'(41x,i1,17x,i1,2x,i5)') j,i,m
               offset(i,j,k) = n+1
               nterms(i,j,k) = m
               write(*,'(a10,4i6)') pname(k),i,j,offset(i,j,k),m
               mcheck = 0
            else if (m .gt. 0) then
               mcheck = mcheck+1
               read(line1,'(79x,a18,a14,a20)') c1,c2,c3
               n = n+1
               if ( mod(n-1,64) .eq. 0 ) then
                  page = page+1
                  if ( page .lt. 494 ) then
                     write(2,20) ff,n,n+63,page
                  else
                     write(2,20) ff,n,n+24,page
                  endif
 20               format(a/6X,'DATA ((COEFFS(I,J),I=1,3),J=',
     &               i5,',',i5,') /   !Page',i5)
               endif
               if ( mod(n,64) .eq. 0 .or.
     &              (ios .lt. 0 .and. k .eq. 8) ) then
                  c4 = '/'
               else
                  c4 = ','
               endif
               write(2,21) c1,c2,c3,c4
 21            format(5x,'&',a18,',',a14,',',a20,a1)
            endif
         enddo
         close(1)
      enddo
      write(2,'(a)') '      END'
      write(*,*) page, n
      if ( ios .gt. 0 .or. n .ne. 31577 ) then
         close(2,status='delete')
         stop ' Error creating VSOP87.TAIL!'
      else
         close(2)
      endif
!
!     Write the file containing the indices & offsets (head)
      open(3,file='vsop87.head',status='unknown',form='formatted')
      write(3,30) n
 30   format(
     &'      BLOCK DATA VSOP87'/
     &'!'/
     &'!     VSOP87D data (cosine terms), all eight planets'/
     &'!'/
     &'!     This file created by make_vsop87.f'/
     &'!'/
     &'!     B. Knapp, 2000-09-01'/
     &'!'/
     &'      implicit none'/
     &'!'/
     &'      integer*4 NDATA'/
     &'      parameter (NDATA=',i5,')'/
     &'      integer*4 OFFSET(0:5,3,8), NTERMS(0:5,3,8), I, J'/
     &'      real*8 COEFFS(3,NDATA)'/
     &'      common /VSOP87_COEFFS/ OFFSET, NTERMS, COEFFS')
      write(3,31) offset
 31   format(
     &    '!'/
     &    '      DATA OFFSET /',
     &23(/'     &',6(I6,',')),
     &  (/'     &',5(I6,','),I6,'/'))
      write(3,32) nterms
 32   format(
     &    '!'/
     &    '      DATA NTERMS /',
     &23(/'     &',6(I6,',')),
     &  (/'     &',5(I6,','),I6,'/'))
      close(3)
!
      end

      SUBROUTINE str_trim(mode,src,dest,n)
C
C     Trim leading and/or trailing white space from a string.
C
C     B. Knapp, 1988-01-28, 2001-04-24 (i*4)
C
C     Input:
C        mode  -- INTEGER*4 variable or constant.
C                 Set mode = 1, 2, or 3 for leading, trailing, or both.
C        src   -- CHARACTER string variable or constant.
C                 String to be stripped.
C
C     Output:
C        dest  -- CHARACTER string variable.  Note: if src is entered
C                 twice in the parameter list, replacing dest, then
C                 dest overwrites src.
C        n     -- INTEGER*4 variable, returning the length (of the
C                 stripped portion) of dest.
C
      IMPLICIT NONE
              
      INTEGER*4 mode,n
      CHARACTER*(*) src,dest
      INTEGER*4 i,j

      i = 1
      j = LEN(src)

      IF (mode.GE.2) THEN
C        Set j = index of last non-blank character in src
   10    IF (ICHAR(src(j:j)).LE.32 .AND. j.GT.1) THEN
            j = j-1
            GOTO 10
         ENDIF
      ENDIF

      IF (MOD(mode,2) .EQ. 1) THEN
C        Set i = index of first non-blank character in src
   20    IF (ICHAR(src(i:i)).LE.32 .AND. i.LT.j) THEN
            i = i+1
            GOTO 20
         ENDIF
      ENDIF

      dest = src(i:j)
      n = j-i+1
      RETURN
      END