;+
; Project     :	Computational Physics
;
; Name        :	ORBPLOT
;
; Purpose     :	Creates plot environment and axes for polar plots of an orbit.
;
; Explanation :	This procedure creates the plot grid with labeled x & y axes
;               with the Sun at the origin.  This procedure makes prettier axes
;               than the IDL AXIS procedure and accounts for orbits with large
;               eccentricities.
;
; Use         :	ORBPLOT, ASEMI, ORBECC, MTITLE [, UNITS]
;
; Inputs      :	ASEMI:  The semimajor axis of the orbit in units in UNITS.
;
;               ORBECC: The eccentricity of the orbit.
;
;               MTITLE: The main title of the plot.
;
; Opt. Inputs :	UNITS:  The units that the x & y axes are in.  [DEFAULT = 'AU']
;                       Note that if a unit other than AU is passed, a
;                       conversion is made internally to AU.  The plots
;                       generated in this routine are always in AU.
;                       Allowed values are: 'AU' (astronomical units), 'km'
;                       (kilometers), 'm' (meters), and 'cm' (centimeters).
;
; Outputs     :	None.
;
; Opt. Outputs:	None.
;
; Keywords    :	AXRANGE:  Four element array containing [!X.CRANGE, !Y.CRANGE].
;                       This is meant to be used for output.  It is ignored on
;                       input.
;
;               HARDCOPY: If set, the plot is being generated for a hardcopy
;                       postscript printer.  When this is set, the lines are
;                       drawn at 3 times the normal thickness.
;
;
;               SCDIST: Scale factor for the distance if AU not used.
;
;               ERRMSG: If defined and passed, then any error messages will
;                       be returned to the user in this parameter rather than
;                       being internally printed.  If no errors are
;                       encountered, then a null string is returned.  In order
;                       to use this feature, the string ERRMSG must be defined
;                       first, e.g.,
;
;                           ERRMSG = ''
;                           ORBPLOT, 5.2, 0.012, 'AU', ERRMSG=ERRMSG
;                           IF ERRMSG(0) NE '' THEN ...
;
; Calls       :	ORBTICKS.
;
; Common      :	None.
;
; Restrictions:	None.
;
; Side effects:	None.
;
; Category    :	Computation physics.
;
; Prev. Hist. :	None.
;
; Written     :	Donald G. Luttermoser, ETSU/Physics, 9 December 2002
;
; Modified    :	Version 1, Donald G. Luttermoser, ETSU/Physics, 9 Dec 2002
;                   Initial program.
;               Version 2, Donald G. Luttermoser, ETSU/Physics, 3 May 2007
;                   Replaced the IDL MESSAGE utility with a PRINT statement.
;               Version 3, Donald G. Luttermoser, ETSU/Physics, 23 Nov 2009
;                   Changed variable array indicies symbols from "()" to "[]".
;
; Version     :	Version 3,  23 November 2009.
;
;-
;
PRO ORBPLOT, ASEMI, ORBECC, MTITLE, UNITS, AXRANGE=AXRANGE, HARDCOPY=HARDCOPY, $
    SCDIST=SCDIST, ERRMSG=ERRMSG
;
EMESS = ''	; Set to non-null string if error is encountered.
;
CASE N_PARAMS() OF
    3:  SUNITS = '!6AU'
    4:  SUNITS = '!6' + STRTRIM(STRING(UNITS), 2)
 ELSE:  BEGIN
        EMESS = 'Syntax:  ORBPLOT, ASEMI, ORBECC, MTITLE [, UNITS, ' + $
            'AXRANGE=AXRANGE]'
        GOTO, HANDLE_ERROR
        END
ENDCASE
;
; Check keywords.
;
IF KEYWORD_SET(HARDCOPY) THEN PASS = 1 ELSE PASS = 0
IF PASS EQ 1 THEN CHSC = 0.6 ELSE CHSC = 1.0
;
; Convert to Astronomical Units if needed.
;
CASE STRUPCASE(SUNITS) OF
    '!6KM': SCDIST = 6.6845E-09
    '!6M':  SCDIST = 6.6845E-12
    '!6CM': SCDIST = 6.6845E-14
    '!6AU': SCDIST = 1.0
    ELSE:   BEGIN
            SCDIST = 0.0
            EMESS = 'Unknown unit requested: ' + SUNITS
            GOTO, HANDLE_ERROR
            END
ENDCASE
SEMIMAJ = ASEMI * SCDIST
;
; Internal parameters.
;
RPERI = SEMIMAJ * (1.0D0 - ORBECC)
RAP = SEMIMAJ * (1.0D0 + ORBECC)
BSEMI = SEMIMAJ * SQRT(1.0D0 - ORBECC^2)
IF ORBECC GE 0.7 THEN SCPERI = 10. ELSE SCPERI = 1.2
SDOT = FINDGEN(16) * 2*!PI / 15.
SUNSIZE = 2.0
XSUN = SUNSIZE * COS(SDOT)
YSUN = SUNSIZE * SIN(SDOT)
USERSYM, XSUN, YSUN, /FILL
;
; Figure out major tick mark locations.
;
ORBTICKS, RAP, XTICKS, SXTICKS, NXTIC, NXMINOR, XMINSEP, ERRMSG=EMESS
IF EMESS NE '' THEN GOTO, HANDLE_ERROR
ORBTICKS, BSEMI, YTICKS, SYTICKS, NYTIC, NYMINOR, YMINSEP, ERRMSG=EMESS
IF EMESS NE '' THEN GOTO, HANDLE_ERROR
;
; Plot axes.
;
PLOT, [-SCPERI*RPERI, 1.2*RAP], [-1.2*BSEMI, 1.2*BSEMI], XSTYLE=4+1, $
    YSTYLE=4+1, /NODATA
AXRANGE = [!X.CRANGE, !Y.CRANGE]
OPLOT, !X.CRANGE, [0, 0], THICK=3*PASS+2
OPLOT, [0, 0], !Y.CRANGE, THICK=3*PASS+2
XYOUTS, 1.18*RAP, 1.18*BSEMI, MTITLE, SIZE=1.3*CHSC, ALIGN=1.0
;
XDIFF = !X.CRANGE[1] - !X.CRANGE[0]
YDIFF = !Y.CRANGE[1] - !Y.CRANGE[0]
;
; Plot tick marks along x axis.
;
FOR I = 0, NXTIC-1 DO BEGIN
;
; Positive x goes out to RAP.
;
    IF XTICKS[I] LE !X.CRANGE[1] THEN BEGIN
        OPLOT, [XTICKS[I], XTICKS[I]], [0., 0.03*YDIFF], THICK=3*PASS+2
        XYOUTS, XTICKS[I], -0.05*YDIFF, SXTICKS[I], SIZE=1.2*CHSC, ALIGN=0.5
    ENDIF
;
; Negative x goes out to RPERI.
;
    IF XTICKS[I] LE ABS(!X.CRANGE[0]) THEN BEGIN
        OPLOT, [-XTICKS[I], -XTICKS[I]], [0., 0.03*YDIFF], THICK=3*PASS+2
        XYOUTS, -XTICKS[I], -0.05*YDIFF, SXTICKS[I], SIZE=1.2*CHSC, ALIGN=0.5
    ENDIF
ENDFOR
;
; Plot tick marks along y axis.
;
FOR I = 0, NYTIC-1 DO BEGIN
;
; +/- y goes out to +/- BSEMI.
;
    IF YTICKS[I] LE !Y.CRANGE[1] THEN BEGIN
        OPLOT, [0., 0.03*XDIFF], [YTICKS[I], YTICKS[I]], THICK=3*PASS+2
        XYOUTS, -0.04*XDIFF, YTICKS[I]-0.01*YDIFF, SYTICKS[I], SIZE=1.2*CHSC, $
            ALIGN=0.5
        OPLOT, [0., 0.03*XDIFF], [-YTICKS[I], -YTICKS[I]], THICK=3*PASS+2
        XYOUTS, -0.04*XDIFF, -YTICKS[I]-0.01*YDIFF, SYTICKS[I], SIZE=1.2*CHSC,$
            ALIGN=0.5
    ENDIF
ENDFOR
;
; Plot the minor tick marks.
;
XPTMIN = XTICKS[0]
WHILE XPTMIN LE !X.CRANGE[1] DO BEGIN
    OPLOT, [XPTMIN, XPTMIN], [0., 0.02*YDIFF], THICK=3*PASS+2
    XPTMIN = XPTMIN + XMINSEP
ENDWHILE
XMTMIN = XTICKS[0]
WHILE XMTMIN LE ABS(!X.CRANGE[0]) DO BEGIN
    OPLOT, [-XMTMIN, -XMTMIN], [0., 0.02*YDIFF], THICK=3*PASS+2
    XMTMIN = XMTMIN + XMINSEP
ENDWHILE
YTMIN = YTICKS[0]
WHILE YTMIN LE !Y.CRANGE[1] DO BEGIN
    OPLOT, [0., 0.015*XDIFF], [YTMIN, YTMIN], THICK=3*PASS+2
    OPLOT, [0., 0.015*XDIFF], [-YTMIN, -YTMIN], THICK=3*PASS+2
    YTMIN = YTMIN + YMINSEP
ENDWHILE
;
; Plot and label Sun.
;
XYOUTS, 0.03*XDIFF, 0.05*YDIFF, '!6Sun', SIZE=1.3*CHSC
OPLOT, [100*RAP, 0], [100*RAP, 0], PSYM=8
;
; Label axes.
;
XYOUTS, 0.2*XDIFF, -0.10*YDIFF, '!6r ('+SUNITS+')', SIZE=1.3*CHSC
;
IF EMESS NE '' THEN GOTO, HANDLE_ERROR
IF N_ELEMENTS(ERRMSG) GT 0 THEN ERRMSG = EMESS
RETURN
;
; Error handling portion of the procedure.
;
HANDLE_ERROR:
IF N_ELEMENTS(ERRMSG) NE 0 THEN ERRMSG = EMESS ELSE PRINT, EMESS
;
RETURN
END