The MultiDrizzle Handbook


5.5 Optimizing the Drizzle Parameters for Your Data

5.5.1 Creating Custom Association Tables

Association tables for dithered, REPEAT-OBS, or CR-SPLIT observations are generated by the calibration pipeline. These tables keep track of the input exposure filenames and the output product filenames. Some types of observations, however, will not have association tables. Others will have multiple tables from different visits which need to be combined into one. In the following discussion, we present the methodology for creating custom association tables, either by merging association tables or creating them from scratch. Users also have the capability to manually edit association tables to include any known image offsets. This can be done using the ttools task tedit in IRAF, where the user adds the columns XOFFSET, YOFFSET and/or ROTATION to the association table. Alternately, an ascii file with the known image offsets may be specified and the association automatically updated. Merging Association Tables

Observing programs which cover a large portion of the sky will generally be made up of multiple pointings. Each of these pointings may be dithered or split for cosmic ray rejection and will possess their own association table. In order for PyDrizzle to produce a single combined product, a composite association table must be built from the separate association tables. Users can easily create this composite by merging the individual tables using the ttools task tmerge with "option = append".

The default product rootname will be taken from the first listed DTH-PROD in the composite association table. This rootname can be modified by the user (with tedit) to suit the desired naming convention. Generally, this rootname should match the rootname of the composite association table. A detailed example of merging association tables is given in step 1 of Example 4 from the ACS Data Handbook, Section 3.5.2. Creating Association Tables from Scratch

In some cases, observations of the same target will not have accompanying association tables, for example observations taken in different visits or dithered observations taken using POS-TARG offsets. The PyRAF task buildasn (Section 5.5.2) has been developed for the dither package to create association tables which may be used for reprocessing with PyDrizzle or MultiDrizzle.

The following illustrates how to use Buildasn (Section 5.5.2) to create an association table from all calibrated input files with the FLT suffix found in a single directory. Within PyRAF and using the IRAF syntax, we type:

  pyraf> buildasn mymosaic suffix='flt.fits'

Alternately, the same result may be obtained using the Python syntax:

  pyraf> iraf.buildasn('mymosaic',suffix='flt.fits')

The association table will have the name `mymosaic_asn.fits', the product name will be `mymosaic_drz.fits', and all files with the FLT suffix found in the working directory will be included in the association. To specify a subset of images within a given directory, the user may specify the `suffix' parameter to be a filelist (`@filename'), a wild-card list (`*8cd*flt*'), or any portion of the filename (`crj' or `f555w').

If user determined offsets are available, buildasn (Section 5.5.2) has the capability of incorporating them into the association table. These offsets (XOFFSET, YOFFSET, and/or ROTATION) are specified by the file listed in the "shiftfile" parameter.This option allows users to fine-tune the final image combination by providing corrections to the header WCS information, if required.

  pyraf> buildasn mymosaic suffix='flt' shiftfile='shift.txt'

MultiDrizzle utilizes the buildasn (Section 5.5.2) task and creates association tables automatically when a file suffix (i.e. `flt.fits') or a file list and/or a shift file is specified in the input parameters.

5.5.2 Buildasn Parameters

The PyRAF task buildasn was written to help users create an association table for data which was not originally considered part of an association. It is a member of the IRAF dither package and requires the following inputs:

Table 5.9: Parameters for Buildasn
Rootname for the ASN table
input file description for association
members, may contain wildcards and "@" list
optional file containing user defined
shifts in order of x,y,rotation
output filename for the ASN table
Print additional information while
building the table?
boolean [False,True] Example Usage

In order to create a new association table, you must first make sure that the correct packages are loaded. In the Python shell (this will work for both PyRAF and a regular Python shell) you can use the asnUtil package in one of three ways:

 python>from glob import glob 
python>from pytools import asnutil

 and then:

Or from the PyRAF shell:

 pyraf>buildasn mymosiac suffix='flt'

both of these command sets build an association table called mymosaic_asn.fits using all the *_flt.fits images that are available in the current directory. If you were to input this asn file into PyDrizzle or MultiDrizzle, then the output product would be called mymosaic_drz.fits. MultiDrizzle can take as input any combination of asn tables, filenames or wildcard selection criteria and will use the asnutil package to generate the association table to feed to PyDrizzle. Shift Files

The shift file is an ASCII file containing the user-computed shifts for a set of images and is used for updating the offsets in the association table. Observers may use either the association (ASN) table delivered by the pipeline or create their own. The shift file can contain shifts for files other than those contained in the ASN table you wish to update, only those entries found in the ASN table will be used to update the table. Thus, a single shift file can be generated for a whole set of observations which are represented by multiple association tables. For more detailed information on how to calculate shifts for your images, please see the section on Alignment (Section 5.5.3).

The shift file uses the following format:

  # units: pixels             
 # frame: input              
 # form: absolute            
 # refimage: <filename>      
 filename1  xshift1 yshift1  [rotation1  [scale1]] 
 filename2  xshift2 yshift2  [rotation2  [scale2]]

The refimage filename requires the specification of the extension from which the WCS information will be read, but only in the case where the refimage is a multi-extension FITS image.

The specification of the filenames for the input images need to include the full filename, with suffix and extension; for example, j94f05bgq_flt.fits.

The entries for rotation and scale for each input image are optional. However, if a scale difference is noted for an image, a value must be given for the rotation as well.

The units for the items specified in the first three lines are:

    units                    Values: pixels (default), arcseconds
   frame                    Values: input (default), output
   form                     Values: absolute (default), delta  

The first three lines specify the units of the shifts; the frame of reference, whether the shifts were computed using the input images (input) or singly drizzled images (output); the form of the shifts, whether they should be added to the header information as additional incremental shifts (delta) or the explicit shift that should be applied in total (absolute). The `#' symbol is required to signify that these are commentary lines in the file. When shifts are given in the `output' frame of reference, it is necessary to specify the name of the reference image used. This reference image contains a copy of the WCS information from the image that was used as a reference to calculate the shifts and is only required when `frame = output'. If shifts were calculated from single drizzled images which have been distortion corrected, rotated or altered in scale then this image must be one of the drizzled files so that the corrections are taken into account when the shifts are added to the original WCS information.


Absolute Shifts

Figure 5.13: Absolute Shifts in the Input and Output Reference Frames

Absolute Shifts in the `input' and `output' frame of reference. The direction of north is indicated by the bold arrows.

Figure 5.13 shows an absolute shift which is the total shift between two images and includes offsets implied by the header WCS. Absolute shifts may be determined by separately drizzling images to their own unique WCS frame. PyDrizzle will use the header information to optimally place the drizzled image in the frame. If absolute offsets are provided in the shift file, PyDrizzle will populate the OFFSET columns in the association and fill the DELTA values with zeroes.

To derive absolute shifts in the `input' frame of reference, each distorted image is drizzled separately, and PyDrizzle uses the unique WCS information from each frame to choose the central RA/Dec, the image orientation, and final image size. No additional rotation or scaling is applied while drizzling. (Alternately, the user may catalog sources in each distorted image and then apply the distortion model to the X/Y positions. However, this requires tasks to transform coordinates from distorted to undistorted space)

Absolute shifts in the `output' frame of reference may be computed by separately drizzling each distorted image, with the same rotation and scaling, onto a unique WCS frame. Target lists in each drizzled image may then be matched to derive the absolute shifts.


Delta Shifts

Figure 5.14: Delta Shifts in the Input and Output Reference Frames

Delta Shifts in the `input' and `output' frame of reference.

Figure 5.14 shows a delta shift defined as the residual shift required to align sources, after applying the offsets implied by the header WCS. Delta shifts may be determined by separately drizzling images onto a common WCS frame with the same central right ascension and declination. This is performed in the `driz_separate' step of MultiDrizzle (Section 5.4) or by setting the PyDrizzle parameter `single=yes'. Object lists derived for each drizzled image may then be matched and fit to derive a single delta shift for each image.

Delta shifts are defined to be in the `input' frame of reference when the distortion-corrected, drizzled image has the same orientation and scale as the distorted image. Shifts will be in the `output' frame of reference if any rotation or scaling was applied while drizzling. When the shifts are given in the `output' frame of reference, PyDrizzle requires the specification of the reference image to properly account for any orientation and scale changes when interpreting the shifts. It is recommended that shifts be specified in the output frame with a reference image to insure the proper interpretation of the shifts.

If delta shifts are provided, those values will be used to populate the DELTA columns and the OFFSET columns will be zero. If shifts are given in both the OFFSET and DELTA columns, the OFFSET column will be applied and the DELTA columns will be ignored. This convention used by buildasn eliminates any ambiguity as to which values are applied to the data. Editing the ASN table directly will then allow the user to further update the offsets using either delta or absolute shifts.

It should be noted that shift files can generally be provided directly to MultiDrizzle without having to edit the association tables directly. In fact, MultiDrizzle allows the user to specify the input files without the need for an association table at all, and it will create one as needed on its own for use in combining the images.

The following is an example of a shift file in units of pixels, measured in the `input' frame of reference, with the form specified as `delta' shifts for an ACS image set. An x and y shift and a rotation (in degrees) is specified for each dataset. If only a simple shift is required, then the rotation column need not be specified.

  #units: pixels
 #frame: input
 #form: delta
 j8c0b1skq_flt.fits  -2.4 -1.2 -0.002
 j8c0b1snq_flt.fits  -4.3 -2.3  0.001
 j8c0b1sqq_flt.fits  -6.9 -3.5  0.003

This is what the output ASN table from buildasn should contain if the above shift file was used as input:

Table 5.10: Output ASN Table from buildasn

5.5.3 Alignment Alignment Error Sources

Users are encouraged to reference the section on HST Pointing Accuracy and Stability (Section 4.1.1) for a discussion on error sources that derive from the stability of the telescope. Additional errors may arise from the nature of the science data itself such as: Visit Alignment

The astrometry of images taken within the same visit and orbit of the telescope are generally limited to the accuracy of the telescope pointing which is controlled by the Fine Guidance Sensors (FGS) and the Guide Star Catalog. If you have external knowledge of the field (for example other ground or space based images that the current dataset can be linked to) then you can improve the inherent absolute astrometry of your image. If your images are populated with sufficient well exposed point sources, then the relative astrometry of each image can be improved through source detection and centroiding. These updates can then be folded back in to the astrometry information in the header of each image and used to combine all the exposures into a single well aligned and drizzled mosaic. Inter-visit alignment

If you have more than one visit for your entire dataset, enough time has probably elapsed that the telescope has used different guide stars for the observations. This could result in alignment errors of many pixels depending on the platescale of the instrument that was used. There are two paths you can take to get the best alignment for the final mosaic. If you cannot improve the inherent astrometry of the images inside a single visit, then you can calculate shifts for all the images at the same time using one image from the entire dataset as a reference, as described in the previous section for visit alignment. You might choose the exposure with the longest observing time (this would provide the best SNR and number of available objects for matching) or some other image based on the requirements of the dataset.

If you have improved the alignment of the images inside each visit, you can then calculate the inter-visit offsets between the combined visit level mosaics by following these proscribed steps: Combining large mosaics or data from multiple programs

The same methods are used to drizzle and align large mosaics and multiple programs as smaller ones. Every effort has been taken to ensure the drizzle algorithms are structured to provide the fastest computation and memory management. However, the user should consider limitations which exist due the size of the data, and the amount of memory available on the processing computer. Processing large datasets on a 32-bit system will be limited by the OS restrictions for addressing memory, so only 2Gb of memory can be used for all the output arrays. This results in a limiting size of about 16000 x 16000 for the combined image, if no context image is generated. This limitation can be avoided by running the MultiDrizzle code on a 64-bit system where the OS can address much larger blocks of memory. User Defined Shifts

The ability of PyDrizzle and MultiDrizzle to properly align images depends on the accuracy of the World Coordinate System (WCS) and any other keywords that may deal with astrometry or distortion information in the header of the image. Unfortunately this does not always provide the most precise alignment, even for images taken within the same visit. Incorrect shifts between exposures can degrade image quality and corrupt cosmic ray rejection when left uncorrected. To resolve this, users can specify additional shifts that can be added to the header calculation. By deriving residual (delta) shifts based on the position of objects in the data after aligning the images based on the WCS information, the user may refine the alignment based on the WCS header information to create a precisely-aligned drizzled product, especially for images taken in multiple-visits. Many observers have developed their own methods of comparing images and computing offsets. The conventions described here should support the majority of users. In practice, the shifts are applied to the data using a shiftfile (Section 5.5.2) which can be incorporated into the association table for the dataset using the buildasn function (Section 5.5.2).

Deriving these shifts generally follows the same set of steps; namely,

Use of IRAF tasks for the fit, such as geomap, assume that all rotations occur about (0.,0.), which by default falls on the corner of each image. However, the drizzle algorithm rotates all input images around the center of the image, not the corner. This requires that the positions for all sources need to be specified relative to the center of the output frame instead of the corner. This can be done by simply subtracting the coordinates of the output frame's center pixel from all the sources position prior to performing the fit. If this is not done, a residual offset will be introduced after applying the computed fit. The task tweakshifts, as described in the next section, demonstrates how to determine the shifts using IRAF based tasks while taking into account the correct center of rotation when performing the fit between source positions.

These shifts, regardless of how they were computed, should describe the offset which needs to be applied to the image to shift the objects so that they align with the same objects in the reference frame. Tweakshifts

Tweakshifts provides an automated interface for computing residual shifts between input exposures being combined using MultiDrizzle or PyDrizzle.The shifts computed by Tweakshifts correspond to pointing differences after applying the WCS information from the input image's headers. Such errors would, for example, be due to errors in guide-star positions when combining observations from different observing visits or from slight offsets introduced upon re-acquiring the guide stars in a slightly different position. This task was written using Python and relies on the Python interface to the IRAF environment, PyRAF, for its operation. As a result, this task can only be run under PyRAF. The primary implementation of tweakshifts involves using existing IRAF tasks such as daofind, geomap, and crossdriz to compute the shifts. Tweakshifts automates the operation of these tasks to produce a shiftfile and reference WCS file which can serve as input to MultiDrizzle. The shiftfile will contain the residual shifts for each input image, while the reference WCS file consists of a FITS header with all the WCS keywords appropriate for the reference image without any image data. The reference WCS will then be used by MultiDrizzle to properly interpret the shifts from the shiftfile for the final alignment of the images.

Tweakshifts supports multiple methods for determining the residual shifts between images: primarily, catalog matching and cross-correlation. Each mode of operation has its own requirements for input file formats in order to properly compute the shifts between the images.


Catalog Matching

A widely utilized method for computing offsets between images consists of identifying sources in each image, matching them up and solving the matched sets of positions to find offsets. This technique requires that the image contain recognizable sources which are point-sources, or similar enough to point-sources so that they can be identified by the software as a source. In addition, there has to be enough overlap between the sources identified from each input exposure to positively identify the same target in each image and allow the production of a sorted matched list of targets. Tweakshifts relies on either the daofind task or Source Extractor as the object identification routine. The first input image gets selected as the reference image for the final product. This results in each of the remaining image's source lists being matched to the reference image source list using xyxymatch producing an output matched list for each image relative to the reference image positions. These matched lists will only contain sources that are found in common between each image and the reference image. The matched lists will then be used as input to geomap to compute the final shift, rotation and scale change for each input image relative to the reference image.

Normally, distortion-free input images would be required as input in order to allow positive identification of sources relative to the reference image. For ACS, this would require the use of MultiDrizzle (or PyDrizzle) to generate distortion-corrected images, such as the singly-drizzled images from MultiDrizzle. However, tweakshifts supports the use of calibrated, distorted images (such as FLT images for ACS, or.c0h images for WFPC2) as input images. The use of distorted input images requires that the undistort parameter be set to yes to turn on distortion correction of the source objects positions from each image. This removes the necessity to actually run MultiDrizzle, PyDrizzle, or drizzle on each input.

Several other parameters act to directly control the operation of the underlying daofind, xyxymatch and geomap tasks to allow the user to fine-tune the results. The "computesig" parameter, for example, sets Tweakshifts so that it will compute a global value for the sigma based on the sky value found in the reference image. This value will then be used in daofind to tune the object-finding algorithm for the images being combined. Since it is a global value, all input exposures need to be normalized to the same units, either count rates or counts or electrons, in order for the object identification to work the same on all input images.

The source lists from each image generally will include cosmic rays as detected sources, sometimes significantly confusing object identification between images. Long-exposure observations often have more cosmic ray events that source objects, so weeding them out in those cases would improve the efficiency of identifying common sources between images. One such method for trimming potentially bad or confusing sources from each source list would be to set a flux limit and only use sources above that limit. The fluxes reported in the default daofind source object lists are given as magnitude values. Thus, setting a limit based on the daofind magnitudes for the reference image as the fluxmax or fluxmin parameters and setting the ascend parameter to yes would allow the source lists to be trimmed of all sources fainter than the provided limit. This new trimmed source list would then be used in matching sources between images and for the final fitting for the shifts.



The use of cross-correlation to determine shifts can be selected by setting the parameter findmode to "cross-corr". This technique allows shifts to be computed between images which consist primarily of large extended sources with few or no point-sources. The algorithm implemented by Tweakshifts relies on running the tasks crossdriz to perform the cross-correlations between the images and then shiftfind to determine the shifts from the cross-correlation images. As with the catalog finding method, several parameters have been provided to directly control the operation of these underlying IRAF tasks. The inputs for this step, however, MUST be distortion-free images which all share the same WCS. These can be generated as the single-drizzle products from MultiDrizzle by turning off all steps past the single drizzle processing.


Available Options

Table 5.11: Available Parameters for Tweakshifts
Default Value
ASN table, ASCII file single image
input shiftfile with initial values
filename of the OUTPUT reference WCS
output shift file
boolean [True, False]
mode for finding shifts
string [catalog, cross-cor]
generate catalog with this task
string [sextractor, daofind]
parameter listing
sextractor parameters
apply distortion correction to input image positions
boolean [False, True]
automatically compute sigma for all inputs
boolean [False, True]
key for selecting idc table
string [idctab,cubic,trauger,none]
remove intermediate files
boolean [True, False]
print extra messages during processing
boolean [False, True]
Coordinate File Description
file containing coordinate filenames for input files
column number for x positions
column number for y positions
column number for flux/magnitude values
maximum value for valid objects
minimum value for valid objects
units of flux values used for sorting
string [counts,cps,mag]
number of brightest objects to keep after sorting
Object Detection Parameters
minimum number of objects acceptable for matching
maximum number of objects to match
the matching algorithm
string [tolerance,trianges]
the matching tolerance in pixels
the fwhm of the psf in scale units
standard deviation of the background in counts
minimum good data value
maximum good data value
threshold in sigma for feature detection
width of convolution kernel in sigma
fitting geometry
string [shift,xyscale,rotate,rscale,rxyscale,general]
surface type
string [checbyshev,legendre,polynomial]
maximum number of rejection iterations
rejection limit in sigma units
Cross-correlation Parameters
reference image for cross correlation
margin to strip down
edge region to taper
pad working area to prevent wraparound affects?
boolean [False, True]
cross correlation peak fwhm
cross correlation peak ellipticity
cross correlation peak position angle
box size used to fit cross correlation peak

5.5.4 Sky Subtraction Introduction

The default behavior of MultiDrizzle in the pipeline is to estimate the sky by means of the statistical distribution of pixels in the image, and subtract that from copies of the input files that are used to create the final drizzled image. The input files that are delivered to the user from the archive are not modified, ie they are not sky subtracted, but the sky values calculated by MultiDrizzle are contained in the header keyword mdrizsky.

For cameras with multiple detectors (such as ACS/WFC, WFPC2, or WFC3), the sky values in each exposure are first measured separately for the different detectors. Then these different values are compared, and the lowest measured sky value is used as the estimate for all the detectors for that exposure. This is based on the premise that for targets consisting of large extended or bright sources, the pixel intensity distribution in one or more of the detectors may be significantly skewed toward the bright end by the target itself, thereby overestimating the sky on that detector. If the other detector is less affected by such a target, then its sky value will be lower, and can therefore also be substituted as the sky value for the detector with the bright source.

The default behavior for MultiDrizzle in the pipeline is to perform iterative sigma-clipping, starting with the full range of pixel intensity values and calculating the standard deviation, then rejecting pixels with values deviating from the mean by more than 4 sigma either way, repeating this for a total of 5 iterations, and then using the median value of the final distribution as the sky value. This gives good results for a wide range of datasets, but occasionally the sky is still slightly overestimated and can be improved by post-pipeline processing.

We recommend sky subtraction to be turned on for broad-band data or other data with significant amounts of sky. This is because the sky background is retained in the science data header files; the pipeline MultiDrizzle calculates a sky background value and subtracts it on-the-fly before creating the pipeline drizzled product, but does not modify the science pixel values in the original science files. Therefore, turning on this step when running MultiDrizzle offline will ensure that the sky background keywords are read from the header and are applied before drizzling the data.

Of course, it should be kept in mind that if the background of a field is this strongly affected by the bright sources within it, then the purpose of the background measurement is primarily to ensure that there are no more varying offset levels in all the exposures prior to combination, due to sky variations from one exposure to the next; the true background may be difficult or impossible to determine in such cases.

This is also why background subtraction is turned off by default when running MultiDrizzle in the pipeline on narrow-band data and UV observations that are "dark" (i.e., have no geocoronal emission in the filter bandpass), since: 1) the sky background through such filters is much lower than through an optical broad-band filter, and (2) such observations are often of extended diffuse emission-line targets whose flux is much higher than the background, thus any automated attempt to measure the background may introduce errors that are larger than the background itself. However, if the user is able to determine an accurate background level in such cases, then the above mechanism may be used to propagate those values directly to MultiDrizzle.

Sky subtraction is generally recommended for optimal flagging and removal of CR's, when the sky background is more than a few electrons. However, some science applications may require the sky to not be removed. Thus, the final drizzle step can be performed with no sky subtraction. If you turn off sky subtraction, you should also set drizzle.pixfrac =1. Otherwise variations in sky between images will add noise to your data. Methodology

The clipped mode is computed for each input chip and scaled to a reference plate scale as an estimate of the sky background. The lowest scaled value for each chip of an observation (file) is then re-scaled to the chip's plate scale and subtracted as the sky value. The primary header of each input image is updated with this value. In lieu of having MultiDrizzle compute the sky value, the user can supply their own sky value as a keyword in the input file header. This keyword name would then be given to MultiDrizzle in the `skyuser' parameter.

Input: Aside from the input parameters, this step requires opening each input image to access the science (SCI) extensions for computing the sky values.

Output: The input file primary headers get updated with the computed sky value, and each input image's SCI array (or copy, if `workinplace' is set to False), gets sky subtracted. Parameter Details

skysub: Turn on or off sky subtraction on the input data.

skywidth: Bin width, in sigma, used to sample the distribution of pixel flux values in order to compute the sky background statistics.

skystat: Statistical method for determining the sky value from the image pixel values. Valid options are:

skylower: Lower limit of usable pixel values for computing the sky. This value should be specified in units of electrons.

skyupper: Upper limit of usable pixel values for computing the sky. This value should be specified in units of electrons.

skyclip: Number of clipping iterations to use when computing the sky value.

skylsigma: Lower clipping limit, in sigma, used when computing the sky value.

skyusigma: Upper clipping limit, in sigma, used when computing the sky value.

skyuser: Name of header keyword which records the sky value already subtracted from the image by the user. Basic Example

Figure 5.15: Histogram Resulting from an Image Dominated by Scattered Light

Histogram of pixel intensity values in the calibrated science data of j8bt07oyq_flt.fits, which is dominated by scattered light from bright stars, introducing a significant bright skew to the distribution of pixel intensities near the background, thereby strongly affecting any automated estimate of the background. In cases such as this, the choice median, mean or mode can significantly affect the value assigned to the sky, and users may wish to try the various options, or carefully examine the statistics of their image before combining them.

5.5.5 Cosmic Ray Rejection

Few HST observing proposals have sufficient time to take a number of exposures at each of several dither positions. Therefore, if dithering is to be of wide-spread use, one must be able to remove cosmic rays from data where repeated images at the same position on the sky are not taken. We have therefore adapted Drizzle to the removal of cosmic rays, and automated this entire process in the task MultiDrizzle (Section 5.4). As the techniques involved in cosmic ray removal are also valuable in characterizing the image fidelity of drizzle, we will discuss them first. Here then is a short description of the method we use for the removal of cosmic rays:

  1. Drizzle each image onto a separate sub-sampled output image using pixfrac = 1.0 where each image has the same WCS as the final output image.
  2. Take the median of the resulting aligned drizzled images. This provides a first estimate of an image free of cosmic rays.
  3. Map the median image back to the input plane of each of the individual images, taking into account the image shifts and geometric distortion. This can be done by interpolating the values of the median image using the IRAF program named blot.
  4. Take the spatial derivative of each of the blotted output images. This derivative image is used in the next step to estimate the degree to which errors in the computed image shift or the blurring effect of taking the median could have distorted the value of the blotted estimate.
  5. Compare each original image with the corresponding blotted image. Where the difference is larger than can be explained by noise statistics, the flattening effect of taking the median, or an error in the shift, the suspect pixel is masked.
  6. Repeat the previous step on pixels adjacent to pixels already masked, using a more stringent comparison criterion.
  7. Finally, drizzle the input images onto a single output image using the pixel masks created in the previous steps. For this final combination, a smaller pixfrac than in the first step will usually be used in order to maximize the resolution of the final image. Image Alignment Requirement

One of the primary complicating factors in accurately determining what pixels are affected by cosmic rays remains the alignment of the images. Any misalignment of a star from one image to the next by more than 0.1 pixel will generally lead to a misidentification of a cosmic ray for the misaligned image. This level of error can result from pointing errors or inaccurate distortion models, both of which must be addressed prior to accurate identification of cosmic rays. The images generated in the first step of this process, the images drizzled to the final output WCS, can be used to verify the alignment and accuracy of the distortion-correction for each of the input images. Any processing which cross-matches sources from one image to the next can be used to compute the offset (shift) between the images. The task Tweakshifts (Section 5.5.3) demonstrates how these images can be used to determine the shifts using either catalog matching of sources or cross-correlation based on IRAF tasks such as daofind, xyxymatch and geomap.

5.5.6 Selecting the Optimal Scale and Pixfrac

Pixels in the original input images are mapped into pixels in the subsampled output image, taking into account shifts and rotations between images and the optical distortion of the camera. However, in order to avoid re-convolving the image with the large pixel "footprint" of the camera, drizzle allows the user to shrink the pixel before it is averaged into the output image. The new shrunken pixels, or "drops", rain down upon the sub-sampled output. The value of an input pixel is averaged into an output pixel with a weight proportional to the area of overlap between the "drop" and the output pixel. If the drop size is sufficiently small not all output pixels will have data added to them from each input image. One must therefore choose a drop size that is small enough to avoid degrading the image, but large enough so that after all images are drizzled the coverage is reasonably uniform. The drop size is controlled by a user-adjustable parameter called pixfrac, which is simply the ratio of the linear size of the drop to the input pixel (before any adjustment due to the geometric distortion of the camera). Thus interlacing is equivalent to drizzling in the limit as pixfrac 0.0, while shift-and-add is equivalent to pixfrac = 1.0. The degree of subsampling of the output is controlled by the user through the scale parameter, s, which is the ratio of the linear size of an output pixel to an input pixel.

This can be formally represented through the following relations. When a pixel from an input image with data value and user defined weight is added to an output image pixel with the value , weight and fractional pixel overlap , the resulting values and weights of that same pixel, and are:

where a factor is introduced to conserve surface intensity, and there and are used to distinguish the input and output pixel indices. In practice, drizzle applies this iterative procedure to the input data, pixel-by-pixel, image-by-image. Thus, after each input image is processed, there is a usable output image and weight, I and W. The final output images, after all input have been processed, can be written as:

In nearly all cases , since very few input pixels overlap a given output pixel. When the dithered positions of the input images map directly onto the centers of the output grid, and pixfrac and scale are chosen so that p is only slightly greater than s, one obtains the full advantages of interlacing: because the power in an output pixel is almost entirely determined by input pixels centered on that output pixel, the convolutions with both p and G effectively drop away. Nonetheless, the small overlap between adjacent drops fills in missing data.

Noise considerations must also be taken into account and users are strongly encouraged to examine the section on Weight Maps and Correlated Noise (Section 3.3). Choosing the drizzle.scale Parameter

The values used for the drizzle.pixfrac and drizzle.scale parameters depends on the number of input images and the size of the shifts. In general, the full-width at half-maximum of a PSF in the final image should be around 2.5 pixels for a well-sampled output image. In most cases choosing an output pixel size equal to 0.6 or 0.5 the input pixel works well with HST instruments. The output pixel size specified with the drizzle.scale parameter is given in arcseconds. This is because some instruments, such as WFPC2, have more than one pixel size associated with them. The drizzle.pixfrac parameter, however, is given as a fraction of the input pixel size. Choosing the drizzle.pixfrac Parameter

For observers who generally have three to four dither pointings with 0.5 pixel shift increments, the drop size should be larger than the output pixel size. (for example, for WFPC2, scale=0.05 arcsec and pixfrac=0.6) - remember, the scale is now given in arcseconds while the pixfrac is relative to the input pixel. This choice of parameters allows some of the "drop" to spill over to adjacent pixels, thus recovering some resolution to the image. For offsets that are close to integer values, it is difficult to recover any resolution, and a large drop size is recommended (like pixfrac=0.8 or 0.9 with output pixel one-half the input).

One test of the choice of drizzle.scale and drizzle.pixfrac is to measure the r.m.s of the weight image away from the edges of the output image (preferably under a region of sky). The r.m.s. of the final drizzle weight image should be less than 20% to 30% of the median (use imstat to get these numbers). A larger r.m.s. suggests that the weights are so variable that the ignoring the weight image when doing photometry will add noticeably to the final photometric errors.

Since the behavior of the final image statistics will depend significantly upon the exact number of offsets, as well as the degree of sub-pixel sampling, observers are encouraged to experiment by re-running drizzle a few times with different values of pixfrac, to see which yields the best results for the data.

5.5.7 Selecting the `Bits' Parameter Data Quality Flags for Bad Pixels: The `Bits' Parameter

The data quality flags which were set during calibration can be used as bit masks when drizzling, and the user may specify which bit values should actually be considered `good' and included in image combination. This is done via the parameters `driz_sep_bits' and `final_bits'. Any pixel flagged otherwise will be considered `bad' and will be given zero weight. If a pixel is flagged as bad but has non-zero weight in any other image, the pixel will be replaced with the value from that image during image combination. Otherwise, the pixel will be replaced with the `fill value'. The default `fill value' used by MultiDrizzle during pipeline processing is INDEF. If the `fill value' is set to INDEF during reprocessing with MultiDrizzle, and if there are no pixels with non-zero weight, the pixel will retain its original value.

The flags for the DQ array are presented in the Data Handbook for each instrument. For ACS, for example, the bits value used during pipeline processing is 96 which is the sum of 32 and 64 (two flag values for warm pixels). Note that these pixels were flagged during calibration as "bad" or "suspect", but may have been corrected in later processing steps. The bits parameter indicates which "suspect" pixels to keep. For ACS, for example, when the total counts in a given pixel have exceeded the full well, the DQ flag is 256. However, testing shows that counts still accumulate in a highly linear manner and that ACS photometry of saturated stars is quite practical if using a gain that samples the full well depth.

The default bits value in the off-line MultiDrizzle software available for reprocessing is zero, but can be reset by the user prior to reprocessing. The value chosen for this parameter is completely up to the user and should be selected based on the specific calibration needs. For ACS/HRC processing, the user may want to add 8 (masked by aperture) to the value so that pixels which lie behind the occulting finger mask are not given zero weight. (In the ACS Data Handbook, Figure 4.7 of Section 4.5.1, bit 8 is not included in the bits parameter prior to drizzling. Pixels behind the occulting finger were therefore given zero weight and replaced with the default fill value=zero.) Of course, the photometry in this region will be inaccurate due to improper flat fielding, so this bit would likely only be set for aesthetic reasons. Using the 'Bits' Parameter to Update the Masks

Before executing any of the 7 processing steps, MultiDrizzle completes several initialization steps which are described in Section 5.4.8. ■During ■this initialization process, two separate mask files are created for each chip.

The first of these files is called the single mask, and every pixel which has been ■flagged in the input image DQ array ■will be assigned a value of zero in the mask. ■■All other pixels are assumed to be good and will be assigned a value of 1.0. ■The single mask image is used in the 'driz_separate' step to create the '*_single_sci.fits' images which are used as input for creating the median image. The user may tell MultiDrizzle to ignore specific DQ flag values by specifying the parameter 'driz_sep_bits.' ■This can be particularly helpful when ■objects have been erroneously flagged as cosmic rays during pipeline processing. ■For ACS data, for example, the user would set ■'driz_sep_bits=4096', and all pixels originally flagged as cosmic rays will be assumed to be 'good', and these pixels will instead be given a value of 1.0 in the single mask.

The second mask image is called the final mask. ■The final mask is used in the 'driz_combine' step to create the final drizzled product. ■The parameter 'final_bits' allows the user to tell MultiDrizzle which DQ flags to ignore for the final image combination. When this parameter is left to the default value of zero, every pixel flagged in the input image DQ array is assigned a value of zero in the mask.

The two mask files are created during the software initialization and subsequently updated when running MultiDrizzle steps 1 and 6 to compute the static mask and the cosmic ray mask. ■In step 1, the static mask is retained in ■memory and is combined with the previous versions of the single mask and final mask to create updated versions of these images. In step 6, the cosmic ray mask is computed, and when the parameter 'driz_cr_corr=yes', this mask ■will be written to a file called '*_crmask.fits'. This mask should be inspected and blinked with the original input image to verify the quality of the rejection. ■The final mask is then updated with those pixels flagged in the cosmic ray mask. ■In ■addition, the DQ array of the original input image is updated with the pixels flagged as cosmic rays.

During final drizzle combination in step 7, the DQ array of the input image is used in combination with the 'final_bits' parameter to produce a new version of the final mask which will be used to create the drizzled product '*_drz.fits'.

Space Telescope Science Institute
Voice: (410) 338-1082