pro IRIS_CHECK_WOBBLE, l2dir, opath, $ disps, subdat, corrdat, $ nocorr = nocorr, noplot = noplot, pixels = pixels, pstop = pstop, bitrate=bitrate, trim=trim ;+ ; ; Examine the wobble in a set of level 2 data, in a number of ways: ; 1) Make low-cadence (~3 minute) SJI movies, which tend to make the ; residual wobble much easier to identify; ; 2) Measure frame-to-frame cross-correlation of selected SJI channels ; and plot the results ; ; INPUTS: ; l2dir - Set to the path to the level2 data to look at ; (e.g. /irisa/data/level2_nrt/2015/05/31/20150531_075615_3820356897/) ; opath - Set to the path where the output files (movies, maybe a PNG ; of a plot) should be written ; ; OUTPUTS: ; disps ; ; KEYWORDS: ; /nocorr - If set, then the correlation is not performed ; /noplot - If set, then the output plots are not made (currently doesn't matter) ; /pstop - Stop in the routine for debugging ; /trim - If set, movie is trimmed to include only area that has data in all frames ; pixels= The size of the output movie, in un-summed IRIS pixels (i.e. set ; pixels=400 to play a movie of a 400x400 region of the image) ;- if n_elements(bitrate) eq 0 then bitrate=2e6 wobble_cadence = 180. f = FILE_SEARCH(CONCAT_DIR(l2dir, 'iris_l2_*_SJI_*.fits')) numf = N_ELEMENTS(f) if f[0] eq '' then begin PRINT, 'IRIS_CHECK_WOBBLE: No SJI files found!' RETURN endif s2832 = STRPOS(f, 'SJI_2832') if MAX(s2832) lt 0 then PRINT, $ 'IRIS_CHECK_WOBBLE: Try to use something with 2832 SJI' ; Read in the first SJI channel and check it out READ_IRIS_L2, f[0], hdr, dat exptime = hdr[0].exptime nframe = N_ELEMENTS(hdr) sumspat = hdr[0].sumspat sumsptrl = hdr[0].sumsptrl duration = MAX(hdr.time) pztxs = LIMITS(hdr.pztx) ; Check that the OBS is reasonable for checking wobble if duration lt 7200d then PRINT, $ 'IRIS_CHECK_WOBBLE: Try to use something with at least 2 hour duration' if exptime le 2.0 then PRINT, $ 'IRIS_CHECK_WOBBLE: Try to use something with exposure time longer than 2 seconds' if (pztxs[1] - pztxs[0]) gt 50 then PRINT, $ 'IRIS_CHECK_WOBBLE: Try to use something that doesnt raster so far' ; Set up color tables/backing store width = 700 height = 700 fps = 20 TVLCT, rr, gg, bb, /get oldpmulti = !p.multi WINDOW, 12, retain = 2, xsize = width, ysize = height !p.multi = 0 ; Make a wobble movie for i = 0, numf - 1 do begin if i gt 0 then READ_IRIS_L2, f[i], hdr, dat nframe = N_ELEMENTS(hdr) sumspat = hdr[0].sumspat sumsptrl = hdr[0].sumsptrl tai0 = ANYTIM2TAI(hdr[0].date_obs) thispath = hdr[0].tdesc1 filebase = CONCAT_DIR(opath, thispath + '_') dustdat = IRIS_DUSTBUSTER(hdr, dat) IRIS_LCT, hdr ; Downsample in time to accentuate the wobble cadence = MEDIAN(DERIV(hdr.time)) cadence_sample = ROUND(wobble_cadence / cadence) > 1 nsub = nframe / cadence_sample index = INDGEN(nsub) * cadence_sample if keyword_set(trim) then begin ; Trim down to only that part of the movie that contains data in all frames dmin = MIN(dustdat, dim = 3) dmask = dmin gt -200 dmx = TOTAL(dmask, 1) dmy = TOTAL(dmask, 2) subx = LIMITS(WHERE(dmx gt (MAX(dmx) * 0.8))) suby = LIMITS(WHERE(dmy gt (MAX(dmy) * 0.8))) endif else begin subx = [0,n_elements(dat[0,*,0])-1] suby = [0,n_elements(dat[*,0,0])-1] endelse subdat = IRIS_INTSCALE(dustdat[suby[0]:suby[1], subx[0]:subx[1],*], hdr) sdsize = SIZE(subdat) if KEYWORD_SET(pixels) then begin ; Find the brightest part of the image and subframe it profx = SMOOTH(TOTAL(TOTAL(subdat, 3), 2), 20) profy = SMOOTH(TOTAL(TOTAL(subdat, 3), 1), 20) centx = WHERE(profx eq MAX(profx)) centy = WHERE(profy eq MAX(profy)) xmin = 0 > (centx - (pixels / 2. / sumsptrl)) ymin = 0 > (centx - (pixels / 2. / sumspat)) xmax = (xmin + (pixels / sumsptrl) ) < (sdsize[1]-1) ymax = (ymin + (pixels / sumspat) ) < (sdsize[2]-1) endif else begin ; Just use the whole image xmin = 0 ymin = 0 xmax = sdsize[1]-1 ymax = sdsize[2]-1 endelse vfile = CONCAT_DIR(opath, thispath + '_wobble.mp4') CD, opath, current = old_dir video = IDLffVideoWrite(FILE_BASENAME(vfile)) vidStream = video.AddVideoStream(width, height, fps, bit_rate=bitrate) for j = 0, nsub - 1 do begin PLOT_IMAGE, subdat[xmin:xmax,ymin:ymax,index[j]], chars = 1.5, scale = [sumsptrl, sumspat], $ xtitle = 'Un-summed pixels [X]', ytitle = 'Un-summed pixels [Y]', $ title = thispath + ' ' + TAI2UTC(tai0 + hdr[index[j]].time, /ccsds) thisframe = TVRD(/true) dummy = video.Put(vidStream, thisframe) endfor video.Cleanup CD, old_dir endfor ; Make correlation-tracking plots (using the last SJI channel, ; which is typically the 2832 if not KEYWORD_SET(nocorr) then begin corrdat = subdat ; De-emphasize the slit so it doesn't confuse the correlation ; corrdat[WHERE(subdat lt (MEDIAN(subdat)*0.6))] = MEDIAN(subdat) disps = TR_GET_DISP(corrdat, /shift) ; Make correlation-corrected movie WINDOW, 14, xsize = width*2, ysize = height, retain = 2 !p.multi = [0, 2, 1] vfile = CONCAT_DIR(opath, thispath + '_correl.mp4') CD, opath, current = old_dir video = IDLffVideoWrite(FILE_BASENAME(vfile)) vidStream = video.AddVideoStream(width*2, height, fps, bit_rate=bitrate) for j = 0, nsub - 1 do begin PLOT_IMAGE, subdat[xmin:xmax,ymin:ymax,index[j]], chars = 1.5, scale = [sumsptrl, sumspat], $ xtitle = 'Un-summed pixels [X]', ytitle = 'Un-summed pixels [Y]', $ title = 'Original ' + TAI2UTC(tai0 + hdr[index[j]].time, /ccsds), $ min = 0, max = 255 PLOT_IMAGE, corrdat[xmin:xmax,ymin:ymax,index[j]], chars = 1.5, scale = [sumsptrl, sumspat], $ xtitle = 'Un-summed pixels [X]', ytitle = 'Un-summed pixels [Y]', $ title = 'Corrected ' + TAI2UTC(tai0 + hdr[index[j]].time, /ccsds), $ min = 0, max = 255 thisframe = TVRD(/true) dummy = video.Put(vidStream, thisframe) endfor video.Cleanup CD, old_dir ; Make plot of correlation !p.multi = 0 LOADCT, 0 WINDOW, 12, xsize = width, ysize = height, retain = 2 PLOT, disps[0,*] * sumsptrl, disps[1,*] * sumspat, chars = 1.5, psym = -4, $ xtitle = 'Un-summed Pixels X', ytitle = 'Un-summed Pixels Y', /iso, $ title = 'Frame-to-frame correlation' thisplot = TVRD(/true) WRITE_PNG, CONCAT_DIR(opath, 'correl.png'), thisplot endif ; Make plots of the expected wobble correction if not KEYWORD_SET(noplot) then begin ; IRIS_COMPARE_WOBBLE endif if KEYWORD_SET(pstop) then STOP TVLCT, rr, gg, bb !p.multi = oldpmulti end