
                        PHOENIX DATA REDUCTIONS
                        -----------------------

The data reductions for this Phoenix data set follow the reductions
outlined in the Cryogenic Spectrometer User Manual August 1995 with some
slight modifications - there are probably other ways to do these reductions
depending on the actual data and what you hope to get from it, so be sure
to experiment with different techniques.  The data are provided by Caty 
Pilachowski from a run at Kitt Peak on the 2.1m with Phoenix in late 
February 1997.  The data were observed at 2.37um.

The data set consists of 20 dark frames, 20 flats, 2 observations of an 
object, and 2 observations of a telluric line star.  The dark frames 
were taken with the instrument set up for flatfield integrations but with 
the calibration lamp off.  The flats were taken at the same integration 
time as the darks but with the lamp on.  The observations of the stars 
were taken one right after the other, at the same integration time, 
but the star is on different parts of the slit for each integration.

An assumption is made during this exercise that the user has a 
mimimal working knowledge of IRAF.  If that is not the case then it is 
recommended that you (a) read through "A Beginner's Guide to IRAF" (available
from iraf.noao.edu in the pub/beguide.ps.Z file), and (b) work through
the "intro" exercise (available from iraf.noao.edu in the iraf/misc directory -
see the exer211.readme file).  This Phoenix exercise was updated for
IRAF Version 2.11, although the tasks used should work in a similar way 
in earlier releases.

We will assume that you are running IRAF from an XGTERM window, that your 
display window is XIMTOOL, and that the Phoenix data set is in a subdirectory
of your IRAF login directory (as discussed in the README file accompanying
the data).  Go to directory from which you plan to run IRAF and type the 
following.  

    % cl

This will load IRAF for you. Notice that the prompt is now "cl" (we will
continue to use the "cl" prompt in our examples below even though this
prompt will change on your screen as you load other IRAF packages.  Now go 
to the directory which contains the Phoenix data and list the files.

    cl> cd phx
    cl> dir

The data are stored in FITS formatted files.  We will unpack the data into 
the native IRAF image format using the task RFITS.  We will restore these 
FITS images to disk using the previous names from when they resided on disk.

    cl> unlearn rfits
    cl> lpar rfits
    cl> rfits dark*,flat*,n20* * tmp old+
    cl> dir
    cl> imheader *.imh
    cl> imhead dark_n2021.imh l+

When you list the directory this time you will see additional files with
the .imh extension - these are the IRAF image files.  The first IMHEADER 
command lists information about each image - mainly the image size,
the pixel type, and the title.  The second use of the command lists all of
the header information for the specified image.  Before proceeding further
let's cover one more general topic - image display.  Check to be sure you 
have an XIMTOOL window running.

    cl> unlearn display
    cl> display flat_n2004.imh 1 fill+

XIMTOOL has lots of options so if you are not familiar with it take a few
moments to experiment a bit.  With the cursor in the image window hold
down the right mouse button and move it up/down and left/right to vary
the contrast and brightness.  You may want to look at one of the dark
frames and one of the object frames as well at this time so you know
what your data look like.  Put theses images into different frames of
XIMTOOL so you can compare them.

    cl> display dark_n2033.imh 2 fill+
    cl> display n2069.imh 3 fill+

You can quickly cycle through the frames with the Control/b and/or Control/f
keys while the cursor is in the XIMTOOL window.  You can also bring up the
XTIMOOL Control panel and auto-blink the frames (press the left-most button
in the series of buttons in the right upper corner of the window).

IRAF's default frame buffer, the memory allocated to hold the image in 
XIMTOOL, is 512 x 512 - the fill+ parameter forced DISPLAY to interpolate 
the image so it fit within these constraints.  We could also choose another 
frame buffer size that better matches our images.  The task GDEVICES lists
the available frame buffer sizes. Choosing a size equal to or larger than
our images means we no longer need the fill+ parameter.  Let's re-display 
our star with this new frame buffer size.

    cl> gdevices
    cl> set stdimage=imt14
    cl> display n2069.imh 1

Notice now that you are only seeing part of the image but you are seeing each
pixel.  The part of the image that is visible is outlined in the small
panner window at the top-right of the IMTOOL display.  You can pan around
by moving the image cursor to an area of the image that you wish to center
and click on the middle mouse button (in either the main image window or in
the panner box).  You can also adjust the size of the XIMTOOL window if you 
wish by pulling on a corner of the window with the cursor while holding
down the left mouse button.

This may be a good time to look at the task IMPLOT, as well, which is a 
general plotting task.  This task allows you to make row or column slices
through an image using the XGTERM graphics window.  

    cl> implot dark_n2033
           [use "c" or "l" to make column or row plots using the current
            cursor position in the XGTERM graphics window - type "?" to see a 
            complete list of the cursor commands available to you in this task 
            - type "q" to quit the task

            let's try a few things:

            notice across the top of the plot the label says you are currently 
            looking a line 512 - put the vertical cursor at 100 using
            the x-axis labels (which are column numbers, right?) - press "c"
            and note now that the label across the top says you are looking
            at column 100 (or thereabouts)

            the numbers along the x-axis are now line numbers - but there
            are virtual column labels on the right vertical axis - place
            the horizontal cursor on 200 and press "c"

            the "l" (line or row) key works in a similar way

            now type ":c 100 200" (Return) - notice now you are looking at 
            the average of these columns

            type ":a 10" (Return) and then ":l 300" (Return) - "a" sets the 
            number of columns or rows to be plotted - by now you probably
            see a pattern, single stroke cursor keys do not require a
            "Return" but : commands do

            place both cursors where the vertical cursor cuts 50 on the x-axis 
            and type "e" - then move both cursors to where 40 on the left 
            vertical axis intersects with 200 on the horizontal axis - type "e"
            again - the "box" that we have defined is now displayed in a 
            zoomed mode

            type :l 300 to get the full plot back again

            try some other things - don't forget that "q" quits cursor mode]

If you would rather make plots from the image display then look at the
task IMEXAMINE which works in tandem with the image displayed in XIMTOOL.

    cl> imexamine n2069.imh
           [put the cursor in the IMTOOL window and see the blinking image
            cursor - press "l" for a line plot - type "?" for a list of
            the many useful commands - type "q" to quit]

(See the introductory tutorial for more help on these useful and
probably necessary tasks.)

At this point we are ready to start the data reductions.


Trimming the data frames
------------------------

As you have probably noticed the lower part of the data frames contain
some bad data.  Let's trim that bad section off - otherwise it causes
problems (not unsurmountable) later in the processing.  We will use
the task CCDPROC in the IMRED/CCDRED package under NOAO.  We will use other 
tasks in this package as well later.  First we need to decide the "trim"
parameters - that section of the image we wish to keep - this should be 
the same for all images! (But trim parameters can vary from run to run 
depending on the instrument setup.)  Let's use the task IMPLOT and one of 
the flat fields to determine this information.

    cl> implot flat_n2004.imh
          [the columns look ok so plot the middle column so we can inspect
           the row data - put the vertical cursor on about column 140 and 
           press "c" - expand the first 30 lines or so with the "e" key (place 
           both cursors just below the lower left boundary of the plot, press 
           "e", then place both cursors at the top of the graph at about line 
           50, press "e' again) - measure the value of the first good line to 
           keep with "C" - I'd say 16 - type "q" to quit]

    cl> imhead flat_n2004.imh

Our images have 256 columns and 1024 rows - we wish to trim off the first
15 rows, so our image section representing usable data will be [1:256,16:1024]
(range of good columns followed by range of good lines).  First let's checkout 
the package parameters for CCDRED (yes, even packages can have parameters that 
affect all tasks in the package).  Let's turn on the "verbose" flag so we can 
see what's going on. Notice the "backup" parameter. The CCDRED package is one 
of the few instances in which IRAF overwrites the input data - this parameter 
allows the input data to be saved in a backup directory before processing.  If 
you set this parameter to a directory name be sure to create the directory 
before you execute any tasks.  Next we need to modify the parameters for 
CCDPROC so they look like those below.

    cl> imred
    cl> ccdred
    cl> epar ccdred
    cl> lpar ccdred
   (pixeltype = "real real")    Output and calculation pixel datatypes
     (verbose = yes)            Print log information to the standard output?
     (logfile = "logfile")      Text log file
    (plotfile = "")             Log metacode plot file
      (backup = "")             Backup directory or prefix
  (instrument = "")             CCD instrument file
      (ssfile = "subsets")      Subset translation file
    (graphics = "stdgraph")     Interactive graphics output device
      (cursor = "")             Graphics cursor input
     (version = "2: October 1987") 
        (mode = "ql")           
      ($nargs = 0)   

    cl> unlearn ccdproc
    cl> epar ccdproc
    cl> lpar ccdproc
       images = ""              List of CCD images to correct
     (ccdtype = "")             CCD image type to correct
   (max_cache = 0)              Maximum image caching memory (in Mbytes)
      (noproc = no)             List processing steps only?\n
      (fixpix = no)             Fix bad CCD lines and columns?
    (overscan = no)             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?
    (illumcor = no)             Apply illumination correction?
   (fringecor = no)             Apply fringe correction?
     (readcor = no)             Convert zero level image to readout correction?
     (scancor = no)             Convert flat field image to scan correction?\n
    (readaxis = "line")         Read out axis (column|line)
     (fixfile = "")             File describing the bad lines and columns
     (biassec = "")             Overscan strip image section
     (trimsec = "[1:256,16:1024]") 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
  (minreplace = 1.)             Minimum flat field value
    (scantype = "shortscan")    Scan type (shortscan|longscan)
       (nscan = 1)              Number of short scan lines\n
 (interactive = no)             Fit overscan interactively?
    (function = "legendre")     Fitting function
       (order = 1)              Number of polynomial terms or spline pieces
      (sample = "*")            Sample points to fit
    (naverage = 1)              Number of sample points to combine
    (niterate = 1)              Number of rejection iterations
  (low_reject = 3.)             Low sigma rejection factor
 (high_reject = 3.)             High sigma rejection factor
        (grow = 0.)             Rejection growing radius
        (mode = "ql")         

Remember that when we run this task it will overwrite our input images so
we need to be certain that this is what we want to do - and that we have
backup copies of our original data.  Try trimming one image first and 
inspect it with IMPLOT.  We will use IMCOPY for this check.

    cl> imcopy flat_n2004.imh[1:256,16:1024] test
    cl> imhead test
    cl> implot test 
    cl> imdelete test

If the resultant test.imh file looks ok, then proceed with trimming all
the data. 

    cl> ccdproc *.imh
    cl> imhead n2078



Combining the dark frames
-------------------------

The dark frames were taken back-to-back and they should not really vary
very much.  If we looked carefully at the image header in the introduction
we would have noticed that each integration was 10 seconds.  Let's just
verify that, looking for the EXPTIME keyword.  Using the task HSELECT will
allow us to list just the header information that we are interested in, 
such as the image name (denoted here by $I) and the exposure time (EXPTIME).
The IMSTATISTICS task lets us look at statistics for the list of images
to help us decide if any of the darks had a problem and should be
discarded.

    cl> imhead dark_n2032.imh l+ | page
    cl> hselect dark*.imh $I,EXPTIME yes
    cl> imstat dark*.imh

The task that we want to use for combining the dark frames is DARKCOMBINE 
in the CCDRED package which we used above.  The tasks in the CCDRED package 
are header driven so we want to be sure we have the correct information in 
the headers, or we have compensated for this in the task parameters.  The 
task CCDLIST will tell us if the package recognizes our header information 
correctly.

You will note that CCDLIST reports back that it thinks our "dark" frames are 
"object" frames, but it does recognize the exposure time correctly.  So we
need to keep these things in mind as we setup the parameters for DARKCOMBINE.

    cl> ccdlist dark_n2040.imh l+
    cl> unlearn darkcombine
    cl> lpar darkcombine

The combining algorithm that we wish to use for DARKCOMBINE is "median" with
a rejection algorithm set to sigclip.  It expects the frames to be recognized
as "dark" frames so we need to change that.  The task also is setup to
do some preprocessing before "combining" but there is no preprocessing for
the Phoenix data, so we need to change that.  After running EPAR on the
task, the parameter list should look like the following.

    cl> epar darkcombine
    cl> lpar darkcombine
        input =                 List of dark images to combine
      (output = "Dark")         Output dark image root name
     (combine = "median")       Type of combine operation
      (reject = "sigclip")      Type of rejection
     (ccdtype = "object")       CCD image type to combine
     (process = no)             Process images before combining?
      (delete = no)             Delete input images after combining?
     (clobber = no)             Clobber existing output image?
       (scale = "exposure")     Image scaling
     (statsec = "")             Image section for computing statistics
        (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)
       (pclip = -0.5)           pclip: Percentile clipping parameter
       (blank = 0.)             Value if there are no pixels
        (mode = "ql")

Now we are ready to combine our darks.  The combined image will be called
Dark.imh.  Let's compare the combined dark with the original dark frames.

    cl> darkcombine dark*.imh
    cl> imstat dark*.imh,Dark.imh
    cl> implot dark_n2035.imh
          [and then overplot the same row from the combined dark frame,
           i.e., :i Dark
                  o 
                 :l 505
                 q]


Combining the flat frames
-------------------------

We need to do a similar operation with the flat frames.  The task we will
use in this case however is called FLATCOMBINE.  We should first inspect
the flat frames for consistency and then check the headers.  Again
these are called "object" frames.

    cl> imstat flat*.imh
    cl> ccdlist flat*.imh
    cl> ccdlist flat_n2014.imh l+
    cl> unlearn flatcombine
    cl> lpar flatcombine

We will again need to modify the ccdtype and process parameters to reflect
our specific needs with the Phoenix data.  The algorithm for combining flat 
data is "median" after using a rejection algorithm called "avsigclip" - see 
the help page for COMBINE for details.  The final parameter setup for 
FLATCOMBINE is listed below.

    cl> phelp combine
    cl> epar flatcombine
    cl> lpar flatcombine
        input =                 List of flat field images to combine
      (output = "Flat")         Output flat field root name
     (combine = "median")       Type of combine operation
      (reject = "avsigclip")    Type of rejection
     (ccdtype = "object")       CCD image type to combine
     (process = no)             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 = "median")       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 = "ql")           

Once we are satisfied with everything we can run the task, and then look
at the results as we did with the dark frames.

    cl> flatcombine flat*.imh
    cl> imstat flat*.imh,Flat
    cl> implot flat_n2013.imh
         [and overplot the same line from the combined Flat.imh frame,
          i.e., :i Flat
                :o 
                :l 505
                q]


Preparing the Flat frame
------------------------

The combined Flat frame needs further processing before it can be used by the 
object frames.  We first need to subtract the dark from it, and then we need
to normalize that dark subtracted Flat (scale to approximately unity). 

We will use the task IMARITH to do the dark subtraction.  It is important
that the Dark and Flat frames have the same integration times.

    cl> hselect Flat,Dark EXPTIME yes
    cl> unlearn imarith
    cl> lpar imarith
    cl> imarith Flat - Dark Flatsub v+
    cl> imstat Flat,Flatsub,Dark

Now we should attempt to normalize Flatsub to some degree.  There are
various techniques for doing this, and the method that one might use will
depend on the data.  We will normalize Flatsub using two of these
techniques as a comparison.  The task APFLATTEN is probably the task
of choice for these data (it models both the profile and the overall
spectrum shape, removing it from the flat field, and leaving only the 
sensitivity variations in the output image) but we have noticed that the 
model is being affected by the bad pixels in the array.  Since we do not have 
a ready solution to this problem we will not use this task for these data.

First let's simply normalize Flatsub by dividing it by a constant value;
we will use the mean value computed for the image.

    cl> imstat Flatsub
    cl> unlearn imarith
    cl> imarith Flatsub / 3600 Norm
    cl> implot Norm
        [make column and row plots to see how things look - notice that we
         still have the overall spectrum shape in the Norm image]

The RESPONSE task is in the TWODSPEC.LONGSLIT package.  This task fits a
function to the sum of the data across the dispersion direction, and then
divides that function into each column of the image producing an output
response image.  An image section can be used to describe only a piece of the
input image for fitting purposes.  This task attempts to remove the spectral 
shape in the dispersion direction.  Edit the parameters until they look like 
those below. We have set the "normalization" parameter to the middle of slit 
to avoid any edge affects due to the tilt - you can verify these values with 
IMPLOT.

    cl> twodspec
    cl> longslit
    cl> implot Flatsub
        :o
        :l 50
        :o 
        :l 1000
        [looking at these plots we notice the affect of the tilt - from
         this plot we wish to determine that region of the spectrum where
         light is continuous - we will use only this part of the spectrum
         to determine the fit along the dispersion - use "C" to print
         out the x-axis values - I think 90 and 176 work well - then
         plot the average of these columns]
        :c 90 176
        q
    cl> unlearn response
    cl> epar response
    cl> lpar response
  calibration = "Flatsub"       Longslit calibration images
 normalizatio = "Flatsub[90:176,*]"       Normalization spectrum images
     response = "Resp"          Response function images
 (interactive = yes)            Fit normalization spectrum interactively?
   (threshold = INDEF)          Response threshold
      (sample = "*")            Sample of points to use in fit
    (naverage = 1)              Number of points in sample averaging
    (function = "spline3")      Fitting function
       (order = 1)              Order of fitting function
  (low_reject = 0.)             Low rejection in sigma of fit
 (high_reject = 0.)             High rejection in sigma of fit
    (niterate = 1)              Number of rejection iterations
        (grow = 0.)             Rejection growing radius
    (graphics = "stdgraph")     Graphics output device
      (cursor = "")             Graphics cursor input
        (mode = "ql")       

    cl> response
        [answer yes to the question about fitting interactively]

        [you are now presented with a plot of the sum of the normalization
         image as specified by our image section with an initial fit and the 
         fitting parameters noted above - we want to fit the overall shape of 
         the curve but not the small bumps and wiggles since these are probably 
         real and we want to remove them from our data in a later step]

        [change the order of the fit and then refit the data]
        :order 5
        f
        q

Now look at this output image with IMPLOT and inspect a few rows and columns. 
Note that along the columns things are fairly flat now, but we still have some
profile variation since we did not attempt to take that out.  This should
be ok - our sky subtraction in a later step should remove these residual
effects from our data.

    cl> implot Resp


The data
--------

We are finally ready to reduce the actual data.  There are two images
in this data set of the same star, TY Vir. These data were taken one right 
after the other with the star in different slit positions in each image.
Since the background is so high and variable in these data we will use
the sky region in each image for the background subtraction in the other
image.  That way the same region of the array is used for sky and star.
We will then divide by our normalized flat, extract the two stellar spectra,
and finally add them.  

Let's look at these images with XIMTOOL - use the Control/f key to move 
between the two displayed images.  Notice the shift of the stellar spectrum.

    cl> display n2069.imh 1
    cl> display n2070.imh 2

Another way to look at the data is with IMPLOT, overplotting the second
spectrum on the first.  Notice the two stellar profiles offset from one
another.

    cl> implot n2069
         [:i n2070
          :o
          :l 505] 

These images should have the same exposure time but let's check that just
to be sure.

    cl> hselect n2069,n2070 exptime yes

Use IMARITH to subtract each frame from the other to do the sky subtraction.  
Then divide each frame by the normalized flatfield, Norm or Resp - let's
use Resp. Look at the resultant frames with IMPLOT.  The star is certainly 
easier to see in these processed frames!  (Let's do the same for the
telluric line star while we are at it.)

    cl> imarith n2069 - n2070 tyvir1ss
    cl> imarith n2070 - n2069 tyvir2ss
    cl> imarith n2075 - n2076 avir1ss
    cl> imarith n2076 - n2075 avir2ss
    cl> imarith tyvir1ss / Resp tyvir1
    cl> imarith tyvir2ss / Resp tyvir2
    cl> imarith avir1ss / Resp avir1
    cl> imarith avir2ss / Resp avir2
    cl> implot tyvir1
         :l 400 450
         [Let's measure the full width (FW) of this profile while we are in 
         IMPLOT - we will use this value in the APALL step when we extract 
         these data to 1-d: expand the plot around the profile of the star
         on the left by marking the lower-left corner of the expansion box
         and press "e" - move to the upper-right corner of the expansion box 
         and press "e" again.  Now mark the left edge of the spectrum and press 
        "C" and note the first number - then move to the right edge and do the 
         same - subtract the two numbers.  11 works for me.]

         :i tyvir2 
         :l 400 450

         [Then do the same thing for the spectrum on the right.  11 works for 
         me here as well.}

Now the next big step is to use the task APALL and extract the star in each
image to a 1-d spectrum - one could also use the sky subtraction option in 
this task if it was felt that a second order sky subtraction was needed, but I
think things look pretty good so we will turn that option off.  APALL is
one of the many tasks in IRAF that has enough parameters to make ice cream
and it is probably worth your time to read about things in detail when you 
have some time.  But for now we will go with the defaults and change things
on the fly as we do things interactively.

    cl> apextract
    cl> unlearn apall
    cl> lpar apall
    cl> apall tyvir1
        
        Step 1: setting up the aperture parameters
          answer "no" to the question about finding apertures automatically
          and "yes" to editing the apertures

          when a cross-cut of the spectrum is displayed put the vertical 
          cursor on the center of the stellar profile (on the left)
          and press "m"

          the profile center and the default extraction apertures appear
          on the plot and on the status line at the bottom of the plot -
          use the vertical cursor to see where the extraction apertures
          fall on the profile

          use the "y" key to adjust the extraction apertures - move the 
          horizontal cursor along the profile until you find a spot for
          which everything above should be extracted and press "y" - notice
          the new extraction apertures on the plot and on the status line -
          do this until you are happy with the aperture settings

          type in the profile width that we computed above by typing 
          :width 11 - this is the number of pixels used to find the center
          of the profile

          check the centering radius - the number of pixels that the 
          center is allowed to wander during the trace - we want to restrict
          this so any poor centers caused by bad pixels are excluded but at 
          the same time be wide enough to compensate for any tilt in the
          spectrum - 10 should work ok - :radius 10

          type :line 550, note that things are not quite centered, type
          "c" to center things up

          there are just tons of cursor options in all the IRAF interactive
          tasks so type ? if you would like to see the list - they vary from
          task to task

          now we are ready to trace the spectrum - type "q"

        Step 2: tracing
          answer "yes" to the next series of questions

          you will now be presented with a plot of the column position
          versus the line position for each central point for the star
          spectrum found during the trace - the parameters for the fit
          to these points are printed at the top of the plot - let's
          modify them by typing :order 3, :niterate 2, followed by an "f"
          to redo the fit with these new parameters - the RMS value is the
          error of the fit in pixel units

          type ? if you would like to see more about your fitting options

          this new fit looks good to me - typing "q" will take you to the
          extraction process

        Step 3: extraction
          answer "yes" to the next series of questions - the parameters
          for the aperture including the trace will be saved in a 
          directory call database and will automatically be recalled if
          you look at this image again

          the extracted spectrum is plotted for review - we have some
          bad pixels in the data, don't we - APALL has an option for
          cleaning the data during extraction - you may want to go back
          and experiment with that option but you will need to know the
          gain and readnoise for your data.

          type "q" to exit

Let's see what was produced from this extraction.  Looking in your directory
you should see a spectrum called tyvir1.ms.imh.  Let's look at it further
with SPLOT after first loading the ONEDSPEC package.  SPLOT has its own set
of cursor options so type "?" to view them.  Try the autoexpand key "a"
by placing the vertical cursor first on the left side and then the right 
side of the window you wish to expand.  Also try the "w" option which has
its own set of cursor commands.

    cl> dir ty*
    cl> onedspec
    cl> splot tyvir1.ms.imh
         [if the scaling is off due to bad pixels first place the vertical 
         cursor on left side of the spectrum to the right of any bad pixels, 
         type "a", then place the cursor on right side of spectrum to the left
         of any bad pixels, type "a" - place the horizontal cursor above the 
         good data and type "wt" - place the horizontal cursor just below
         the good data and type "wb"]

         [let's go back to the full spectrum by typing "c" and then try to
          clean out the one or two bad pixel values.  Expand the region around 
          the bad pixels using the "a" key - then using both cursors mark
          the first good pixel left of the first bad pixel - type "x" - move
          the cursors to the right side of the bad pixel and type "x"
          again - the two good pixels should now be joined.  Do the same
          for any other bad pixels.  Now type "c" - does that look better?]

         [you can save the revised spectrum with the "i" key if you would
          like - then exit with "q".]

Do the same processing for tyvir2.imh, avir1, and avir2 - extracting to 1-d.
Be sure to answer "no" to the question about finding the apertures; set 
:width 11, :radius 10 while you are in the edit mode.  Note that these 
settings are remembered for the second and third image.  Also remember to
set the tracing parameters for the first image (:order 3, :niterate 2) -
they will also be remembered for the next two stars.

    cl> apall tyvir2.imh,avir1,avir2

Now use SPLOT to clean the bad pixels from the avir1.ms image and save it
as a new image called Dispim.  We will use this new image to compute our 
wavelength correction.

    cl> splot avir1.ms
        [use the "x" key to interpolate over the bad pixels after expanding
         the plot with the "a" key]

        i
        Dispim
        q


The dispersion correction
-------------------------

The next step is to compute a wavelength calibration using our Dispim spectrum
and the task IDENTIFY. For this particular spectrum we will use methane lines 
selected from the ATMOS line list to be relatively free of blends, and not use 
the water vapor lines.  We will apply this calibration to the extracted 
data using the task DISPCOR.  Print out the file "phx.eps" which contains a 
plot of the Dispim spectrum with the methane line wavenumbers indicated for
those lines in the spectrum that we will use to compute the wavelength 
solution.

First inspect the Dispim image with SPLOT to determine the full width of
the lines that will be used. I first "zoomed" the plot with the "z" key,
and then used "." and "," to move up and down the spectrum while using
the "spacebar" to help measure the line widths (place both cursors on one
side of the line and press "spacebar" - notice the x position - then do the
same for the other side of the line).  I computed a value of 12 - this is not 
quite the full width but close enough for these very strong lines - so 12 
pixels for each line will be used to compute the line center.

Let's modify the parameters for IDENTIFY until they look like those below.
We want to be sure to set the "coordlist" parameter to the null ("")
character since we will be entering our coordinates from the terminal.  We
will be using wavenumbers for our units ("units=wn").  The "ftype"
parameter is set to "absorption" and "fwidth" is set to the value we
determined above, 12.

    cl> splot Dispim
    cl> unlearn identify
    cl> epar identify
    cl> lpar identify
       images = "Dispim"        Images containing features to be identified
        crval =                 Approximate coordinate (at reference pixel)
        cdelt =                 Approximate dispersion
     (section = "middle line")  Section to apply to two dimensional images
    (database = "database")     Database in which to record feature data
   (coordlist = "")             User coordinate list
       (units = "wn")           Coordinate units
        (nsum = "10")           Number of lines/columns/bands to sum in 2D imag
       (match = -3.)            Coordinate list matching limit
 (maxfeatures = 50)             Maximum number of features for automatic identi
      (zwidth = 100.)           Zoom graph width in user units
       (ftype = "absorption")   Feature type
      (fwidth = 12.)            Feature width in pixels
     (cradius = 5.)             Centering radius in pixels
   (threshold = 0.)             Feature threshold for centering
      (minsep = 2.)             Minimum pixel separation
    (function = "spline3")      Coordinate function
       (order = 1)              Order of coordinate function
      (sample = "*")            Coordinate sample regions
    (niterate = 0)              Rejection iterations
  (low_reject = 3.)             Lower rejection sigma
 (high_reject = 3.)             Upper rejection sigma
        (grow = 0.)             Rejection growing radius
   (autowrite = no)             Automatically write to database
    (graphics = "stdgraph")     Graphics output device
      (cursor = "")             Graphics cursor input
     (aidpars = "")             Automatic identification algorithm parameters
        (mode = "ql")           

    cl> identify

You will be presented with a plot of Dispim.  Find the first line on your
chart, place the vertical cursor on that line, press "m", and then type
in the wavenumber value.  Do this for each identified line in the chart.
Once you have identified the lines press "f" to fit a function.  You are
then presented with a new plot showing the residuals of the fit.  Look
at the fitting parameters at the top of the plot - the RMS error is in
wavenumbers.  This should be about 0.005 for these data.  You can vary the
fitting parameters if you wish.

       :func cheby         # change the function
       :order 3            # set the order to a quadratic
       f                   # redo the fit
       l                   # look at the fit in another way with the linear
                             term removed
       q                   # return to the IDENTIFY plot

Notice now that the scale has changed from pixel coordinates to wavenumbers.
You can save this solution by pressing "q" and answering the question
appropriately.  The solution will now be in the database directory with an
"id" prepended to the image name.  The image Dispim is now known to IRAF as 
a "reference" image - defined by the existence of "REFSPEC1= 'Dispim'" in its
image header that was put there by the IDENTIFY task.

    cl> dir database
    cl> page database/idDispim
    cl> imhead Dispim l+

We now need to identify the reference image for each image that we wish
to wavenumber calibrate - this is done with the task REFSPECTRA.  All this
task really does is add a keyword to the header of the input image that says
just what it does in the Dispim header (this can be more complicated but
for our case things are quite simple).  The "reference" parameter is set to
our complete image list, but there is only one reference image in that list
and that is Dispim.  If your own data set has many reference images then
look at the help page for REFSPECTRA to see the various options for
determining reference images and their weights for the data.

    cl> unlearn refspectra
    cl> lpar refspectra
    cl> refspectra *.ms.imh select=nearest sort="" group=""
    cl> hselect *.ms.imh $I,refspec1 yes

The next, and final, step for this exercise is to apply the solution that
was saved in the database file to your extracted 1-d data.  This is
done with the task DISPCOR.  DISPCOR will look in the header of each input
image to determine what to use for the calibration - that's the refspec1
keyword in our image headers.

    cl> unlearn dispcor
    cl> lpar dispcor
    cl> dispcor *.ms.imh *%.ms%d%.imh
    cl> imhead avir1d.imh l+

The data are now calibrated - the pixel values have been interpolated so that
there is a constant delta wavenumber - as set by the DISPCOR parameter
"linearize=yes".  Check things out with SPLOT.

If you wish to remove telluric lines from your spectra then look at the
task TELLURIC in the ONEDSPEC package - we reduced Alpha Virgo for that
purpose.  If you wish to use data in your own software you can convert
the IRAF images to FITS files with WFITS or write them out as text files
with WSPECTEXT.


This is a very brief introduction to the reductions of your Phoenix data.
Before you proceed with your own data reductions look at the help pages
for all the tasks involved so you understand your options - and in some
cases there are certainly lots of those.  Good luck!


References
----------

ATMOS line list:  Molecular line parameters for the atmospheric trace
  molecule spectroscopy experiment, Brown et al., Applied Optics 26, 
  December 1987

A Beginner's Guide to Using IRAF (IRAF Version 2.10), Jeannette Barnes, 
  August 1993 (ftp to iraf.noao.edu, cd pub, binary, get beguide.ps.Z)

Cryogenic Spectrometer User Manual August 1995 (from the KPNO Support Office)

Help pages for IRAF tasks (online with IRAF)


jvb 10/10/97
(with many helpful comments from Lloyd Wallace, NOAO)
