PRO SPLOT
; IDL main program to read and plot DMSP electron spectra saved by DMSPSPEC.PAS
; D P Steele, 12 Sept 1990
; Plot of mean + envelope added 21 Sept 1990 by D P Steele
; Listing of whole file name, and variable energy flux format, added 24 Sept 90
; Plotting of 1-sigma limits, and option not to plot, added 26 April 1992
; Converted to a procedure to faciilitate calling from within a main
; program, DPS, 27 Aug 1992
;
F20 = FLTARR( 20)	; variable type for spectra and energy steps
E = F20
Struct = { Spectrum, UT : 120000, EFlux : 1.0, AvgE : 1.0, Fluxes : F20 }
Str = ''
Filename = Str
Lab0 = Str
Lab1 = Str
N = 1
PRINT, 'Enter the spectrum file name'
READ, Filename
OPENR, 1, Filename
READF, 1, Lab0
READF, 1, N
Spec = REPLICATE( { Spectrum }, N)
Spectra = FLTARR( 20, N)	; to hold spectra only
READF, 1, Lab1
READF, 1, E
READF, 1, Lab1
READF, 1, Lab1
UT = 120000	; initialize variables to set their types
EFlux = 1.0
AvgE = 1.0
NFlux = 1.0
J = F20
MaxEFlux = 0.0
MaxFlux = 1e5
K = 0
K2 = 0
FOR I = 0, N - 1 DO BEGIN
    READF, 1, UT, EFlux, AvgE, NFlux, J, $
              FORMAT = '(I10,3F12.5,3(/6F12.5)/2F12.5)'
    IF (EFlux GT 0) THEN BEGIN		; accept only non-zero spectra
        Spec( K).UT = UT		; hhmmss format
        EFlux = EFlux / !Pi		; erg/cm2 s sr
        Spec( K).EFlux = EFlux
        MaxEFlux = MAX( [MaxEFlux, EFlux])	; highest energy flux yet
        Spec( K).AvgE = AvgE		; keV
        Spec( K).Fluxes = J		; el/cm2 s sr eV
        MaxFlux = MAX( [MaxFlux, MAX( J)])	; highest spectrum point yet
        K = K + 1
        IF Max( J) GT 0.0 THEN BEGIN
            Spectra( *, K2) = J		; same units as above
            K2 = K2 + 1
        ENDIF
    ENDIF
ENDFOR
CLOSE, 1
Sup = F20 & Inf = F20 & SDUp = F20 & SDDown = F20 
Spectra = Spectra( *, 0:K2-1)
MeanEFlux = Spec( 0).EFlux
MeanNFlux = Spec( 0).EFlux / Spec( 0).AvgE	; strange units - don't worry
MeanSpec = Spec( 0).Fluxes
FOR I = 1, K - 1 DO BEGIN
    MeanEFlux = MeanEFlux + Spec( I).EFlux
    MeanNFlux = MeanNFlux + Spec( I).EFlux / Spec( I).AvgE
    MeanSpec = MeanSpec + Spec( I).Fluxes
ENDFOR
MeanAvgE = MeanEFlux / MeanNFlux
MeanEFlux = MeanEFlux / K
MeanSpec = MeanSpec / K
FOR I = 0, 19 DO BEGIN
    Sup( I) = MAX( Spectra( I, *)) & Inf( I) = MIN( Spectra( I, *))
    SD = StDev( Spectra( I, *))
    SDUp( I) = MeanSpec( I) + SD & SDDown( I) = MeanSpec( I) - SD
ENDFOR
; Offer user option to PRINT results rather than PLOT them
Resp=''
PRINT,'There are'+STRCOMPRESS(K2)+' spectra; do you want to plot them?'
READ,Resp
Resp=STRLOWCASE(Resp)
Plot = ((Resp EQ 'y') OR (Resp EQ 'Y'))
IF Plot THEN BEGIN
    !X.TITLE = 'Energy ( eV )'
    !Y.TITLE = 'Flux ( (cm!u2!n s sr eV)!u-1!n )'
    !P.TITLE = STRUPCASE( STRMID( Lab0, 27, STRLEN( Lab0) - 27))
    PRINT,'Present plot title is:'
    PRINT,!P.TITLE
    resp=''
    PRINT,'Do you want to enter a new one? (Y/[N])'
    READ,resp
    IF STRUPCASE(resp) EQ 'Y' THEN BEGIN
	ntitle=''
	PRINT,'Enter the new title now.'
	READ,ntitle
	!P.TITLE=ntitle
    ENDIF
    !X.RANGE = [1,1e5]
    !Y.RANGE = [MaxFlux/1.1E5,MaxFlux]
    DivFact=1.0+!P.MULTI(1)/2.
    !P.CHARSIZE = 1.0/DivFact
    LegSize = 1.0
    IF (!D.NAME EQ 'PS') THEN $
        LegSize = !D.X_SIZE/FLOAT(!D.Y_SIZE*DivFact) $	; PostScript output
        ELSE LegSize = 1.0/DivFact
    LStyle = INDGEN( 5)
    LStyle( 4) = 5		; avoid one line type on TEK
    Lx = 1.5
    Rx = 3
    Cx = 4
ENDIF
IF MaxEFlux GT 1E-5 THEN EFFormat = 'F7.5'
IF MaxEFlux GT 1E-2 THEN EFFormat = 'F7.4'
IF MaxEFlux GT 1E-1 THEN EFFormat = 'F7.3'
IF Plot THEN BEGIN
    FOR J = 0, K - 1, 5 DO BEGIN	; limit of 5 different linestyles
        PLOT_OO, E, Spec( J).Fluxes > 0.1
        YLo = 10^!Y.CRANGE( 0)
        XYOUTS, Lx, 64*YLo, 'Curve    UT      EFlux    AvgE', SIZE = LegSize
        Ty = 32*YLo
        PLOTS, [Lx, Rx], [Ty, Ty]
        Legend = STRING( Spec( J).UT, Spec( J).EFlux, Spec( J).AvgE, $
                         FORMAT = '(I6, 2X, ' + EFFormat + ', 2X, F5.2)')
        XYOUTS, Cx, Ty, Legend, SIZE = LegSize
        Lims = [4, K - 1 - J]
        FOR I = 1, MIN( Lims) DO BEGIN
            OPLOT, E, Spec( I + J).Fluxes > 0.1, Line = LStyle( I)
            Ty = Ty / 2.0
            PLOTS, [Lx, Rx], [Ty, Ty], Linestyle = LStyle( I)
            Legend = STRING( Spec( I + J).UT, Spec( I + J).EFlux, $
                             Spec( I + J).AvgE, $
                             FORMAT = '(I6, 2X, ' + EFFormat + ', 2X, F5.2)')
            XYOUTS, Cx, Ty, Legend, SIZE = LegSize
        ENDFOR
        STOP	; for plotting
    ENDFOR
    PLOT_OO, E, MeanSpec > 0.1
    OPLOT, E, Sup > 0.1, LINESTYLE = 1
    OPLOT, E, Inf > 0.1, LINESTYLE = 1
ENDIF
Title = 'Mean spectrum with envelope of ' + STRCOMPRESS( K) + ' spectra'
IF Plot THEN XYOUTS, Lx, 32*YLo, Title, SIZE = LegSize ELSE PRINT,Title
Title = 'from ' + STRCOMPRESS( Spec( 0).UT) + ' to ' + $
        STRCOMPRESS( Spec( K - 1).UT) + ' UT'
IF Plot THEN XYOUTS, Lx, 16*YLo, Title, SIZE = LegSize ELSE PRINT, Title
Title = STRING( MeanEFlux, ' ergs/(cm!u2!n s sr)  ', MeanAvgE, ' keV', $
                FORMAT = '(' + EFFormat + ', A22, F5.2, A4)')
IF Plot THEN XYOUTS, Lx, 8*YLo, Title, SIZE = LegSize ELSE PRINT, Title
IF Plot THEN BEGIN
    STOP
    PLOT_OO, E, MeanSpec > 0.1
    OPLOT, E, SDUp > 0.1, LINESTYLE = 1
    OPLOT, E, SDDown > 0.1, LINESTYLE = 1
ENDIF
spint, E, SDUp, UpJE, UpAvgE
IF UpAvgE GT 35. THEN UpAvgE = UpAvgE/1000.	; convert to keV
spint, E, SDDown, DownJE, DownAvgE
IF DownAvgE GT 35. THEN DownAvgE = DownAvgE/1000.	; convert to keV
TitlePl = 'Mean spectrum with !9+!x 1 !4r!x limits of ' + $
          STRCOMPRESS( K, /REMOVE_ALL) + ' spectra'
TitlePr = 'Mean spectrum with +/- 1 sigma limits of ' + $
          STRCOMPRESS( K, /REMOVE_ALL) + ' spectra'
IF Plot THEN XYOUTS, Lx, 32*YLo, TitlePl, SIZE = LegSize ELSE PRINT, TitlePr
Title = 'from ' + STRCOMPRESS( Spec( 0).UT) + ' to ' + $
        STRCOMPRESS( Spec( K - 1).UT) + ' UT'
IF Plot THEN XYOUTS, Lx, 16*YLo, Title, SIZE = LegSize ELSE PRINT, Title
Title = STRING( UpJE, ' ergs/(cm!u2!n s sr)  ', UpAvgE, ' keV', $
                FORMAT = '(' + EFFormat + ', A22, F5.2, A4)')
IF Plot THEN XYOUTS, Lx, 8*YLo, Title, SIZE = LegSize ELSE PRINT, Title
Title = STRING( DownJE, ' ergs/(cm!u2!n s sr)  ', DownAvgE, ' keV', $
                FORMAT = '(' + EFFormat + ', A22, F5.2, A4)')
IF Plot THEN XYOUTS, Lx, 4*YLo, Title, SIZE = LegSize ELSE PRINT, Title
RETURN
END