C+
C NAME:
C	smei_orbit_info
C PURPOSE:
C	Defines origin for counting Coriolis orbits
C CALLING SEQUENCE:
c	subroutine smei_orbit_info(n0, t0, p0, dp0)
C OUTPUTS:
C	n0		integer
C	t0(2)		integer
C	p0		double precision
C	dp0		double precision
C CALLS:
C	Time2YDoy
C PROCEDURE:
C	Orbit 1.0 corresponds to 2002/12/31 23:56:41.900 UT
C MODIFICATION HISTORY:
C	JUN-2004, Paul Hick (UCSD/CASS; pphick@ucsd.edu)
C-
c	    integer		n0
c	    integer		t0(2)
c	    double precision	p0
c	    double precision	dp0

c	n0 = 1

c	t0(1) = 0
c	t0(2) = 86201900		! 86400000-198100
c	call Time2YDoy(1,t0,2003,t0)	! Orbit 1.0 is 2002/12/31 23:56:41.900

c	p0  = 0.0705612268519d0		! Nominal period in days
c	dp0 = 0.084d0			! Correction per orbit in msec

c	return
c	end
C+
C NAME:
C	smei_orbit_time
C PURPOSE:
C	Returns time at which Coriolis was in specified orbit
C CALLING SEQUENCE:
c	subroutine smei_orbit_time(orbit, tt)
C INPUTS:
C	orbit		double precision
C					Coriolis orbit number
C OUTPUTS:
C	tt(2)		integer		time at which Coriolis was at specified orbit
C CALLS:
C	smei_orbit_info, Time2YDoy, Time2Add
C PROCEDURE:
C	Orbit 1 starts at t0 = 2002/12/31 23:56:41.900.
C	The orbital period at this time was 0.0705612268519d0 days,
C	and has been decreasing since with 0.084 msec per orbit.
C	(0.42s reduction in SMEI period per 365 days).
C
C	The effective orbital period for orbit n is p0+n*dp0
C	t(0) = t0
C	t(1) = t(  0)+p0+1*dp0
C	t(2) = t(  1)+p0+2*dp0
C	...
C	t(n) = t(n-1)+p0+n*dp0 = t0+n*p0+0.5*n*(n+1)*dp0
C MODIFICATION HISTORY:
C	JUN-2004, Paul Hick (UCSD/CASS; pphick@ucsd.edu)
C-
c	    double precision	orbit
c	    integer		tt(2)

c	integer			t0(2)
c	double precision	p0
c	double precision	dp0
c	double precision	dt
c	double precision	ds

c	call smei_orbit_info(n0,t0,p0,dp0)

c	ds = orbit-n0
c	dt = ds*p0			! Days since t0
c	ds = -ds*(ds+1d0)*0.5d0*dp0	! Additional milli-seconds since t0

c	tt(1) = int(dt)
c	tt(2) = nint( (dt-dble(tt(1)))*86400d3+ds )

c	call Time2Add(t0,tt,tt)

c	return
c	end
C+
C NAME:
C	smei_orbit_period
C PURPOSE:
C	Coriolis orbital period
C CALLING SEQUENCE:
c	double precision function smei_orbit_period(orbit)
C INPUTS:
C	orbit		double precision
C					Coriolis orbit number
C OUTPUTS:
C	smei_orbit_period
C			double precision
C					orbital period in days
C CALLS:
C	smei_orbit_info
C PROCEDURE:
C	See href=smei_orbit_time=
C MODIFICATION HISTORY:
C	JUN-2004, Paul Hick (UCSD/CASS; pphick@ucsd.edu)
C-
c	    double precision	orbit

c	integer			t0(2)
c	double precision	p0
c	double precision	dp0

c	call smei_orbit_info(n0,t0,p0,dp0)

c	smei_orbit_period = p0-((orbit-dble(n0))*dp0)/86400d3

c	return
c	end
C+
C NAME:
C	smei_orbit
C PURPOSE:
C	Coriolis orbit number at specified time
C CALLING SEQUENCE:
c	double precision function smei_orbit(tt)
C INPUTS:
C	tt(2)		integer			time
C OUTPUTS:
C	smei_orbit	double precision	Coriolis orbit number
C CALLS:
C	Time2YDoy, Time2Delta, Time2Day8, smei_orbit_period, smei_orbit_time
C	smei_orbit_info
C PROCEDURE:
C	See href=smei_orbit_time=
C MODIFICATION HISTORY:
C	JUN-2004, Paul Hick (UCSD/CASS; pphick@ucsd.edu)
C-
c	    integer		tt(2)

c	integer			uu(2)
c	double precision	p0
c	double precision	dp0
c	double precision	days
c	double precision	dorbit
c	double precision	smei_orbit_period

c	call smei_orbit_info(n0, uu, p0, dp0)
c	call Time2Delta(tt,uu,uu)	! Time difference 

	! First guess at orbit number using the nominal orbital period
	! This will underestimate the number of orbits elapsed since
	! the orbital period decreases over time

c	call Time2Day8(0,uu,days)

	! Approximate # orbits elapsed

c	dorbit = dble(n0)
c	smei_orbit = dorbit+days/smei_orbit_period(dorbit)
c	smei_orbit = dorbit+days/smei_orbit_period((dorbit+smei_orbit)*0.5d0)

	! The time uu at which the underestimated number of orbits
	! are elapsed is earlier than the input time tt

c	call smei_orbit_time(smei_orbit, uu)	! Exact time for 'orbit'
c	call Time2Delta(tt,uu,uu)		! Difference in days
c	call Time2Day8(0,uu,dorbit)		! dorbit is positive

	! Correct the approximate number of elapsed orbits
	! using the instantaneous orbital period

c	smei_orbit = smei_orbit+dorbit/smei_orbit_period(smei_orbit)

c	return
c	end