IRAF TUTORIALS:
EXTRACTION OF STELLAR LONG SLIT SPECTRA USING DOSLIT

Abstract

This tutorial and exercise illustrates flat fielding of long slit stellar spectra with response and the use of the reduction task doslit. As an exercise, artificial long slit stellar data will be created and used. The reader may then execute the steps illustrated here as well as explore the options. As a tutorial this document may be read alone and the sample output may be viewed. The version you are looking at has the graphical output as part of the document making it self-contained and more easily read as a tutorial but also slower to load. A version of this document which uses optional links to a Postscript viewer may be selected here.

Long slit spectra require creating a special flat field with the task response. which contains the normalized flat field response. After flat fielding the long slit data is reduced using the task doslit. This reduction task is specialized for background sky subtraction, extraction, wavelength calibration, and flux calibration of long slit stellar spectra. It is a command language script which collects and combines the functions and parameters of many general purpose tasks to provide a single complete data reduction path. The task provides a degree of guidance, automation, and record keeping necessary when dealing with many observations.


Introduction

This tutorial and exercise presents a simple, though complete, use of the response task for creating flat fields appropriate for long slit spectra and the special stellar slit spectra reduction task doslit. This task is a complex script that guides the user through the various steps of background sky subtraction, order extraction, wavelength calibration, and flux calibration of stellar long slit spectra. In addition to organizing and guiding the user through these steps the task attempts to provide additional record keeping such that only images and steps not yet performed will be done. This exercise assumes some knowledge of how to use IRAF such as starting the IRAF session, moving about directories, and editing parameters.

There is on-line help for these tasks as well as a typeset user's guide for doslit. These are essentially identical and differ only in visual appearance and where the parameters are described. To read the on-line help type

kp> phelp response
[output]
kp> phelp doslit
[output]
A printed copy of this help can be produced with the commands
kp> help response,doslit | lprint
where the lprint command will print to the default printer which may be changed using the hidden parameters. Finally, the typeset Postscript version of the reduction task may be obtained and read by clicking here

Another useful reference for the reduction of long slit spectra which covers doing the various steps separately or with doslit is the document A User's Guide to Reducing Slit Spectra with IRAF by P. Massey, F. Valdes, and J. Barnes.

Getting Started

Start an IRAF session and go to a directory in which to do this exercise; you could use your home directory if you wish. Now load the artdata, longslit and kpnoslit packages and unlearn all previous parameters as shown below (with the output not shown).

cl> noao
[list of tasks]
no> artdata
[list of tasks]
ar> longslit
[list of tasks]
ar> imred
[list of tasks]
im> kpnoslit
[list of tasks]
kp> unlearn artdata longslit kpnoslit
The artdata package is used to create the sample data we will use, the longslit package contains the special flat field normalization task, and the kpnoslit package includes the special processing task that will extract and calibrate the stellar long slit spectra.

To create some sample long slit data type the following commands. If the data already exists then nothing will happen.

kp> demos mkdoslit
Creating example longslit in image demoarc1 ...
Creating example longslit in image demoobj1 ...
Creating example longslit in image demostd1 ...
Creating example longslit in image demoarc2 ...
kp> mkexample longslit demoflat oseed=4 nseed=1
Creating example longslit in image demoflat ...

In version V2.11 the flat field will also be created by the first command but in previous versions the flat field has to be created with the mkexamples task. The data we have created consists of a flat field image demoflat, two arc images demoarc1 and demoarc2 taken before and after the object exposures, a bright standard star image demostd1, and an object image demoobj1. These images have already been processed to subtract the overscan bias level, trim the overscan region, and apply a zero level calibration.

Let's begin by looking at one of the images. We can use an image display, which requires that you have a display server running on your console (SAOimage or XImtool) and that your session has been setup to display to that window, and/or the plotting task implot.

kp> display demoobj1 1
[The image is shown rotated to save space]

kp> implot demoobj1

In the plotting task you exit by typing 'q' as is true for most graphics tasks. You may also wish to look at the other images and list out a header with imhead. You will see that this data has a spectrum running down the approximate center of the slit, that there are a number of night sky lines, and the slit profile drops off near the edge of the image. Sometimes the slit will not cover the entire image format and other times it will. Also the exercise images are only 100 by 512 pixels which is smaller than typical data. However, this sample data has all the basic characteristics of stellar long slit spectra.

Creating a Flat Field with RESPONSE

One cannot simply divide by the observed flat field image because the flat field image includes the shape of the flat field spectrum. Applying the flat field directly will distort the object spectra and, in particular, the pixel values will no longer be related to the actual observed signal. This means that it would not be possible to use algorithms that depend on using the detector gain and readout noise to estimate the uncertainties of the pixels from the pixel values. These algorithms are for optimally weighted extraction and cosmic ray detection.

Thus, the algorithm we need models the shape of the observed flat field spectrum and normalize it away. This is done by the response task.

The first step is to determine the region of useful data. We can do this with various image plotting tasks. Let's use implot.

kp> implot demoflat

We see the slit transmission falls off fairly strongly below pixel 20 and above pixel 80. So we will use the data between these two points. We could trim the images first with imcopy or ccdproc (which is the recommended way). But in this exercise will do it when we apply the flat field to the data.

We now run the response task to create a flat field for the long slit data.

lo> response
Longslit calibration images: demoflat
Normalization spectrum images: demoflat[20:80,*]
Response function images: Flat
Fit the normalization spectrum for demoflat interactively (yes):

The task takes the normalization spectrum and sums up each point across the dispersion; in this case the columns at each line. We used an image section to limit this to the region of the slit with good signal. The summed data is then plotted and a function is fit to this one dimenisonal spectrum. The function that is fit should represent the basic flat field source spectrum but not any flat field response effects. In particular, the function should not fit any fringing which may occur and which would appear as oscillations in the spectrum.

This function fitting step uses a standard IRAF tool called icfit for interactive curve fitting. Within this mode you can type '?' to get a list of the commands. Typically one would use only :function and :order to change the function and order of the fit, and 'f' to cause the fit to be redone after changing any fitting parameters. Feel free to experiment and when you are done and have an acceptable fit exit with 'q'.

The function representing the basic flat field source spectrum is then divided into each column (since the dispersion runs along columns) to form the normalized flat field. We can look at the flat field with implot.

kp> implot Flat

[output after typing c at column 50]

Note that the values are near one and that a plot along the central column shows little of the shape of the flat field spectrum.

The last step is to apply this flat field to the data. This could be done as part of the basic CCD processing done by ccdproc. In this simple exercise we will use imarith to divide the images by the normalized flat field. We will also trim the images by using image sections to remove the first and last columns.

kp> imarith demo*[20:80,*] / Flat[20:80,*] demo*
This takes the middle 61 columns and all lines and divides by the corresponding part of the flat field image and then replace the original data by the flat fielded data.

Running DOSLIT

In this section we will run doslit on the sample stellar long slit data. Since doslit combines many different operations the discussion below is broken up into the logical steps. The discussion will show the various prompts and interactions with commentary added. Sample output is also given.

In general the one dimensional spectra produced will have the same image name as the original image with the extension .ms. However for arc calibration images you will notice funny names will be used for the extracted spectra. The names are a concatenation of the object spectrum and the arc spectrum. This is because objects in the slit may be in different places and the ideal calibration uses arc data from the same pixels as the objects. So for each object the arc or arcs are extracted covering the same region of the slit occupied by the object and the output arcs are given the concatenated name. In generally you will not need to be concerned with this.

Parameters

To begin we set the task parameters. Below are the suggested parameters for this exercise.

kp> epar doslit
                                   I R A F
                    Image Reduction and Analysis Facility
PACKAGE = kpnoslit
   TASK = doslit

objects =             demoobj1  List of object spectra
(arcs   =    demoarc1,demoarc2) List of arc spectra
(arctabl=                     ) Arc assignment table (optional)
(standar=             demostd1) List of standard star spectra

(readnoi=              rdnoise) Read out noise sigma (photons)
(gain   =                 gain) Photon gain (photons/data number)
(dispaxi=          )_.dispaxis) Dispersion axis (1=along lines, 2=along columns)
(width  =                   5.) Width of profiles (pixels)

(dispcor=                  yes) Dispersion correct spectra?
(extcor =                  yes) Extinction correct spectra?
(fluxcal=                  yes) Flux calibrate spectra?
(resize =                  yes) Automatically resize apertures?
(clean  =                  yes) Detect and replace bad pixels?
(splot  =                  yes) Plot the final spectrum?
(redo   =                   no) Redo operations if previously done?
(update =                   no) Update spectra if cal data changes?
(quicklo=                   no) Minimally interactive quick-look?
(batch  =                   no) Extract objects in batch?
(listonl=                   no) List steps but don't process?

(sparams=                     ) Algorithm parameters
(mode   =                   ql)
The first set of parameters sets the names of the various types of images and the second set of parameters sets some properties of the image. (These parameters are for V2.10-V2.10.2 and in later versions the dispaxis parameter is gone and a datamax parameter is present). The next set of parameters are the processing steps to perform. In this exercise we select all the steps. Normally, doslit will only do each operation on a data image if it has not been done previously, as determined by looking at the images in the directory. The redo parameter can be set to force repeating previous processing.

There is a long set of algorithm parameters contained in the parameter set with the name sparams. In this exercise we will use the default parameters. Very few of these parameters generally need to be modified. If you are interested you can examine these parameters with epar or by typing :e when the cursor is over the sparams parameter.

A useful package parameter to set to give a more descriptive output of what happens during the processing is shown below.

kp> verbose=yes

Defining the Apertures

The first thing the task does is to find the position(s) of the stars in the images. For each image to be processed the task finds the brighest peak in the average of the 10 center lines (see apfind). A region around the peak containing most of the profile is then set using a resizing algorithm (see apresize). You then have the option to review the location, aperture width, and background regions. This should be checked for misdentifications, improper aperture size, and possible problems with contaminating objects.

kp> doslit
List of object spectra (demoobj1):
Searching aperture database ...
Finding apertures ...
Nov 21 14:18: FIND - 1 apertures found for demoobj1
Resizing apertures ...
Nov 21 14:18: APRESIZE  - 1 apertures resized for demoobj1 (-3.97, 3.51)
Edit apertures for demoobj1?  (yes):

The apedit aperture editing mode has many commands, which you can list with the '?' key, and allows you to do many things. For this data one would correct misidentifications and adjust aperture widths. The window can be zoomed in various ways such as with 'X' or 'w' followed by 'x' (see gtools). Since we are going to do background subtraction we should check the positioning of the background regions. Type 'b' on the central order to see the background regions.

The default regions are rather narrow in this case and there are no contaminating objects. The background region graphical interface uses the icfit graphical fitting mode which you have encountered previously and will encounter later. You can also use the 'w' windowing commands from gtools. To change the background regions initialize them with 't' and set them with the cursor by typing pairs of 's' (for sample regions). Or you can set them explicitly with

:sample -20.5:-7.5,7.5:20.5
[Now type 'f']

Note that the positions of the sample regions are relative to the order center and they will shift to follow the order position. Also shown is an example fit to the data in the sample regions. The background subtraction (done during the spectrum extraction) will fit the data in the background regions at each line. There are many options for the fitting. The default shown computes the median in each region and then fits a constant function (which is equivalent to the average of the two medians). You exit the background examination with 'q'. When you are satisfied with the apertures and background regions exit this stage with 'q'.

Fit traced positions for demoobj1 interactively?  (yes):
Tracing apertures...
Fit curve to aperture 1 of demoobj1 interactively?  (yes):

Nov 21 14:18: TRACE - 1 apertures traced in demoobj1.
Nov 21 14:18: DATABASE - 1 apertures for demoobj1 written to database

The position of the spectrum is traced, as described for the aptrace task, from the central lines to the rest of the image. This is done in steps of 10 lines with an average of 10 lines. Once the positions at this subsample of lines have been measured, a function is fit to the positions. This may be done interactively allowing you to adjust the fitting function and order of the function. The interactive graphical step again uses the icfit tool. For trace fitting typically one would use only :order to change the order of the fit, and 'f' to cause the fit to be redone after changing any fitting parameter. Feel free to experiment and when you are done and have an acceptable fit (one which fits the measured positions to fractions of a pixel except for noise points caused by the influence of cosmic rays) exit with 'q'.

These steps of finding, setting the background regions, and tracing the spectrum are repeated for each image. You have the choice of how interactive to make this. In this exercise the only other image is the standard star image. Repeat setting the background regions and checking the fit.

Searching aperture database ...
Finding apertures ...
Nov 21 14:18: FIND - 1 apertures found for demostd1
Resizing apertures ...
Nov 21 14:18: APRESIZE  - 1 apertures resized for demostd1 (-4.03, 3.26)
Edit apertures for demostd1?  (yes):
b
:sample -20.5:-7.5,7.5:20.5
f
q
q
Fit traced positions for demostd1 interactively?  (yes):
Tracing apertures ...
Fit curve to aperture 1 of demostd1 interactively  (yes): n
Nov 21 14:18: TRACE - 1 apertures traced in demostd1.
Nov 21 14:18: DATABASE - 1 apertures for demostd1 written to database

Dispersion Calibration

The doslit task now extracts, as described for the task apsum, the first helium-neon-argon arc spectrum specified in the arc spectra list. This will be the reference arc against which all arcs will be automatically reidentified using the same set of arc lines and dispersion function parameters. The extracted master arc is displayed (with the identify task which is quite powerful) and you must now identify a few lines set the basic dispersion scale against which other lines from the line list (here linelists$idhenear.dat) may be found. Mark the line at pixel position 131 by moving the cursor to this point and typeing 'm'. Respond to the query for the line wavlength with 5015 (if there is no query you did not have the cursor close enough to the line). Now mark the line at pixel position 500 and give the wavelength as 7281. The sample graph below shows how the arc identifications should look at this point.

Extract arc reference image demoarc1
Searching aperture database ...
Finding apertures ...
Nov 21 14:18: FIND - 1 apertures found for demoarc1
Nov 21 14:18: DATABASE - 1 apertures for demoarc1 written to database
Extracting apertures ...
Nov 21 14:18: EXTRACT - Aperture 1 from demoarc1 --> demoarc1.ms
Determine dispersion solution for demoarc1

After two or more lines have been marked (two in this example) we we can have the task find an initial dispersion solution and find additional lines from the line list. Type 'l' and you will see new lines added and a wavelength scale plotted.

Now let us look at the dispersion function by typing 'f'.

This is again the icfit mode and the same commands for adjusting the function, order, and other parameters apply. Note that two lines were misidentified and one has been automatically rejected from the fit as being too deviant (marked with the diamond symbol). So these lines don't affect any other arcs we may reidentify, use the cursor to point to them and type 'd'.

For dispersion function fitting a useful way to look at the fit is to subtract a linear component defined by the fitted wavelengths of the first and last pixels in the spectrum. This is done by typing 'l'.

Here we see a basically quadratic dependence (typical of gratings) which is well fit by the default cubic spline of one piece. You can experiment with the fit and when done exit with 'q'. You will again be in the line marking mode but displayed with the fitted dispersion function scale. For basic usage you can now exit this step with 'q'.

The next step is to define a linear dispersion scale to which all extracted spectra will be sampled. You can set the starting wavelength, the ending wavelength, the wavelength per pixel, and the number of pixels in the final spectrum. You are prompted with the defaults based on the dispersion function just fit. You answer yes to the query if you want to change the defaults.

demoarc1.ms.imh: ap = 1, w1 = 4204.1..., w2 = 7355.4..., dw = 6.16..., nw = 512
Change wavelength coordinate assignments? (yes|no|NO): n

Flux Calibration

Flux calibration consists of extracting all the standard star observations and comparing the observed fluxes against tabulated values for the star as given in the directory specified by the package caldir directory (in this case the Massey standards in onedstds$spec50cal/). The ratios of the observed flux over some bandpass to the tabulate value over the same bandpass are fit by a "sensitivity" function which ultimately provides the calibration.

First each standard star is extracted, dispersion calibrated, and the standard star name defined. The previously defined aperture definitions are read and the standard spectrum is extracted. Then the nearest two arc spectra in time are found from the list of arc spectra (in this case we only have two) and assigned weights based on the times. Each of these arcs is then extracted using the same apertures as the standard star. Here is where the concatenated image names appear. New dispersion functions are fit to the arcs using the same set of lines and function parameters as the master arc. This is done automatically by the reidentify task. The statistics of the new fit are printed and you have the option to look at the spectrum, line identifications, and fit, interactively. What is usually done is to monitor the RMS value and only review those in which the RMS is significantly larger than the orignal master dispersion solution. In this example you may look at the new solutions but they are generally fine. Finally the two arc dispersion solutions are interpolated with the weights and used to resample the standard star spectrum to the previously defined wavelength scale.

Extract standard star spectrum demostd1
Searching aperture database ...
Nov 21 14:19: DATABASE  - 1 apertures read for demostd1 from database
Extracting apertures ...
Nov 21 14:19: EXTRACT - Aperture 1 from demostd1 --> demostd1.ms
Assign arc spectra for demostd1
[demostd1] refspec1='demoarc1 0.403'
[demostd1] refspec2='demoarc2 0.597'
Extract and reidentify arc spectrum demoarc1
Searching aperture database ...
Nov 21 14:19: DATABASE  - 1 apertures read for demostd1 from database
Nov 21 14:19: DATABASE - 1 apertures for demoarc1 written to database
Extracting apertures ...
Nov 21 14:19: EXTRACT - Aperture 1 from demoarc1 --> demostd1demoarc1.ms

REIDENTIFY: NOAO/IRAF V2.10EXPORT valdes@puppis Mon 14:19:07 21-Nov-94
  Reference image = demoarc1.ms, New image = demostd1demoarc1.ms, Refit = yes
          Image Data    Found     Fit Pix Shift  User Shift  Z Shift      RMS
demostd1demoarc1.ms - Ap 1 48/48   48/48    3.24E-4     6.92E-4  2.25E-7    0.27
3
Fit dispersion function interactively? (no|yes|NO|YES) (yes): n
Extract and reidentify arc spectrum demoarc2
Searching aperture database ...
Nov 21 14:19: DATABASE  - 1 apertures read for demostd1 from database
Nov 21 14:19: DATABASE - 1 apertures for demoarc2 written to database
Extracting apertures ...
Nov 21 14:19: EXTRACT - Aperture 1 from demoarc2 --> demostd1demoarc2.ms

REIDENTIFY: NOAO/IRAF V2.10EXPORT valdes@puppis Mon 14:19:12 21-Nov-94
  Reference image = demoarc1.ms, New image = demostd1demoarc2.ms, Refit = yes
          Image Data    Found     Fit Pix Shift  User Shift  Z Shift      RMS
demostd1demoarc2.ms - Ap 1 48/48   48/48    0.00302      0.0174  3.64E-6    0.28
1
Fit dispersion function interactively? (no|yes|NO|YES) (no):
Dispersion correct demostd1
demostd1.ms.imh: ap = 1, w1 = 4204.107, w2 = 7355.405, dw = 6.166924, nw = 512

The query for the standard star name prints out the image title to help you remember which standard star was observed. You then give the standard star name using one of the names known by the database of standard stars. You can try typing '?' to get a list of names. In this example we respond with hz44. The task used is standard. The ratios of the measured standard star fluxes to the tabulated fluxes are written to the file std.

Compile standard star fluxes for demostd1
demostd1.ms.imh[1]: Example artificial long slit image
Star name in calibration list: hz44

The task sensfunc is used to fit the sensitivity function. This task has many options for fitting, deleting points or stars, determining residual extinctions from the default extinction function, and displaying. In this tutorial/exercise we only have a single star so there isn't too much to do. You can experiment with this task and quit with 'q' when you are done. You have the option to look at the fits or use the last set of fitting parameters.

Compute sensitivity function
Fit aperture 1 interactively? (no|yes|NO|YES) (no|yes|NO|YES) (yes):

Sensitivity function for all apertures --> sens
Flux and/or extinction calibrate standard stars
[demostd1.ms.imh][1]: Example artificial long slit image
  Extinction correction applied
  Flux calibration applied
Standard stars:
Splot spectrum? (no|yes|NO|YES): n

The output of the sensitivity function fitting is a sensitivity function spectrum which have the name sens. The standard star spectra are then flux calibrated by applying an extinction correction based on the air mass of the observation and the extinction function given by the extinction package parameter (which is onedstds$kpnoextinct.dat in this demonstration) and the sensitivity function just derived. You are then given the opportunity to review the flux calibrated standard star spectra with the general spectrum plotting task splot.

Extraction of the Object Spectra

Having now set the dispersion calibration and the flux calibration, doslit will extract and calibrate each object spectrum with little interaction required. The objects are extracted with the apertures defined earlier. Then one or more arcs from the arc list are assigned to the object. The arc spectra are extracted using the apertures for the object. The arc lines are automatically reidentified and a new dispersion function is fit based on the same lines and function parameters used with the master arc. The object spectrum is then resampled to linear dispersion. Finally the object spectrum is extinction and flux calibrated.

Extract object spectrum demoobj1
Searching aperture database ...
Nov 21 14:20: DATABASE  - 1 apertures read for demoobj1 from database
Extracting apertures ...
Nov 21 14:20: EXTRACT - Aperture 1 from demoobj1 --> demoobj1.ms
Assign arc spectra for demoobj1
[demoobj1] refspec1='demoarc1 0.403'
[demoobj1] refspec2='demoarc2 0.597'
Extract and reidentify arc spectrum demoarc1
Searching aperture database ...
Nov 21 14:20: DATABASE  - 1 apertures read for demoobj1 from database
Nov 21 14:20: DATABASE - 1 apertures for demoarc1 written to database
Extracting apertures ...
Nov 21 14:20: EXTRACT - Aperture 1 from demoarc1 --> demoobj1demoarc1.ms

REIDENTIFY: NOAO/IRAF V2.10EXPORT valdes@puppis Mon 14:20:22 21-Nov-94
  Reference image = demoarc1.ms, New image = demoobj1demoarc1.ms, Refit = yes
          Image Data    Found     Fit Pix Shift  User Shift  Z Shift      RMS
demoobj1demoarc1.ms - Ap 1 48/48   48/48   -1.31E-4    -0.00208  -1.7E-7    0.27
2
Fit dispersion function interactively? (no|yes|NO|YES) (no):
Extract and reidentify arc spectrum demoarc2
Searching aperture database ...
Nov 21 14:20: DATABASE  - 1 apertures read for demoobj1 from database
Nov 21 14:20: DATABASE - 1 apertures for demoarc2 written to database
Extracting apertures ...
Nov 21 14:20: EXTRACT - Aperture 1 from demoarc2 --> demoobj1demoarc2.ms

REIDENTIFY: NOAO/IRAF V2.10EXPORT valdes@puppis Mon 14:20:48 21-Nov-94
  Reference image = demoarc1.ms, New image = demoobj1demoarc2.ms, Refit = yes
          Image Data    Found     Fit Pix Shift  User Shift  Z Shift      RMS
demoobj1demoarc2.ms - Ap 1 48/48   48/48    0.00239      0.0135  2.95E-6    0.27
9
Fit dispersion function interactively? (no|yes|NO|YES) (no):
Dispersion correct demoobj1
demoobj1.ms.imh: ap = 1, w1 = 4204.107, w2 = 7355.405, dw = 6.166924, nw = 512
Extinction correct demoobj1
Flux calibrate demoobj1
[demoobj1.ms.imh][1]: Example artificial long slit image
  Extinction correction applied
  Flux calibration applied

One of the options when you run doslit is to allow you to look at each spectrum after it has been extracted using the commonly used task splot. You may answer yes, no or the capitalized forms to suppress other queries. Exit with 'q' as usual.

demoobj1.ms.imh:
Splot spectrum? (no|yes|NO|YES) (no): y

The extraction and calibration of all the objects, in this case only one, continues until done. That concludes doslit.

Additional Extractions

Unless the redo parameter is set to yes, any new object spectra specified in additional executions will be processed largely automatically Also any spectra which have already been processed will be ignored. For instance try

kp> doslit
List of object spectra (demoobj1): 

That is all for this tutorial exercise. You may run through it again, look at the logfile and standard star file std with page, and look at any of the extracted spectra and sensitivity function with splot.

References

The two main tasks illustrated here, response and doslit have on-line help pages. There is also a typeset Postscript version of doslit.

There is a general user's guide to long slit reductions, A User's Guide to Reducing Slit Spectra with IRAF by P. Massey, F. Valdes, and J. Barnes.

There are on-line help pages for the various tasks used by the doslit script. These tasks are: icfit, gtools, response, apfind, apresize, apedit, aptrace, apsum, identify, reidentify, standard, sensfunc, splot

There are some other IRAF exercises, including one on long slit data reductions using real data. These are obtained from the anonymous ftp account at iraf.noao.edu in the directory iraf/misc. You may retrieve the V2.10.2 version or V2.10.3 version here.


Francisco Valdes
November, 1994