From 99f9a51010da210445e6af97d70cfde08abbb1f4 Mon Sep 17 00:00:00 2001 From: Ben Sunnquist Date: Wed, 24 Jul 2024 14:50:12 -0400 Subject: [PATCH 1/9] Adding code to remove wisps from NIRCam data --- .../nircam_wisp_subtraction.ipynb | 462 +++++++++++++++ .../NIRCam_wisp_subtraction/requirements.txt | 2 + .../NIRCam_wisp_subtraction/subtract_wisp.py | 559 ++++++++++++++++++ 3 files changed, 1023 insertions(+) create mode 100644 notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb create mode 100644 notebooks/NIRCam/NIRCam_wisp_subtraction/requirements.txt create mode 100644 notebooks/NIRCam/NIRCam_wisp_subtraction/subtract_wisp.py diff --git a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb new file mode 100644 index 000000000..965839714 --- /dev/null +++ b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb @@ -0,0 +1,462 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "00417c2b-3e15-4ebd-bd0a-43307a5e1708", + "metadata": {}, + "source": [ + "# NIRCam Wisp Removal\n", + "---\n", + "**Author**: Ben Sunnquist (bsunnquist@stsci.edu)| **Latest Update**: 24 July 2024" + ] + }, + { + "cell_type": "markdown", + "id": "3f29f4b4-f889-4cbb-89be-4dc5d8c7e5d1", + "metadata": {}, + "source": [ + "\n", + "## Table of Contents\n", + "* [Introduction](#intro)\n", + "* [Imports](#imports)\n", + "* [Data](#data)\n", + "* [Run the wisp subtraction code](#wisp_sub)\n", + "* [Run the wisp subtraction code using custom settings](#wisp_sub_custom)\n", + " * [New wisp scaling settings](#ex1)\n", + " * [Data quality flagging](#ex2)\n", + " * [Other misc settings](#ex3)\n", + "* [Using custom wisp templates](#custom_templates)\n", + "* [Running the code in pieces](#pieces)" + ] + }, + { + "cell_type": "markdown", + "id": "338e0783-ee8e-4291-9574-755faaa73610", + "metadata": {}, + "source": [ + "\n", + "## Introduction" + ] + }, + { + "cell_type": "markdown", + "id": "fa9f29bf-62d3-42c4-984e-4b4b2b16d88d", + "metadata": {}, + "source": [ + "This notebook demonstrates how to remove [wisps](https://jwst-docs.stsci.edu/known-issues-with-jwst-data/nircam-known-issues/nircam-scattered-light-artifacts#NIRCamScatteredLightArtifacts-wispsWisps) from NIRCam imaging data. Wisps are a scattered light feature affecting detectors A3, A4, B3, and B4. For a given filter, wisps appear in the same detector location with only their brightness varying between exposures; therefore, they can be removed from science data by scaling and subtracting a wisp template (i.e. a median combination of all wisp appearances).\n", + "\n", + "Wisp templates used by this notebook can be downloaded from the version3 folder in the [NIRCam wisp template Box folder](https://stsci.box.com/s/1bymvf1lkrqbdn9rnkluzqk30e8o2bne). For this notebook, only the F200W templates are needed.\n", + "\n", + "This notebook uses the `subtract_wisp.py` code to scale and subtract the wisps. That code can be used by itself within python, and is preferred if calibrating a large number of files in parallel, but this notebook will be used to demonstrate the various parameters available to optimize wisp removal. For each notebook cell, we'll also show the corresponding command to run the equivalent in python." + ] + }, + { + "cell_type": "markdown", + "id": "0a824094-c56b-44b3-80e6-9466a6f4f57a", + "metadata": {}, + "source": [ + "\n", + "## Imports" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "fa348563-e007-4135-9481-7fed3e6c09e6", + "metadata": {}, + "outputs": [], + "source": [ + "# How I made this environment:\n", + "# conda create -n jwst1140 python=3.11 notebook\n", + "# source activate jwst1140\n", + "# pip install jwst==1.14.0\n", + "# pip install astroquery\n", + "\n", + "from subtract_wisp import make_segmap, process_file, process_files, subtract_wisp\n", + "\n", + "from astropy.io import fits\n", + "from astroquery.mast import Mast, Observations\n", + "import glob\n", + "import matplotlib\n", + "%matplotlib inline\n" + ] + }, + { + "cell_type": "markdown", + "id": "cb679760-0f63-4377-b99a-73c690000ff7", + "metadata": {}, + "source": [ + "\n", + "## Data" + ] + }, + { + "cell_type": "markdown", + "id": "38e477be-b001-49be-9b8f-7fb503686199", + "metadata": {}, + "source": [ + "We'll use F200W imaging data from program 1063 observation 6 taken on detectors A3, A4, B3, and B4. These are images of sparse fields that are affected by bright wisps. We'll also use the corresponding longwave images for source detection." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "499b9d8d-1f8a-4258-925d-3b60d44c7de7", + "metadata": {}, + "outputs": [], + "source": [ + "# Download the relevant data from MAST\n", + "\n", + "params = {\"columns\": \"dataURI, filename, exp_type\",\n", + " \"filters\": [{\"paramName\": \"program\", \"values\": ['1063']},\n", + " {\"paramName\": \"observtn\", \"values\": ['6']},\n", + " {\"paramName\": \"exposure\", \"values\": ['00005']},\n", + " {\"paramName\": \"visit\", \"values\": ['004']},\n", + " {\"paramName\": \"detector\", \"values\": ['NRCA3','NRCA4','NRCB3','NRCB4','NRCALONG','NRCBLONG']},\n", + " {\"paramName\": \"productLevel\", \"values\": ['2b']}]}\n", + "t = Mast().service_request('Mast.Jwst.Filtered.Nircam', params)\n", + "for row in t:\n", + " if '_cal' in row['filename']: # only want cal files\n", + " result = Observations().download_file(row['dataURI'], cache=False)\n" + ] + }, + { + "cell_type": "markdown", + "id": "4ee4c38e-d86d-4096-833e-9e92a87cb78a", + "metadata": {}, + "source": [ + "\n", + "## Run the wisp subtraction code" + ] + }, + { + "cell_type": "markdown", + "id": "32af21b1-b980-433c-907a-578d34629339", + "metadata": {}, + "source": [ + "In this section, we show how to run the wisp subtraction code on all of the shortwave images using the default parameters. The exception is we set `show_plot=True` to allow the diagnostic plots to display in the notebook.\n", + "\n", + "The basic process the code uses is the following:\n", + "\n", + "1. Generate a source mask using the corresponding longwave image, since the longwave isn't affected by wisps. Ignore these pixels going forward. Occurs in the `make_segmap()` function.\n", + "2. Load the relevant wisp template. Occurs in the `process_file()` function.\n", + "3. Apply a series of scale factors to the template. For each scale factor, multiply it onto the wisp template and subtract the result from the input science data. Record the residual noise within the wisp region for each of the scale factors tested. A correction to residual 1/f (horizontal) noise is applied during this step to help with the fitting. Fit a polynomial to the residuals, and choose the scale factor with the lowest noise based on this fit. Occurs in the `subtract_wisp()` function.\n", + "4. Multiply the chosen wisp factor onto the wisp template to generate the final wisp model. Occurs in the `subtract_wisp()` function.\n", + "5. Subtract the final wisp model from the input data, and write out the wisp-subtracted data, the wisp model, and diagnostic plots summarizing the results. Occurs in the `subtract_wisp()` function." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "c196dbd7-902a-42b1-9f85-5402b5220fcb", + "metadata": {}, + "outputs": [], + "source": [ + "files = glob.glob('*_cal.fits')\n", + "files = [f for f in files if 'long' not in f] # only want shortwave files\n", + "for file in files:\n", + " results = process_file(file, show_plot=True)\n" + ] + }, + { + "cell_type": "markdown", + "id": "9418c448-5a14-41ac-b14e-8fad1011fc49", + "metadata": {}, + "source": [ + "Python equivalent: `python subtract_wisp.py --files *_cal.fits`" + ] + }, + { + "cell_type": "markdown", + "id": "3b16a76b-75cf-4861-81f0-6b249d0533f7", + "metadata": {}, + "source": [ + "As seen above, the code generates the following:\n", + "* A wisp-corrected version of the input file, with the same name as the input file besides a _wisp suffix added (note: the suffix used here can be changed using the `suffix` argument, which can also be set to an empty string to overwrite the input file).\n", + "* The wisp model subtracted from the input image, with the same name as the input file besides a _wisp_model suffix added.\n", + "* A diagnostic plot showing the original image, the wisp-corrected image, the wisp model subtracted, and a plot showing all of the wisp scale factors tested along with their corresponding residual noise. The factor with the lowest noise is the one applied to the wisp template to generate the wisp model.\n", + "\n", + "Note: calling `process_files()` rather than `process_file()` in the cell above will process the files in parallel. This will save time, but the plots won't display inline in the notebook." + ] + }, + { + "cell_type": "markdown", + "id": "eb65ba67-65c1-4623-ae03-a3267f0e2d37", + "metadata": {}, + "source": [ + "\n", + "## Run the wisp subtraction code using custom settings" + ] + }, + { + "cell_type": "markdown", + "id": "3c9cd47f-7d44-4817-b2a4-e754db0e0888", + "metadata": {}, + "source": [ + "In this section, we'll run the wisp subtraction code using custom settings that may be helpful to optimize the wisp removal. The different settings shown here can be used in combination with those from other examples or as a standalone." + ] + }, + { + "cell_type": "markdown", + "id": "996ff57a-cd4e-44e5-8f74-9fe21dd90023", + "metadata": {}, + "source": [ + "\n", + "### Example 1: New wisp scaling method" + ] + }, + { + "cell_type": "markdown", + "id": "752133eb-2271-4e22-b1c1-8a605f8ad807", + "metadata": {}, + "source": [ + "By default, the wisp factor is decided based on that which results in the lowest residual noise. However, the scaling can instead be based on minimizing the overall signal level differences inside and outside the wisp-affected detector regions by setting `scale_method='median'`. For this scaling method, it's recommended to also turn off the polynomial fitting on the residuals (`poly_degree=0`) to simply choose the factor that results in the lowest difference.\n", + "\n", + "During the scaling, residual 1/f noise is corrected by default (`correct_rows=True`); however, for datasets that have strong e.g. amp offsets or odd-even column residuals, this same correction can be applied in the vertical direction as well (`correct_cols=True`). Note that these corrections aren't subtracted in the output file, they are only used to help the wisp scaling procedure.\n", + "\n", + "Below, we perform a new run using these settings." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "67c3d9ae-c62e-4eb7-a550-dc181183c159", + "metadata": {}, + "outputs": [], + "source": [ + "results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits', scale_method='median', poly_degree=0, \n", + " correct_cols=True, show_plot=True)\n" + ] + }, + { + "cell_type": "markdown", + "id": "0a3f8b44-659c-40b0-9ebc-fc52657ecfd8", + "metadata": {}, + "source": [ + "Python equivalent: `python subtract_wisp.py --files jw01063006004_02101_00005_nrcb4_cal.fits --scale_method median --poly_degree 0 --correct_cols`" + ] + }, + { + "cell_type": "markdown", + "id": "e69a261c-00ed-4b78-b5a0-fa1a38ee736d", + "metadata": {}, + "source": [ + "Sometimes, e.g. for crowded fields, there may not be enough empty regions to scale the wisp properly. In this case, simply subtracting the wisp template with no scaling may provide the best results (`scale_wisp=False`, or similarly, in python `--no-scale_wisp`)." + ] + }, + { + "cell_type": "markdown", + "id": "2bbf490a-45f4-4475-b9dd-aa5c89408126", + "metadata": {}, + "source": [ + "\n", + "### Example 2: Data quality flagging" + ] + }, + { + "cell_type": "markdown", + "id": "e1817da3-aa13-471c-aa3f-ff2af64e0f10", + "metadata": {}, + "source": [ + "Another option is to flag any pixels with wisp signals above a certain threshold in the image's data quality array, but not subtract the wisp itself from the science data; this can be accomplished by setting `flag_wisp_thresh` to the chosen wisp signal threshold (units are assumed to be the same as the input file) and `sub_wisp=False`. The data quality value used here is set by the `dq_val` parameter.\n", + "\n", + "The default value of `dq_val` is 1, which is equivalent to the JWST Pipeline DO_NOT_USE flag, and it is ignored by default during image3 drizzling. If you would prefer the pixels to be flagged, but not ignored during that step, another option is using `dq_val=1073741824`, which is equivalent to the JWST Pipeine OTHER_BAD_PIXEL flag. More details on the various data quality flags can be found [here](https://jwst-pipeline.readthedocs.io/en/latest/jwst/references_general/references_general.html#data-quality-flags)." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "e69fcf4f-a289-48f4-88a3-21609e87b770", + "metadata": {}, + "outputs": [], + "source": [ + "results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits', flag_wisp_thresh=0.03, dq_val=1073741824, sub_wisp=False,\n", + " show_plot=True)\n" + ] + }, + { + "cell_type": "markdown", + "id": "98429059-9038-40a3-be3c-ca5edb81b0d8", + "metadata": {}, + "source": [ + "Python equivalent: `python subtract_wisp.py --files jw01063006004_02101_00005_nrcb4_cal.fits --flag_wisp_thresh 0.03 --dq_val 1073741824 --no-sub_wisp`\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "e5736e6e-1e40-4e60-8490-c5740a94d13d", + "metadata": {}, + "outputs": [], + "source": [ + "# Check that the data quality array in this file was updated appropriately\n", + "\n", + "dq = fits.getdata('jw01063006004_02101_00005_nrcb4_cal_wisp.fits', 'DQ')\n", + "dq = (dq&1073741824!=0).astype(int) # only want to see pixels flagged as OTHER_BAD_PIXEL, i.e. the dq_val used above\n", + "plt.imshow(dq, cmap='gray', origin='lower', vmin=0, vmax=0.1)\n" + ] + }, + { + "cell_type": "markdown", + "id": "275c048c-f642-46c0-9f0b-68b9bcda0d53", + "metadata": {}, + "source": [ + "\n", + "### Example 3: Other misc settings" + ] + }, + { + "cell_type": "markdown", + "id": "6fbfb5c2-0b40-4fa4-8890-3fdb1a66e284", + "metadata": {}, + "source": [ + "This example will show a variety of other settings that may be useful either used alone or in combination with any of the other settings shown above. The full list of optional settings can be seen in the `make_segmap()` and `subtract_wisp()` docstrings.\n", + "\n", + "In this example, we'll do the source-finding on the input image itself (rather than the corresponding longwave image) by setting `seg_from_lw=False`. We'll increase the source detection `sigma` to help avoid flagging the wisp itself, and we'll save the the resulting segmentation mask by setting `save_segmap=True` for inspection. The segmentation mask filename is the same as the input file, just with the _seg suffix added.\n", + "\n", + "We'll also only subtract wisp signals above 0.01 MJy/sr (`min_wisp=0.01`). To save a little time, we'll only test wisp scale factors between 0.5 and 1.5 in 0.05 steps." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "8eaf492c-bff4-48e9-9220-e14e56e43e9d", + "metadata": {}, + "outputs": [], + "source": [ + "results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits', seg_from_lw=False, sigma=1.5, save_segmap=True, \n", + " min_wisp=0.01, factor_min=0.5, factor_max=1.5, factor_step=0.05, show_plot=True)\n" + ] + }, + { + "cell_type": "markdown", + "id": "8195caf8-ca91-497a-b904-e54d4a990e8d", + "metadata": {}, + "source": [ + "Python equivalent: `python subtract_wisp.py --files jw01063006004_02101_00005_nrcb4_cal.fits --no-seg_from_lw --sigma 1.5 --save_segmap --min_wisp 0.01 --factor_min 0.5 --factor_max 1.5 --factor_step 0.05`\n" + ] + }, + { + "cell_type": "markdown", + "id": "793c22b0-9ad0-4692-91bb-d96792696ab1", + "metadata": {}, + "source": [ + "As seen above, these custom settings resulted in an oversubtraction. This is likely due to poor source finding in the wisp region (e.g. the wisp itself being flagged as a source), and why source finding using the longwave image is recommended. We can confirm this by inspecting the generated segmentation map." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "acaaaaf0-9285-42ee-a6b8-3bf84f4cc826", + "metadata": {}, + "outputs": [], + "source": [ + "data = fits.getdata('jw01063006004_02101_00005_nrcb4_cal.fits')\n", + "segmap = fits.getdata('jw01063006004_02101_00005_nrcb4_cal_seg.fits')\n", + "\n", + "fig, axes = plt.subplots(1, 2, figsize=(20,10))\n", + "axes[0].imshow(data, origin='lower', cmap='gray', vmin=0.18, vmax=0.3)\n", + "axes[1].imshow(segmap, origin='lower', cmap='gray', vmin=0, vmax=0.1)\n" + ] + }, + { + "cell_type": "markdown", + "id": "b57fa256-94df-4577-a8d7-404e2824a1f4", + "metadata": {}, + "source": [ + "\n", + "## Using custom wisp templates" + ] + }, + { + "cell_type": "markdown", + "id": "a01cb24d-13af-4745-b89b-3ccc7d0e1dc4", + "metadata": {}, + "source": [ + "This section shows how you can incorporate custom wisp templates while still using the wisp subtraction procedure contained in the main code. This may be useful if you have your own wisp templates or want to apply a manual scale factor to the existing templates.\n", + "\n", + "In this example, we'll use the existing template, but multiply it by our own scale factor and turn the scaling procedure off (`scale_wisp=False`) so the code simply applies the wisp model as given. Since the scaling procedure is off, the source finding is irrelevant, and so we can skip the segmap creation function `make_segmap()` and run `subtract_wisp()` directly. " + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "17dbbd90-a406-416e-a70d-2fc45b0a7b38", + "metadata": {}, + "outputs": [], + "source": [ + "# Get the wisp template and apply a manual scale factor\n", + "wisp_data = fits.getdata('WISP_NRCB4_F200W_CLEAR.fits') * 1.05\n", + "\n", + "# Process the file with the custom wisp template\n", + "results = subtract_wisp('jw01063006004_02101_00005_nrcb4_cal.fits', wisp_data, scale_wisp=False, show_plot=True)\n" + ] + }, + { + "cell_type": "markdown", + "id": "5aedb011-4c13-49ce-bd18-c5d60417d0fc", + "metadata": {}, + "source": [ + "\n", + "## Runing the code in pieces" + ] + }, + { + "cell_type": "markdown", + "id": "a712be10-d672-4d67-be18-b90a94393b34", + "metadata": {}, + "source": [ + "In this section, we'll show how to run the full wisp removal process in pieces. This may be useful to provide more flexibility or when finding the optimal settings to use for your data. We'll use a variety of custom settings - see the function docstrings for more details." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "aaf41481-d1f1-4edc-9b8a-9e2d0f06527f", + "metadata": {}, + "outputs": [], + "source": [ + "# The file to process\n", + "f = 'jw01063006004_02101_00005_nrcb4_cal.fits'\n", + "\n", + "# The wisp template to use. We'll use a \"custom\" template again.\n", + "wisp_data = fits.getdata('WISP_NRCB4_F200W_CLEAR.fits') * 1.05\n", + "\n", + "# Generate the segmentation map\n", + "segmap_data = make_segmap(f, sigma=1.0, npixels=8, dilate_segmap=10)\n", + "\n", + "# Scale and subtract the wisp template\n", + "results = subtract_wisp(f, wisp_data=wisp_data, segmap_data=segmap_data, show_plot=True)\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "69682d25-1b97-4b00-801b-36acd502bc81", + "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" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/notebooks/NIRCam/NIRCam_wisp_subtraction/requirements.txt b/notebooks/NIRCam/NIRCam_wisp_subtraction/requirements.txt new file mode 100644 index 000000000..f6365f22e --- /dev/null +++ b/notebooks/NIRCam/NIRCam_wisp_subtraction/requirements.txt @@ -0,0 +1,2 @@ +astroquery==0.4.7 +jwst==1.14.0 \ No newline at end of file diff --git a/notebooks/NIRCam/NIRCam_wisp_subtraction/subtract_wisp.py b/notebooks/NIRCam/NIRCam_wisp_subtraction/subtract_wisp.py new file mode 100644 index 000000000..e1f3c0d54 --- /dev/null +++ b/notebooks/NIRCam/NIRCam_wisp_subtraction/subtract_wisp.py @@ -0,0 +1,559 @@ +#! /usr/bin/env python + +""" +Scales and subtracts wisp templates from full-frame NIRCam images. + +Authors +------- + Ben Sunnquist, 2024 + +Use +--- + To run from the command line using all default arguments (assumes input files + and wisp templates are in the current working directory): + >>> python subtract_wisp.py + + To run on all NRCB4 files with no wisp scaling or segmap creation, and flagging + wisp values above 0.01 in the input images' data quality arrays: + >>> python subtract_wisp.py --files ./*nrcb4_cal.fits --no-scale_wisp --no-create_segmap --flag_wisp_thresh 0.01 + + To run the same example above within python e.g. a Jupyter notebook: + >>> from subtract_wisp import process_files + >>> process_files(files='./*nrcb4_cal.fits', create_segmap=False, scale_wisp=False, flag_wisp_thresh=0.01) + +Notes +----- + Wisp templates designed to work with this code are available at the following STScI Box link: + https://stsci.app.box.com/s/1bymvf1lkrqbdn9rnkluzqk30e8o2bne + + This code was developed using the JWST Pipeline v1.14.0 environment: + https://github.com/spacetelescope/jwst + > conda create -n python=3.11 + > conda activate + > pip install jwst==1.14.0 + Other environments may work fine, but be sure to use Python 3.9+. +""" + +import argparse +from functools import partial +import multiprocessing +import os +import warnings +warnings.filterwarnings('ignore', message="Input data contains invalid values*") # nan values expected throughout code +warnings.filterwarnings('ignore', message="All-NaN slice encountered*") +warnings.filterwarnings('ignore', message="'obsfix' made the change*") # from astropy wcs during lw segmap blotting +warnings.filterwarnings('ignore', message="'datfix' made the change*") # from astropy wcs during lw segmap blotting + +from astropy.convolution import convolve, Gaussian2DKernel +from astropy.io import fits +from astropy.stats import median_absolute_deviation, sigma_clipped_stats +from astropy.wcs import WCS +import matplotlib +matplotlib.use('Agg') +import matplotlib.pyplot as plt +import numpy as np +from photutils.segmentation import detect_sources, detect_threshold +from scipy.ndimage import binary_dilation, generate_binary_structure + +# ----------------------------------------------------------------------------- + +def make_segmap(f, seg_from_lw=True, sigma=0.8, npixels=10, dilate_segmap=5, save_segmap=False): + """ + Make a segmentation map for the input file. + + Parameters + ---------- + f : str + The filename to make a segmentation map for, e.g. + jw01063006004_02101_00005_nrca3_cal.fits. Source-finding is only performed + on cal files, so the cal file must be present if a rate file is input here. + + seg_from_lw : bool + Option to perform the source-finding on the corresponding + longwave image. This is generally preferred since the longwave is not + affected by wisps. + + sigma: float + The source detection sigma threshold. + + npixels: int + The number of connected pixels above sigma to include as a source. + + dilate_segmap: int + The number of pixels to grow the segmap outwards by. + + save_segmap : bool + Option to save the generated segmentation map. + + Returns + ------- + segmap : numpy.ndarray + The generated segmentation map. + + Outputs + ------- + {f}_seg.fits : fits image + The generated segmentation map. + """ + + print('Making segmap for {}'.format(f)) + outfile = f.replace('.fits', '_seg.fits') + + # Get the input data; always source-find on the cal image + detector = os.path.basename(f).split('_')[-2].lower() + f_sw = f.replace('_rate.fits', '_cal.fits') + if (seg_from_lw) & ('long' not in detector): + if 'a' in detector: + f_lw = f.replace(detector.lower(), 'nrcalong').replace('_rate.fits', '_cal.fits') + if 'b' in detector: + f_lw = f.replace(detector.lower(), 'nrcblong').replace('_rate.fits', '_cal.fits') + data = fits.getdata(f_lw, 'SCI') + dq = fits.getdata(f_lw, 'DQ') + else: + data = fits.getdata(f_sw, 'SCI') + dq = fits.getdata(f_sw, 'DQ') + + # Make the segmentation map + threshold = detect_threshold(data, sigma) + g = Gaussian2DKernel(x_stddev=3) + data_conv = convolve(data, g, mask=dq&1!=0) # Smooth input image before detecting sources + seg = detect_sources(data_conv, threshold, npixels=npixels, mask=dq&1!=0) # avoid bad pixels as sources + segmap_data = seg.data + segmap_data[segmap_data!=0] = 1 + + # Dilate the segmap outwards + if dilate_segmap != 0: + segmap_data = binary_dilation(segmap_data, iterations=dilate_segmap, structure=generate_binary_structure(2, 2)) + segmap_data[segmap_data!=0] = 1 + + # Blot LW segmap back onto SW detector space + if (seg_from_lw) & ('long' not in detector): + segmap_tmp = np.zeros(segmap_data.shape).astype(int) + wcs = WCS(fits.getheader(f_lw, 'SCI')) # lw cal wcs + seg_y, seg_x = np.where(segmap_data!=0) + sky_coords = wcs.pixel_to_world(seg_x, seg_y) + wcs = WCS(fits.getheader(f_sw, 'SCI')) # sw cal wcs + coords = wcs.world_to_pixel(sky_coords) + for i in np.arange(len(coords[0])): + y,x = int(coords[1][i]), int(coords[0][i]) + if (y<2048) & (x<2048) & (y>=0) & (x>=0): + segmap_tmp[y,x] = 1 + # Dilate to compensate for y,x rounding due to different lw/sw pixel scales + segmap_data = binary_dilation(segmap_tmp, iterations=1, structure=generate_binary_structure(2, 2)).astype(int) + + # Save the segmap + segmap_data = segmap_data.astype(int) + if save_segmap: + fits.writeto(outfile, segmap_data, overwrite=True) + + return segmap_data + +# ----------------------------------------------------------------------------- + +def process_file(f, wisp_dir='./', create_segmap=True, seg_from_lw=True, sigma=0.8, npixels=10, dilate_segmap=5, + save_segmap=False, sub_wisp=True, gauss_smooth_wisp=False, gauss_stddev=3.0, scale_wisp=True, + scale_method='mad', poly_degree=5, factor_min=0.0, factor_max=2.0, factor_step=0.01, min_wisp=None, + flag_wisp_thresh=None, dq_val=1, correct_rows=True, correct_cols=False, save_data=True, save_model=True, + plot=True, show_plot=False, suffix='_wisp'): + """ + The main processing function. Combines the segmap creation and wisp scaling/subtraction steps together. + + See make_segmap() and subtract_wisp() docstrings for more details on all of the input arguments. + + Parameters + ---------- + f : str + The input filename to subtract the wisp from. Must be a full-frame NIRCam image, e.g. + jw01063006004_02101_00005_nrca3_cal.fits. + + wisp_dir : str + The directory containing the wisp templates. The templates are assumed to have the + form WISP_{DETECTOR}_{FILTER}_{PUPIL}.fits. + """ + + # Get the relevant wisp template + print('Processing {}'.format(f)) + header = fits.getheader(f) + det, fltr, pupil = header['DETECTOR'], header['FILTER'], header['PUPIL'] + file_type = f.split('_')[-1].replace('.fits', '').upper() # CAL or RATE + wisp_data = fits.getdata(os.path.join(wisp_dir, 'WISP_{}_{}_{}.fits'.format(det, fltr, pupil)), file_type) + + # Make the segmentation map + if create_segmap: + segmap_data = make_segmap(f, seg_from_lw=seg_from_lw, sigma=sigma, npixels=npixels, + dilate_segmap=dilate_segmap, save_segmap=save_segmap) + else: + segmap_data = np.zeros(wisp_data.shape).astype(int) + + # Scale and subtract wisp template + results = subtract_wisp(f, wisp_data=wisp_data, segmap_data=segmap_data, sub_wisp=sub_wisp, + gauss_smooth_wisp=gauss_smooth_wisp, gauss_stddev=gauss_stddev, + scale_wisp=scale_wisp, scale_method=scale_method, poly_degree=poly_degree, + factor_min=factor_min, factor_max=factor_max, factor_step=factor_step, + min_wisp=min_wisp, flag_wisp_thresh=flag_wisp_thresh, dq_val=dq_val, + correct_rows=correct_rows, correct_cols=correct_cols, save_data=save_data, + save_model=save_model, plot=plot, show_plot=show_plot, suffix=suffix) + print('Processing complete for {}'.format(f)) + +# ----------------------------------------------------------------------------- + +def process_files(files, nproc=6, **kwargs): + """"Wrapper around the process_file() function to allow for multiprocessing.""" + + # Remove any files that are not in a detector impacted by wisps + files = [f for f in files if any(substring in f for substring in ['nrca3', 'nrca4', 'nrcb3', 'nrcb4'])] + print('Found {} relevant input files.'.format(len(files))) + + # Proess the files + process_file_partial = partial(process_file, **kwargs) + p = multiprocessing.Pool(nproc) + results = p.map(process_file_partial, files) + p.close() + p.join() + +# ----------------------------------------------------------------------------- + +def subtract_wisp(f, wisp_data, segmap_data=None, sub_wisp=True, gauss_smooth_wisp=False, gauss_stddev=3.0, scale_wisp=True, + scale_method='mad', poly_degree=5, factor_min=0.0, factor_max=2.0, factor_step=0.01, min_wisp=None, + flag_wisp_thresh=None, dq_val=1, correct_rows=True, correct_cols=False, save_data=True, save_model=True, plot=True, + show_plot=False, suffix='_wisp'): + """Scales and subtracts a wisp template from the input file. + + Parameters + ---------- + f : str + The file to subtract the wisp from. Must be a full-frame NIRCam image, e.g. + jw01063006004_02101_00005_nrca3_cal.fits. + + wisp_data : numpy.ndarray + The wisp template data. Must be the same shape as input file data. + + segmap_data : numpy.ndarray + A mask of the sources in the input file. Must be the same shape as input file data. + + sub_wisp : bool + Option to subtract the wisp template from the input file. Generally this is always True + except for cases where you only want to flag high wisp values in the input file's data + quality array, i.e. using flag_wisp_thresh. + + gauss_smooth_wisp: bool + Option to smooth the wisp template with a Gaussian filter before scaling/subtracting it. + + gauss_stddev: float + The standard deviation of the Gaussian filter to apply to the wisp if gauss_smooth_wisp is True. + + scale_wisp : bool + Option to scale the wisp templates before subtracting it from the input file. + + scale_method : str + The method used to scale the wisp template. Options are 'mad' and 'median'. The former + chooses the scale factor by minimizing the median absolute deviation of the wisp-corrected + image. The latter chooses the scale factor by minimizing the overall signal level + difference between the wisp-affected region and the remaining detector region. + + poly_degree : int + The degree of the fitting polynomial applied to the residuals of the chosen scale_method. + The minimum of this polynomial fit decides the wisp scale factor. Set to 0 to ignore, + which is generally recommended when scale_method is 'median'. + + factor_min: float + The lowest factor to scale the wisp template by. + + factor_max: float + The highest factor to scale the wisp template by. + + factor_step: float + The step distance between factor_min and factor_max to try scaling the wisp template by. + Scale factors from factor_min to factor_max in iterations equal to the step distance + specified here will all be tested. + + min_wisp : float + The minimum wisp value to perform wisp subtraction. Everything below this value will be set + to zero in the wisp template before subtracting from the input file. The units of this + value should match the units in the input file. Set to None to ignore. + + flag_wisp_thresh : float + Flag values in the wisp template above this threshold in the input file's data quality + array. The units of this value should match the units in the input file. Set to None to ignore. + The DQ value used here is set by the dq_val parameter. + + dq_val : int + The data quality value applied to the input file's DQ array to pixels above flag_wisp_thresh. + The default of 1, i.e. DO_NOT_USE in the JWST Pipeline, results in these pixels being ignored + in the JWST Pipeline image3 drizzling step by default. 1073741824, i.e. OTHER_BAD_PIXEL, is another + option to use if you want these pixels flagged in the input file's DQ array, but still used in the + default drizzling step. Only relevant if flag_wisp_thresh is True. + + correct_rows : bool + Option to correct horizontal noise in the input file (e.g. 1/f noise residuals). This won't be + subtracted from the input file; it is only used to help the wisp scaling step. + + correct_cols : bool + Option to correct vertical noise in the input file (e.g. residual odd-even column effects + and amp offsets). This noise is generally much lower than the horizontal noise in correct_rows. + This won't be subtracted from the input file; it is only used to help the wisp scaling step. + + save_data : bool + Option to save the wisp-subtracted data. The filename will be {f}{suffix}.fits. + + save_model : bool + Option to save the wisp model subtracted fom the input file. The filename will + be {f}{suffix}_model.fits. + + plot : bool + Option to save a diagnostic plot showing the original data, wisp-subtracted data, final wisp model, + and a plot of the various wisp scale factors tested with their corresponding residual noise. + The filename will be {f}{suffix}_plots.png. + + show_plot : bool + Option to show the diagnostic plot in e.g. jupyter notebook (note you'll have to set + %matplotlib inline in the notebook as well for this to work). + + suffix : str + The suffix added to all of the output products. If left blank, the original input file will be + overwritten if save_data is True. + + Returns + ------- + new_data : numpy.ndarray + The wisp-subtracted data. + + wisp_model : numpy.ndarray + The wisp model subtracted from the input file. + + factors : numpy.ndarray + The wisp scale factors tested. + + residuals : numpy.ndarray + The residuals that resulted when applying each of the wisp scale factors tested. The exact + meaning of the residuals depends on the chosen scale_method. + + factor : float + The chosen scale factor applied to the wisp template to create wisp_model. + + Outputs + ------- + {f}{suffix}.fits : fits image + The wisp-subtracted input file. Only generated if save_data is True. + + {f}{suffix}_model.fits : fits image + The wisp model that was subtracted from the input file. Only generated if save_model is True. + + {f}{suffix}_plots.png : png image + A diagnostic plot showing the original data, wisp-subtracted data, final wisp model, and + a plot of the various wisp scale factors tested with their corresponding residuals. + Only generated if plot is True. + """ + + # Get the input data + print('Applying wisp template to {}'.format(f)) + h = fits.open(f) + data = h['SCI'].data + dq = h['DQ'].data + if segmap_data is None: + print('Warning: No segmap data given for {}. Assuming no sources.'.format(f)) + segmap_data = np.zeros(data.shape) + + # Smooth wisp template + if gauss_smooth_wisp: + g = Gaussian2DKernel(x_stddev=gauss_stddev) + wisp_data = convolve(wisp_data, g) + + # Scale the wisp template + if scale_wisp: + # Make a mask of the brightest wisp region + mean, med, stddev = sigma_clipped_stats(wisp_data) + wisp_mask = (wisp_data > med + 2 * stddev).astype(int) + + # Make versions of the original data and wisp model where only good + # pixels in the wisp region are unmasked. + data_masked = np.copy(data) + wisp_data_masked = np.copy(wisp_data) + data_masked[(dq&1!=0) | (segmap_data!=0) | (wisp_mask==0)] = np.nan + wisp_data_masked[(dq&1!=0) | (segmap_data!=0) |(wisp_mask==0)] = np.nan + + # Make a version of the original data where only good pixels outside + # the wisp region are unmasked. + data_masked_ff = np.copy(data) + data_masked_ff[(dq&1!=0) | (segmap_data!=0) | (wisp_mask!=0)] = np.nan + + # Correct median-collapsed row/column offsets, representing the 1/f residuals + # and odd-even column residuals and amp offsets, respectively. + med = np.nanmedian(data_masked_ff) + if correct_rows: + collapsed_rows = np.nanmedian(data_masked_ff - med, axis=1) + else: + collapsed_rows = np.zeros(2048) + if correct_cols: + collapsed_cols = np.nanmedian(data_masked_ff - med, axis=0) + else: + collapsed_cols = np.zeros(2048) + correction_image = np.tile(collapsed_cols, (2048, 1)) + \ + np.swapaxes(np.tile(collapsed_rows, (2048, 1)), 0, 1) + data_masked = data_masked - correction_image + data_masked_ff = data_masked_ff - correction_image + med = np.nanmedian(data_masked_ff) + + # Scale wisp template and record residuals + factors = np.arange(factor_min, factor_max, factor_step) + residuals = [] + for factor in factors: + wisp_model = wisp_data_masked * factor + new_data = data_masked - wisp_model + if scale_method == 'mad': + residual = median_absolute_deviation(new_data, ignore_nan=True) + if scale_method == 'median': + residual = abs(med - np.nanmedian(new_data)) + residuals.append(residual) + residuals = np.array(residuals) + + # Choose the wisp template with the lowest noise + if poly_degree != 0: + fn = np.poly1d(np.polyfit(factors, residuals, poly_degree)) # smooth results + factor = factors[np.argmin(fn(factors))] + else: + factor = factors[np.argmin(residuals)] + wisp_model = wisp_data * factor + else: + # Don't scale the wisp data + factors, residuals = np.ones(1), np.zeros(1) # dummy data + factor = 1.0 + wisp_model = wisp_data * factor + + # Only subtract wisp values above the specified threshold + if min_wisp is not None: + wisp_model[wisp_modelflag_wisp_thresh)] += dq_val + else: + new_dq = dq + + # Subtract the final wisp template + if sub_wisp: + new_data = data - wisp_model + else: + new_data = data + + # Save the wisp-subtracted data and model + if save_data: + h['SCI'].data = new_data.astype('float32') + h['DQ'].data = new_dq + h.writeto(f.replace('.fits', '{}.fits'.format(suffix)), overwrite=True) + if save_model: + fits.writeto(f.replace('.fits', '{}_model.fits'.format(suffix)), wisp_model, overwrite=True) + h.close() + + # Make diagnostic plots + if plot: + fig, axes = plt.subplots(1, 4, figsize=(60,10)) + for ax in axes: + ax.tick_params(axis='both', which='major', labelsize=20) + # Plot original image + mean, med, stddev = sigma_clipped_stats(data) + vmin, vmax = med-3*stddev, med+3*stddev + im = axes[0].imshow(data, origin='lower', vmin=vmin, vmax=vmax, cmap='gray') + cbar = fig.colorbar(im, ax=axes[0]) + cbar.ax.tick_params(labelsize=20) + axes[0].set_title('Original', fontsize=35) + # Plot wisp-corrected image + im = axes[1].imshow(new_data, origin='lower', vmin=vmin, vmax=vmax, cmap='gray') + cbar = fig.colorbar(im, ax=axes[1]) + cbar.ax.tick_params(labelsize=20) + axes[1].set_title('Wisp Corrected', fontsize=35) + # Plot wisp model + mean, med, stddev = sigma_clipped_stats(wisp_model) + im = axes[2].imshow(wisp_model, origin='lower', vmin=0, vmax=med+5*stddev, cmap='gray') + cbar = fig.colorbar(im, ax=axes[2]) + cbar.ax.tick_params(labelsize=20) + axes[2].set_title('Wisp Model', fontsize=35) + # Plot scatterplot of wisp scale factors tested vs residuals + if (scale_wisp) & (poly_degree != 0): + axes[3].plot(factors, fn(factors), color='black') + axes[3].scatter(factors, residuals) + axes[3].axvline(factor, color='black') + y_label = 'Residual Median Absolute Deviation' if scale_method == 'mad' else 'Residual Median Difference' + axes[3].set_ylabel(y_label, fontsize=25) + axes[3].set_xlabel('Wisp Factor', fontsize=25) + axes[3].set_title('Chosen Wisp Factor: {:.3f}'.format(factor), fontsize=35) + axes[3].grid(ls='--', alpha=0.5) + fig.savefig(f.replace('.fits', '{}_plots.png'.format(suffix)), dpi=100, bbox_inches='tight') + if not show_plot: + plt.close() + + return new_data, wisp_model, factors, residuals, factor + +# ----------------------------------------------------------------------------- +# ----------------------------------------------------------------------------- + +def parse_args(): + """ + Parses command line arguments. + + See make_segmap() and subtract_wisp() docstrings for more details on all of the input arguments. + + Returns + ------- + args : object + Contains the input arguments. + """ + + # Make the help strings + files_help = 'The files to subtract the wisp templates from. Wildcards are supported, e.g. ./data/*nrca3*_cal.fits' + nproc_help = 'The number of processes to use during multiprocessing.' + wisp_dir_help = 'The directory containing the wisp templates. The templates are assumed to have the form WISP_{DETECTOR}_{FILTER}_{PUPIL}.fits.' + create_segmap_help = 'Option to make a source segmentation map to help with scaling the wisp template.' + + # Add the potential arguments + parser = argparse.ArgumentParser() + parser.add_argument('--files', dest='files', action='store', nargs='+', type=str, required=False, help=files_help, default='./*_cal.fits') + parser.add_argument('--nproc', dest='nproc', action='store', type=int, required=False, help=nproc_help, default=6) + parser.add_argument('--wisp_dir', dest='wisp_dir', action='store', type=str, required=False, help=wisp_dir_help, default='./') + parser.add_argument('--create_segmap', dest='create_segmap', action=argparse.BooleanOptionalAction, required=False, help=create_segmap_help, default=True) + + # Add arguments for make_segmap() + parser.add_argument('--seg_from_lw', dest='seg_from_lw', action=argparse.BooleanOptionalAction, required=False, default=True) + parser.add_argument('--sigma', dest='sigma', action='store', type=float, required=False, default=0.8) + parser.add_argument('--npixels', dest='npixels', action='store', type=int, required=False, default=10) + parser.add_argument('--dilate_segmap', dest='dilate_segmap', action='store', type=int, required=False, default=5) + parser.add_argument('--save_segmap', dest='save_segmap', action=argparse.BooleanOptionalAction, required=False, default=False) + + # Add arguments for subtract_wisp() + parser.add_argument('--sub_wisp', dest='sub_wisp', action=argparse.BooleanOptionalAction, required=False, default=True) + parser.add_argument('--gauss_smooth_wisp', dest='gauss_smooth_wisp', action=argparse.BooleanOptionalAction, required=False, default=False) + parser.add_argument('--gauss_stddev', dest='gauss_stddev', action='store', type=float, required=False, default=3.0) + parser.add_argument('--scale_wisp', dest='scale_wisp', action=argparse.BooleanOptionalAction, required=False, default=True) + parser.add_argument('--scale_method', dest='scale_method', action='store', type=str, required=False, default='mad') + parser.add_argument('--poly_degree', dest='poly_degree', action='store', type=int, required=False, default=5) + parser.add_argument('--factor_min', dest='factor_min', action='store', type=float, required=False, default=0.0) + parser.add_argument('--factor_max', dest='factor_max', action='store', type=float, required=False, default=2.0) + parser.add_argument('--factor_step', dest='factor_step', action='store', type=float, required=False, default=0.01) + parser.add_argument('--min_wisp', dest='min_wisp', action='store', type=float, required=False, default=None) + parser.add_argument('--flag_wisp_thresh', dest='flag_wisp_thresh', action='store', type=float, required=False, default=None) + parser.add_argument('--dq_val', dest='dq_val', action='store', type=int, required=False, default=1) + parser.add_argument('--correct_rows', dest='correct_rows', action=argparse.BooleanOptionalAction, required=False, default=True) + parser.add_argument('--correct_cols', dest='correct_cols', action=argparse.BooleanOptionalAction, required=False, default=False) + parser.add_argument('--save_data', dest='save_data', action=argparse.BooleanOptionalAction, required=False, default=True) + parser.add_argument('--save_model', dest='save_model', action=argparse.BooleanOptionalAction, required=False, default=True) + parser.add_argument('--plot', dest='plot', action=argparse.BooleanOptionalAction, required=False, default=True) + parser.add_argument('--show_plot', dest='show_plot', action=argparse.BooleanOptionalAction, required=False, default=False) + parser.add_argument('--suffix', dest='suffix', action='store', type=str, required=False, default='_wisp') + + # Get the arguments + args = parser.parse_args() + + return args + +# ----------------------------------------------------------------------------- +# ----------------------------------------------------------------------------- + +if __name__ == '__main__': + + # Get the command line arguments + args = parse_args() + + # Process the input files + results = process_files(**vars(args)) + print('subtract_wisp.py complete.') From 301c0a8b753a797229ec7d4c73aa562bc4465215 Mon Sep 17 00:00:00 2001 From: Ben Sunnquist Date: Thu, 25 Jul 2024 11:49:19 -0400 Subject: [PATCH 2/9] pep8 updates --- .../nircam_wisp_subtraction.ipynb | 32 ++++----- .../NIRCam_wisp_subtraction/subtract_wisp.py | 67 +++++++++++-------- 2 files changed, 55 insertions(+), 44 deletions(-) diff --git a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb index 965839714..d750c5579 100644 --- a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb +++ b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb @@ -45,7 +45,7 @@ "source": [ "This notebook demonstrates how to remove [wisps](https://jwst-docs.stsci.edu/known-issues-with-jwst-data/nircam-known-issues/nircam-scattered-light-artifacts#NIRCamScatteredLightArtifacts-wispsWisps) from NIRCam imaging data. Wisps are a scattered light feature affecting detectors A3, A4, B3, and B4. For a given filter, wisps appear in the same detector location with only their brightness varying between exposures; therefore, they can be removed from science data by scaling and subtracting a wisp template (i.e. a median combination of all wisp appearances).\n", "\n", - "Wisp templates used by this notebook can be downloaded from the version3 folder in the [NIRCam wisp template Box folder](https://stsci.box.com/s/1bymvf1lkrqbdn9rnkluzqk30e8o2bne). For this notebook, only the F200W templates are needed.\n", + "Wisp templates used by this notebook are available in the version3 folder in the [NIRCam wisp template Box folder](https://stsci.box.com/s/1bymvf1lkrqbdn9rnkluzqk30e8o2bne). For this notebook, only the F200W templates are needed (i.e. WISP_NRCA3_F200W_CLEAR.fits, WISP_NRCA4_F200W_CLEAR.fits, WISP_NRCB3_F200W_CLEAR.fits and WISP_NRCB4_F200W_CLEAR.fits).\n", "\n", "This notebook uses the `subtract_wisp.py` code to scale and subtract the wisps. That code can be used by itself within python, and is preferred if calibrating a large number of files in parallel, but this notebook will be used to demonstrate the various parameters available to optimize wisp removal. For each notebook cell, we'll also show the corresponding command to run the equivalent in python." ] @@ -72,13 +72,13 @@ "# pip install jwst==1.14.0\n", "# pip install astroquery\n", "\n", - "from subtract_wisp import make_segmap, process_file, process_files, subtract_wisp\n", + "from subtract_wisp import make_segmap, process_file, subtract_wisp\n", "\n", "from astropy.io import fits\n", "from astroquery.mast import Mast, Observations\n", "import glob\n", - "import matplotlib\n", - "%matplotlib inline\n" + "import matplotlib.pyplot as plt\n", + "%matplotlib inline" ] }, { @@ -112,12 +112,12 @@ " {\"paramName\": \"observtn\", \"values\": ['6']},\n", " {\"paramName\": \"exposure\", \"values\": ['00005']},\n", " {\"paramName\": \"visit\", \"values\": ['004']},\n", - " {\"paramName\": \"detector\", \"values\": ['NRCA3','NRCA4','NRCB3','NRCB4','NRCALONG','NRCBLONG']},\n", + " {\"paramName\": \"detector\", \"values\": ['NRCA3', 'NRCA4', 'NRCB3', 'NRCB4', 'NRCALONG', 'NRCBLONG']},\n", " {\"paramName\": \"productLevel\", \"values\": ['2b']}]}\n", "t = Mast().service_request('Mast.Jwst.Filtered.Nircam', params)\n", "for row in t:\n", " if '_cal' in row['filename']: # only want cal files\n", - " result = Observations().download_file(row['dataURI'], cache=False)\n" + " result = Observations().download_file(row['dataURI'], cache=False)" ] }, { @@ -155,7 +155,7 @@ "files = glob.glob('*_cal.fits')\n", "files = [f for f in files if 'long' not in f] # only want shortwave files\n", "for file in files:\n", - " results = process_file(file, show_plot=True)\n" + " results = process_file(file, show_plot=True)" ] }, { @@ -225,7 +225,7 @@ "outputs": [], "source": [ "results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits', scale_method='median', poly_degree=0, \n", - " correct_cols=True, show_plot=True)\n" + " correct_cols=True, show_plot=True)" ] }, { @@ -271,7 +271,7 @@ "outputs": [], "source": [ "results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits', flag_wisp_thresh=0.03, dq_val=1073741824, sub_wisp=False,\n", - " show_plot=True)\n" + " show_plot=True)" ] }, { @@ -292,8 +292,8 @@ "# Check that the data quality array in this file was updated appropriately\n", "\n", "dq = fits.getdata('jw01063006004_02101_00005_nrcb4_cal_wisp.fits', 'DQ')\n", - "dq = (dq&1073741824!=0).astype(int) # only want to see pixels flagged as OTHER_BAD_PIXEL, i.e. the dq_val used above\n", - "plt.imshow(dq, cmap='gray', origin='lower', vmin=0, vmax=0.1)\n" + "dq = (dq & 1073741824 != 0).astype(int) # only want to see pixels flagged as OTHER_BAD_PIXEL, i.e. the dq_val used above\n", + "plt.imshow(dq, cmap='gray', origin='lower', vmin=0, vmax=0.1)" ] }, { @@ -325,7 +325,7 @@ "outputs": [], "source": [ "results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits', seg_from_lw=False, sigma=1.5, save_segmap=True, \n", - " min_wisp=0.01, factor_min=0.5, factor_max=1.5, factor_step=0.05, show_plot=True)\n" + " min_wisp=0.01, factor_min=0.5, factor_max=1.5, factor_step=0.05, show_plot=True)" ] }, { @@ -354,9 +354,9 @@ "data = fits.getdata('jw01063006004_02101_00005_nrcb4_cal.fits')\n", "segmap = fits.getdata('jw01063006004_02101_00005_nrcb4_cal_seg.fits')\n", "\n", - "fig, axes = plt.subplots(1, 2, figsize=(20,10))\n", + "fig, axes = plt.subplots(1, 2, figsize=(20, 10))\n", "axes[0].imshow(data, origin='lower', cmap='gray', vmin=0.18, vmax=0.3)\n", - "axes[1].imshow(segmap, origin='lower', cmap='gray', vmin=0, vmax=0.1)\n" + "axes[1].imshow(segmap, origin='lower', cmap='gray', vmin=0, vmax=0.1)" ] }, { @@ -389,7 +389,7 @@ "wisp_data = fits.getdata('WISP_NRCB4_F200W_CLEAR.fits') * 1.05\n", "\n", "# Process the file with the custom wisp template\n", - "results = subtract_wisp('jw01063006004_02101_00005_nrcb4_cal.fits', wisp_data, scale_wisp=False, show_plot=True)\n" + "results = subtract_wisp('jw01063006004_02101_00005_nrcb4_cal.fits', wisp_data, scale_wisp=False, show_plot=True)" ] }, { @@ -426,7 +426,7 @@ "segmap_data = make_segmap(f, sigma=1.0, npixels=8, dilate_segmap=10)\n", "\n", "# Scale and subtract the wisp template\n", - "results = subtract_wisp(f, wisp_data=wisp_data, segmap_data=segmap_data, show_plot=True)\n" + "results = subtract_wisp(f, wisp_data=wisp_data, segmap_data=segmap_data, show_plot=True)" ] }, { diff --git a/notebooks/NIRCam/NIRCam_wisp_subtraction/subtract_wisp.py b/notebooks/NIRCam/NIRCam_wisp_subtraction/subtract_wisp.py index e1f3c0d54..7d93839d0 100644 --- a/notebooks/NIRCam/NIRCam_wisp_subtraction/subtract_wisp.py +++ b/notebooks/NIRCam/NIRCam_wisp_subtraction/subtract_wisp.py @@ -39,10 +39,6 @@ import multiprocessing import os import warnings -warnings.filterwarnings('ignore', message="Input data contains invalid values*") # nan values expected throughout code -warnings.filterwarnings('ignore', message="All-NaN slice encountered*") -warnings.filterwarnings('ignore', message="'obsfix' made the change*") # from astropy wcs during lw segmap blotting -warnings.filterwarnings('ignore', message="'datfix' made the change*") # from astropy wcs during lw segmap blotting from astropy.convolution import convolve, Gaussian2DKernel from astropy.io import fits @@ -55,8 +51,15 @@ from photutils.segmentation import detect_sources, detect_threshold from scipy.ndimage import binary_dilation, generate_binary_structure +warnings.filterwarnings('ignore', message="Input data contains invalid values*") # nan values expected throughout code +warnings.filterwarnings('ignore', message="All-NaN slice encountered*") +warnings.filterwarnings('ignore', message="'obsfix' made the change*") # from astropy wcs during lw segmap blotting +warnings.filterwarnings('ignore', message="'datfix' made the change*") # from astropy wcs during lw segmap blotting + + # ----------------------------------------------------------------------------- + def make_segmap(f, seg_from_lw=True, sigma=0.8, npixels=10, dilate_segmap=5, save_segmap=False): """ Make a segmentation map for the input file. @@ -116,15 +119,15 @@ def make_segmap(f, seg_from_lw=True, sigma=0.8, npixels=10, dilate_segmap=5, sav # Make the segmentation map threshold = detect_threshold(data, sigma) g = Gaussian2DKernel(x_stddev=3) - data_conv = convolve(data, g, mask=dq&1!=0) # Smooth input image before detecting sources - seg = detect_sources(data_conv, threshold, npixels=npixels, mask=dq&1!=0) # avoid bad pixels as sources + data_conv = convolve(data, g, mask=dq & 1 != 0) # Smooth input image before detecting sources + seg = detect_sources(data_conv, threshold, npixels=npixels, mask=dq & 1 != 0) # avoid bad pixels as sources segmap_data = seg.data - segmap_data[segmap_data!=0] = 1 + segmap_data[segmap_data != 0] = 1 # Dilate the segmap outwards if dilate_segmap != 0: segmap_data = binary_dilation(segmap_data, iterations=dilate_segmap, structure=generate_binary_structure(2, 2)) - segmap_data[segmap_data!=0] = 1 + segmap_data[segmap_data != 0] = 1 # Blot LW segmap back onto SW detector space if (seg_from_lw) & ('long' not in detector): @@ -135,9 +138,9 @@ def make_segmap(f, seg_from_lw=True, sigma=0.8, npixels=10, dilate_segmap=5, sav wcs = WCS(fits.getheader(f_sw, 'SCI')) # sw cal wcs coords = wcs.world_to_pixel(sky_coords) for i in np.arange(len(coords[0])): - y,x = int(coords[1][i]), int(coords[0][i]) - if (y<2048) & (x<2048) & (y>=0) & (x>=0): - segmap_tmp[y,x] = 1 + y, x = int(coords[1][i]), int(coords[0][i]) + if (y < 2048) & (x < 2048) & (y >= 0) & (x >= 0): + segmap_tmp[y, x] = 1 # Dilate to compensate for y,x rounding due to different lw/sw pixel scales segmap_data = binary_dilation(segmap_tmp, iterations=1, structure=generate_binary_structure(2, 2)).astype(int) @@ -148,8 +151,10 @@ def make_segmap(f, seg_from_lw=True, sigma=0.8, npixels=10, dilate_segmap=5, sav return segmap_data + # ----------------------------------------------------------------------------- + def process_file(f, wisp_dir='./', create_segmap=True, seg_from_lw=True, sigma=0.8, npixels=10, dilate_segmap=5, save_segmap=False, sub_wisp=True, gauss_smooth_wisp=False, gauss_stddev=3.0, scale_wisp=True, scale_method='mad', poly_degree=5, factor_min=0.0, factor_max=2.0, factor_step=0.01, min_wisp=None, @@ -186,17 +191,18 @@ def process_file(f, wisp_dir='./', create_segmap=True, seg_from_lw=True, sigma=0 segmap_data = np.zeros(wisp_data.shape).astype(int) # Scale and subtract wisp template - results = subtract_wisp(f, wisp_data=wisp_data, segmap_data=segmap_data, sub_wisp=sub_wisp, - gauss_smooth_wisp=gauss_smooth_wisp, gauss_stddev=gauss_stddev, - scale_wisp=scale_wisp, scale_method=scale_method, poly_degree=poly_degree, - factor_min=factor_min, factor_max=factor_max, factor_step=factor_step, - min_wisp=min_wisp, flag_wisp_thresh=flag_wisp_thresh, dq_val=dq_val, - correct_rows=correct_rows, correct_cols=correct_cols, save_data=save_data, - save_model=save_model, plot=plot, show_plot=show_plot, suffix=suffix) + _ = subtract_wisp(f, wisp_data=wisp_data, segmap_data=segmap_data, sub_wisp=sub_wisp, + gauss_smooth_wisp=gauss_smooth_wisp, gauss_stddev=gauss_stddev, + scale_wisp=scale_wisp, scale_method=scale_method, poly_degree=poly_degree, + factor_min=factor_min, factor_max=factor_max, factor_step=factor_step, + min_wisp=min_wisp, flag_wisp_thresh=flag_wisp_thresh, dq_val=dq_val, + correct_rows=correct_rows, correct_cols=correct_cols, save_data=save_data, + save_model=save_model, plot=plot, show_plot=show_plot, suffix=suffix) print('Processing complete for {}'.format(f)) # ----------------------------------------------------------------------------- + def process_files(files, nproc=6, **kwargs): """"Wrapper around the process_file() function to allow for multiprocessing.""" @@ -207,12 +213,14 @@ def process_files(files, nproc=6, **kwargs): # Proess the files process_file_partial = partial(process_file, **kwargs) p = multiprocessing.Pool(nproc) - results = p.map(process_file_partial, files) + _ = p.map(process_file_partial, files) p.close() p.join() + # ----------------------------------------------------------------------------- + def subtract_wisp(f, wisp_data, segmap_data=None, sub_wisp=True, gauss_smooth_wisp=False, gauss_stddev=3.0, scale_wisp=True, scale_method='mad', poly_degree=5, factor_min=0.0, factor_max=2.0, factor_step=0.01, min_wisp=None, flag_wisp_thresh=None, dq_val=1, correct_rows=True, correct_cols=False, save_data=True, save_model=True, plot=True, @@ -369,13 +377,13 @@ def subtract_wisp(f, wisp_data, segmap_data=None, sub_wisp=True, gauss_smooth_wi # pixels in the wisp region are unmasked. data_masked = np.copy(data) wisp_data_masked = np.copy(wisp_data) - data_masked[(dq&1!=0) | (segmap_data!=0) | (wisp_mask==0)] = np.nan - wisp_data_masked[(dq&1!=0) | (segmap_data!=0) |(wisp_mask==0)] = np.nan + data_masked[(dq & 1 != 0) | (segmap_data != 0) | (wisp_mask == 0)] = np.nan + wisp_data_masked[(dq & 1 != 0) | (segmap_data != 0) |(wisp_mask == 0)] = np.nan # Make a version of the original data where only good pixels outside # the wisp region are unmasked. data_masked_ff = np.copy(data) - data_masked_ff[(dq&1!=0) | (segmap_data!=0) | (wisp_mask!=0)] = np.nan + data_masked_ff[(dq & 1 != 0) | (segmap_data != 0) | (wisp_mask != 0)] = np.nan # Correct median-collapsed row/column offsets, representing the 1/f residuals # and odd-even column residuals and amp offsets, respectively. @@ -388,8 +396,7 @@ def subtract_wisp(f, wisp_data, segmap_data=None, sub_wisp=True, gauss_smooth_wi collapsed_cols = np.nanmedian(data_masked_ff - med, axis=0) else: collapsed_cols = np.zeros(2048) - correction_image = np.tile(collapsed_cols, (2048, 1)) + \ - np.swapaxes(np.tile(collapsed_rows, (2048, 1)), 0, 1) + correction_image = np.tile(collapsed_cols, (2048, 1)) + np.swapaxes(np.tile(collapsed_rows, (2048, 1)), 0, 1) data_masked = data_masked - correction_image data_masked_ff = data_masked_ff - correction_image med = np.nanmedian(data_masked_ff) @@ -422,12 +429,12 @@ def subtract_wisp(f, wisp_data, segmap_data=None, sub_wisp=True, gauss_smooth_wi # Only subtract wisp values above the specified threshold if min_wisp is not None: - wisp_model[wisp_modelflag_wisp_thresh)] += dq_val + new_dq[(dq & dq_val == 0) & (wisp_model > flag_wisp_thresh)] += dq_val else: new_dq = dq @@ -448,7 +455,7 @@ def subtract_wisp(f, wisp_data, segmap_data=None, sub_wisp=True, gauss_smooth_wi # Make diagnostic plots if plot: - fig, axes = plt.subplots(1, 4, figsize=(60,10)) + fig, axes = plt.subplots(1, 4, figsize=(60, 10)) for ax in axes: ax.tick_params(axis='both', which='major', labelsize=20) # Plot original image @@ -485,9 +492,11 @@ def subtract_wisp(f, wisp_data, segmap_data=None, sub_wisp=True, gauss_smooth_wi return new_data, wisp_model, factors, residuals, factor + # ----------------------------------------------------------------------------- # ----------------------------------------------------------------------------- + def parse_args(): """ Parses command line arguments. @@ -546,11 +555,13 @@ def parse_args(): return args + # ----------------------------------------------------------------------------- # ----------------------------------------------------------------------------- + if __name__ == '__main__': - + # Get the command line arguments args = parse_args() From 111794d1b21912af86ccd0bc293241db600ad63d Mon Sep 17 00:00:00 2001 From: Ben Sunnquist Date: Thu, 25 Jul 2024 11:53:53 -0400 Subject: [PATCH 3/9] pep8 updates --- notebooks/NIRCam/NIRCam_wisp_subtraction/subtract_wisp.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/notebooks/NIRCam/NIRCam_wisp_subtraction/subtract_wisp.py b/notebooks/NIRCam/NIRCam_wisp_subtraction/subtract_wisp.py index 7d93839d0..a5840d31b 100644 --- a/notebooks/NIRCam/NIRCam_wisp_subtraction/subtract_wisp.py +++ b/notebooks/NIRCam/NIRCam_wisp_subtraction/subtract_wisp.py @@ -45,11 +45,11 @@ from astropy.stats import median_absolute_deviation, sigma_clipped_stats from astropy.wcs import WCS import matplotlib -matplotlib.use('Agg') import matplotlib.pyplot as plt import numpy as np from photutils.segmentation import detect_sources, detect_threshold from scipy.ndimage import binary_dilation, generate_binary_structure +matplotlib.use('Agg') warnings.filterwarnings('ignore', message="Input data contains invalid values*") # nan values expected throughout code warnings.filterwarnings('ignore', message="All-NaN slice encountered*") @@ -133,7 +133,7 @@ def make_segmap(f, seg_from_lw=True, sigma=0.8, npixels=10, dilate_segmap=5, sav if (seg_from_lw) & ('long' not in detector): segmap_tmp = np.zeros(segmap_data.shape).astype(int) wcs = WCS(fits.getheader(f_lw, 'SCI')) # lw cal wcs - seg_y, seg_x = np.where(segmap_data!=0) + seg_y, seg_x = np.where(segmap_data != 0) sky_coords = wcs.pixel_to_world(seg_x, seg_y) wcs = WCS(fits.getheader(f_sw, 'SCI')) # sw cal wcs coords = wcs.world_to_pixel(sky_coords) @@ -378,7 +378,7 @@ def subtract_wisp(f, wisp_data, segmap_data=None, sub_wisp=True, gauss_smooth_wi data_masked = np.copy(data) wisp_data_masked = np.copy(wisp_data) data_masked[(dq & 1 != 0) | (segmap_data != 0) | (wisp_mask == 0)] = np.nan - wisp_data_masked[(dq & 1 != 0) | (segmap_data != 0) |(wisp_mask == 0)] = np.nan + wisp_data_masked[(dq & 1 != 0) | (segmap_data != 0) | (wisp_mask == 0)] = np.nan # Make a version of the original data where only good pixels outside # the wisp region are unmasked. From e821c8ace411fed4da8e9fcf4a48bf762f978ce1 Mon Sep 17 00:00:00 2001 From: Ben Sunnquist Date: Tue, 13 Aug 2024 10:36:04 -0400 Subject: [PATCH 4/9] Added clarification when using custom templates --- .../NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb index d750c5579..254b7c26e 100644 --- a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb +++ b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb @@ -375,7 +375,7 @@ "source": [ "This section shows how you can incorporate custom wisp templates while still using the wisp subtraction procedure contained in the main code. This may be useful if you have your own wisp templates or want to apply a manual scale factor to the existing templates.\n", "\n", - "In this example, we'll use the existing template, but multiply it by our own scale factor and turn the scaling procedure off (`scale_wisp=False`) so the code simply applies the wisp model as given. Since the scaling procedure is off, the source finding is irrelevant, and so we can skip the segmap creation function `make_segmap()` and run `subtract_wisp()` directly. " + "In this example, we'll use the existing template, but multiply it by our own scale factor and turn the scaling procedure off (`scale_wisp=False`) so the code simply applies the wisp model as given. Note that turning off the scaling procedure is not required when using a custom template; users can both provide a custom template and still keep the scaling procedure on. Since the scaling procedure is off in this example, the source finding is irrelevant, and so we can skip the segmap creation function `make_segmap()` and run `subtract_wisp()` directly." ] }, { From 6883f6dea033f89c9399a8183a61102e75fd96b4 Mon Sep 17 00:00:00 2001 From: Ben Sunnquist Date: Fri, 13 Sep 2024 13:40:58 -0400 Subject: [PATCH 5/9] downloading wisp templates from jdat box folder --- .../nircam_wisp_subtraction.ipynb | 29 +++++++++++++++++-- 1 file changed, 27 insertions(+), 2 deletions(-) diff --git a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb index 254b7c26e..d7be2f141 100644 --- a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb +++ b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb @@ -45,7 +45,7 @@ "source": [ "This notebook demonstrates how to remove [wisps](https://jwst-docs.stsci.edu/known-issues-with-jwst-data/nircam-known-issues/nircam-scattered-light-artifacts#NIRCamScatteredLightArtifacts-wispsWisps) from NIRCam imaging data. Wisps are a scattered light feature affecting detectors A3, A4, B3, and B4. For a given filter, wisps appear in the same detector location with only their brightness varying between exposures; therefore, they can be removed from science data by scaling and subtracting a wisp template (i.e. a median combination of all wisp appearances).\n", "\n", - "Wisp templates used by this notebook are available in the version3 folder in the [NIRCam wisp template Box folder](https://stsci.box.com/s/1bymvf1lkrqbdn9rnkluzqk30e8o2bne). For this notebook, only the F200W templates are needed (i.e. WISP_NRCA3_F200W_CLEAR.fits, WISP_NRCA4_F200W_CLEAR.fits, WISP_NRCB3_F200W_CLEAR.fits and WISP_NRCB4_F200W_CLEAR.fits).\n", + "The wisp templates used in this notebook will be downloaded in the [Data](#data) section, but they are also available in the \"version3\" folder within the [NIRCam wisp template Box folder](https://stsci.box.com/s/1bymvf1lkrqbdn9rnkluzqk30e8o2bne). Future updates to these templates will be added to this same Box folder, and users are encouraged to check this folder regularly to ensure they have the latest version.\n", "\n", "This notebook uses the `subtract_wisp.py` code to scale and subtract the wisps. That code can be used by itself within python, and is preferred if calibrating a large number of files in parallel, but this notebook will be used to demonstrate the various parameters available to optimize wisp removal. For each notebook cell, we'll also show the corresponding command to run the equivalent in python." ] @@ -78,7 +78,8 @@ "from astroquery.mast import Mast, Observations\n", "import glob\n", "import matplotlib.pyplot as plt\n", - "%matplotlib inline" + "%matplotlib inline\n", + "import urllib" ] }, { @@ -120,6 +121,30 @@ " result = Observations().download_file(row['dataURI'], cache=False)" ] }, + { + "cell_type": "markdown", + "id": "2faab36d-49a8-4877-b4eb-4693ce665de1", + "metadata": {}, + "source": [ + "Next, we'll download the wisp templates. These templates are equivalent to those in the \"version3\" folder within the [NIRCam wisp template Box folder](https://stsci.box.com/s/1bymvf1lkrqbdn9rnkluzqk30e8o2bne). Users are encouraged to check this folder regularly for any updates. For this notebook, only the F200W templates are needed (i.e. WISP_NRCA3_F200W_CLEAR.fits, WISP_NRCA4_F200W_CLEAR.fits, WISP_NRCB3_F200W_CLEAR.fits and WISP_NRCB4_F200W_CLEAR.fits)." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "641cd018-8760-47fb-beb2-a4f0d9a11061", + "metadata": {}, + "outputs": [], + "source": [ + "# Download the wisp templates\n", + "\n", + "boxlink = 'https://data.science.stsci.edu/redirect/JWST/jwst-data_analysis_tools/nircam_wisp_templates/'\n", + "\n", + "for detector in ['NRCA3', 'NRCA4', 'NRCB3', 'NRCB4']:\n", + " boxfile = os.path.join(boxlink, 'WISP_{}_F200W_CLEAR.fits'.format(detector))\n", + " urllib.request.urlretrieve(boxfile, os.path.basename(boxfile))" + ] + }, { "cell_type": "markdown", "id": "4ee4c38e-d86d-4096-833e-9e92a87cb78a", From f42bd073947529ee8db1b3d1ef7a8f880475ba11 Mon Sep 17 00:00:00 2001 From: Ben Sunnquist Date: Fri, 13 Sep 2024 13:45:43 -0400 Subject: [PATCH 6/9] forgot to import os --- .../NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb | 1 + 1 file changed, 1 insertion(+) diff --git a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb index d7be2f141..2ff5b4b50 100644 --- a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb +++ b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb @@ -79,6 +79,7 @@ "import glob\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline\n", + "import os\n", "import urllib" ] }, From a572aaaf38ca0190c73041bbd5d8864ce1ad8065 Mon Sep 17 00:00:00 2001 From: Ben Sunnquist Date: Fri, 13 Sep 2024 15:11:14 -0400 Subject: [PATCH 7/9] implemented Rossy technical review --- .../nircam_wisp_subtraction.ipynb | 256 +++++++++++------- 1 file changed, 165 insertions(+), 91 deletions(-) diff --git a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb index 2ff5b4b50..5416b5365 100644 --- a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb +++ b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb @@ -7,47 +7,51 @@ "source": [ "# NIRCam Wisp Removal\n", "---\n", - "**Author**: Ben Sunnquist (bsunnquist@stsci.edu)| **Latest Update**: 24 July 2024" - ] - }, - { - "cell_type": "markdown", - "id": "3f29f4b4-f889-4cbb-89be-4dc5d8c7e5d1", - "metadata": {}, - "source": [ - "\n", - "## Table of Contents\n", - "* [Introduction](#intro)\n", - "* [Imports](#imports)\n", - "* [Data](#data)\n", - "* [Run the wisp subtraction code](#wisp_sub)\n", - "* [Run the wisp subtraction code using custom settings](#wisp_sub_custom)\n", - " * [New wisp scaling settings](#ex1)\n", - " * [Data quality flagging](#ex2)\n", - " * [Other misc settings](#ex3)\n", - "* [Using custom wisp templates](#custom_templates)\n", - "* [Running the code in pieces](#pieces)" + "**Author**: Ben Sunnquist (bsunnquist@stsci.edu)
\n", + "**Latest Update**: 24 July 2024
\n", + "**Use case**: NIRCam Imaging detectors A3, A4, B3, and B4.
\n", + "**Data**: PID 1063 Obs 6 Imaging Flats
\n", + "**Test Pipeline Version**: 1.15.1" ] }, { "cell_type": "markdown", "id": "338e0783-ee8e-4291-9574-755faaa73610", - "metadata": {}, + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, "source": [ "\n", - "## Introduction" + "## Introduction\n", + "\n", + "This notebook demonstrates how to remove \n", + "[wisps](https://jwst-docs.stsci.edu/known-issues-with-jwst-data/nircam-known-issues/nircam-scattered-light-artifacts#NIRCamScatteredLightArtifacts-wispsWisps) from NIRCam imaging data. Wisps are a scattered light feature affecting detectors A3, A4, B3, and B4. For a given filter, wisps appear in the same detector location with only their brightness varying between exposures; therefore, they can be removed from science data by scaling and subtracting a wisp template (i.e. a median combination of all wisp appearances).\n", + "\n", + "The wisp templates used in this notebook will be downloaded in the [Data](#data) section, but they are also available in the \"version3\" folder within the [NIRCam wisp template Box folder](https://stsci.box.com/s/1bymvf1lkrqbdn9rnkluzqk30e8o2bne). Future updates to these templates will be added to this same Box folder, and users are encouraged to check this folder regularly to ensure they have the latest version.\n", + "\n", + "This notebook uses the `subtract_wisp.py` code to scale and subtract the wisps. That code can be used by itself within python, and is preferred if calibrating a large number of files in parallel, but this notebook will be used to demonstrate the various parameters available to optimize wisp removal. For each notebook cell, we'll also show the corresponding command to run the equivalent in python." ] }, { "cell_type": "markdown", - "id": "fa9f29bf-62d3-42c4-984e-4b4b2b16d88d", + "id": "3f29f4b4-f889-4cbb-89be-4dc5d8c7e5d1", "metadata": {}, "source": [ - "This notebook demonstrates how to remove [wisps](https://jwst-docs.stsci.edu/known-issues-with-jwst-data/nircam-known-issues/nircam-scattered-light-artifacts#NIRCamScatteredLightArtifacts-wispsWisps) from NIRCam imaging data. Wisps are a scattered light feature affecting detectors A3, A4, B3, and B4. For a given filter, wisps appear in the same detector location with only their brightness varying between exposures; therefore, they can be removed from science data by scaling and subtracting a wisp template (i.e. a median combination of all wisp appearances).\n", - "\n", - "The wisp templates used in this notebook will be downloaded in the [Data](#data) section, but they are also available in the \"version3\" folder within the [NIRCam wisp template Box folder](https://stsci.box.com/s/1bymvf1lkrqbdn9rnkluzqk30e8o2bne). Future updates to these templates will be added to this same Box folder, and users are encouraged to check this folder regularly to ensure they have the latest version.\n", - "\n", - "This notebook uses the `subtract_wisp.py` code to scale and subtract the wisps. That code can be used by itself within python, and is preferred if calibrating a large number of files in parallel, but this notebook will be used to demonstrate the various parameters available to optimize wisp removal. For each notebook cell, we'll also show the corresponding command to run the equivalent in python." + "\n", + "## Table of Contents\n", + "* [Imports](#imports)\n", + "* [Data](#data)\n", + "* [Run the wisp subtraction code](#Run-the-wisp-subtraction-code)\n", + "* [Run the wisp subtraction code using custom settings](#Run-the-wisp-subtraction-code-using-custom-settings)\n", + " * [New wisp scaling settings](#Example-1:-New-wisp-scaling-method)\n", + " * [Data quality flagging](#Example-2:-Data-quality-flagging)\n", + " * [Other misc settings](#Example-3:-Other-misc-settings)\n", + "* [Using custom wisp templates](#Using-custom-wisp-templates)\n", + "* [Running the code in pieces](#Runing-the-code-in-pieces)" ] }, { @@ -56,7 +60,18 @@ "metadata": {}, "source": [ "\n", - "## Imports" + "## Imports\n", + "\n", + "How to make an environment to run this notebook:\n", + "\n", + "\n", + "conda create -n jwst1140 python=3.11 notebook\n", + "source activate jwst1140\n", + "pip install -r requirements.txt\n", + "\n", + "\n", + "
\n", + "The version of the pipeline used by this notebook can be updated in the `requirments.txt` file before installation." ] }, { @@ -66,30 +81,18 @@ "metadata": {}, "outputs": [], "source": [ - "# How I made this environment:\n", - "# conda create -n jwst1140 python=3.11 notebook\n", - "# source activate jwst1140\n", - "# pip install jwst==1.14.0\n", - "# pip install astroquery\n", - "\n", + "# load functions from the subtract_wisp.py code \n", "from subtract_wisp import make_segmap, process_file, subtract_wisp\n", "\n", + "# load important packages\n", "from astropy.io import fits\n", "from astroquery.mast import Mast, Observations\n", "import glob\n", - "import matplotlib.pyplot as plt\n", - "%matplotlib inline\n", "import os\n", - "import urllib" - ] - }, - { - "cell_type": "markdown", - "id": "cb679760-0f63-4377-b99a-73c690000ff7", - "metadata": {}, - "source": [ - "\n", - "## Data" + "import urllib\n", + "import zipfile\n", + "import matplotlib.pyplot as plt\n", + "%matplotlib inline" ] }, { @@ -97,7 +100,20 @@ "id": "38e477be-b001-49be-9b8f-7fb503686199", "metadata": {}, "source": [ - "We'll use F200W imaging data from program 1063 observation 6 taken on detectors A3, A4, B3, and B4. These are images of sparse fields that are affected by bright wisps. We'll also use the corresponding longwave images for source detection." + "## Data\n", + "\n", + "We'll use F200W imaging data from program 1063 observation 6 taken on\n", + "detectors A3, A4, B3, and B4. These are 6 images of sparse fields that\n", + "are affected by bright wisps. We'll also use the corresponding\n", + "longwave images for source detection since the longwave is not impacted by wisps.\n", + "\n", + "Download the relevant data from MAST using\n", + "[astroquery](https://astroquery.readthedocs.io/en/latest/mast/mast.html).
\n", + "\n", + "
\n", + "Warning: If this notebook is halted during this step the downloaded\n", + "file may be incomplete, and cause crashes later on!\n", + "
" ] }, { @@ -107,8 +123,7 @@ "metadata": {}, "outputs": [], "source": [ - "# Download the relevant data from MAST\n", - "\n", + "# Select the data relevant data from MAST\n", "params = {\"columns\": \"dataURI, filename, exp_type\",\n", " \"filters\": [{\"paramName\": \"program\", \"values\": ['1063']},\n", " {\"paramName\": \"observtn\", \"values\": ['6']},\n", @@ -124,21 +139,25 @@ }, { "cell_type": "markdown", - "id": "2faab36d-49a8-4877-b4eb-4693ce665de1", + "id": "ff808e6b-ddbf-4e58-ac0b-26b1bf24606a", "metadata": {}, "source": [ - "Next, we'll download the wisp templates. These templates are equivalent to those in the \"version3\" folder within the [NIRCam wisp template Box folder](https://stsci.box.com/s/1bymvf1lkrqbdn9rnkluzqk30e8o2bne). Users are encouraged to check this folder regularly for any updates. For this notebook, only the F200W templates are needed (i.e. WISP_NRCA3_F200W_CLEAR.fits, WISP_NRCA4_F200W_CLEAR.fits, WISP_NRCB3_F200W_CLEAR.fits and WISP_NRCB4_F200W_CLEAR.fits)." + "Next, we'll download the wisp templates. These templates are equivalent to those in the \"version3\" folder within the [NIRCam wisp template Box folder](https://stsci.box.com/s/1bymvf1lkrqbdn9rnkluzqk30e8o2bne). Users are encouraged to check this folder regularly for any updates. For this notebook, only the F200W templates are needed:\n", + "\n", + "* WISP_NRCA3_F200W_CLEAR.fits\n", + "* WISP_NRCA4_F200W_CLEAR.fits\n", + "* WISP_NRCB3_F200W_CLEAR.fits\n", + "* WISP_NRCB4_F200W_CLEAR.fits" ] }, { "cell_type": "code", "execution_count": null, - "id": "641cd018-8760-47fb-beb2-a4f0d9a11061", + "id": "ef409ef8-280c-4156-bd78-f1e06e5820e9", "metadata": {}, "outputs": [], "source": [ "# Download the wisp templates\n", - "\n", "boxlink = 'https://data.science.stsci.edu/redirect/JWST/jwst-data_analysis_tools/nircam_wisp_templates/'\n", "\n", "for detector in ['NRCA3', 'NRCA4', 'NRCB3', 'NRCB4']:\n", @@ -147,12 +166,19 @@ ] }, { - "cell_type": "markdown", - "id": "4ee4c38e-d86d-4096-833e-9e92a87cb78a", + "cell_type": "code", + "execution_count": null, + "id": "ad7f9846-e6b2-4d94-9d96-94057487c980", "metadata": {}, + "outputs": [], "source": [ - "\n", - "## Run the wisp subtraction code" + "# Check if template files are downloaded\n", + "template_files = ['WISP_NRCA3_F200W_CLEAR.fits', 'WISP_NRCA4_F200W_CLEAR.fits',\n", + " 'WISP_NRCB3_F200W_CLEAR.fits', 'WISP_NRCB4_F200W_CLEAR.fits']\n", + "\n", + "for tem_file in template_files:\n", + " if not os.path.isfile(tem_file):\n", + " print(f'file {tem_file} does not exist, please downoad it from Box')" ] }, { @@ -160,12 +186,14 @@ "id": "32af21b1-b980-433c-907a-578d34629339", "metadata": {}, "source": [ + "## Run the wisp subtraction code\n", + "\n", "In this section, we show how to run the wisp subtraction code on all of the shortwave images using the default parameters. The exception is we set `show_plot=True` to allow the diagnostic plots to display in the notebook.\n", "\n", - "The basic process the code uses is the following:\n", + "The main processing function `process_file()` combines the segmap creation and wisp scaling/subtraction steps together. The basic process the code uses is the following:\n", "\n", "1. Generate a source mask using the corresponding longwave image, since the longwave isn't affected by wisps. Ignore these pixels going forward. Occurs in the `make_segmap()` function.\n", - "2. Load the relevant wisp template. Occurs in the `process_file()` function.\n", + "2. Load the relevant wisp template.\n", "3. Apply a series of scale factors to the template. For each scale factor, multiply it onto the wisp template and subtract the result from the input science data. Record the residual noise within the wisp region for each of the scale factors tested. A correction to residual 1/f (horizontal) noise is applied during this step to help with the fitting. Fit a polynomial to the residuals, and choose the scale factor with the lowest noise based on this fit. Occurs in the `subtract_wisp()` function.\n", "4. Multiply the chosen wisp factor onto the wisp template to generate the final wisp model. Occurs in the `subtract_wisp()` function.\n", "5. Subtract the final wisp model from the input data, and write out the wisp-subtracted data, the wisp model, and diagnostic plots summarizing the results. Occurs in the `subtract_wisp()` function." @@ -174,35 +202,51 @@ { "cell_type": "code", "execution_count": null, - "id": "c196dbd7-902a-42b1-9f85-5402b5220fcb", - "metadata": {}, + "id": "24f58c3a-1f7b-4fdc-aa7f-fbe824c2358f", + "metadata": { + "scrolled": true + }, "outputs": [], "source": [ + "# Collect all the science data files \n", "files = glob.glob('*_cal.fits')\n", - "files = [f for f in files if 'long' not in f] # only want shortwave files\n", - "for file in files:\n", - " results = process_file(file, show_plot=True)" + "\n", + "files = [cal_file for cal_file in files if 'long' not in cal_file] # only want shortwave files" ] }, { "cell_type": "markdown", - "id": "9418c448-5a14-41ac-b14e-8fad1011fc49", + "id": "3b16a76b-75cf-4861-81f0-6b249d0533f7", "metadata": {}, "source": [ - "Python equivalent: `python subtract_wisp.py --files *_cal.fits`" + "Run the `process_file()` function in each dataset. The code generates the following:\n", + "\n", + "* A wisp-corrected version of the input file, with the same name as the input file besides a _wisp suffix added (note: the suffix used here can be changed using the `suffix` argument, which can also be set to an empty string to overwrite the input file).\n", + "* The wisp model subtracted from the input image, with the same name as the input file besides a _wisp_model suffix added.\n", + "* A diagnostic plot showing the original image, the wisp-corrected image, the wisp model subtracted, and a plot showing all of the wisp scale factors tested along with their corresponding residual noise. The factor with the lowest noise is the one applied to the wisp template to generate the wisp model.\n", + "\n", + "Note: calling `process_files()` rather than `process_file()` will process the files in parallel. This will save time, but the plots won't display inline in the notebook." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "c196dbd7-902a-42b1-9f85-5402b5220fcb", + "metadata": {}, + "outputs": [], + "source": [ + "# Remove wisps from the science data\n", + "for file in files:\n", + " results = process_file(file, show_plot=True)" ] }, { "cell_type": "markdown", - "id": "3b16a76b-75cf-4861-81f0-6b249d0533f7", + "id": "7e5c7f06-a38f-4a37-85d1-5a4cdbc24dc5", "metadata": {}, "source": [ - "As seen above, the code generates the following:\n", - "* A wisp-corrected version of the input file, with the same name as the input file besides a _wisp suffix added (note: the suffix used here can be changed using the `suffix` argument, which can also be set to an empty string to overwrite the input file).\n", - "* The wisp model subtracted from the input image, with the same name as the input file besides a _wisp_model suffix added.\n", - "* A diagnostic plot showing the original image, the wisp-corrected image, the wisp model subtracted, and a plot showing all of the wisp scale factors tested along with their corresponding residual noise. The factor with the lowest noise is the one applied to the wisp template to generate the wisp model.\n", - "\n", - "Note: calling `process_files()` rather than `process_file()` in the cell above will process the files in parallel. This will save time, but the plots won't display inline in the notebook." + "Python equivalent: \n", + "`python subtract_wisp.py --files *_cal.fits`" ] }, { @@ -250,13 +294,15 @@ "metadata": {}, "outputs": [], "source": [ - "results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits', scale_method='median', poly_degree=0, \n", + "# Remove wisps using the scaling method\n", + "results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits',\n", + " scale_method='median', poly_degree=0, \n", " correct_cols=True, show_plot=True)" ] }, { "cell_type": "markdown", - "id": "0a3f8b44-659c-40b0-9ebc-fc52657ecfd8", + "id": "e5f0466b-afed-4df9-a82f-bb6cb39eb9f0", "metadata": {}, "source": [ "Python equivalent: `python subtract_wisp.py --files jw01063006004_02101_00005_nrcb4_cal.fits --scale_method median --poly_degree 0 --correct_cols`" @@ -284,7 +330,7 @@ "id": "e1817da3-aa13-471c-aa3f-ff2af64e0f10", "metadata": {}, "source": [ - "Another option is to flag any pixels with wisp signals above a certain threshold in the image's data quality array, but not subtract the wisp itself from the science data; this can be accomplished by setting `flag_wisp_thresh` to the chosen wisp signal threshold (units are assumed to be the same as the input file) and `sub_wisp=False`. The data quality value used here is set by the `dq_val` parameter.\n", + "Another option is to flag any pixels with wisp signals above a certain threshold in the image's data quality array, but not subtract the wisp itself from the science data. This may be useful if you want to ignore these pixels when combining several images into a mosaic or flag them while performing photometry. This can be accomplished by setting `flag_wisp_thresh` to the chosen wisp signal threshold (units are assumed to be the same as the input file) and `sub_wisp=False`. The data quality value used here is set by the `dq_val` parameter.\n", "\n", "The default value of `dq_val` is 1, which is equivalent to the JWST Pipeline DO_NOT_USE flag, and it is ignored by default during image3 drizzling. If you would prefer the pixels to be flagged, but not ignored during that step, another option is using `dq_val=1073741824`, which is equivalent to the JWST Pipeine OTHER_BAD_PIXEL flag. More details on the various data quality flags can be found [here](https://jwst-pipeline.readthedocs.io/en/latest/jwst/references_general/references_general.html#data-quality-flags)." ] @@ -296,8 +342,10 @@ "metadata": {}, "outputs": [], "source": [ - "results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits', flag_wisp_thresh=0.03, dq_val=1073741824, sub_wisp=False,\n", - " show_plot=True)" + "# Do not remove wisps, just flag them\n", + "results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits',\n", + " flag_wisp_thresh=0.03, dq_val=1073741824,\n", + " sub_wisp=False, show_plot=True)" ] }, { @@ -308,6 +356,14 @@ "Python equivalent: `python subtract_wisp.py --files jw01063006004_02101_00005_nrcb4_cal.fits --flag_wisp_thresh 0.03 --dq_val 1073741824 --no-sub_wisp`\n" ] }, + { + "cell_type": "markdown", + "id": "20506cbb-6f6a-41a0-b6cf-4492b7ce4859", + "metadata": {}, + "source": [ + "Check that the data quality array in this file was updated appropriately" + ] + }, { "cell_type": "code", "execution_count": null, @@ -315,10 +371,11 @@ "metadata": {}, "outputs": [], "source": [ - "# Check that the data quality array in this file was updated appropriately\n", - "\n", + "# Check data quality (DQ) array \n", "dq = fits.getdata('jw01063006004_02101_00005_nrcb4_cal_wisp.fits', 'DQ')\n", - "dq = (dq & 1073741824 != 0).astype(int) # only want to see pixels flagged as OTHER_BAD_PIXEL, i.e. the dq_val used above\n", + "\n", + "# Show pixels flagged as OTHER_BAD_PIXEL\n", + "dq = (dq & 1073741824 != 0).astype(int) # only want to see pixels with dq_val 1073741824\n", "plt.imshow(dq, cmap='gray', origin='lower', vmin=0, vmax=0.1)" ] }, @@ -336,11 +393,22 @@ "id": "6fbfb5c2-0b40-4fa4-8890-3fdb1a66e284", "metadata": {}, "source": [ - "This example will show a variety of other settings that may be useful either used alone or in combination with any of the other settings shown above. The full list of optional settings can be seen in the `make_segmap()` and `subtract_wisp()` docstrings.\n", + "This example will show a variety of other settings that may be\n", + "useful either used alone or in combination with any of the\n", + "other settings shown above. The full list of optional settings\n", + "can be seen in the `make_segmap()` and `subtract_wisp()` docstrings.\n", "\n", - "In this example, we'll do the source-finding on the input image itself (rather than the corresponding longwave image) by setting `seg_from_lw=False`. We'll increase the source detection `sigma` to help avoid flagging the wisp itself, and we'll save the the resulting segmentation mask by setting `save_segmap=True` for inspection. The segmentation mask filename is the same as the input file, just with the _seg suffix added.\n", + "In this example, we'll do the source-finding on the input image\n", + "itself (rather than the corresponding longwave image) by\n", + "setting `seg_from_lw=False`. We'll increase the source\n", + "detection `sigma` to help avoid flagging the wisp itself, and\n", + "we'll save the resulting segmentation mask by setting\n", + "`save_segmap=True` for inspection. The segmentation mask\n", + "filename is the same as the input file, just with the\n", + "_seg suffix added.\n", "\n", - "We'll also only subtract wisp signals above 0.01 MJy/sr (`min_wisp=0.01`). To save a little time, we'll only test wisp scale factors between 0.5 and 1.5 in 0.05 steps." + "We'll also only subtract wisp signals above 0.01 MJy/sr\n", + "(`min_wisp=0.01`). To save a little time, we'll only test wispscale factors between 0.5 and 1.5 in 0.05 steps." ] }, { @@ -350,8 +418,13 @@ "metadata": {}, "outputs": [], "source": [ - "results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits', seg_from_lw=False, sigma=1.5, save_segmap=True, \n", - " min_wisp=0.01, factor_min=0.5, factor_max=1.5, factor_step=0.05, show_plot=True)" + "# Get full list of optional settings for make_segmap() using help(make_segmap)\n", + "#\n", + "# Use special settings\n", + "results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits',\n", + " seg_from_lw=False, sigma=1.5, save_segmap=True, \n", + " min_wisp=0.01, factor_min=0.5, factor_max=1.5,\n", + " factor_step=0.05, show_plot=True)" ] }, { @@ -380,6 +453,7 @@ "data = fits.getdata('jw01063006004_02101_00005_nrcb4_cal.fits')\n", "segmap = fits.getdata('jw01063006004_02101_00005_nrcb4_cal_seg.fits')\n", "\n", + "# Show the images\n", "fig, axes = plt.subplots(1, 2, figsize=(20, 10))\n", "axes[0].imshow(data, origin='lower', cmap='gray', vmin=0.18, vmax=0.3)\n", "axes[1].imshow(segmap, origin='lower', cmap='gray', vmin=0, vmax=0.1)" @@ -443,22 +517,22 @@ "outputs": [], "source": [ "# The file to process\n", - "f = 'jw01063006004_02101_00005_nrcb4_cal.fits'\n", + "cal_file = 'jw01063006004_02101_00005_nrcb4_cal.fits'\n", "\n", "# The wisp template to use. We'll use a \"custom\" template again.\n", "wisp_data = fits.getdata('WISP_NRCB4_F200W_CLEAR.fits') * 1.05\n", "\n", "# Generate the segmentation map\n", - "segmap_data = make_segmap(f, sigma=1.0, npixels=8, dilate_segmap=10)\n", + "segmap_data = make_segmap(cal_file, sigma=1.0, npixels=8, dilate_segmap=10)\n", "\n", "# Scale and subtract the wisp template\n", - "results = subtract_wisp(f, wisp_data=wisp_data, segmap_data=segmap_data, show_plot=True)" + "results = subtract_wisp(cal_file, wisp_data=wisp_data, segmap_data=segmap_data, show_plot=True)" ] }, { "cell_type": "code", "execution_count": null, - "id": "69682d25-1b97-4b00-801b-36acd502bc81", + "id": "669c194b-4144-45ce-b5c0-cf6a349853de", "metadata": {}, "outputs": [], "source": [] From 0504b7ea2c17739d830c9d7bc913f94cbf19e1ee Mon Sep 17 00:00:00 2001 From: Ben Sunnquist Date: Fri, 13 Sep 2024 15:14:21 -0400 Subject: [PATCH 8/9] removed unnecessary import --- .../NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb | 1 - 1 file changed, 1 deletion(-) diff --git a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb index 5416b5365..ecf27d95a 100644 --- a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb +++ b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb @@ -90,7 +90,6 @@ "import glob\n", "import os\n", "import urllib\n", - "import zipfile\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline" ] From 1f09e8e5454304556cfcd6c88d619bf1f4672e1b Mon Sep 17 00:00:00 2001 From: Ben Sunnquist Date: Fri, 13 Sep 2024 16:56:40 -0400 Subject: [PATCH 9/9] updating requirements to jwst 1.15.1 and removing empty cell --- .../NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb | 8 -------- notebooks/NIRCam/NIRCam_wisp_subtraction/requirements.txt | 2 +- 2 files changed, 1 insertion(+), 9 deletions(-) diff --git a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb index ecf27d95a..573968867 100644 --- a/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb +++ b/notebooks/NIRCam/NIRCam_wisp_subtraction/nircam_wisp_subtraction.ipynb @@ -527,14 +527,6 @@ "# Scale and subtract the wisp template\n", "results = subtract_wisp(cal_file, wisp_data=wisp_data, segmap_data=segmap_data, show_plot=True)" ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "669c194b-4144-45ce-b5c0-cf6a349853de", - "metadata": {}, - "outputs": [], - "source": [] } ], "metadata": { diff --git a/notebooks/NIRCam/NIRCam_wisp_subtraction/requirements.txt b/notebooks/NIRCam/NIRCam_wisp_subtraction/requirements.txt index f6365f22e..161385ab4 100644 --- a/notebooks/NIRCam/NIRCam_wisp_subtraction/requirements.txt +++ b/notebooks/NIRCam/NIRCam_wisp_subtraction/requirements.txt @@ -1,2 +1,2 @@ astroquery==0.4.7 -jwst==1.14.0 \ No newline at end of file +jwst==1.15.1 \ No newline at end of file