pro eit_proton_summary, time0, time1, $
wave=wave, full=full, $ ; pass to > eit_files
eitfiles=eitfiles, $ ; optional filelist
cadence=cadence, $ ; optional cadence (mins_
subfov=subfov, threshold=threshold, $ ; pass to > cosmic_stat
nodisplay=nodisplay, www=www, movie=movie, $ ; display & www summary options
quiet=quiet, cstats=cstats, index=index, $
mwindow=mwindow
;+
;
; Name: eit_proton_summary
;
; Purpose: "standard" processing pipeline/summary plot for proton events
;
; Calls:
; plot_goes, plot_goesp, eit_files, eit_prep, cosmic_stat
; time_window, evt_grid, plus the usual ssw suspects.
;
; Input Parameters:
; time0, time1 - time range of interest
; wave - desired wavelength (default = 195)
; full - if set, use the FULL DISK data base
; subfov - desired box for cosmic ray counts - per cosmic_stat.pro
; [default = [0,0,1024,50]
; threshold - per cosmic_stat.pro
; nodisplay - if set, dont make the summary plots (just return cstat)
; mwindow - time window (minutes) +/- around EIT Peak for movie
; cstats (OUTPUT) - the affected area statistics (%pixels affected)
; index (OUTPUT) - index (header) vector
;-
;
if not keyword_set(wave) then wave='195'
full=keyword_set(full)
nodisplay=keyword_set(nodisplay)
display=1-nodisplay
movie=keyword_set(movie)
www=(keyword_set(www) or movie) and display
quiet=keyword_set(quiet)
; -------- find the files of interest
case 1 of
data_chk(eitfiles,/string): efiles=eitfiles ; user supplied list
n_params() eq 2: begin
if full then begin
efiles=eit_files(time0,time1,full=full, wave=wave)
sscnt=([n_elements(efiles),0])(efiles(0) eq '')
endif else begin
efiles=eit_files(time0,time1)
read_eit, efiles, header
ss=where(gt_tagval(header,/naxis1) eq 1024 and $
gt_tagval(header,/naxis2) eq 1024 and $
gt_tagval(header,/nmissb) le 1 and $
gt_tagval(header,/object,missing='') eq "full FOV", sscnt)
endelse
endcase
else: begin
box_message,['Must supply time range or EITFILES list',$
'IDL> eit_proton_summary,time0,time1 [OPTIONS]', $
' -OR-',$
'IDL> eit_proton_summary,eitfiles=FILELIST [OPTIONS]']
return
endcase
endcase
if sscnt eq 0 then begin
box_message,"No eit files match desired criteria, returning..."
return
endif
if keyword_set(cadence) then begin
ss=grid_data(file2time(efiles),minutes=cadence)
endif else if n_element(ss) eq 0 then ss=lindgen(n_elements(efiles))
read_eit, efiles(ss), index ; read headers
eit_prep, efiles(ss), pindex, pdata, /save_zero ; prep data
pindex=index
delvarx,data ; memory cleanup
if not keyword_set(subfov) then subfov=[0,0,data_chk(pdata,/nx),50]
if not keyword_set(threshold) then threshold=20
if not quiet then box_message,$
'Using subfield boxsize [x0,y0,nx,ny] >> [' + arr2str(subfov,/trim) + ']'
cstats=cosmic_stat(pindex, pdata, threshold, /percent, $
peak_spectra=peak_spectra, subfov=subfov)
psave=!p.multi
time_window,index, t0, t1
deit=deriv_arr([cstats(0),cstats])
epeakd=(where(deit eq max(deit)))(0)
epeak=(where(cstats eq max(cstats)))(0)
epeak50=where(cstats ge .5*cstats(epeak)) ; .5 max
epeakd=epeak50(0)
rd_gxd,t0,t1,/goes8,/five_min,gxd
gpeak=(where(gxd.lo eq max(gxd.lo)))(0)
dgoes=deriv_arr([gxd(0).lo,gxd.lo])
gpeakd=(where(dgoes eq max(dgoes)))(0)
dgp2ep=ssw_deltat(pindex(epeak),ref=gxd(gpeak),/hour)
dgr2er=ssw_deltat(pindex(epeakd), ref=gxd(gpeakd),/hour)
if display then begin
wdef,xx,zbuffer=www, 850,850
!p.multi=[0,1,3]
linecolors
plot_goes, t0, t1, background=11, color=7, xstyle=5, charsize=2.5, $
ymargin=[4,2],xmargin=[8,0]
evt_grid,pindex(epeakd),color=5,linestyle=0, thick=2, label='EIT Peak Event Half Max'+$
anytim(pindex(epeakd),/vms,/trunc),labpos=.95,align=0
plot_goesp, t0, t1, xstyle=5, /nowindow, charsize=2.5, /log, $
ymargin=[2,0], xmargin=[8,0], /proton_only
utplot, index, cstats, ytitle='Affected Area (%Pixels)', $
timerange=[t0,t1],xstyle=9, psym=2, color=5, charsize=2.5, $
ymargin=[4,2], xmargin=[8,0], title='EIT CCD Affected Area'
endif
!p.multi=psave
if www then begin
root='eit_protons_' + time2file(t0)
gdir=concat_dir('$path_http','protons')
if not file_exist(gdir) then gdir=get_logenv('path_http')
gname=concat_dir(gdir,root+'.gif') ; summary plot
hdoc=concat_dir(gdir,root+'.html') ; top http
zbuff2file,gname ; Z->gif
movie_html=''
cfile=root+'_cadence'
sfile=root+'_peak_spectra.gif'
sname=concat_dir(gdir,sfile)
pfile=root+'_peak_events.gif'
pname=concat_dir(gdir,pfile)
genxfile=concat_dir(gdir,root+'.genx')
; --- make peak event image ------
wdef,xx,im=pdata(*,*,epeak),/zbuffer
ct2rgb,3,gamma=.8
tvscl,safe_log10(pdata(*,*,epeak)<3000,/byte)
zbuff2file,pname
ptf=str_replace(pname,'.gif','_thumb.gif')
pthumb=mkthumb(ingif=pname,outfile=ptf,nx=200)
; --- plot Peak spectra ---------
wdef,xx,512,512,/zbuffer
linecolors
plot, peak_spectra, psym=10, background=11,color=7,$
title='Energy Spectra for Image @' + anytim(pindex(epeak),/vms,/trunc), $
ytitle='Frequency (#Pixels)', xtitle='CCD DN'
zbuff2file,sname
stf=str_replace(sname,'.gif','_thumb.gif')
sthumb=mkthumb(ingif=sname,outfile=stf,nx=200)
thumbhtml=thumbnail_table_html([pname,sname],[pthumb,sthumb])
;make a cadence plot to append to movie -----------
cname=concat_dir(gdir,cfile+'.gif')
wdef,xx,/zbuffer,data_chk(pdata,/nx),data_chk(pdata,/nx)/3
loadct,3
!p.multi=0
utplot, index, cstats, ytitle='EIT CCD Affected Area (%Pixels)', $
timerange=[t0,t1],xstyle=1, psym=2, color=0, charsize=1, $
ymargin=[4,2], xmargin=[8,0], background=!p.color-10, $
mtitle='EIT Affected Area'
evt_grid,index,ticklen=.02,tickpos=.95,color=5,thick=2,linestyle=0
evt_grid,index(epeakd),label='EIT Peak Rise (half max)'+anytim(index(epeakd),/vms,/trunc),$
align=0,color=50,thick=2, tickpos=.99, labpos=.88,labcolor=20,/arrow
evt_grid,index(epeakd),color=20
zbuff2file,cname
eit_fill_cube, pdata
if movie then begin
if n_elements(mwindow) ne 2 then mwindow=[-1,3]
time_window,[anytim(gxd(gpeakd),/vms), $
anytim(pindex(epeak),/vms)],mt0, mt1, hour=mwindow
mss=sel_timrange(pindex,mt0,mt1)
pdata=cube_edit(pdata,/data,ss=mss)
pindex=pindex(mss)
box_message,'Cleaning data via <ssw_unspike_cube>'
cleandata=confac(temporary(ssw_unspike_cube(pindex,pdata)),.5)
pdata=confac(temporary(pdata),.5)
mdata=[pdata,cleandata]
mreadfits_fixup,pindex,mdata
sdata=safe_log10(pindex,mdata,/byte)
event_movie2,pindex, data=sdata, moviedata ,/reverse
special_movie,pindex,moviedata,/inc_times,html=movie_html,$
movie_dir=gdir,thumbsize=128, /no_htmldoc, movie_name=root,$
table=3
endif
html_doc,hdoc,/header
file_append,hdoc,'<h3> EIT Proton Event Summary</h3>'
file_append,hdoc,['Created by <em>eit_proton_summary.pro</em> ', $
' at ' + systime(),'<p>']
header=['Start Time', 'Stop Time',$
'GOES XRay Peak Time',$
'Maximum EIT Affected Area',$
'dT(GOES XR Peak to EIT Peak)', $
'dt(GOES XR Rise to EIT Rise)', $
'Median Threshold','CCD Sub Field [x0,y0,nx,ny]']
tdata=[anytim([time0,time1],/ecs,/trunc),$
anytim(gxd(gpeak),/ecs,/trunc), $
string(max(cstats),format='(f5.1)') +'%', $
string(dgp2ep,format='(f5.1)') +' (hr)',$
string(dgr2er,format='(f5.1)') +' (hr)',$
string(threshold,format='(F5.1)'), $
'[ ' + arr2str(subfov,/trim) + ' ]']
file_append,hdoc,strtab2html([[header],[tdata]],/row0)
gtab=transpose(['GOES X-Ray / PROTONS / EIT Affected Area', $
'<img src="' + root + '.gif">' ])
file_append,hdoc,strtab2html(gtab,/row0)
file_append,hdoc,['<p>',$
'<h3> Image and Spectra at time of EIT event peak: '+ $
anytim(index(epeak),/ecs) + '</h3>',$
thumbhtml]
file_append,hdoc,['<p>',$
'<h3>Event Movie - Original and "cleaned" version</h3>',$
movie_html]
html_doc,hdoc,/trailer
savegen,file=genxfile,root,time0,time1,threshold,subfov,index,pindex,cstats
endif
return
end
Last revised: - Wed May 9 21:45:20 2007- F. Auchère