pro fg_get_flat, index, flat_index, flat_image $ , flatdir=flatdir, deltaday=deltaday $ , flatfile=flatfile $ , version=progver, name=prognam $ , retrn=retrn ;+ ; NAME: ; FG_GET_FLAT ; ; PURPOSE: ; Return an FG flat field image appropriate for the index ; structure. This routine is called from FG_FLATFIELD by default ; when a user-provided flat image is not available. ; ; CATEGORY: ; FITS processing ; ; SAMPLE CALLING SEQUENCE: ; fg_get_flat, index, flat_index, flat_image $ ; [, flatdir=flatdir, deltaday=deltaday, flatfile=flatfile $ ; , version=version, name=name ] ; ; INPUTS: ; INDEX - the SSW index structure ; ; OUTPUTS: ; FLAT_INDEX - the SSW header index structure for the flat field ; FLAT_IMAGE - 2D full FOV flat field image without ROI extraction ; ; OPTIONAL INPUT KEYWORD PARAMETERS: ; DARKDIR - Alternate directory for processed dark current fits images, ; default = concat_dir('$SOT_SSWDB','flats') ; ; OUTPUT KEYWORD PARAMETERS: ; DELTADAY - Time difference (day) between the image and the flat ; FLATFILE - Filename of the flat field FITS file ; VERSION - The current version number of the program ; NAME - The name of this routine ; RETRN - Flag for branching on unusual return ; ; PROCEDURE: ; ; MODIFICATION HISTORY: ;V0.1 16-Apr-2007 Created by YK ;V0.2 24-Apr-2007 Changed variable named 'SORT' to ST to avoid ; internal routine name. TEB. ;v1.0 made sum and bin checking more robust. 24-Oct-2007. TEB. ;v1.1 Corrected "img_waveid le 13" to "le 23". TEB 25-Oct-2007. ;v1.1.1 Return deltaday=-1 with unity flatfield. TEB 2-Nov-2007. ;v1.1.2 Minor bug fix on line 91. TEB 31-Jan-2008. ;v1.2 "No Ca II 3968 flat" switch removed. TEB 4-Feb-2009. ; ;- ;----------------------------------------------------------------------------- ; prognam = 'FG_GET_FLAT.PRO' progver = 'V1.2' ; default directory if N_ELEMENTS(flatdir) eq 0 then flatdir = concat_dir(getenv('SOT_SSWDB_FG_CAL'),'flats') ; get size, wavelength & summing img_nx0 = gt_tagval(index[0],/naxis1) img_ny0 = gt_tagval(index[0],/naxis2) img_waveid = gt_tagval(index[0], /waveid) img_psum = gt_tagval(index[0], /campsum) img_ssum = gt_tagval(index[0], /camssum) img_binx = gt_tagval(index[0], /fgbinx) img_biny = gt_tagval(index[0], /fgbiny) img_t_fgccd = gt_tagval(index[0], /t_fgccd) img_t_fgceb = gt_tagval(index[0], /t_fgceb) ; wavelength as string wave_bfi = ['3883', '3968', '4305', '4504', '5550', '6684'] wave_nfi = ['5172', '5250', '5576', '6302', '6563', '5896'] set_file: if img_waveid ge 1 and img_waveid le 6 then $ wave = 'bfi_'+wave_bfi[img_waveid-1] $ ; BFI else if img_waveid ge 8 and img_waveid le 23 then $ wave = 'nfi_'+wave_nfi[(img_waveid-8) mod 10] $ else begin MESSAGE,/CONT,'The wavelength is not supported now. Returning...' flat_index = -1 flat_image = -1 retrn = 1 deltaday = -1 RETURN endelse ; summing sum = '1x1' if (img_ssum gt 1) or (img_psum gt 1) then sum = STRING(img_psum,img_ssum, form='(I1,"x",I1)') if (img_binx gt 1) or (img_biny gt 1) then sum = STRING(img_binx,img_biny, form='(I1,"x",I1)') ;if img_waveid eq 3 then sum='1x1' ;use the 1x1 G-band flat only. ; return file list file = 'flat_'+wave+'_'+sum+'_'+'*.fits' filelist = FILE_SEARCH(flatdir, file, count=nfile) if nfile le 0 then begin if (img_ssum gt 1) or (img_psum gt 1) then begin ;summed flat doesn't exist img_ssum = 1 img_psum = 1 GOTO, set_file end if (img_binx gt 1) or (img_biny gt 1) then begin ;same as above img_binx = 1 img_biny = 1 GOTO, set_file end ; If not found, return array of 1s. MESSAGE, /CONT, 'No flat file found. Returning unity array...' flat_index = index flat_image = REPLICATE(1.0,img_nx0, img_ny0) flatfile = 'No flatfield found in database. Using unity array.' deltaday=-1 RETURN endif ; Check time difference ; read header to get DATE_OBS mreadfits, filelist, flat_index, /nod utc_obs = str2utc(index[0].date_obs) utc_flat = str2utc(flat_index[0].date_obs) deltaday = (utc_flat.mjd-utc_obs.mjd) + (utc_flat.time - utc_obs.time)/8.64e7 ; use minimum DT deltaday = MIN(ABS(deltaday),imin) flatfile = filelist[imin] ; read header again mreadfits, flatfile, flat_index, /nodata ; Read flat file case img_waveid of 6: begin ; Red continuum n_grid = gt_tagval(flat_index, 'N_TGRID') flat_t_fgccd = fltarr(n_grid) for i=0, n_grid-1 do begin index_ex = fitshead2struct(headfits(flatfile,ext=i+1)) flat_t_fgccd[i] = gt_tagval(index_ex[0], 'T_FGCCD') endfor ; find closest and 2nd closest temperatures st = SORT(abs(img_t_fgccd-flat_t_fgccd)) dt0 = img_t_fgccd - flat_t_fgccd[st[0]] dt1 = img_t_fgccd - flat_t_fgccd[st[1]] ; Get flat image if dt0*dt1 ge 0 then begin ; Return the closest one if the temperature is out of the range flat_image = mrdfits(flatfile, st[0]+1, /fscale, /silent) ; Modify output header flat_index.t_fgccd = flat_t_fgccd[st[0]] endif else begin ; Linear interolation if the temperature is within the range flat0 = mrdfits(flatfile, st[0]+1, /fscale, /silent) flat1 = mrdfits(flatfile, st[1]+1, /fscale, /silent) flat_image = (-flat0*dt1 + flat1*dt0)/(dt0-dt1) flat_image = flat_image>0<2 flat_index.t_fgccd = img_t_fgccd endelse ; Update NAXIS* keyword for image flat_index.naxis = 2 flat_index = add_tag(flat_index,index_ex.naxis1,'NAXIS1',index='NAXIS') flat_index = add_tag(flat_index,index_ex.naxis2,'NAXIS2',index='NAXIS1') end Else: begin flat_image = mrdfits(flatfile, /fscale, /silent) end end RETURN END