Viewing contents of file '../idllib/ghrs/pro/calhrs.pro'
;**************************************************************************
;+
;*NAME:
; calhrs
;
;*PURPOSE:
; Reduce (calibrate) GHRS science data.
;
;*PARAMETERS:
; CALLING SEQUENCE:
; calhrs,name,params,readouts,wave,flux,eps,err,headers,log
;
; INPUTS:
; name - name of input file without qualifier
; or the id number of the file.
; OPTIONAL INPUTS:
; params - parameter description.
; 0 - use all defaults
; 1 - interative selection of parameters
; string - parameter defintion of the form
; 'parname1=value1,parname2=value2,,,,'
; string array - each element of the array has
; form 'parname=value'
; string of form '$filename' where filename is the
; name of the file containing one parameter
; per line in the form;
; parname=value
; Any parameter not specified is set to its default.
; Defaults for the parameters can be found by using interactive
; selection params=1 or examinining the default (text) file,
; ZDEF:CALHRS.DEF.
;
; The following parameters are availble.
;
; OUTFILE Write results to output science file
; (0-no 1-yes)
; INDENT Output file identification (string, output file
; name is <input file name><indent>.SCI
; ERRORS propagate statistical errors (0-no 1-yes)
; LOG write processing log (0-no, 1-yes)
; BCK_MED Background, median filter width
; BCK_MEAN Background, mean filter width
; SKY_MED Sky mean filter width
; SKY_MEAN Sky median filter width
; DEFAULTS name of text file containing default
; reference file names
; DCOFF name of dispersion coefficient table
; used to compute offset from DCTAB file.
; The following parameters control the processing steps that
; are preformed and the intermediate output written. Each
; parameter is a 2 digit integer where the first digit specifies
; the value for accum. mode and the second digit specifies the
; value for direct-downlink mode. Values for the digits are
; assigned as follows:
; 0 - do not perform the step
; 1 - perform step, don't write output from the step
; 2 - perform step, write data
; 3 - perform step, write data and epsilon
; 4 - perform step, write data and errors
; 5 - perform step, write data, epsilon, and errors
; 6 - perform step, write data, epsilon, and errors,
; and do not perform any following steps.
;
; MAP Perform mapping function
; DQI data quality initialization
; EXP Convert to count rates
; PPC Perform paired pulse correction
; DIO Diode non-uniformity correction
; DOP Correct for doppler compensation
; PHC Correct for photocathode granularity
; VIG Correct for vignetting
; MER merge substep bins
; BCK subtract background
; ADC Apply dispersion coeffients
; THM Correct for thermal motion
; IAC apply incidence angle correction
; RIP perform echelle ripple correction
; ABS absolute sensitivity correction
; HEL helocentric wavelengths
; AIR convert to air wavelengths
;
; The following parameters control which reference file to
; use. There default value is DEF specifying the default
; reference file. Default reference file names are stored in
; the text file ZCAL:DEFAULTS.TXT
;
; LMAPTAB line mapping function coefficient table
; SMAPTAB sample mapping function coefficient table
; CCTAB carrousel calibration table
; DCTAB dispersion coefficient table
; IATAB incidence angle table (SSA to LSA)
; SCLTAB offsets to spectral calibration lamp apertures
; GRATTAB echelle grating parameter table
; RIPTAB echelle ripple coefficient table
; DQIFILE data quality initialization file
; DIOFILE diode response file name
; PHCFILE photocathode response file name
; VIGFILE vignetting response file name
; ABSFILE absolute sensitivity file name
; WAVFILE wavelengths for ABSFILE
; BLEMTAB contains photocathode blemish locations
;
; readouts - vector of scalar giving the readout numbers to
; process. If set to 0, the last readout is processed.
; If not supplied, or set to 'ALL', all readouts are
; processed.
;
; OPTIONAL OUTPUTS:
; wave - two dimensional wavelength array. wave(*,i) contain the
; wavelengths for readout i+1.
; flux - two dimensional flux array (same form as wave)
; eps - two dimensional data quality array
; err - two dimensional statistical error array
; headers - two dimensional data header arrays. Size=(128xNREADOUTS)
; log - processing log (string array in internal FITS header format)
;
; OPERATIONAL NOTES:
; If you want to like to watch the progress of the reduction set
; !DUMP to 1 or greater. Text output will then be printed as
; controlled by !TEXTOUT.
;
; EXAMPLES:
; Reduce observation 20000 with all defaults, results will be
; placed in H20000.SCI and H20000.PLH (processing log).
; CALHRS,20000
;
; Reduce the same observation but stop after background subtraction,
; and use MYRESP.HHH for the diode response reference file.
; CALHRS,20000,'BCK=66,DIOFILE=MYFILE'
;
; Reduce observation file Z01Y0201.SCI with interactive parameter
; selection. return wavelengths and flux in idl variables.
; CALHRS,'Z01Y0201',1,'ALL',WAVE,FLUX
;
; Restore observation 22090, readout 10, with no file output and
; results returned in IDL variables.
; CALHRS,22090,'OUTFILE=0,LOG=0',10,WAVE,FLUX,EPS,ERR,IH,LOG
;
; Process 50 observations, specify text output, and place the
; text output for all 50 observations in file TODAYS.PRT
; !textout=3 Specify file output
; textopen,'today' Open textfile
; !textout=5 Inhibit CALHRS from opening
; new textfile for each obs.
; for i=20000,20049 do calhrs,i Reduce them
; !textout=3 This allows textclose to close
; the file
; textclose
;
; HISTORY:
; version 1, D. Lindler April 1989
; version 1.1, D. Lindler Feb 6. 1991
; added DCOFF input parameter and FPSPLIT disp. coef.
; offset computation in CALHRS_WAV. implemented
; disp. coef. thermal motion. and cubic dispersion
; coef term.
; version 1.2 D. Lindler Mar 1991 Added background scale factors
; BACK_A,...,BACK_D and polynomial fit to compute
; a smooth background. Switched order of paired
; pulse and diode response correction. Enhanced
; epsilon computation
; 9-SEP-1991 JKF/ACC - added COADD_DDL keyword
; parameter to coadd direct downlink obs.
; prior to calibration.
; 5-NOV-1991 JKF/ACC - changed COADD_DDL logical per DJL.
; Version 1.4 10-Jun-1992 DJL/ACC - added photocathode blemish flagging.
; Allows use of global coefficients for computing
; dispersion coefficients for arbituary carrousel
; positions.
; version 1.5 added 2 temperature thermal model and time model.
; version 1.6 added ability to select default reference files by
; time, grating mode, detector, and aperture
; version 1.7 added GIMP and Echelle scattered light defaults
; table (DJL, Feb 24, 1993)
;-
;-----------------------------------------------------------------------------
pro calhrs,name,params,readouts,wave,flux,eps,err,headers,log
VERSION = 1.7 ;version number
on_error,2
if n_params(0) lt 1 then $
message,'Calling Sequence: ' + $
'calhrs,name,params,readouts,wave,flux,eps,err,headers,log'
;
; get default parameters
;
if n_params(0) lt 2 then params=0 ;all defaults
if n_params(0) lt 3 then readouts='ALL' ;all readouts
getpar,params,'calhrs',par
;
; open input file and get observation mode info ------------------------------
; grating is a number 0-7 0 (not grating)
; sdtype is 'ACC' -accumulation mode
; 'DDL' - direct downlink
; 'MAP' - field maps
; 'NONE' - no science data in file
hrs_open,name,infile,'read','',log
filename=strtrim(string(byte(infile(10:99))))
hrs_mode,infile,detector,grating,sdtype
if sdtype eq 'NONE' then begin
print,'CALHRS-- ERROR No Science data found in file'+filename
goto,done
retall
endif
if sdtype eq 'MAP' then begin
print,'CALHRS-- Input observation contains field maps'
print,'CALHRS-- No processing performed'
write_log=0
goto,done
endif
hist='CALHRS version'+string(version,'(F5.2)')+' '+!stime+' '+ $
'ID = '+strtrim(name,2)
if !dump gt 0 then begin
textopen,'calhrs'
printf,!textunit,' '
end
sxaddhist,hist,log
if !dump gt 0 then printf,!textunit,hist
;
; read UDL's and SHP's
;
hrs_read,infile,'shp',ihshp
hrs_read,infile,'udl',ihudl,udls
if !err lt 0 then goto,done
;
; determine target aperture
;
sclamp = udls(2)
binid1 = udls(25)
if sclamp gt 0 then begin
aperture = 'SC'+strtrim(sclamp,2)
end else begin
if binid1 eq 2 then aperture = 'LSA' else aperture = 'SSA'
end
;
; determine observation date
;
obsdate = string(byte(ihshp(54:65),0,23))
;
; get default reference files for ones not supplied by user
;
calhrs_def,par(4),detector,grating,obsdate,aperture,par,log
;
; extract information from keyword parameter vector
;
outflag=fix(par(0))
ident=strtrim(par(1),2)
error_flag=fix(par(2))
write_log=fix(par(3))
fwidths=intarr(6)
for i=0,5 do fwidths(i)=fix(par(i+11))
dc_offtab = strtrim(par(20),2)
bck_scale = par(81:84)
pflags=intarr(17)
for i=0,16 do pflags(i)=fix(par(30+i))
coadd_ddl = fix( par(85) ) ; coadd DDL's flag
gimp = [float(par(9)),float(par(10))]
;
; Direct Downlink Obs.
;
if sdtype eq 'DDL' then pflags=pflags mod 10 else pflags=pflags/10
last=where(pflags gt 5) ;last step
if !err gt 0 then begin ;do not do anything past it
last=last(0)
if last lt 16 then pflags(last+1:16)=0
endif
;
; read data and open output data set ----------------------------------------
;
;
; Determine science records to read for direct down-link data
; ...ingore the first 2 and last 1 (Hysterisis and pass deflections)
;
recnums = readouts
if sdtype eq 'DDL' then begin
rectypes=infile(100:n_elements(infile)-1)
number=long(total(rectypes eq 10)) ;number of reads
sz = size(recnums) & ndim = sz(0)
if ndim gt 0 then begin ;vector supplied?
good = where((recnums gt 2) and (recnums lt number),ngood)
if ngood lt n_elements(recnums) then begin
message,' Error: You asked for an invalid' + $
' direct downlink readout number',/continue
message,' Readout numbers should be from 3 to the'+$
' number of readouts minus 1',/continue
message,' First 2 and last readout are not science data'
endif
end else begin ;scalar supplied
if strupcase(strtrim(recnums)) eq 'ALL' then $
recnums = indgen(number-3) + 3
end
end
hrs_read,infile,'raw',ih,data,recnums
if !err lt 0 then goto,done
;
; DDL mode, user can request to have readouts coadded before calibration...
; i.e. coadd_ddl=1 - coadd DDL readouts before calibration. This
; feature helps with virtual memory required in later steps since
; there can be thousands of readouts in DDL mode.
;
if sdtype eq 'DDL' then begin
if coadd_ddl then begin
data = sum(data,1)
ih = ih(*,0)
ih(46) = number-3
end
end
data=float(data)
if !prelaunch then begin
epsilon=data*0 ;no raw data quality vectors
end else begin
hrs_read,infile,'raw/eps',iheps,epsilon,readouts
if !err lt 0 then epsilon = data*0
fill = where(epsilon eq 16) & n=!err ;fill data
if n gt 0 then epsilon(fill) = 800
bad = where(epsilon eq 1) & n=!err ;reed-solomon error
if n gt 0 then epsilon(bad) = 100
end
if sdtype eq 'DDL' then $
if coadd_ddl then epsilon = epsilon(*,0)
if error_flag then errors=sqrt(data)>1 else errors=0
if outflag then hrs_open,name,outfile,'W',ident else pflags=pflags<1
;
; Perform reduction ----------------------------------------------------------
;
; extract some needed info from the unique data logs
;
calhrs_udl,ihudl,ih,udls,u1,u2 ;match sdps to udls
calhrs_xdef,udls,u1,u2,ih,log ;determine init. x-deflections
;
; determine number of bins to merge
;
calhrs_sort,udls,ih,first_bin,nbins,nfound,nmerge,log
if sdtype eq 'DDL' then nmerge=0 ;don't merge direct downlink
;
; mapping function with optional GIMP correction
;
if pflags(0) gt 0 then begin
calhrs_map,par(50:51),ih,log
if gimp(detector-1) ne 0.0 then $
calhrs_gimp,gimp,ihudl,udls,u1,u2,sdtype,ih,log
endif
calhrs_out,outfile,'raw',pflags(0),ih,data,epsilon,errors
;
; data quality initialization and blemish flagging
;
if (pflags(1) gt 0) and (strupcase(strtrim(par(60))) ne 'NONE') then $
calhrs_dqi,par(60),udls,epsilon,log
if (pflags(1) gt 0) and (pflags(0) gt 0) and $
(strupcase(strtrim(par(66))) ne 'NONE') then $
calhrs_blem,par(66),ih,epsilon,log
calhrs_out,outfile,'raw',pflags(1),ih,data,epsilon,errors
;
; conversion to count rates
;
if pflags(2) gt 0 then calhrs_exp,ih,data,errors,log
calhrs_out,outfile,'exp',pflags(2),ih,data,epsilon,errors
;
; diode response
;
if pflags(3) gt 0 then $
calhrs_dio,par(61),udls,ih,data,epsilon,errors,log
calhrs_out,outfile,'dio',pflags(3),ih,data,epsilon,errors
;
; paired pulse correction
;
if pflags(4) gt 0 then calhrs_ppc,ih,data,epsilon,errors,log
calhrs_out,outfile,'ppc',pflags(4),ih,data,epsilon,errors
;
; doppler correction
;
if pflags(5) gt 0 then begin
print,'Doppler compensation processing not yet implemented'
; calhrs_dop,ih,ihshp,ihudl,udls,first_bin, $
; nbins,nfound,dop_offsets,log
end else begin
dop_offsets=0
end
;
; photocathode response
;
if pflags(6) gt 0 then $
print,'Photocathode granularity removal not yet available'
; calhrs_phc,par(62),first_bin,dop_offsets, $
; ih,data,epsilon,errors,log
;;; calhrs_out,outfile,'phc',pflags(6),ih,data,epsilon,errors
;
; vignetting response
;
if pflags(7) gt 0 then $
calhrs_vig,par(63),first_bin,dop_offsets, $
ih,data,epsilon,errors,log
calhrs_out,outfile,'vig',pflags(7),ih,data,epsilon,errors
;
; merge spectral data
;
if pflags(8) eq 0 then begin
nmerge=0
n=size(ih) & n=n_elements(ih)/s(1)
first_bin=indgen(n) & nfound=replicate(1,n)
endif
calhrs_mer,ih,data,epsilon,errors, $
first_bin,nbins,nfound,nmerge, $
headers,flux,eps,err,log
calhrs_out,outfile,'mer',pflags(8),headers,flux,eps,err
;
; wavelength assignments and thermal motion correction
;
if pflags(9) gt 0 then $
calhrs_wav,pflags(10),par(52:54),dc_offtab, $
headers,wave,log
calhrs_out,outfile,'adc',pflags(9)>pflags(10),headers,wave
;
; incidence angle correction
;
if pflags(11) gt 0 then $
calhrs_iac,par(55),headers,wave,log
calhrs_out,outfile,'iac',pflags(11),headers,wave
;
; background substraction
;
if pflags(12) gt 0 then calhrs_bck,ih,data,epsilon,fwidths, $
par(54),bck_scale,first_bin,nfound,nmerge,headers,flux, $
ihb,back,epsb,log
calhrs_out,outfile,'bck',pflags(12),ihb,back,epsb,errb
calhrs_out,outfile,'net',pflags(12),headers,flux,eps,err
;
; echelle ripple
;
if pflags(13) gt 0 then calhrs_rip,par(57:58),headers,flux,err,log
calhrs_out,outfile,'rip',pflags(13),headers,flux,eps,err
;
; absolute sensitivity
;
if pflags(14) gt 0 then $
calhrs_abs,par(64:65),headers,wave,flux,err,log
calhrs_out,outfile,'abs',pflags(14),headers,flux,eps,err
;
; heliocentric correction
;
if pflags(15) gt 0 then calhrs_hel,headers,wave,log
if !err gt 0 then calhrs_out,outfile,'hel',pflags(15),headers,wave
;
; vacumm to air
;
if pflags(16) gt 0 then calhrs_air,wave,headers,log
calhrs_out,outfile,'air',pflags(16),headers,wave
;
; done
;
done:
textclose
if write_log then $
sxhwrite,strtrim(byte(infile(10:99)))+ident+'.plh',log
hrs_close,infile
if n_elements(outfile) gt 100 then hrs_close,outfile
return
end