Reduction of Images with PyRAF¶

This tutorial will use observations from program GS-2006B-Q-18 (PI: J. Arias): narrow- and broad-band imaging of three adjacent regions in M8 (or NGC 6523, the Lagoon Nebula). The science goals likely included studying the ionization and small-scale structures in this star forming region.

The observations were obtained in Queue mode during several nights between 2006 Sep. 10 and Oct. 10. The science images were obtained in 5 filters: Ha, HaC, SII, r, and i. (See the description of the filters.) The exposures were obtained with \(2\times2\) binning (or about \(0.16\times0.16\) arcsec, which still samples the PSF very well). Contemporaneous dayCal exposures were obtained including Twilight flats in all the filters, but not including photometric standards.

For other tutorials, see the following links:

Retrieve & Organize Data¶

The first step is to retrieve the data from the Gemini Observatory Archive (see Archive Searches). You may search the GOA yourself, or instead just cut-and-paste the direct URL:

# M8 imaging data:
https://archive.gemini.edu/searchform//sr=600/GS-2006B-Q-18/notengineering/GMOS-S/imaging/NotFail/present/canonical

After retrieving the science data, click the Load Associated Calibrations tab on the search results page and download the associated bias and flat-field exposures. Unpack all of them in a subdirectory of your working directory named /raw. Be sure to uncompress the files. See Retrieving Data for details.

Processing Preparation¶

Reference Files¶

Two of the required MasterCals will be created in this tutorial.

Bias Residual (created below)
Flat-field for each filter (created below)
Bad Pixel Mask (download: bpm_gmos-s_EEV_v1_2x2_img_MEF.fits.bz2)

After downloading, uncompress the BPM in your work directory for use in the reductions.

Software¶

You must create an observing log database of the data in the ./raw subdirectory. Download: obslog.py to the ./raw subdirectory, and execute it from the unix prompt.

python obslog.py obsLog.sqlite3

See Creating an Observing Log for details.

Also retrieve the python file selection module, which includes template SQL statements for selecting files, and functions for specifying metadata on which to perform selections.

Download: fileSelect.py

Place this module in your work directory; it is used by the reduction script (below). You can perform all of the processing steps for this tutorial by downloading the Imaging Tutorial python script.

Download: gmos_img_proc.py

Place the script in the work directory, and execute the reduction script either from the Unix prompt,

python gmos_img_proc.py

or from an active PyRAF session:

from pyraf import iraf
from pyraf.iraf import gemini, gemtools, gmos
import fileSelect as fs
import gmos_img_proc

You may find it useful to download the script to follow this tutorial in detail, and use it as the basis for reducing other imaging observations.

Caution

This script may use a large amount of available storage if Google Drive is connected to your desktop. Disable synching for the working directory before running the script or a hidden directory named ‘.tmp.driveupload’, which is in the same directory as ‘vm_transfer.’, can take up 200+ GB of space.

Building MasterCals¶

The next steps will create the necessary MasterCal reference files that are used to calibrate the science files. Files are selected by matching specific exposure metadata in the observing log database (see GMOS Processing Keywords). Within the PyRAF session, first create the Bias Residual MasterCal:

# Observing log database
dbFile='./raw/obsLog.sqlite3'

# Create the query dictionary of essential parameter=value pairs.
qd = {'use_me':1,
      'Instrument':'GMOS-S','CcdBin':'2 2','RoI':'Full','Object':'M8-%',
      'DateObs':'2006-09-01:2006-10-30'
      }
# Now prepare to call the IRAF processing task.
# Set the task parameters.
gmos.gbias.unlearn()
biasFlags = {
    'logfile':'biasLog.txt','rawpath':'./raw/','fl_vardq':'yes',
    'verbose':'no'
}
# The following SQL generates the list of files to process.
SQL = fs.createQuery('bias', qd)
biasFiles = fs.fileListQuery(dbFile, SQL, qd)

# The str.join() function is needed to transform a python list into a
# comma-separated list of files that IRAF can understand.
if len(biasFiles) > 1:
    gmos.gbias(','.join(str(x) for x in biasFiles), 'MCbias.fits',
        **biasFlags)

Now create the Flat-field MasterCal files from the twilight images.

# Select flats obtained contemporaneously with the observations.
qd.update({'DateObs':'2006-09-10:2006-10-10'})
# Set the task parameters.
gmos.giflat.unlearn()
flatFlags = {
     'fl_scale':'yes','sctype':'mean','fl_vardq':'yes',
     'rawpath':'./raw/','logfile':'giflatLog.txt','verbose':'no'
    }
filters = ['Ha', 'HaC', 'SII', 'r', 'i']
for f in filters:
    # Select filter name using a substring of the full designation.
    qd["Filter2"] = f + '_G%'
    mcName = 'MCflat_%s.fits' % (f)
    flatFiles = fs.fileListQuery(dbFile, fs.createQuery('twiFlat', qd), qd)
    if len(flatFiles) > 0:
        gmos.giflat(','.join(str(x) for x in flatFiles), mcName,
            bias='MCbias', **flatFlags)

Science Processing¶

The gireduce task has more than 50 parameters; the table below lists the defaults for the processing flag keywords—i.e., the keywords with logical values to indicate whether to perform an operation. For the most part you can use the default parameter values; exceptions are noted explicitly in the code blocks below.

**gireduce Processing Flag Defaults**¶
Flag	Default	Description
`fl_addmdf`	No	Append MDF extension? Not applicable to imaging.
`fl_bias`	Yes	Subtract bias residual?
`fl_dark`	No	Subtract scaled dark image?
`fl_flat`	Yes	Apply flat-field correction?
`fl_mult`	Yes	Multiply by the CCD gains?
`fl_inter`	No	Fit overscan levels interactively?
`fl_over`	Yes	Perform overscan correction?
`fl_qecorr`	No	Apply QE correction?
`fl_trim`	Yes	Trim overscan region?
`fl_vardq`	No	Propagate VAR and DQ?

Turning now to the science reductions, we first remove the restriction on the range of dates and set some task parameters:

# Set task parameters.
# Employ the imaging Static BPM for this set of detectors.
gmos.gireduce.unlearn()
sciFlags = {
    'fl_over':'yes','fl_trim':'yes','fl_bias':'yes','fl_dark':'no',
    'fl_flat':'yes','logfile':'gireduceLog.txt','rawpath':'./raw/',
    'fl_vardq':'yes','bpm':'bpm_gmos-s_EEV_v1_2x2_img_MEF.fits','verbose':'no'
    }
gemtools.gemextn.unlearn()    # disarms a bug in gmosaic
gmos.gmosaic.unlearn()
mosaicFlags = {
    'fl_paste':'no','fl_fixpix':'no','fl_clean':'yes','geointer':'nearest',
    'logfile':'gmosaicLog.txt','fl_vardq':'yes','fl_fulldq':'yes','verbose':'no'
    }

After gireduce, the next step is to mosaic the multiple image extensions in each file into one extension. The gmosaic task will transform the image extensions, taking into account the relative offsets and orientations of the CCDs, and generate a single (interpolated) image.

Since the gmosaic task does not support file lists especially well, we perform that step in a loop immediately following basic reductions.

Caution

It is important to use the gmosaic task to combine the extensions, as information about the relative position and orientation of the CCDs in the FPA is hard-coded within the task. There are other options beyond imcoadd (noted below) to combine exposures once they are mosaiced.

qd["DateObs"] = '*'
prefix = 'rg'
for f in filters:
    qd["Filter2"] = f + '_G%'
    flatFile = 'MCflat_' + f + '.fits'
    sciFiles = fs.fileListQuery(dbFile, fs.createQuery('sciImg', qd), qd)
    if len(sciFiles) > 0:
        gmos.gireduce (','.join(str(x) for x in sciFiles), bias='MCbias',
                       flat1=flatFile, **sciFlags)
        for file in sciFiles:
            gmos.gmosaic (prefix+file, **mosaicFlags)

Examine the output files (now prefixed with mrg) to assess data quality, and adjust the processing or the lists of input files if necessary.

Warning

Though one might ordinarily use a file list of all science files and set each filter in gireduce with flat1=,... flat4=, this will likely result in an unhelpful error message if the set of input images includes more than 4 filters (as it does in this tutorial). So in the above code we process the exposures in a given filter separately.

Fringe Correction¶

A fringe pattern may be apparent in the extreme red (i-band) for some CCDs. This pattern is ordinarily derived from many low-background science images (with very few extended targets) by masking the targets and combining many dithered (preferably, non-overlapping) images that have been flat-fielded. The gmos.gifringe task would be used to construct the fringe frame, and the gmos.girmfinge task would scale and remove the pattern from science images. Since this science program obtained images of a single, extremely extended emission line region, the genuine astrophysical intensity variations completely overwhelm the fringe signature. The purist may wish to derive a fringe frame from another science program (but with the same CCDs), obtained under similar conditions of sky background. Such an approach is beyond the scope of this tutorial, however.

This marks the end of calibration processing for GMOS images; what to do next depends on your specific science goals. Other tools in the gemini package, or other third-party tools may be essential for your analysis.

Advanced Processing¶

Image Stacking¶

If there is more then one overlapping exposure of the same field taken with the same filter you may wish to stack them to, e.g., create a deep image for source identification or, as in this case, to create a mosaic of offset exposures to cover an extended target. This can be accomplished in part with the gemtools.imcoadd task, which handles DQ masking and outlier rejection. It does have some limitations, however:

It takes a lot of patience and trial-and-error tweaking of parameters to get good results
There is little control over sky background
The output image is no bigger than the first (reference) image, rather than the union of the image footprints

The best procedure for imcoadd stacking depends upon the details of the program and the observing sequence. In this program, M8 was observed in 3 marginally overlapping positions, and several dithered exposures were obtained at each position. The complex background requires that the dithered images at each position be combined first; tiling the images from the 3 positions requires subsequent processing.

# Use primarily the default task parameters.
gemtools.imcoadd.unlearn()
coaddFlags = {
    'fwhm':3,'datamax':6.e4,'geointer':'nearest','logfile':'imcoaddLog.txt'
    }
targets = ['M8-1', 'M8-2', 'M8-3']
prefix = 'mrg'
for f in filters:
    qd['Filter2'] = f + '_G%'
    for t in targets:
        qd['Object'] = t + '%'
        outImage = t + '_' + f + '.fits'
        coAddFiles = fs.fileListQuery(dbFile, fs.createQuery('sciImg', qd),
                                      qd)
        gemtools.imcoadd(','.join(prefix+str(x) for x in coAddFiles),
                 outimage=outImage, **coaddFlags)

Warning

Images with a complex background, such as narrow-band images of a region of extended nebulosity, may not combine well. This is in large part because imcoadd uses a very simple characterization of the background, and confuses it with genuine astrophysical background. Unfortunately the statsec parameter (where the sky is sampled) is specified in image pixels, rather than world (Ra, Dec) coordinates, so that regions of low astrophysical background cannot be specified.

Inspect the output images for quality; tweak the task parameters as necessary for improved results. These co-added images may be tiled after the background levels have been brought into agreement.

Photometric Calibration¶

Photometric calibration for images is beyond the scope of this Cookbook, and there no standard way of recording the calibration within the images. However the mechanism in its simplest form is straightforward: Use a photometry program such as SExtractor to measure instrumental magnitudes in each passband for detected stars, and determine the photometric zero-point for each image. Complications arise for very crowded fields (where stellar PSFs overlap), and PSF shapes that vary over the FoV.

If you are interested in deriving photometric zero-points and color terms, you may find it helpful to consult the paper by Jorgensen 2009, [CP] where these quantities are derived for the broad-band SDSS filters over a two-year period early in the life of GMOS.