From f172eb77fb2720411f0e3d5935c321f88f6665ea Mon Sep 17 00:00:00 2001 From: ojustino Date: Fri, 21 May 2021 14:51:14 -0400 Subject: [PATCH 1/7] Wrote notebook on synthetic MOSviz image creation --- .../concepts/generate_mosviz_photometry.ipynb | 480 ++++++++++++++++++ 1 file changed, 480 insertions(+) create mode 100644 notebooks/concepts/generate_mosviz_photometry.ipynb diff --git a/notebooks/concepts/generate_mosviz_photometry.ipynb b/notebooks/concepts/generate_mosviz_photometry.ipynb new file mode 100644 index 0000000000..53037a6321 --- /dev/null +++ b/notebooks/concepts/generate_mosviz_photometry.ipynb @@ -0,0 +1,480 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "e5e8735c-0e9a-4533-a4e2-f2b6b9ece9f9", + "metadata": {}, + "source": [ + "# Synthetic image creation for MOSviz pipeline data\n", + "\n", + "**Motiviation**: The synthetic dataset we've received from the JWST data pipeline team for use in MOSviz contains simulated spectra from NIRSpec, but no simulated photometry from NIRCam. We'd like to have test imagery to display in MOSviz alongside the 2D and 1D spectra we already possess. We've learned that the pipeline team has no plans to produce any, so we attempt to piece together our own.\n", + "\n", + "**Goal**: Populate an image of background noise with properly-scaled galaxy cutouts sourced from an HST image and placed at their respective WCS locations. These galaxies' real spectra don't necessarily correspond with those in our dataset, but at this point we care more about the veneer of having photometry to match with our spectra.\n", + "\n", + "**Execution**: We pull our galaxy cutouts and catalog information from the Hubble Deep Field image. ~~ASTRODEEP's image of the [Abell 2774 Parallel](http://astrodeep.u-strasbg.fr/ff/?img=JH140?cm=grayscale) | [MACS J0416.1-2403 Parallel](http://astrodeep.u-strasbg.fr/ff/?ffid=FF_M0416PAR&id=1264&cm=grayscale)~~. We sought to use [Artifactory](https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth) to obtain a range of RA/Dec over which to project our synthetic image, but that information was absent. Instead, we place the image over a manually chosen RA/Dec range and place the cutouts in randomly selected locations within that field of view.\n", + "\n", + "**Issues**:\n", + "- We wanted to scrape RA/Dec information the data pipeline products to get a range of coordinates over which to scale our synthetic image, but it appears that the pipeline's data products don't have `\"TARG_RA\"` or `\"TARG_DEC\"` keywords in their headers.\n", + " - The data products also don't appear to have WCS information. We don't strictly need it to achieve this notebook's goals, but it would be convenient to have.\n", + "- There doesn't appear to be an observation with level 2 data in `jwst-pipeline/truth/test_nirspec_mos_spec2` and level 3 data in `jwst-pipeline/truth/test_nirspec_mos_spec3`. All observations are either level 2 only or level 3 only.\n", + "- _(Resolved)_ A good number of the cutouts from the first couple of field images we tested had intrustions from other galaxies due to crowding. We settled on the Hubble Deep Field as a good source image, but if we hadn't, we may have considered using galaxies modeled with Sersic profiles to get cleaner cutouts to inject into our synthetic image.\n", + "\n", + "*Writers: Robel Geda and O. Justin Otor*" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "aba137f6-9bfa-4371-b927-21f30bf9b7fd", + "metadata": {}, + "outputs": [], + "source": [ + "from astropy.io import fits\n", + "from astropy.nddata import block_reduce, Cutout2D\n", + "from astropy.stats import sigma_clipped_stats, sigma_clip\n", + "from astropy.table import Table, join\n", + "from astropy.wcs import WCS\n", + "from glob import glob\n", + "\n", + "import matplotlib.pyplot as plt\n", + "import numpy as np" + ] + }, + { + "cell_type": "markdown", + "id": "fffd44de-ee89-4a46-8c75-198a94cf124d", + "metadata": {}, + "source": [ + "### Generate galaxy cutouts" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "e172b883-f625-4a41-b8a5-a2d13d41b396", + "metadata": {}, + "outputs": [], + "source": [ + "# download the image containing sources to be cut out later\n", + "image_fits = fits.open('https://archive.stsci.edu/pub/hlsp/hdf/v2/mosaics/x4096/f814_mosaic_v2.fits')\n", + "image_header = image_fits[0].header\n", + "image_data = image_fits[0].data\n", + "\n", + "image_data.shape" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "354bce9f-dc43-434b-a992-b6a636eda576", + "metadata": {}, + "outputs": [], + "source": [ + "# download those sources' locations in the image and sort them by brightness\n", + "source_info1 = Table.read('https://archive.stsci.edu/pub/hlsp/hdf/wfpc_hdfn_v2catalog/HDFN_wfpc_v2generic.cat',\n", + " format='ascii')\n", + "source_info2 = Table.read('https://archive.stsci.edu/pub/hlsp/hdf/wfpc_hdfn_v2catalog/HDFN_f814_v2.cat',\n", + " format='ascii')\n", + "sources = table_join(source_info1, source_info2)\n", + "\n", + "# confirm that both tables contain the same objects in the same order\n", + "(source_info1['NUMBER'] == source_info2['NUMBER']).sum() == len(source_info1) == len(source_info2)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "dedfe0f0-b1e4-49e7-b51a-bdf7f069fa06", + "metadata": {}, + "outputs": [], + "source": [ + "# sort sources by flux contained within 71.1 pixel diameter of source (11)\n", + "sources.sort('FLUX_APER_11', reverse=True)\n", + "#sources.sort('KRON_RADIUS', reverse=True)\n", + "#sources.sort('FWHM_IMAGE')\n", + "\n", + "# filter out likely stars and sources with negative flux\n", + "sources = sources[(sources['CLASS_STAR'] < .5)\n", + " & (sources['FLUX_APER_8'] > 0)]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "3efabf0f-ddbb-4900-9ba8-528998b4bb69", + "metadata": {}, + "outputs": [], + "source": [ + "sources[:5]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "f5dcc4be-4e10-45c0-a084-4d224abd6ba1", + "metadata": {}, + "outputs": [], + "source": [ + "# convert the sources' WCS locations to in-image pixel values\n", + "image_wcs = WCS(image_fits[0].header)\n", + "sources_x, sources_y = image_wcs.world_to_pixel_values(sources['ALPHA_J2000'],\n", + " sources['DELTA_J2000'])" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "34f1f230-88da-44c3-803d-bb405fe6f6e7", + "metadata": { + "scrolled": true, + "tags": [] + }, + "outputs": [], + "source": [ + "# save a list of good cutouts for later use\n", + "cutout_list = []\n", + "first_source = 0\n", + "catalog_size = 20\n", + "downsample_factor = 2\n", + "patch_length = 100\n", + "\n", + "for x, y in list(zip(sources_x, sources_y))[first_source:]:\n", + " # use pixel locations to cut a source from the image\n", + " cutout = Cutout2D(image_data, (x, y),\n", + " patch_length * downsample_factor).data\n", + " \n", + " # bin by downsample_factor to increase field of view\n", + " cutout = block_reduce(cutout, downsample_factor)\n", + " \n", + " # skip any cutouts that extend past the image border\n", + " if ( np.all(cutout[-1] <= 0) or np.all(cutout[0] <= 0)\n", + " or np.all(cutout[:,-1] <= 0) or np.all(cutout[:,0] <= 0) ):\n", + " continue\n", + " \n", + " # save and plot the new cutout\n", + " cutout_list.append(cutout)\n", + " \n", + " plt.imshow(cutout, vmin=-1e-5, vmax=image_data.std(),\n", + " origin='lower', cmap='bone')\n", + " plt.show()\n", + " \n", + " if len(cutout_list) == catalog_size:\n", + " break" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "5b0a6505-3bdd-4f8d-8c2c-e2b7b489d078", + "metadata": {}, + "outputs": [], + "source": [ + "clipped_mean, clipped_median, clipped_stddev = sigma_clipped_stats(image_data, sigma=3.)" + ] + }, + { + "cell_type": "markdown", + "id": "53891064-621c-4f6f-91da-3a24741ad632", + "metadata": {}, + "source": [ + "### Extract destination RA/Dec from spectra files" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "9436aa5d-54c9-4a1e-a983-7f8658920be6", + "metadata": {}, + "outputs": [], + "source": [ + "# save level 3 spectra FITS header information\n", + "# (is there a way to do so using the URL?)\n", + "# 'https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth%2Ftest_nirspec_mos_spec3%2Fjw00626-o030_s00000_nirspec_f170lp-g235m_s2d.fits'\n", + "x1d_header = fits.getheader('/Users/jotor/Downloads/jw00626-o030_s00000_nirspec_f170lp-g235m_x1d.fits')\n", + "s2d_header = fits.getheader('/Users/jotor/Downloads/jw00626-o030_s00000_nirspec_f170lp-g235m_s2d.fits')" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "bb2858ac-e8e9-4c80-a2ac-72947ed7195b", + "metadata": {}, + "outputs": [], + "source": [ + "for k, v in x1d_header.items():\n", + " if 'or' in k.lower():\n", + " print(k, v)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "96dda8bd-5010-45ef-a937-b99ee54b8091", + "metadata": {}, + "outputs": [], + "source": [ + "(x1d_header['TARG_RA'], x1d_header['TARG_DEC'],\n", + " s2d_header['TARG_RA'], s2d_header['TARG_DEC'])" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "a4ecfbe4-bc85-4d94-b660-e50d0bf621ba", + "metadata": {}, + "outputs": [], + "source": [ + "# no WCS data... not what we want\n", + "WCS(x1d_header)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "879040fb-eb58-4e4c-8173-a289b21a61e2", + "metadata": {}, + "outputs": [], + "source": [ + "# search for RA/Dec information from Artifactory observation files\n", + "x1d_header_list = [fits.getheader(file) for file in glob('/Users/jotor/Downloads/jw00626*x1d.fits')]\n", + "\n", + "# all the same... not what we want\n", + "ras, decs = np.array([[h['TARG_RA'], h['TARG_DEC']] for h in x1d_header_list]).T\n", + "ras, decs" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "9931cde4-5f77-4a40-96e5-32d252bec47c", + "metadata": {}, + "outputs": [], + "source": [ + "# since all information is the same, we randomly generate our own RAs/decs\n", + "# (ranges based on NIRSpec MSA's on-sky projection size of 3.6x3.4 arcmins)\n", + "np.random.seed(19)\n", + "ras = np.random.uniform(0, 1/15, catalog_size)\n", + "decs = np.random.uniform(-1/30, 1/30, catalog_size)" + ] + }, + { + "cell_type": "markdown", + "id": "c4b5da62-c424-4603-adc3-4c550291d13a", + "metadata": {}, + "source": [ + "### Create synthetic image" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "6d1f9e28-a5fa-4091-bea6-72c88b48e564", + "metadata": {}, + "outputs": [], + "source": [ + "# create synthetic image onto which cutouts will be pasted\n", + "synth_img_size = 1000\n", + "synth_image = np.zeros((synth_img_size, synth_img_size))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "cb856895-9b87-44c8-8c4b-2cbab45edb76", + "metadata": {}, + "outputs": [], + "source": [ + "# add noise\n", + "synth_image += np.random.normal(loc=clipped_mean, scale=clipped_stddev*8,\n", + " size=synth_image.shape)\n", + "# synth_image += np.random.normal(loc=image_data.mean(), scale=image_data.std(),\n", + "# size=synth_image.shape)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "3b52414b-453c-43eb-984b-e856c284cb95", + "metadata": {}, + "outputs": [], + "source": [ + "plt.imshow(synth_image, cmap='bone')\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "id": "674017bb-2cc1-437a-bf3a-449f306fa506", + "metadata": {}, + "source": [ + "### Fill out new WCS object for `synth_image`" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "1cd4ee37-d41f-4ddd-960c-447a72b1faa2", + "metadata": {}, + "outputs": [], + "source": [ + "synth_wcs = WCS(naxis=2)\n", + "synth_wcs" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "0d8aa070-1b6f-4a0a-9975-58bf5634ff86", + "metadata": { + "tags": [] + }, + "outputs": [], + "source": [ + "ra_bounds = np.array([ras.max(), ras.min()])\n", + "dec_bounds = np.array([decs.max(), decs.min()])" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "80f9d82c-2a63-42f9-8cdc-2dc5d979062e", + "metadata": {}, + "outputs": [], + "source": [ + "delta_ra = ras.max() - ras.min()\n", + "delta_dec = decs.max() - decs.min()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "d70d1d2a-8c9c-40b7-b31f-3a72377279df", + "metadata": {}, + "outputs": [], + "source": [ + "# set minimum FOV as maximum range in RA or dec\n", + "if delta_ra > delta_dec:\n", + " min_image_fov = abs(delta_ra * np.cos(np.pi/180 * dec_bounds.sum()/2))\n", + "else:\n", + " min_image_fov = delta_dec\n", + " \n", + "min_image_fov" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "308eb86b-5a67-4fbe-bf31-04a08642537a", + "metadata": {}, + "outputs": [], + "source": [ + "# scale this FOV by pixels\n", + "pix_scale = min_image_fov / synth_img_size\n", + "\n", + "# add a buffer to the borders\n", + "pix_scale *= 1.5\n", + "pix_scale" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "09dadad8-902d-4119-a317-c10064968597", + "metadata": { + "tags": [] + }, + "outputs": [], + "source": [ + "synth_wcs.wcs.ctype = ['RA---TAN', 'DEC--TAN']\n", + "\n", + "# match value of center pixel of detector to value of FOV's central coordinate in the sky\n", + "synth_wcs.wcs.crpix = [synth_img_size / 2, synth_img_size / 2]\n", + "synth_wcs.wcs.crval = [ra_bounds.sum() / 2, dec_bounds.sum() / 2]\n", + "\n", + "# distance (in sky coordinates) traversed by one pixel length in each dimension\n", + "synth_wcs.wcs.cdelt = [-pix_scale, pix_scale]\n", + "\n", + "synth_wcs" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "8f92545b-5c5d-4ab3-a9d6-ec61e2d329a6", + "metadata": {}, + "outputs": [], + "source": [ + "# convert source RAs/decs from real coordinates to pixels \n", + "ras_pix, decs_pix = np.round(synth_wcs.world_to_pixel_values(ras, decs)).astype(int)\n", + "ras_pix, decs_pix" + ] + }, + { + "cell_type": "markdown", + "id": "5792df4f-c0d6-49a5-be72-a7c5ab39d9b4", + "metadata": {}, + "source": [ + "### Populate `synth_image` with the cutouts" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "c587e262-15fc-4fec-b5d3-b9d173afc286", + "metadata": {}, + "outputs": [], + "source": [ + "cutout_half_delta = patch_length // 2\n", + "\n", + "for i in range(len(ras_pix)):\n", + " synth_image[ras_pix[i] - cutout_half_delta : ras_pix[i] + cutout_half_delta,\n", + " decs_pix[i] - cutout_half_delta : decs_pix[i] + cutout_half_delta] += cutout_list[i]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "9bccee92-2930-4116-b599-0269fdf1d2dc", + "metadata": {}, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(figsize=(10,10))\n", + "ax.imshow(synth_image, vmin=0, vmax=synth_image.std()*3, origin='lower', cmap='bone')\n", + "#plt.imshow(synth_image, vmin=0, vmax=synth_image.mean()*3, origin='lower', cmap='bone')\n", + "plt.show()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "337dc041-618e-4f5d-bcb2-845cacbfac31", + "metadata": {}, + "outputs": [], + "source": [ + "fits.writeto('synthetic_HDF_more.fits', synth_image,\n", + " header=synth_wcs.to_header(), overwrite=True)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "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.8.8" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} From e93abb1385b1137575bfe84ee7ab03e9984e0c04 Mon Sep 17 00:00:00 2001 From: ojustino Date: Tue, 25 May 2021 12:39:36 -0400 Subject: [PATCH 2/7] Cleaned up synthetic image creation notebook --- .../concepts/generate_mosviz_photometry.ipynb | 215 ++++++++++++------ 1 file changed, 142 insertions(+), 73 deletions(-) diff --git a/notebooks/concepts/generate_mosviz_photometry.ipynb b/notebooks/concepts/generate_mosviz_photometry.ipynb index 53037a6321..a30c6386d0 100644 --- a/notebooks/concepts/generate_mosviz_photometry.ipynb +++ b/notebooks/concepts/generate_mosviz_photometry.ipynb @@ -2,24 +2,39 @@ "cells": [ { "cell_type": "markdown", - "id": "e5e8735c-0e9a-4533-a4e2-f2b6b9ece9f9", + "id": "01609adc-98ac-41d2-8030-8eb545ca2fb0", "metadata": {}, "source": [ "# Synthetic image creation for MOSviz pipeline data\n", "\n", - "**Motiviation**: The synthetic dataset we've received from the JWST data pipeline team for use in MOSviz contains simulated spectra from NIRSpec, but no simulated photometry from NIRCam. We'd like to have test imagery to display in MOSviz alongside the 2D and 1D spectra we already possess. We've learned that the pipeline team has no plans to produce any, so we attempt to piece together our own.\n", + "**Motiviation**: The synthetic dataset we currently possess from the JWST data pipeline team for use in MOSviz contains simulated spectra from NIRSpec but no associated NIRCam photometry. We would like to have test imagery to display in MOSviz alongside these 2D and 1D spectra. We are unsure whether the pipeline team has plans to produce any.\n", "\n", - "**Goal**: Populate an image of background noise with properly-scaled galaxy cutouts sourced from an HST image and placed at their respective WCS locations. These galaxies' real spectra don't necessarily correspond with those in our dataset, but at this point we care more about the veneer of having photometry to match with our spectra.\n", + "**Goal**: Populate an image of background noise with properly-scaled galaxy cutouts sourced from a Hubble Space Telescope image and placed at their analogous locations in the new image. These galaxies' real spectra do not necessarily correspond with those in our dataset, but we care more about the veneer of having photometry to match with our spectra at this point.\n", "\n", "**Execution**: We pull our galaxy cutouts and catalog information from the Hubble Deep Field image. ~~ASTRODEEP's image of the [Abell 2774 Parallel](http://astrodeep.u-strasbg.fr/ff/?img=JH140?cm=grayscale) | [MACS J0416.1-2403 Parallel](http://astrodeep.u-strasbg.fr/ff/?ffid=FF_M0416PAR&id=1264&cm=grayscale)~~. We sought to use [Artifactory](https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth) to obtain a range of RA/Dec over which to project our synthetic image, but that information was absent. Instead, we place the image over a manually chosen RA/Dec range and place the cutouts in randomly selected locations within that field of view.\n", "\n", "**Issues**:\n", - "- We wanted to scrape RA/Dec information the data pipeline products to get a range of coordinates over which to scale our synthetic image, but it appears that the pipeline's data products don't have `\"TARG_RA\"` or `\"TARG_DEC\"` keywords in their headers.\n", + "- We wanted to scrape RA/Dec information the data pipeline products to get a range of coordinates over which to scale our synthetic image, but it appears that the pipeline's data products lack `\"TARG_RA\"` or `\"TARG_DEC\"` keywords in their headers.\n", " - The data products also don't appear to have WCS information. We don't strictly need it to achieve this notebook's goals, but it would be convenient to have.\n", - "- There doesn't appear to be an observation with level 2 data in `jwst-pipeline/truth/test_nirspec_mos_spec2` and level 3 data in `jwst-pipeline/truth/test_nirspec_mos_spec3`. All observations are either level 2 only or level 3 only.\n", - "- _(Resolved)_ A good number of the cutouts from the first couple of field images we tested had intrustions from other galaxies due to crowding. We settled on the Hubble Deep Field as a good source image, but if we hadn't, we may have considered using galaxies modeled with Sersic profiles to get cleaner cutouts to inject into our synthetic image.\n", + "- There does not appear to be an observation with level 2 data in `jwst-pipeline/truth/test_nirspec_mos_spec2` and level 3 data in `jwst-pipeline/truth/test_nirspec_mos_spec3`. All observations are either level 2 only or level 3 only.\n", + "- _(Resolved)_ A good number of the cutouts from the first couple of field images we tested had intrustions from other galaxies due to crowding. We settled on the Hubble Deep Field as a good source image, but had we not, we may have considered using galaxies modeled with Sersic profiles to get cleaner cutouts to inject into our synthetic image." + ] + }, + { + "cell_type": "markdown", + "id": "57df882f-e15d-4f16-931c-9f495893c0e1", + "metadata": {}, + "source": [ + "### Import packages\n", "\n", - "*Writers: Robel Geda and O. Justin Otor*" + "- We use `astropy.io.fits` to read in existing FITS files and write a new one with the synthetic image.\n", + "- The objects from `astropy.nddata` help with creating cutouts once we've identified galaxies we'd like to take from the field image.\n", + "- The methods from `astropy.stats` work with image data that's clipped to within a certain number of deviations from the mean.\n", + "- The objects from `astropy.table` help with reading an modifiying tabular data.\n", + "- `astropy.wcs.WCS` creates a World Coordinate System object that's useful for transforming on-sky data from one patch to another.\n", + "- `glob.glob` lists local files that match a given pattern.\n", + "- We use `matplotlib.pyplot` to preview the field image, the cutouts, and various stages of our synthetic image.\n", + "- We use `numpy` to facilitate several specialized mathematical and array-related operations." ] }, { @@ -45,7 +60,9 @@ "id": "fffd44de-ee89-4a46-8c75-198a94cf124d", "metadata": {}, "source": [ - "### Generate galaxy cutouts" + "### Generate galaxy cutouts\n", + "\n", + "The galaxy cutouts come from the Hubble Deep Field image. To retrieve them, we download the image itself and its associated catalog data, search the catalog for the brightest galaxies, then find those galaxies in the image." ] }, { @@ -70,15 +87,17 @@ "metadata": {}, "outputs": [], "source": [ - "# download those sources' locations in the image and sort them by brightness\n", + "# download sources' location and flux information\n", "source_info1 = Table.read('https://archive.stsci.edu/pub/hlsp/hdf/wfpc_hdfn_v2catalog/HDFN_wfpc_v2generic.cat',\n", " format='ascii')\n", "source_info2 = Table.read('https://archive.stsci.edu/pub/hlsp/hdf/wfpc_hdfn_v2catalog/HDFN_f814_v2.cat',\n", " format='ascii')\n", - "sources = table_join(source_info1, source_info2)\n", + "sources = join(source_info1, source_info2)\n", "\n", - "# confirm that both tables contain the same objects in the same order\n", - "(source_info1['NUMBER'] == source_info2['NUMBER']).sum() == len(source_info1) == len(source_info2)" + "# confirm that both tables contain the same objects in the same order (True)\n", + "( (source_info1['NUMBER'] == source_info2['NUMBER']).sum()\n", + " == len(source_info1)\n", + " == len(source_info2) )" ] }, { @@ -88,10 +107,8 @@ "metadata": {}, "outputs": [], "source": [ - "# sort sources by flux contained within 71.1 pixel diameter of source (11)\n", + "# sort sources by flux within 71.1 pixel diameter of source, or aperture 11\n", "sources.sort('FLUX_APER_11', reverse=True)\n", - "#sources.sort('KRON_RADIUS', reverse=True)\n", - "#sources.sort('FWHM_IMAGE')\n", "\n", "# filter out likely stars and sources with negative flux\n", "sources = sources[(sources['CLASS_STAR'] < .5)\n", @@ -121,6 +138,14 @@ " sources['DELTA_J2000'])" ] }, + { + "cell_type": "markdown", + "id": "6b4f1421-e3c3-4112-a475-dc1b7ad07146", + "metadata": {}, + "source": [ + "Note: Depending on the value of `catalog_size`, the following cell can produce a lot of output. Right-click the cell and select \"Enable Scrolling for Outputs\" to expand it or \"Disable Scrolling for Outputs\" to condense it." + ] + }, { "cell_type": "code", "execution_count": null, @@ -162,6 +187,14 @@ " break" ] }, + { + "cell_type": "markdown", + "id": "c3504d50-99cf-4d73-9678-1eff2ca43676", + "metadata": {}, + "source": [ + "We also save image statistics calculated from pixels within a chosen number of standard deviations from the image's mean intensity. Some of them may be useful in creating the synthetic image later on." + ] + }, { "cell_type": "code", "execution_count": null, @@ -169,7 +202,8 @@ "metadata": {}, "outputs": [], "source": [ - "clipped_mean, clipped_median, clipped_stddev = sigma_clipped_stats(image_data, sigma=3.)" + "clipped_mean, clipped_median, clipped_stddev = sigma_clipped_stats(image_data,\n", + " sigma=3.)" ] }, { @@ -177,7 +211,11 @@ "id": "53891064-621c-4f6f-91da-3a24741ad632", "metadata": {}, "source": [ - "### Extract destination RA/Dec from spectra files" + "### Extract destination RA/Dec from spectra files\n", + "\n", + "Note that the following cells only run if on a connection with access to the STScI VPN. They use copies of JWST pipeline files from May 19, 2020.\n", + "\n", + "To run them elsewhere and/or use more current files, visit the [Artifactory](https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth%2Ftest_nirspec_mos_spec3), ensure the `jwst-pipeline/dev/truth/test_nirspec_mos_spec3/` folder is selected in the directory tree, and download all files whose names begin with `jw00626-o030` and end with `nirspec_f170lp-g235m_x1d.fits`. (The original version of the notebook uses six of these Level 3 data products.) Then, change `filepath` in the next cell to the location where you saved the downloaded files. Remember the trailing forward slash." ] }, { @@ -187,45 +225,38 @@ "metadata": {}, "outputs": [], "source": [ - "# save level 3 spectra FITS header information\n", - "# (is there a way to do so using the URL?)\n", - "# 'https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth%2Ftest_nirspec_mos_spec3%2Fjw00626-o030_s00000_nirspec_f170lp-g235m_s2d.fits'\n", - "x1d_header = fits.getheader('/Users/jotor/Downloads/jw00626-o030_s00000_nirspec_f170lp-g235m_x1d.fits')\n", - "s2d_header = fits.getheader('/Users/jotor/Downloads/jw00626-o030_s00000_nirspec_f170lp-g235m_s2d.fits')" + "# view level 3 spectra FITS header information\n", + "filepath = '/user/jotor/jwst-pipeline-lvl3/'\n", + "x1d_header = fits.getheader(filepath + 'jw00626-o030_s00000_nirspec_f170lp-g235m_x1d.fits')\n", + "#s2d_header = fits.getheader(filepath + 'jw00626-o030_s00000_nirspec_f170lp-g235m_s2d.fits')" ] }, { "cell_type": "code", "execution_count": null, - "id": "bb2858ac-e8e9-4c80-a2ac-72947ed7195b", + "id": "96dda8bd-5010-45ef-a937-b99ee54b8091", "metadata": {}, "outputs": [], "source": [ - "for k, v in x1d_header.items():\n", - " if 'or' in k.lower():\n", - " print(k, v)" + "x1d_header['TARG_RA'], x1d_header['TARG_DEC']" ] }, { "cell_type": "code", "execution_count": null, - "id": "96dda8bd-5010-45ef-a937-b99ee54b8091", + "id": "ade024d2-9788-4572-9304-af69c30e86da", "metadata": {}, "outputs": [], "source": [ - "(x1d_header['TARG_RA'], x1d_header['TARG_DEC'],\n", - " s2d_header['TARG_RA'], s2d_header['TARG_DEC'])" + "WCS(x1d_header)" ] }, { - "cell_type": "code", - "execution_count": null, - "id": "a4ecfbe4-bc85-4d94-b660-e50d0bf621ba", + "cell_type": "markdown", + "id": "809be77d-b099-40a3-8a85-42a7a2cf3244", "metadata": {}, - "outputs": [], "source": [ - "# no WCS data... not what we want\n", - "WCS(x1d_header)" + "Notice that this header lacks WCS information. Additionally, examining all of this observation's file headers reveals that they all have the same RA/Dec, which is not what we expect for its different \"pointings.\"" ] }, { @@ -236,13 +267,21 @@ "outputs": [], "source": [ "# search for RA/Dec information from Artifactory observation files\n", - "x1d_header_list = [fits.getheader(file) for file in glob('/Users/jotor/Downloads/jw00626*x1d.fits')]\n", + "x1d_header_list = [fits.getheader(file)\n", + " for file in glob(filepath + 'jw00626*x1d.fits')]\n", "\n", - "# all the same... not what we want\n", "ras, decs = np.array([[h['TARG_RA'], h['TARG_DEC']] for h in x1d_header_list]).T\n", "ras, decs" ] }, + { + "cell_type": "markdown", + "id": "debb3ccc-88ec-4e2d-bf50-69e8361f7d4b", + "metadata": {}, + "source": [ + "Since the headers lack the information we seek, we randomly generate our sources' RA/Dec information in predetermined patch of sky. We take the patch size from the size of NIRSpec Micro-Shutter Assembly (MSA) -- about 3.6'x3.4'-- approximated to a square field of view." + ] + }, { "cell_type": "code", "execution_count": null, @@ -250,8 +289,6 @@ "metadata": {}, "outputs": [], "source": [ - "# since all information is the same, we randomly generate our own RAs/decs\n", - "# (ranges based on NIRSpec MSA's on-sky projection size of 3.6x3.4 arcmins)\n", "np.random.seed(19)\n", "ras = np.random.uniform(0, 1/15, catalog_size)\n", "decs = np.random.uniform(-1/30, 1/30, catalog_size)" @@ -262,7 +299,9 @@ "id": "c4b5da62-c424-4603-adc3-4c550291d13a", "metadata": {}, "source": [ - "### Create synthetic image" + "### Create synthetic image\n", + "\n", + "We initialize a `numpy` array and fill it with normally-distributed background noise based on some of the clipped image statistics that were calculated earlier." ] }, { @@ -272,7 +311,6 @@ "metadata": {}, "outputs": [], "source": [ - "# create synthetic image onto which cutouts will be pasted\n", "synth_img_size = 1000\n", "synth_image = np.zeros((synth_img_size, synth_img_size))" ] @@ -286,9 +324,7 @@ "source": [ "# add noise\n", "synth_image += np.random.normal(loc=clipped_mean, scale=clipped_stddev*8,\n", - " size=synth_image.shape)\n", - "# synth_image += np.random.normal(loc=image_data.mean(), scale=image_data.std(),\n", - "# size=synth_image.shape)" + " size=synth_image.shape)" ] }, { @@ -298,7 +334,8 @@ "metadata": {}, "outputs": [], "source": [ - "plt.imshow(synth_image, cmap='bone')\n", + "imshow_params = {'cmap': 'bone', 'origin': 'lower'}\n", + "plt.imshow(synth_image, **imshow_params)\n", "plt.show()" ] }, @@ -307,7 +344,9 @@ "id": "674017bb-2cc1-437a-bf3a-449f306fa506", "metadata": {}, "source": [ - "### Fill out new WCS object for `synth_image`" + "### Fill out new WCS object for `synth_image`\n", + "\n", + "Creating a `WCS` object for `synth_image` allows it to be mathematically transformed into a projection on the sky. That projection can then be compared to other FITS images with their own `WCS` information." ] }, { @@ -322,27 +361,26 @@ ] }, { - "cell_type": "code", - "execution_count": null, - "id": "0d8aa070-1b6f-4a0a-9975-58bf5634ff86", - "metadata": { - "tags": [] - }, - "outputs": [], + "cell_type": "markdown", + "id": "c1a2876e-963b-4e62-8dd1-d68fb4b4e0a7", + "metadata": {}, "source": [ - "ra_bounds = np.array([ras.max(), ras.min()])\n", - "dec_bounds = np.array([decs.max(), decs.min()])" + "The next step is to calculate field of view information for `synth_image`." ] }, { "cell_type": "code", "execution_count": null, - "id": "80f9d82c-2a63-42f9-8cdc-2dc5d979062e", + "id": "068af407-9a8a-429c-9f7e-03f4f0a0b867", "metadata": {}, "outputs": [], "source": [ - "delta_ra = ras.max() - ras.min()\n", - "delta_dec = decs.max() - decs.min()" + "# find the range of sources in RA and dec\n", + "ra_bounds = np.array([ras.max(), ras.min()])\n", + "dec_bounds = np.array([decs.max(), decs.min()])\n", + "\n", + "delta_ra = np.ptp(ras)\n", + "delta_dec = np.ptp(decs)" ] }, { @@ -352,9 +390,9 @@ "metadata": {}, "outputs": [], "source": [ - "# set minimum FOV as maximum range in RA or dec\n", + "# save the maximum span in coordinates, RA or dec\n", "if delta_ra > delta_dec:\n", - " min_image_fov = abs(delta_ra * np.cos(np.pi/180 * dec_bounds.sum()/2))\n", + " min_image_fov = abs(delta_ra * np.cos(np.pi / 180 * dec_bounds.sum() / 2))\n", "else:\n", " min_image_fov = delta_dec\n", " \n", @@ -368,14 +406,22 @@ "metadata": {}, "outputs": [], "source": [ - "# scale this FOV by pixels\n", + "# scale this field of view (FOV) by pixels\n", "pix_scale = min_image_fov / synth_img_size\n", "\n", - "# add a buffer to the borders\n", + "# add a buffer to the FOV's borders to prevent clipping sources\n", "pix_scale *= 1.5\n", "pix_scale" ] }, + { + "cell_type": "markdown", + "id": "d28eb2e7-d3a1-4190-978a-8be1b30d3a3c", + "metadata": {}, + "source": [ + "With those calculations done, the `WCS` object is ready to be filled out." + ] + }, { "cell_type": "code", "execution_count": null, @@ -397,10 +443,20 @@ "synth_wcs" ] }, + { + "cell_type": "markdown", + "id": "5792df4f-c0d6-49a5-be72-a7c5ab39d9b4", + "metadata": {}, + "source": [ + "### Populate `synth_image` with the cutouts\n", + "\n", + "Finally, we convert the cutouts' coordinates to pixels and add them into `synth_image` to complete the creation of the mock image." + ] + }, { "cell_type": "code", "execution_count": null, - "id": "8f92545b-5c5d-4ab3-a9d6-ec61e2d329a6", + "id": "34c92d83-8e78-4cc5-8df3-728cb51d6567", "metadata": {}, "outputs": [], "source": [ @@ -409,14 +465,6 @@ "ras_pix, decs_pix" ] }, - { - "cell_type": "markdown", - "id": "5792df4f-c0d6-49a5-be72-a7c5ab39d9b4", - "metadata": {}, - "source": [ - "### Populate `synth_image` with the cutouts" - ] - }, { "cell_type": "code", "execution_count": null, @@ -439,11 +487,19 @@ "outputs": [], "source": [ "fig, ax = plt.subplots(figsize=(10,10))\n", - "ax.imshow(synth_image, vmin=0, vmax=synth_image.std()*3, origin='lower', cmap='bone')\n", - "#plt.imshow(synth_image, vmin=0, vmax=synth_image.mean()*3, origin='lower', cmap='bone')\n", + "ax.imshow(synth_image, vmin=0, vmax=synth_image.std()*3, **imshow_params)\n", + "\n", "plt.show()" ] }, + { + "cell_type": "markdown", + "id": "eb2a0045-74b6-4ac2-85a8-c535e1a64800", + "metadata": {}, + "source": [ + "Save the image, overwriting the previous file if it already exists in the current directory." + ] + }, { "cell_type": "code", "execution_count": null, @@ -454,6 +510,19 @@ "fits.writeto('synthetic_HDF_more.fits', synth_image,\n", " header=synth_wcs.to_header(), overwrite=True)" ] + }, + { + "cell_type": "markdown", + "id": "fe8817e2-51f4-4417-84ce-659c9449e87a", + "metadata": {}, + "source": [ + "

\n", + " Authors: Robel Geda and O. Justin Otor \n", + " \"Space\n", + "

\n", + "\n", + "-----" + ] } ], "metadata": { From 3a76e26bb364a7a0774e93e7bff5e476e9fe6193 Mon Sep 17 00:00:00 2001 From: ojustino Date: Fri, 28 May 2021 18:50:54 -0400 Subject: [PATCH 3/7] Incorporated photutils source detection --- .../concepts/generate_mosviz_photometry.ipynb | 428 +++++++++++++----- 1 file changed, 314 insertions(+), 114 deletions(-) diff --git a/notebooks/concepts/generate_mosviz_photometry.ipynb b/notebooks/concepts/generate_mosviz_photometry.ipynb index a30c6386d0..8ef3dd93f6 100644 --- a/notebooks/concepts/generate_mosviz_photometry.ipynb +++ b/notebooks/concepts/generate_mosviz_photometry.ipynb @@ -2,12 +2,12 @@ "cells": [ { "cell_type": "markdown", - "id": "01609adc-98ac-41d2-8030-8eb545ca2fb0", + "id": "9263c0c7-e830-4609-8217-0942a7f373fc", "metadata": {}, "source": [ - "# Synthetic image creation for MOSviz pipeline data\n", + "# Synthetic image creation for Mosviz pipeline data\n", "\n", - "**Motiviation**: The synthetic dataset we currently possess from the JWST data pipeline team for use in MOSviz contains simulated spectra from NIRSpec but no associated NIRCam photometry. We would like to have test imagery to display in MOSviz alongside these 2D and 1D spectra. We are unsure whether the pipeline team has plans to produce any.\n", + "**Motiviation**: The synthetic dataset we currently possess from the JWST data pipeline team for use in Mosviz contains simulated spectra from NIRSpec but no associated NIRCam photometry. We would like to have test imagery to display in Mosviz alongside these 2D and 1D spectra. We are unsure whether the pipeline team has plans to produce any.\n", "\n", "**Goal**: Populate an image of background noise with properly-scaled galaxy cutouts sourced from a Hubble Space Telescope image and placed at their analogous locations in the new image. These galaxies' real spectra do not necessarily correspond with those in our dataset, but we care more about the veneer of having photometry to match with our spectra at this point.\n", "\n", @@ -22,17 +22,19 @@ }, { "cell_type": "markdown", - "id": "57df882f-e15d-4f16-931c-9f495893c0e1", + "id": "bfede422-125f-4347-9f62-67d2dd934041", "metadata": {}, "source": [ "### Import packages\n", "\n", + "- `astropy.convolution.Gaussian2DKernel` smooths two-dimensional data by convolving it with two Gaussian functions.\n", "- We use `astropy.io.fits` to read in existing FITS files and write a new one with the synthetic image.\n", - "- The objects from `astropy.nddata` help with creating cutouts once we've identified galaxies we'd like to take from the field image.\n", - "- The methods from `astropy.stats` work with image data that's clipped to within a certain number of deviations from the mean.\n", + "- The objects from `astropy.nddata` assist in creating cutouts once we have identified galaxies we'd like to take from the field image.\n", + "- The `sigma_*` methods from `astropy.stats` work with image data that's clipped to within a certain number of deviations from the mean. `gaussian_fwhm_to_sigma` is a float for converting full width at half maximum values to standard deviations from the image mean.\n", "- The objects from `astropy.table` help with reading an modifiying tabular data.\n", - "- `astropy.wcs.WCS` creates a World Coordinate System object that's useful for transforming on-sky data from one patch to another.\n", + "- `astropy.wcs.WCS` creates a World Coordinate System object that is useful for transforming on-sky data from one patch to another.\n", "- `glob.glob` lists local files that match a given pattern.\n", + "- The objects from `photutils.segmentation` do the work of analytically finding bright sources, separating those that overlap, and creating a catalog of the resulting information.\n", "- We use `matplotlib.pyplot` to preview the field image, the cutouts, and various stages of our synthetic image.\n", "- We use `numpy` to facilitate several specialized mathematical and array-related operations." ] @@ -40,16 +42,22 @@ { "cell_type": "code", "execution_count": null, - "id": "aba137f6-9bfa-4371-b927-21f30bf9b7fd", - "metadata": {}, + "id": "d4c4daca-ca20-45a0-b42c-212f9dde03b5", + "metadata": { + "tags": [] + }, "outputs": [], "source": [ + "from astropy.convolution import Gaussian2DKernel\n", "from astropy.io import fits\n", "from astropy.nddata import block_reduce, Cutout2D\n", - "from astropy.stats import sigma_clipped_stats, sigma_clip\n", + "from astropy.stats import (gaussian_fwhm_to_sigma,\n", + " sigma_clipped_stats, sigma_clip)\n", "from astropy.table import Table, join\n", "from astropy.wcs import WCS\n", "from glob import glob\n", + "from photutils.segmentation import (detect_sources, deblend_sources,\n", + " SourceCatalog)\n", "\n", "import matplotlib.pyplot as plt\n", "import numpy as np" @@ -57,22 +65,21 @@ }, { "cell_type": "markdown", - "id": "fffd44de-ee89-4a46-8c75-198a94cf124d", + "id": "7f13f717-4236-49e0-b229-b0a1c04d7ad4", "metadata": {}, "source": [ - "### Generate galaxy cutouts\n", + "### Download an HST galaxy field image\n", "\n", - "The galaxy cutouts come from the Hubble Deep Field image. To retrieve them, we download the image itself and its associated catalog data, search the catalog for the brightest galaxies, then find those galaxies in the image." + "The galaxy cutouts come from the Hubble Deep Field image and save its header and flux data to separate variables." ] }, { "cell_type": "code", "execution_count": null, - "id": "e172b883-f625-4a41-b8a5-a2d13d41b396", + "id": "5ea5fefe-deda-4ef5-ad08-32316e45ac41", "metadata": {}, "outputs": [], "source": [ - "# download the image containing sources to be cut out later\n", "image_fits = fits.open('https://archive.stsci.edu/pub/hlsp/hdf/v2/mosaics/x4096/f814_mosaic_v2.fits')\n", "image_header = image_fits[0].header\n", "image_data = image_fits[0].data\n", @@ -83,145 +90,329 @@ { "cell_type": "code", "execution_count": null, - "id": "354bce9f-dc43-434b-a992-b6a636eda576", + "id": "a26b5846-4280-40e9-abb1-99372e2c7f49", + "metadata": {}, + "outputs": [], + "source": [ + "imshow_params = {'cmap': 'bone', 'origin': 'lower'}\n", + "plt.imshow(image_data, vmin=0, vmax=image_data.mean()*3, **imshow_params)" + ] + }, + { + "cell_type": "markdown", + "id": "1936e3dd-31b6-44fb-b84c-460307b6a808", + "metadata": {}, + "source": [ + "We also save image statistics calculated from pixels within a chosen number of standard deviations from the image's mean intensity. Some of them are useful in creating the synthetic image later on." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "6892fa59-cdcf-453b-9ce2-ead0675d7afb", + "metadata": {}, + "outputs": [], + "source": [ + "clipped_mean, clipped_median, clipped_stddev = sigma_clipped_stats(image_data,\n", + " sigma=3.)" + ] + }, + { + "cell_type": "markdown", + "id": "fb0f606e-93fa-44e1-a1e1-969ade968f53", + "metadata": {}, + "source": [ + "### Generate analytical galaxy cutouts\n", + "\n", + "Whether or not the image comes with a catalog, we can use the methods imported from `photutils` to make our own `SourceCatalog`. We follow a workflow modified from `photuils`' [documentation on segmentation](https://photutils.readthedocs.io/en/stable/segmentation.html)." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "8614cc7f-442f-414b-aa7e-084b68bd72c6", + "metadata": { + "scrolled": true, + "tags": [] + }, + "outputs": [], + "source": [ + "npixels = 10 # minimum source length in pixels\n", + "sigma = 3 * gaussian_fwhm_to_sigma # from a FWHM of 3\n", + "\n", + "# define a kernel for smoothing image data in the cutouts\n", + "kernel = Gaussian2DKernel(sigma, x_size=3, y_size=3)\n", + "kernel.normalize()\n", + "#kernel = None\n", + "\n", + "# take first pass at identifying sources\n", + "segments = detect_sources(data=image_data, threshold=8*clipped_stddev,\n", + " npixels=npixels**2, filter_kernel=kernel, mask=None)\n", + "\n", + "# filter out any sources that are too close to the image border\n", + "segments.remove_border_labels(200)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "fa737746-d202-49ff-bdb0-3156dc3dab68", + "metadata": { + "scrolled": true, + "tags": [] + }, + "outputs": [], + "source": [ + "# separate overlapping sources into distinct entries\n", + "segm_deblend = deblend_sources(image_data, segments, npixels=npixels**2,\n", + " filter_kernel=kernel, nlevels=32,\n", + " contrast=1e-1)\n", + "\n", + "# create an astropy Table of source information, sorted by area\n", + "sources = SourceCatalog(image_data, segm_deblend).to_table()\n", + "sources.sort('area', reverse=True)" + ] + }, + { + "cell_type": "markdown", + "id": "b15e29de-269f-423d-b4e8-43c283d29766", + "metadata": {}, + "source": [ + "We can also add a column to `sources` that contains RA/Dec information by creating a `WCS` object and using the field image's header in converting pixel locations to coordinates." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "0cc50110-95e0-4d6d-98e9-c27694303126", + "metadata": {}, + "outputs": [], + "source": [ + "image_wcs = WCS(image_header)\n", + "\n", + "sources.add_columns(\n", + " image_wcs.pixel_to_world_values(sources['xcentroid'], sources['ycentroid']),\n", + " names=['RA', 'Dec'])" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "45af517f-360f-46a7-9305-a9ee5984f994", + "metadata": {}, + "outputs": [], + "source": [ + "sources[:5]" + ] + }, + { + "cell_type": "markdown", + "id": "5f2ade07-e354-499d-b27d-fd383e3353e3", + "metadata": {}, + "source": [ + "At this point, we are ready to ready to create and save cutouts using the source locations found earlier.\n", + "\n", + "_Note: Depending on the value of `catalog_size`, the following cell can produce a lot of output. Right-click the cell and select \"Enable Scrolling for Outputs\" to expand it or \"Disable Scrolling for Outputs\" to condense it._" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "92a73a1b-f834-4d4c-b566-fc5238945d16", + "metadata": { + "scrolled": true, + "tags": [] + }, + "outputs": [], + "source": [ + "# define parameters for cutout loop\n", + "all_cutouts = []\n", + "catalog_size = 20 # how many sources to include in the cutout list\n", + "cutout_length = 100 # the maximum length of a cutout in pixels\n", + "\n", + "for i in sources['label']:\n", + " # save and plot the new cutout\n", + " segm_source = segm_deblend.segments[i - 1]\n", + " # (sources[\"label\"] is 1-indexed, so subtract 1 for matching segm index)\n", + "\n", + " cutout = segm_source.make_cutout(image_data, masked_array=True)\n", + " plt.imshow(cutout, vmin=0, vmax=image_data.std()*3, **imshow_params)\n", + " plt.show()\n", + " \n", + " all_cutouts.append(cutout)\n", + " if len(all_cutouts) == catalog_size:\n", + " break\n", + " \n", + "# indices 9 and 11 for HDF are stars..." + ] + }, + { + "cell_type": "markdown", + "id": "f3cf9ca9-b443-4a10-8d18-b7210b2bfe3d", + "metadata": {}, + "source": [ + "___" + ] + }, + { + "cell_type": "markdown", + "id": "ace709b6-aba7-4ec1-8733-99acdebc4dda", + "metadata": {}, + "source": [ + "___" + ] + }, + { + "cell_type": "markdown", + "id": "87691179-7fd9-497c-9260-8f24c8b31d41", + "metadata": {}, + "source": [ + "### _[Generate galaxy cutouts from catalog file(s)]_\n", + "\n", + "It is also possible to build a list of cutouts from an already-existing source catalog; this section takes a brief excursion to show how it is done.\n", + "\n", + "Depending on the catalog, doing it this way can give you more detailed information about the resulting sources, but the process is less generalizable from one telescope's/instrument's/detector's images to another's. This path may also be less desirable for synthetic image creation because the cutouts will be square patches and thus look a little more, well, synthetic." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "b2f3f95b-2d95-413d-9080-79651e8f03ac", "metadata": {}, "outputs": [], "source": [ "# download sources' location and flux information\n", - "source_info1 = Table.read('https://archive.stsci.edu/pub/hlsp/hdf/wfpc_hdfn_v2catalog/HDFN_wfpc_v2generic.cat',\n", - " format='ascii')\n", - "source_info2 = Table.read('https://archive.stsci.edu/pub/hlsp/hdf/wfpc_hdfn_v2catalog/HDFN_f814_v2.cat',\n", - " format='ascii')\n", - "sources = join(source_info1, source_info2)\n", + "_source_info1 = Table.read('https://archive.stsci.edu/pub/hlsp/hdf/wfpc_hdfn_v2catalog/HDFN_wfpc_v2generic.cat',\n", + " format='ascii')\n", + "_source_info2 = Table.read('https://archive.stsci.edu/pub/hlsp/hdf/wfpc_hdfn_v2catalog/HDFN_f814_v2.cat',\n", + " format='ascii')\n", + "_sources = join(_source_info1, _source_info2)\n", "\n", "# confirm that both tables contain the same objects in the same order (True)\n", - "( (source_info1['NUMBER'] == source_info2['NUMBER']).sum()\n", - " == len(source_info1)\n", - " == len(source_info2) )" + "( (_source_info1['NUMBER'] == _source_info2['NUMBER']).sum()\n", + " == len(_source_info1)\n", + " == len(_source_info2) )" ] }, { "cell_type": "code", "execution_count": null, - "id": "dedfe0f0-b1e4-49e7-b51a-bdf7f069fa06", + "id": "d1fab869-ffa6-4ffb-8fa6-4199cc6e929d", "metadata": {}, "outputs": [], "source": [ "# sort sources by flux within 71.1 pixel diameter of source, or aperture 11\n", - "sources.sort('FLUX_APER_11', reverse=True)\n", + "_sources.sort('FLUX_APER_11', reverse=True)\n", "\n", "# filter out likely stars and sources with negative flux\n", - "sources = sources[(sources['CLASS_STAR'] < .5)\n", - " & (sources['FLUX_APER_8'] > 0)]" + "_sources = _sources[(_sources['CLASS_STAR'] < .5)\n", + " & (_sources['FLUX_APER_8'] > 0)]" ] }, { "cell_type": "code", "execution_count": null, - "id": "3efabf0f-ddbb-4900-9ba8-528998b4bb69", + "id": "2413bcaf-c66a-4040-ade9-a04de1495e6b", "metadata": {}, "outputs": [], "source": [ - "sources[:5]" + "_sources[:5]" ] }, { "cell_type": "code", "execution_count": null, - "id": "f5dcc4be-4e10-45c0-a084-4d224abd6ba1", + "id": "01e5a612-e0f2-48fc-8465-253bd1cc3a9a", "metadata": {}, "outputs": [], "source": [ "# convert the sources' WCS locations to in-image pixel values\n", - "image_wcs = WCS(image_fits[0].header)\n", - "sources_x, sources_y = image_wcs.world_to_pixel_values(sources['ALPHA_J2000'],\n", - " sources['DELTA_J2000'])" + "_image_wcs = WCS(image_header)\n", + "_sources_x, _sources_y = image_wcs.world_to_pixel_values(_sources['ALPHA_J2000'],\n", + " _sources['DELTA_J2000'])" ] }, { "cell_type": "markdown", - "id": "6b4f1421-e3c3-4112-a475-dc1b7ad07146", + "id": "4c5e6b73-584c-4bc3-9775-1636190ff1e3", "metadata": {}, "source": [ - "Note: Depending on the value of `catalog_size`, the following cell can produce a lot of output. Right-click the cell and select \"Enable Scrolling for Outputs\" to expand it or \"Disable Scrolling for Outputs\" to condense it." + "_Note: Depending on the value of `_catalog_size`, the following cell can produce a lot of output. Right-click the cell and select \"Enable Scrolling for Outputs\" to expand it or \"Disable Scrolling for Outputs\" to condense it._" ] }, { "cell_type": "code", "execution_count": null, - "id": "34f1f230-88da-44c3-803d-bb405fe6f6e7", + "id": "1764e645-d1ba-4933-b27b-cbe7422e7116", "metadata": { "scrolled": true, "tags": [] }, "outputs": [], "source": [ - "# save a list of good cutouts for later use\n", - "cutout_list = []\n", - "first_source = 0\n", - "catalog_size = 20\n", - "downsample_factor = 2\n", - "patch_length = 100\n", + "# define parameters for cutout loop\n", + "_cutout_list = []\n", + "_first_source = 0 # which source to cut out first (in descending order by flux)\n", + "_catalog_size = 20 # how many sources to include in the cutout list\n", + "_patch_length = 100 # the length of a cutout in pixels\n", + "_downsample_factor = 2\n", "\n", - "for x, y in list(zip(sources_x, sources_y))[first_source:]:\n", + "# save a list of good cutouts for later use\n", + "for x, y in list(zip(_sources_x, _sources_y))[_first_source:]:\n", " # use pixel locations to cut a source from the image\n", - " cutout = Cutout2D(image_data, (x, y),\n", - " patch_length * downsample_factor).data\n", + " _cutout = Cutout2D(image_data, (x, y),\n", + " _patch_length * _downsample_factor).data\n", " \n", " # bin by downsample_factor to increase field of view\n", - " cutout = block_reduce(cutout, downsample_factor)\n", + " cutout = block_reduce(_cutout, _downsample_factor)\n", " \n", " # skip any cutouts that extend past the image border\n", - " if ( np.all(cutout[-1] <= 0) or np.all(cutout[0] <= 0)\n", - " or np.all(cutout[:,-1] <= 0) or np.all(cutout[:,0] <= 0) ):\n", + " if ( np.all(_cutout[-1] <= 0) or np.all(_cutout[0] <= 0)\n", + " or np.all(_cutout[:,-1] <= 0) or np.all(_cutout[:,0] <= 0) ):\n", " continue\n", " \n", " # save and plot the new cutout\n", - " cutout_list.append(cutout)\n", + " _cutout_list.append(_cutout)\n", " \n", - " plt.imshow(cutout, vmin=-1e-5, vmax=image_data.std(),\n", - " origin='lower', cmap='bone')\n", + " plt.imshow(_cutout, vmin=-1e-5, vmax=image_data.std(), **imshow_params)\n", " plt.show()\n", " \n", - " if len(cutout_list) == catalog_size:\n", + " if len(_cutout_list) == _catalog_size:\n", " break" ] }, { "cell_type": "markdown", - "id": "c3504d50-99cf-4d73-9678-1eff2ca43676", + "id": "775f3dc2-2a42-43c4-967c-44c5f4d7cb6d", "metadata": {}, "source": [ - "We also save image statistics calculated from pixels within a chosen number of standard deviations from the image's mean intensity. Some of them may be useful in creating the synthetic image later on." + "___" ] }, { - "cell_type": "code", - "execution_count": null, - "id": "5b0a6505-3bdd-4f8d-8c2c-e2b7b489d078", + "cell_type": "markdown", + "id": "5a6b90a4-e2db-4994-9436-580bd72245b2", "metadata": {}, - "outputs": [], "source": [ - "clipped_mean, clipped_median, clipped_stddev = sigma_clipped_stats(image_data,\n", - " sigma=3.)" + "___" ] }, { "cell_type": "markdown", - "id": "53891064-621c-4f6f-91da-3a24741ad632", + "id": "746288d9-3430-4032-898e-156cd0e39d36", "metadata": {}, "source": [ "### Extract destination RA/Dec from spectra files\n", "\n", "Note that the following cells only run if on a connection with access to the STScI VPN. They use copies of JWST pipeline files from May 19, 2020.\n", "\n", - "To run them elsewhere and/or use more current files, visit the [Artifactory](https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth%2Ftest_nirspec_mos_spec3), ensure the `jwst-pipeline/dev/truth/test_nirspec_mos_spec3/` folder is selected in the directory tree, and download all files whose names begin with `jw00626-o030` and end with `nirspec_f170lp-g235m_x1d.fits`. (The original version of the notebook uses six of these Level 3 data products.) Then, change `filepath` in the next cell to the location where you saved the downloaded files. Remember the trailing forward slash." + "To run them elsewhere and/or use more current files, visit the [Artifactory](https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth%2Ftest_nirspec_mos_spec3), ensure the `jwst-pipeline/dev/truth/test_nirspec_mos_spec3/` folder is the one selected in the directory tree, and download all files whose names begin with `jw00626-o030` and end with `nirspec_f170lp-g235m_x1d.fits`. (The original version of the notebook uses six of these Level 3 data products.) Then, change `filepath` in the next cell to the location where you saved the downloaded files. Remember the trailing forward slash." ] }, { "cell_type": "code", "execution_count": null, - "id": "9436aa5d-54c9-4a1e-a983-7f8658920be6", + "id": "3c54fde5-ad04-43f1-8acc-ccc33dc2c5b1", "metadata": {}, "outputs": [], "source": [ @@ -234,7 +425,7 @@ { "cell_type": "code", "execution_count": null, - "id": "96dda8bd-5010-45ef-a937-b99ee54b8091", + "id": "63908a5c-83d1-4cba-a516-cd0a8b28f76c", "metadata": {}, "outputs": [], "source": [ @@ -244,7 +435,7 @@ { "cell_type": "code", "execution_count": null, - "id": "ade024d2-9788-4572-9304-af69c30e86da", + "id": "e94df0cd-dbf0-4fcf-96f0-e3f3d446229b", "metadata": {}, "outputs": [], "source": [ @@ -253,7 +444,7 @@ }, { "cell_type": "markdown", - "id": "809be77d-b099-40a3-8a85-42a7a2cf3244", + "id": "a225379a-d1b7-403c-9421-02dbf3386ced", "metadata": {}, "source": [ "Notice that this header lacks WCS information. Additionally, examining all of this observation's file headers reveals that they all have the same RA/Dec, which is not what we expect for its different \"pointings.\"" @@ -262,13 +453,12 @@ { "cell_type": "code", "execution_count": null, - "id": "879040fb-eb58-4e4c-8173-a289b21a61e2", + "id": "15610ff2-8103-4b79-9957-d537c126e7ec", "metadata": {}, "outputs": [], "source": [ "# search for RA/Dec information from Artifactory observation files\n", - "x1d_header_list = [fits.getheader(file)\n", - " for file in glob(filepath + 'jw00626*x1d.fits')]\n", + "x1d_header_list = [fits.getheader(file) for file in glob(filepath + 'jw00626*x1d.fits')]\n", "\n", "ras, decs = np.array([[h['TARG_RA'], h['TARG_DEC']] for h in x1d_header_list]).T\n", "ras, decs" @@ -276,7 +466,7 @@ }, { "cell_type": "markdown", - "id": "debb3ccc-88ec-4e2d-bf50-69e8361f7d4b", + "id": "e238ac42-f7ea-408d-a2e5-7199ba315989", "metadata": {}, "source": [ "Since the headers lack the information we seek, we randomly generate our sources' RA/Dec information in predetermined patch of sky. We take the patch size from the size of NIRSpec Micro-Shutter Assembly (MSA) -- about 3.6'x3.4'-- approximated to a square field of view." @@ -285,10 +475,11 @@ { "cell_type": "code", "execution_count": null, - "id": "9931cde4-5f77-4a40-96e5-32d252bec47c", + "id": "82a845e5-8c3d-4b97-b3a9-e265b9d8d63d", "metadata": {}, "outputs": [], "source": [ + "# (ranges based on NIRSpec MSA's on-sky projection size of 3.6x3.4 arcmins)\n", "np.random.seed(19)\n", "ras = np.random.uniform(0, 1/15, catalog_size)\n", "decs = np.random.uniform(-1/30, 1/30, catalog_size)" @@ -296,7 +487,7 @@ }, { "cell_type": "markdown", - "id": "c4b5da62-c424-4603-adc3-4c550291d13a", + "id": "c9fd855a-f16e-4ee3-bd80-bb68023238c3", "metadata": {}, "source": [ "### Create synthetic image\n", @@ -307,10 +498,13 @@ { "cell_type": "code", "execution_count": null, - "id": "6d1f9e28-a5fa-4091-bea6-72c88b48e564", - "metadata": {}, + "id": "ba5ff697-c9b0-4218-a1d0-95d2fc8ff8e9", + "metadata": { + "tags": [] + }, "outputs": [], "source": [ + "# create synthetic image onto which cutouts will be pasted\n", "synth_img_size = 1000\n", "synth_image = np.zeros((synth_img_size, synth_img_size))" ] @@ -318,7 +512,7 @@ { "cell_type": "code", "execution_count": null, - "id": "cb856895-9b87-44c8-8c4b-2cbab45edb76", + "id": "7ceac16a-8bfc-4e8e-8713-2dfe6251f60e", "metadata": {}, "outputs": [], "source": [ @@ -330,39 +524,37 @@ { "cell_type": "code", "execution_count": null, - "id": "3b52414b-453c-43eb-984b-e856c284cb95", + "id": "5a40c665-935a-40f4-8a85-15fd943285a0", "metadata": {}, "outputs": [], "source": [ - "imshow_params = {'cmap': 'bone', 'origin': 'lower'}\n", "plt.imshow(synth_image, **imshow_params)\n", "plt.show()" ] }, { "cell_type": "markdown", - "id": "674017bb-2cc1-437a-bf3a-449f306fa506", + "id": "502d2402-572c-4892-81a4-692e59dbb208", "metadata": {}, "source": [ "### Fill out new WCS object for `synth_image`\n", "\n", - "Creating a `WCS` object for `synth_image` allows it to be mathematically transformed into a projection on the sky. That projection can then be compared to other FITS images with their own `WCS` information." + "Creating a `WCS` object for `synth_image` allows it to be mathematically transformed into a projection on the sky. That projection can then be compared to other FITS images with their own WCS information." ] }, { "cell_type": "code", "execution_count": null, - "id": "1cd4ee37-d41f-4ddd-960c-447a72b1faa2", + "id": "0f22d48b-67ba-4cba-9f41-2d0d24bba9e2", "metadata": {}, "outputs": [], "source": [ - "synth_wcs = WCS(naxis=2)\n", - "synth_wcs" + "synth_wcs = WCS(naxis=2)" ] }, { "cell_type": "markdown", - "id": "c1a2876e-963b-4e62-8dd1-d68fb4b4e0a7", + "id": "f8609331-a34e-4b60-a227-2bae558cb2f3", "metadata": {}, "source": [ "The next step is to calculate field of view information for `synth_image`." @@ -371,8 +563,10 @@ { "cell_type": "code", "execution_count": null, - "id": "068af407-9a8a-429c-9f7e-03f4f0a0b867", - "metadata": {}, + "id": "c5eaf810-6ed5-4471-9fed-8c8eed8e6f2d", + "metadata": { + "tags": [] + }, "outputs": [], "source": [ "# find the range of sources in RA and dec\n", @@ -386,7 +580,7 @@ { "cell_type": "code", "execution_count": null, - "id": "d70d1d2a-8c9c-40b7-b31f-3a72377279df", + "id": "adee6b70-b098-4efe-b551-0bbd212c3e36", "metadata": {}, "outputs": [], "source": [ @@ -402,7 +596,7 @@ { "cell_type": "code", "execution_count": null, - "id": "308eb86b-5a67-4fbe-bf31-04a08642537a", + "id": "46aabc3b-7368-40db-9c2b-9e92343a60fd", "metadata": {}, "outputs": [], "source": [ @@ -416,7 +610,7 @@ }, { "cell_type": "markdown", - "id": "d28eb2e7-d3a1-4190-978a-8be1b30d3a3c", + "id": "3366b9b1-1944-42be-828f-bc609a3aeb6b", "metadata": {}, "source": [ "With those calculations done, the `WCS` object is ready to be filled out." @@ -425,7 +619,7 @@ { "cell_type": "code", "execution_count": null, - "id": "09dadad8-902d-4119-a317-c10064968597", + "id": "85ffc3e6-0688-4a88-b3bc-73ac42ef8fae", "metadata": { "tags": [] }, @@ -443,20 +637,10 @@ "synth_wcs" ] }, - { - "cell_type": "markdown", - "id": "5792df4f-c0d6-49a5-be72-a7c5ab39d9b4", - "metadata": {}, - "source": [ - "### Populate `synth_image` with the cutouts\n", - "\n", - "Finally, we convert the cutouts' coordinates to pixels and add them into `synth_image` to complete the creation of the mock image." - ] - }, { "cell_type": "code", "execution_count": null, - "id": "34c92d83-8e78-4cc5-8df3-728cb51d6567", + "id": "a50df507-e367-46c1-911b-bc3b1bff7e6d", "metadata": {}, "outputs": [], "source": [ @@ -465,36 +649,52 @@ "ras_pix, decs_pix" ] }, + { + "cell_type": "markdown", + "id": "b3e474ba-5a5d-431f-95d3-6b99e621ef7b", + "metadata": {}, + "source": [ + "### Populate `synth_image` with the cutouts\n", + "\n", + "Finally, we add the cutouts into `synth_image` to complete its creation." + ] + }, { "cell_type": "code", "execution_count": null, - "id": "c587e262-15fc-4fec-b5d3-b9d173afc286", + "id": "4162eedc-5388-48f9-a637-d3b5447508c4", "metadata": {}, "outputs": [], "source": [ - "cutout_half_delta = patch_length // 2\n", - "\n", - "for i in range(len(ras_pix)):\n", - " synth_image[ras_pix[i] - cutout_half_delta : ras_pix[i] + cutout_half_delta,\n", - " decs_pix[i] - cutout_half_delta : decs_pix[i] + cutout_half_delta] += cutout_list[i]" + "for i, co in enumerate(all_cutouts):\n", + " # first, calculate cutout's distances from center to each edge, accounting\n", + " # for odd numbers by using np.floor to shift that axis' values by 0.5\n", + " cutout_ranges = np.array([np.floor([-n/2, n/2]) for n in co.shape],\n", + " dtype='int')\n", + " \n", + " synth_image[\n", + " ras_pix[i] + cutout_ranges[0,0] : ras_pix[i] + cutout_ranges[0,1],\n", + " decs_pix[i] + cutout_ranges[1,0] : decs_pix[i] + cutout_ranges[1,1]\n", + " ] += co" ] }, { "cell_type": "code", "execution_count": null, - "id": "9bccee92-2930-4116-b599-0269fdf1d2dc", + "id": "24ae2b48-ef41-4e45-9412-85dc091c70d3", "metadata": {}, "outputs": [], "source": [ "fig, ax = plt.subplots(figsize=(10,10))\n", - "ax.imshow(synth_image, vmin=0, vmax=synth_image.std()*3, **imshow_params)\n", - "\n", + "ax.imshow(synth_image, vmin=0, vmax=synth_image.std(), **imshow_params)\n", + "#ax.imshow(synth_image, vmin=0, vmax=synth_image.mean()*3, **imshow_params)\n", + "ax.tick_params(labelsize=14)\n", "plt.show()" ] }, { "cell_type": "markdown", - "id": "eb2a0045-74b6-4ac2-85a8-c535e1a64800", + "id": "8b518b39-3d22-42c4-a632-4c17801900a5", "metadata": {}, "source": [ "Save the image, overwriting the previous file if it already exists in the current directory." @@ -503,7 +703,7 @@ { "cell_type": "code", "execution_count": null, - "id": "337dc041-618e-4f5d-bcb2-845cacbfac31", + "id": "9b4f35a3-45a8-492e-98f8-039d5831a29f", "metadata": {}, "outputs": [], "source": [ @@ -513,11 +713,11 @@ }, { "cell_type": "markdown", - "id": "fe8817e2-51f4-4417-84ce-659c9449e87a", + "id": "608047a0-706e-4437-a0b4-8f699cac8f82", "metadata": {}, "source": [ "

\n", - " Authors: Robel Geda and O. Justin Otor \n", + " Authors: O. Justin Otor and Robel Geda \n", " \"Space\n", "

\n", "\n", From 7b992543bbd25c9d829d8fe2d839d92eda5d3034 Mon Sep 17 00:00:00 2001 From: ojustino Date: Tue, 8 Jun 2021 14:01:36 -0400 Subject: [PATCH 4/7] Added steps to make cutouts from synthetic image --- .../concepts/generate_mosviz_photometry.ipynb | 180 +++++++++++++----- 1 file changed, 135 insertions(+), 45 deletions(-) diff --git a/notebooks/concepts/generate_mosviz_photometry.ipynb b/notebooks/concepts/generate_mosviz_photometry.ipynb index 8ef3dd93f6..71a39bbcc2 100644 --- a/notebooks/concepts/generate_mosviz_photometry.ipynb +++ b/notebooks/concepts/generate_mosviz_photometry.ipynb @@ -5,18 +5,18 @@ "id": "9263c0c7-e830-4609-8217-0942a7f373fc", "metadata": {}, "source": [ - "# Synthetic image creation for Mosviz pipeline data\n", + "# Synthetic image and cutout creation for Mosviz pipeline data\n", "\n", - "**Motiviation**: The synthetic dataset we currently possess from the JWST data pipeline team for use in Mosviz contains simulated spectra from NIRSpec but no associated NIRCam photometry. We would like to have test imagery to display in Mosviz alongside these 2D and 1D spectra. We are unsure whether the pipeline team has plans to produce any.\n", + "**Motiviation**: The synthetic dataset we currently possess from the JWST data pipeline team for use in Mosviz contains simulated spectra from NIRSpec but no associated NIRCam photometry. We would like to have test imagery to display in Mosviz alongside these 2D and 1D spectra since we are not aware that the pipeline team has plans to produce any.\n", "\n", - "**Goal**: Populate an image of background noise with properly-scaled galaxy cutouts sourced from a Hubble Space Telescope image and placed at their analogous locations in the new image. These galaxies' real spectra do not necessarily correspond with those in our dataset, but we care more about the veneer of having photometry to match with our spectra at this point.\n", + "**Goal**: Populate an image of background noise with properly-scaled galaxy cutouts sourced from a Hubble Space Telescope image and placed at their analogous locations in the new image. These galaxies' real spectra do not necessarily correspond with those in our dataset, but we care more about the veneer of having photometry to match with our spectra in Mosviz at this point.\n", "\n", - "**Execution**: We pull our galaxy cutouts and catalog information from the Hubble Deep Field image. ~~ASTRODEEP's image of the [Abell 2774 Parallel](http://astrodeep.u-strasbg.fr/ff/?img=JH140?cm=grayscale) | [MACS J0416.1-2403 Parallel](http://astrodeep.u-strasbg.fr/ff/?ffid=FF_M0416PAR&id=1264&cm=grayscale)~~. We sought to use [Artifactory](https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth) to obtain a range of RA/Dec over which to project our synthetic image, but that information was absent. Instead, we place the image over a manually chosen RA/Dec range and place the cutouts in randomly selected locations within that field of view.\n", + "**Execution**: We pull our galaxy cutouts and catalog information from the Hubble Deep Field image. _(Formerly used ASTRODEEP's image of the [Abell 2774 Parallel](http://astrodeep.u-strasbg.fr/ff/?img=JH140?cm=grayscale) and [MACS J0416.1-2403 Parallel](http://astrodeep.u-strasbg.fr/ff/?ffid=FF_M0416PAR&id=1264&cm=grayscale).)_ We use `photutils` to select bright sources from the image for use in our own. We sought to use [Artifactory](https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth) to obtain a range of RA/Dec over which to project our synthetic image, but that information was absent. Instead, we place the image over a manually chosen RA/Dec range and place the cutouts in randomly selected locations within that field of view.\n", "\n", "**Issues**:\n", - "- We wanted to scrape RA/Dec information the data pipeline products to get a range of coordinates over which to scale our synthetic image, but it appears that the pipeline's data products lack `\"TARG_RA\"` or `\"TARG_DEC\"` keywords in their headers.\n", - " - The data products also don't appear to have WCS information. We don't strictly need it to achieve this notebook's goals, but it would be convenient to have.\n", - "- There does not appear to be an observation with level 2 data in `jwst-pipeline/truth/test_nirspec_mos_spec2` and level 3 data in `jwst-pipeline/truth/test_nirspec_mos_spec3`. All observations are either level 2 only or level 3 only.\n", + "- We wanted to scrape RA/Dec information from the data pipeline products to get a range of coordinates over which to scale our synthetic image, but it appears that the pipeline's data products lack `\"TARG_RA\"` or `\"TARG_DEC\"` keywords in their headers.\n", + " - The data products also do not appear to have WCS information. Though it is not strictly needed to achieve this notebook's goals, it would be convenient to have.\n", + "- There does not appear to be an observation with Level 2 data in `jwst-pipeline/truth/test_nirspec_mos_spec2` and Level 3 data in `jwst-pipeline/truth/test_nirspec_mos_spec3`. All observations are either Level 2 only or Level 3 only.\n", "- _(Resolved)_ A good number of the cutouts from the first couple of field images we tested had intrustions from other galaxies due to crowding. We settled on the Hubble Deep Field as a good source image, but had we not, we may have considered using galaxies modeled with Sersic profiles to get cleaner cutouts to inject into our synthetic image." ] }, @@ -27,10 +27,12 @@ "source": [ "### Import packages\n", "\n", - "- `astropy.convolution.Gaussian2DKernel` smooths two-dimensional data by convolving it with two Gaussian functions.\n", + "- We use `astropy.units` to convert between degrees and arcseconds.\n", + "- `astropy.coordinates.SkyCoord` helps merge lists of right ascensions (RA) and declinations (Dec) into one object.\n", + "- `astropy.convolution.Gaussian2DKernel` smooths two-dimensional data by performing convolutions with two Gaussian functions.\n", "- We use `astropy.io.fits` to read in existing FITS files and write a new one with the synthetic image.\n", - "- The objects from `astropy.nddata` assist in creating cutouts once we have identified galaxies we'd like to take from the field image.\n", - "- The `sigma_*` methods from `astropy.stats` work with image data that's clipped to within a certain number of deviations from the mean. `gaussian_fwhm_to_sigma` is a float for converting full width at half maximum values to standard deviations from the image mean.\n", + "- The objects from `astropy.nddata` assist in creating cutouts once we have identified galaxies to take from the field image and also with re-importing our synthetic image.\n", + "- The `sigma_*` methods from `astropy.stats` work with image data that is clipped to within a certain number of deviations from the mean. `gaussian_fwhm_to_sigma` is a float for converting full width at half maximum values to standard deviations from the image mean.\n", "- The objects from `astropy.table` help with reading an modifiying tabular data.\n", "- `astropy.wcs.WCS` creates a World Coordinate System object that is useful for transforming on-sky data from one patch to another.\n", "- `glob.glob` lists local files that match a given pattern.\n", @@ -48,9 +50,11 @@ }, "outputs": [], "source": [ + "from astropy import units as u\n", + "from astropy.coordinates import SkyCoord\n", "from astropy.convolution import Gaussian2DKernel\n", "from astropy.io import fits\n", - "from astropy.nddata import block_reduce, Cutout2D\n", + "from astropy.nddata import block_reduce, Cutout2D, CCDData\n", "from astropy.stats import (gaussian_fwhm_to_sigma,\n", " sigma_clipped_stats, sigma_clip)\n", "from astropy.table import Table, join\n", @@ -70,7 +74,7 @@ "source": [ "### Download an HST galaxy field image\n", "\n", - "The galaxy cutouts come from the Hubble Deep Field image and save its header and flux data to separate variables." + "The galaxy cutouts come from the Hubble Deep Field image, whose header and flux data we save as separate variables." ] }, { @@ -212,13 +216,13 @@ "source": [ "At this point, we are ready to ready to create and save cutouts using the source locations found earlier.\n", "\n", - "_Note: Depending on the value of `catalog_size`, the following cell can produce a lot of output. Right-click the cell and select \"Enable Scrolling for Outputs\" to expand it or \"Disable Scrolling for Outputs\" to condense it._" + "_Note: Depending on the value of `catalog_size`, the following cell can produce a lot of output. Right-click the cell and select \"Disable Scrolling for Outputs\" to expand it or \"Enable Scrolling for Outputs\" to condense it._" ] }, { "cell_type": "code", "execution_count": null, - "id": "92a73a1b-f834-4d4c-b566-fc5238945d16", + "id": "edf1b0a4-1132-42e8-b6fd-12436f2355b5", "metadata": { "scrolled": true, "tags": [] @@ -227,23 +231,25 @@ "source": [ "# define parameters for cutout loop\n", "all_cutouts = []\n", - "catalog_size = 20 # how many sources to include in the cutout list\n", + "catalog_size = 33 # how many sources to include in the cutout list\n", "cutout_length = 100 # the maximum length of a cutout in pixels\n", + "downsample_factor = 4 # scale sources down to proper size for the new image\n", "\n", "for i in sources['label']:\n", " # save and plot the new cutout\n", " segm_source = segm_deblend.segments[i - 1]\n", " # (sources[\"label\"] is 1-indexed, so subtract 1 for matching segm index)\n", "\n", + " # implement the downsample\n", " cutout = segm_source.make_cutout(image_data, masked_array=True)\n", + " cutout = block_reduce(cutout, downsample_factor) / downsample_factor**2\n", + " \n", " plt.imshow(cutout, vmin=0, vmax=image_data.std()*3, **imshow_params)\n", " plt.show()\n", " \n", " all_cutouts.append(cutout)\n", " if len(all_cutouts) == catalog_size:\n", - " break\n", - " \n", - "# indices 9 and 11 for HDF are stars..." + " break" ] }, { @@ -267,7 +273,7 @@ "id": "87691179-7fd9-497c-9260-8f24c8b31d41", "metadata": {}, "source": [ - "### _[Generate galaxy cutouts from catalog file(s)]_\n", + "### _[An alternate cutout generation process using catalog file(s)]_\n", "\n", "It is also possible to build a list of cutouts from an already-existing source catalog; this section takes a brief excursion to show how it is done.\n", "\n", @@ -337,7 +343,7 @@ "id": "4c5e6b73-584c-4bc3-9775-1636190ff1e3", "metadata": {}, "source": [ - "_Note: Depending on the value of `_catalog_size`, the following cell can produce a lot of output. Right-click the cell and select \"Enable Scrolling for Outputs\" to expand it or \"Disable Scrolling for Outputs\" to condense it._" + "_Note: Depending on the value of `_catalog_size`, the following cell can produce a lot of output. Right-click the cell and select \"Disable Scrolling for Outputs\" to expand it or \"Enable Scrolling for Outputs\" to condense it._" ] }, { @@ -364,7 +370,7 @@ " _patch_length * _downsample_factor).data\n", " \n", " # bin by downsample_factor to increase field of view\n", - " cutout = block_reduce(_cutout, _downsample_factor)\n", + " _cutout = block_reduce(_cutout, _downsample_factor)\n", " \n", " # skip any cutouts that extend past the image border\n", " if ( np.all(_cutout[-1] <= 0) or np.all(_cutout[0] <= 0)\n", @@ -402,11 +408,11 @@ "id": "746288d9-3430-4032-898e-156cd0e39d36", "metadata": {}, "source": [ - "### Extract destination RA/Dec from spectra files\n", + "### (Try to) Extract destination RA/Dec from spectra files\n", "\n", - "Note that the following cells only run if on a connection with access to the STScI VPN. They use copies of JWST pipeline files from May 19, 2020.\n", + "Note that this section's cells only run on a connection with access to the STScI VPN. They use copies of JWST pipeline files from May 19, 2021.\n", "\n", - "To run them elsewhere and/or use more current files, visit the [Artifactory](https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth%2Ftest_nirspec_mos_spec3), ensure the `jwst-pipeline/dev/truth/test_nirspec_mos_spec3/` folder is the one selected in the directory tree, and download all files whose names begin with `jw00626-o030` and end with `nirspec_f170lp-g235m_x1d.fits`. (The original version of the notebook uses six of these Level 3 data products.) Then, change `filepath` in the next cell to the location where you saved the downloaded files. Remember the trailing forward slash." + "To run them elsewhere and/or use more current files, visit the [Artifactory](https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth%2Ftest_nirspec_mos_spec3), ensure the `jwst-pipeline/dev/truth/test_nirspec_mos_spec3/` folder is the one selected in the directory tree, and download all files whose names begin with `jw00626-o030` and end with `nirspec_f170lp-g235m_x1d.fits`. (The original version of this notebook uses six of these Level 3 data products.) Then, change `filepath` in the next cell to the location where you saved the downloaded files. Remember the trailing forward slash." ] }, { @@ -417,7 +423,7 @@ "outputs": [], "source": [ "# view level 3 spectra FITS header information\n", - "filepath = '/user/jotor/jwst-pipeline-lvl3/'\n", + "filepath = '/user/jotor/jwst-pipeline-lvl3-copy/'\n", "x1d_header = fits.getheader(filepath + 'jw00626-o030_s00000_nirspec_f170lp-g235m_x1d.fits')\n", "#s2d_header = fits.getheader(filepath + 'jw00626-o030_s00000_nirspec_f170lp-g235m_s2d.fits')" ] @@ -447,7 +453,7 @@ "id": "a225379a-d1b7-403c-9421-02dbf3386ced", "metadata": {}, "source": [ - "Notice that this header lacks WCS information. Additionally, examining all of this observation's file headers reveals that they all have the same RA/Dec, which is not what we expect for its different \"pointings.\"" + "Notice that the `WCS` object above contains its default settings; this is because the spectrum's header lacks WCS information. Additionally, examining all of this observation's file headers reveals that they all have the same RA/Dec, which is not what we expect for its different \"pointings.\"" ] }, { @@ -458,24 +464,27 @@ "outputs": [], "source": [ "# search for RA/Dec information from Artifactory observation files\n", - "x1d_header_list = [fits.getheader(file) for file in glob(filepath + 'jw00626*x1d.fits')]\n", + "x1d_header_list = [fits.getheader(file)\n", + " for file in glob(filepath + 'jw00626*x1d.fits')]\n", "\n", - "ras, decs = np.array([[h['TARG_RA'], h['TARG_DEC']] for h in x1d_header_list]).T\n", - "ras, decs" + "_ras, _decs = np.array([[h['TARG_RA'], h['TARG_DEC']] for h in x1d_header_list]).T\n", + "_ras, _decs" ] }, { "cell_type": "markdown", - "id": "e238ac42-f7ea-408d-a2e5-7199ba315989", + "id": "c24cff3e-8455-4fac-9e00-7a3a537bb9ea", "metadata": {}, "source": [ - "Since the headers lack the information we seek, we randomly generate our sources' RA/Dec information in predetermined patch of sky. We take the patch size from the size of NIRSpec Micro-Shutter Assembly (MSA) -- about 3.6'x3.4'-- approximated to a square field of view." + "### Create synthetic image\n", + "\n", + "Since the headers lack the information we seek, we randomly generate our sources' RA/Dec information in a predetermined patch of sky. We choose the patch size by approximating NIRSpec's 3.6'x3.4' Micro-Shutter Assembly (MSA) to a square field of view." ] }, { "cell_type": "code", "execution_count": null, - "id": "82a845e5-8c3d-4b97-b3a9-e265b9d8d63d", + "id": "21585c3f-30d1-401a-9958-d50cc17d4c17", "metadata": {}, "outputs": [], "source": [ @@ -487,12 +496,10 @@ }, { "cell_type": "markdown", - "id": "c9fd855a-f16e-4ee3-bd80-bb68023238c3", + "id": "cf39ba39-0264-4dfd-87c7-1a5df076e251", "metadata": {}, "source": [ - "### Create synthetic image\n", - "\n", - "We initialize a `numpy` array and fill it with normally-distributed background noise based on some of the clipped image statistics that were calculated earlier." + "We initialize a `numpy` array and fill it with a normally-distributed background noise level based on some of the clipped image statistics that were calculated earlier." ] }, { @@ -604,8 +611,8 @@ "pix_scale = min_image_fov / synth_img_size\n", "\n", "# add a buffer to the FOV's borders to prevent clipping sources\n", - "pix_scale *= 1.5\n", - "pix_scale" + "pix_scale *= 1.2\n", + "(pix_scale * u.deg).to(u.arcsecond)" ] }, { @@ -656,7 +663,7 @@ "source": [ "### Populate `synth_image` with the cutouts\n", "\n", - "Finally, we add the cutouts into `synth_image` to complete its creation." + "Now is the time to add the cutouts into `synth_image` and complete its creation." ] }, { @@ -672,9 +679,10 @@ " cutout_ranges = np.array([np.floor([-n/2, n/2]) for n in co.shape],\n", " dtype='int')\n", " \n", + " # remember to flip RA/dec to dec/RA when indexing image array\n", " synth_image[\n", - " ras_pix[i] + cutout_ranges[0,0] : ras_pix[i] + cutout_ranges[0,1],\n", - " decs_pix[i] + cutout_ranges[1,0] : decs_pix[i] + cutout_ranges[1,1]\n", + " decs_pix[i] + cutout_ranges[0,0] : decs_pix[i] + cutout_ranges[0,1],\n", + " ras_pix[i] + cutout_ranges[1,0] : ras_pix[i] + cutout_ranges[1,1]\n", " ] += co" ] }, @@ -697,7 +705,9 @@ "id": "8b518b39-3d22-42c4-a632-4c17801900a5", "metadata": {}, "source": [ - "Save the image, overwriting the previous file if it already exists in the current directory." + "We save the image to a local location. The default is this notebook's current directory; you can change that by adjusting `savepath` to your preferred path and remembering the trailing slash.\n", + "\n", + "This cell will overwrite a same-named file if it already exists in that location." ] }, { @@ -707,15 +717,95 @@ "metadata": {}, "outputs": [], "source": [ - "fits.writeto('synthetic_HDF_more.fits', synth_image,\n", + "savepath = ''\n", + "synth_file = 'synthetic_mosviz_image.fits'\n", + "\n", + "fits.writeto(savepath + synth_file, synth_image,\n", " header=synth_wcs.to_header(), overwrite=True)" ] }, { "cell_type": "markdown", - "id": "608047a0-706e-4437-a0b4-8f699cac8f82", + "id": "12fb131f-a580-4624-99e1-d1f1443cc15c", "metadata": {}, "source": [ + "### Harvest cutouts from synthetic image\n", + "\n", + "Finally, we can create image cutouts from the new synthetic image to associate with the pipeline's synthetic spectral data in [Mosviz](https://jdaviz.readthedocs.io/en/latest/mosviz/index.html). We begin by reopening the synthetic image with `astropy`'s `CCDData` class in order to handle the WCS information in its header more smoothly.\n", + "\n", + "Then, we use `Cutout2D` to create square cutouts of each source in the image by using coordinate information previously taken from the synthetic sprectra files." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "1eb13a31-7342-4876-98f5-8d249c439b38", + "metadata": {}, + "outputs": [], + "source": [ + "ccd_synth_img = CCDData.read(savepath + synth_file, unit='electron/s')" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "07ff36c5-8345-4700-8c23-ef7c4a69eb21", + "metadata": {}, + "outputs": [], + "source": [ + "synth_cutout_size = 50 # pixels\n", + "synth_cutout_list = [Cutout2D(ccd_synth_img, SkyCoord(r, d), synth_cutout_size)\n", + " for r, d in zip(ras * u.deg, decs * u.deg)]\n", + "\n", + "plt.imshow(synth_cutout_list[0].data, **imshow_params,\n", + " vmin=0, vmax=ccd_synth_img.data.std())" + ] + }, + { + "cell_type": "markdown", + "id": "c8ec6814-e350-400b-b314-fdf9ac4e05f2", + "metadata": {}, + "source": [ + "These cutout images are also ready for saving after adding some extra information to their headers. It may be useful to create a new folder and save them there by modifying the `savepath2` variable below.\n", + "\n", + "_Note: Depending on the length of `synth_cutout_list`, the following cell can produce a lot of output. Right-click the cell and select \"Disable Scrolling for Outputs\" to expand it or \"Enable Scrolling for Outputs\" to condense it._" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "cd58fa14-739c-423d-af3b-8cbb11c80ba2", + "metadata": { + "scrolled": true, + "tags": [] + }, + "outputs": [], + "source": [ + "savepath2 = ''\n", + "\n", + "for i, cut in enumerate(synth_cutout_list):\n", + " sdt = cut.data\n", + " shd = cut.wcs.to_header()\n", + " \n", + " # add metadata for proper reading in Mosviz\n", + " shd['OBJ_RA'] = (ras[i], 'Source right ascension')\n", + " shd['OBJ_DEC'] = (decs[i], 'Source declination')\n", + " shd['OBJ_ROT'] = (0, '')\n", + " \n", + " fits.writeto(savepath2 + f\"/user/jotor/mosviz-cutouts-copy/synth_cutout{src}.fits\",\n", + " sdt, header=shd, overwrite=True)\n", + " \n", + " plt.imshow(sdt, **imshow_params, vmin=0, vmax=ccd_synth_img.data.std())\n", + " plt.show()" + ] + }, + { + "cell_type": "markdown", + "id": "6b4c5cb2-f304-471e-b04b-5e2fa6baec17", + "metadata": {}, + "source": [ + "These cutout files and the spectral files we examined earlier can serve as the data behind a workflow in [Mosviz](https://jdaviz.readthedocs.io/en/latest/mosviz/index.html). Please see other notebooks in this directory for more information on how to use that tool.\n", + "\n", "

\n", " Authors: O. Justin Otor and Robel Geda \n", " \"Space\n", @@ -741,7 +831,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.8.8" + "version": "3.8.10" } }, "nbformat": 4, From 381a97f6b5f44f349c17ef3198a5f7d503d56e2e Mon Sep 17 00:00:00 2001 From: ojustino Date: Tue, 15 Jun 2021 17:35:00 -0400 Subject: [PATCH 5/7] Implemented Jupyter magic; fixed saving of cutouts --- .../concepts/generate_mosviz_photometry.ipynb | 42 +++++++++---------- 1 file changed, 20 insertions(+), 22 deletions(-) diff --git a/notebooks/concepts/generate_mosviz_photometry.ipynb b/notebooks/concepts/generate_mosviz_photometry.ipynb index 71a39bbcc2..b4870a267e 100644 --- a/notebooks/concepts/generate_mosviz_photometry.ipynb +++ b/notebooks/concepts/generate_mosviz_photometry.ipynb @@ -38,7 +38,8 @@ "- `glob.glob` lists local files that match a given pattern.\n", "- The objects from `photutils.segmentation` do the work of analytically finding bright sources, separating those that overlap, and creating a catalog of the resulting information.\n", "- We use `matplotlib.pyplot` to preview the field image, the cutouts, and various stages of our synthetic image.\n", - "- We use `numpy` to facilitate several specialized mathematical and array-related operations." + "- We use `numpy` to facilitate several specialized mathematical and array-related operations.\n", + "- `os` helps us create a new directory for the final cutouts from the synthetic image." ] }, { @@ -64,7 +65,8 @@ " SourceCatalog)\n", "\n", "import matplotlib.pyplot as plt\n", - "import numpy as np" + "import numpy as np\n", + "import os" ] }, { @@ -252,14 +254,6 @@ " break" ] }, - { - "cell_type": "markdown", - "id": "f3cf9ca9-b443-4a10-8d18-b7210b2bfe3d", - "metadata": {}, - "source": [ - "___" - ] - }, { "cell_type": "markdown", "id": "ace709b6-aba7-4ec1-8733-99acdebc4dda", @@ -275,7 +269,7 @@ "source": [ "### _[An alternate cutout generation process using catalog file(s)]_\n", "\n", - "It is also possible to build a list of cutouts from an already-existing source catalog; this section takes a brief excursion to show how it is done.\n", + "It is also possible to build a list of cutouts from an already-existing source catalog; this section serves as a reference for how it is done. Its cells are not necessary for the completion of the notebook and do not run due to the `%%script false --no-raise-error` statements at the top of each one.\n", "\n", "Depending on the catalog, doing it this way can give you more detailed information about the resulting sources, but the process is less generalizable from one telescope's/instrument's/detector's images to another's. This path may also be less desirable for synthetic image creation because the cutouts will be square patches and thus look a little more, well, synthetic." ] @@ -287,6 +281,8 @@ "metadata": {}, "outputs": [], "source": [ + "%%script false --no-raise-error\n", + "\n", "# download sources' location and flux information\n", "_source_info1 = Table.read('https://archive.stsci.edu/pub/hlsp/hdf/wfpc_hdfn_v2catalog/HDFN_wfpc_v2generic.cat',\n", " format='ascii')\n", @@ -307,6 +303,8 @@ "metadata": {}, "outputs": [], "source": [ + "%%script false --no-raise-error\n", + "\n", "# sort sources by flux within 71.1 pixel diameter of source, or aperture 11\n", "_sources.sort('FLUX_APER_11', reverse=True)\n", "\n", @@ -322,6 +320,7 @@ "metadata": {}, "outputs": [], "source": [ + "%%script false --no-raise-error\n", "_sources[:5]" ] }, @@ -332,6 +331,8 @@ "metadata": {}, "outputs": [], "source": [ + "%%script false --no-raise-error\n", + "\n", "# convert the sources' WCS locations to in-image pixel values\n", "_image_wcs = WCS(image_header)\n", "_sources_x, _sources_y = image_wcs.world_to_pixel_values(_sources['ALPHA_J2000'],\n", @@ -356,6 +357,8 @@ }, "outputs": [], "source": [ + "%%script false --no-raise-error\n", + "\n", "# define parameters for cutout loop\n", "_cutout_list = []\n", "_first_source = 0 # which source to cut out first (in descending order by flux)\n", @@ -395,14 +398,6 @@ "___" ] }, - { - "cell_type": "markdown", - "id": "5a6b90a4-e2db-4994-9436-580bd72245b2", - "metadata": {}, - "source": [ - "___" - ] - }, { "cell_type": "markdown", "id": "746288d9-3430-4032-898e-156cd0e39d36", @@ -717,7 +712,7 @@ "metadata": {}, "outputs": [], "source": [ - "savepath = ''\n", + "savepath = './'\n", "synth_file = 'synthetic_mosviz_image.fits'\n", "\n", "fits.writeto(savepath + synth_file, synth_image,\n", @@ -781,7 +776,10 @@ }, "outputs": [], "source": [ - "savepath2 = ''\n", + "savepath2 = './mosviz-cutouts/'\n", + "\n", + "if not os.path.exists(savepath2):\n", + " os.makedirs(savepath2)\n", "\n", "for i, cut in enumerate(synth_cutout_list):\n", " sdt = cut.data\n", @@ -792,7 +790,7 @@ " shd['OBJ_DEC'] = (decs[i], 'Source declination')\n", " shd['OBJ_ROT'] = (0, '')\n", " \n", - " fits.writeto(savepath2 + f\"/user/jotor/mosviz-cutouts-copy/synth_cutout{src}.fits\",\n", + " fits.writeto(savepath2 + f\"synth_cutout{i}.fits\",\n", " sdt, header=shd, overwrite=True)\n", " \n", " plt.imshow(sdt, **imshow_params, vmin=0, vmax=ccd_synth_img.data.std())\n", From a207acbb5e6f07dda6c2aec1c548d3ad63c1ef74 Mon Sep 17 00:00:00 2001 From: ojustino Date: Wed, 16 Jun 2021 18:20:03 -0400 Subject: [PATCH 6/7] Substituted TARG_RA/DEC for SRCRA/DEC --- .../concepts/generate_mosviz_photometry.ipynb | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/notebooks/concepts/generate_mosviz_photometry.ipynb b/notebooks/concepts/generate_mosviz_photometry.ipynb index b4870a267e..3799813478 100644 --- a/notebooks/concepts/generate_mosviz_photometry.ipynb +++ b/notebooks/concepts/generate_mosviz_photometry.ipynb @@ -267,7 +267,7 @@ "id": "87691179-7fd9-497c-9260-8f24c8b31d41", "metadata": {}, "source": [ - "### _[An alternate cutout generation process using catalog file(s)]_\n", + "### _[An alternate, catalog file-based cutout generation process]_\n", "\n", "It is also possible to build a list of cutouts from an already-existing source catalog; this section serves as a reference for how it is done. Its cells are not necessary for the completion of the notebook and do not run due to the `%%script false --no-raise-error` statements at the top of each one.\n", "\n", @@ -407,7 +407,7 @@ "\n", "Note that this section's cells only run on a connection with access to the STScI VPN. They use copies of JWST pipeline files from May 19, 2021.\n", "\n", - "To run them elsewhere and/or use more current files, visit the [Artifactory](https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth%2Ftest_nirspec_mos_spec3), ensure the `jwst-pipeline/dev/truth/test_nirspec_mos_spec3/` folder is the one selected in the directory tree, and download all files whose names begin with `jw00626-o030` and end with `nirspec_f170lp-g235m_x1d.fits`. (The original version of this notebook uses six of these Level 3 data products.) Then, change `filepath` in the next cell to the location where you saved the downloaded files. Remember the trailing forward slash." + "To run them elsewhere and/or use more current files, visit the [Artifactory](https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth%2Ftest_nirspec_mos_spec3), ensure the `jwst-pipeline/dev/truth/test_nirspec_mos_spec3/` folder is the one selected in the directory tree, and download all files whose names begin with `jw00626-o030` and end with `nirspec_f170lp-g235m_x1d.fits` or `nirspec_f170lp-g235m_s2d.fits`. (The original version of this notebook uses six of these Level 3 data products.) Then, change `filepath` in the next cell to the location where you saved the downloaded files, remembering the trailing forward slash." ] }, { @@ -419,8 +419,8 @@ "source": [ "# view level 3 spectra FITS header information\n", "filepath = '/user/jotor/jwst-pipeline-lvl3-copy/'\n", - "x1d_header = fits.getheader(filepath + 'jw00626-o030_s00000_nirspec_f170lp-g235m_x1d.fits')\n", - "#s2d_header = fits.getheader(filepath + 'jw00626-o030_s00000_nirspec_f170lp-g235m_s2d.fits')" + "x1d_header = fits.getheader(filepath + 'jw00626-o030_s00000_nirspec_f170lp-g235m_x1d.fits', ext=1)\n", + "#s2d_header = fits.getheader(filepath + 'jw00626-o030_s00000_nirspec_f170lp-g235m_s2d.fits'', ext=1)" ] }, { @@ -430,7 +430,7 @@ "metadata": {}, "outputs": [], "source": [ - "x1d_header['TARG_RA'], x1d_header['TARG_DEC']" + "x1d_header['SRCRA'], x1d_header['SRCDEC']" ] }, { @@ -448,7 +448,7 @@ "id": "a225379a-d1b7-403c-9421-02dbf3386ced", "metadata": {}, "source": [ - "Notice that the `WCS` object above contains its default settings; this is because the spectrum's header lacks WCS information. Additionally, examining all of this observation's file headers reveals that they all have the same RA/Dec, which is not what we expect for its different \"pointings.\"" + "This header's listed source RA/Dec of (0, 0) seems suspicious. Examining all of this observation's file headers reveals that it is indeed not in line with the other pointings." ] }, { @@ -459,10 +459,10 @@ "outputs": [], "source": [ "# search for RA/Dec information from Artifactory observation files\n", - "x1d_header_list = [fits.getheader(file)\n", + "x1d_header_list = [fits.getheader(file, ext=1)\n", " for file in glob(filepath + 'jw00626*x1d.fits')]\n", "\n", - "_ras, _decs = np.array([[h['TARG_RA'], h['TARG_DEC']] for h in x1d_header_list]).T\n", + "_ras, _decs = np.array([[h['SRCRA'], h['SRCDEC']] for h in x1d_header_list]).T\n", "_ras, _decs" ] }, @@ -473,7 +473,7 @@ "source": [ "### Create synthetic image\n", "\n", - "Since the headers lack the information we seek, we randomly generate our sources' RA/Dec information in a predetermined patch of sky. We choose the patch size by approximating NIRSpec's 3.6'x3.4' Micro-Shutter Assembly (MSA) to a square field of view." + "To get a more realistic field of view, we randomly generate our sources' RA/Dec information in a predetermined patch of sky. We choose the patch size by approximating NIRSpec's 3.6'x3.4' Micro-Shutter Assembly (MSA) to a square field of view." ] }, { From b2643c1c2fb09a90db66ffa30125056a193b3aba Mon Sep 17 00:00:00 2001 From: ojustino Date: Wed, 23 Jun 2021 14:06:06 -0400 Subject: [PATCH 7/7] Made edits in response to review --- .../concepts/generate_mosviz_photometry.ipynb | 30 +++++++++++++------ 1 file changed, 21 insertions(+), 9 deletions(-) diff --git a/notebooks/concepts/generate_mosviz_photometry.ipynb b/notebooks/concepts/generate_mosviz_photometry.ipynb index 3799813478..dddd6bd18a 100644 --- a/notebooks/concepts/generate_mosviz_photometry.ipynb +++ b/notebooks/concepts/generate_mosviz_photometry.ipynb @@ -7,7 +7,7 @@ "source": [ "# Synthetic image and cutout creation for Mosviz pipeline data\n", "\n", - "**Motiviation**: The synthetic dataset we currently possess from the JWST data pipeline team for use in Mosviz contains simulated spectra from NIRSpec but no associated NIRCam photometry. We would like to have test imagery to display in Mosviz alongside these 2D and 1D spectra since we are not aware that the pipeline team has plans to produce any.\n", + "**Motivation**: The synthetic dataset we currently possess from the JWST data pipeline team for use in Mosviz contains simulated spectra from NIRSpec but no associated NIRCam photometry. We would like to have test imagery to display in Mosviz alongside these 2D and 1D spectra since we are not aware that the pipeline team has plans to produce any.\n", "\n", "**Goal**: Populate an image of background noise with properly-scaled galaxy cutouts sourced from a Hubble Space Telescope image and placed at their analogous locations in the new image. These galaxies' real spectra do not necessarily correspond with those in our dataset, but we care more about the veneer of having photometry to match with our spectra in Mosviz at this point.\n", "\n", @@ -17,7 +17,7 @@ "- We wanted to scrape RA/Dec information from the data pipeline products to get a range of coordinates over which to scale our synthetic image, but it appears that the pipeline's data products lack `\"TARG_RA\"` or `\"TARG_DEC\"` keywords in their headers.\n", " - The data products also do not appear to have WCS information. Though it is not strictly needed to achieve this notebook's goals, it would be convenient to have.\n", "- There does not appear to be an observation with Level 2 data in `jwst-pipeline/truth/test_nirspec_mos_spec2` and Level 3 data in `jwst-pipeline/truth/test_nirspec_mos_spec3`. All observations are either Level 2 only or Level 3 only.\n", - "- _(Resolved)_ A good number of the cutouts from the first couple of field images we tested had intrustions from other galaxies due to crowding. We settled on the Hubble Deep Field as a good source image, but had we not, we may have considered using galaxies modeled with Sersic profiles to get cleaner cutouts to inject into our synthetic image." + "- _(Resolved)_ A good number of the cutouts from the first couple of field images we tested had intrusions from other galaxies due to crowding. We settled on the Hubble Deep Field as a good source image, but had we not, we may have considered using galaxies modeled with Sersic profiles to get cleaner cutouts to inject into our synthetic image." ] }, { @@ -39,7 +39,8 @@ "- The objects from `photutils.segmentation` do the work of analytically finding bright sources, separating those that overlap, and creating a catalog of the resulting information.\n", "- We use `matplotlib.pyplot` to preview the field image, the cutouts, and various stages of our synthetic image.\n", "- We use `numpy` to facilitate several specialized mathematical and array-related operations.\n", - "- `os` helps us create a new directory for the final cutouts from the synthetic image." + "- `os` helps us create a new directory for the final cutouts from the synthetic image.\n", + "- `%matplotlib inline` is notebook \"magic\" that helps display in-notebook plots consistently for different users." ] }, { @@ -51,6 +52,9 @@ }, "outputs": [], "source": [ + "import os\n", + "from glob import glob\n", + "\n", "from astropy import units as u\n", "from astropy.coordinates import SkyCoord\n", "from astropy.convolution import Gaussian2DKernel\n", @@ -60,13 +64,21 @@ " sigma_clipped_stats, sigma_clip)\n", "from astropy.table import Table, join\n", "from astropy.wcs import WCS\n", - "from glob import glob\n", "from photutils.segmentation import (detect_sources, deblend_sources,\n", " SourceCatalog)\n", "\n", "import matplotlib.pyplot as plt\n", - "import numpy as np\n", - "import os" + "import numpy as np" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "a451ac82-45ce-4ed9-a0fc-b122ceb128cc", + "metadata": {}, + "outputs": [], + "source": [ + "%matplotlib inline" ] }, { @@ -407,7 +419,7 @@ "\n", "Note that this section's cells only run on a connection with access to the STScI VPN. They use copies of JWST pipeline files from May 19, 2021.\n", "\n", - "To run them elsewhere and/or use more current files, visit the [Artifactory](https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth%2Ftest_nirspec_mos_spec3), ensure the `jwst-pipeline/dev/truth/test_nirspec_mos_spec3/` folder is the one selected in the directory tree, and download all files whose names begin with `jw00626-o030` and end with `nirspec_f170lp-g235m_x1d.fits` or `nirspec_f170lp-g235m_s2d.fits`. (The original version of this notebook uses six of these Level 3 data products.) Then, change `filepath` in the next cell to the location where you saved the downloaded files, remembering the trailing forward slash." + "To run them elsewhere and/or use more current files, visit the [Artifactory](https://bytesalad.stsci.edu/ui/repos/tree/General/jwst-pipeline%2Fdev%2Ftruth%2Ftest_nirspec_mos_spec3), ensure the `jwst-pipeline/dev/truth/test_nirspec_mos_spec3/` folder is the one selected in the directory tree, and download all files whose names begin with `jw00626-o030` and end with `nirspec_f170lp-g235m_x1d.fits` or `nirspec_f170lp-g235m_s2d.fits`. (The original version of this notebook uses six of these Level 3 data products.) Then, change `filepath` in the next cell to the location where you saved the downloaded files, remembering the trailing slash." ] }, { @@ -420,7 +432,7 @@ "# view level 3 spectra FITS header information\n", "filepath = '/user/jotor/jwst-pipeline-lvl3-copy/'\n", "x1d_header = fits.getheader(filepath + 'jw00626-o030_s00000_nirspec_f170lp-g235m_x1d.fits', ext=1)\n", - "#s2d_header = fits.getheader(filepath + 'jw00626-o030_s00000_nirspec_f170lp-g235m_s2d.fits'', ext=1)" + "s2d_header = fits.getheader(filepath + 'jw00626-o030_s00000_nirspec_f170lp-g235m_s2d.fits', ext=1)" ] }, { @@ -430,7 +442,7 @@ "metadata": {}, "outputs": [], "source": [ - "x1d_header['SRCRA'], x1d_header['SRCDEC']" + "(x1d_header['SRCRA'], x1d_header['SRCDEC']), (s2d_header['SRCRA'], s2d_header['SRCDEC'])" ] }, {