;+
; NAME:
;	big_orbit
; PURPOSE:
;	Get coordinates arrays for the orbit of one of the
;	bodies supported by the href=big_eph= ephemeris.
; CATEGORY:
;	gen/idl/ephem
; CALLING SEQUENCE:
	FUNCTION big_orbit, TT		, $
		body		= body		, $
		precess 	= precess	, $
		to_sphere	= to_sphere , $
		degrees 	= degrees	, $
		silent		= silent	, $
		_extra		= _extra
; INPUTS:
;	TT				array[1]; type: standard time structure
; OPTIONAL INPUT PARAMETERS:
;	body=body		scalar; type: string; default: jpl_body(/earth)
;	/precess		passed to big_eph
;					precesses the orbital locations from
;					J2000 to epoch TT
;	/to_sphere		by default rectangular coordinates are returned
;					If /to_sphere is set these are converted
;					to spherical coordinates.
;	/degrees		if set, angles returned when /to_sphere is
;					set are in degrees; default is radians.
;	_extra=_extra	passed to href=CvSky=
;					this should be one of the to_* keywords
;					allowed for CvSky and determines the
;					reference frame of the orbital coordinates
; OUTPUTS:
;	Result		array[3,*]	type: float
;					coordinate arrays for orbit of 'body'
;					The unit of length is AU; if /to_sphere
;					is set the units for the angles depends
;					on the setting of /degrees.
;					If the body is not recognized then -1 is
;					returned	
; INCLUDE:
	@compile_opt.pro			; On error, return to caller
; CALLS:
;	InitVar, jpl_body, usno_body, mpc_body, TimeSet, TimeOp, gridgen
;	big_eph, CvSky, mpc_eph_range
; PROCEDURE:
;	The orbit is set up by using NewcombSun on an array of times
; 	covering one Julian year centered on TT.
;	Then CvSky is used to transform to heliographic coordinates.
; MODIFICATION HISTORY:
;	JUN-1997, Paul Hick (UCSD/CASS)
;	OCT-2002, Paul Hick (UCSD/CASS)
;		Fixed bug in orbit calculation by adding /heliocentric
;		keyword to NewcombSun. The problem only was visible
;		for viewing locations extremely close (within a few
;		earth radii) of Earth.
;	JUN-2004, Paul Hick (UCSD/CASS)
;		Switched to JPL ephemeris instead of using NewcombSun.
;		This allowed adding orbits for other planets and Moon.
;		Renamed from RemoteView_EarthOrbit to RemoteView_PlanetOrbit
;	OCT-2006, Paul Hick (UCSD/CASS)
;		Now using big_eph instead of jpl_eph. Allows adding of
;		partial orbits for comets and asteroids.
;		Renamed to RemoteView_BodyOrbit
;	OCT-2006, Paul Hick (UCSD/CASS; pphick@ucsd.edu)
;		Better documentation. Renamed to big_orbit.
;-

InitVar, body, jpl_body(/earth,/string)
InitVar, to_sphere, /key

; T covers one Julian year in 365 steps. The number of steps must
; be uneven to make sure that TT itself is on the list.

; Orbital period in Julian years

CASE 1 OF
jpl_body (body		  ) NE -1	: P = ([0.0,0.2408,0.6152,1.0,1.8809,11.862,29.458,84.01,164.79,248.54,0.074743])[jpl_body(body)]
usno_body(body,/silent) NE ''	: P =  1.0
body EQ big_body('Ulysses')	 	: P =  1.0
body EQ big_body('Stereo A')	: P =  1.0
body EQ big_body('Stereo B')	: P =  1.0
body EQ big_body('Messenger')	: P =  1.0
mpc_body (body		  ) NE ''	: P = -1.0 
ELSE: BEGIN
	IF silent LE 0 THEN message, /info, "ignoring unknown body '"+strcompress(body,/rem)+"'"
	P = 0.0
END 
ENDCASE

IF P NE 0.0 THEN BEGIN

	IF P GT 0.0 THEN BEGIN

		range = P*!earth.julianyr/2*[-1,1]
		T = TimeOp(/add, TT, TimeSet(day=gridgen(1001, range=range)))

	ENDIF ELSE BEGIN

		; For a comet (P=-1) the ephemeris file (an html file downloaded
		; from MPC) only covers a restricted time range. The time array
		; set up here covers the whole period. No orbit locations are
		; returned if the input time TT is outside the range covered by
		; the comet ephemeris.

		range = mpc_eph_range(body,/silent)
		T = TimeOp(/subtract, range, TT, TimeUnit(/day))
		IF T[0] LE 0 AND T[1] GE 0 THEN T = gridgen(1001,range=range) ELSE destroyvar, T

	ENDELSE

	IF IsTime(T) THEN BEGIN

		; Set up the orbit by getting locations along the orbit in
		; a sidereal coordinate frame (equatorial) at all times TT
		; (precess if necessary)

		Orbit = big_eph(T, body=body, precess=precess, /to_equatorial,	$
			to_sphere=to_sphere, degrees=degrees, silent=silent)
		Orbit = reform(Orbit, /overwrite)

		; Convert from equatorial coordinates to coordinates
		; specified in _extra at single input time TT.
		; Note that it is NOT a good idea to use big_eph to go
		; straight to the coordinate frame specified in _extra.

		Orbit = CvSky(TT, from_equatorial=Orbit, _extra=_extra,		$
			rectangular=1-to_sphere, degrees=degrees, silent=silent)

	ENDIF ELSE BEGIN

		Orbit = -1

	ENDELSE

ENDIF ELSE BEGIN

	Orbit = -1

ENDELSE

RETURN, Orbit  &  END