DP02_15_Survey_Property_Maps.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "749b0ddf",
   "metadata": {},
   "source": [
    "# Visualizing Survey Property Maps\n",
    "\n",
    "<img align=\"left\" src = https://project.lsst.org/sites/default/files/Rubin-O-Logo_0.png width=250 style=\"padding: 10px\"> \n",
    "<br>\n",
    "Contact author: Eli Rykoff <br>\n",
    "Last verified to run: 2024-12-17 <br>\n",
    "LSST Science Pipelines version: Weekly 2024_50 <br>\n",
    "Container Size: large <br>\n",
    "Targeted learning level: intermediate <br>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9da1a210-d858-42fe-8591-570965b8be1a",
   "metadata": {},
   "source": [
    "**Description:** Demonstrate tools to visualize survey property maps."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "80a0baf5-51ad-40ec-8991-060a7b27c289",
   "metadata": {},
   "source": [
    "**Skills:** Load and visualize survey property maps using healsparse and skyproj."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "393da88f-7978-4920-aa4a-a9830df6eed9",
   "metadata": {},
   "source": [
    "**LSST Data Products:** Survey property maps."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5c67fab9-136a-4adc-bb42-142b91ab69dd",
   "metadata": {},
   "source": [
    "**Packages:** healsparse, skyproj, lsst.daf.butler"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8f72b27f",
   "metadata": {},
   "source": [
    "**Credit:**\n",
    "This notebook was originally developed by Eli Rykoff with contributions from Alex Drlica-Wagner and Melissa Graham."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "28e91cbf-ab7f-4e26-9276-b00299d6065e",
   "metadata": {},
   "source": [
    "**Get Support:**\n",
    "Find DP0-related documentation and resources at <a href=\"https://dp0.lsst.io\">dp0.lsst.io</a>. Questions are welcome as new topics in the <a href=\"https://community.lsst.org/c/support/dp0\">Support - Data Preview 0 Category</a> of the Rubin Community Forum. Rubin staff will respond to all questions posted there."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cfc73be0",
   "metadata": {},
   "source": [
    "## 1. Introduction\n",
    "\n",
    "This notebook will teach the user how to load and visualize survey property maps generated by the Rubin Science Pipelines. Data products are accessed through the Butler, and the user is expected to be familiar with the content of the introductory Butler tutorial. It introduces two new packages, `healsparse` and `skyproj`, for working with survey property maps."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "dc36f107",
   "metadata": {},
   "source": [
    "### 1.1 Package Imports\n",
    "\n",
    "Import general python packages and the Butler from the science pipelines.\n",
    "\n",
    "Import two additional packages for working with the survey property maps.\n",
    "\n",
    "The `healsparse` package provides utilities for reading and manipulating sparse healpix maps.\n",
    "More information can be found in the documentation \"[HealSparse: A sparse implementation of HEALPix](https://healsparse.readthedocs.io/en/latest/)\".\n",
    "\n",
    "The `skyproj` package provides utilities for visualizing both sparse and dense HEALPix maps, as described in the documentation \"[SkyProj: Sky Projections with matplotlib and PROJ](https://skyproj.readthedocs.io/en/latest/)\"."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cddc1458",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# general python packages\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from astropy.visualization import ZScaleInterval, LinearStretch, ImageNormalize\n",
    "from astropy.wcs import WCS\n",
    "\n",
    "# packages for working with sparse healpix maps\n",
    "import healsparse as hsp\n",
    "import skyproj\n",
    "\n",
    "# LSST packages\n",
    "from lsst.daf.butler import Butler\n",
    "import lsst.geom as geom\n",
    "\n",
    "# allow interactive plots\n",
    "%matplotlib widget\n",
    "\n",
    "# default plot style is accessible\n",
    "plt.style.use('tableau-colorblind10')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ec51ac0b",
   "metadata": {
    "tags": []
   },
   "source": [
    "## 2. Access Survey Property Maps\n",
    "\n",
    "Survey property maps are created as part of the LSST Science Pipelines.\n",
    "They take the form of sparse HEALPix maps, where the survey property at each spatial pixel is identified by a pixel number/pixel value pair.\n",
    "\n",
    "Start by creating an instance of the Butler and using it to access these maps for DP0."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "81b8cd59-1ba3-4eaa-846f-6478ed0c3cf5",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "config = 'dp02'\n",
    "collections = '2.2i/runs/DP0.2'\n",
    "butler = Butler(config, collections=collections)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8ae887d0-8c38-4460-91e0-004bf7ea198b",
   "metadata": {},
   "source": [
    "Determine which property maps are available for the survey."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "32153567-355b-46fa-840d-6db79763b5c1",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "for dtype in sorted(butler.registry.queryDatasetTypes(expression=\"*consolidated_map*\")):\n",
    "    print(dtype.name)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "703779fe-9d88-4dd5-a5b0-470670996542",
   "metadata": {},
   "source": [
    "Each of these products represents a healsparse map containing the value of an individual survey property.\n",
    "\n",
    "The meaning of these types is:\n",
    "\n",
    "* `deepCoadd_dcr_ddec_consolidated_map_weighted_mean`: Average effect of differential chromatic refraction (DCR) in declination direction\n",
    "* `deepCoadd_dcr_dra_consolidated_map_weighted_mean`: Average effect of differential chromatic refraction (DCR) in right ascension direction\n",
    "* `deepCoadd_dcr_e1_consolidated_map_weighted_mean`: Average effect of differential chromatic refraction (DCR) on psf e1\n",
    "* `deepCoadd_dcr_e2_consolidated_map_weighted_mean`: Average effect of differential chromatic refraction (DCR) on psf e2\n",
    "* `deepCoadd_exposure_time_consolidated_map_sum`: Total exposure time (seconds)\n",
    "* `deepCoadd_psf_e1_consolidated_map_weighted_mean`: Weighted mean of psf e1 of images input to coadd\n",
    "* `deepCoadd_psf_e2_consolidated_map_weighted_mean`: Weighted mean of psf e2 of images input to coadd\n",
    "* `deepCoadd_psf_maglim_consolidated_map_weighted_mean`: PSF Flux 5-sigma magnitude limit (AB)\n",
    "* `deepCoadd_psf_size_consolidated_map_weighted_mean`: Weighted mean of psf size of images input to coadd (pixels)\n",
    "* `deepCoadd_sky_background_consolidated_map_weighted_mean`: Weighted mean of sky background of images input to coadd (ADU)\n",
    "* `deepCoadd_sky_noise_consolidated_map_weighted_mean`: Weighted mean of sky noise of images input to coadd (ADU)\n",
    "\n",
    "\n",
    "Note that the DCR maps are proportionality maps; that is, the expected effect will be proportional to the value in the map with an arbitrary/empirically derived constant of proportionality.\n",
    "\n",
    "Read a map by specifying the map name and a band."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5497b232-0d72-43fb-9999-cd756bd039bc",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "hspmap = butler.get('deepCoadd_psf_maglim_consolidated_map_weighted_mean', band='i')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "72b36c44-f582-4158-bd95-09a21a73e814",
   "metadata": {},
   "source": [
    "## 3. Manipulating Survey Property Maps\n",
    "\n",
    "The survey property maps are provided as HealSparseMap objects.\n",
    "\n",
    "We provide a few very brief examples here."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a0275ed9-e997-464b-b39c-e259b3a4b765",
   "metadata": {},
   "source": [
    "To conserve memory, HealSparse uses a dual-map approach, where a low-resolution full-sky “coverage map” is combined with a high resolution map containing the pixel data where it is available. It is easy to find the resolution of these maps."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3f719f7a-9aa2-4630-966c-95faaa714069",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "print(hspmap)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c76d13d6-ca16-4c86-b1d3-1ee1047de108",
   "metadata": {},
   "source": [
    "Each pixel of the healsparse map corresponds to a small region of the sky.\n",
    "The value of the map corresponds to the value of the survey property at that location.\n",
    "\n",
    "To access the survey property value at a specific location or set of locations, query for the map value using the `get_values_pos` functionality."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a4b69492-fd6c-4e67-861c-ec9d8d064705",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "print(hspmap.get_values_pos(60, -37))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a9503446-0376-4de8-8919-f13cc1063a60",
   "metadata": {},
   "source": [
    "Query for the map value at an array of locations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d4283aad-0dff-4f21-a600-a0b07cef0c4a",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "ra = np.linspace(59.5, 60.5, 5)\n",
    "print('RA: ', ra)\n",
    "dec = np.linspace(-37.5, -36.5, 5)\n",
    "print('Dec: ', dec)\n",
    "for d in dec:\n",
    "    print(hspmap.get_values_pos(ra, d))\n",
    "del ra, dec"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "37d4efcc-27c0-4fe1-bf71-8825bb0872df",
   "metadata": {},
   "source": [
    "If you ask for the value of the map outside of the region where it is defined, you will get a sentinel value."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3798eef5-e153-475b-acd4-da03a215ffc5",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "print(hspmap.get_values_pos(180, 0))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "063d743a-7c65-49b8-8521-b559aea71c0b",
   "metadata": {},
   "source": [
    "## 4. Visualizing Survey Property Maps\n",
    "\n",
    "Now that we know how to access the values of the healsparse map, we can put together our own simple visualization by creating a grid of RA, Dec values and asking for the map values.\n",
    "We can then plot these values with matplotlib.\n",
    "Note that if you pan or zoom, the map does not update in resolution or coverage area."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5e98a50c-b02f-4a79-8edb-0cb4afb58eed",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "ra = np.linspace(59.5, 60.5, 250)\n",
    "dec = np.linspace(-37.5, -36.5, 250)\n",
    "x, y = np.meshgrid(ra, dec)\n",
    "values = hspmap.get_values_pos(x, y)\n",
    "\n",
    "fig = plt.figure(figsize=(8, 5))\n",
    "plt.pcolormesh(x, y, values)\n",
    "plt.xlabel(\"Right Ascension (deg)\")\n",
    "plt.ylabel(\"Declination (deg)\")\n",
    "plt.colorbar(label=\"PSF Maglim (i-band)\")\n",
    "plt.show()\n",
    "\n",
    "del fig, ra, dec, x, y, values"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "782a6975-7ca7-4254-8b78-08182c9772d2",
   "metadata": {},
   "source": [
    "> Figure 1: The DP0.2 survey property map for PSF magnitude limit in the $i$-band, for a small sky region.\n",
    "\n",
    "The `skyproj` package provides much more advanced plotting capabilities.\n",
    "Here we will demonstrate some basic use-cases for dealing with survey property maps.\n",
    "\n",
    "The use of `McBrydeSkyproj` with `lon_0=65.0` creates a visualization using the McBryde-Thomas Flat Polar Quartic projection, centered at 65 deg longitude, which is appropriate for the DP0.2 footprint."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6ea07e2d-8b5d-4d91-8aeb-230bcd46eef5",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(8, 5))\n",
    "sp = skyproj.McBrydeSkyproj(ax=ax, lon_0=65.0)\n",
    "sp.draw_hspmap(hspmap)\n",
    "sp.draw_colorbar(label='PSF Maglim (i-band)')\n",
    "plt.show()\n",
    "\n",
    "del fig, ax, sp"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8270eda5-2e19-4eca-ac46-941a86e67f3d",
   "metadata": {},
   "source": [
    "> Figure 2: The same survey property map data as in Figure 2, but visualized with the McBryde-Thomas projection over the whole area of DP0.2.\n",
    "\n",
    "Notice that the edges of the survey are pulling the color scale and making it hard to see variations in the main survey area. We can create another visualization with a new colormap range to emphasize these smaller veriations.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3b937aa2-ed9b-43ad-87e7-bb6e3d950caa",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(8, 5))\n",
    "sp = skyproj.McBrydeSkyproj(ax=ax, lon_0=65.0)\n",
    "sp.draw_hspmap(hspmap, vmin=26.0, vmax=26.3)\n",
    "sp.draw_colorbar(label='PSF Maglim (i-band)')\n",
    "plt.show()\n",
    "\n",
    "del fig, ax, sp"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "31593b77-7567-4522-ade8-1a15b6ef3399",
   "metadata": {},
   "source": [
    "> Figure 3: The same data as in Figure 2, but with a different scaling to emphasize small variations across the field.\n",
    "\n",
    "The above maps are interactive thanks to our use of `%matplotlib widget` at the beginning of the notebook.  You can use the matplotlib \"Zoom to rectangle\" tool to draw a box and zoom in.  As you move the mouse around the map, the bottom label will show the longitude, latitude, and map value.  When you zoom in the default is that the colorbar will autoscale according to the data in the box.  If you want to turn this off you can specify `sp.set_autorescale(False)`\n",
    "\n",
    "You can also zoom on a specific longitude/latitude range, as shown below, with the vmin/vmax changed to make it look better."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "59007020-1594-435b-96c4-cf48bd810522",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(figsize=(8, 5))\n",
    "sp = skyproj.McBrydeSkyproj(ax=ax, lon_0=65.0)\n",
    "sp.draw_hspmap(hspmap, vmin=26.0, vmax=26.3, lon_range=[55, 60], lat_range=[-40, -35])\n",
    "sp.draw_colorbar(label='PSF Maglim (i-band)')\n",
    "plt.show()\n",
    "\n",
    "del fig, ax, sp"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fab0ed38-2de3-42d3-bb72-916f3f799c06",
   "metadata": {},
   "source": [
    "> Figure 4: The same data and scaling as in Figure 3, but zoomed in to a smaller region of sky.\n",
    "\n",
    "## 5. Visualizing multiple maps, side-by-side with a deepCoadd image\n",
    "\n",
    "We can also create a 2x2 grid with the deep coadded image alongside the survey property maps of the the magnitude limit, PSF size, and sky background.\n",
    "\n",
    "### 5.1. Retrieve an i-band deepCoadd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ff0d1e48-2aa3-48f7-8d5d-e886d8abcb96",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "ra_deg = 55.745834\n",
    "dec_deg = -32.269167\n",
    "spherePoint = geom.SpherePoint(ra_deg*geom.degrees, dec_deg*geom.degrees)\n",
    "skymap = butler.get('skyMap')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "35063a3c-fd31-4960-8db7-8b01db88513e",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "tract = skymap.findTract(spherePoint)\n",
    "patch = tract.findPatch(spherePoint)\n",
    "tract = tract.tract_id\n",
    "patch = patch.getSequentialIndex()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "29ae74a8-37a3-4d51-8eb0-8c1aebe55c18",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "dataId = {'band': 'i', 'tract': tract, 'patch': patch}\n",
    "deepCoadd = butler.get('deepCoadd', dataId=dataId)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "84a799e4-b818-4b73-83db-2964f723397f",
   "metadata": {},
   "source": [
    "Use the properties of the deepCoadd to obtain the WCS (world coordinate system), the bounding box, and the corners."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7fc324af-dcb5-49d3-aff1-8360d6e4bf60",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "deepCoadd_wcs = deepCoadd.wcs\n",
    "deepCoadd_bbox = deepCoadd.getBBox()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6acc81e0-1821-454a-bc6e-761311c0760e",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "deepCoadd_extent = (deepCoadd.getBBox().beginX, deepCoadd.getBBox().endX,\n",
    "                    deepCoadd.getBBox().beginY, deepCoadd.getBBox().endY)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d2216ee0-f368-4031-88f7-b6de516a88ac",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "corners = [deepCoadd_wcs.pixelToSky(deepCoadd_bbox.beginX, deepCoadd_bbox.beginY),\n",
    "           deepCoadd_wcs.pixelToSky(deepCoadd_bbox.beginX, deepCoadd_bbox.endY),\n",
    "           deepCoadd_wcs.pixelToSky(deepCoadd_bbox.endX, deepCoadd_bbox.endY),\n",
    "           deepCoadd_wcs.pixelToSky(deepCoadd_bbox.endX, deepCoadd_bbox.beginY)]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "43135e2a-4f0d-4636-8a5a-0c66a1b5601b",
   "metadata": {},
   "source": [
    "Option to display the deepCoadd on its own."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "02937c63-61cb-4859-af64-1f729aad2f16",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# temp = deepCoadd.image.array.flatten()\n",
    "# norm = ImageNormalize(temp, interval=ZScaleInterval(), stretch=LinearStretch())\n",
    "\n",
    "# fig = plt.figure()\n",
    "# plt.subplot(projection=WCS(deepCoadd.getWcs().getFitsMetadata()))\n",
    "# im = plt.imshow(deepCoadd.image.array, cmap='gray', norm=norm,\n",
    "#                 extent=deepCoadd_extent, origin='lower')\n",
    "# plt.grid(color='white', ls='solid')\n",
    "# plt.xlabel('Right Ascension')\n",
    "# plt.ylabel('Declination')\n",
    "# plt.show()\n",
    "\n",
    "# del fig, im, temp, norm"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fa6560d7-db7b-4f68-9c6a-c38051122f46",
   "metadata": {},
   "source": [
    "### 5.2. Retrieve the PSF size and sky background survey property maps for the i-band"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b5baed88-8ef6-4061-a2d1-f7913ea539ad",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "hspmap_psf = butler.get('deepCoadd_psf_size_consolidated_map_weighted_mean', band='i')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "43baf0b7-b632-48e5-9ccc-75ef75b73b04",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "hspmap_sky = butler.get('deepCoadd_sky_background_consolidated_map_weighted_mean', band='i')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "29588cab-fb48-4a90-8453-1c9c9ba4cb1a",
   "metadata": {},
   "source": [
    "Define the linear space to use for the grid of survey property map data, which matches the extent of the deepCoadd."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8f9fdba7-75f1-45c0-9bd0-eb9312b4e635",
   "metadata": {},
   "outputs": [],
   "source": [
    "for corner in corners:\n",
    "    print('%6.3f %5.2f' % (corner.getRa().asDegrees(), corner.getDec().asDegrees()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0d02779e-d100-447b-9274-3c7492defc78",
   "metadata": {},
   "outputs": [],
   "source": [
    "ra = np.linspace(55.514, 55.790, 150)\n",
    "dec = np.linspace(-32.09, -32.32, 150)\n",
    "x, y = np.meshgrid(ra, dec)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f3f7128a-61c5-4eba-89e4-a464482fc28c",
   "metadata": {},
   "source": [
    "### 5.3. Create the multi-panel figure\n",
    "\n",
    "It's challenging to overlay a projected RA, Dec grid in a subplot.\n",
    "\n",
    "Instead, prepare to substitute pixels for sky coordinates as tick labels for the image."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "81afd11f-268e-49b7-b043-f1128462fe16",
   "metadata": {},
   "outputs": [],
   "source": [
    "xtick_locs = np.linspace(deepCoadd_extent[0], deepCoadd_extent[1], 5)\n",
    "ytick_locs = np.linspace(deepCoadd_extent[2], deepCoadd_extent[3], 5)\n",
    "xtick_lbls = []\n",
    "ytick_lbls = []\n",
    "for xt, yt in zip(xtick_locs, ytick_locs):\n",
    "    temp = deepCoadd_wcs.pixelToSky(xt, yt)\n",
    "    xtick_lbls.append(str(np.round(temp.getRa().asDegrees(), 2)))\n",
    "    ytick_lbls.append(str(np.round(temp.getDec().asDegrees(), 2)))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "babe4851-176e-432c-8bfb-945b6c0e319b",
   "metadata": {},
   "source": [
    "Option to print the tick locations and labels."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ae68bbc0-6dff-48d7-b9e8-11a291d8deaf",
   "metadata": {},
   "outputs": [],
   "source": [
    "# print(xtick_locs)\n",
    "# print(xtick_lbls)\n",
    "# print(ytick_locs)\n",
    "# print(ytick_lbls)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5f6e8381-53af-4e07-be60-5d689be663a9",
   "metadata": {},
   "source": [
    "Show the deepCoadd and the three survey property maps as a 4-subplot figure."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b94ffd86-659e-494a-917c-cee9b4dc2ca4",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(2, 2, figsize=(9, 6))\n",
    "\n",
    "temp = deepCoadd.image.array.flatten()\n",
    "norm = ImageNormalize(temp, interval=ZScaleInterval(), stretch=LinearStretch())\n",
    "im = ax[0, 0].imshow(deepCoadd.image.array, cmap='gray', norm=norm,\n",
    "                     extent=deepCoadd_extent, origin='lower')\n",
    "fig.colorbar(im, ax=ax[0, 0], label=\"Pixel Value (i-band)\")\n",
    "ax[0, 0].set_xticks(xtick_locs, xtick_lbls)\n",
    "ax[0, 0].set_yticks(ytick_locs, ytick_lbls)\n",
    "ax[0, 0].axis('tight')\n",
    "\n",
    "values = hspmap.get_values_pos(x, y)\n",
    "pcm = ax[0, 1].pcolormesh(x, y, values, cmap='viridis')\n",
    "fig.colorbar(pcm, ax=ax[0, 1], label=\"PSF Mag Limit (i-band)\")\n",
    "ax[0, 1].axis('tight')\n",
    "ax[0, 1].invert_xaxis()\n",
    "del pcm, values\n",
    "\n",
    "values = hspmap_psf.get_values_pos(x, y)\n",
    "pcm = ax[1, 0].pcolormesh(x, y, values, cmap='plasma')\n",
    "fig.colorbar(pcm, ax=ax[1, 0], label=\"PSF Size (i-band)\")\n",
    "ax[1, 0].axis('tight')\n",
    "ax[1, 0].invert_xaxis()\n",
    "del pcm, values\n",
    "\n",
    "values = hspmap_sky.get_values_pos(x, y)\n",
    "pcm = ax[1, 1].pcolormesh(x, y, values, cmap='cividis')\n",
    "fig.colorbar(pcm, ax=ax[1, 1], label=\"Sky Background (i-band)\")\n",
    "ax[1, 1].axis('tight')\n",
    "ax[1, 1].invert_xaxis()\n",
    "del pcm, values\n",
    "\n",
    "plt.tight_layout()\n",
    "\n",
    "del fig, ax, temp, norm"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0411a19e-0fae-403c-9b6b-2f80c5a0b503",
   "metadata": {},
   "source": [
    "> Figure 5: Four different views on DP0.2: the $i$-band deepCoadd image pixel data (grayscale, upper left);\n",
    "> the $i$-band PSF magnitude limit (viridis, upper right); the $i$-band PSF size (plasma, lower left);\n",
    "> and the $i$-band sky background (cividis, lower right)."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "LSST",
   "language": "python",
   "name": "lsst"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.9"
  },
  "toc-autonumbering": false
 },
 "nbformat": 4,
 "nbformat_minor": 5
}