Aligning HST Mosaics

Note: The notebook in this repository 'Initialization.ipynb' goes over many of the basic concepts such as the setup of the environment/package installation and should be read first if you are new to HST images, DrizzlePac, or Astroquery.

Introduction

This notebook demonstrates how to align and drizzle mosaicked tiles of the Eagle Nebula (M16) obtained with WFC3 with both UVIS and IR detectors. It is based on the example highlighted in the following WFC3 technical report: ISR 2015-09: Combining WFC3 Mosaics of M16 with DrizzlePac and highlights special features in DrizzlePac to improve mosaics.

In prior alignment tutorials, building up an aligned set of tiles required an iterative approach. Now, mosaic alignment can be achieved in a single step by building up an expanded reference catalog ‘on-the-fly’. New sky matching options make it easier to produce seamless mosaics, which can be challenging for extended sources with little or no blank sky.

In [1]:
from astroquery.mast import Observations
from ccdproc import ImageFileCollection
from astropy.table import Table
from astropy.io import fits
from astropy.io import ascii
from astropy.visualization import ZScaleInterval
from IPython.display import Image
import glob
import matplotlib.pyplot as plt
import numpy as np
import os
import shutil
from drizzlepac import tweakreg
from drizzlepac import astrodrizzle
The following task in the stsci.skypac package can be run with TEAL:
                                    skymatch                                    

ModuleNotFoundErrorTraceback (most recent call last)
ModuleNotFoundError: No module named 'numpy.core._multiarray_umath'
 Coordinate transformation and image resampling library NOT found!

 Please check the installation of this package to insure C code was built successfully.

ImportErrorTraceback (most recent call last)
/opt/conda/envs/notebooks_env/lib/python3.6/site-packages/drizzlepac/ablot.py in <module>
     22 try:
---> 23     from . import cdriz
     24 except ImportError:

ImportError: numpy.core.multiarray failed to import

During handling of the above exception, another exception occurred:

ImportErrorTraceback (most recent call last)
<ipython-input-1-a076c3b1ab48> in <module>
     11 import os
     12 import shutil
---> 13 from drizzlepac import tweakreg
     14 from drizzlepac import astrodrizzle

/opt/conda/envs/notebooks_env/lib/python3.6/site-packages/drizzlepac/__init__.py in <module>
     19 from .version import *
     20 
---> 21 from . import ablot
     22 from . import adrizzle
     23 from . import astrodrizzle

/opt/conda/envs/notebooks_env/lib/python3.6/site-packages/drizzlepac/ablot.py in <module>
     26     print('\n Coordinate transformation and image resampling library NOT found!')
     27     print('\n Please check the installation of this package to insure C code was built successfully.')
---> 28     raise ImportError
     29 
     30 from .version import *

ImportError: 

1. Observations

Mosaics of the Eagle Nebula were acquired by HST GO/DD program 13926 in September 2014 for HST's 25th Annivery. A 2x2 tile mosaic with the IR detector (~4 arcmin across) was observed in the F110W and F160W filters.

A slightly larger 2x2 mosaic with the UVIS detector (~5 arcmin across) was observed with the F502N, F657N, and F673N filters. Small dithers between exposures in a given tile will fill in the UVIS chip gap and allow for the rejection of cosmic rays and detector artifacts. More detail on the observing strategy may be found in the Phase II file.

Two additional UVIS tiles overlap the central portion of 2x2 mosaic in order to have very high signal-to-noise in the Eagle's pillars. These two visits (09,10) were not included in this example for brevity. The data used in this notebook example is also limited to the IR/F160W filter (visits 01-04) and UVIS/F657N filters (visits 05-08), shown in the diagrams below.

 IR Mosaic           UVIS Mosaic         UVIS(overlap)
 ____ ____            ____ ____            ____ 
|    |    |          |    |    |          |    |   
| 02 | 01 |          | 06 | 05 |          | 09 | 
|____|____|          |____|____|          |____|
|    |    |          |    |    |          |    |  
| 04 | 03 |          | 08 | 07 |          | 10 |
|____|____|          |____|____|          |____|

2. MAST Download

In the steps below, the calibrated IR data (*_flt.fits) and the calibrated, CTE-corrected UVIS data (*_flc.fits) are retrieved from MAST and placed in the same directory as this notebook, along with the associated telemetry and engineering files (*_spt.fits).

In [ ]:
# Retrieve the IR/F160W calibrated FLT and SPT data products  
science_list = Observations.query_criteria(proposal_id='13926', filters='F160W')
Observations.download_products(science_list['obsid'], mrp_only=False, download_dir='./science',
                               productSubGroupDescription=['FLT', 'SPT'])

science_files = glob.glob(os.path.join(os.curdir, 'science', 'mastDownload', 'HST', '*', '*fits'))
for im in science_files:
    root = im.split('/')[-1]
    os.rename(im, './' + root)
shutil.rmtree('science/')
In [ ]:
# Obtain the UVIS/F657N calibrated FLC and SPT data products from visits 05-08 
science_list = Observations.query_criteria(proposal_id='13926', filters='F657N', obs_id='ICK90[5678]*')
Observations.download_products(science_list['obsid'], mrp_only=False, download_dir='./science',
                               productSubGroupDescription=['FLC', 'SPT'])

science_files = glob.glob(os.path.join(os.curdir, 'science', 'mastDownload', 'HST', '*', '*fits'))
for im in science_files:
    root = im.split('/')[-1]
    os.rename(im, './' + root)
shutil.rmtree('science/')

3. Dithers

3a. IR detector

IR exposures were obtained in Visits 01-04. (The visit ID is found in the 5th and 6th character of the filename). Each visit (mosaic tile) consists of a pair of exposures using the WFC3-IR-DITHER-BLOB dither of 7.2" along the y-axis (pattern_orient=90 degrees). This dither can be seen when comparing the POSTARG2 keyword between pairs of exposures in a given visit in the table below.

Pairs of IR exposures making up each visit are referred to as v01a and v01b in this notebook. The first four images listed in the table below are associated with v01a and the last four with v01b.

In [ ]:
collect_ir = ImageFileCollection('./', glob_include="*flt.fits", ext=0,
                                 keywords=["asn_id", "detector", "filter", "nsamp",
                                           "exptime", "postarg1", "postarg2"])

ir_table = collect_ir.summary
ir_table['exptime'].format = '7.1f'
ir_table['postarg1'].format = '7.2f'
ir_table['postarg2'].format = '7.2f'
ir_table

3b. UVIS detector

UVIS exposures were acquired in Visits 05-08. Each UVIS visit (tile) consists of a set of 3 dithered exposures using the WFC3-UVIS-MOSAIC-LINE pattern, with an offset ~12" along a 65 degree diagonal. This dither can be seen in the POSTARG1, POSTARG2 offsets which are ~5" in X and ~10" in Y between exposures in a given visit.

Sets of three exposures making up each UVIS visit are referred to as v05a, v05b, v05c in this notebook.

In [ ]:
collect_uvis = ImageFileCollection('./', glob_include="*flc.fits", ext=0,
                                   keywords=["asn_id", "detector", "filter", "exptime", "postarg1", "postarg2"])
    
uvis_table = collect_uvis.summary
uvis_table['exptime'].format = '7.1f'
uvis_table['postarg1'].format = '7.2f'
uvis_table['postarg2'].format = '7.2f'
uvis_table

4. TweakReg for Mosaics

Before combining observations with AstroDrizzle, the WCS keywords in the header of each input frame should be aligned to sub-pixel accuracy. This may be achieved with TweakReg, which allows users to align sets of images to one another or to an external astrometric reference frame. TweakReg has been enhanced to support the alignment of observations that cover a large area on the sky. Making use of the expand_refcat parameter, TweakReg will build up an expanded reference catalog on the sky to be used for alignment. When set to 'True', TweakReg selects two images from the input list with the largest overlap on the sky, generates source catalogs for each image, and computes a fit (shift, rotation, and/or scale change) from the matched source list.

Next, the algorithm computes the area of overlap of each of these two images with the rest of the input images, and the one with the largest total overlap on the sky is selected as the reference image. Sources from the second image that have not been matched to the reference image catalog are considered good new sources and are added to the reference catalog. In this way, the reference catalog keeps expanding with each new matched image. With a large (expanded) reference catalog it is therefore possible to align images that had no direct overlap with the starting image.

4a. IR Alignment

For this large multi-filter dataset, the user should carefully consider which observations to align and combine first. These will serve as a reference image for aligning additional filters. The broadband IR images of M16 contain a large number of stars distributed uniformly over the field of view. The UVIS frames, on the other hand, are largely devoid of point sources and full of cosmic-rays which can trip up TweakReg when trying to compute a fit. Even though the IR detector has a smaller footprint on the sky and the IR PSF is more undersampled, the high density of stars makes it a better anchor for aligning the UVIS tiles.

Of the two IR filters, F160W has the largest number of point sources and therefore makes a good choice for the reference image. To generate source lists for matching, the TweakReg parameter conv_width should be set to approximately twice the FWHM of the PSF, ~2.5 pixels for IR observations and ~3.5 pixels for UVIS observations. TweakReg will automatically compute the standard deviation of the sky background (skysigma), so the number of sources in each catalog may be controlled simply by changing the ‘threshold’ parameter.

In this example, TweakReg is run in ‘non-interactive’ mode (interactive='False') so that the astrometric fit residuals and vectors diagrams will be saved as png files in the user’s local directory for inspection. Once the parameters have been fine-tuned and the fit looks adequate, users may run TweakReg a second time (see below) to update the image header WCS keywords by setting the parameter updatehdr to True.

In [ ]:
tweakreg.TweakReg('*flt.fits',
                  imagefindcfg={'threshold': 50, 'conv_width': 2.5},
                  expand_refcat=True,
                  enforce_user_order=False,
                  shiftfile=True,
                  outshifts='shift160_flt.txt',
                  searchrad=2.0,
                  ylimit=0.3,
                  updatehdr=True,
                  reusename=True,
                  wcsname='IR_FLT',
                  interactive=False)

4b. Inspect the shift file to verify the pointing residuals

With the threshold parameter set to 50 sigma, TweakReg finds ~500 objects per FLC image, matches the individual catalogs, and computes residual shifts between exposures. These offsets (given in pixels at the native IR scale=0.1283”/pixel) are recorded in an output “shift file” which is shown below. Note that the fourth row of the shift file corresponds to the reference image ‘ick902neq_flt.fits’ which was automatically selected by TweakReg. The computed offsets reflect updates to the header WCS required to correct for small pointing errors.

In [ ]:
shift_table = Table.read('shift160_flt.txt',
                         format='ascii.no_header', 
                         names=['file', 'dx', 'dy', 'rot', 'scale', 'xrms', 'yrms'])

formats = ['.2f', '.2f', '.3f', '.5f', '.2f', '.2f']
for i, col in enumerate(shift_table.colnames[1:]):
    shift_table[col].format = formats[i]
shift_table

The expected pointing accuracy for various observing scenerios is summarized in the DrizzlePac Handbook Appendix B. Exposures making up visit-level drizzled products are typically aligned to 2-5 milliarcsecond (mas) accuracy with fine-lock on 2 guide stars. The shift file above confirms this, with offsets (dx,dy) between pairs of exposures in the same visit visit averaging ~0.05 IR pixels (6 mas). For different visits using the same set of guide stars, offsets of ~50-100 mas (0.6-1.2 IR pixels are expected. For visits with different sets of guide stars, the pointing accuracy is typically 0.2-0.5 arcseconds.

The cell below shows how to check the image header for the 'Dominant' and 'Secondary' guide stars used. Since the M16 tiles (visits) each used different guide star pairs, the relatively large offsets (>1 arcsec) required for visits 03 and 04 to match the visit 02 reference image are not surprising, though are larger than expected. Note that the TweakReg searchrad parameter was increased to 2.0" from the default value of 1.0" to allow TweakReg to find the correct fit for these visits.

In [ ]:
collect_spt = ImageFileCollection('./', glob_include="ick9*0_spt.fits", ext=0,
                                  keywords=["asn_id", "config", "dgestar", "sgestar"])
table_spt = collect_spt.summary
table_spt

With enforce_user_order='True', the FLC files were aligned in the order shown below. The exposure v02b was automatically selected as the reference and then v02a was selected as the image with the most overlap on the sky. The two input catalogs had 601 and 574 sources, respectively, and this gave 483 matches for the two exposures, after sigma-clipping. The reference catalog was then expanded by adding 73 new objects from the dithered exposure v02a to the matched catalog for a total of 674 sources.

TweakReg next moved over to tile 04 and to align exposure v04a to the expanded catalog. In the overlap region between tiles, 63 matches were found along the upper edge of v04a. The reference catalog was expanded once again, adding the 411 sources from v04a to give a total of 1087 sources. Now when v04b is aligned, the number of matches is 388, since the expanded catalog now includes v04a. TweakReg continues in this way until all input frames are aligned to the expanded catalog, and the number of matches for each exposure is listed in the ascii table below.

The table below lists visit/exposure ID, the filename, and the number of matched sources for the 8 IR exposures. ID FILE N_Match v02b ick902neq_flt - v02a ick902n9q_flt 483 v04a ick904obq_flt 63 edge v04b ick904ogq_flt 388 v01a ick901hzq_flt 44 edge v01b ick901i7q_flt 454 v03a ick903n4q_flt 57 edge v03b ick903ncq_flt 408

4c. Inspect the IR fits

To verify that TweakReg obtained a good fit between matched source catalogs, it is useful to inspect the results before updating the image header WCS. Below sources matched with the reference frame (v02b) are overplotted on the first input image (v02a). It can be useful to check that TweakReg locked onto stars and not hot pixels or other detector artifacts before proceeding. Next, the vector residuals plot is displayed and checked for any systematics. Finally, the 4-panel plot of fit residuals: dx, dy vs X and Y is inspected to verify that the residuals cluster around zero.

In [ ]:
# v02a matches
plt.figure(figsize = (20, 7))
data = fits.open('ick902n9q_flt.fits')['SCI', 1].data
zscale = ZScaleInterval()
z1, z2 = zscale.get_limits(data)
plt.imshow(data, cmap='Greys', origin='lower', vmin=z1, vmax=z2)
match_tab = ascii.read('ick902n9q_flt_catalog_fit.match')  # load match file in astropy table
x_coord, y_coord = match_tab['col11'], match_tab['col12']
plt.scatter(x_coord, y_coord, s=30, edgecolor='r', facecolor='None')
plt.ylim(0, 1014)
plt.xlim(0, 1014)
plt.title('Match: v02a to v02b(Ref)', fontsize=20)
In [ ]:
# v02a vector residuals
Image(filename='vector_ick902n9q_flt.png', width=500, height=300)
In [ ]:
# v02a fit residuals
Image(filename='residuals_ick902n9q_flt.png', width=500, height=300)
In [ ]:
# v04a matches (Edge overlap region only)
plt.figure(figsize = (20, 7))
data = fits.open('ick904obq_flt.fits')['SCI', 1].data
z1, z2 = zscale.get_limits(data)
plt.imshow(data, cmap='Greys', origin='lower', vmin=z1, vmax=z2)
match_tab = ascii.read('ick904obq_flt_catalog_fit.match')  # load match file in astropy table
x_coord, y_coord = match_tab['col11'], match_tab['col12']
plt.scatter(x_coord, y_coord, s=30, edgecolor='r', facecolor='None')
plt.ylim(0, 1014)
plt.xlim(0, 1014)
plt.title('Match: v04a to v02b(Ref)', fontsize=20)
In [ ]:
# v04a vector residuals (Edge gives narrow range of Y-values)
Image(filename='vector_ick904obq_flt.png', width=500, height=300)
In [ ]:
# v04a fit residuals (Edge gives narrow range of Y-values)
Image(filename='residuals_ick904obq_flt.png', width=500, height=300)
In [ ]:
# v04b matches (Full image since v04a already aligned)
plt.figure(figsize = (20, 7))
data = fits.open('ick904ogq_flt.fits')['SCI', 1].data
z1, z2 = zscale.get_limits(data)
plt.imshow(data, cmap='Greys', origin='lower', vmin=z1, vmax=z2)
match_tab = ascii.read('ick904ogq_flt_catalog_fit.match')  # load match file in astropy table
x_coord, y_coord = match_tab['col11'], match_tab['col12']
plt.scatter(x_coord, y_coord, s=30, edgecolor='r', facecolor='None')
plt.ylim(0, 1014)
plt.xlim(0, 1014)
plt.title('Match: v04b to v02b(Ref)', fontsize=20)
In [ ]:
# v04b vector residuals
Image(filename='vector_ick904ogq_flt.png', width=500, height=300)
In [ ]:
# v04b fit residuals
Image(filename='residuals_ick904ogq_flt.png', width=500, height=300)
In [ ]:
# v01a matches (Edge overlap region only)
plt.figure(figsize = (20, 7))
data = fits.open('ick901hzq_flt.fits')['SCI',1].data
z1, z2 = zscale.get_limits(data)
plt.imshow(data, cmap='Greys', origin = 'lower', vmin = z1, vmax = z2)
match_tab = ascii.read('ick901hzq_flt_catalog_fit.match')                 #load match file in astropy table
x_coord, y_coord = match_tab['col11'], match_tab['col12']
plt.scatter(x_coord, y_coord, s = 30, edgecolor = 'r', facecolor = 'None')
plt.ylim(0,1014)
plt.xlim(0,1014)
plt.title('Match: v01a to v02b(Ref)', fontsize=20)
In [ ]:
# v01a vector residuals (Left edge gives narrow range of X-values)
Image(filename='vector_ick901hzq_flt.png', width=500, height=300)
In [ ]:
# v01a fit residuals (Left edge gives narrow range of X-values)
Image(filename='residuals_ick901hzq_flt.png', width=500, height=300)
In [ ]:
# v01b matches (Full image since v01a already aligned)

plt.figure(figsize = (20,7))
data = fits.open('ick901i7q_flt.fits')['SCI',1].data
z1, z2 = zscale.get_limits(data)
plt.imshow(data, cmap='Greys', origin = 'lower', vmin = z1, vmax = z2)
match_tab = ascii.read('ick901i7q_flt_catalog_fit.match')                 #load match file in astropy table
x_coord, y_coord = match_tab['col11'], match_tab['col12']
plt.scatter(x_coord, y_coord, s = 30, edgecolor = 'r', facecolor = 'None')
plt.ylim(0,1014)
plt.xlim(0,1014)
plt.title('Match: v01b to v02b(Ref)', fontsize=20)
In [ ]:
# v01b vector residuals
Image(filename='vector_ick901i7q_flt.png', width=500, height=300)
In [ ]:
# v01b fit residuals
Image(filename='residuals_ick901i7q_flt.png', width=500, height=300)
In [ ]:
# v03a matches (Top and left edge overlap regions only)
plt.figure(figsize = (20, 7))
data = fits.open('ick903n4q_flt.fits')['SCI', 1].data
z1, z2 = zscale.get_limits(data)
plt.imshow(data, cmap='Greys', origin='lower', vmin=z1, vmax=z2)
match_tab = ascii.read('ick903n4q_flt_catalog_fit.match')  # load match file in astropy table
x_coord, y_coord = match_tab['col11'], match_tab['col12']
plt.scatter(x_coord, y_coord, s=30, edgecolor='r', facecolor='None')
plt.ylim(0, 1014)
plt.xlim(0, 1014)
plt.title('Match: v03a to v02b(Ref)', fontsize=20)
In [ ]:
# v03a vector residuals (Top and left edges)
Image(filename='vector_ick903n4q_flt.png', width=500, height=300)
In [ ]:
# v03a fit residuals
Image(filename='residuals_ick903n4q_flt.png', width=500, height=300)
In [ ]:
# v03b matches (Full image since v03a already aligned)
plt.figure(figsize = (20, 7))
data = fits.open('ick903ncq_flt.fits')['SCI', 1].data
z1, z2 = zscale.get_limits(data)
plt.imshow(data, cmap='Greys', origin='lower', vmin=z1, vmax=z2)
match_tab = ascii.read('ick903ncq_flt_catalog_fit.match')                 #load match file in astropy table
x_coord, y_coord = match_tab['col11'], match_tab['col12']
plt.scatter(x_coord, y_coord, s=30, edgecolor='r', facecolor='None')
plt.ylim(0, 1014)
plt.xlim(0, 1014)
plt.title('Match: v03b to v02b(Ref)', fontsize=20)
In [ ]:
# v03b vector residuals
Image(filename='vector_ick903ncq_flt.png', width=500, height=300)
In [ ]:
# v03b fit residuals
Image(filename='residuals_ick903ncq_flt.png', width=500, height=300)

The dx,dy residuals for all IR exposures are all clustered around dx,dy=0 and have an RMS less than 0.1 pixels, indicating a good fit. If the alignment needs to be fine-tuned, changing the parameters threshold, sigma, and searchrad may help TweakReg to lock onto an accurate solution.

5. Mosaicking Features in AstroDrizzle

AstroDrizzle now makes it easier for users to match the sky background when tiling together large mosaics. In prior versions of the software, the sky background was based on clipped statistics in each image separately. The sky background was measured for each chip and the lowest sky value (in electrons/arcsec^2) among all of the chips was adopted. For observations of sparse fields, this approach generally works well. However, when large extended objects fill the detector, there is no true 'blank sky' and the background value will be an overestimate. Additionally, when extended targets are observed as mosaics (e.g. with large dithers), the 'scene' can change significantly between exposures and bias the background estimate.

An error in determining the sky background may in turn impact the cosmic ray rejection, and if severe enough, the resulting photometry. Additionally, by not properly matching the sky background before combining frames, correlated noise will be added to the final drizzled products when differences in the background levels are significant. Until now, the recommended workaround has been for users to give AstroDrizzle an ASCII file (skyfile) containing user-defined background values.

AstroDrizzle now features several new options for computing the sky. One of these, skymethod='match', is useful for “equalizing” the sky background across large mosaics. This method computes differences in sky values using only pixels in common between images. The sky values will then be set relative to the value computed for the input frame with the lowest sky value for which the MDRIZSKY keyword will be set to 0. In this way, the sky background is not removed, but instead equalized before the data are combined. For more details on the sky matching functions used by AstroDrizzle, see the following webpage.

5a. Drizzle the IR/F160W Mosaic

Now AstroDrizzle can be used to combine the full set of F160W frames. In this case the final orientation has been set to -35 degrees so that the pillars will be oriented vertically. Note that users must first set the parameter final_wcs='True' in order to turn on parameters in AstroDrizzle’s step 7a: Custom WCS for Final Output. For these observations, the IR scale (0.08”/pixel) is chosen to be exactly twice that for the UVIS mosaics (0.04”/pixel) by setting final_scale=0.08, and the drizzled images have been oversized slightly to match the sky area on the sky covered by the UVIS. The sky background may be equalized across mosaic tiles by setting the parameter skymethod='match'.

The parameter final_bits defines which DQ flags in the FLT image to treat as good. All other pixels with non-zero DQ values will be assumed to be bad and rejected from the final mosaic. For IR data, these two parameters are typically set to 64+512 in the pipeline, corresponding to warm pixels and IR blobs. This program included a blob dither, however, and so the 512 flag may be removed from the list of good DQ values such that these pixels will be replaced with non-flagged pixels from the accompanying dithered pair.

For IR data, cosmic-rays have already rejected via the 'up-the-ramp' fitting by calwf3, so Steps 3, 4, 5, 6 in AstroDrizzle have been turned off when combining the FLT exposures. See this reference for more information about the various input parameters to AstroDrizzle.

In [ ]:
astrodrizzle.AstroDrizzle('*flt.fits',
                          output='f160w',
                          preserve=False,
                          clean=True,
                          build=False,
                          context=False,
                          skymethod='match', 
                          driz_separate=False,
                          median=False,
                          blot=False,
                          driz_cr=False,
                          final_bits='64',
                          final_wcs=True,
                          final_scale=0.08,
                          final_rot=-35,
                          final_ra=274.721587,
                          final_dec=-13.841549,
                          final_outnx=4000,
                          final_outny=4500)

5b. Display the combined DRZ science and weight images

In [ ]:
sci = fits.getdata('f160w_drz_sci.fits')
fig = plt.figure(figsize=(14, 14))
plt.imshow(sci, vmin=1, vmax=6, cmap='Greys_r', origin='lower')
In [ ]:
wht = fits.getdata('f160w_drz_wht.fits')
fig = plt.figure(figsize=(14, 14))
plt.imshow(wht, vmin=0, vmax=1700, cmap='Greys_r', origin='lower')

6. Align the UVIS FLC frames to the IR mosaic

The F160W drizzled mosaic defines the reference frame for aligning the UVIS filters. In ISR 2015-09, the UVIS visit-level DRC frames were aligned directly to the IR reference image, and this approach was chosen because the long UVIS exposures contain numerous cosmic-rays and relatively few point sources. The Drizzlepac task TweakBack was then used to propagate the updated WCS from the drizzled image header back to the individual FLC input frames making up each association prior to drizzling.

In this notebook, the FLC frames may be aligned directly to the IR mosaic by making use of a new parameter in TweakReg which allows for specific flags in the DQ array of the FLC frames to be used or ignored. The imagefindpars parameter dqbits may be prepended with ‘~’ to the string value to indicate which DQ flags to consider as "bad" pixels. For example, when deriving source catalogs, Tweakreg will ignore any pixels flagged as cosmic-ray flags in the MAST visit-level drizzled data products when dqbits is set to ~4096. This dramatically cuts down the number of false detections due to cosmic rays in the input FLC science arrays. More details on imagefindpars options may be found on the following webpage.

In this example, the threshold value was manually adjusted to get ~50 matches per UVIS exposure. Note that setting the threshold to a very low value does not necessarily translate to a better solution, since all sources are weighted equally when computing fits to match catalogs. This is especially relevant for UVIS data where CTE tails can shift the centroid position slightly along the readout direction for faint sources and potentially bias the fit.

In [ ]:
tweakreg.TweakReg('*_flc.fits',
                  enforce_user_order=False,
                  imagefindcfg={'threshold': 200, 'conv_width': 3.5, 'dqbits': ~4096},
                  refimage='f160w_drz_sci.fits', 
                  refimagefindcfg={'threshold': 50, 'conv_width': 2.5},
                  shiftfile=True,
                  outshifts='shift657_flc.txt',
                  searchrad=5.0,
                  ylimit=0.6,
                  updatehdr=True,
                  wcsname='UVIS_FLC',
                  reusename=True,
                  interactive=False)

6a. Inspect the shift file to verify the pointing residuals

With a threshold of 200 sigma, TweakReg generates catalogs with several hundred objects per UVIS image and it matches 20-70 objects to the IR reference catalog. The computed offsets written to the shift file are reported at the scale of the reference image (0.08”/pixel for the drizzled IR mosaic). Thus a fit rms ~0.20 pixels (below) at the IR scale is equivalent to an rms ~0.10 pixels at the native UVIS scale (0.04"/pixel).

Because each visit was acquired using a unique pair of Guide Stars, the three FLC exposures making up each visit should have roughly similar residual corrections to the WCS, which can be seen in the table below. Note that residual offsets for visit 05 are much larger than those found for the other three visits at nearly 3".

In [ ]:
shift_table=Table.read('shift657_flc.txt', format='ascii.no_header',
                       names=['file', 'dx', 'dy', 'rot', 'scale', 'xrms', 'yrms'])

formats = ['.2f', '.2f', '.3f', '.5f', '.2f', '.2f']
for i, col in enumerate(shift_table.colnames[1:]):
    shift_table[col].format = formats[i]
shift_table

6b. Inspect the UVIS fits

In [ ]:
#v06a matched sources (chip 1)
plt.figure(figsize=(20, 10))
chip1_data = fits.open('ick906l5q_flc.fits')['SCI', 2].data
z1, z2 = zscale.get_limits(chip1_data)
plt.imshow(chip1_data, cmap='Greys', origin='lower', vmin=z1, vmax=z2)
match_tab = ascii.read('ick906l5q_flc_catalog_fit.match')  # load match file in astropy table
match_tab_chip1 = match_tab[match_tab['col15'] == 2]  # filter table for sources on chip 1 (on ext 4)
x_cord, y_cord = match_tab_chip1['col11'], match_tab_chip1['col12']
plt.scatter(x_cord, y_cord, s=50, edgecolor='r', facecolor='None', label='Matched Sources, Chip 1')
plt.ylim(0, 2051)
plt.xlim(0, 4096)
plt.legend(loc='best', fontsize=20)
In [ ]:
#v06a matched sources (chip 2)
plt.figure(figsize=(20, 10))
chip1_data = fits.open('ick906l5q_flc.fits')['SCI', 1].data
z1, z2 = zscale.get_limits(chip1_data)
plt.imshow(chip1_data, cmap='Greys', origin='lower', vmin=z1, vmax=z2)
match_tab = ascii.read('ick906l5q_flc_catalog_fit.match')  # load match file in astropy table
match_tab_chip2 = match_tab[match_tab['col15'] == 1]  # filter table for sources on chip 2 (on ext 1)
x_cord, y_cord = match_tab_chip2['col11'], match_tab_chip2['col12']
plt.scatter(x_cord, y_cord, s=50, edgecolor='r', facecolor='None', label='Matched Sources, Chip 2')
plt.ylim(0, 2051)
plt.xlim(0, 4096)
plt.legend(loc='best', fontsize=20)
In [ ]:
# v06a vector residuals
Image(filename='vector_ick906l5q_flc.png', width=600, height=300)
In [ ]:
# v06a fit residuals
Image(filename='residuals_ick906l5q_flc.png', width=600, height=300)

7. Drizzling the F657N mosaic

The full set of 12 FLC frames have now been aligned and may be drizzled to a mosaic exactly half the scale (0.04"/pixel) of the original IR mosaic (0.08"/pixel). The same output WCS paramters: final_rot,final_ra, and final_dec values are used, but the final_outnx and final_outny are now twice the size of the IR mosaic at 8000x9000 pixels.

To save time for this notebook, the visit-level cosmic ray flags are assumed to be adequate, Steps 3, 4, 5, 6 are turned off, and the parameter resetbits is set to '0' to avoid wiping out the 4096 flags in the FLC data quality arrays. Also, the sky background levels have been provided via the skyfile parameter, but users would typically set skymethod='match' to compute the sky background levels.

To further improve cosmic-ray rejection in the chip gap, AstroDrizzle may alternatively be run with all steps turned on as shown in the text below and with resetbits=4096 to update the DQ flags.

astrodrizzle.AstroDrizzle('*flc.fits',
                          output='f657n_improved',
                          preserve=False,
                          clean=False,
                          build=False,
                          context=False,
                          resetbits=4096,
                          skymethod='match',
                          combine_type='minmed',
                          final_bits='64,16',
                          final_wcs=True,
                          final_scale=0.04,
                          final_rot=-35,
                          final_ra=274.721587,
                          final_dec=-13.841549,
                          final_outnx=8000,
                          final_outny=9000)
In [ ]:
astrodrizzle.AstroDrizzle('*flc.fits',
                          output='f657n',
                          preserve=False,
                          clean=False,
                          build=False,
                          context=False,
                          resetbits=0,
                          skyfile='skyfile.txt',
                          driz_separate=False,
                          median=False,
                          blot=False,
                          driz_cr=False,
                          final_bits='64,16',
                          final_wcs=True,
                          final_scale=0.04,
                          final_rot=-35,
                          final_ra=274.721587,
                          final_dec=-13.841549,
                          final_outnx=8000,
                          final_outny=9000)

7a. Display the combined DRC science and weight images

In [ ]:
sci = fits.getdata('f657n_drc_sci.fits')
fig = plt.figure(figsize=(14, 14))
plt.imshow(sci, vmin=0, vmax=1, cmap='Greys_r', origin='lower')
In [ ]:
sci = fits.getdata('f657n_drc_wht.fits')
fig = plt.figure(figsize=(14, 14))
plt.imshow(sci, vmin=0, vmax=10000, cmap='Greys_r', origin='lower')

Summary

This notebook provides the methodology for creating UVIS and IR mosaics of M16, as well as recommendations for key parameters. As such, it is relevant for users combining multi-visit observations from any HST imaging program, whether mosaics or single pointings. In this example, the IR F160W mosaic was chosen to define the reference frame, since it contained the largest number of point sources, despite the pixels being undersampled. Typically, the detector with the best resolution and PSF sampling would be selected as the reference, but in this case the UVIS catalogs were too sparse.

When building up an expanded reference catalog on the sky, TweakReg was allowed to select the order in which tiles were aligned, since it made little impact on the results. For aligning mosaics with more tiles, however, users are recommended to experiment with changing the order in which images are aligned and seeing how this changes the astrometric residuals. Generally, users are advised to start from the center of the mosaic and work their way out to avoid propogation of errors across the mosaic.

Alternatively, mosaic tiles can be aligned directly to an absolute reference catalog such as GAIA, when a sufficient number of stars are available in the input frames for alignment. For more detail on this methodology, see the notebook 'aligning_to_catalogs.ipynb'.

About this Notebook

Author: J. Mack, STScI WFC3 Team  
Updated: December 17, 2018