pro eit_prep,filename,outheader,outimage, $ ; ;----------------------Calibration, processing options no_calibrate = no_calibrate, cosmic=cosmic,$ surround=surround,fill_block=fill_block,$ no_index=no_index,$ ; ;----------------------Output controls outfits=outfits,nodata=nodata,outdir=outdir,$ ; ;----------------------Misc. verbose=verbose,save_zero=save_zero, $ nrl=nrl,n_block=n_block,data=data, $ image_no=image_no, noprompt = noprompt,no_roll=no_roll, _extra=extra ;+ ; NAME: ; EIT_PREP ; ; PURPOSE: ; Process an EIT image. ; ; Processing routine for preparing Level-1 EIT data. ; Presently this version reads in the named FITS ; image and produces a background subtracted, degridded, ; flat-fielded, and degradation corrected output array or FITS file. ; In the future this routine will also vignette correct. ; The input may also be the output index structure and data cube as ; produced by read_eit. ; ; CATEGORY: ; FITS processing ; ; CALLING SEQUENCE: ; EIT_PREP,filename,[outheader,outimage] [data=data] [,/outfits] ; [,/verbose] [,/cosmic] [,/no_calibrate] [,/no_index] ; [,/save_zero] [,/nrl] [n_block = n_block] ; [,/surround] [,/nodata] [,outdir = outdir] ; [,image_no = image_no] ; ; INPUTS: ; filename - FITS file name or the output of eit_catrd ; may also be index structure (e.g. output of read_eit) ; ; ; OPTIONAL INPUT KEYWORD PARAMETERS: ; outfits - set if want to produce an output FITS file ; (or set string filename) ; outdir - Destination directory for prepped FITS files ; (if outdir is set, /outfits is assumed) ; nrl - set if input file was processed at NRL ; surround - set to replace missing block data with the average of ; surrounding blocks. ; image_no - Image number(s) (or 'ALL') to process in multiple image LZ file. ; fill_block - intarr(32,32) which is used to replace missing data blocks ; save_zero- set to retain images that have only 0 counts. ; verbose - set for a few more messages ; data - data array(cube) corresponding to index structure ; nodata - set if don't want output array (auto set if lt 3 params) ; cosmic - set if want to remove cosmic rays - may remove small ; real features ; no_calibrate - set if want "raw" returned images, i.e. only background subtracted ; no_index - set if don't want to use recommended Index correlation ; for final response scaling ; no_roll - do not correct for SOHO roll (mainly for testing) ; ; OPTIONAL OUTPUTS: ; outheader- FITS header or output index structure ; outimage - Processed EIT image. Default data type is Float ; ; OPTIONAL OUTPUT KEYWORD PARAMETERS: ; n_block - Number of blocks repaired ; ; COMMON BLOCKS: none. ; ; SIDE EFFECTS: ; If /outfits is set, then a FITS file is written. ; ; RESTRICTIONS: ; 1) Processes only one 3-Dim file at a time, either a single image ; vector, or the entire file. ; ; 2) Does not process a mixture of 2D and 3D files ; ; 3) If filename is a vector, there is no checks to make sure ; that all the images are the same size. User must be careful ; to check this. ; ; PROCEDURE: ; The FITS file is read in. Any missing blocks are replaced ; with the present dark noise value. The dark noise is ; subtracted off the entire image. The image is degridded. ; The image is flat fielded and corrected for degradation. ; The image is optionally returned. ; A new FITS file with an updated header is written out ; with the format eit_l1_yyyymmdd_hhmmss.fits ; ; Note: the De-gridding and Flat-fielding are multiplicative ; operations with an AVERAGE value over the full field ; of 1.0 ; ; MODIFICATION HISTORY: ; Written by: J. Newmark Date. July 1996 ; 1996 Sep 24 - J. Newmark - add flat fielding ; 1996 Oct 8 - J. Newmark - allow NRL format files ; 1996 Oct 16, J. R. Lemen, Changed calling arguments. ; 1996 Dec 09 - J. Newmark - add usage of image_no throughout. ; 1997 Feb 10 - J. Newmark - allow input if index structure and data ; cube ; 1997 Mar 09 - J. Newmark - write out FITS files from indx/data ; structures ; 1997 Mar 10 - J. Newmark - fix output FITS file names. ; 1997 Apr 09 - J. Newmark - account for new exptime ; 1997 Apr 12 - J. Newmark - account for correction of Al+1 to ClEAR ; 1997 Aug 14 - J. Newmark - implemented Handy requests of outdir, use ; of nodata ; 1997 Aug 18 - J. Newmark - updated conversion of struct->fitshead ; fix up code for /surround ; clean up code for output file name ; 1997 Oct 06 - J. Newmark - added preliminary response normalization ; 1997 Oct 30 - J. Newmark - update header for LZ Partial Fov images ; 1997 Nov 10 - J. Newmark - add new grids after July 1997 ; ; Version 2.0 ; 1998 Mar 27 - J. Newmark - support ALL or vector images in single 3D ; file, speed up 3D processing, ; clean up history, pathway checks, ; update documentation, fix n_block bug ; added adjust keyword ; 1998 Jun 3 - D. Zarro (SAC/GSFC) ; - sprinkled liberally with calls to temporary ; 1998 Jul 22 - J. Newmark - add phot_norm keyword ; 1998 Jul 23 - J. Newmark - add cosmic keyword ; 1998 Dec 08 - J. Newmark - change adjust keyword ; 2000 Jan 05 - J. Newmark - new response normalization, deleted ; phot_norm keyword ; 2000 Aug 04 - J. Newmark - NOPROMPT keyword for batch processing ; large data sets, cleaned up memory managment ; 2001 Mar 15 - F. Auchere - added time varying dark current ; ;Version 3.0 ; 2001 May 18 - J. Newmark/F. Auchere - new response normalization ; 2001 Jun 21 - N. Rich - Add NRL check at (near) line 332 ; 2001 Aug 15 - Zarro (EITI/GSFC) - add optional FITS out file name in ; OUTFITS ; ;Version 4.0 ; 2001 Dec 10 - J. Newmark - changed default options, ; new calibration, new default options lead ; to default output type as float ; 2002 Jan 12 - F. Auchere - modified call to eit_dark to handle ; different offset in binned images ; 2003 Apr 17 - D. Biesecker - modified call to is_fits to correctly ; handle collapsed directory path names ; 2003 Jul 17 - C.A. Young - added 180 degree roll correction ; 2003 Sep 13 - C.A. Young - added fix for multifile list ; 2004 Feb 4 - D.M. Zarro - Update header for rotation correction ; 2008 March 25 - Zarro (ADNET) ; - fixed a bug in roll correction that only ; occurred when an array of index structures is input ; 2008 August 20 - F. Auchere - no exptime normalize if exptime ; is 0 ; 2010 January 10 - Zarro (ADNET) ; - increased internal memory limit to 1GB ; - replaced calls to is_fits with valid_fits ; 2013 August 5 - D.M. fecit ; - Since valid_fits returns only a single status, ; - added file_search call to determine if file exists ; - (prompted by cases of catalog entries with no ; - corresponding files in archive) ; 2013 August 23 - N.Rich ; - Previous mod breaks this for /NRL; omit eit_file2path in that case. ; 2014 March 1 - Zarro (ADNET) ; - Merged Zarro and D.M. Fecit mods. ; - Use SWITCH to check for all likely file ; locations. ;- utc_date_19960327 = anytim2utc('1996/03/27') utc_date_19960323 = anytim2utc('1996/03/23') utc_obe_good = anytim2utc('1996/07/18') t_obe = utc_obe_good.mjd + 1.d-3*utc_obe_good.time/86400. flag_3d = 0 use_fits = 1 filter_factor = [0.49, 0.49, 0.33, 0.29] wavelength = [171, 195, 284, 304] if is_string(outfits) then outfile=outfits if n_elements(outdir) then outfits = 1 else outdir = '' yes_data = ~keyword_set(nodata) or n_params() eq 3 t0 = systime(1) ; Keep track of running time if datatype(filename) eq 'STC' then begin if keyword_set(data) then begin use_fits = 0 image_no = 0 endif else begin message, /info, 'No Data supplied, reading FITS files...' filename = filename.filename endelse endif nfiles = n_elements(filename) vv = fltarr(nfiles) n_block = intarr(nfiles) if use_fits then begin ;-- use SWITCH to check for all likely input file locations and formats switch 1 of 1: begin files=filename existing_file = file_search(files[0], count = n_existing_files) if n_existing_files gt 0 then break end 2: begin files=eit_file2path(filename,/check_eit_environment) existing_file = file_search(files[0], count = n_existing_files) if n_existing_files gt 0 then break end 3: begin files=eit_file2path(filename,/check_eit_environment,/collapse) existing_file = file_search(files[0], count = n_existing_files) if n_existing_files gt 0 then break end else: begin message, '$EIT_PREP-E-NO_FILE, no file found with path ' + filename[0] + '.',/info return end endswitch ;-- if we got here, then we located the file. Now check it's a valid FITS file stat = valid_fits(files[0]) if stat eq 0 then begin message,'The file[s] do not appear to be FITS format.',/info return endif if yes_data then begin noprompt = keyword_set(noprompt) or (get_logenv('ssw_batch') ne '') prompt = 1 - noprompt tsize = total(file_stat(files,/size))/1.e6 if tsize ge 1000 and prompt then begin print,'You have requested '+string(tsize)+' MB of data' print,'You can select the output to be only FITS files if this '+$ 'is too large for your machine"s memory.' ans='' read,ans,prompt='Enter C (continue) or F (FITS only): ' if strupcase(ans) eq 'F' then begin yes_data = 0 outfits = 1 endif endif endif if keyword_set(image_no) and nfiles gt 1 then begin message,'This program only works on 1 3-Dim file at a time.',/info return endif endif i = 0 repeat begin ; MAIN LOOP: begins if ((i+1) mod 10) eq 0 then message,/info, $ string('Working on',i+1,'/',strtrim(nfiles,2), $ ' ETC=',(nfiles-i+1)*(systime(1)-t0)/60./float(i+1), $ ' mins.',form='(a,i4,3a,f5.1,a)') if use_fits then begin thisfile = files[i] if n_elements(image_no) eq 0 then begin image = readfits(thisfile,header,silent=~keyword_set(verbose)) image_no = 0 endif else begin if flag_3d then image_no = images3d[i] header = headfits(thisfile) if n_elements(image_no) gt 1 then begin nfiles = n_elements(image_no) files = replicate(files,nfiles) images3d = image_no image_no = images3d[0] flag_3d = 1 vv = fltarr(nfiles) n_block = intarr(nfiles) endif else if strupcase(image_no) eq 'ALL' then begin nfiles = eit_fxpar(header,'naxis3') files = replicate(files,nfiles) images3d = indgen(nfiles) image_no = 0 flag_3d = 1 vv = fltarr(nfiles) n_block = intarr(nfiles) endif if image_no eq 0 then begin image = readfits(thisfile,silent=~keyword_set(verbose)) endif else begin image = readfits(thisfile,silent=~keyword_set(verbose),$ nslice = image_no) endelse endelse if keyword_set(nrl) then header = nrl2eit(header) endif else begin header = filename[i] thisfile = eit_file2path(header.filename) image = data[*,*,i] endelse ; set final keyword for degrid routine if strpos(strlowcase(thisfile), 'efz') ge 0 or strtrim(eit_fxpar(header,$ 'DATASRC'),2) eq 'LZ file' then final = 1 else final = 0 date_obs = eit_fxpar(header,'date_obs',image_no=image_no) ;== Added 15 March 2001 F. Auchere utc_date_obs = anytim2utc(date_obs) hist = eit_fxpar(header,'history') IF keyword_set(NRL) THEN version = 3.0 ELSE version = float(strmid(hist(0),8,3)) ;==Add NRL check, 21Jun01, nbr ; Fix missing blocks, dark subtract: if keyword_set(fill_block) and keyword_set(surround) then $ message,'Cannot set both fill_bock and surround' if keyword_set(cosmic) then begin image = cosmicr(image) tagval = 'Replaced Cosmic Rays with Median Filtering' if use_fits then fxaddpar, header,'HISTORY', tagval else $ header = rep_tag_value(header,[header.history,tagval],'HISTORY') endif bin = round(eit_fxpar(header, 'CDELT1')/eit_pixsize()) case 1 of keyword_set(surround): begin if not keyword_set(no_degrid) then $ image = raise_missing_blocks(image, repl=eit_dark(date_obs,bin=bin),n_block=nblock,/no_copy) $ else image = raise_missing_blocks(image, /surround,n_block=nblock,/no_copy) ;== Modified 15 March 2001 F. Auchere tagval = 'Replaced missing blocks with SURROUND' end keyword_set(fill_block) : begin image =raise_missing_blocks(image,fill_block=fill_block,n_block=nblock,/no_copy) tagval = 'Replaced missing blocks with FILL_BLOCK' end else: begin image = raise_missing_blocks(image, repl=eit_dark(date_obs, bin=bin),n_block=nblock,/no_copy) ;== Modified 15 March 2001 12 jan 2002 F. Auchere tagval = 'Replaced missing blocks with EIT_DARK' end endcase n_block(i) = nblock if use_fits then fxaddpar, header,'HISTORY', tagval else $ header = rep_tag_value(header,[header.history,tagval],'HISTORY') image = temporary(image) - eit_dark(date_obs, bin=bin) ;== Modified 15 March 2001, 12-jan 2002 F. Auchere tagval = 'Subtracted dark current, EIT_DARK' if use_fits then fxaddpar, header,'HISTORY', tagval else $ header = rep_tag_value(header,[header.history,tagval],'HISTORY') ; calibrate - degrid, flat field, normalize by exposure time, time varying ; response correction if ~keyword_set(no_calibrate) then begin image = eit_degrid(image, header, final = final, image_no = image_no,$ /no_copy) if keyword_set(surround) then $ image = raise_missing_blocks(image, /surround,n_block=nblock,/no_copy) tagval = 'Degridded image with EIT_DEGRID' if use_fits then fxaddpar, header,'HISTORY', tagval else $ header = rep_tag_value(header,[header.history,tagval],'HISTORY') image = eit_flat(image, header, image_no = image_no,/no_copy) tagval = 'Flat fielded image with EIT_FLAT' if use_fits then fxaddpar, header,'HISTORY', tagval else $ header = rep_tag_value(header,[header.history,tagval],'HISTORY') if version(0) ge 2.0 then shut_time = 0 else begin t_obs = utc_date_obs.mjd + 1.d-3*utc_date_obs.time/86400 shut_time = eit_fxpar(header,'SHUTTER CLOSE TIME',image_no=image_no) if shut_time le 0 or t_obs lt t_obe then shut_time = 2.1 endelse exptime=eit_fxpar(header,'exptime',image_no=image_no)+shut_time if exptime gt 0.0 then image = temporary(image) / exptime tagval = 'Exposure normalized (per sec)' if use_fits then fxaddpar, header,'HISTORY', tagval else $ header = rep_tag_value(header,[header.history,tagval],'HISTORY') filter = strupcase(strtrim(eit_fxpar(header,'filter',image_no=image_no),2)) wave = eit_fxpar(header,'wavelnth',image_no=image_no) i_wave = where(wavelength eq wave) if (not final) and (utc_date_obs.mjd lt utc_date_19960327.mjd) and $ (utc_date_obs.mjd ne utc_date_19960323.mjd) then begin if filter eq 'CLEAR' then image = temporary(image) / filter_factor(i_wave[0]) endif else if filter eq 'AL +1' then image = temporary(image)/filter_factor(i_wave[0]) if n_elements(leak_284) ne 0 and filter eq 'CLEAR' and wave eq 284 then $ image = temporary(image)/leak_284 tagval = 'Flux normalized to CLEAR Filter' if use_fits then fxaddpar, header,'HISTORY', tagval else $ header = rep_tag_value(header,[header.history,tagval],'HISTORY') wave = eit_fxpar(header,'wavelnth',image_no=image_no) image = temporary(image) / eit_norm_response(utc_date_obs,wave,header,$ no_index=no_index) tagval = 'Flux normalized for time response changes' if use_fits then fxaddpar, header,'HISTORY', tagval else $ header = rep_tag_value(header,[header.history,tagval],'HISTORY') endif else image = fix(round(temporary(image))) ; ; FITS header updates, needed for LZ partial field of view. if use_fits and final and eit_fxpar(header,'naxis3') gt 0 then begin fxaddpar, header, 'WAVELNTH', eit_fxpar(header,'wavelnth',image_no=image_no) fxaddpar, header, 'CFTEMP', eit_fxpar(header,'cftemp',image_no=image_no) fxaddpar, header, 'DATE_OBS', eit_fxpar(header,'date_obs',image_no=image_no) fxaddpar, header, 'FILTER', eit_fxpar(header,'filter',image_no=image_no) fxaddpar, header, 'EXPTIME', eit_fxpar(header,'exptime',image_no=image_no) up_corr = where(strpos(header,'COMMENT CORREC') ne -1) if up_corr[0] ne -1 then header[up_corr] = 'COMMENT '+$ 'CORRECTED DATE_OBS = '+ anytim2utc(eit_fxpar(header,$ 'corrected date_obs',image_no=image_no),/ccsds) endif ; Pointing adjustment - 3/30/98 JSN if version lt 2.6 and use_fits then begin fxaddpar, header, 'CDELT1', eit_pixsize() fxaddpar, header, 'CDELT2', eit_pixsize() xy = eit_point(utc_date_obs) fxaddpar, header, 'CRPIX1', xy[0] fxaddpar, header, 'CRPIX2', xy[1] endif ; Set up the output file name: dat_ob=anytim2cal(eit_fxpar(header,'date_obs',image_no=image_no),form=8) if is_string(outfile) then outname=outfile else $ outname = 'eit_l1_'+strmid(dat_ob,0,8)+'_'+strmid(dat_ob,8,6)+'.fits' outname = concat_dir(outdir,outname) ; Write the FITS file if keyword_set(outfits) then begin if use_fits then writefits, outname, image, header else begin hdr = struct2fitshead(header) writefits,outname,image,hdr delvarx,hdr endelse endif ; Concatenate the headers and images if i eq 0 then outheader = header else begin if use_fits then begin szo = size(outheader) szh = size(header) case 1 of szh[1] lt szo[1]: header = [header,' '] szh[1] gt szo[1]: header = header[0:szo[1]-1] else: endcase outheader = [[outheader],[header]] endif else begin out_str = outheader[0] outheader = [str_copy_tags(out_str,outheader), str_copy_tags(out_str,header)] endelse endelse if i eq 0 and yes_data then begin ss = size(image) if keyword_set(no_calibrate) then outimage = intarr(ss[1],ss[2],nfiles) else $ outimage = fltarr(ss[1],ss[2],nfiles) endif vv[i] = max(image) if yes_data then outimage[0,0,i] = temporary(image) i = i + 1 endrep until i eq nfiles ; MAIN LOOP ends if ~keyword_set(save_zero) then begin ; Reject blank images jj = where(vv ne 0,ncount) if ncount eq 0 then begin message,/info,'All images are zero' outimage = 0 outheader = '' endif else begin if yes_data and ncount ne nfiles then outimage = outimage[*,*,jj] if use_fits then outheader = outheader[*,jj] else $ outheader = outheader[jj] n_block = n_block[jj] if n_elements(jj) eq 1 then n_block = n_block[0] endelse endif ;-- correct for 180 degree roll ;-- check in case OUTHEADER is a structure. If so, we handle it ; differently. if ~keyword_set(no_roll) then begin in_struct=have_tag(outheader,'SC_ROLL',rindex) for i=0,nfiles-1 do begin if in_struct then begin header=outheader[i] sc_roll=header.(rindex) endif else begin header=outheader[*,i] sc_roll = fxpar(header, 'SC_ROLL') endelse if abs(fix(sc_roll)) eq 180 then begin outimage[*,*,i]=rotate(outimage[*,*,i],2) if in_struct then outheader[i]=rot_fits_head(header) else $ outheader[*,i]=rot_fits_head(header) endif endfor endif message,/info,string('took',systime(1)-t0,' seconds',format='(a,f5.1,a)') end