program xspacecraft_pv
!
!     Exercise driver for spacecraft_pv
!
!     B. Knapp, 1999-08-18, 2000-05-01
!
      implicit none
!
      integer status, i, id, pathlen
      character*(1024) path
      real*8 jd, pv(6), dpv(6), f, r, dr, v, dv, re
      parameter (re=6378.135d0)
!
!     Subroutines
      external spacecraft_pv
!
!     Get the spacecraft id & test date
 1    write(*,'(/a$)') ' Spacecraft ID, Julian Date? '
      read(*,*) id, jd
!
!     Call spacecraft_pv
      path = ''
      pathlen = 0
      call spacecraft_pv( id, pathlen, path, jd, pv, dpv, f, status )
      write(*,*) ' Interpolation fraction: ', f
      write(*,*) ' Status: ', status
      write(*,*) ''
!
      do i=1,3
         write(*,'(2f24.16)') pv(i)*re, dpv(i)*re
      enddo
      do i=4,6
         write(*,'(2f24.16)') pv(i)*re/60.d0, dpv(i)*re/60.d0
      enddo
      write(*,*) ''

!     Compute distance from center of earth and orbital speed
      r = 0.d0
      dr = 0.d0
      v = 0.d0
      dv = 0.d0
      do i=1,3
         r = r + pv(i)**2
         dr = dr + (pv(i)*dpv(i))**2
         v = v + pv(i+3)**2
         dv = dv + (pv(i+3)*dpv(i+3))**2
      enddo
      r = sqrt(r)
      dr = sqrt(dr)/r
      v = sqrt(v)
      dv = sqrt(dv)/v

      write(*,40) '  r,  x,  y,  z:',  r, (pv(i),i=1,3)
      write(*,40) ' dr, dx, dy, dz:', dr, (dpv(i),i=1,3)
      write(*,40) '  v, vx, vy, vz:',  v, (pv(i),i=4,6)
      write(*,40) ' dv,dvx,dvy,dvz:', dv, (dpv(i),i=4,6)
 40   format(a,1p,4e16.8)
!
!
!     Now, repeat, but convert to km, km/s
      r = 0.d0
      dr = 0.d0
      v = 0.d0
      dv = 0.d0
      do i=1,3
         pv(i) = pv(i)*re
         dpv(i) = dpv(i)*re
         r = r + pv(i)**2
         dr = dr + (pv(i)*dpv(i))**2
         pv(i+3) = pv(i+3)*re/60.d0
         dpv(i+3) = dpv(i+3)*re/60.d0
         v = v + pv(i+3)**2
         dv = dv + (pv(i+3)*dpv(i+3))**2
      enddo
      r = sqrt(r)
      dr = sqrt(dr)/r
      v = sqrt(v)
      dv = sqrt(dv)/v

      write(*,*) ' '
      write(*,40) '  r,  x,  y,  z:',  r, (pv(i),i=1,3)
      write(*,40) ' dr, dx, dy, dz:', dr, (dpv(i),i=1,3)
      write(*,40) '  v, vx, vy, vz:',  v, (pv(i),i=4,6)
      write(*,40) ' dv,dvx,dvy,dvz:', dv, (dpv(i),i=4,6)
!
      goto 1
!
      end