FUNCTION sdo_eis_overlay_movie, sdolist, eis_obj, individ=individ, $ scale_factor=scale_factor, tgap=tgap, $ tstart=tstart, tend=tend, $ datmax=datmax, goes_tstart=goes_tstart, $ goes_tend=goes_tend, write_png=write_png, $ png_start=png_start, png_dir=png_dir, $ count=count, eis_xoff=eis_xoff, $ eis_yoff=eis_yoff, contlist=contlist, $ light_curve=light_curve, log=log, $ tref=tref, mov_y=mov_y, pddx=pddx, $ write_jpeg=write_jpeg, jpeg_dir=jpeg_dir, $ remove_saturated=remove_saturated ;+ ; NAME ; ; SDO_EIS_OVERLAY_MOVIE ; ; EXPLANATION ; ; This routine takes a list of SDO FITS files, reads them, and ; creates a 3D image array showing the SDO images. By also ; specifying an EIS data object from the same time period as the ; SDO data, the position of the EIS slit is over-plotted on the ; SDO images. ; ; Note that, due to the different time formats used by AIA and ; HMI, this routine converts all times to a common format using ; the ANYTIM suite of routines. The times displayed in the movie ; are UTC format (not TAI). In addition, all times represent the ; midpoint of the exposures, and not the start time. ; ; INPUTS ; ; SDOLIST A list of SDO filenames. These will be read and ; plotted to make a 3D movie array. Note that if EIS_OBJ ; is also specified, then only those elements of SDOLIST ; that bracket the EIS observation will be used. ; ; OPTIONAL INPUTS ; ; EIS_OBJ An EIS data object. Time and pointing information from ; the object will be used to display the position of the ; EIS slit on the SDO images. ; ; TGAP By default, the routine begins the SDO movie 30s ; before the EIS raster starts, and finishes 30s after ; the EIS raster finishes. This time gap can be adjusted ; by setting TGAP. This can be useful if you intend to ; concatenate several movies. ; ; MOV_Y The movie will have image frames that are MOV_Y pixels ; high. The default value is 370. Note that the X-size ; of the frames is computed by the routine based on the ; aspect ratio of the AIA images. ; ; DATMAX Normally, the image data are automatically scaled in ; intensity (see SCALE_FACTOR and INDIVID keywords), but ; by specifying DATMAX the value for the maximum ; intensity value can be specified. This can be useful ; for concatenating several movie arrays. If DATMAX is ; undefined in the call, then DATMAX will contain the ; maximum intensity value when the routine exits. ; ; PNG_START If /WRITE_PNG is set, then the output images are ; labelled image_0000.png, image_0001.png, etc. By ; specifying PNG_START, the indexing can be made to ; start from the value PNG_START instead. This is useful ; when concatenating several movie arrays. ; ; EIS_XOFF This performs a X-shift of the EIS slit location by ; EIS_XOFF arcseconds. ; ; EIS_YOFF This performs a Y-shift of the EIS slit location by ; EIS_XOFF arcseconds. ; ; TREF A time in a standard format accepted by the ANYTIM ; routines. See discussion in PROGRAMMING NOTES section ; below for how to use this. ; ; SCALE_FACTOR This specifies how the image intensity is ; scaled. The minimum intensity is always 0. If the ; maximum intensity is MAX, then the routine will ; plot intensities between 0 and ; MAX*SCALE_FACTOR. The default value is 0.5. If ; /LOG is set, then SCALE_FACTOR is ignored. ; ; PDDX The output plot consists of two panels, and the ; routine automatically sizes the panels for ; display. Sometimes however, the sizing may result in ; some of the text being outside of the displayed ; window. The parameter PDDX allows the arrangement of ; the windows to be tweaked. It sets the spacing between ; between the two panels as a fraction of the plot ; width. The default value is 0.10. Setting values of ; between 0.10 and 0.20 should help with the text ; display. ; ; TSTART This allows the user to manually specify the start ; time of the image sequence. Should be given a time ; format accepted by the anytim2tai.pro routine. ; ; TEND This allows the user to manually specify the end ; time of the image sequence. Should be given a time ; format accepted by the anytim2tai.pro routine. ; ; GOES_TSTART This allows the start time of the GOES plot to be ; manually adjusted. Should be in a format accepted by ; the anytim2tai.pro routine. ; ; GOES_TEND This allows the start time of the GOES plot to be ; manually adjusted. Should be in a format accepted by ; the anytim2tai.pro routine. ; ; KEYWORDS ; ; INDIVID By default the intensities in the images are scaled ; relative to the complete time sequence of images. If ; /INDIVID is set then each image is scaled separately. ; ; WRITE_PNG If set, then the individual image frames are written ; as png files in the current working directory (or ; PNG_DIR). The images are labeled image_0000.png, ; image_0001.png, etc, although see the PNG_START ; keyword. ; ; WRITE_JPEG If set, then the individual image frames are written ; as jpeg files in the current working directory (or ; JPEG_DIR). The images are labeled image_0000.jpg, ; image_0001.jpg, etc, although see the PNG_START ; keyword. ; ; JPEG_DIR Directory to which jpeg files (see WRITE_JPEG) will be ; sent. If not set, then will be written to the current ; working directory. ; ; LOG The images are plotted on log intensity scale. Note ; that scale_factor is ignored in this case. ; ; REMOVE_SATURATED This attempts to remove image frames that are ; badly saturated. The routine counts the number ; of saturated pixels (values >16383), and a ; frame is deemed to be saturated if there are ; more than 300 such pixels. ; ; OPTIONAL OUTPUTS ; ; LIGHT_CURVE An IDL structure that contains a light curve for ; the image sequence. The tags are: ; ; .time Time since 1st exposure (seconds). ; .intensity Average intensity of image. ; .t_init Time in CCSDS format of first exposure. ; ; The light curve can be plotted with ; ; IDL> utplot,light_curve.time,light_curve.intensity,light_curve.t_init ; ; OUTPUT ; ; A 3D byte array of size (NX,NY,N) where NX is computed by the ; routine based on the aspect ratio of the AIA images, NY is the ; same as MOV_Y (see above), and N is the number of AIA image ; frames. For playing the movie it is recommended to use: ; ; IDL> xmovie, movie, rel=0.7, resiz=2 ; ; Trying playing with the REL and RESIZ keywords to adjust how big ; the movie appears on your screen. ; ; PROGRAMMING NOTES ; ; TREF keyword ; ------------ ; In order to produce a movie that doesn't show little ; pixel jumps during the sequence, the maps are plotted with an ; x-range and y-range that are different from those implicit in ; the map structure. A reference image is identified from which ; the x-range, y-range, xcen and ycen are taken. For all other ; images, the reference xcen and ycen are differentially rotated to this ; reference image, and then a new x-range and y-range are computed ; with the same widths as the reference image but with a center ; defined by the new xcen and ycen. The reference image is assumed ; to be the first image in the sequence, but if the time TREF is ; specified then the nearest image to TREF will be used. For long, ; multi-day sequences it is best to set TREF to be near the middle ; of the sequence. ; ; Time range ; ---------- ; The period of time covered by the output movie will vary as ; follows: ; (1) if tstart and/or tend are specified then these ; over-ride everything else. ; (2) if an EIS object is specified then the movie will span the ; period TGAP seconds either side of the EIS raster. ; (3) otherwise the movie will contain all frames of the input ; SDOLIST sequence. ; ; CALLS ; ; READ_SDO, EIS_GET_EXPOSURE_INFO ; ; HISTORY ; ; Ver.1, 21-Apr-2012, Peter Young ; Ver.2, 26-Apr-2012, Peter Young ; now restricts movie to images that bracket the EIS ; observation. ; Ver.3, 18-May-2012, Peter Young ; modified to work with magnetograms; added /noshell and ; /uncomp_delete keywords to read_sdo. ; Ver.4, 24-May-2012, Peter Young ; added TREF= and LIGHT_CURVE= keywords; modified the format of ; the EIS time ; Ver.5, 25-May-2012, Peter Young ; updated header and now correctly handles the different EIS ; slit widths. ; Ver.6, 30-May-2012, Peter Young ; fixed bug with last update ; Ver.7, 1-Jun-2012, Peter Young ; implemented MOV_Y= optional input and added PDDX= optional ; input. ; Ver.8, 8-Jun-2012, Peter Young ; removed time-stamp from GOES plot. ; Ver.9, 4-Aug-2012, Peter Young ; now only plots the low energy GOES channel, and ysty=3 is ; now used for the GOES plot. ; Ver.10, 7-Aug-2012, Peter Young ; added TSTART and TEND keywords; modified how time sequence of ; images is created. ; Ver.11, 26-Sep-2012, Peter Young ; added /REMOVE_SATURATED. ; Ver.12, 4-Jan-2012, Peter Young ; no longer normalize exposures to set min value to zero. ;- IF n_elements(scale_factor) EQ 0 THEN scale_factor=0.50 IF n_elements(tgap) EQ 0 THEN tgap=30.0 IF n_params() LT 1 THEN BEGIN print,'' print,'Use: IDL> movie=sdo_eis_overlay_movie( sdofilelist [, eis_obj, /individ, scale_factor=,' print,' tgap=, datmax=, goes_tstart=, goes_tend=, /write_png, png_start=,' print,' png_dir=, count=, eis_xoff=, eis_yoff=, mov_y=, tstart=, tend=] )' print,'Default values:' print,' scale_factor = '+trim(string(format='(f8.3)',scale_factor)) print,' tgap = '+trim(string(format='(f8.1)',tgap))+' s' print,'' return,-1 ENDIF IF keyword_set(log) THEN scale_factor=1.0 ; ; Check if EIS object has been specified ; IF datatype(eis_obj) EQ 'OBJ' THEN BEGIN ti1=eis_obj->getti_1() ti2=eis_obj->getti_2() timid=(ti2-ti1)/2. + ti1 eis_tai=eis_obj->ti2tai(timid) ; eis_start_tai=anytim2tai(eis_obj->getinfo('DATE_OBS')) eis_end_tai=anytim2tai(eis_obj->getinfo('DATE_END')) ; xpos=eis_obj->getxpos() ypos=eis_obj->getypos(195.12) ycen=eis_obj->getinfo('ycen') ; ycen=eis_obj->getycen(195.12) fovy=eis_obj->getfovy() ; slit_ind=eis_obj->getinfo('SLIT_IND') CASE slit_ind OF 0: slit_wid=1.0 1: slit_wid=266.0 2: slit_wid=2.0 3: slit_wid=40.0 ENDCASE ; eis_get_exposure_info, eis_obj, nexp=nexp, exp_times=exposures, tot_exposures=tot_exposures, $ npos=npos IF nexp EQ 1 THEN BEGIN eis_exp_time=exposures[0] ENDIF ELSE BEGIN print,"%SDO_EIS_OVERLAY_MOVIE: I haven't coded this case yet !!" return,-1 ENDELSE ENDIF IF n_elements(png_start) EQ 0 THEN png_start=0 IF n_elements(png_dir) EQ 0 THEN png_dir='' IF n_elements(eis_xoff) EQ 0 THEN eis_xoff=0.0 IF n_elements(eis_yoff) EQ 0 THEN eis_yoff=0.0 n=n_elements(sdolist) ; ; The following works out the time interval for which the movie is ; required, which can vary depending on whether an EIS object is ; specified, or the TSTART and TEND keywords are ; given - see 'Programming Notes' section in header for more details. ; ; Once the time interval is computed, the SDO data are read. ; read_sdo,sdolist,index,/uncomp_delete,/noshell t_obs_tai=anytim2tai(index.t_obs) ; IF n_elements(tstart) EQ 0 THEN BEGIN IF n_elements(eis_start_tai) EQ 0 THEN BEGIN tstart_tai=t_obs_tai[0] ENDIF ELSE BEGIN tstart_tai=eis_start_tai-tgap ENDELSE ENDIF ELSE BEGIN tstart_tai=anytim2tai(tstart) ENDELSE ; IF n_elements(tend) EQ 0 THEN BEGIN IF n_elements(eis_end_tai) EQ 0 THEN BEGIN tend_tai=t_obs_tai[-1] ENDIF ELSE BEGIN tend_tai=eis_end_tai+tgap ENDELSE ENDIF ELSE BEGIN tend_tai=anytim2tai(tend) ENDELSE ; k=where(t_obs_tai GE tstart_tai AND t_obs_tai LE tend_tai,n) ; IF n EQ 0 THEN BEGIN print,'%SDO_EIS_OVERLAY_MOVIE: there are no SDO files in the specified time interval. ' print,' Please check your inputs. Returning...' return,-1 ENDIF ; read_sdo,sdolist[k],index,data,/uncomp_delete,/noshell IF keyword_set(remove_saturated) THEN BEGIN ind=-1 FOR i=0,n-1 DO BEGIN img=data[*,*,i] k=where(img GE 16383,nk) IF nk LT 300 THEN ind=[ind,i] ENDFOR ind=ind[1:*] index=index[ind] data=data[*,*,ind] ; IF n_elements(ind) NE n THEN BEGIN print,'**Removed '+trim(n-n_elements(ind))+' saturated images.**' n=n_elements(index) ENDIF ENDIF ;; IF n_elements(goes_tstart) NE 0 OR n_elements(goes_tend) NE 0 THEN BEGIN ;; IF n_elements(goes_tstart) EQ 0 THEN tstart_tai=anytim2tai(index[0].t_obs) $ ;; ELSE tstart_tai=anytim2tai(goes_tstart) ;; IF n_elements(goes_tend) EQ 0 THEN tend_tai=anytim2tai(index[-1].t_obs)$ ;; ELSE tend_tai=anytim2tai(goes_tend) ;; ; ;; t_obs_tai=anytim2tai(index.t_obs) ;; k=where(t_obs_tai GE tstart_tai AND t_obs_tai LE tend_tai,n) ;; IF n EQ 0 THEN BEGIN ;; print,'%SDO_EIS_OVERLAY_MOVIE: there are no SDO files for time period specified by ' ;; print,' goes_tstart and goes_tend. Returning...' ;; return,-1 ;; ENDIF ;; read_sdo,sdolist[k],index,data, /uncomp_delete, /noshell ;; eis_swtch=0 ;; ENDIF ELSE BEGIN ;; read_sdo,sdolist,index,data, /uncomp_delete, /noshell ;; eis_swtch=1 ;; ENDELSE ;; The following restricts the SDO images to only those that occur ;; within +/- TGAP seconds of the EIS raster. ;; IF eis_swtch EQ 1 THEN BEGIN ;; t_obs_tai=anytim2tai(index.t_obs) ;; IF n_elements(eis_start_tai) NE 0 THEN BEGIN ;; k=where(t_obs_tai GE eis_start_tai-tgap AND t_obs_tai LE eis_end_tai+tgap,nk) ;; IF nk GT 0 THEN BEGIN ;; index=index[k] ;; data=data[*,*,k] ;; n=nk ;; ENDIF ELSE BEGIN ;; print,'%SDO_EIS_OVERLAY_MOVIE: there are no SDO images that match the EIS raster. Returning...' ;; return,-1 ;; ENDELSE ;; ENDIF ;; ENDIF ; ; I want to keep AIA images and HMI separate from each other as I need ; special code for magnetograms. ; inst=strlowcase(strmid(index.instrume,0,3)) chck=where(inst EQ 'aia',n_aia) chck=where(inst EQ 'hmi',n_hmi) ; IF n_aia NE n AND n_hmi NE n THEN BEGIN print,'%SDO_EIS_OVERLAY_MOVIE: the file list contains both HMI and AIA images. Please separate them before' print,' calling this routine.' print,' Returning...' return,-1 ENDIF ; ; Flag magnetogram sequences with keyword 'magnetogram'. Bail out if ; the HMI sequence contains other types of data. ; IF n_hmi GT 0 THEN BEGIN k=where(index.content EQ 'MAGNETOGRAM',nk) IF nk NE 0 AND nk NE n THEN BEGIN print,'%SDO_EIS_OVERLAY_MOVIE: the file list contains mixed HMI filter images, including '+trim(nk)+' magnetograms.' print,' Please process the magnetograms separately.' print,' Returning...' data=0 return,-1 ENDIF IF nk NE 0 THEN magnetogram=1 ELSE magnetogram=0 ENDIF count=n siz=size(data) nx=siz[1] ny=siz[2] ; ; These processing steps don't apply to magnetograms. ; IF NOT keyword_set(magnetogram) THEN BEGIN ; ; Since the exposure time can be variable, it's necessary to ; scale the individual frames by the exposure time. ; newdata=fltarr(nx,ny,siz[3]) FOR i=0,n-1 DO newdata[*,*,i]=data[*,*,i]/index[i].exptime data=0 data=temporary(newdata) ; ; Set minimum of image to be zero. ; PRY, 4-Jan-2012, I've commented the following out as some ; rogue exposures can have large negative values. ; ; IF NOT keyword_set(log) THEN data=data-min(data) ENDIF ; ; Create a light curve from the data and store in structure. The ; intensity is just the average of image intensity over the field of ; view. (Note that light curve is created after the correction has ; been made for different exposure times.) ; light_curve={time: dblarr(n), intensity: fltarr(n), t_init: ''} light_curve.time=anytim2tai(index.t_obs) - anytim2tai(index[0].t_obs) light_curve.intensity=average(data,[1,2]) light_curve.t_init=anytim2utc(index[0].t_obs,/ccsds) ; ; Set the maximum intensity value to be used for plotting. ; IF n_elements(datmax) EQ 0 THEN datmax=max(abs(data))*scale_factor ; ; T_OBS is the observation time of the AIA image. ; sdo_tai=anytim2tai(index.t_obs) sdo_ccsds=anytim2utc(index.t_obs,/ccsds) ; ; Create a GOES object and set the start and end times suitable for ; the displayed SDO images. ; IF n_elements(goes_tstart) EQ 0 THEN goes_tstart=index[0].t_obs IF n_elements(goes_tend) EQ 0 THEN goes_tend=index[n-1].t_obs ; goes_tstart_tai=anytim2tai(goes_tstart) ; g=ogoes() g->set,tstart=goes_tstart, tend=goes_tend,/sdac,mode=1 ; ; These are parameters that specify how the two plots are arranged in ; the movie. See the POSITION= keyword in the later plot commands. ; px0=0.00 px1=0.03 pdx1=0.55 pdx2=(1-px0-px1-pdx1) IF n_elements(pddx) EQ 0 THEN pddx=0.10 py0=0.07 py1=0.15 pdy=(1-py0-py1) pddy=0.05 ; ; I'm setting the movie array to be mov_y pixels tall in the ; Y-direction. I can compute how big the array has to be in the ; X-direction using AIA image array size and the plot parameters. ; IF n_elements(mov_y) EQ 0 THEN mov_y=370 ELSE mov_y=fix(mov_y) mov_x=fix(float(nx)/float(ny)*(pdy-pddy)/(pdx1-pddx)*mov_y) movie=bytarr(mov_x,mov_y,n) ; ; The following code computes how the (xcen,ycen) of the first image ; changes during the sequence. This is used to set the xrange and ; yrange of the plot images (see call to plot_map later). The idea is ; get rid of the little pixel jumps in the movie sequence that occur ; when each image is plotted with plot_map. ; ; The routine can be quite slow for long image sequences. I ; implemented a /QUICK keyword to see if it would speed things up, but ; it didn't seem to so I disabled it. ; t_obs=index.t_obs t_obs_tai=anytim2tai(index.t_obs) ;IF NOT keyword_set(quick) THEN BEGIN IF n_elements(tref) EQ 0 THEN BEGIN iref=0 tref_tai=t_obs_tai[iref] ENDIF ELSE BEGIN tref_tai=anytim2tai(tref) getmin=min(abs(t_obs_tai-tref_tai),iref) tref_tai=t_obs_tai[iref] ENDELSE xc_0=index[iref].xcen yc_0=index[iref].ycen t_obs_0=index[iref].t_obs lonlat_0=xy2lonlat([xc_0,yc_0],t_obs_0) ddays=(anytim2tai(t_obs)-anytim2tai(t_obs_0))/86400. drot=diff_rot(ddays,lonlat_0[1]) lonlat=fltarr(2,n) lonlat[0,*]=lonlat_0[0]+drot lonlat[1,*]=lonlat_0[1] xy_c=lonlat2xy(lonlat,t_obs) ;ENDIF all_min=fltarr(n) FOR i=0,n-1 DO all_min[i]=min(data[*,*,i]) med_min=median(all_min) FOR i=0,n-1 DO BEGIN set_plot,'z' device,set_res=[mov_x,mov_y] !p.multi=[0,2,1] img=data[*,*,i] ; IF keyword_set(magnetogram) THEN BEGIN dmin=-datmax dmax=datmax t_exp='' ENDIF ELSE BEGIN IF keyword_set(individ) THEN BEGIN IF keyword_set(log) THEN BEGIN img=alog10(img > 1) dmax=max(img) img=dmax-img dmin=dmax-dmax*scale_factor ENDIF ELSE BEGIN img=max(img)-img dmax=max(img) dmin=max(img)-max(img)*scale_factor ENDELSE ; ENDIF ELSE BEGIN IF keyword_set(log) THEN BEGIN img=alog10(img > 1) dmax=alog10(datmax) img=dmax-img dmax=alog10(datmax)-alog10(med_min) dmin=0. ENDIF ELSE BEGIN img=imgplot,showclass=1,xsty=1,legend=0, $ position=[px0+pdx1+pddx,py0+pddy,px0+pdx1+pdx2,1-py1], $ /no_timestamp,/low,ysty=3 outplot,(sdo_tai[i]-goes_tstart_tai)*[1,1],[1e-10,1e10] im=tvrd() im=bytscl(temporary(im)) movie[*,*,i]=im IF keyword_set(write_png) THEN BEGIN pngfile='image_'+strpad(trim(i+png_start),4,fill='0')+'.png' IF png_dir NE '' THEN pngfile=concat_dir(png_dir,pngfile) write_png,pngfile,im ENDIF IF keyword_set(write_jpeg) THEN BEGIN jpgfile='image_'+strpad(trim(i+png_start),4,fill='0')+'.jpg' IF jpeg_dir NE '' THEN jpgfile=concat_dir(jpeg_dir,jpgfile) write_jpeg,jpgfile,im ENDIF !p.multi=0 device,/close set_plot,'x' ENDFOR ; ; Clean up objects ; obj_destroy,g return,movie END