;+ function sdc_verify,index,data,verbose=verbose,progress=progress, $ correct_col_blem=correct_col_blem,version=version ;NAME: ; SDC_VERIFY ;PURPOSE: ; Check the SDC database for bad dark frames ;CATEGORY: ;CALLING SEQUENCE: ; quality = sdc_verify(index,data) ;INPUTS: ; index = SDC index ; data = SDC data ;OPTIONAL INPUT PARAMETERS: ;KEYWORD PARAMETERS ; /correct_col_blem = attempt to fix column blemish problems by ; interpolating over the bad areas. WARNING: if this ; keyword is set, the data will be modified if any areas ; of the images appear to have bleed or blemishes. NOte ; that the quality structure will have the values for the ; *uncorrected* data so that data which did not need to ; be corrected will be used in favor of data which did ; need to be corrected. ; /progress = print out progress reports as the checks are underway. ; Like /verbose, but with fewer diagnostics and no ; images. ; /verbose = spit out lots of diagnostics and display the good and bad ; images. ; version = Outputs with the version number ;OUTPUTS: ; quality = index type structure for images ; quality.badness = measure of how bad the image is. Allows ranking. ; 0 means "perfect", higher numbers are worse. ; 0-20 is probably OK; ; 20-100 is marginal; ; >=100 is probably bad; ; >= 1000. do not use the image. ;COMMON BLOCKS: ;SIDE EFFECTS: ;RESTRICTIONS: ;PROCEDURE: ; ;MODIFICATION HISTORY: ; T. Metcalf 1998-Sep-21 ; 1998-Sep-29 TRM Added check for exposure mode and BLS mode ; 1998-Dec-16 TRM Major modification to the tests ; 1999-May-24 TRM Added minimum 8ms exposure to gt_dpe call ; 1999-Jun-15 TRM Added column blemish correction ; 2002-Dec-02 TRM Small change to "ignore vertical lines" so that ; streaks are picked up more readily. ; 2002-Dec-03 TRM Version 1.0. Added version number. ; 2002-Dec-03 TRM Version 1.1. Added TFMS check and percent ; check. Increased pixcut from 0 to 30. ; 2002-Dec-12 TRM Version 1.2. Added check for nx=1024 images ;- version = fix(1.2 * 1000.) qualitystr = {sdc_qcheck, $ gen:{sdc_qcheck_gen,day:0,time:0L,rbm_status:0b}, $ sxt:{sdc_qcheck_sxt,imgparam:0b,explevmode:0b,pfi_ffi:0b,periph:0b, $ corner_cmd:intarr(2),expdur:0}, $ date_obs:'XXXXXXXXXXXXXXXXXX', $ ; Date/time of image reason:'XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX', $ ; Was an error detected? nsaturated: 0, $ ; # saturated pixels nbleed: 0, $ ; # contiguous columns with bleed nbleed_tot: 0, $ ; # columns with bleed (not nec. contiguous) nhot: 0, $ ; # hot pixels nstreak: 0, $ ; # pixels in a "streak" bls: 0, $ ; 1 if BLS mode is on, 0 if OK mode: 0, $ ; 1 if NOT dark mode, 0 if OK experror: 0., $ ; Percentage exposure duration is off tfms:0., $ ; Try to avoid early orbit darks stripimage:0., $ ; Avoid nx=1024 images percent:0., $ ; Combination of percentd and percentover badness: 0. $ ; # degree of "badness" for image } quality = replicate(qualitystr,n_elements(index)) if keyword_set(verbose) then begin max_xsize = max(index.sxt.shape_cmd[0]) max_ysize = max(index.sxt.shape_cmd[1]) endif for image=0L,n_elements(index)-1 do begin xsize = index[image].sxt.shape_cmd[0] ysize = index[image].sxt.shape_cmd[1] d0 = data[0:xsize-1,0:ysize-1,image] & i0 = index[image] if (datatype(i0,2)) NE 8 then message,'Index must be a structure' res=gt_res(i0) ;-------------- Check quality flags ------------------ percent = (255.0 - index[image].sxt.percentd) + index[image].sxt.percentover if index[image].sxt.shape_cmd[0] EQ 1024 AND $ index[image].sxt.shape_cmd[1] LT 512 then stripimage=1 else stripimage=0 ;-------------- Check time in orbit ------------------ tfms = (sxt_uvf_info(index[image],/tfms))[0] if tfms GE 10.0 then tfms = 0.0 $ else if tfms GE 4.5 then tfms = (10.0-tfms) $ else tfms = (10.0-tfms)^2 - 24.75 ; the 24.75 makes the transition at 4.5 smooth ;-------------- Check for exposure problems ---------- ; Exposure difference in percent check_exp = 100.*abs(gt_expdur(i0) - (gt_dpe(i0,/conv)>8.)) / $ (gt_dpe(i0,/conv)>8.) ;-------------- Check for Saturation ----------------- ; Saturation expected for long exposures so only look at shorter exps dpe_hi_cutoffs = [21, 23, 19] ; should match values in sel_dc_image dpe_lo_cutoffs = [13, 13, 13] satarray = sxt_satpix(i0,d0) satarray[*,0:4] = 0 isat = where((2.*(satarray)-shift(satarray,1,0)) eq 1,nsat) ; Saturation only counts as an error for short exposures ; The 2% cutoff below needs to be calibrated if gt_comp(i0) EQ 0 then d1 = sxt_decomp(d0) else d1=fix(d0) if gt_dpe(i0) GE dpe_hi_cutoffs[res] and $ nsat GT 0 and nsat LT n_elements(d0)*0.02 then begin ; If saturatiion it may be OK, so remove it before the dstreek call d1_desat = rm_feature(temporary(d1),isat) ; bilinear interpolation nsat = 0 d1 = temporary(d1_desat) endif ;--------------- Check for Streaks (e.g. cosmic ray hits) -------- dhot = d1 d1 = median(temporary(d1),3) ; Get rid of salt and pepper noise sd1 = stdev(d1) mn1 = mean(d1) sdlimit = 4.0*(1.0+1.0*(gt_expdur(i0)/30208.)) pixcut = 30 vnorm = total(d1,2)/ysize d2 = d1 for i=0,xsize-1 do d2[i,*]=d2[i,*]/vnorm[i] d2 = ((d2-mean(d2))/stdev(d2)) GE sdlimit d2[*,0:pixcut/(2^res)] = 0 d2[0:pixcut/(2^res),*] = 0 d2[*,ysize-1-pixcut/(2^res):*] = 0 d2[xsize-1-pixcut/(2^res):*,*] = 0 ; Ignore vertical lines ; TRM 2002-Dec-02 ;vertical = where( d2 AND $ ; ( (shift(d2,+0,+1) AND $ ; shift(d2,+0,+2) AND $ ; shift(d2,+0,+3) ) $ ; OR $ ; (shift(d2,+0,-1) AND $ ; shift(d2,+0,-2) AND $ ; shift(d2,+0,-3))),nvertical) vertical = where(d2 AND ((shift(d2,+0,+1) AND shift(d2,+0,-1)) OR $ (shift(d2,+0,+1) AND shift(d2,+0,+2)) OR $ (shift(d2,+0,-1) AND shift(d2,+0,-2))),nvertical) ;vertical = where( $ ; d2 AND (shift(d2,+0,+1) OR shift(d2,+0,-1)) $ ;AND $ ; ;(shift(d2,+0,+2) OR shift(d2,+0,-2)) AND $ ; ;(shift(d2,+0,+3) OR shift(d2,+0,-3)) $ ;AND $ ; ;(shift(d2,+0,+4) OR shift(d2,+0,-4)) AND $ ; ;(shift(d2,+0,+5) OR shift(d2,+0,-5)) AND $ ; ;(shift(d2,+0,+6) OR shift(d2,+0,-6)) $ ; ,nvertical) if nvertical GT 0 then d2[vertical] = 0 ; Look for a roughly linear structure 3 pixels or more long streek = where( $ (d2 AND shift(d2,+0,+1) AND $ (shift(d2,+0,+2) OR $ shift(d2,+1,+2) OR $ shift(d2,-1,+2))) OR $ (d2 AND shift(d2,+0,-1) AND $ (shift(d2,+0,-2) OR $ shift(d2,+1,-2) OR $ shift(d2,-1,-2))) OR $ (d2 AND shift(d2,+1,+0) AND $ (shift(d2,+2,+0) OR $ shift(d2,+2,+1) OR $ shift(d2,+2,-1))) OR $ (d2 AND shift(d2,+1,+1) AND $ (shift(d2,+1,+2) OR $ shift(d2,+2,+1) OR $ shift(d2,+2,+2))) OR $ (d2 AND shift(d2,+1,-1) AND $ (shift(d2,+2,-1) OR $ shift(d2,+1,-2) OR $ shift(d2,+1,-2))) OR $ (d2 AND shift(d2,-1,+0) AND $ (shift(d2,-2,+0) OR $ shift(d2,-2,+1) OR $ shift(d2,-2,-1))) OR $ (d2 AND shift(d2,-1,+1) AND $ (shift(d2,-1,+2) OR $ shift(d2,-2,+1) OR $ shift(d2,-2,+1))) OR $ (d2 AND shift(d2,-1,-1) AND $ (shift(d2,-1,-2) OR $ shift(d2,-2,-1) OR $ shift(d2,-2,-2))), nbadpixels) d3 = d1 ;--------------- Check for hot pixels ------------------------ ;dhot = sxt_hotpix_fudge(i0,dhot,status=hpstat) rmeandhot = smooth(dhot,25/(2^res),/edge_truncate) ; running mean sdhot = stdev(dhot) dd = (dhot[25/(2^res):xsize-25/(2^res),25/(2^res):ysize-25/(2^res)] - $ rmeandhot[25/(2^res):xsize-25/(2^res),25/(2^res):ysize-25/(2^res)])/sdhot sdlimit = 50.*(1.0+1.0*(gt_expdur(i0)/30208.)) checkhot = where(dd GT sdlimit,nhot) nhot = nhot; + (hpstat GE 1) ; --------------- Check for bleed and blemishes ---------------- ; Get rid of salt and pepper noise d5 = median(fix(d3),5) ; This is very slow.. Was 7 changed to 5 md5 = mean(d5) sd5 = stdev(d5) ; Cut out the edges d5 = d5[*,25/(2^res):ysize-26/(2^res)] ; Compress along the y-dimesion, to give a 1-d vector d5 = total(d5,2) rmeand5 = smooth(d5,25/(2^res),/edge_truncate) ; running mean xstart = 25/(2^res) d5 = d5(xstart:xsize-26/(2^res)) rmeand5 = rmeand5(xstart:xsize-26/(2^res)) ;d5 = abs(d5-rmeand5)/stdev(d5) d5 = (d5-rmeand5)/stdev(d5) sdlimit = 0.8*(1.0+1.0*(gt_expdur(i0)/30208.)) d6 = d5 GT sdlimit d6[0] = 0 & d6[n_elements(d6)-1]=0 check2 = where(d6,nbleed_tot) ; # number columns check2 = where(d6 and $ shift(d6,1) and $ shift(d6,2),nbleed) ; # 3 contiguous columns if keyword_set(correct_col_blem) and nbleed GT 1 then begin ; Attempt to interpolate over the bad areas c = intarr(xsize,ysize) & c[*] = 0 colss = check2+xstart ; Since there must be 3 contiguous columns but the above where ; statement only picks up the last, add the two columns before ; the tagged column. This would have to ; change if the number of columns required for ; a "hit" were to change. colss = [min(colss)-2,min(colss)-1,colss]>0 c[colss,*] = 1 ; Find out where in the column the problem is ss = where(d3 GT (md5+sd5) AND c NE 0,nss) ; interpolate the problem away if nss GT 0 then begin message,/info,'Interpolating SDC column blemish: '+fmt_tim(i0) dcorrect = rm_feature(data[0:xsize-1,0:ysize-1,image],ss) ;erase ;tvscl,data[0:xsize-1,0:ysize-1,image] & tvscl,dcorrect,0,ysize ;xyouts,xsize+20,ysize/2,'Uncorrected',/device,charsize=2 ;xyouts,xsize+20,3*ysize/2,'Corrected',/device,charsize=2 data[0:xsize-1,0:ysize-1,image] = dcorrect endif endif ; Save string describing errors bad_reason = '' if nbadpixels GT 1 then bad_reason = bad_reason + 'Streak ' if nsat GT 0 then bad_reason = bad_reason + 'Saturation ' if gt_expmode(i0) NE 1 then bad_reason = bad_reason + 'Expmode ' if mask(i0.sxt.pfi_ffi,3,3) NE 0 then $ bad_reason = bad_reason + 'BLS ' if nhot GT 1 then bad_reason = bad_reason + 'Hot ' if nbleed GT 1 then bad_reason = bad_reason + 'Bleed ' if check_exp GT 10 then bad_reason = bad_reason + 'Exposure ' if tfms GT 5.5 then bad_reason = bad_reason + 'TFMS ' if stripimage then bad_reason = bad_reason + 'Strip_Image ' if percent GT 10. then bad_reason = bad_reason + 'Bad_Data ' quality[image].date_obs = fmt_tim(i0) quality[image].reason = bad_reason quality[image].nsaturated = nsat quality[image].nbleed = nbleed quality[image].nbleed_tot = nbleed_tot quality[image].nstreak = nbadpixels quality[image].nhot = nhot quality[image].bls = mask(i0.sxt.pfi_ffi,3,3) NE 0 quality[image].mode = gt_expmode(i0) NE 1 quality[image].experror = check_exp quality[image].tfms = tfms quality[image].stripimage = stripimage quality[image].percent = percent ; The badness is a measure of how bad each image is. The weights ; are supposed to make each bad "event" have a score of about 100. ; BLS and non-dark mode should eliminate the image, so give them 1000. quality[image].badness = 8.*quality[image].nsaturated + $ 20.*quality[image].nbleed + $ 10.*quality[image].nstreak + $ 20.*quality[image].nhot + $ 1000.*quality[image].bls + $ 1000.*quality[image].mode + $ 1.*quality[image].experror^2 + $ 10.*quality[image].tfms + $ 9999.*quality[image].stripimage + $ 1000.*quality[image].percent quality[image].gen.day = i0.gen.day quality[image].gen.time = i0.gen.time quality[image].gen.rbm_status = i0.gen.rbm_status quality[image].sxt.corner_cmd = i0.sxt.corner_cmd quality[image].sxt.periph = i0.sxt.periph quality[image].sxt.imgparam = i0.sxt.imgparam quality[image].sxt.explevmode = i0.sxt.explevmode quality[image].sxt.pfi_ffi = i0.sxt.pfi_ffi quality[image].sxt.expdur = i0.sxt.expdur if nbadpixels GT 1 OR $ ; Streak (like a cosmic ray hit) nsat GT 0 OR $ ; Saturation gt_expmode(i0) NE 1 OR $ ; In Dark mode mask(i0.sxt.pfi_ffi,3,3) NE 0 OR $ ; No BLS mode nhot GT 1 OR $ check_exp GT 10. OR $ nbleed GT 1 then begin ; Bleed ; found error ; If verbose mode save data for display if keyword_set(verbose) then begin bigdata = intarr(max_xsize,max_ysize) bigdata[0,0] = d0 if n_elements(bad_data) EQ 0 then begin iibad = image bad_index = index[image] bad_data = bigdata endif else begin iibad = [iibad,image] bad_index = [bad_index,index[image]] bad_data = [[[bad_data]],[[bigdata]]] endelse endif if keyword_set(progress) or keyword_set(verbose) then $ print,get_info2(i0,/noninteractive),' ',bad_reason endif else begin if keyword_set(verbose) then begin bigdata = intarr(max_xsize,max_ysize) bigdata[0,0] = d0 if n_elements(good_data) EQ 0 then begin good_index = index[image] good_data = bigdata endif else begin good_index = [good_index,index[image]] good_data = [[[good_data]],[[bigdata]]] endelse endif endelse endfor if n_elements(good_data) GT 0 and keyword_set(verbose) then begin print,'Displaying GOOD images ... ',n_elements(good_index) xstepper,good_data,get_info(good_index) endif if n_elements(bad_data) GT 0 and keyword_set(verbose) then begin print,'Displaying BAD images ... ',n_elements(bad_index) info = get_info(bad_index) for i=0,n_elements(info)-1 do $ info[i]=info[i]+' '+quality[iibad[i]].reason xstepper,bad_data,info endif return,quality end