pro get_dc_image, struct_in, index, data, imap_out, match_qual=match_qual, qstop=qstop, qdebug=qdebug, nodata=nodata, $ times_input=times_input, res_input=res_input, dpe_input=dpe_input, temp_input=temp_input, $ interpolate=interpolate, dc_interpolate=dc_interpolate, qerr=qerr, q2debug=q2debug, float=float, $ unc_data=unc_data, $ orbit_correct=orbit_correct ; ; ;+ ;NAME: ; get_dc_image ;PURPOSE: ; Given an index or roadmap structure, match and ; extract the closest dark current image for that exposure ;CALLING SEQUENCE: ; get_dc_image, index(3), dc_index, dc_data ; get_dc_image, roadmap, dc_index, dc_data, imap ; ; get_dc_image, xxx, dc_index, dc_data, times=times, res=res, dpe=dpe ; get_dc_image, xxx, dc_index, dc_data, times='1-sep-92 1:00', res=0, dpe=15 ; ; unc_data = 1 ;YOU HAVE TO Pass in a defined variable to get uncertainties ; get_dc_image, index(3), dc_index, dc_data, unc_data=unc_data ;INPUT: ; struct_in - An index or roadmap structure. If the input ; is an array, then data will be a matrix and it ; is necessary to use the output variable "imap" ; to match up which dark frame goes with which input ; image. This is done so that unnessary duplication ; of the same image is not done. ;OUTPUT: ; index - The dark current index structure ; data - The dark current image DECOMPRESSED ; imap - A list of indicies telling which image in "data" to ; use for each "struct_in" image. The output is an ; array the same length as "struct_in" with values that ; go from 0 to the number of output images minus one. ;KEYWORD OUTPUT: ; match_qual - The quality of the match. It is a structure with ; the following fields. ; .dpe = # of DPE levels away from the actual ; 0 is perfect match, positive means DC ; is longer exposure than requested. ; .time = the number of days from the requested day. ; positive is for DC images taken past ; the time of the input image ;OPTIONAL KEYWORD INPUT: ; times - MANUAL METHOD: A list of input times (any of 3 standard ; formats) ; res - MANUAL METHOD: The resolution required (0=full res, ; 1=half res, 2=quarter res) ; dpe - MANUAL METHOD: The DPE level (exposure level) ; nodata - If set, just return the DC_INDEX and IMAP, do not read ; the data. ; dc_interpolate - If set, then perform interpolation on the dark frame ; images to get a better backround subtraction. Deafult is to ; take the image closest in exposure. ; interpolate - Same as /DC_INTERPOLATE ; orbit_correct - Use the orbital correction (must also specify the ; interpolate keyword for this to work). ;OPTIONAL KEYWORD INPUT/OUTPUT: ; unc_data - The uncertainty in the dark current, due to ; decompression and interpolation. Note that this will only ; be calculated if some non-zero variable is passed in using ; this keyword. ;METHOD: ; The program will try to make an exact match to the exposure ; level passed. If it cannot match, it will take the next exposure ; larger which is available. ; ; If the images requested for extraction match the images that are ; passed in, then the data file is not read again. ; ; GET_DC_IMAGE can be called manually by passing "times", "res" ; and "dpe" using the keywords. All three need to be passed, however ; it is possible to pass an array of times, but have res or dpe ; be a scalar. ;HISTORY: ; Written 20-Apr-92 by M.Morrison ; 21-Apr-92 (MDM) - General mods/additions ; 1-May-92 (MDM) - Fixed trouble with "bad_match" option ; 15-Oct-92 (MDM) - Added /NODATA option ; 6-May-93 (MDM) - Modified to work on VMS (call FILE_LIST) ; 10-Jun-93 (MDM) - Slight modification to the /NODATA option ; 8-Jul-93 (MDM) - Modified to check DC_ROADMAP variable instead of ; DC_INFIL variable to determine if need to read SDC ; roadmap ; 16-Aug-93 (MDM) - Expanded search criteria such that for full resolution ; PFI and ORs, it looks to get the FFI image which ; will most likely overlap the PFI/OR (this is because ; the FFI FR images are taken in 1024x512 strips). ; 18-Aug-93 (MDM) - Changed to use DIR_SXT_SDC ; 24-Aug-93 (MDM) - Adjusted criteria for FR search so that it works on ; older data when there was only one FR dark image per set ;V3.00 24-Aug-93 (MDM) - Extracted much of the selection code into the routine ; SEL_DC_IMAGE (allowed cooled and warm to go through ; one routine) ; - Added TEMP_INPUT option ; - Added /INTERPOLATE option ; 25-Aug-93 (MDM) - Use EXT_SUBIMG3 ; - Modification since SEL_DC_IMAGE now returns the ; MK_DSET_STR structure (speeds up access to data file) ;V3.01 1-Sep-93 (MDM) - Corrected code for case where SEL_DC_IMAGE does not ; return a proper image ;V3.02 7-Sep-93 (MDM) - Corrected error for when it is called for long time ; spans. It was using the same image over the whole ; span when the characteristics matched. ;V3.03 16-Dec-93 (MDM) - Modified SEL_DC_IMAGE to not chose dark images on ; the "other side" of 13-Nov-92 ; - Corrected code for case where no matches are found ;V3.04 14-Jan-94 (MDM) - Added code to not crash if no SDC datasets are available ; - Added /DC_INTERPOLATE option ; - Added QERR ;V3.05 3-Jun-94 (MDM) - Modification to not crash when no SDC datasets are ; available AND the /INTERPOLATE option is set. ; 21-Sep-94 (MDM) - Updated documention header ;V3.06 30-Nov-94 (MDM) - Updated SEL_DC_IMAGE logic quite a bit ; - Logic changed so that the first exposure shorter and ; longer than the input exposure is selected during ; /DC_INTERPOLATE ; - Modification to not perform interpolation if the ; DPE of the exposure found matches identically. ; Did it by checking MATCH_QUAL.DPE - easy! ;V3.07 1-Dec-94 (MDM) - Changed the logic on printing messages about ; insufficient DC images available ; 1-Dec-94 (MDM) - Added /FLOAT option ; 15-Dec-94 (MDM) - Modified to have /INTERPOLATE option work with ; input specified manually ;V3.08 15-Dec-94 (MDM) - Modified to have the DC_INDEX DPE/FilterA/ExpDur ; updated accordingly if interpolation was used. ;V3.09 7-Feb-95 (MDM) - Updated to use GT_EXPDUR with the /SHUTDUR when ; finding out exposure factor when doing interpolation. ; 8-feb-95 (jmm) - Added UNC_DATA keyword, for uncertainty calculation ;V3.10 27-Feb-95 (MDM) - Put McTiernan changes on-line ; 2-Mar-95 (SLF) - protect agains udata0 not defined (why??) ; 2-Mar-95 (MDM) - Define UDATA0 for case where a scalar dark current ; value is used. ;V3.11 2-Feb-98 (TRM) - Added tfactor variable to account for the orbital ; variation of the dark current. ; 23-Mar-99 (TRM) - Added tfms to set of criteria which determines if ; images can use the same dark frame, when ; /orbit_correct is set. This was a bug fix. ; Changed the critera to a string from long so that ; tfms would fit. ; 02-Apr-1999 (TRM) Added week id into the val test so that new images ; are selected when a week bounday is crossed. ; 11-Jan-2001 (LWA) [message, 'Interpolation of dark images ; enabled', /info] commented out. ; 3-Dec_2001 (LWA) Commented out "perfect match" print statement. ;- ; common get_dc_image, qprint1, last_time ; progverno = 3.11*1000 ; quncert = keyword_set(unc_data) qinterp = keyword_set(interpolate) or keyword_set(dc_interpolate) if (getenv('ys_dark_interp') ne '') then qinterp = 1 if (n_elements(qprint1) eq 0) then begin qprint1 = 1 end else begin del = int2secarr(!stime, last_time)/60 qprint1 = (del gt 10) ;it's been over 10 minutes end if (qinterp and qprint1) then begin ;only want the message printed once per hour ; message, 'Interpolation of dark images enabled', /info last_time = !stime qprint1 = 1 end ; sel_dc_image, struct_in, dc_infil, dset, imap, match_qual=match_qual, qstop=qstop, qdebug=qdebug, qmanual=qmanual, $ times_input=times_input, res_input=res_input, dpe_input=dpe_input, temp_input=temp_input, q2debug=q2debug ; ;----- imap is the map of source dark frames to use, imap_out is the map of what is output since the same frame can be ; used multiple times if it is an OR which moves, but same exposure ; if (keyword_set(dc_infil)) then begin ;if it found something rd_xda, dc_infil, dset, index, data, dc_roadmap, /nodata if (keyword_set(nodata)) then return ;just return the dc_index qerr = 0 end else begin qerr = 1 end nin = n_elements(imap) ;number of input images specified ; if (qmanual) then begin val = res_input + dpe_input*64 nx = 1024./(2^min(res_input)) ny = nx end else begin xcorn = fix(gt_corner(struct_in, /lower, /x)/4) ;0-255 ycorn = fix(gt_corner(struct_in, /lower, /y)/4) res = gt_res(struct_in) dpe = gt_dpe(struct_in) ;week = tim2weekid(struct_in) week = strarr(n_elements(struct_in)) for iweek=0,n_elements(week)-1 do week(iweek) = tim2weekid(struct_in(iweek)) if (not qinterp and not qerr) then dpe = gt_dpe(index(imap)) ;for the case where we are using the closest image, we don't ;want to save the same dark image twice just becuase there are two ;exposures of similar duration ;; ;do all FR in order of DPE, then all HR, then all QR ;; val = xcorn + ycorn*256L + dpe*256L*256L + res*256L*256L*64 + imap*256L*256L*256L ;; ; low 8-bits 0:7 second 8-bits 8:15 third 8-bits high part of 3rd fourth 8-bits ;; ;^^ Added +imap 7-Sep-93 ;; if keyword_set(orbit_correct) then begin ;; ; 4 bits is not enough! Re think this scheme! ;; val = val + bytscl(tfms,min=0.,max=70.,top=63)*256L*256L*256L*64L ;; endif ; TRM made val a string to avoid bit overflow problems 1999-Mar-23 ; TRM added week id string so that the dark frame returned will be for ; the same week as the data. ;do all FR in order of DPE, then all HR, then all QR val = string(xcorn,format='(i5)') + $ string(ycorn,format='(i5)') + $ string(dpe,format='(i2)') + $ string(res,format='(i1)') + $ string(imap,format='(i5)') + $ week if keyword_set(orbit_correct) then begin tfms = sxt_uvf_info(struct_in,/tfms) ; 0.01 minute granularity essentially ; guarantees that all val's will be unique val = val + string(fix(tfms*100.+0.5),format='(i5)') endif ispaces = where(byte(val) EQ 32b,nspaces) ; Convert leading spaces to zeroes so that the sort will work if nspaces GT 0 then begin vtemp = byte(val) vtemp(ispaces) = 48b val = string(vtemp) endif nx = max(gt_shape(struct_in, /x)) ny = max(gt_shape(struct_in, /y)) end if (n_elements(val) eq 1) then ss_val = 0 else ss_val = uniq(val, sort(val)) ;in terms of datasets of input images uval = val( ss_val ) imap_out = intarr(nin) for i=0,n_elements(uval)-1 do imap_out(where(val eq uval(i))) = i ;gotta be a better way nout = n_elements(uval) ;number of uniq dark images if (keyword_set(float)) then data = fltarr(nx, ny, nout) $ else data = intarr(nx, ny, nout) if(quncert) then unc_data = data ;jmm, 8-feb-95 if (not qerr) then begin index = replicate(index(0), nout) his_index, /enable his_index, index end if (keyword_set(qdebug)) then print, 'qmanual,qinterp',qmanual,qinterp if (qinterp) then begin sel_dc_image, struct_in, dc_infil1, dset1, imap1, match_qual=match_qual1, $ times_input=times_input, $ res_input=res_input, dpe_input=dpe_input, temp_input=temp_input, /short, $ qdebug=qdebug, q2debug=q2debug sel_dc_image, struct_in, dc_infil2, dset2, imap2, match_qual=match_qual2, $ times_input=times_input, $ res_input=res_input, dpe_input=dpe_input, temp_input=temp_input, /long, $ qdebug=qdebug, q2debug=q2debug, mmdpe=match_qual1.dpeval ;MDM added 30-Nov-94 end ;------------------------------------------------------------------------------ ; Loop through each dark image to be read ;------------------------------------------------------------------------------ qprint2 = 1 for i=0,nout-1 do begin idset = ss_val(i) ;the dataset of the uniq input "struct_in" image that is being processed qextract = 0 if (not qmanual) then begin ;a index or roadmap structure was passed in qffi = gt_pfi_ffi(struct_in(idset), /true_ffi) ;0 for PFI,OR,FFI2PFI - 1 for true FFI qextract = (qffi eq 0) end qavoid = 0 if (qinterp) then begin if (idset(0) eq -1) then qavoid = 1 ;MDM added conditional 3-Jun-94 if (not qavoid) then qavoid = (imap1(idset) eq -1) or (imap2(idset) eq -1) ;added 1-Sep-93 if (qavoid) then message, 'Insufficient DC images for interpolation. Using single image', /info if (not qavoid) then begin ; TRM add orbit_correct to guarantee that the orbital correction ; is applied. qavoid = (match_qual(idset).dpe eq 0 AND $ NOT keyword_set(orbit_correct)) ;MDM added 30-Nov-94 if (qavoid) then begin qprint2 = 0 ; message, '/DC_INTERPOLATE but found perfect match. No interpolation needed', /info end else begin if (not qavoid) then qavoid = (match_qual1(idset).dpeval eq match_qual2(idset).dpeval) ;MDM 30-Nov-94 if (qavoid and not keyword_set(orbit_correct)) then message, 'Trouble with DC images for interpolation. Using single image', /info end end end if (qinterp and (not qavoid)) then begin ;---------------------------------------- Interpolation case dc_dset1 = dset1(imap1(idset)) ;imap converts between input dataset number and dark image dataset number dc_dset2 = dset2(imap2(idset)) rd_xda, dc_infil1, dc_dset1, index1, data1, dc_roadmap ;short exposure rd_xda, dc_infil2, dc_dset2, index2, data2, dc_roadmap ;long exposure his_index, index1 if (qextract) then begin data1 = ext_subimg3(index1, temporary(data1), ref_image=struct_in(idset), /update_index, /use_avg, qdebug=qdebug) data2 = ext_subimg3(index2, temporary(data2), ref_image=struct_in(idset), /update_index, /use_avg, qdebug=qdebug) end if(quncert) then begin ;jmm, 8-feb-95 data1 = sxt_decomp(temporary(data1), udata1) data2 = sxt_decomp(temporary(data2), udata2) end else begin data1 = sxt_decomp(temporary(data1)) data2 = sxt_decomp(temporary(data2)) end struct_in0 = struct_in(idset) exp1 = gt_expdur(index1, /shutdur) ;short exp2 = gt_expdur(index2, /shutdur) ;long if (qmanual) then dpe0 = dpe_input $ else dpe0 = gt_dpe(struct_in0) exp0 = gt_dpe(dpe0, /conv) tags = tag_names(struct_in0) if (tags(0) eq 'GEN') then exp0 = gt_expdur(struct_in0, /shutdur) ; ; tfactor added by TRM 1998-Feb-02, Correction for orbital variation if keyword_set(orbit_correct) then begin tfactor = dc_orbit_correct(index2,struct_in0,[0]) tfactor = tfactor(0) ;;print,'Dark: ',fmt_tim(index2),sxt_uvf_info(index2,/tfms) endif if (n_elements(tfactor) eq 0) then tfactor = 1.0 ;what is this? can't remember factor = (exp0-exp1) / (exp2-exp1) factor = factor(0) ;make scalar if (factor eq 0) then factor = tfactor else factor = factor*tfactor if (keyword_set(float)) then begin ;jmm, 11-feb-95 data0 = (data2-data1)*factor + data1 if (quncert) then udata0 = sqrt(((1.-factor)*udata1)^2+(factor*udata2)^2) end else begin data0 = fix( (data2-data1)*factor +0.5) + data1 if (quncert) then udata0 = fix(0.5+sqrt(((1.-factor)*udata1)^2+(factor*udata2)^2)) ;The sqrt must be done as floats, 11-feb-95,jmm end time_dark = [gt_time(index1), gt_time(index2)] day_dark = [gt_day(index1), gt_day(index2)] if (qdebug and not qmanual) then begin print, 'Input Image Time: ', fmt_tim(struct_in0) print, 'Short Exposure: ', fmt_tim(index1), ' Long Exposure: ', fmt_tim(index2) print, 'Short Exposure: ', exp1, ' Long Exposure: ', exp2 print, 'Input Expsoure: ', exp0 print, 'Factor = ', factor print, 'TFactor = ', tfactor end index0 = index1 update_dpe, index0, dpe0 end else begin ;---------------------------------------- Simple closest image case if (imap(idset) eq -1) then begin message, 'No dark data images available. Using scalar value', /info struct_in0 = struct_in(idset) data0 = intarr(nx, ny) + get_res_offset(struct_in0) if (quncert) then dummy = sxt_decomp(data0, udata0) ;MDM added 2-Mar-95 index = struct_in(idset) ;hope that it was an index that was passed in (?) index0 = struct_in(idset) ;hope that it was an index that was passed in (?) time_dark = get_res_offset(struct_in0) day_dark = 0 end else begin ; TRM 1999-Mar-22 if (keyword_set(orbit_correct)) then $ message,/info,'/orbit_correct set, but only one exposure. No orbital correction' dc_dset = dset(imap(idset)) ;imap converts between input dataset number and dark image dataset number rd_xda, dc_infil, dc_dset, index0, data0, dc_roadmap his_index, index0 if (qextract) then data0 = ext_subimg3(index0, temporary(data0), ref_image=struct_in(idset), $ /update_index, /use_avg, qdebug=qdebug) if (quncert) then data0 = sxt_decomp(temporary(data0), udata0) $ ;jmm, 8-feb-95 else data0 = sxt_decomp(temporary(data0)) time_dark = gt_time(index0) day_dark = gt_day(index0) end end data(0,0,i) = temporary(data0) if (quncert) and keyword_set(udata0) then unc_data(0, 0, i) = temporary(udata0) ;jmm, 2-17-95 index_copy = str_copy_tags(index(0), index0) ;copy the "index0" into the output array structure type index(i) = index_copy his_index, index, i, 'q_dark_sub', progverno his_index, index, i, 'time_dark', time_dark his_index, index, i, 'day_dark', day_dark end ; if (keyword_set(qstop)) then stop end