;+
; Project     :	COURSES
;
; Name        :	HW3415ION
;
; Purpose     :	Make a plot for Problem 8.10c in the Carroll & Ostlie textbook.
;
; Explanation :	This procedure makes a plot for Problem 8.10c in the Carroll &
;		Ostlie textbook used in the Astrophysics (ASTR-3415) course.
;
; Use         :	HW3415ION (no arguments)
;
; Inputs      :	None.
;
; Opt. Inputs :	None.
;
; Outputs     :	None.
;
; Opt. Outputs:	None.
;
; Keywords    :	ERRMSG:	If defined and passed, then any error messages will 
;			be returned to the user in this parameter rather than
;			being printed internally by this procedure.  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 = ''
;				HW3415ION, ERRMSG=ERRMSG
;				IF ERRMSG(0) NE '' THEN ...
;
; Calls       :	None.
;
; Common      :	None.
;
; Restrictions:	None.
;
; Side effects:	None.
;
; Category    :	Saha equation, ionization, homework.
;
; Prev. Hist. :	None.
;
; Written     :	Donald G. Luttermoser, ETSU/Physics, 20 March 2003
;
; Modified    :	Version 1, Donald G. Luttermoser, ETSU/Physics, 20 March 2003
;			Initial program.
;
; Version     :	Version 1,  20 March 2003.
;
;-
;
PRO HW3415ION, ERRMSG=ERRMSG
;
ON_ERROR, 2   ; Return to caller if error is encountered.
EMESSAGE = ''  ; Set to non-null string if error is encountered.
;
; Check input.
;
IF N_PARAMS() NE 0 THEN BEGIN
	EMESSAGE = 'Syntax:  HW3415ION (no arguments)'
	GOTO, HANDLE_ERROR
ENDIF
;
; ************************ Global Internal Parameters ************************
;
TERM = !D.NAME  &  AFONT = [-1, 0]  &  PASS = 0
;
; Physical constants and parameters for helium.
;
ME = 9.1093897D-28  &  KB = 1.380658D-16  &  EV = 1.60217733D-12
H = 6.6260755D-27
CHI = [24.6D0, 54.4D0, 0.0D0]  &  ZPART = [1.0D0, 2.0D0, 1.0D0]
;
TEMP = FINDGEN(1000) * 20 + 5000.  &  PELEC = 2.0D2
;
; Ratio of first ionization stage to neutral stage.
;
N2TO1 = ((2. * KB * TEMP * ZPART(1)) / (PELEC * ZPART(0))) * $
	((2. * !DPI * ME * KB * TEMP) / H^2)^(3./2.) * $
	EXP(-CHI(0)*EV/(KB*TEMP))
;
; Ratio of second ionization stage to first ionization stage.
;
N3TO2 = ((2. * KB * TEMP * ZPART(2)) / (PELEC * ZPART(1))) * $
	((2. * !DPI * ME * KB * TEMP) / H^2)^(3./2.) * $
	EXP(-CHI(1)*EV/(KB*TEMP))
;
; Ratio of second ionization stage to total species density.
;
N2TOTOT = 1.0D0 / (1.0D0 + N3TO2 + 1.0D0/N2TO1)
;
; Plot characteristics.
;
!P.CHARSIZE = 1.5  &  !X.CHARSIZE = 1.0  &  !Y.CHARSIZE = 1.0
!P.THICK = 2.4  &  !X.THICK = 2.4  &  !Y.THICK = 2.4
;
NEXTPLOT:
;
PLOT, TEMP/1000, N2TOTOT, THICK=3*(2*PASS+1), XTITLE='!6Temperature (1000 K)',$
	YTITLE='!6N!DII!N/N!Dtotal!N', XSTYLE=1, YSTYLE=1, XRANGE=[5,25], $
	YRANGE=[-0.02,1.02], FONT=AFONT(PASS)
;
IF PASS EQ 0 THEN BEGIN
	SCREEN, 'PS', /ENCAPS, /SQUARE
	DEVICE, /SCHOOLBOOK, /BOLD, FONT_INDEX=17
	DEVICE, /SYMBOL, FONT_INDEX=7
	DEVICE, /BKMAN, /DEMI, FONT_INDEX=6
	DEVICE, /HELVETICA, /BOLD, FONT_INDEX=5
	PASS = 1
	GOTO, NEXTPLOT
ENDIF ELSE BEGIN
	PASS = 0
	DEVICE, /CLOSE_FILE
	SCREEN, TERM
	SPAWN, 'mv idl.eps hw3415ion.eps'
	PRINT, 'Plot saved in file:  hw3415ion.eps'
ENDELSE
;
IF EMESSAGE NE '' THEN GOTO, HANDLE_ERROR
IF N_ELEMENTS(ERRMSG) GT 0 THEN ERRMSG=EMESSAGE
RETURN
;
; Error handling portion of the procedure.
;
HANDLE_ERROR:
;
IF N_ELEMENTS(ERRMSG) EQ 0 THEN PRINT, EMESSAGE
ERRMSG = EMESSAGE
RETURN
;
END