;+
; NAME:
;   ocpluto_pos
;
; PURPOSE:
;   Calculates position of Pluto, Charon, barycenter, and photocenter
;
; DESCRIPTION:
;   Calculates position of Pluto, Charon, barycenter, and photocenter
;
; CALLING SEQUENCE:
;
; INPUTS:
;   et - the ET at which you want to know the coordinates (may be an
;        array).  Units are seconds of TDB (barycentric dynamical
;        time) after J2000.
;
; KEYWORD INPUT
;   lonlatalt - lon, lat, altitude for topocentric
;               ephemerides. Default geocentric
; KEYWORD OUTPUT:
;   pluto - pluto body center. pluto[0,*] = ra, etc
;   Charon - Charon body center
;   barycener - barycenter
;   photocenter - system photocenter
;
; RESTRICTIONS
;   For geocentric observations, ltime and position can be off by
;    Re/c = 21 ms.
; 
;  REVISON HISTORY:
;   2006 Jan 02 Leslie Young.  Based on single_eph
;   2006 May 11 LAY Added ltime
;-

pro oc_pluto_pos, et, lonlatalt = lonlatalt, $
                  fluxrat_cp = fluxrat_cp, $
                  pluto=pluto, charon=charon, $
                  barycenter=barycenter, photocenter=photocenter
    
    net = n_elements(et)

    plutoeph = oc_ephicy_et(999, 399, et, neph, $
                    lonlatalt=lonlatalt)
    pluto = dblarr(2,net)
    pluto[0,*] = plutoeph.ra
    pluto[1,*] = plutoeph.dec
    charoneph = oc_ephicy_et(901, 399, et, neph, $
                    lonlatalt=lonlatalt)
    charon = dblarr(2,net)
    charon[0,*] = charoneph.ra
    charon[1,*] = charoneph.dec
    barycentereph = oc_ephicy_et(9, 399, et, neph, $
                    lonlatalt=lonlatalt)
    barycenter = dblarr(2,net)
    barycenter[0,*] = barycentereph.ra
    barycenter[1,*] = barycentereph.dec
    
    if not keyword_set(fluxrat_cp) then fluxrat_cp = 0.18
    
    photocenter = dblarr(2,net)
    photocenter[0,*] = fluxrat_cp*charon[0,*] + (1-fluxrat_cp)*pluto[0,*]
    photocenter[1,*] = fluxrat_cp*charon[1,*] + (1-fluxrat_cp)*pluto[1,*]
    
end