;;;read in counts vs. MLAT @12 MLT data file
;;;do Gaussian fit and calculate the fehm.

;openr, lun,'./CHAMP_360km_AACGM_peaks_30%_normalized@12MLT.dat',/get_lun
;openr, lun,'./combine_CHAMPGRACE_AACGM_peaks_30%_normalized@12MLT.dat',/get_lun
openr, lun,'../../GRACE/GRACE_470km_maxden_30%/GRACE_470km_AACGM_peaks_30%_normalized@12MLT.dat',/get_lun
mlat = fltarr(200)
counts = fltarr(200)
m = 0.0
c = 0.0
npt = 0
while (not eof(lun)) do begin
   readf,lun,m,c
   mlat(npt) = m
   counts(npt) = c
   npt = npt+1
endwhile
mlat = mlat(0:npt-1)
counts = counts(0:npt-1)

set_plot,'ps'
device,filename = 'fwhm.ps'
device,/color
loadct, 39

fp = where((mlat lt 0.) and (mlat gt -85.))
mlat1 = mlat(fp)
counts1 = counts(fp)

yfit = GaussFit(mlat1, counts1, coeff, NTERMS=4)
fwhm = 2 * SQRT(2 * ALOG(2)) * coeff(2)

plot, mlat, counts,xr =[-90, 0],xstyle = 1,xticks = 9,thick = 2, $
      xtitle='MLAT (SH)', ytitle='counts' 
oplot, mlat1, yfit, color = 250, linestyle = 2, thick = 3
;xyouts, -50,0.02,'Maximum:'+string(coeff(0),format='(F0.2)')
xyouts, -50,0.015,'Center:'+string(coeff(1),format='(F0.2)')
xyouts, -50,0.010,'FWHM:'+string(fwhm,format='(F0.2)')

fp = where((mlat gt 0.) and (mlat lt 85.))
mlat2 = mlat(fp)
counts2 = counts(fp)

yfit = GaussFit(mlat2, counts2, coeff, NTERMS=4)
fwhm = 2 * SQRT(2 * ALOG(2)) * coeff(2)

plot, mlat, counts, xr=[0, 90], xstyle = 1,xticks = 9,thick = 2, $
       xtitle='MLAT (NH)', ytitle='counts'
oplot, mlat2, yfit, color = 250, linestyle = 2, thick = 3
;xyouts, 40,0.02,'Maximum:'+string(coeff(0),format='(F0.2)')
xyouts, 40,0.015,'Center:'+string(coeff(1),format='(F0.2)')
xyouts, 40,0.010,'FWHM:'+string(fwhm,format='(F0.2)')
device,/close
set_plot,'x'
end