bHROS Data Reduction

Elements of bHROS data processing that will be incorporated into the Gemini IRAF package are not yet available.

In the meantime, the following processing recipe was developed by Michael Murphy, Antonio Hales and the bHROS Demonstration Science team and appears to work well in the cases tested so far. As this is still a work in progress, users should be critical of the recipe steps and results.

"The most difficult issue in the data reduction is the fact that lines of constant wavelength across an echelle order are not parallel to the CCD rows. Furthermore, this "tilt" of the arc lines is different from order to order. It appears that the best way to solve this is to use an arc lamp exposure to map the tilt of the arc lines and transform the image to one in which the arc lines run parallel to the rows. For 1D spectra, this can be done using a combination of the "fitcoords" and "transform" packages in IRAF. But this cannot easily be applied as it stands to 2D echelle spectra (i.e. when several orders are present and a 2D wavelength solution is required). If one simply extracted the orders as they are, a loss in resolution and wavelength precision would result. 

[We] have attempted to get around the tilted arc-lines problem by creating several (14 for the single-fibre mode, up to 10 for the dual-slicer mode) "sub-apertures" across each echelle order. A trace is defined by tracing along a single slice (of, for example, the 7 slices in each order of a "large-slicer" exposure). The trace is then translated to each sub-aperture and the data are extracted. This seems to work reasonably well. [The] attached recipe of suggested IRAF commands applies this algorithm to single-fibre exposures (... currently working on reducing dual-fibre data). One needs to make a different directory for each science exposure and make the scripts trace_create.pl and trace_verify.csh available within that directory. One also needs the attached hcopy.cl script for IRAF. The data so produced can be combined into a single 1D vacuum-heliocentric spectrum using UVES_popler, a code [Michael Murphy has] written to combine many overlapping (or not) echelle orders and exposures of VLT/UVES data. This code has been modified to read in the IRAF data products. The code and an incomplete manual are available at http://www.ast.cam.ac.uk/~mim/UVES_popler.html . The code will also be modified to perform the sky subtraction for dual-fibre data.

(This recipe is also available as a text file, reduction.txt, and the various perl, shell and IRAF scripts named above are available in gzip-compressed tarfile hrosDataRed.tar.gz)

# Suggested IRAF recipe

# Assumptions:
# badpix.dat object.fits thar.fits flat.fits bias.fits all exist
# We know the spatial binning, in this case 1x1
# We know the trimsection for the data, in this case [1:2048,1:4608]
# We know whether this is single- and dual-aperture data
# We know the number of orders, and which ones, we are going to extract
# We know which slice we are going to trace on
# We know how wide the orders are in pixels

observatory set gemini-south
# Setting package information to echelle defautls; exectue ":q", ":q"
setinstrument echelle
# Set observ=gemini-south dispaxis=2 databas=database
epar echelle

# Make sure CCDSEC, DETSEC and DATASEC all agree
hcopy input1=object.fits[1] input2=object.fits[1] keywor1=DATASEC keywor2=CCDSEC keytype=string
hcopy input1=object.fits[1] input2=object.fits[1] keywor1=DATASEC keywor2=DETSEC keytype=string
hcopy input1=flat.fits[1] input2=flat.fits[1] keywor1=DATASEC keywor2=CCDSEC keytype=string
hcopy input1=flat.fits[1] input2=flat.fits[1] keywor1=DATASEC keywor2=DETSEC keytype=string
hcopy input1=thar.fits[1] input2=thar.fits[1] keywor1=DATASEC keywor2=CCDSEC keytype=string
hcopy input1=thar.fits[1] input2=thar.fits[1] keywor1=DATASEC keywor2=DETSEC keytype=string
hcopy input1=bias.fits[0] input2=bias.fits[0] keywor1=DATASEC keywor2=CCDSEC keytype=string
hcopy input1=bias.fits[0] input2=bias.fits[0] keywor1=DATASEC keywor2=DETSEC keytype=string

# Bias subtraction
!/bin/ls object.fits thar.fits flat.fits | sed 's/.fits/.fits[1]/' > in1
!/bin/ls object.fits thar.fits flat.fits | sed 's/.fits/_b.fits/' > in2
ccdproc images=@in1 output=@in2 ccdtype=' ' fixpix+ oversca+ trim+ zerocor+ darkcor- flatcor- readaxis=line fixfile=badpix.dat biassec=image trimsec=[1:2048,1:4608] zero=bias.fits interac+

# Copy over relevant header keywords
hcopy input1=object.fits[0] input2=object_b.fits keywor1=OBSEPOCH keywor2=' ' keytype=real
hcopy input1=object.fits[0] input2=object_b.fits keywor1=OBSEPOCH keywor2=' ' keytype=real
hcopy input1=object.fits[0] input2=object_b.fits keywor1=DATE-OBS keywor2=' ' keytype=string
hcopy input1=object.fits[0] input2=object_b.fits keywor1=UTSTART keywor2=' ' keytype=string
hcopy input1=object.fits[0] input2=object_b.fits keywor1=UTSTART keywor2=' ' keytype=string
hedit images=object_b.fits fields=LAT_OBS value=-30.228333 add+ addonly- delete- verify- show+ update+
hedit images=object_b.fits fields=LONG_OBS value=70.723333 add+ addonly- delete- verify- show+ update+
hedit images=object_b.fits fields=ALT_OBS value=2737.0 add+ addonly- delete- verify- show+ update+
hcopy input1=object.fits[0] input2=object_b.fits keywor1=RA keywor2=' ' keytype=string
hcopy input1=object.fits[0] input2=object_b.fits keywor1=RA keywor2=' ' keytype=string
hcopy input1=object.fits[0] input2=object_b.fits keywor1=DEC keywor2=' ' keytype=string
hcopy input1=object.fits[0] input2=object_b.fits keywor1=DEC keywor2=' ' keytype=string
hcopy input1=object.fits[0] input2=object_b.fits keywor1=EQUINOX keywor2=' ' keytype=real
hcopy input1=object.fits[0] input2=object_b.fits keywor1=EXPTIME keywor2=' ' keytype=real

# Define flats for each sub-aperture
!ln -s flat_b.fits trace.fits
!ln -s flat_b.fits order_trace.fits

# Define an aperture covering the slice with the most coverage for
#  tracing each order and do the trace
apfind trace.fits apertur=9 referen=' ' interac+ find+ recente+ resize+ edit+ nsum=30 nfind=9 minsep=200.0 maxsep=400.0 order=increasing
aptrace trace.fits interac+ find- recente- resize- edit- trace+ fittrac+ nsum=30 step=30 nlost=10 functio=chebyshev order=4 niterat=4

# Redefine the apertures for tracing the entire orders and apply the
#  above trace. Note that you must try to include an extra two or three
#  pixels either side of the order.
apfind order_trace.fits referen=trace.fits apertur=9 interac+ find- recente- resize- edit+ nsum=30 nfind=9 order=increasing

# Flatten the flatfield.
apflatten flat_b.fits flat_bn.fits referen=order_trace.fits interac+ find- recente- resize- edit- trace- fittrac- nsum=30 thresho=20.0 readnoi=rdnoise gain=gain functio=chebyshev order=11 niterat=3

!ln -s flat_b.fits flat_b_01.fits
!ln -s flat_b.fits flat_b_02.fits
!ln -s flat_b.fits flat_b_03.fits
!ln -s flat_b.fits flat_b_04.fits
!ln -s flat_b.fits flat_b_05.fits
!ln -s flat_b.fits flat_b_06.fits
!ln -s flat_b.fits flat_b_07.fits
!ln -s flat_b.fits flat_b_08.fits
!ln -s flat_b.fits flat_b_09.fits
!ln -s flat_b.fits flat_b_10.fits
!ln -s flat_b.fits flat_b_11.fits
!ln -s flat_b.fits flat_b_12.fits
!ln -s flat_b.fits flat_b_13.fits
!ln -s flat_b.fits flat_b_14.fits

# Create shifted traces for sub-apertures: TO BE EXECUTED IN SHELL, NOT CL
cd database; ../trace_create.pl trace flat_b 5 7 201 1 ; cd ..

# Flatfield the object and thar files
!echo "object_b.fits\nthar_b.fits" > in1
!sed 's/_b/_bn/' in1 > in2
ccdproc images=@in1 output=@in2 fixpix+ oversca- trim- zerocor- flatcor+ fixfile=badpix.dat flat=flat_bn.fits

# Divide flattened flat into original flat to obtain blaze function
!ln -s flat_b.fits blaze_b.fits
ccdproc images=blaze_b.fits output=blaze_bn.fits fixpix+ oversca- trim- zerocor- flatcor+ fixfile=badpix.dat flat=flat_bn.fits

# Fit scattered light distribution in object and blaze exposure
# Set functio=chebyshev order=2
epar apscat1
# Set functio=chebyshev order=5
epar apscat2
apscatter input=object_bn.fits output=object_bns.fits referen=order_trace.fits interac+ find- recente- resize- edit- trace- fittrac- subtrac+ smooth+ fitscat+ fitsmoo+ nsum=50 buffer=2
apscatter input=blaze_bn.fits output=blaze_bns.fits referen=order_trace.fits interac+ find- recente- resize- edit- trace- fittrac- subtrac+ smooth+ fitscat+ fitsmoo+ nsum=50 buffer=2

!ln -s object_bns.fits object_bns_01.fits
!ln -s object_bns.fits object_bns_02.fits
!ln -s object_bns.fits object_bns_03.fits
!ln -s object_bns.fits object_bns_04.fits
!ln -s object_bns.fits object_bns_05.fits
!ln -s object_bns.fits object_bns_06.fits
!ln -s object_bns.fits object_bns_07.fits
!ln -s object_bns.fits object_bns_08.fits
!ln -s object_bns.fits object_bns_09.fits
!ln -s object_bns.fits object_bns_10.fits
!ln -s object_bns.fits object_bns_11.fits
!ln -s object_bns.fits object_bns_12.fits
!ln -s object_bns.fits object_bns_13.fits
!ln -s object_bns.fits object_bns_14.fits
!ln -s blaze_bns.fits blaze_bns_01.fits
!ln -s blaze_bns.fits blaze_bns_02.fits
!ln -s blaze_bns.fits blaze_bns_03.fits
!ln -s blaze_bns.fits blaze_bns_04.fits
!ln -s blaze_bns.fits blaze_bns_05.fits
!ln -s blaze_bns.fits blaze_bns_06.fits
!ln -s blaze_bns.fits blaze_bns_07.fits
!ln -s blaze_bns.fits blaze_bns_08.fits
!ln -s blaze_bns.fits blaze_bns_09.fits
!ln -s blaze_bns.fits blaze_bns_10.fits
!ln -s blaze_bns.fits blaze_bns_11.fits
!ln -s blaze_bns.fits blaze_bns_12.fits
!ln -s blaze_bns.fits blaze_bns_13.fits
!ln -s blaze_bns.fits blaze_bns_14.fits

# Extract object spectra from each subaperture
!/bin/ls object_bns_*.fits > in1
!sed 's/_bns/_bnse/' in1 > in2
!/bin/ls flat_b_*.fits > in3
apall input=@in1 output=@in2 format=echelle referen=@in3 interac- find- recente- resize- edit- trace- fittrac- extract+ extras+ review- nsum=5 backgro=none weights=variance clean+ saturat=30000.0 readnoi=rdnoise gain=gain lsigma=4.0 usigma=4.0 nsubaps=1

# Extract blaze spectra from each subaperture
!/bin/ls blaze_bns_*.fits > in1
!sed 's/_bns/_bnse/' in1 > in2
!/bin/ls flat_b_*.fits > in3
apall input=@in1 output=@in2 format=echelle referen=@in3 interac- find- recente- resize- edit- trace- fittrac- extract+ extras- review- nsum=5 backgro=none weights=variance clean+ saturat=30000.0 readnoi=rdnoise gain=gain lsigma=4.0 usigma=4.0 nsubaps=1

# Normalize the blaze spectra by a fitting a constant value
!/bin/ls blaze_bnse_*.fits > in1
continuum input=@in1 output=@in1 type=ratio replace- wavesca- logscal- overrid- listonl- interac- naverag=1 function=chebyshev order=1 low_rej=3.0 high_re=3.0 niterat=5 grow=0.0 markrej-

# Split the flux and error arrays of object spectra
!/bin/ls object_bnse_*.fits > in1
!sed 's/object_bnse/flx/' in1 > in2
!sed 's/object_bnse/sig/' in1 > in3
scopy input=@in1 output=@in2 bands=1 renumbe- offset=0 merge- rebin-
scopy input=@in1 output=@in3 bands=3 renumbe- offset=0 merge- rebin-

# Blaze-correct flux and sigma spectra
!/bin/ls flx_*.fits > in1
!/bin/ls sig_*.fits > in2
!/bin/ls blaze_bnse_*.fits > in3
sarith input1=@in1 op=/ input2=@in3 output=@in1 rebin- clobber+
sarith input1=@in2 op=/ input2=@in3 output=@in2 rebin- clobber+

noao.imred.echelle
# Set nsum=9 extras- t_step=30 t_funct=chebyshev t_order=4 t_niter=5
#  weights=variance thresho=40.0 coordli=linelists$thar.dat match=0.02
#  fwidth=5.0 cradius=9.0 i_funct=chebyshev i_xorde=4 i_yorde=4
#  i_niter=5 select=match sort=none group=none refit+ lineari- log-
#  flux+
epar sparams

!ln -s thar_bn.fits thar_bn_01.fits
!ln -s thar_bn.fits thar_bn_02.fits
!ln -s thar_bn.fits thar_bn_03.fits
!ln -s thar_bn.fits thar_bn_04.fits
!ln -s thar_bn.fits thar_bn_05.fits
!ln -s thar_bn.fits thar_bn_06.fits
!ln -s thar_bn.fits thar_bn_07.fits
!ln -s thar_bn.fits thar_bn_08.fits
!ln -s thar_bn.fits thar_bn_09.fits
!ln -s thar_bn.fits thar_bn_10.fits
!ln -s thar_bn.fits thar_bn_11.fits
!ln -s thar_bn.fits thar_bn_12.fits
!ln -s thar_bn.fits thar_bn_13.fits
!ln -s thar_bn.fits thar_bn_14.fits

# Extract and identify ThAr spectra of subaperture with best coverage
apall input=thar_bn_09.fits output=wav_09.fits format=echelle referen=flat_b_09.fits interac- find- recente- resize- edit- trace- fittrac- extract+ extras+ review- nsum=5 backgro=none weights=variance clean+ saturat=30000.0 readnoi=rdnoise gain=gain lsigma=3.0 usigma=3.0 nsubaps=1
ecidentify images=wav_09.fits databas=database coordli=linelists$thar.dat match=0.02 maxfeat=100 zwidth=10.0 ftype=emission fwidth=5.0 cradius=9.0 thresho=80.0 minsep=4.0 functio=chebyshev xorder=3 yorder=3 niterat=5 lowreje=3.0 highrej=3.0 autowri-

# Extract and identify ThAr spectra for other apertures using previous
#  solution as first guess.
!sed 's/wav_09/wav_01/' database/ecwav_09 > database/ecwav_01
!sed 's/wav_09/wav_02/' database/ecwav_09 > database/ecwav_02
!sed 's/wav_09/wav_03/' database/ecwav_09 > database/ecwav_03
!sed 's/wav_09/wav_04/' database/ecwav_09 > database/ecwav_04
!sed 's/wav_09/wav_05/' database/ecwav_09 > database/ecwav_05
!sed 's/wav_09/wav_06/' database/ecwav_09 > database/ecwav_06
!sed 's/wav_09/wav_07/' database/ecwav_09 > database/ecwav_07
!sed 's/wav_09/wav_08/' database/ecwav_09 > database/ecwav_08
!sed 's/wav_09/wav_10/' database/ecwav_09 > database/ecwav_10
!sed 's/wav_09/wav_11/' database/ecwav_09 > database/ecwav_11
!sed 's/wav_09/wav_12/' database/ecwav_09 > database/ecwav_12
!sed 's/wav_09/wav_13/' database/ecwav_09 > database/ecwav_13
!sed 's/wav_09/wav_14/' database/ecwav_09 > database/ecwav_14
!/bin/ls thar_bn_*.fits | grep -v "09" > in1
!sed 's/thar_bn/wav/' in1 > in2
!/bin/ls flat_b_*.fits | grep -v "09" > in3
apall input=@in1 output=@in2 format=echelle referen=@in3 interac- find- recente- resize- edit- trace- fittrac- extract+ extras+ review- nsum=5 backgro=none weights=variance clean+ saturat=30000.0 readnoi=rdnoise gain=gain lsigma=3.0 usigma=3.0 nsubaps=1
ecidentify images=@in2 databas=database coordli=linelists$thar.dat match=0.02 maxfeat=100 zwidth=10.0 ftype=emission fwidth=5.0 cradius=9.0 thresho=80.0 minsep=4.0 functio=chebyshev xorder=5 yorder=5 niterat=5 lowreje=3.0 highrej=3.0 autowri-

# Attach wavelength solutions to corresponding flux and sigma spectra
!/bin/ls flx_*.fits > in1
!/bin/ls sig_*.fits > in2
!/bin/ls wav_*.fits > in3
refspectra input=@in1 referen=@in3 refaps=' ' ignorea+ select=match sort=none group=none time=no overrid+ confirm- assign+
dispcor input=@in1 output=@in1 lineari=no databas=database log=no flux=yes blank=0.0 samedis- global- ignorea- confirm- listonl- verbose-
refspectra input=@in2 referen=@in3 refaps=' ' ignorea+ select=match sort=none group=none time=no overrid+ confirm- assign+
dispcor input=@in2 output=@in2 lineari=no databas=database log=no flux=yes blank=0.0 samedis- global- ignorea- confirm- listonl- verbose-

# UVES_popler hd140283_red.list -filetype 1 -disp 0.7 -ordsig 13.0 -ordsignbr 8.0 -ordsigzero 11.0 -clipsig 9.0 -scalclip 9.0 -pctlred 0.05 -rsigredl 7.0 -rsigredu 9.0 -vcred 3500.0

Last update August 30, 2006; Steven Margheim