program testpos
      implicit none
      save

c..this program exercises the main position subroutine

c..declare
      character*60     string
      character*40     aname
      character*12     otype(3)
      character*3      amon(12)

      integer          n,iyear,imonth,iday,icoord
      double precision rhour,glon,glat,height,tjd,
     1                 xequ,yequ,zequ,vxequ,vyequ,vzequ,
     2                 ra,dec,rah,ram,ras,decd,decm,decs,
     3                 xecl,yecl,zecl,vxecl,vyecl,vzecl,
     4                 lam,bet,lamd,lamm,lams,betd,betm,bets,
     5                 xalt,yalt,zalt,vxalt,vyalt,vzalt,
     6                 alt,azi,altd,altm,alts,azid,azim,azis,
     7                 dist,drdt,dsun,dsundt,
     8                 semi_maj,qper,xecen,lan,xinc,aop,xmanom,
     9                 period,tjdp

      integer ii

c..for picking which ephemeris to use
      integer          which_eph
      common /dop/     which_eph


c..data statements
      data       amon /'jan' , 'feb' , 'mar' , 'apr' , 'may' , 'jun' , 
     1                 'jul' , 'aug' , 'sep' , 'oct' , 'nov' , 'dec' /

      data      otype /'geocentric', 'heliocentric','topographic'/  


c..popular format statements
 02   format(1x,a,a,/,
     1       1x,a,i2.2,a,i4,' +',0pf5.2' hrs',t40, a,0pf13.4,/,
     2       1x,a,a,/,
     3       1x,a,1pe11.3,t27,a,1pe11.3,/)

 03   format(1x,a,/,
     1       1x,a,1p3e16.8,/,
     2       1x,a,1p3e16.8,/,
     3       1x,a,1p2e16.8,/,
     4       1x,a,i4.2,':',i2.2,':',0pf7.4,
     5        ' ',i4.2,':',i2.2,':',0pf7.4,/)

 04   format(1x,a,1pe15.8,t35,a,1pe15.8,/,
     1       1x,a,1pe15.8,t35,a,1pe15.8,/)

 05   format(1x,a,/,
     1     1x,a,1pe15.8,/,
     2     1x,a,1pe15.8,/,
     3     1x,a,1pe15.8,/,
     4     1x,a,1pe15.8,/,
     5     1x,a,1pe15.8,/,
     6     1x,a,1pe15.8,/,
     7     1x,a,1pe15.8,/,
     8     1x,a,1pe15.8,/,
     9     1x,a,0p1f14.5,/)

 06   format(1x,a,/,
     1       1x,0pf12.7,/
     2       1x,0pf12.7,/
     3       1x,0pf12.7,/
     4       1x,0pf12.7,/)




c..initialize the name database
      which_eph = 405
      call init_names(1)


c..get the input vector
c..n      = body number
c..iyear  = year
c..imonth = month
c..iday   = day
c..rhour  = ut hours past midnight of iday
c..icoord = orgin type  1=geocentric 2=heliocentric 3=topographic
c..glon   = geographic longitude in decimal degrees (+west, -east of greenwich)
c..glat   = geographic latitude in decimal degrees (+north, -south of equator)
c..height = height above mean sea level in meters

 10   write(6,*) 
     1  'give body iyear imonth iday rhour icoord glon glat height=>'
      read(5,*)  n,iyear,imonth,iday,rhour,icoord,glon,glat,height


c..get the julian day number
      call juldat(iyear,imonth,iday,rhour,tjd)


c      ii = 1
c      do n=1,11


c..get the state vectors
      call position(
     1              n,tjd,icoord,glon,glat,height,
     2              aname,
     3              xequ,yequ,zequ,vxequ,vyequ,vzequ,
     4              ra,dec,rah,ram,ras,decd,decm,decs,
     5              xecl,yecl,zecl,vxecl,vyecl,vzecl,
     6              lam,bet,lamd,lamm,lams,betd,betm,bets,
     7              xalt,yalt,zalt,vxalt,vyalt,vzalt,
     8              alt,azi,altd,altm,alts,azid,azim,azis,
     9              dist,drdt,dsun,dsundt,
     &              semi_maj,qper,xecen,lan,xinc,aop,xmanom,period,tjdp)


c      if (ii .eq. 1) goto 100

c..say what we got
      write(6,*) 
      write(6,02) 
     1   'body      = ',aname,
     2   'date      = ',iday,amon(imonth),iyear,rhour,
     3   'julian day= ',tjd ,
     4   'origin    = ',otype(icoord),
     5   'longitude =',glon,
     6   'latitude  =',glat

      write(6,03)
     1   'equatorial coordinates:',
     2   'x y z    =',xequ,yequ,zequ,
     3   'vx vy vz =',vxequ,vyequ,vzequ,
     4   'ra dec   =',ra,dec,
     5   'ra dec   = ',int(rah),int(abs(ram)),abs(ras),
     6                 int(decd),int(abs(decm)),abs(decs)

      write(6,03)
     1   'ecliptic coordinates:',
     2   'x y z    =',xecl,yecl,zecl,
     3   'vx vy vz =',vxecl,vyecl,vzecl,
     4   'lam bet  =',lam,bet,
     5   'lam bet  = ',int(lamd),int(abs(lamm)),abs(lams),
     6                 int(betd),int(abs(betm)),abs(bets)

      if (icoord .eq. 3) then
       write(6,03)
     1   'altitude-azimuth coordinates:',
     2   'x y z    =',xalt,yalt,zalt,
     3   'vx vy vz =',vxalt,vyalt,vzalt,
     4   'alt azi  =',alt,azi,
     5   'alt azi  = ',int(altd),int(abs(altm)),abs(alts),
     6                 int(azid),int(abs(azim)),abs(azis)
      end if


      write(6,04)
     1   'dist from earth =',dist,' dr/dt    = ',drdt, 
     2   'dist from sun   =',dsun,' drsun/dt = ',dsundt

      string = 'osculating elements for a heliocentric orbit:'
      if (n.eq.2) string = 'osculating elements for a geocentric orbit:'

      write(6,05) string,
     1   'semi_major axis in au         =',semi_maj,
     2   'perihelion distance in au     =',qper,
     3   'eccentricity                  =',xecen,
     4   'inclination in degrees        =',xinc,
     5   'longitude of ascending node   =',lan,
     6   'angle of perhelion in degrees =',aop,
     7   'mean anomoly in degrees       =',xmanom,
     8   'period in days                =',period,
     9   'julian date of perihelion     =',tjdp
 

c 100  continue
c      write(6,06) aname,15.0d0*ra,dec,lam,bet
c
c      enddo


c..and go back for another input vector
      goto 10
      end




      include 'position_routine.f'
      include 'glue.f'