c*******************************************************************************
c$ Component_name:
c	TABLE3
c$ Abstract:
c	Sample program to demonstrate the use of FKEP_SOLVEA() for locating
c	satellite resonances.
c$ Keywords:
c	KEPLER, ORBITAL_MOTION
c	FORTRAN, SAMPLE_PROGRAM
c$ Declarations:
c	none
c$ Inputs:
c	none
c$ Outputs:
c	none
c$ Returns:
c	none
c$ Detailed_description:
c	This simple test program illustrates the use of the KEPLER library
c	routines to solve for the location of resonances.  When run, it should
c	reproduce Table 3 of Lissauer and Cuzzi (1982), Resonances in Saturn's
c	Rings, Astron. J. 87, 1051-1058.  This table lists the Lindblad
c	resonance locations of 1980S28 in units of Saturn radii.
c$ External_references:
c	FKEP_SETPLANET(), FKEP_SETERROR(), FKEP_SOLVEA()
c$ Examples:
c	The output of this program should be as follows:
c	   2/  1     1.45287
c	   3/  2     1.74835
c	   4/  3     1.88816
c	   5/  4     1.96982
c	   6/  5     2.02341
c	   7/  6     2.06128
c	   8/  7     2.08947
c	   9/  8     2.11127
c	  10/  9     2.12864
c	  20/ 19     2.20597
c	  30/ 29     2.23147
c	  40/ 39     2.24416
c	  50/ 49     2.25176
c$ Error_handling:
c	none
c$ Limitations:
c	none
c$ Author_and_institution:
c	Mark R. Showalter
c	NASA/Ames Research Center
c$ Version_and_date:
c	1991 June 18
c$ Change_history:
c	none
c*******************************************************************************

	program TABLE3

	integer		m
	real*8		rate, a
	real*8		fKep_SolveA
	real*8		pi
	data		pi/ 3.141592653589793 /

c Parameters of Saturn's gravitational field
	real*8		Rs, GM, Js(3)
	data		Rs/ 60330. /
	data		GM/ 3.792914d7 /
	data		Js/ 0.0162991, -0.0009167, 0.0000813 /

c Parameters of 1980S28
	real*8		period
	data		period/ 51993.4 /

c***********************************************************************
c The equation to be solved is
c	m * omega(1980S28) = m * omega(a) - kappa(a)
c for specified values of m.
c***********************************************************************

c Define the planetary gravity field
	call FKep_SetPlanet( Rs, GM, Js, 3 )

c Reduce the accuracy sought, since we only need 6 digits
	call FKep_SetError( 1.d-6 )

c Determine mean motion of 1980S28
	rate = 2.*pi / period

c Solve for the first 9 resonances
	do m = 2, 10
		a = FKep_SolveA( m*rate, m, -1, 0 )
		write(*,10)  m, m-1, a/Rs
	end do

c Solve for resonances up to order 50 in steps of 10.
	do m = 20, 50, 10
		a = FKep_SolveA( m*rate, m, -1, 0 )
		write(*,10)  m, m-1, a/Rs
	end do

10	format( 1x, i3, '/', i3, f12.5 )

	stop
	end

c***********************************************************************