diff --git a/_config.yml b/_config.yml index 99e39594f..79a5b28bc 100644 --- a/_config.yml +++ b/_config.yml @@ -48,5 +48,4 @@ html: use_issues_button: true use_repository_button: true # Exclude notebooks that have as-yet unresolved build failures (see tickets SPB-1153SPB-1160, SPB-1168) -exclude_patterns: [notebooks/DrizzlePac/drizzle_wfpc2/drizzle_wfpc2.ipynb, - notebooks/WFC3/dash/dash.ipynb] \ No newline at end of file +exclude_patterns: [notebooks/DrizzlePac/drizzle_wfpc2/drizzle_wfpc2.ipynb] \ No newline at end of file diff --git a/notebooks/WFC3/dash/README.md b/notebooks/WFC3/dash/README.md deleted file mode 100644 index 2b9ef992f..000000000 --- a/notebooks/WFC3/dash/README.md +++ /dev/null @@ -1,23 +0,0 @@ -This notebook is the first in a new Drift And SHift (DASH) pipeline workflow developed to ease the process of reducing DASH data. The pipeline is customizable, able to be changed according to scientific goals of the user, and this first tutorial walks the user from data download to a finished product ready for science analysis. - -Installation ------------- -`wfc3_dash` requires `lacosmicx` to be installed in the environment from [WFC3 Library's](https://github.com/spacetelescope/WFC3Library) installation instructions. `lacosmicx` implements the Laplacian cosmic ray detection algorithm and is needed to reduce dash data. After activating the virtual environment, follow the steps below to install the package: - -1. Install `cython` to allow the installation of `lacosmicx`: -``` ->>> conda install cython -``` -2. Go to https://github.com/cmccully/lacosmicx and clone the repository in the current directory: -``` ->>> git clone https://github.com/cmccully/lacosmicx.git -``` -3. Change directory to the `lacosmicx` repository: -``` ->>> cd lacosmicx -``` -4. Install the `lacosmicx` package: -``` ->>> python setup.py develop -``` -5. **The package is installed; now you will be able to complete the WFC3 Dash notebook.** diff --git a/notebooks/WFC3/dash/dash.ipynb b/notebooks/WFC3/dash/dash.ipynb deleted file mode 100644 index 56f4bb411..000000000 --- a/notebooks/WFC3/dash/dash.ipynb +++ /dev/null @@ -1,817 +0,0 @@ -{ - "cells": [ - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "\n", - "# How to use `wfc3_dash` on DASH data\n", - "\n", - "***\n", - "## Learning Goals:\n", - "By the end of this tutorial, you will:\n", - "- Create difference files, association tables, and segmentation maps using `wfc3_dash`.\n", - "- Subtract background and fix cosmic rays from newly generated FLTs.\n", - "- Align reads to each other for a final product.\n", - "\n", - "## Table of Contents:\n", - "[Introduction](#introduction)
\n", - "[1. Imports](#imports)
\n", - "[2. Download relevant data](#downloads)
\n", - "[3. Run DASH](#DASH)
\n", - "- [3.1 Create DashData object](#object)
\n", - "- [3.2 Create difference files from reads](#diff_files)
\n", - "- [3.3 Create an association table](#asn_table)
\n", - "- [3.4 Create a segmentation map](#seg_map)
\n", - "- [3.5 Subtract background from the difference files](#subtract_ext)
\n", - "- [3.6 Fix cosmic rays](#cosmic_rays)
\n", - "- [3.7 Align reads to each other](#align_each_other)
\n", - "\n", - "[4. Plot original IMA vs. DASH pipeline science result](#plot)
\n", - "[5. Conclusions](#conclusions)
\n", - "[Additional Resources](#resources)
\n", - "[About the Notebook](#about)
\n", - "[Citations](#cite)
" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - " \n", - "## Introduction\n", - "\n", - "This notebook is the first in a new Drift And SHift (DASH) pipeline workflow developed to ease the process of reducing DASH data. The pipeline is customizable, able to be changed according to scientific goals of the user, and this first tutorial walks the user from data download to a finished product ready for science analysis. For more information, see [Momcheva et. al 2016](https://arxiv.org/pdf/1603.00465.pdf) and [WFC3 ISR 2021-01](https://www.stsci.edu/files/live/sites/www/files/home/hst/instrumentation/wfc3/documentation/instrument-science-reports-isrs/_documents/2021/2021-02.pdf)." - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "## 1. Imports\n", - "This notebook assumes you have created the virtual environment in [WFC3 Library's](https://github.com/spacetelescope/WFC3Library) installation instructions.\n", - "\n", - "We import:\n", - "- *os* for setting environment variables\n", - "- *glob* for querying through directories\n", - "- *numpy* for handling array functions\n", - "- *matplotlib.pyplot* for plotting data\n", - "- *astropy* for astronomy related functions\n", - "- *astroquery.Observations* for downloading data\n", - "- *drizzlepac.astrodrizzle* for combining images\n", - "- *reduce_dash* for reducing DASH data" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "import os\n", - "from glob import glob\n", - "\n", - "import numpy as np\n", - "import matplotlib.pyplot as plt \n", - "\n", - "from astropy.io import fits \n", - "from astropy.table import Table\n", - "from astropy.io import ascii\n", - "from astroquery.mast import Observations\n", - "\n", - "from reduce_dash import DashData\n", - "%matplotlib inline" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "## 2. Download relevant data\n", - "\n", - "Retrieve the table of observations associated with 15238." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "obsTable = Observations.query_criteria(proposal_id=['15238'], obs_id=['IDNM0J030'])" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Get the full list of products associated to the table and restrict the list to IMA files." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "product_list = Observations.get_product_list(obsTable)\n", - "BM = (product_list['productSubGroupDescription'] == 'IMA') \n", - "product_list = product_list[BM]\n", - "\n", - "product_list.show_in_notebook(display_length=5)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Choose a single exposure file to work on. In this example, we choose the first exposure. To create usable data, you will have to follow this work flow on all individual IMA files in your dataset." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "myID = product_list['obsID'][0:1]" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Download the IMA and FLT files for that exposure. The standard pipeline-FLT will be used for comparison with the detrended final product." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "download = Observations.download_products(myID, mrp_only=False, productSubGroupDescription=['IMA', 'FLT'])" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Display the results of the download operation." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "download" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Read the files that were just downloaded locally. In addition, have the path be just the rootname, i.e. without the file extension." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "scrolled": true - }, - "outputs": [], - "source": [ - "localpathtofile = download['Local Path'][1][:-8]\n", - "localpathtofile\n", - "\n", - "original_ima = fits.open(localpathtofile+'ima.fits')\n", - "original_flt = fits.open(localpathtofile+'flt.fits')\n", - "original_ima.info()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Print the number of samples and plot the individual reads of the IMA file.\n", - "\n", - "**Note: the individual 'SCI' extensions are stored in reverse order, with 'SCI', 1 corresponding to the last read.**" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "nsamp = original_ima[0].header['NSAMP']\n", - "print('NSAMP', nsamp)\n", - "\n", - "fig, axarr = plt.subplots((nsamp+3)//4, 4, figsize=(9, 3*((nsamp+3)//4)))\n", - "\n", - "for i in range(1, 4*((nsamp+3)//4)+1):\n", - "\n", - " row = (i-1)//4\n", - " col = (i-1) % 4\n", - " if (i <= nsamp):\n", - " immed = np.nanmedian(original_ima['SCI', i].data)\n", - " stdev = np.nanstd(original_ima['SCI', i].data)\n", - " axarr[row, col].imshow(original_ima['SCI', i].data, clim=[immed-.3*stdev, immed+.5*stdev], cmap='Greys', origin='lower')\n", - " axarr[row, col].set_title('SCI '+str(i))\n", - " axarr[row, col].set_xticks([]) \n", - " axarr[row, col].set_yticks([]) \n", - " else:\n", - " fig.delaxes(axarr[row, col])\n", - "\n", - "fig.tight_layout()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "## Query CRDS for reference files \n", - "\n", - "Before running `reduce_dash`, we need to set some environment variables for several subsequent calibration tasks.\n", - "\n", - "We will point to a subdirectory called `crds_cache/` using the IREF environment variable. The IREF variable is used for WFC3 reference files. Other instruments use other variables, e.g., JREF for ACS." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "os.environ['CRDS_SERVER_URL'] = 'https://hst-crds.stsci.edu'\n", - "os.environ['CRDS_SERVER'] = 'https://hst-crds.stsci.edu'\n", - "os.environ['CRDS_PATH'] = './crds_cache'\n", - "os.environ['iref'] = 'iref/'\n", - "if not os.path.exists('iref'):\n", - " os.mkdir('iref')" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "The code block below will query CRDS for the best reference files currently available for these datasets and update the header keywords to point to these new files. We will use the Python package `os` to run terminal commands. In the terminal, the line would be:\n", - "\n", - "```crds bestrefs --files [filename] --sync-references=1 --update-bestrefs```\n", - "\n", - "...where 'filename' is the name of your fits file." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "ima_files = glob('*_ima.fits') \n", - "\n", - "for file in ima_files:\n", - " command_line_input = 'crds bestrefs --files {:} --sync-references=1 --update-bestrefs'.format(file)\n", - " os.system(command_line_input)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "## 3. Run DASH\n", - "Run the DASH pipeline for a single exposure. This procedure showcases the capabilities and customization options of the DASH pipeline.\n", - "\n", - "**Note: the following will only work if you are using the notebooks inside of the Notebook directory. `wfc3_dash` submodule will be properly packaged and installed within the `wfc3_tools` module sometime in the future.**\n", - "\n", - "If you move the notebooks and want to use them elsewhere, you can still provide a `temp_path` to the dash codes and remove the comments below. " - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "# import sys\n", - "\n", - "# module_path = os.path.abspath(os.path.join('..'))\n", - "# if module_path not in sys.path:\n", - "# sys.path.append(module_path)\n", - "\n", - "# tmp_path = \".../wfc3_dash/wfc3_dash\"\n", - "# if tmp_path not in sys.path:\n", - "# sys.path.append(tmp_path)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "### 3.1 Create DashData object\n", - "\n", - "We use the both IMA and FLT extensions of our local image to create a DashData object." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "myDash = DashData(localpathtofile+'ima.fits', flt_file_name=localpathtofile+'flt.fits')\n", - "print(myDash.root)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "### 3.2 Create difference files from reads\n", - "\n", - "A difference (diff) file contains the counts accumulated between two reads of the IMA file. The diff files are written to disk in a directory named `./diff` under the current working directory (cwd). In creating diff files from the readouts of the IMA, the first difference, between the 1-st and 0-th read is ignored because of its very short exposure time of 2.9 seconds, resulting in a noisy image. \n", - "\n", - "In order to create the best possible results, the `split_ima()` method uses the `bestrefs` function from [CRDS](https://hst-crds.stsci.edu/) to ensure all reference files are up to date and available." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "myDash.split_ima()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Print the number of diff files and plot the diff files." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "ndiff = len(myDash.diff_files_list)\n", - "print('Number of diff files', ndiff)\n", - "\n", - "if ndiff > 4: \n", - " fig,axarr = plt.subplots((ndiff+3)//4, 4, figsize=(9, 3*((ndiff+3)//4)))\n", - "\n", - " for i in range(4*((ndiff+3)//4)):\n", - "\n", - " row = (i)//4\n", - " col = (i) % 4\n", - " if (i < ndiff):\n", - " diff_i = fits.open(myDash.diff_files_list[i]+'_diff.fits')\n", - " immed = np.nanmedian(diff_i['SCI'].data)\n", - " stdev = np.nanstd(diff_i['SCI'].data)\n", - " axarr[row,col].imshow(diff_i['SCI'].data,clim=[immed-.3*stdev,immed+.5*stdev],cmap='Greys',origin='lower')\n", - " axarr[row,col].set_title('Diff:'+str(i+1))\n", - " axarr[row,col].set_xticks([]) \n", - " axarr[row,col].set_yticks([]) \n", - " else:\n", - " fig.delaxes(axarr[row,col])\n", - "else:\n", - " fig, axarr = plt.subplots(1, ndiff, figsize=(15, 15))\n", - " for i in range(ndiff):\n", - " immed = np.nanmedian(diff_i['SCI'].data)\n", - " stdev = np.nanstd(diff_i['SCI'].data)\n", - " diff_i = fits.open(myDash.diff_files_list[i]+'_diff.fits')\n", - " axarr[i].imshow(diff_i['SCI'].data, clim=[immed-.3*stdev, immed+.5*stdev], cmap='Greys', origin='lower')\n", - " axarr[i].set_title('Diff:'+str(i+1))\n", - " axarr[i].set_xticks([]) \n", - " axarr[i].set_yticks([]) \n", - "\n", - "fig.tight_layout()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "### 3.3 Create an association table\n", - "\n", - "This file mimics a typical association file for dithered exposures, which is used by AstroDrizzle to align and stack multiple exposures taken at the same sky position with small dithers. \n", - "\n", - "We exploit the fact that a WFC3/IR exposure taken under gyro control can be effectively split into individual pseudo-exposures (the diff images we created in [Section 3.2](#diff_files)). From there, AstroDrizzle can treat such pseudo-expsoures as individual dithers, and combine them into a single exposure." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "myDash.make_pointing_asn()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Show the content of the asn file, which was created in `./diff`." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "asn_filename = 'diff/{}_asn.fits'.format(myDash.root)\n", - "asn_table = Table(fits.getdata(asn_filename, ext=1))\n", - "asn_table.show_in_notebook()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "### 3.4 Create a segmentation map\n", - "\n", - "Create segmentation map from original FLT image to assist with background subtraction and fixing of cosmic ray flags using `create_seg_map()`. This method makes a directory called `./segmentation_maps`, which holds the outputs." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "myDash.create_seg_map()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Plot segmentation map." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "scrolled": true - }, - "outputs": [], - "source": [ - "rootname = myDash.root\n", - "segmap_name = ('segmentation_maps/'+ rootname + '_seg.fits')\n", - "segmap = fits.getdata(segmap_name)\n", - "\n", - "fig = plt.figure(figsize=(6, 8))\n", - "plt.title(segmap_name)\n", - "plt.imshow(segmap, origin='lower', vmin=0, vmax=1, cmap='Greys_r')" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Print and read source list." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "sourcelist_name = ('segmentation_maps/' + rootname + '_source_list.dat')\n", - "sourcelist = ascii.read(sourcelist_name)\n", - "print(sourcelist)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Let's create a segmentation map and source list from the difference files. We need to make source lists from our difference files created from the IMA so that `TweakReg` can better align these difference files to catalogs, each other, etc.\n", - "\n", - "First, generate a list of difference files that contain the full path name." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "diffpath = os.path.dirname(os.path.abspath('diff/{}_*_diff.fits'.format(rootname)))\n", - "cat_images = sorted([os.path.basename(x) for x in glob('diff/{}_*_diff.fits'.format(rootname))])\n", - "sc_diff_files = [diffpath + '/' + s for s in cat_images]" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Then, create a difference segmentation map using `diff_seg_map()` and the diff files." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "myDash.diff_seg_map(cat_images=sc_diff_files)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Plot the segmentation map from a diffrence files." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "segmap_name = ('segmentation_maps/' + rootname + '_01_diff_seg.fits')\n", - "segmap = fits.getdata(segmap_name)\n", - "fig = plt.figure(figsize=(6, 8))\n", - "plt.title(segmap_name)\n", - "plt.imshow(segmap, origin='lower', vmin=0, vmax=1, cmap='Greys_r')" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "### 3.5 Subtract background from difference files\n", - "\n", - "Subtract background from the individual reads taken from the original IMA file using the DRZ and SEG images produced in the background subtraction of the original FLT. \n", - "\n", - "By default, `subtract_background_reads()` will subtract the background and write it to the header. By setting `subtract=False`, the background will not be subtracted and will only be written to the header. In addition, setting `reset_stars_dq=True` will reset cosmic rays within objects to 0 since the centers of the stars are flagged." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "myDash.subtract_background_reads()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "### 3.6 Fix cosmic rays\n", - "\n", - "Now, we can use `fix_cosmic_rays()` to reset cosmic rays within the segmentation maps of objects and use [L.A.Cosmic](https://arxiv.org/pdf/astro-ph/0108003.pdf) to find them again." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "myDash.fix_cosmic_rays()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "### 3.7 Align reads to each other\n", - "Align reads from the IMA to one another by aligning each difference file to the first diff file.\n", - "\n", - "Listed below are all the parameters available to `myDash.align()`. `align()` uses `TweakReg` to update the WCS information in the headers of the diff files, then drizzles the images together using `AstroDrizzle`. There are more parameters available to users when working with `TweakReg` and `AstroDrizzle` that could be an integral part of the workflow for users of DASH. The example below lists the default values set for every input:\n", - "\n", - "``myDash.align(self, subtract_background = True, \n", - " align_method = None, \n", - " ref_catalog = None, \n", - " create_diff_source_lists = True,\n", - " updatehdr = True, \n", - " updatewcs = True, \n", - " wcsname = 'DASH', \n", - " threshold = 50., \n", - " cw = 3.5, \n", - " searchrad = 20., \n", - " astrodriz = True, \n", - " cat_file = 'catalogs/diff_catfile.cat',\n", - " drz_output = None, \n", - " move_files = False)``\n", - " \n", - "Refer to documentation to customize parameters for [TweakReg](https://drizzlepac.readthedocs.io/en/latest/tweakreg.html) and [AstroDrizzle](https://drizzlepac.readthedocs.io/en/latest/astrodrizzle.html). \n", - "\n", - "Note: the error `UnboundLocalError: local variable 'sig' referenced before assignment` can be solved by lowering threshold parameter." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "myDash.align(updatehdr=False, updatewcs=True, astrodriz=False)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Print the shifts file to analyze how well the alignment went. Do not update header until shifts, as seen in the `xrms` and `yrms` columns, are satisfactory. Further information about the outputs in the shift file and what is 'satisfactory' can be found in the [Drizzlepac Handbook](https://hst-docs.stsci.edu/drizzpac)." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "shift_file = glob('shifts/shifts_*.txt')\n", - "shift_file_name = shift_file[0]\n", - "\n", - "\n", - "shift_tab = Table.read(shift_file_name,\n", - " format='ascii.no_header',\n", - " names=['file','dx','dy','rot','scale','xrms','yrms'])\n", - "\n", - "formats = ['.2f', '.2f', '.3f', '.5f', '.2f', '.2f']\n", - "for i, col in enumerate(shift_tab.colnames[1:]):\n", - " shift_tab[col].format = formats[i]\n", - "shift_tab" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "Let's align our images with a threshold of 20, and update the headers and WCS information." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "myDash.align(threshold = 20.)" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "## 4. Plot original IMA and DASH pipeline science result\n", - "\n", - "Plot the final DRZ image and compare to the original IMA." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": { - "scrolled": true - }, - "outputs": [], - "source": [ - "sci_name = myDash.root + '_drz_sci.fits'\n", - "og_flt_name = 'mastDownload/HST/' + myDash.root + '/' + myDash.root + '_ima.fits'\n", - "sci = fits.getdata(sci_name)\n", - "og_flt = fits.getdata(og_flt_name)\n", - "\n", - "fig = plt.figure(figsize=(9, 4))\n", - "ax1 = fig.add_subplot(1,2,2)\n", - "ax2 = fig.add_subplot(1,2,1)\n", - "\n", - "ax1.set_title('DASH Pipeline Reduced Science File')\n", - "ax2.set_title('Original IMA (not reduced using pipeline)')\n", - "\n", - "ax1.set_xlim(-10,1120)\n", - "ax2.set_xlim(-10,1120)\n", - "\n", - "ax1.set_ylim(-10,1050)\n", - "ax2.set_ylim(-10,1050)\n", - "\n", - "ax1.imshow(sci, vmin=0, vmax=40, cmap='Greys_r', origin='lower')\n", - "ax2.imshow(og_flt, vmin=0, vmax=40, cmap='Greys_r', origin='lower')\n", - "plt.show()" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "\n", - "## 5. Conclusions\n", - "\n", - "Thank you for walking through this notebook. Now using WFC3 data, you should be more familiar with:\n", - "\n", - "- Creating difference files, association tables, and segmentation maps using `wfc3_dash`.\n", - "- Subtracting background and fixing cosmic rays from newly generated FLTs.\n", - "- Aligning reads to each other for a final product.\n", - "\n", - "Congratulations, you have completed the notebook!" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "## \n", - "## Additional Resources\n", - "Below are some additional resources that may be helpful. Please send any questions through the [HST Helpdesk](https://stsci.service-now.com/hst).\n", - "\n", - "- [WFC3 Website](https://www.stsci.edu/hst/instrumentation/wfc3)\n", - "- [WFC3 Instrument Handbook](https://hst-docs.stsci.edu/wfc3ihb)\n", - "- [WFC3 Data Handbook](https://hst-docs.stsci.edu/wfc3dhb)\n", - " - see sections 9.5.4 for reference to this notebook\n", - " \n", - "\n", - "## About this Notebook\n", - "\n", - "**Authors:** Catherine Martlin; WFC3 Instrument Team\n", - "\n", - "**Updated on:** 2023-02-28\n", - "\n", - "\n", - "## Citations\n", - "\n", - "If you use `numpy`, `matplotlib`, `astropy`, `astroquery`, `photutils`, or `drizzlepac` for published research, please cite the authors. Follow these links for more information about citing the libraries below:\n", - "\n", - "* [Citing `numpy`](https://numpy.org/citing-numpy/)\n", - "* [Citing `matplotlib`](https://matplotlib.org/stable/users/project/citing.html)\n", - "* [Citing `astropy`](https://www.astropy.org/acknowledging.html)\n", - "* [Citing `astroquery`](https://astroquery.readthedocs.io/en/latest/license.html)\n", - "* [Citing `photutils`](https://photutils.readthedocs.io/en/stable/citation.html)\n", - "* [Citing `drizzlepac`](https://drizzlepac.readthedocs.io/en/latest/LICENSE.html)\n", - "\n", - "***\n", - "[Top of Page](#title)\n", - "\"Space " - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [] - } - ], - "metadata": { - "kernelspec": { - "display_name": "Python 3 (ipykernel)", - "language": "python", - "name": "python3" - }, - "language_info": { - "codemirror_mode": { - "name": "ipython", - "version": 3 - }, - "file_extension": ".py", - "mimetype": "text/x-python", - "name": "python", - "nbconvert_exporter": "python", - "pygments_lexer": "ipython3", - "version": "3.11.8" - }, - "toc": { - "base_numbering": 1, - "nav_menu": {}, - "number_sections": false, - "sideBar": true, - "skip_h1_title": false, - "title_cell": "Table of Contents", - "title_sidebar": "Contents", - "toc_cell": false, - "toc_position": {}, - "toc_section_display": true, - "toc_window_display": true - } - }, - "nbformat": 4, - "nbformat_minor": 4 -} diff --git a/notebooks/WFC3/dash/pre-requirements.txt b/notebooks/WFC3/dash/pre-requirements.txt deleted file mode 100644 index 43c5f747a..000000000 --- a/notebooks/WFC3/dash/pre-requirements.txt +++ /dev/null @@ -1,2 +0,0 @@ -cython -numpy==1.23.4 diff --git a/notebooks/WFC3/dash/reduce_dash.py b/notebooks/WFC3/dash/reduce_dash.py deleted file mode 100644 index 796014599..000000000 --- a/notebooks/WFC3/dash/reduce_dash.py +++ /dev/null @@ -1,1047 +0,0 @@ -#! /usr/bin/env python - -""" Reduces DASH/IR data. - -This script serves as a pipeline to reduce DASH/IR data. The -FLT/IMA files available from MAST are inappropiate to use to -perform science upon (Momcheva et al. 2017). This script will -create an aligned image that is the combination of the -differences between adjacent reads. This file is appropriate -to use for science work. - - -The following output products are created by this script from -each IMA file: - -1. Individual fits files of differences between IMA reads with - thier own headers. -2. A master ASN file of those read fits files for a single IMA. -3. An aligned version of the IMA read fits files. - -Authors -------- - Rosalia O'Brien 2019 - Iva Momcheva 2018 - Catherine Martlin 2018, 2023 - Mario Gennaro 2018 - -Use ---- - This script is intended to be executed via the command line - as such: - :: - - python reduce_dash.py [-f|--file] - - ``-f --file`` - The IMA file name/path. - -Notes ------ - - Lacosmicx (used to identify cosmic rays): https://github.com/cmccully/lacosmicx - -References ----------- - - - Ivelina G. Momcheva, Pieter G. van Dokkum, Arjen van der Wel, - Gabriel B. Brammer, John MacKenty, Erica J. Nelson, Joel Leja, - Adam Muzzin, and Marijn Franx. 2017 January. - doi:10.1088/1538-3873/129/971/015004 - -""" -from glob import glob -import os - -from astropy.convolution import convolve -from astropy.convolution import Gaussian2DKernel -from astropy.io import ascii -from astropy.io import fits -from astropy.stats import gaussian_fwhm_to_sigma -from astropy.table import Table -from drizzlepac import tweakreg -from drizzlepac import astrodrizzle -import lacosmicx # for fix_cosmic_rays -import numpy as np -from stsci.tools import teal -import stwcs -from photutils import detect_sources # Needed for create_seg_map -from photutils import detect_threshold # Needed for create_seg_map -from photutils.segmentation import SourceCatalog - -from utils import get_flat -from utils import get_IDCtable - - -class DashData(object): - def __init__(self, file_name=None, flt_file_name=None): - """ - The init method performs a series of tests to make sure that the file - fed to the DashData class is a valid IMA file with units in e/s. Will - also add root attribute to DashData class. - - Paramaters - ---------- - self : object - DashData object created from an individual IMA file. - file_name : fits file - Name of IMA file that is fed to DashData class. - flt_file_name : fits file - Name of FLT file that is fed to DashData class. - - Outputs - ------- - N/A - """ - - if "_ima.fits" not in file_name: - raise Exception("Input needs to be an IMA file.") - else: - self.file_name = file_name - # First test whether the file exists - try: - self.ima_file = fits.open(self.file_name) - ### If it does - test its content - - ###### Test whether it is a wfc3/ir image - if ( - (self.ima_file[0].header["INSTRUME"].strip() == "WFC3") - & (self.ima_file[0].header["DETECTOR"].strip() == "IR") - ) == False: - raise Exception( - "This observation was not performed with WFC3/IR, instrument is set to: {}, and detector is set to: {}".format( - self.ima_file[0].header["INSTRUME"], - self.ima_file[0].header["DETECTOR"], - ) - ) - - ###### Test whether it has more than one science extension (FLTs have only 1) - nsci = [ext.header["EXTNAME "] for ext in self.ima_file[1:]].count( - "SCI" - ) - if nsci == 1: - raise Exception( - "This file has only one science extension, it cannot be a WFC3/IR ima" - ) - - ###### Test that this file is NOT a RAW - calib_keys = [ - "DQICORR", - "ZSIGCORR", - "ZOFFCORR", - "DARKCORR", - "BLEVCORR", - "NLINCORR", - "FLATCORR", - "CRCORR", - "UNITCORR", - "PHOTCORR", - "RPTCORR", - "DRIZCORR", - ] - performed = 0 - for ck in calib_keys: - if self.ima_file[0].header[ck] == "COMPLETE": - performed = performed + 1 - if performed == 0: - raise Exception("This file looks like a RAW file") - - ###### Test that the units of the individual images are e/s - - bu = self.ima_file[1].header["BUNIT"] - if bu != "ELECTRONS/S": - uc = self.ima_file[0].header["UNITCORR"] - fc = self.ima_file[0].header["FLATCORR"] - raise Exception( - 'BUNIT in the "SCI" extensions of this file is set to "{}", but should be set to "ELECTRONS/S"\n' - 'This is a consequence of UNITCORR set to "{}" and FLATCORR set to "{}".\n' - 'Please rerun calwf3 on this file after setting both UNITCORR and FLATCORR to "PERFORM" in the 0-th extension header'.format( - bu, uc, fc - ) - ) - - except IOError: - print("Cannot read file.") - - self.flt_file_name = flt_file_name - self.root = self.file_name.split("/")[-1].split("_ima")[0] - - def align( - self, - subtract_background=True, - align_method=None, - ref_catalog=None, - create_diff_source_lists=True, - updatehdr=True, - updatewcs=True, - wcsname="DASH", - threshold=50.0, - cw=3.5, - searchrad=20.0, - astrodriz=True, - cat_file="catalogs/diff_catfile.cat", - drz_output=None, - move_files=False, - ): - - """ - Aligns new FLT's to reference catalog. - - Parameters - ---------- - self : object - DashData object created from an individual IMA file. - subtract_background : bool, optional - If True, runs subtract_background_reads functions. - align_method : string, optional - Defines alignment method to be used. Default is None (input files - will align to each other). - ref_catalog : cat file, optional - Defines reference image that will be referenced for CATALOG - alignment method. - create_diff_source_lists : bool, optional - Specifies whether or not to create a segmentation image and source - list from the difference files. - updatehdr : bool, optional - Specifies whether to update the headers after aligning during - TweakReg. Default is True. - updatewcs : bool, optional - Specifies whether to update the WCS after aligning during TweakReg. - Default is True. - wcsname : str, optional (as long as name you choose doesn't already exist) - Specifies name of WCS. Default is 'DASH'. - threshold : float, optional - The object detection threshold above the local background in units of - sigma used by TweakReg. - cw : float, optional - conv_width or the convolution kernel width in pixels used by TweakReg. - searchrad : float, optional - Radius (in pixels) that TweakReg will search around sources to find - matches. Default is 20 pixels. - astrodriz : bool, optional - Specifies whether to drizzle images together using Astrodrizzle. - Default is True. - cat_file : str, optional - Name of catfile to be used to align sources in TweakReg. Default is - the catfile created by setting create_diff_source_lists to True, - catalogs/diff_catfile.cat - drz_output : str, optional - Name of output file after drizzling using AstroDrizzle. Default is - the root name of the original IMA. - move_files : bool, optional - If True, move files from alignment steps to folders. - - - Outputs - ------- - Shifts file : txt file - File containing shifts during TweakReg. - WCS Shifts file : fits file - File containing WCS shifts during TweakReg. - Residual plot : png - Plot showing residuals from TweakReg alignment. - Vector plot : png - Plot showing vector poitnings from TweakReg alignment. - Aligned FLT's : fits file - Same diff files created in split_ima that have been aligned using TweakReg. - Drizzled Image : fits file - Setting astrodriz to True will output a drizzled image form a single - exposure (produced only if astrodriz is set to True). - """ - - # Set name for output drizzled image to the rootname of the original IMA if it is not specified - if drz_output is None: - drz_output = self.root - - input_images = sorted(glob("diff/{}_*_diff.fits".format(self.root))) - - if not os.path.exists("shifts"): - os.mkdir("shifts") - outshifts = "shifts/shifts_{}.txt".format(self.root) - outwcs = "shifts/shifts_{}_wcs.fits".format(self.root) - - if subtract_background: - self.subtract_background_reads() - - # Align images to a catalog - if align_method == "CATALOG": - if ref_catalog is not None: - - ##Create source list and segmentation maps based on difference files - if create_diff_source_lists is True: - diffpath = os.path.dirname( - os.path.abspath("diff/{}_*_diff.fits".format(self.root)) - ) - cat_images = sorted( - [ - os.path.basename(x) - for x in glob("diff/{}_*_diff.fits".format(self.root)) - ] - ) - - sc_diff_files = [diffpath + "/" + s for s in cat_images] - - self.diff_seg_map(cat_images=sc_diff_files) - - teal.unlearn("tweakreg") - teal.unlearn("imagefindpars") - - tweakreg.TweakReg( - input_images, - refcat=ref_catalog, - catfile=cat_file, - xcol=2, - ycol=3, - updatehdr=updatehdr, - updatewcs=updatewcs, - wcsname=wcsname, - verbose=True, - imagefindcfg={"threshold": threshold, "conv_width": cw}, - searchrad=searchrad, - searchunits="pixels", - shiftfile=True, - outshifts=outshifts, - outwcs=outwcs, - interactive=False, - fitgeometry="rscale", - minobj=5, - ) - - else: - - raise Exception( - "Need to specify reference catalog and reference image." - ) - - # Align images to the first image - else: - - teal.unlearn("tweakreg") - teal.unlearn("imagefindpars") - - tweakreg.TweakReg( - input_images, - catfile=cat_file, - xcol=2, - ycol=3, - updatehdr=updatehdr, - updatewcs=updatewcs, - wcsname=wcsname, - verbose=True, - imagefindcfg={"threshold": threshold, "conv_width": cw}, - searchrad=searchrad, - searchunits="pixels", - shiftfile=True, - outshifts=outshifts, - outwcs=outwcs, - interactive=False, - fitgeometry="rscale", - minobj=5, - ) - - pass - - # Drizzle the images together - if astrodriz is True: - - # Do not have drizzle take 256 flags into account - no_tfs = 2, 4, 8, 16, 32, 64, 128, 512, 2048, 4096, 8192, 16384 - - astrodrizzle.AstroDrizzle( - input_images, - output=drz_output, - clean=False, - final_pixfrac=1.0, - context=False, - resetbits=0, - preserve=False, - driz_cr_snr="8.0 5.0", - driz_cr_scale="2.5 0.7", - driz_sep_bits=no_tfs, - final_bits=no_tfs, - num_cores=1, - ) # added num cores = 1 - - if move_files is True: - self.move_files() - - def create_seg_map(self): - """ - Creates segmentation map, from original FLT file, that is used in - background subtraction and to fix cosmic rays. - - Parameters - ---------- - self : object - DashData object created from an individual IMA file. - - Output - ------ - Segmentation Image : fits file - Segmentation map. - Source List : .dat file - List of sources and their properties. - """ - - flt = fits.open(self.flt_file_name) - data = flt[1].data - - threshold = detect_threshold(data, nsigma=3.0) - - sigma = 3.0 * gaussian_fwhm_to_sigma # FWHM = 3. - kernel = Gaussian2DKernel(sigma, x_size=3, y_size=3) - kernel.normalize() - convolved_data = convolve(data, kernel) - segm = detect_sources(convolved_data, threshold, npixels=10) - - hdu = fits.PrimaryHDU(segm.data) - if not os.path.exists("segmentation_maps"): - os.mkdir("segmentation_maps") - hdu.writeto(("segmentation_maps/{}_seg.fits").format(self.root), overwrite=True) - - # Create source list - cat = SourceCatalog(data, segm) - - tbl = cat.to_table() - tbl["xcentroid"].info.format = ".2f" - tbl["ycentroid"].info.format = ".2f" - # tbl['cxx'].info.format = '.2f' - # tbl['cxy'].info.format = '.2f' - # tbl['cyy'].info.format = '.2f' - - ascii.write( - tbl, - "segmentation_maps/{}_source_list.dat".format(self.root), - overwrite=True, - ) - - def diff_seg_map( - self, cat_images=None, remove_column_names=True, nsigma=1.0, sig=5.0, npixels=5 - ): - """ - Creates segmentation image and source list from difference files. - - Parameters - ---------- - self : object - DashData object created from an individual IMA file. - cat_images : list, str - List of difference files with full path name. - remove_column_names : bool - Specifies whether to remove the header from the source lists so - TweakReg can read them. - nsigma : float - The number of standard deviations per pixel above the background - for which to consider a pixel as possibly being part of a source. - Used by photutils.detect_threshold. - sig : float - Factor used to adjust gaussian_fwhm_to_sigma. - npixels : float - The number of connected pixels, each greater than threshold, that - an object must have to be detected. npixels must be a positive . - - Outputs - ------- - Segmentation Image : fits file - Segmentation map. - Source List : .dat file - List of sources and their properties. - """ - - input_images = sorted(glob("diff/{}_*_diff.fits".format(self.root))) - - for index, diff in enumerate(input_images, start=1): - - diff = fits.open(diff) - data = diff[1].data - - threshold = detect_threshold(data, nsigma=nsigma) - - sigma = sig * gaussian_fwhm_to_sigma - kernel = Gaussian2DKernel(sigma, x_size=sig, y_size=sig) - kernel.normalize() - convolved_data = convolve(data, kernel) - segm = detect_sources(convolved_data, threshold, npixels=npixels) - - hdu = fits.PrimaryHDU(segm.data) - if not os.path.exists("segmentation_maps"): - os.mkdir("segmentation_maps") - hdu.writeto( - ("segmentation_maps/{}_{:02d}_diff_seg.fits").format(self.root, index), - overwrite=True, - ) - - # Create source list - cat = SourceCatalog(data, segm) - - tbl = cat.to_table() - tbl["xcentroid"].info.format = ".2f" - tbl["ycentroid"].info.format = ".2f" - # tbl['cxx'].info.format = '.2f' - # tbl['cxy'].info.format = '.2f' - # tbl['cyy'].info.format = '.2f' - - ascii.write( - tbl, - "segmentation_maps/{}_{:02d}_diff_source_list.dat".format( - self.root, index - ), - overwrite=True, - ) - - if remove_column_names is True: - # Remove headers to source lists so tweakreg can read them - n = 1 - nfirstlines = [] - - with open( - "segmentation_maps/{}_{:02d}_diff_source_list.dat".format( - self.root, index - ) - ) as f, open("temp_sl.dat", "w") as out: - for x in range(n): - nfirstlines.append(next(f)) - for line in f: - out.write(line) - - os.remove( - "segmentation_maps/{}_{:02d}_diff_source_list.dat".format( - self.root, index - ) - ) - os.rename( - "temp_sl.dat", - "segmentation_maps/{}_{:02d}_diff_source_list.dat".format( - self.root, index - ), - ) - - if cat_images is not None: - x = np.array(cat_images) - y = np.array( - sorted( - glob( - ("segmentation_maps/{}_*_diff_source_list.dat").format( - self.root - ) - ) - ) - ) - catdata = Table([x, y], names=["Diff File", "Source List"]) - # Save catfile in catalogs folder - if not os.path.exists("catalogs"): - os.mkdir("catalogs") - ascii.write(catdata, "catalogs/diff_catfile.cat", overwrite=True) - else: - raise Exception( - "Need to input list of difference files in order to make source list. List should include full path." - ) - - def fix_cosmic_rays(self, rm_custom=False, flag=None, **lacosmic_param): - """ - Resets cosmic rays within the seg maps of objects and uses L.A.Cosmic - to find them again. - - Parameters - ---------- - self : object - DashData object created from an individual IMA file. - rm_custom : bool - Specifies whether or not the user would like to remove custom flags - within the boundaries of sources, as defined by the segmentation map - created from the original FLT. - flag : int - Specifies flag the user would like the remove within the boundaries - of sources. - lacosmic_param : dic - Dictionary of the L.A.Cosmic parameters that users may want to specify. - If not set, then presets are used. - - Output - ------ - Fixed for cosmic rays diff files : fits - Same diff files created in split_ima that have now been corrected - for cosmic ray errors. - """ - - asn_exposures = sorted(glob("diff/" + self.root + "_*_diff.fits")) - - seg = fits.open("segmentation_maps/{}_seg.fits".format(self.root)) - seg_data = np.cast[np.float32](seg[0].data) - - flt_full = fits.open(self.flt_file_name) - flt_full_wcs = stwcs.wcsutil.HSTWCS(flt_full, ext=1) - - EXPTIME = flt_full[0].header["EXPTIME"] - - if lacosmic_param: - gain = lacosmic_param["gain"] - readnoise = lacosmic_param["readnoise"] - objlim = lacosmic_param["objlim"] - pssl = lacosmic_param["pssl"] - verbose = lacosmic_param["verbose"] - else: - gain = 1.0 - readnoise = 20.0 - objlim = 15.0 - pssl = 0.0 - verbose = True - # Have lacosmicx locate cosmic rays - crmask, clean = lacosmicx.lacosmicx( - flt_full[1].data, - gain=gain, - readnoise=readnoise, - objlim=objlim, - pssl=pssl, - verbose=verbose, - ) - - yi, xi = np.indices((1014, 1014)) - - # Remove all 4096 flags within the boundaries of objects - for exp in asn_exposures: - flt = fits.open(exp, mode="update") - flagged_stars = ((flt["DQ"].data & 4096) > 0) & (seg_data > 0) - flt["DQ"].data[flagged_stars] -= 4096 - new_cr = ( - (crmask == 1) - & ((flt["DQ"].data & 4096) == 0) - & ((seg_data == 0) | ((seg_data > 0) & (flt["SCI"].data < 1.0))) - & (xi > 915) - & (yi < 295) - ) - flt["DQ"].data[new_cr] += 4096 - flt.flush() - - # Remove custom flags - if rm_custom is True: - if flag is not None: - for exp in asn_exposures: - flt = fits.open(exp, mode="update") - flagged_stars = ((flt["DQ"].data & flag) > 0) & (seg_data > 0) - flt["DQ"].data[flagged_stars] -= flag - new_cr = ( - (crmask == 1) - & ((flt["DQ"].data & flag) == 0) - & ((seg_data == 0) | ((seg_data > 0) & (flt["SCI"].data < 1.0))) - & (xi > 915) - & (yi < 295) - ) - flt["DQ"].data[new_cr] += flag - flt.flush() - else: - raise Exception("Must specify which flags to remove.") - - def make_pointing_asn(self): - """ - Makes a new association table for the reads extracted from a given IMA. - - Parameters - ---------- - self : object - DashData object created from an individual IMA file. - - Outputs - ---------- - ASN files : fits file - Fits file formatted as an association file that holds - the names of the difference files created by split_ima - and the root name of the individual IMA file. - """ - - asn_filename = "diff/{}_asn.fits".format(self.root) - asn_list = self.diff_files_list.copy() - asn_list.append(self.root) - - # Create Primary HDU: - hdr = fits.Header() - hdr["FILENAME"] = asn_filename - hdr["FILETYPE"] = "ASN_TABLE" - hdr["ASN_ID"] = self.root - hdr["ASN_TABLE"] = asn_filename - hdr[ - "COMMENT" - ] = "This association table is for the read differences for the IMA." - primary_hdu = fits.PrimaryHDU(header=hdr) - - # Create the information in the asn file - num_mem = len(asn_list) - - asn_mem_names = np.array(asn_list) - asn_mem_types = np.full(num_mem, "EXP-DTH", dtype=np.chararray) - asn_mem_types[-1] = "PROD-DTH" - asn_mem_prsnt = np.ones(num_mem, dtype=np.bool_) - asn_mem_prsnt[-1] = 0 - - hdu_data = fits.BinTableHDU().from_columns( - [ - fits.Column(name="MEMNAME", format="14A", array=asn_mem_names), - fits.Column(name="MEMTYPE", format="14A", array=asn_mem_types), - fits.Column(name="MEMPRSNT", format="L", array=asn_mem_prsnt), - ] - ) - - # Create the final asn file - hdu = fits.HDUList([primary_hdu, hdu_data]) - - if "EXTEND" not in hdu[0].header.keys(): - hdu[0].header.update("EXTEND", True, after="NAXIS") - - hdu.writeto(asn_filename, overwrite=True) - - # Create property of the object that is the asn filename. - self.asn_filename = asn_filename - - def move_files(self): - """ - Moves files from alignment steps to folders. - - Parameters - ---------- - self : object - DashData object created from an individual IMA file. - - Outputs - ------- - residuals : dir - Directory containing all residual and vector plots. - histograms : dir - Directory containing all histogram png's. - misC_tweakreg_files : dir - Directory containing all misc. files. - """ - - # Move all residual plots into residuals folder - if not os.path.exists("residuals"): - os.mkdir("residuals") - - residuals = sorted(glob("residuals_{}_*_diff.png".format(self.root))) - for image in residuals: - os.rename(image, "residuals/" + image) - - vectors = sorted(glob("vector_{}_*_diff.png".format(self.root))) - for image in vectors: - os.rename(image, "residuals/" + image) - - # Move all histograms plots into histograms folder - if not os.path.exists("histograms"): - os.mkdir("histograms") - - graphs = sorted(glob("hist2d_{}_*_diff.png".format(self.root))) - for image in graphs: - os.rename(image, "histograms/" + image) - - # Move misc TweakReg files to misc_tweakreg folder - if not os.path.exists("misc_tweakreg_files"): - os.mkdir("misc_tweakreg_files") - - files = sorted(glob("{}_*_diff*".format(self.root))) - for image in files: - os.rename(image, "misc_tweakreg_files/" + image) - fits = sorted(glob("{}_*_mask*".format(self.root))) - for image in fits: - os.rename(image, "misc_tweakreg_files/" + image) - - def split_ima(self): - """ - Will create individual files for the difference between - adjacent reads of a IMA file. Will also add more attributes - to the DashData object. - - Parameters - ---------- - self : object - DashData object created from an individual IMA file. - - Outputs - ---------- - N files : fits files - Fits files of the difference between adjacent IMA reads. - - """ - FLAT = fits.open(get_flat(self.file_name)) - IDCtable = fits.open(get_IDCtable(self.file_name)) - - NSAMP = self.ima_file[0].header["NSAMP"] - shape = self.ima_file["SCI", 1].shape - - cube = np.zeros((NSAMP, shape[0], shape[1]), dtype="float32") - dq = np.zeros((NSAMP, shape[0], shape[1]), dtype=np.int32) - time = np.zeros(NSAMP) - - for i in range(NSAMP): - cube[NSAMP - 1 - i, :, :] = ( - self.ima_file["SCI", i + 1].data - * self.ima_file["TIME", i + 1].header["PIXVALUE"] - ) - dq[NSAMP - 1 - i, :, :] = self.ima_file["DQ", i + 1].data - time[NSAMP - 1 - i] = self.ima_file["TIME", i + 1].header["PIXVALUE"] - - diff = np.diff(cube, axis=0) - self.diff = diff[1:] - dt = np.diff(time) - self.dt = dt[1:] - - self.dq = dq[1:] - - self.readnoise_2D = np.zeros((1024, 1024), dtype="float32") - self.readnoise_2D[512:, 0:512] += self.ima_file[0].header["READNSEA"] - self.readnoise_2D[0:512, 0:512] += self.ima_file[0].header["READNSEB"] - self.readnoise_2D[0:512, 512:] += self.ima_file[0].header["READNSEC"] - self.readnoise_2D[512:, 512:] += self.ima_file[0].header["READNSED"] - self.readnoise_2D = self.readnoise_2D**2 - - self.diff_files_list = [] - for j in range(1, NSAMP - 1): - - hdu0 = fits.PrimaryHDU(header=self.ima_file[0].header) - hdu0.header["EXPTIME"] = dt[j] - hdu0.header["IMA2FLT"] = (1, "FLT {} extracted from IMA file".format(j)) - hdu0.header["NEXTEND"] = 5 - hdu0.header["OBSMODE"] = "ACCUM" - hdu0.header["NSAMP"] = 1 - - science_data = diff[j, :, :] / dt[j] - hdu1 = fits.ImageHDU( - data=science_data[5:-5, 5:-5], - header=self.ima_file["SCI", NSAMP - j - 1].header, - name="SCI", - ) - hdu1.header["EXTVER"] = 1 - hdu1.header["ROOTNAME"] = "{}_{:02d}".format(self.root, j) - - var = 2 * self.readnoise_2D + science_data * FLAT["SCI"].data * dt[j] - err = np.sqrt(var) / dt[j] - - hdu2 = fits.ImageHDU( - data=err[5:-5, 5:-5], - header=self.ima_file["ERR", NSAMP - j - 1].header, - name="ERR", - ) - hdu2.header["EXTVER"] = 1 - - hdu3 = fits.ImageHDU( - data=dq[j + 1][5:-5, 5:-5], - header=self.ima_file["DQ", NSAMP - j - 1].header, - name="DQ", - ) - hdu3.header["EXTVER"] = 1 - - # Turn the 8192 cosmic ray flag to the standard 4096 - hdu3.data[(hdu3.data & 8192) > 0] -= 4096 - # remove the 32 flag, these are not consistently bad - hdu3.data[(hdu3.data & 32) > 0] -= 32 - - hdu4 = fits.ImageHDU( - data=(np.zeros((1014, 1014), dtype=np.int16) + 1), - header=self.ima_file["SAMP", NSAMP - j - 1].header, - name="SAMP", - ) - hdu5 = fits.ImageHDU( - np.zeros((1014, 1014)) + dt[j], - header=self.ima_file["TIME", NSAMP - j - 1].header, - name="TIME", - ) - hdu4.header["EXTVER"] = 1 - hdu5.header["EXTVER"] = 1 - - hdu = fits.HDUList([hdu0, hdu1, hdu2, hdu3, hdu4, hdu5]) - print("Writing {}_{:02d}_diff.fits".format(self.root, j)) - - self.hdu = hdu - - if not os.path.exists("diff"): - os.mkdir("diff") - - hdu.writeto("diff/{}_{:02d}_diff.fits".format(self.root, j), overwrite=True) - - self.diff_files_list.append("diff/{}_{:02d}".format(self.root, j)) - - def subtract_background_reads(self, subtract=True, reset_stars_dq=False): - """ - Performs median background subtraction for each individual difference file. - Uses the DRZ and SEG images produced in FLT background subtraction. - - Parameters - ---------- - self : object - DashData object created from an individual IMA file. - subtract : bool - Does not subtract the background by default, but still writes it to - the header. - Set to True to subtract background. - reset_stars_dq : bool - Set to True to reset cosmic rays within objects to 0 because the - centers of stars are flagged. - - Outputs - ------- - Background Subtracted N files : fits files - Fits files of the difference between adjacent IMA reads that have - been background subtracted. - """ - - seg = fits.open("segmentation_maps/{}_seg.fits".format(self.root)) - seg_data = np.cast[np.float32](seg[0].data) - - yi, xi = np.indices((1014, 1014)) - - self.bg_models = [] - - for ii, exp in enumerate(self.diff_files_list): - - diff = fits.open("{}_diff.fits".format(exp), mode="update") - diff_wcs = stwcs.wcsutil.HSTWCS(diff, ext=1) - - mask = ( - (seg_data == 0) - & (diff["DQ"].data == 0) - & (diff[1].data > -1) - & (xi > 10) - & (yi > 10) - & (xi < 1004) - & (yi < 1004) - ) - mask &= diff[1].data < 5 * np.median(diff[1].data[mask]) - data_range = np.percentile(diff[1].data[mask], [2.5, 97.5]) - mask &= (diff[1].data >= data_range[0]) & (diff[1].data <= data_range[1]) - data_range = np.percentile(diff[2].data[mask], [0.05, 99.5]) - mask &= (diff[2].data >= data_range[0]) & (diff[2].data <= data_range[1]) - - sky_level = np.median(diff[1].data[mask]) - model = diff[1].data * 0.0 + sky_level - - diff[1].header["MDRIZSKY"] = sky_level - if not subtract: - if "BG_SUB" not in (key for key in diff[1].header.keys()): - flt[1].header["BG_SUB"] = "No" - else: - if "BG_SUB" not in (key for key in diff[1].header.keys()): - diff[1].data -= model - diff[1].header["BG_SUB"] = "Yes" - else: - print( - "Keyword BG_SUB already set to {}. Skipping background subtraction.".format( - diff[1].header["BG_SUB"] - ) - ) - - if reset_stars_dq: - flagged_stars = ((diff["DQ"].data & 4096) > 0) & (blotted_seg > 0) - diff["DQ"].data[flagged_stars] -= 4096 - - diff.flush() - print("Background subtraction, {}_diff.fits: {}".format(exp, sky_level)) - - -def main( - ima_file_name=None, - flt_file_name=None, - align_method=None, - ref_catalog=None, - drz_output=None, - subtract_background=False, - wcsname="DASH", - threshold=50.0, - cw=3.5, - updatehdr=True, - updatewcs=True, - searchrad=20.0, - astrodriz=True, - cat_file="catalogs/diff_catfile.cat", -): - - """ - Runs entire DashData pipeline under a single function. - - Parameters - ---------- - ima_file_name : str - File name of ima file. - flt_file_name : str - File name of FLT file. - align_method : str, optional - Method to align difference files using TweakReg. Default is None, which - aligns reads to the first read. - Setting align_method equal to 'CATALOG' will align the reads to a catalog. - ref_catalog : str, optional - Catalog to be aligned to if using CATALOG align method. - drz_output : str, optional - Name to call final drizzled output image ( + '_drz_sci.fits' ). - Default is root_name + '_drz_sci.fits'. - subtract_background : bool - Determines whether background is subtracted or not during align - function. - Default is False since during this main function, - the background is subtracted separately. - wcsname : str - Name for WCS during TweakReg. - threshold : float - The object detection threshold above the local background in units of - sigma used by TweakReg. - cw : float - conv_width or the convolution kernel width in pixels used by TweakReg. - updatehdr : bool - Determines whether to update headers during TweakReg. Default is True. - updatewcs : bool - Determines whether to update WCS information during TweakReg. Default - is True. - searchrad : float - The search radius for a match used by TweakReg. - astrodriz : bool - Determines whether or not to run astrodrizzle. Default is True. - cat_file : str, optional - Name of catfile to be used to align sources in TweakReg. Default is - the catfile created by setting create_diff_source_lists to True, - catalogs/diff_catfile.cat - - Outputs - ------- - N files : fits files - Fits files of the difference between adjacent IMA reads. - N files : fits files - Fits files of the difference between adjacent IMA reads. - Shifts file : txt file - File containing shifts during TweakReg. - WCS Shifts file : fits filr - File containing WCS shifts during TweakReg. - Drizzled science image : fits - Drizzled science image from one exposure reduced using the DASH - pipeline. - Weighted science image : fits - Weighted drizzled science image from one exposure reduced using the - DASH pipeline. - """ - - myDash = DashData(ima_file_name, flt_file_name) - - myDash.split_ima() - - myDash.create_seg_map() - - diffpath = os.path.dirname( - os.path.abspath("diff/{}_*_diff.fits".format(myDash.root)) - ) - cat_images = sorted( - [os.path.basename(x) for x in glob("diff/{}_*_diff.fits".format(myDash.root))] - ) - sc_diff_files = [diffpath + "/" + s for s in cat_images] - myDash.diff_seg_map(cat_images=sc_diff_files) - - myDash.subtract_background_reads() - - myDash.fix_cosmic_rays() - - myDash.align( - align_method=align_method, - ref_catalog=ref_catalog, - drz_output=drz_output, - subtract_background=subtract_background, - wcsname=wcsname, - threshold=threshold, - cw=cw, - updatehdr=updatehdr, - updatewcs=updatewcs, - cat_file=cat_file, - searchrad=searchrad, - astrodriz=astrodriz, - ) diff --git a/notebooks/WFC3/dash/requirements.txt b/notebooks/WFC3/dash/requirements.txt deleted file mode 100644 index 81a6a7d8a..000000000 --- a/notebooks/WFC3/dash/requirements.txt +++ /dev/null @@ -1,12 +0,0 @@ -astropy==5.3.3 -astroquery==0.4.6 -drizzlepac==3.5.1 -git+https://github.com/cmccully/lacosmicx.git -matplotlib==3.7.0 -photutils==1.6.0 -stsci.image==2.3.5 -stsci.imagestats==1.6.3 -stsci.skypac==1.0.9 -stsci.stimage==0.2.6 -stsci.tools==4.0.1 -stwcs==1.7.2 diff --git a/notebooks/WFC3/dash/utils.py b/notebooks/WFC3/dash/utils.py deleted file mode 100644 index f228d7c50..000000000 --- a/notebooks/WFC3/dash/utils.py +++ /dev/null @@ -1,80 +0,0 @@ -#! /usr/bin/env python - -""" Tools to help with reducing DASH/IR data. - -Authors -------- - Rosalia O'Brien 2019 - Catherine Martlin 2018/2019 - Iva Momcheva 2018 - Mario Gennaro 2018 - -""" -import os - -from astropy.io import fits -from urllib.request import urlretrieve - - -def get_flat(file_name): - ''' - Will check if user has proper reference file directories - and files. Will also return flat field file appropriate for - the input file. - - Parameters - ---------- - file_name : string - File name of input IMA. - - Returns - ---------- - reffile_name : string - File name of flat field for that file. - - ''' - os.environ['iref'] = 'iref/' - if not os.path.exists('iref'): - os.mkdir('iref') - - base_url = 'https://hst-crds.stsci.edu/unchecked_get/references/hst/' - - with fits.open(file_name) as fitsfile: - reffile_name = fitsfile[0].header['PFLTFILE'].replace('$', '/') - if not os.path.exists(reffile_name): - urlretrieve(base_url + os.path.basename(reffile_name), reffile_name) - - return reffile_name - - -def get_IDCtable(file_name): - ''' - Will check if user has proper reference file directories - and files. Will also return Instrument Distortion Calibration - reference file appropriate for the input file. - - Parameters - ---------- - file_name : string - File name of input IMA. - - Returns - ---------- - reffile_name : string - File name of Instrument Distortion Calibration reference file - for that file. - - ''' - - os.environ['iref'] = 'iref/' - if not os.path.exists('iref'): - os.mkdir('iref') - - base_url = 'https://hst-crds.stsci.edu/unchecked_get/references/hst/' - - with fits.open(file_name) as fitsfile: - reffile_name = fitsfile[0].header['IDCTAB'].replace('$', '/') - if not os.path.exists(reffile_name): - urlretrieve(base_url + os.path.basename(reffile_name), reffile_name) - - return reffile_name