Reducing NORRIS Spectrograph Data
Todd Small, Stephen Odewahn & Roy Gal

This document contains a detailed example run-through of data reduction for a typical set of images obtained on a single setup of the NORRIS multifiber spectrograph. The data reduction pathway was developed by Todd Small. Detailed documentation of the step-by-step process was created by Stephen Odewahn with help from Todd Small and Roy Gal. This web page was adapted from that documentation by Roy Gal. If you have difficulty with any of the steps outlined here, you can try to contact us:
Roy Galrrg at pha.jhu.edu
Todd Small tas at astro.caltech.edu

We use standard IRAF routines, along with some IRAF scripts written by Todd Small. These scripts can be retrieved in a single gzipped tar file here.

The IRAF setup for this reduction requires some additons to the generic login.cl. A file with these additions, called norris.cl, is included in the above tar file. You will also need to define the location of the special reduction scripts in your login.cl:

scrpto = /home/galatea/sco/iraf_scripts/tas  (for instance)
Finally, we recommend you use xgterm when running IRAF, as it has a nicer graphics window.

We use the following file naming conventions:
.msMultispec format image
.scatImage corrected for scattered light
.respFiber-to-fiber response correction image

Norris images have the spectra dispersed in the Y-direction (spectra dispersed along a column). The fiber spectra are stacked in the X-direction (each spectrum intersects a line).


We are ready to go! Start up IRAF. Set up the instrument-specific parameters and the ccd reduction parameters:

cl> epar setinstrument
   instrument = "norris"        Instrument ID (type ? for a list)
        query = "q"             Instrument ID (type q to quit)
        (site = "")             Site ID
   (directory = "scrpto$")      Instrument directory
      (review = yes)            Review instrument parameters?
        (mode = "ql")           


cl> epar ccdred

PACKAGE = imred
TASK = ccdred

(pixelty=                 real) Output and calculation pixel datatypes
(verbose=                  yes) Print log information to the standard output?
(logfile=              logfile) Text log file
(plotfil=             plotfile) Log metacode plot file
(backup =                     ) Backup directory or prefix
(instrum=    scrpto$norris.dat) CCD instrument file
(ssfile =                     ) Subset translation file
(graphic=             stdgraph) Interactive graphics output device
(cursor =                     ) Graphics cursor input
(version=      2: October 1987)
(mode   =                   ql)
Now we need to define which types of data are in each of our images (i.e. flat, arc, etc.). Make a file called list.1 with the following contents:
  dome flat
  object
  arc
  internal flat

Define and execute the settype task. You will have to enter the image type for each file:

cl> task settype = scrpto$settype.cl
cl>  settype
List of images: *.fits 
This allows us to list our data files using ccdlist properly. Below we show the files used in our example run:
cl> ccdlist *.fits
s101n.fits[2048,2042][real][comp][][OT]:FE+AR ARC
s102n.fits[2048,2042][real][object][][OT]:479_1_3
s103n.fits[2048,2042][real][comp][][OT]:FE+AR
s104n.fits[2048,2042][real][object][][OT]:479_1_3
s105n.fits[2048,2042][real][comp][][OT]:FE+AR
s106n.fits[2048,2042][real][flat][][OT]:UNDITH IFF
s107n.fits[2048,2042][real][object][][OT]:479_1_3
s108n.fits[2048,2042][real][comp][][OT]:FE+AR
s109n.fits[2048,2042][real][object][][OT]:479_1_3
s110n.fits[2048,2042][real][comp][][OT]:FE+AR
s111n.fits[2048,2042][real][comp][][OT]:FE+AR
s123n.fits[2048,2042][real][flat][][OT]:479_1_3 UNDITH DFF
s124n.fits[2048,2042][real][flat][][OT]:479_1_3 UNDITH DFF
s125n.fits[2048,2042][real][flat][][OT]:479_1_3 UNDITH DFF
s126n.fits[2048,2042][real][flat][][OT]:479_1_3 UNDITH DFF
s127n.fits[2048,2042][real][flat][][OT]:479_1_3 UNDITH DFF
s128n.fits[2048,2042][real][flat][][OT]:479_1_3 UNDITH DFF
s129n.fits[2048,2042][real][flat][][OT]:479_1_3 UNDITH DFF
s130n.fits[2048,2042][real][flat][][OT]:479_1_3 UNDITH DFF
s131n.fits[2048,2042][real][flat][][OT]:479_1_3 UNDITH DFF
s132n.fits[2048,2042][real][flat][][OT]:479_1_3 UNDITH DFF

Now the true reduction process begins!

Step 1: bias subtract using overscan section, and trim image using ccdproc:

Use implot to determine trim and bias regions:

cl> implot s130n
You want to average lines (1000:1100 is good) and select the column region for a good bias measure, ususally [2060:2095,*]. The following commands in the plot window are handy:
:l 900 1100     line average
:c 900 1100    column average
"e" to reset BLC,TRC
"r" to redraw
"space" to read cursor position

Now process the data with ccdproc:

cl> epar ccdproc 
PACKAGE = ccdred
   TASK = ccdproc

images  =               *.fits  List of CCD images to correct
(ccdtype=                     ) CCD image type to correct
(max_cac=                    0) Maximum image caching memory (in Mbytes)
(noproc =                   no) List processing steps only?
(fixpix =                   no) Fix bad CCD lines and columns?
(oversca=                  yes) Apply overscan strip correction?
(trim   =                  yes) Trim the image?
(zerocor=                   no) Apply zero level correction?
(darkcor=                   no) Apply dark count correction?
(flatcor=                   no) Apply flat field correction?
(illumco=                   no) Apply illumination correction?
(fringec=                   no) Apply fringe correction?
(readcor=                   no) Convert zero level image to readout correction?
(scancor=                   no) Convert flat field image to scan correction?
(readaxi=                 line) Read out axis (column|line)
(fixfile=                     ) File describing the bad lines and columns
(biassec=        [2060:2095,*]) Overscan strip image section
(trimsec=      [1:2048,4:2045]) Trim data section
(zero   =                     ) Zero level calibration image
(dark   =                     ) Dark count calibration image
(flat   =                     ) Flat field images
(illum  =                     ) Illumination correction images
(fringe =                     ) Fringe correction images
(minrepl=                   1.) Minimum flat field value
(scantyp=            shortscan) Scan type (shortscan|longscan)
(nscan  =                    1) Number of short scan lines
(interac=                  yes) Fit overscan interactively?
(functio=              spline3) Fitting function
(order  =                    3) Number of polynomial terms or spline pieces
(sample =                    *) Sample points to fit
(naverag=                    1) Number of sample points to combine
(niterat=                    5) Number of rejection iterations
(low_rej=                   3.) Low sigma rejection factor
(high_re=                   3.) High sigma rejection factor
(grow   =                   1.) Rejection growing radius
(mode   =                   ql)
Fit the first few overscans interactively (yes) to make sure the fit is correct. Then enter "NO" to do the rest automatically.

Step 2: Combine the flats using flatcombine:

cl> epar flatcombine
PACKAGE = ccdred
   TASK = flatcombine

input   =  s12*.fits,s13*.fits  List of flat field images to combine
(output =             DomeFlat) Output flat field root name
(combine=               median) Type of combine operation
(reject =            avsigclip) Type of rejection
(ccdtype=                 flat) CCD image type to combine
(process=                  yes) Process images before combining?
(subsets=                  yes) Combine images by subset parameter?
(delete =                   no) Delete input images after combining?
(clobber=                   no) Clobber existing output image?
(scale  =                 mode) Image scaling
(statsec=                     ) Image section for computing statistics
(nlow   =                    1) minmax: Number of low pixels to reject
(nhigh  =                    1) minmax: Number of high pixels to reject
(nkeep  =                    1) Minimum to keep (pos) or maximum to reject (neg)
(mclip  =                  yes) Use median in sigma clipping algorithms?
(lsigma =                   3.) Lower sigma clipping factor
(hsigma =                   3.) Upper sigma clipping factor
(rdnoise=                   0.) ccdclip: CCD readout noise (electrons)
(gain   =                   1.) ccdclip: CCD gain (electrons/DN)
(snoise =                   0.) ccdclip: Sensitivity noise (fraction)
(pclip  =                 -0.5) pclip: Percentile clipping parameter
(blank  =                   1.) Value if there are no pixels
(mode   =                    h)
This produces a high S/N 2-D image of the dome flat (note that we excluded the internal dome flat = s106n.fits (UNDITH IFF).

Step 3: Prepare aperture ID's:

Use the awk script ap_id2.awk to create an aperture ID file for the IRAF task dofibers from the raw autofid2 output.
NOTE: The path to the awk script ap_id2.awk is hard-coded into the script called "ap_id2", so be sure the correct path is there.
NOTE: The ap_id2 and ap_id2.awk files must be executable (you may have to chmod them if they are not).

Get the appropriate autofid file (479_1_3.dat) NOTE: Sometimes you may find a whole set of the *.dat files. These correspond to different hour angles, but according to TAS, any of these may be for running ap_id2.

Execute the script:

cl> !/home/galatea/sco/iraf_scripts/tas/ap_id2 479_1_3.dat 479_1_3.ap
Step 4: Assign and trace apertures using aptrace:

We use the internal dome flat (s106n.fits) to do the trace. Set up the spectroscopy reduction parameters in specred:

cl> epar specred

PACKAGE = imred
   TASK = specred

(extinct= onedstds$kpnoextinct.dat) Extinction file
(caldir = onedstds$spec16redcal/) Standard star calibration directory
(observa=          observatory) Observatory of data
(interp =                poly5) Interpolation type
(dispaxi=                    2) Image axis for 2D/3D images
(nsum   =                    1) Number of lines/columns/bands to sum for 2D/3D images
(databas=             database) Database
(verbose=                   no) Verbose output?
(logfile=              logfile) Log file
(plotfil=                     ) Plot file
(records=                     ) Record number extensions
(version= SPECRED V3: April 1992)
(mode   =                    h)
($nargs =                    0)


cl> epar apdefault

PACKAGE = specred
   TASK = apdefault

(lower  =                  -5.) Lower aperture limit relative to center
(upper  =                   5.) Upper aperture limit relative to center
(apidtab=           479_1_3.ap) Aperture ID table
(b_funct=            chebyshev) Background function
(b_order=                    1) Background function order
(b_sampl=          -10:-6,6:10) Background sample regions
(b_naver=                   -3) Background average or median
(b_niter=                    0) Background rejection iterations
(b_low_r=                   3.) Background lower rejection sigma
(b_high_=                   3.) Background upper rejection sigma
(b_grow =                   0.) Background rejection growing radius
(mode   =                    h)

(use :wq to write the parameters in the uparm file) 
Now set up and execute aptrace:
cl> epar aptrace 
PACKAGE = specred
   TASK = aptrace

input   =                s106n  List of input images to trace
(apertur=                     ) Apertures
(referen=                     ) List of reference images
(interac=                  yes) Run task interactively?
(find   =                  yes) Find apertures?
(recente=                  yes) Recenter apertures?
(resize =                  yes) Resize apertures?
(edit   =                  yes) Edit apertures?
(trace  =                  yes) Trace apertures?
(fittrac=                  yes) Fit the traced points interactively?
(line   =                INDEF) Starting dispersion line
(nsum   =                   10) Number of dispersion lines to sum
(step   =                   10) Tracing step
(nlost  =                    3) Number of consecutive times profile is lost befgiving up
(functio=             legendre) Trace fitting function
(order  =                    4) Trace fitting function order
(sample =                    *) Trace sample regions
(naverag=                    1) Trace average or median
(niterat=                    5) Trace rejection iterations
(low_rej=                   3.) Trace lower rejection sigma
(high_re=                   3.) Trace upper rejection sigma
(grow   =                   0.) Trace rejection growing radius
(mode   =                    h)

   Find apertures for s214n.fits?  (yes): yes
   Number of apertures to be found automatically: 150
   Resize apertures for s214n.fits?  (yes): yes
   Edit apertures for s214n.fits?  (yes): yes
You will be confronted with a plot of the numbered apertures.
Zoom in on an aperture that we recognize. Use "w" followed by "e" at the corners of a desired region to zoom in on it.
The apertures on either side of the central blank space are 131 (left) and 43 (right).
ID one of them by placing the cursor on the peak and typing "o", followed by aperture number. The other fibers should get numbered.
Use "w" followed by "r" or "l" to move through the plot.
If there is a blank space where there should be a fiber (often the case for fiber 14 on the right):
  1. Use the "d" key to delete aperture to the right.
  2. Use "n" to force the aperture that was missing.
  3. Use "m" to mark the right-hand aperture, then "o" to assign the correct aperture number (73).
  4. Use "z" to resize the active aperture.
  5. Use "d" to delete the forced (but correctly numbered) aperture for the missing fiber.
You might also want to add a new aperture (using "n") on the far left, where there is usually half a spectrum.

IMPORTANT NOTE: It is crucial that you make sure there are no apertures left out of the ordering, and that no apertures are accidentally double-counted. A failure to do so at this point will screw everything up downstream since the aperture numbering will be wrong!!!

You can use :.snap to print the plot.
To quit use "q". Respond "yes" to tracing apertures.
Be patient, and eventually you'll be asked to fit traces interactively.
-Use "o #" to change the fit order.
-Use "f" to refit.
-Use "q" to quit the trace of current aperture.

Step 5: Subtract scattered light using apscatter:

sp> epar apscatter

PACKAGE = specred
   TASK = apscatter

input   =             DomeFlat  List of input images to subtract scattered light
output  =        DomeFlat.scat  List of output corrected images
(apertur=                     ) Apertures
(scatter=                     ) List of scattered light images (optional)
(referen=                s106n) List of aperture reference images
(interac=                  yes) Run task interactively?
(find   =                  yes) Find apertures?
(recente=                  yes) Recenter apertures?
(resize =                   no) Resize apertures?
(edit   =                   no) Edit apertures?
(trace  =                   no) Trace apertures?
(fittrac=                   no) Fit the traced points interactively?
(subtrac=                  yes) Subtract scattered light?
(smooth =                  yes) Smooth scattered light along the dispersion?
(fitscat=                  yes) Fit scattered light interactively?
(fitsmoo=                  yes) Smooth the scattered light interactively?
(line   =                INDEF) Dispersion line
(nsum   =                   10) Number of dispersion lines to sum or median
(buffer =                   1.) Buffer distance from apertures
(apscat1=                     ) Fitting parameters across the dispersion
(apscat2=                     ) Fitting parameters along the dispersion
(mode   =                    h)
You will see plots of the background values from the areas between apertures (along columns=spatial direction). You should interactively fit some of the lines (especially around line=2000):
-Use the "line ###" command to specify the line number (after a "q").
-You can use "s" to set fitting areas. Use a "t" to delete these zones. -Change the order with "o #". You might go as high as 20. -After hitting q twice, you are asked about smoothing the fit. Say "yes".
You will then be presented with plots of the scattered light along the dispersion direction (lines). Examine this at a few different columns. Hit "q" twice to finish.

The end result is a dome flat with the scattered light component removed.

Step 6: Remove scattered light from object frames with apscatter:

sp> epar apscatter

input   =          s102n,s105n  List of input images to subtract scattered light
output  =  s102n.scat,s105n.scat  List of output corrected images
(apertur=                     ) Apertures
(scatter=                     ) List of scattered light images (optional)
(referen=                s106n) List of aperture reference images
(interac=                  yes) Run task interactively?
(find   =                  yes) Find apertures?
(recente=                  yes) Recenter apertures?
(resize =                   no) Resize apertures?
(edit   =                   no) Edit apertures?
(trace  =                   no) Trace apertures?
(fittrac=                   no) Fit the traced points interactively?
(subtrac=                  yes) Subtract scattered light?
(smooth =                  yes) Smooth scattered light along the dispersion?
(fitscat=                  yes) Fit scattered light interactively?
(fitsmoo=                  yes) Smooth the scattered light interactively?
(line   =                INDEF) Dispersion line
(nsum   =                   10) Number of dispersion lines to sum or median
(buffer =                   1.) Buffer distance from apertures
(apscat1=                     ) Fitting parameters across the dispersion
(apscat2=                     ) Fitting parameters along the dispersion
(mode   =                    h)
Check each fit. Change fitting ranges with "s" and the fit order with ":o" to get better fits. Two fits are made (one for X, one for Y).

Step 7: Extract apertures in the domeflat using apall:

sp> epar apall

PACKAGE = specred
   TASK = apall
    
input   =        DomeFlat.scat) List of input images
(output =     DomeFlat.scat.ms) List of output spectra
(apertur=                     ) Apertures
(format =            multispec) Extracted spectra format
(referen=                s106n) List of aperture reference images
(profile=                     ) List of aperture profile images
(interac=                   no) Run task interactively?
(find   =                   no) Find apertures?
(recente=                  yes) Recenter apertures?
(resize =                   no) Resize apertures?
(edit   =                   no) Edit apertures?
(trace  =                   no) Trace apertures?
(fittrac=                   no) Fit the traced points interactively?
(extract=                  yes) Extract spectra?
(extras =                   no) Extract sky, sigma, etc.?
(review =                   no) Review extractions?
(line   =                INDEF) Dispersion line
(nsum   =                   10) Number of dispersion lines to sum or median

                                # DEFAULT APERTURE PARAMETERS

(lower  =                  -5.) Lower aperture limit relative to center
(upper  =                   5.) Upper aperture limit relative to center
(apidtab=                     ) Aperture ID table (optional)

                                # DEFAULT BACKGROUND PARAMETERS

(b_funct=            chebyshev) Background function
(b_order=                    1) Background function order
(b_sampl=          -10:-6,6:10) Background sample regions
(b_naver=                   -3) Background average or median
(b_niter=                    0) Background rejection iterations
(b_low_r=                   3.) Background lower rejection sigma
(b_high_=                   3.) Background upper rejection sigma
(b_grow =                   0.) Background rejection growing radius

                                # APERTURE CENTERING PARAMETERS

(width  =                   5.) Profile centering width
(radius =                  10.) Profile centering radius
(thresho=                   0.) Detection threshold for profile centering

                                # AUTOMATIC FINDING AND ORDERING PARAMETERS

nfind   =                       Number of apertures to be found automatically
(minsep =                   5.) Minimum separation between spectra
(maxsep =                1000.) Maximum separation between spectra
(order  =           increasing) Order of apertures

                                # RECENTERING PARAMETERS

(aprecen=                     ) Apertures for recentering calculation
(npeaks =                INDEF) Select brightest peaks
(shift  =                  yes) Use average shift instead of recentering?

                                # RESIZING PARAMETERS

(llimit =                INDEF) Lower aperture limit relative to center
(ulimit =                INDEF) Upper aperture limit relative to center
(ylevel =                  0.1) Fraction of peak or intensity for automatic widt
(peak   =                  yes) Is ylevel a fraction of the peak?
(bkg    =                  yes) Subtract background in automatic width?
(r_grow =                   0.) Grow limits by this factor
(avglimi=                   no) Average limits over all apertures?

                                # TRACING PARAMETERS

(t_nsum =                   10) Number of dispersion lines to sum
(t_step =                   10) Tracing step
(t_nlost=                    3) Number of consecutive times profile is lost befo
(t_funct=             legendre) Trace fitting function
(t_order=                    2) Trace fitting function order
(t_sampl=                    *) Trace sample regions
(t_naver=                    1) Trace average or median
(t_niter=                    0) Trace rejection iterations
(t_low_r=                   3.) Trace lower rejection sigma
(t_high_=                   3.) Trace upper rejection sigma
(t_grow =                   0.) Trace rejection growing radius

                                # EXTRACTION PARAMETERS

(backgro=                 none) Background to subtract
(skybox =                    1) Box car smoothing length for sky
(weights=             variance) Extraction weights (none|variance)
(pfit   =                fit1d) Profile fitting type (fit1d|fit2d)
(clean  =                   no) Detect and replace bad pixels?
(saturat=               28000.) Saturation level
(readnoi=                  6.3) Read out noise sigma (photons)
(gain   =                  1.6) Photon gain (photons/data number)
(lsigma =                   3.) Lower rejection threshold
(usigma =                   3.) Upper rejection threshold
(nsubaps=                    1) Number of subapertures per aperture
(mode   =                    h)
The resultant extracted spectra will be used to make throughput corrections for each aperture. This will be a 1d wavelength dependent throughput function for each fiber.

Step 8: Extract arc spectra using the same apall parameters.

Step 9: Find apertures that have extremely low throughput with lapertures. These will be excluded when calculating the mean of all the fiber throughputs.

sp> task lapertures = scrpto$lapertures.e
sp> lapertures DomeFlat.scat.ms.fits >junk
sp> !nawk '{print $2}' junk | statstrm h c a s r 
      Count    Average    Std Dev    Minimum    Maximum
        148    17049.7    5909.23   -6.01271    25812.9
sp> !nawk '$2/17049.7 < 0.04 {print $0}' junk
Sourcing skicat.init...
155 -5.1936121
108 -6.0127134
67 -4.6270208
19 -4.5255575
5 -4.9722753
82 -4.6096549
So, the bad fibers are: 155,108,67,19,5,82

Step 10: Make a temporary dome flat image with bad fibers deleted:

sp> epar scopy
   PACKAGE = specred
   TASK = scopy

input   =     DomeFlat.scat.ms  List of input spectra
output  =             temp1.ms  List of output spectra
(w1     =                INDEF) Starting wavelength
(w2     =                INDEF) Ending wavelength
(apertur=  !155,108,67,19,5,82) List of input apertures or columns/lines
(bands  =                     ) List of input bands or lines/bands
(beams  =                     ) List of beams or echelle orders
(apmodul=                    0) Input aperture modulus (0=none)
(format =            multispec) Output spectra format
(renumbe=                   no) Renumber output apertures?
(offset =                    0) Output aperture number offset
(clobber=                   no) Modify existing output images?
(merge  =                   no) Merge with existing output images?
(rebin  =                  yes) Rebin to exact wavelength region?
(verbose=                  yes) Print operations?
(mode   =                   ql)
The "!" in the aperture list indicates that these apertures are NOT copied.

Step 11: Generate a median domeflat from all good fibers using blkavg:

sp> epar blkavg
PACKAGE = imgeom
   TASK = blkavg

input   =             temp1.ms  Input images
output  =        DomeFlat_mean  Output block averaged images
(option =              average) Type of operation
b1      =                    1  blocking factor in dimension 1 (x or column)
b2      =                  160  blocking factor in dimension 2 (y or line)
b3      =                    1  blocking factor in dimension 3 (z or band)
b4      =                    1  blocking factor in dimension 4
b5      =                    1  blocking factor in dimension 5
b6      =                    1  blocking factor in dimension 6
b7      =                    1  blocking factor in dimension 7
(mode   =                    h)
Step 12: Generate normalized fiber-to-fiber thorughput variations with sarith:
sl> epar sarith 
PACKAGE = specred
   TASK = sarith

input1  =             temp1.ms  List of input spectra
op      =                    /  Operation
input2  =        DomeFlat_mean  List of input spectra or constants
output  =     DomeFlat.resp.ms  List of output spectra
(w1     =                INDEF) Starting wavelength
(w2     =                INDEF) Ending wavelength
(apertur=                     ) List of input apertures or columns/lines
(bands  =                     ) List of input bands or lines/bands
(beams  =                     ) List of input beams or echelle orders
(apmodul=                    0) Input aperture modulus (0=none)
(reverse=                   no) Reverse order of operands in binary operation?
(ignorea=                  yes) Ignore second operand aperture numbers?
(format =            multispec) Output spectral format
(renumbe=                   no) Renumber output apertures?
(offset =                    0) Output aperture number offset
(clobber=                   no) Modify existing output images?
(merge  =                   no) Merge with existing output images?
(rebin  =                  yes) Rebin to exact wavelength region?
(errval =                   0.) Arithmetic error replacement value
(verbose=                  yes) Print operations?
(mode   =                    h)
Step 13: Smooth the fiber responses with fit1d:

sl> epar fit1d

PACKAGE = imfit
   TASK = fit1d

input   =     DomeFlat.resp.ms  Images to be fit
output  =     DomeFlat.resp.ms  Output images
(axis   =                    1) Axis to be fit
type    =                  fit  Type of output (fit, difference, ratio)
(interac=                  yes) Set fitting parameters interactively?
(sample =                    *) Sample points to use in fit
(naverag=                    1) Number of points in sample averaging
(functio=              spline3) Fitting function
(order  =                   20) Order of fitting function
(low_rej=                   3.) Low rejection in sigma of fit
(high_re=                   3.) High rejection in sigma of fit
(niterat=                    5) Number of rejection iterations
(grow   =                   0.) Rejection growing radius in pixels
(graphic=             stdgraph) Graphics output device
(cursor =                     ) Graphics cursor input
(mode   =                    h)
Do a few interactively, then batch the rest (by entering a return).

Step 14: Extract the object spectra with apall:

sl> epar apall
PACKAGE = specred
   TASK = apall
    
input   =            s102n.scat  List of input images
(output =        s102n.scat.ms) List of output spectra
(apertur=                     ) Apertures
(format =            multispec) Extracted spectra format
(referen=                s106n) List of aperture reference images
(profile=                     ) List of aperture profile images
(interac=                   no) Run task interactively?
(find   =                   no) Find apertures?
(recente=                  yes) Recenter apertures?
(resize =                   no) Resize apertures?
(edit   =                   no) Edit apertures?
(trace  =                   no) Trace apertures?
(fittrac=                   no) Fit the traced points interactively?
(extract=                  yes) Extract spectra?
(extras =                   no) Extract sky, sigma, etc.?
(review =                   no) Review extractions?
(line   =                INDEF) Dispersion line
(nsum   =                   10) Number of dispersion lines to sum or median

                                # DEFAULT APERTURE PARAMETERS

(lower  =                  -5.) Lower aperture limit relative to center
(upper  =                   5.) Upper aperture limit relative to center
(apidtab=                     ) Aperture ID table (optional)

                                # DEFAULT BACKGROUND PARAMETERS

(b_funct=            chebyshev) Background function
(b_order=                    1) Background function order
(b_sampl=          -10:-6,6:10) Background sample regions
(b_naver=                   -3) Background average or median
(b_niter=                    0) Background rejection iterations
(b_low_r=                   3.) Background lower rejection sigma
(b_high_=                   3.) Background upper rejection sigma
(b_grow =                   0.) Background rejection growing radius

                                # APERTURE CENTERING PARAMETERS

(width  =                   5.) Profile centering width
(radius =                  10.) Profile centering radius
(thresho=                   0.) Detection threshold for profile centering

                                # AUTOMATIC FINDING AND ORDERING PARAMETERS

nfind   =                       Number of apertures to be found automatically
(minsep =                   5.) Minimum separation between spectra
(maxsep =                1000.) Maximum separation between spectra
(order  =           increasing) Order of apertures

                                # RECENTERING PARAMETERS

(aprecen=                     ) Apertures for recentering calculation
(npeaks =                INDEF) Select brightest peaks
(shift  =                  yes) Use average shift instead of recentering?

                                # RESIZING PARAMETERS

(llimit =                INDEF) Lower aperture limit relative to center
(ulimit =                INDEF) Upper aperture limit relative to center
(ylevel =                  0.1) Fraction of peak or intensity for automatic widt
(peak   =                  yes) Is ylevel a fraction of the peak?
(bkg    =                  yes) Subtract background in automatic width?
(r_grow =                   0.) Grow limits by this factor
(avglimi=                   no) Average limits over all apertures?

                                # TRACING PARAMETERS

(t_nsum =                   10) Number of dispersion lines to sum
(t_step =                   10) Tracing step
(t_nlost=                    3) Number of consecutive times profile is lost befo
(t_funct=             legendre) Trace fitting function
(t_order=                    2) Trace fitting function order
(t_sampl=                    *) Trace sample regions
(t_naver=                    1) Trace average or median
(t_niter=                    0) Trace rejection iterations
(t_low_r=                   3.) Trace lower rejection sigma
(t_high_=                   3.) Trace upper rejection sigma
(t_grow =                   0.) Trace rejection growing radius

                                # EXTRACTION PARAMETERS

(backgro=                 none) Background to subtract
(skybox =                    1) Box car smoothing length for sky
(weights=             variance) Extraction weights (none|variance)
(pfit   =                fit1d) Profile fitting type (fit1d|fit2d)
(clean  =                  yes) Detect and replace bad pixels?
(saturat=               28000.) Saturation level
(readnoi=                  6.3) Read out noise sigma (photons)
(gain   =                  1.6) Photon gain (photons/data number)
(lsigma =                   3.) Lower rejection threshold
(usigma =                   3.) Upper rejection threshold
(nsubaps=                    1) Number of subapertures per aperture
(mode   =                    h)
NOTE: extracting the ~150 norris spectra takes about 6-7 minutes on an unloaded Ultra1.

Step 15: Test the fiber-to-fiber throughput correction:

Be sure that you have deleted the temp1.ms from a previous reduction step (using imdel).

IMPORTANT NOTE: It is crucial to execute flpr here! Because a file called temp1.ms was made before, there seems to be some leftover knowledge of it, which causes sarith to barf. We believe this problem may be Solaris-specific.

The sarith/scopy/specplot steps are only a test to see if the response correction was successful. After checking the result of each image (using specplot temp2.ms) you can delete temp1.ms and temp2.ms.

sl> epar sarith
PACKAGE = specred
   TASK = sarith

input1  =         s102n.scat.ms  List of input spectra
op      =                    /  Operation
input2  =      DomeFlat.resp.ms  List of input spectra or constants
output  =              temp1.ms  List of output spectra
(w1     =                INDEF) Starting wavelength
(w2     =                INDEF) Ending wavelength
(apertur=                     ) List of input apertures or columns/lines
(bands  =                     ) List of input bands or lines/bands
(beams  =                     ) List of input beams or echelle orders
(apmodul=                    0) Input aperture modulus (0=none)
(reverse=                   no) Reverse order of operands in binary operation?
(ignorea=                   no) Ignore second operand aperture numbers?
(format =            multispec) Output spectral format
(renumbe=                   no) Renumber output apertures?
(offset =                    0) Output aperture number offset
(clobber=                   no) Modify existing output images?
(merge  =                   no) Merge with existing output images?
(rebin  =                  yes) Rebin to exact wavelength region?
(errval =                   0.) Arithmetic error replacement value
(verbose=                  yes) Print operations?
(mode   =                    h)
Copy out the sky fibers (beam=0) to see if they all have the same level:
sl> epar scopy

PACKAGE = specred
   TASK = scopy

input   =             temp1.ms  List of input spectra
output  =             temp2.ms  List of output spectra
(w1     =                INDEF) Starting wavelength
(w2     =                INDEF) Ending wavelength
(apertur=                     ) List of input apertures or columns/lines
(bands  =                     ) List of input bands or lines/bands
(beams  =                    0) List of beams or echelle orders
(apmodul=                    0) Input aperture modulus (0=none)
(format =            multispec) Output spectra format
(renumbe=                   no) Renumber output apertures?
(offset =                    0) Output aperture number offset
(clobber=                   no) Modify existing output images?
(merge  =                   no) Merge with existing output images?
(rebin  =                  yes) Rebin to exact wavelength region?
(verbose=                  yes) Print operations?
(mode   =                   ql)
Inspect the spectra and make sure they all have the same flux:
cl> epar specplot

PACKAGE = specred
   TASK = specplot

spectra =              temp2.ms  List of spectra to plot
(apertur=                     ) Apertures to plot
(bands  =                    1) Bands of 3D images to plot
(autolay=                   no) Use automatic layout algorithm?
(autosca=                   no) Scale to common mean for automatic layout?
(fractio=                   1.) Fraction of automatic minimum separation step
(units  =                     ) Coordinate units
(scale  =                   1.) Intensity scale (value, @file, keyword)
(offset =                   0.) Intensity offset (value, @file, keyword)
(step   =                   0.) Default separation step
(ptype  =                    1) Plotting type
(labels =                 user) Type of labels
(ulabels=                     ) User labels (file)
(xlpos  =                 1.02) X label position (fraction of range)
(ylpos  =                   0.) Y label position (fraction of mean)
(sysid  =                  yes) Include system banner and step value?
(yscale =                  yes) Draw Y axis scale?
(title  =                     ) Plot title
(xlabel =                     ) X axis label
(ylabel =                     ) Y axis label
(xmin   =                INDEF) X axis left limit
(xmax   =                INDEF) X axis right limit
(ymin   =                INDEF) Y axis bottom limit
(ymax   =                INDEF) Y axis top limit
(logfile=                     ) Logfile
(graphic=             stdgraph) Graphics output device
(cursor =                     ) Cursor input
(mode   =                   ql)
Step 16: Apply the fiber-to-fiber throughput correction:

Now that we have checked that the response correction works, we apply Dome_response to all object frames for real. This is exactly the same as step 15, except that we now know it works and we do these one at a time, overwriting the images (by setting clobber=yes):

PACKAGE = specred
   TASK = sarith

input1  =         s102n.scat.ms List of input spectra
op      =                    /  Operation
input2  =      DomeFlat.resp.ms List of input spectra or constants
output  =        s102n.scat.ms  List of output spectra
(w1     =                INDEF) Starting wavelength
(w2     =                INDEF) Ending wavelength
(apertur=                     ) List of input apertures or columns/lines
(bands  =                     ) List of input bands or lines/bands
(beams  =                     ) List of input beams or echelle orders
(apmodul=                    0) Input aperture modulus (0=none)

(reverse=                   no) Reverse order of operands in binary operation?
(ignorea=                   no) Ignore second operand aperture numbers?

(format =            multispec) Output spectral format
(renumbe=                   no) Renumber output apertures?
(offset =                    0) Output aperture number offset
(clobber=                  yes) Modify existing output images?
(merge  =                   no) Merge with existing output images?
(rebin  =                  yes) Rebin to exact wavelength region?
(errval =                   0.) Arithmetic error replacement value
(verbose=                  yes) Print operations?
(mode   =                    h)
Step 17: Identify arc lines in one fiber of one arc with identify:
cl> epar identify 

PACKAGE = specred
   TASK = identify

images  =             s105n.ms  Images containing features to be identified
(section=          middle line) Section to apply to two dimensional images
(databas=             database) Database in which to record feature data
(coordli=     scrpto$fe+ar.arc) User coordinate list
(units  =                     ) Coordinate units
(nsum   =                   10) Number of lines/columns/bands to sum in 2D image
(match  =                  -3.) Coordinate list matching limit
(maxfeat=                   50) Maximum number of features for automatic identif
(zwidth =                 100.) Zoom graph width in user units
(ftype  =             emission) Feature type
(fwidth =                   4.) Feature width in pixels
(cradius=                   5.) Centering radius in pixels
(thresho=                   0.) Feature threshold for centering
(minsep =                   2.) Minimum pixel separation
(functio=              spline3) Coordinate function
(order  =                    1) Order of coordinate function
(sample =                    *) Coordinate sample regions
(niterat=                    0) Rejection iterations
(low_rej=                   3.) Lower rejection sigma
(high_re=                   3.) Upper rejection sigma
(grow   =                   0.) Rejection growing radius
(autowri=                   no) Automatically write to database
(graphic=             stdgraph) Graphics output device
(cursor =                     ) Graphics cursor input
crval   =                       Approximate coordinate (at reference pixel)
cdelt   =                       Approximate dispersion
(aidpars=                     ) Automatic identification algorithm parameters
(mode   =                   ql)
Visually id a line and hit "m". Enter the wavelenth in Angstroms. Use "r" to redraw. ID 4 or 5 lines. Use "f" to fit. Use "q" in fitting window to get back to arc line map. Use "l" to add lines from the linelist. Use "d" in either the arc line map window or the fit window to delete unwanted lines.
To put labels on the tic-mark positions of lines use ":label coord" in the arc line plot.
A good fit should have 30 or so lines with an RMS around 0.10.
Hit "q" to quit.

Step 18: Identify and fit arc lines in other fibers with reid:

cl> epar reidentify  

PACKAGE = specred
   TASK = reidentify

referenc=             s105n.ms  Reference image
images  =             s105n.ms  Images to be reidentified
(interac=                  yes) Interactive fitting?
(section=          middle line) Section to apply to two dimensional images
(newaps =                  yes) Reidentify apertures in images not in reference?
(overrid=                   no) Override previous solutions?
(refit  =                  yes) Refit coordinate function?
(trace  =                  yes) Trace reference image?
(step   =                    1) Step in lines/columns/bands for tracing an image
(nsum   =                   10) Number of lines/columns/bands to sum
(shift  =                   0.) Shift to add to reference features (INDEF to sea
(search =                   0.) Search radius
(nlost  =                   15) Maximum number of features which may be lost
(cradius=                   5.) Centering radius
(thresho=                   0.) Feature threshold for centering
(addfeat=                   no) Add features from a line list?
(coordli=     scrpto$fe+ar.arc) User coordinate list
(match  =                  -3.) Coordinate list matching limit
(maxfeat=                   50) Maximum number of features for automatic identif
(minsep =                   2.) Minimum pixel separation
(databas=             database) Database
(logfile=              logfile) List of log files
(plotfil=                     ) Plot file for residuals
(verbose=                   no) Verbose output?
(graphic=             stdgraph) Graphics output device
(cursor =                     ) Graphics cursor input
answer  =                  yes  Fit dispersion function interactively?
crval   =                       Approximate coordinate (at reference pixel)
cdelt   =                       Approximate dispersion
(aidpars=                     ) Automatic identification algorithm parameters
(mode   =                   ql)
Inspect the RMS values for each fiber. If <= 0.10 (or so) say "no" to interactive fit. Otherwise, say "y". Hit "f" to check fit and "d" to delete bad points."q" to get out of fitting and see arc plot; "m" to mark lines.

Step 19: Identify and fit arc lines in other images using reid:

cl> epar reidentify  

PACKAGE = specred
   TASK = reidentify

referenc=             s105n.ms  Reference image
images  =             s109n.ms  Images to be reidentified
Step 20: Apply the dispersion solutions to the object frames using refspec and dispcor:

Make listings of object and arc frames:

cl> hselect *.ms.fits $I "imagetyp = 'object'" >obj.list
cl> hselect *.ms.fits $I "imagetyp = 'arc'" >arc.list      
Assign reference arcs to object frames:
cl> epar refspectra

PACKAGE = specred
   TASK = refspectra

input   =            @obj.list  List of input spectra
(referen=            @arc.list) List of reference spectra
(apertur=                     ) Input aperture selection list
(refaps =                     ) Reference aperture selection list
(ignorea=                  yes) Ignore input and reference apertures?
(select =               interp) Selection method for reference spectra
(sort   =                frame) Sort key
(group  =                     ) Group key
(time   =                   no) Is sort key a time?
(timewra=                  17.) Time wrap point for time sorting
(overrid=                   no) Override previous assignments?
(confirm=                  yes) Confirm reference spectrum assignments?
(assign =                  yes) Assign the reference spectra to the input spectr
(logfile=       STDOUT,logfile) List of logfiles
(verbose=                   no) Verbose log output?
answer  =                       Accept assignment?
(mode   =                   ql)
Apply the dispersion correction:
cl> epar dispcor 

PACKAGE = specred
   TASK = dispcor

input   =            @obj.list  List of input spectra
output  =            @obj.list  List of output spectra
(lineari=                  yes) Linearize (interpolate) spectra?
(databas=             database) Dispersion solution database
(table  =                     ) Wavelength table for apertures
(w1     =                INDEF) Starting wavelength
(w2     =                INDEF) Ending wavelength
(dw     =                INDEF) Wavelength interval per pixel
(nw     =                INDEF) Number of output pixels
(log    =                   no) Logarithmic wavelength scale?
(flux   =                  yes) Conserve flux?
(samedis=                   no) Same dispersion in all apertures?
(global =                   no) Apply global defaults?
(ignorea=                   no) Ignore apertures?
(confirm=                   no) Confirm dispersion coordinates?
(listonl=                   no) List the dispersion coordinates only?
(verbose=                  yes) Print linear dispersion assignments?
(logfile=                     ) Log file
(mode   =                   ql)
Step 21: Produce a high S/N sky spectrum from sky fibers using scopy:
NOTE: You will need to perform steps 21 and 22 for each object frame!

sl> epar scopy
PACKAGE = specred
   TASK = scopy

input   =        s102n.scat.ms  List of input spectra
output  =                  sky  List of output spectra
(w1     =                INDEF) Starting wavelength
(w2     =                INDEF) Ending wavelength
(apertur=                     ) List of input apertures or columns/lines
(bands  =                     ) List of input bands or lines/bands
(beams  =                    0) List of beams or echelle orders
(apmodul=                    0) Input aperture modulus (0=none)
(format =             onedspec) Output spectra format
(renumbe=                   no) Renumber output apertures?
(offset =                    0) Output aperture number offset
(clobber=                   no) Modify existing output images?
(merge  =                   no) Merge with existing output images?
(rebin  =                  yes) Rebin to exact wavelength region?
(verbose=                  yes) Print operations?
(mode   =                   ql)
Use specplot to examine sky spectra:
sl> epar specplot
PACKAGE = specred
   TASK = specplot

spectra =            sky*.fits  List of spectra to plot
(apertur=                     ) Apertures to plot
(bands  =                     ) Bands of 3D images to plot
(autolay=                   no) Use automatic layout algorithm?
(autosca=                   no) Scale to common mean for automatic layout?
(fractio=                   1.) Fraction of automatic minimum separation step
(units  =                     ) Coordinate units
(scale  =                   1.) Intensity scale (value, @file, keyword)
(offset =                   0.) Intensity offset (value, @file, keyword)
(step   =                   0.) Default separation step
(ptype  =                    1) Plotting type
(labels =                 user) Type of labels
(ulabels=                     ) User labels (file)
(xlpos  =                 1.02) X label position (fraction of range)
(ylpos  =                   0.) Y label position (fraction of mean)
(sysid  =                  yes) Include system banner and step value?
(yscale =                  yes) Draw Y axis scale?
(title  =                     ) Plot title
(xlabel =                     ) X axis label
(ylabel =                     ) Y axis label
(xmin   =                INDEF) X axis left limit
(xmax   =                INDEF) X axis right limit
(ymin   =                INDEF) Y axis bottom limit
(ymax   =                INDEF) Y axis top limit
(logfile=                     ) Logfile
(graphic=             stdgraph) Graphics output device
(cursor =                     ) Cursor input
(mode   =                   ql)
Delete any fibers that look out of line.
Use scombine to make final sky spectra:
cl> epar scombine 
PACKAGE = specred
   TASK = scombine

input   =            sky*.fits  List of input spectra
output  =       s102n_sky.fits  List of output spectra
(noutput=                     ) List of output number combined spectra
(logfile=               STDOUT) Log file
(apertur=                     ) Apertures to combine
(group  =                  all) Grouping option
(combine=              average) Type of combine operation
(reject =              sigclip) Type of rejection
(first  =                   no) Use first spectrum for dispersion?
(w1     =                INDEF) Starting wavelength of output spectra
(w2     =                INDEF) Ending wavelength of output spectra
(dw     =                INDEF) Wavelength increment of output spectra
(nw     =                INDEF) Length of output spectra
(log    =                   no) Logarithmic increments?
(scale  =                 none) Image scaling
(zero   =                 none) Image zero point offset
(weight =                 none) Image weights
(sample =                     ) Wavelength sample regions for statistics
(lthresh=                INDEF) Lower threshold
(hthresh=                INDEF) Upper threshold
(nlow   =                    1) minmax: Number of low pixels to reject
(nhigh  =                    1) minmax: Number of high pixels to reject
(nkeep  =                    1) Minimum to keep (pos) or maximum to reject (neg)
(mclip  =                  yes) Use median in sigma clipping algorithms?
(lsigma =                   3.) Lower sigma clipping factor
(hsigma =                   3.) Upper sigma clipping factor
(rdnoise=                   0.) ccdclip: CCD readout noise (electrons)
(gain   =                   1.) ccdclip: CCD gain (electrons/DN)
(snoise =                   0.) ccdclip: Sensitivity noise (fraction)
(sigscal=                  0.1) Tolerance for sigma clipping scaling corrections
(pclip  =                 -0.5) pclip: Percentile clipping parameter
(grow   =                    0) Radius (pixels) for 1D neighbor rejection
(blank  =                   0.) Value if there are no pixels
(mode   =                   ql)
Check the output with splot. If the mean sky is good, delete the individual skys (using imdel sky*.fits).

Step 22: Subtract the sky from each fiber:

sl> epar sarith
PACKAGE = specred
   TASK = sarith

input1  =   s102n.scat.ms.fits  List of input spectra
op      =                    -  Operation
input2  =       s102n_sky.fits  List of input spectra or constants
output  =   s102n.scat.ms.fits  List of output spectra
(w1     =                INDEF) Starting wavelength
(w2     =                INDEF) Ending wavelength
(apertur=                     ) List of input apertures or columns/lines
(bands  =                     ) List of input bands or lines/bands
(beams  =                     ) List of input beams or echelle orders
(apmodul=                    0) Input aperture modulus (0=none)
(reverse=                   no) Reverse order of operands in binary operation?
(ignorea=                  yes) Ignore second operand aperture numbers?
(format =            multispec) Output spectral format
(renumbe=                   no) Renumber output apertures?
(offset =                    0) Output aperture number offset
(clobber=                  yes) Modify existing output images?
(merge  =                   no) Merge with existing output images?
(rebin  =                  yes) Rebin to exact wavelength region?
(errval =                   0.) Arithmetic error replacement value
(verbose=                  yes) Print operations?
(mode   =                    h)
Step 23: Combine all object frames for a given setup with scombine:
sl> epar scombine
PACKAGE = specred
   TASK = scombine

input   =      s102n.scat.ms,s104n.scat.ms  List of input spectra
output  =           poss314.ms  List of output spectra
(noutput=                     ) List of output number combined spectra
(logfile=               STDOUT) Log file

(apertur=                     ) Apertures to combine
(group  =            apertures) Grouping option
(combine=              average) Type of combine operation
(reject =               minmax) Type of rejection

(first  =                   no) Use first spectrum for dispersion?
(w1     =                INDEF) Starting wavelength of output spectra
(w2     =                INDEF) Ending wavelength of output spectra
(dw     =                INDEF) Wavelength increment of output spectra
(nw     =                INDEF) Length of output spectra
(log    =                   no) Logarithmic increments?
(scale  =                 mode) Image scaling
(zero   =                 none) Image zero point offset
(weight =                 mode) Image weights
(sample =                     ) Wavelength sample regions for statistics
(lthresh=                INDEF) Lower threshold
(hthresh=                INDEF) Upper threshold
(nlow   =                    0) minmax: Number of low pixels to reject
(nhigh  =                    1) minmax: Number of high pixels to reject
(nkeep  =                   -1) Minimum to keep (pos) or maximum to reject (neg)
(mclip  =                  yes) Use median in sigma clipping algorithms?
(lsigma =                   3.) Lower sigma clipping factor
(hsigma =                   3.) Upper sigma clipping factor
(rdnoise=                   0.) ccdclip: CCD readout noise (electrons)
(gain   =                   1.) ccdclip: CCD gain (electrons/DN)
(snoise =                   0.) ccdclip: Sensitivity noise (fraction)
(sigscal=                  0.1) Tolerance for sigma clipping scaling corrections
(pclip  =                 -0.5) pclip: Percentile clipping parameter
(grow   =                    0) Radius (pixels) for 1D neighbor rejection
(blank  =                   0.) Value if there are no pixels
(mode   =                   ql)

Believe it or not, you are done! You can use scopy to copy out individual apertures; splot to display your spectra.

Roy Gal
11/10/99