exercises/SEDs_in_Firefly.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "1d69721b",
   "metadata": {},
   "source": [
    "# Using Firefly visualization tools in Python to vet SEDs\n",
    "\n",
    "## Learning Goals"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "18f7a018",
   "metadata": {},
   "source": [
    "By the end of this tutorial, the following goals can be accomplished:\n",
    "* Select a position and categories of data for search of possible nearby sources\n",
    "* Pull photometry for different bands via TAP queries on chosen IRSA catalogs\n",
    "* Add data columns from result to table, including energy density values (converted from magnitudes)\n",
    "* Create a simple SED plot using compiled photometric data\n",
    "* Display matching image cutouts from IRSA data sets for selected position"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "90d5623c",
   "metadata": {},
   "source": [
    "## Introduction"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7d1f386a",
   "metadata": {},
   "source": [
    "The purpose of this notebook is to show how to make and vet Spectral Energy Distributions (SEDs) using Firefly. All examples use the python package `pyvo` to search for catalogs and images from IRSA."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3e70fcbb",
   "metadata": {},
   "source": [
    "Detailed specifications for this notebook come from a paper led by Luisa Rebull on the selection of possible Young Stellar Objects (YSOs) in IC 417 (https://doi.org/10.3847/1538-3881/ace32f)."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3dbd3b9b",
   "metadata": {},
   "source": [
    "In this case, the SEDs generated by photometry data (from several different survey catalogs) near a given position can be compared with corresponding images to verify the location of a point source that could also be considered as a YSO if the graph shape is sound."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e29cb285",
   "metadata": {},
   "source": [
    "## Imports"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "df4f32de",
   "metadata": {},
   "source": [
    "Prior to starting Jupyterlab, you should specify the Firefly server by setting the `FIREFLY_URL` environment variable (use the value https://irsa.ipac.caltech.edu/irsaviewer if you're not sure), or by one of the other ways listed on the [jupyter_firefly_extension page](https://github.com/Caltech-IPAC/jupyter_firefly_extensions)."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f162af5f",
   "metadata": {},
   "source": [
    "The packages needed for this notebook are in the `requirements.txt` file. They can be installed with `pip install -r requirements.txt` before you start your Jupyterlab session."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "25d938ff",
   "metadata": {},
   "source": [
    "* _astropy.units_ - for specifying and manipulating units\n",
    "* _astropy.constants_ - for constants\n",
    "* _astropy.coordinates.Skycoord_ - for representing coordinates\n",
    "* _astropy.table.QTable_ - for tables with units\n",
    "* _firefly_client.FireflyClient_ - Python API to Firefly for displaying tables, images and charts\n",
    "* _numpy_ - for working with arrays\n",
    "* _pyvo_ - for queries to Virtual Observatory services at the archives"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "67e316b9",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Uncomment the next line to install dependencies if needed.\n",
    "# !pip install numpy astropy pyvo firefly_client"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a3b6620f",
   "metadata": {},
   "source": [
    "Call the imports below:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3a4a6a0f",
   "metadata": {},
   "outputs": [],
   "source": [
    "import astropy.units as u\n",
    "from astropy import constants as const\n",
    "from astropy.coordinates import SkyCoord\n",
    "from astropy.table import QTable\n",
    "from astropy.utils.data import download_file\n",
    "from firefly_client import FireflyClient\n",
    "import numpy as np\n",
    "import pyvo"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6c216159",
   "metadata": {},
   "source": [
    "## Choose Target as Example"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ffefec11",
   "metadata": {},
   "source": [
    "From [Figure 10](https://iopscience.iop.org/article/10.3847/1538-3881/ace32f#ajace32ff10) of the referenced paper, pick the source in the upper right corner (J052736.37+344940.6):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "30345e44",
   "metadata": {},
   "outputs": [],
   "source": [
    "target = SkyCoord(ra=\"05h27m36.37s\", dec=\"+34d49m40.6s\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "63530051",
   "metadata": {},
   "source": [
    "## Get Photometry from IRSA"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "455da0b9",
   "metadata": {},
   "source": [
    "The available data at IRSA include the following catalogs (with matching radius at 1 arcsec except where noted by reference paper):\n",
    "* Gaia DR3\n",
    "* 2MASS\n",
    "* Spitzer-IRAC (GLIMPSE? SEIP Source List?)\n",
    "* Spitzer-MIPS-24 (GLIMPSE? SEIP Source List?)\n",
    "* AllWISE\n",
    "* CatWISE\n",
    "* unWISE\n",
    "* MSX\n",
    "* Akari\n",
    "* Herschel-PACS (70um and 160um)\n",
    "\n",
    "The [NITARP Wiki page](https://vmcoolwiki.ipac.caltech.edu/index.php/Central_wavelengths_and_zero_points) on central wavelengths may be helpful for finding zero-point information for specific catalogs, as well as the more general [astroquery service](https://astroquery.readthedocs.io/en/latest/svo_fps/svo_fps.html) on SVO filter data."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "669b5740",
   "metadata": {},
   "source": [
    "Whlie this notebook will focus on a subset of the above list (2MASS, Gaia DR3, and WISE data), note that there is a section where additional data can be included as well."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "324fbc3c",
   "metadata": {},
   "source": [
    "The variable `phot_tbl` is the table containing the photometry for the target and all the derived quantities needed to make SED plots:\n",
    "* Right ascension (in degrees)\n",
    "* Declination (in degrees)\n",
    "* Wavelength for the bandpass (in microns)\n",
    "* Photometry measurement\n",
    "* Photometry uncertainty\n",
    "* Is the measurement an upper limit? (Note: 1 means yes, 0 means no)\n",
    "* Zeropoint for the bandpass in units of Janskys (needed for converting photometric measurements to energy densities)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d658cfac",
   "metadata": {},
   "source": [
    "Set the variables for the photmetry table and each quantity in that data structure:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2f64f482",
   "metadata": {},
   "outputs": [],
   "source": [
    "phot_tbl = QTable(names=['ra', 'dec', 'wavelength', 'mag', 'magerr', 'upper_limit', 'zeropoint'],\n",
    "                  units=[u.deg, u.deg, u.micron, None, None, None, u.Jy])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ecbbf496",
   "metadata": {},
   "source": [
    "We query catalog data using IRSA's TAP (VO Table Access Protocol) service. It's sufficient to set up the service once, at the beginning of the session:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8ba7b2aa",
   "metadata": {},
   "outputs": [],
   "source": [
    "irsa_tap = pyvo.dal.TAPService(\"https://irsa.ipac.caltech.edu/TAP\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "619e4dce",
   "metadata": {},
   "source": [
    "### Retrieve 2MASS photometry"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cde3f711",
   "metadata": {},
   "source": [
    "To start, define a function to retrieve 2MASS photometry via a TAP query to IRSA that includes upper limits (note that only the `ph_qual` column on photometric quality of source is used and not the `rd_flg` column on \"read flag\"):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "14d8331f",
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_2mass_phot(target, radius=1*u.arcsec):\n",
    "    \"\"\"Get photometry from 2MASS\n",
    "\n",
    "    Parameters:\n",
    "    -----------\n",
    "    target: astropy.coordinates.SkyCoord\n",
    "          coordinates of the target(s)\n",
    "\n",
    "    radius: astropy.units.Quantity\n",
    "          matching radius, default 1*u.arcsec\n",
    "\n",
    "    Returns:\n",
    "    --------\n",
    "    rows: list\n",
    "         lists of photometry points to add to table\n",
    "    \"\"\"\n",
    "\n",
    "    query = f\"\"\"\n",
    "    SELECT ra,dec,j_m,j_msigcom,h_m,h_msigcom,k_m,k_msigcom,ph_qual,rd_flg\n",
    "    FROM fp_psc\n",
    "    WHERE CONTAINS(POINT('ICRS', ra, dec),\n",
    "    CIRCLE('ICRS', {target.ra.deg}, {target.dec.deg},\n",
    "    {radius.to('deg').value}))=1\n",
    "    \"\"\"\n",
    "\n",
    "    results = irsa_tap.search(query)\n",
    "\n",
    "    rowj = [results[\"ra\"]*u.deg,\n",
    "            results[\"dec\"]*u.deg,\n",
    "            1.235*u.micron,\n",
    "            results[\"j_m\"][0],\n",
    "            results[\"j_msigcom\"][0],\n",
    "            1 if results[\"ph_qual\"] == \"U\" else 0,\n",
    "            1594*u.Jy]\n",
    "\n",
    "    rowh = [results[\"ra\"]*u.deg,\n",
    "            results[\"dec\"]*u.deg,\n",
    "            1.662*u.micron,\n",
    "            results[\"h_m\"][0],\n",
    "            results[\"h_msigcom\"][0],\n",
    "            1 if results[\"ph_qual\"] == \"U\" else 0,\n",
    "            1024*u.Jy]\n",
    "\n",
    "    rowks = [results[\"ra\"]*u.deg,\n",
    "            results[\"dec\"]*u.deg,\n",
    "            2.159*u.micron,\n",
    "            results[\"k_m\"][0],\n",
    "            results[\"k_msigcom\"][0],\n",
    "            1 if results[\"ph_qual\"] == \"U\" else 0,\n",
    "            666.7*u.Jy]\n",
    "\n",
    "    return([rowj, rowh, rowks])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0c93b83b",
   "metadata": {},
   "source": [
    "Consider that while the search radius is currently set at 1 degree from the position, that value can be adjusted by changing what the number on the `radius` line is.\n",
    "\n",
    "Now gather the rows of 2MASS data and add each of them to the photometry table:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5e2eac75",
   "metadata": {},
   "outputs": [],
   "source": [
    "rows = get_2mass_phot(target)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c55365db",
   "metadata": {},
   "outputs": [],
   "source": [
    "for row in rows:\n",
    "    phot_tbl.add_row(row)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "180db683",
   "metadata": {},
   "source": [
    "### Retrieve AllWISE Photometry"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b1854735",
   "metadata": {},
   "source": [
    "According to the paper referenced above, the detections from the AllWISE catalog were retained if the data quality flags were A, B, or C; if the data quality flag was Z, then the data were provisionally retained with a very large error bar, 30% larger than what appears in the catalog."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e6fa648f",
   "metadata": {},
   "source": [
    "Based on that information, define a function to release AllWISE photometry via a TAP query to IRSA that includes the needed data quality flags and retains upper limits:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d1973798",
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_allwise_phot(target, radius=1*u.arcsec):\n",
    "    \"\"\"Get photometry from AllWISE\n",
    "\n",
    "    Parameters:\n",
    "    -----------\n",
    "    target: astropy.coordinates.SkyCoord\n",
    "          coordinates of the target(s)\n",
    "\n",
    "    radius: astropy.units.Quantity\n",
    "          matching radius, defaults to 1*u.arcsec\n",
    "\n",
    "    Returns:\n",
    "    --------\n",
    "    rows: list\n",
    "         lists of photometry points to add to table\n",
    "    \"\"\"\n",
    "\n",
    "    query = f\"\"\"\n",
    "    SELECT ra,dec,w1mpro,w1sigmpro,w2mpro,w2sigmpro,w3mpro,\n",
    "          w3sigmpro,w4mpro,w4sigmpro,ph_qual\n",
    "    FROM allwise_p3as_psd\n",
    "    WHERE CONTAINS(POINT('ICRS', ra, dec),\n",
    "    CIRCLE('ICRS', {target.ra.deg}, {target.dec.deg},\n",
    "    {radius.to('deg').value}))=1\n",
    "    \"\"\"\n",
    "\n",
    "    results = irsa_tap.search(query)\n",
    "\n",
    "    roww1 = [results[\"ra\"]*u.deg,\n",
    "             results[\"dec\"]*u.deg,\n",
    "             3.4*u.micron,\n",
    "             results[\"w1mpro\"][0],\n",
    "             results[\"w1sigmpro\"][0],\n",
    "             1 if results[\"ph_qual\"] == \"U\" else 0,\n",
    "             309.54*u.Jy]\n",
    "\n",
    "    roww2 = [results[\"ra\"]*u.deg,\n",
    "             results[\"dec\"]*u.deg,\n",
    "             4.6*u.micron,\n",
    "             results[\"w2mpro\"][0],\n",
    "             results[\"w2sigmpro\"][0],\n",
    "             1 if results[\"ph_qual\"] == \"U\" else 0,\n",
    "             171.79*u.Jy]\n",
    "\n",
    "    roww3 = [results[\"ra\"]*u.deg,\n",
    "             results[\"dec\"]*u.deg,\n",
    "             12*u.micron,\n",
    "             results[\"w3mpro\"][0],\n",
    "             results[\"w3sigmpro\"][0],\n",
    "             1 if results[\"ph_qual\"] == \"U\" else 0,\n",
    "             31.676*u.Jy]\n",
    "\n",
    "    roww4 = [results[\"ra\"]*u.deg,\n",
    "             results[\"dec\"]*u.deg,\n",
    "             22*u.micron,\n",
    "             results[\"w4mpro\"][0],\n",
    "             results[\"w4sigmpro\"][0],\n",
    "             1 if results[\"ph_qual\"] == \"U\" else 0,\n",
    "             8.3635*u.Jy]\n",
    "\n",
    "    return([roww1, roww2, roww3, roww4])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "33cf6bbe",
   "metadata": {},
   "outputs": [],
   "source": [
    "rows = get_allwise_phot(target)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6a73505e",
   "metadata": {},
   "outputs": [],
   "source": [
    "for row in rows:\n",
    "    phot_tbl.add_row(row)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d8cae3ac",
   "metadata": {},
   "source": [
    "### Retrieve Gaia DR3 Photometry"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "32b98869",
   "metadata": {},
   "source": [
    "Define a function to retrieve Gaia DR3 photometry via a TAP query to IRSA (not considering upper limits as Gaia DR3 ony provides \"detections or not\" info so actual limits are not given):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "76e740f2",
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_gaia_phot(target, radius=1*u.arcsec):\n",
    "    \"\"\"Get photometry from Gaia DR3\n",
    "\n",
    "    Parameters:\n",
    "    -----------\n",
    "    target: astropy.coordinates.SkyCoord\n",
    "          coordinates of the target(s)\n",
    "\n",
    "    radius: astropy.units.Quantity\n",
    "          matching radius, default is 1*u.arcsec\n",
    "\n",
    "    Returns:\n",
    "    --------\n",
    "    rows: list\n",
    "         lists of photometry points to add to table\n",
    "    \"\"\"\n",
    "\n",
    "    query = f\"\"\"\n",
    "    SELECT designation,ra,ra_error,dec,dec_error,ra_dec_corr,phot_g_mean_flux_over_error,phot_g_mean_mag,\n",
    "       phot_bp_mean_flux_over_error,phot_bp_mean_mag,phot_rp_mean_flux_over_error,phot_rp_mean_mag\n",
    "    FROM gaia_dr3_source\n",
    "    WHERE CONTAINS(POINT('ICRS', ra, dec),\n",
    "    CIRCLE('ICRS', {target.ra.deg}, {target.dec.deg},\n",
    "    {radius.to('deg').value}))=1\n",
    "    \"\"\"\n",
    "\n",
    "    results = irsa_tap.search(query)\n",
    "\n",
    "    rowg = [results[\"ra\"]*u.deg,\n",
    "            results[\"dec\"]*u.deg,\n",
    "            0.622*u.micron,\n",
    "            results[\"phot_g_mean_mag\"][0],\n",
    "            2.5/(np.log(10)*results[\"phot_g_mean_flux_over_error\"][0]),\n",
    "            0,\n",
    "            3228.75*u.Jy]\n",
    "\n",
    "    rowbp = [results[\"ra\"]*u.deg,\n",
    "            results[\"dec\"]*u.deg,\n",
    "            0.511*u.micron,\n",
    "            results[\"phot_bp_mean_mag\"][0],\n",
    "            2.5/(np.log(10)*results[\"phot_bp_mean_flux_over_error\"][0]),\n",
    "            0,\n",
    "            3552.01*u.Jy]\n",
    "\n",
    "    rowrp = [results[\"ra\"]*u.deg,\n",
    "            results[\"dec\"]*u.deg,\n",
    "            0.7769*u.micron,\n",
    "            results[\"phot_rp_mean_mag\"][0],\n",
    "            2.5/(np.log(10)*results[\"phot_rp_mean_flux_over_error\"][0]),\n",
    "            0,\n",
    "            2554.95*u.Jy]\n",
    "\n",
    "    return([rowg, rowbp, rowrp])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3346fd6e",
   "metadata": {},
   "outputs": [],
   "source": [
    "rows = get_gaia_phot(target)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "54928c40",
   "metadata": {},
   "outputs": [],
   "source": [
    "for row in rows:\n",
    "    phot_tbl.add_row(row)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "57d0a0fb",
   "metadata": {},
   "source": [
    "### Retrieve CatWISE Photometry"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "485001fb",
   "metadata": {},
   "source": [
    "Define a function to retrieve CatWISE photometry via a TAP query to IRSA that includes upper limits (note that CatWISE2020 Catalog is being used rather than CatWISE Preliminary Catalog):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "54976ae7",
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_catwise_phot(target, radius=1*u.arcsec):\n",
    "    \"\"\"Get photometry from CatWISE\n",
    "\n",
    "    Parameters:\n",
    "    -----------\n",
    "    target: astropy.coordinates.SkyCoord\n",
    "          coordinates of the target(s)\n",
    "\n",
    "    radius: astropy.units.Quantity\n",
    "          matching radius, defaults to 1*u.arcsec\n",
    "\n",
    "    Returns:\n",
    "    --------\n",
    "    rows: list\n",
    "         lists of photometry points to add to table\n",
    "    \"\"\"\n",
    "\n",
    "    query = f\"\"\"\n",
    "    SELECT ra,dec,w1mpro,w1sigmpro,w2mpro,w2sigmpro\n",
    "    FROM catwise_2020\n",
    "    WHERE CONTAINS(POINT('ICRS', ra, dec),\n",
    "    CIRCLE('ICRS', {target.ra.deg}, {target.dec.deg},\n",
    "    {radius.to('deg').value}))=1\n",
    "    \"\"\"\n",
    "\n",
    "    results = irsa_tap.search(query)\n",
    "\n",
    "    roww1 = [results[\"ra\"]*u.deg,\n",
    "             results[\"dec\"]*u.deg,\n",
    "             3.4*u.micron,\n",
    "             results[\"w1mpro\"][0],\n",
    "             results[\"w1sigmpro\"][0],\n",
    "             1 if results[\"w1sigmpro\"] is None else 0,\n",
    "             309.54*u.Jy]\n",
    "\n",
    "    roww2 = [results[\"ra\"]*u.deg,\n",
    "             results[\"dec\"]*u.deg,\n",
    "             4.6*u.micron,\n",
    "             results[\"w2mpro\"][0],\n",
    "             results[\"w2sigmpro\"][0],\n",
    "             1 if results[\"w2sigmpro\"] is None else 0,\n",
    "             171.79*u.Jy]\n",
    "\n",
    "    return([roww1, roww2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bf0af77e",
   "metadata": {},
   "outputs": [],
   "source": [
    "rows = get_catwise_phot(target)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "287f94d9",
   "metadata": {},
   "outputs": [],
   "source": [
    "for row in rows:\n",
    "    phot_tbl.add_row(row)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "38207b7b",
   "metadata": {},
   "source": [
    "### Retrieve unWISE Photometry"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "865a12bb",
   "metadata": {},
   "source": [
    "Define a function to retrieve unWISE photometry via a TAP query to IRSA(note that flux units are given in nano-mag, zero-point is carried over from AllWISE, and no flag for upper limits is given):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b2fa6abd",
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_unwise_phot(target, radius=1*u.arcsec):\n",
    "    \"\"\"Get photometry from CatWISE\n",
    "\n",
    "    Parameters:\n",
    "    -----------\n",
    "    target: astropy.coordinates.SkyCoord\n",
    "          coordinates of the target(s)\n",
    "\n",
    "    radius: astropy.units.Quantity\n",
    "          matching radius, defaults to 1*u.arcsec\n",
    "\n",
    "    Returns:\n",
    "    --------\n",
    "    rows: list\n",
    "         lists of photometry points to add to table\n",
    "    \"\"\"\n",
    "\n",
    "    query = f\"\"\"\n",
    "    SELECT ra,dec,flux_1,dflux_1,flux_2,dflux_2,fluxlbs_1\n",
    "    FROM unwise_2019\n",
    "    WHERE CONTAINS(POINT('ICRS', ra, dec),\n",
    "    CIRCLE('ICRS', {target.ra.deg}, {target.dec.deg},\n",
    "    {radius.to('deg').value}))=1\n",
    "    \"\"\"\n",
    "\n",
    "    results = irsa_tap.search(query)\n",
    "\n",
    "    roww1 = [results[\"ra\"]*u.deg,\n",
    "             results[\"dec\"]*u.deg,\n",
    "             3.4*u.micron,\n",
    "             pow(10,9)*results[\"flux_1\"][0],\n",
    "             pow(10,9)*results[\"dflux_1\"][0],\n",
    "             0,\n",
    "             309.54*u.Jy]\n",
    "\n",
    "    roww2 = [results[\"ra\"]*u.deg,\n",
    "             results[\"dec\"]*u.deg,\n",
    "             4.6*u.micron,\n",
    "             pow(10,9)*results[\"flux_2\"][0],\n",
    "             pow(10,9)*results[\"dflux_2\"][0],\n",
    "             0,\n",
    "             171.79*u.Jy]\n",
    "\n",
    "    return([roww1, roww2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b2161db3",
   "metadata": {},
   "outputs": [],
   "source": [
    "rows = get_unwise_phot(target)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c3e7a3c5",
   "metadata": {},
   "outputs": [],
   "source": [
    "for row in rows:\n",
    "    phot_tbl.add_row(row)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "37e848c3",
   "metadata": {},
   "source": [
    "### Retrieve Optional Photometry"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "91293ec6",
   "metadata": {},
   "source": [
    "Insert more data retrievals here, following the above formats.\n",
    "\n",
    "Note: After more columns are added in the conversion to spectral energy density, the `phot_tbl.add_row` functions will thrown an exception that number of values doesn't match the number of columns."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c1875c30",
   "metadata": {},
   "outputs": [],
   "source": [
    "# SEIP\n",
    "# GLIMPSE\n",
    "\n",
    "# MSX, AKARI, Herschel PACS\n",
    "\n",
    "# Pan-STARRS, UKIDSS?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "87016f0b",
   "metadata": {},
   "source": [
    "## Convert Magnitudes to Energy Densities"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6edd1c3d",
   "metadata": {},
   "source": [
    "Use the zeropoint information embedded by the catalog retrieval functions to compute energy densities:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e788f946",
   "metadata": {},
   "outputs": [],
   "source": [
    "factor = (const.c*phot_tbl[\"wavelength\"]**-1).to(u.Hz)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d8f41171",
   "metadata": {},
   "outputs": [],
   "source": [
    "phot_tbl[\"lam_flam\"] = ((factor*phot_tbl[\"zeropoint\"] *\n",
    "                            10**(-0.4*phot_tbl[\"mag\"]))\n",
    "                       ).to(u.erg*u.s**-1*u.cm**-2)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ac8bfa91",
   "metadata": {},
   "source": [
    "As the above formula was for the photometry measurements given, compute the photometry uncertainty using a similar formula:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3ae93a14",
   "metadata": {},
   "outputs": [],
   "source": [
    "phot_tbl[\"lam_flam_err\"] = phot_tbl[\"lam_flam\"]*phot_tbl[\"magerr\"]*np.log(10)/2.5"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5ad9a91a",
   "metadata": {},
   "source": [
    "Check the energy density values for each converted value below:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d48d1ef5",
   "metadata": {},
   "outputs": [],
   "source": [
    "phot_tbl[\"lam_flam\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cd7c8a7d",
   "metadata": {},
   "outputs": [],
   "source": [
    "phot_tbl[\"lam_flam_err\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6a2ecb4b",
   "metadata": {},
   "source": [
    "Run the following cell to view the full photometry table and its values for each column (ex. zeropoint, upper limit):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3798991e",
   "metadata": {},
   "outputs": [],
   "source": [
    "phot_tbl"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "703ab293",
   "metadata": {},
   "source": [
    "## Upload to Firefly"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "525d5af1",
   "metadata": {},
   "source": [
    "There are two ways to initialize a Firefly client from Python, depending on whether you're running the notebook in JupyterLab or not. Assuming you have `jupyter-firefly-extensions` set up in your environment as explained [here](https://github.com/Caltech-IPAC/jupyter_firefly_extensions/blob/master/README.md), you can use `make_lab_client()` in JupyterLab, which will open the Firefly viewer in a new tab within the Lab. Otherwise, you can use `make_client()` in a Jupyter Notebook (or even a Python shell), which will open the Firefly viewer in a new web browser tab.\n",
    "\n",
    "You also need a Firefly server to communicate with your Firefly Python client. In this notebook, we use a public Firefly server: the IRSA Viewer (https://irsa.ipac.caltech.edu/irsaviewer). However, you can also run a local Firefly server via a [Firefly Docker image](https://hub.docker.com/r/ipac/firefly) and access it at `http://localhost:8080/firefly`. The URL of the Firefly server is read by both `make_client()` and `make_lab_client()` through the environment variable `FIREFLY_URL`. However, `make_client()` also allows you to pass the URL directly as the `url` parameter."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "93e547bd-846a-42aa-ad2d-ee8d7fc0f52f",
   "metadata": {},
   "source": [
    "**Exercise**: Use either the `FireflyClient.make_lab_client()` method or the `FireflyClient.make_client()` method below:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7e7a446d-1044-45bf-8834-f627f4e549e0",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Complete this line of code to initialize the client\n",
    "fc = "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "19b90092-9e15-410f-8e15-653d1b4469fc",
   "metadata": {},
   "source": [
    "**Exercise**: Run the `reinit_viewer` method to clear the state of the display, if needed."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "93bba22a-5174-4704-9977-ad505de02cef",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Fill out code to clear the display"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7926b90d",
   "metadata": {},
   "source": [
    "**Exercise**: In the event that there are problems with the tab opened above, run the `fc.display_url()` command to obtain a web link that can be opened in a browser directly:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "535454de",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Fill out code to display the web link"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "72d193f5",
   "metadata": {},
   "source": [
    "Write the photometry table to a FITS file:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1cd9a0f5",
   "metadata": {},
   "outputs": [],
   "source": [
    "tblpath = \"./phot_tbl.fits\"\n",
    "phot_tbl.write(tblpath, format=\"fits\", overwrite=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8d62357b",
   "metadata": {},
   "source": [
    "Upload the table (FITS file) to Firefly, using `fc.upload_file()`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "59d445a6",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Complete this line of code\n",
    "tval = "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3ab2e7b0",
   "metadata": {},
   "source": [
    "**Exercise**: Show the table in Firefly, using the `show_table` method.\n",
    "Set the table ID to \"photometry\" and give it a title."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6d29d2c7",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Fill out this code cell to show the table"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "55f296ae",
   "metadata": {},
   "source": [
    "Specify the columns and markers for the SED plot:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b4ac26cc",
   "metadata": {},
   "outputs": [],
   "source": [
    "trace1 = dict(tbl_id=\"photometry\", x=\"tables::wavelength\", y=\"tables::lam_flam\", mode=\"markers\",\n",
    "            type=\"scatter\", error_y=dict(array=\"tables::lam_flam_err\"),\n",
    "            marker=dict(size=4))\n",
    "trace_data = [trace1]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "75f5ddf0",
   "metadata": {},
   "source": [
    "Specify the title and axes for the SED plot:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3d3fd945",
   "metadata": {},
   "outputs": [],
   "source": [
    "layout_s = dict(title=\"SED\", xaxis=dict(title=\"Wavelength (microns)\", type=\"log\"),\n",
    "                yaxis=dict(title=\"Energy density (erg/s/cm2)\", type=\"log\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d60a1b09",
   "metadata": {},
   "source": [
    "**Exercise**: Display the SED plot in Firefly, using the `show_chart` method, and\n",
    "specifying the data and layout given above:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6e71ee0c",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Fill out code to display the SED plot"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2ebcfed9",
   "metadata": {},
   "source": [
    "## Retrieve and Display Images"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "304f30fa",
   "metadata": {},
   "source": [
    "Use queries directed at Simple Image Access services to retrieve images:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "492bbd56",
   "metadata": {},
   "outputs": [],
   "source": [
    "irsa_sia = pyvo.dal.SIAService(\"https://irsa.ipac.caltech.edu/cgi-bin/2MASS/IM/nph-im_sia?type=at&ds=asky&\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d4772f4d",
   "metadata": {},
   "source": [
    "### Get 2MASS Images"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5239c06f",
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_2mass_images(target, search_size=1*u.arcsec):\n",
    "    \"\"\"Retrieve 2MASS images\n",
    "\n",
    "    Parameters:\n",
    "    -----------\n",
    "    target: astropy.coordinates.SkyCoord\n",
    "        coordinates of the target(s)\n",
    "\n",
    "    search_size: astropy.units.Quantity\n",
    "        matching radius, defaults to 1*u.arcsec\n",
    "\n",
    "    Returns:\n",
    "    --------\n",
    "    j_filename: `str`\n",
    "         path to J-band image in FITS format\n",
    "\n",
    "    h_filename: `str`\n",
    "         path to H-band image in FITS format\n",
    "\n",
    "    k_filename: `str`\n",
    "         path to K-band image in FITS format\n",
    "    \"\"\"\n",
    "    im_table = irsa_sia.search(pos=target, size=search_size)\n",
    "    url_j = url_h = url_k = None\n",
    "    # Get URL for the first instance of each band\n",
    "    for i in range(len(im_table)):\n",
    "        if url_j is None and im_table[i][\"band\"] == \"J\":\n",
    "            url_j = im_table[i].getdataurl()\n",
    "        if url_h is None and im_table[i][\"band\"] == \"H\":\n",
    "            url_h = im_table[i].getdataurl()\n",
    "        if url_k is None and im_table[i][\"band\"] == \"K\":\n",
    "            url_k = im_table[i].getdataurl()\n",
    "    j_filename = download_file(url_j, cache=True)\n",
    "    h_filename = download_file(url_h, cache=True)\n",
    "    k_filename = download_file(url_k, cache=True)\n",
    "    return (j_filename, h_filename, k_filename)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b76912d1",
   "metadata": {},
   "outputs": [],
   "source": [
    "j_fname, h_fname, k_fname = get_2mass_images(target)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3c730894-6e8b-4f05-a725-18144dcb00ad",
   "metadata": {},
   "source": [
    "**Exercise**: Use the `upload_file` and `show_fits` methods to upload and display the H and K band images, similar to how the J band image is uploaded and displayed.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "df0c82c4",
   "metadata": {},
   "outputs": [],
   "source": [
    "fc.show_fits(fc.upload_file(j_fname), plot_id=\"2MASS J\", title=\"2MASS J\")\n",
    "# Fill out commands to upload the H and K band images\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "23337691-146e-40d4-b270-e744f574d4eb",
   "metadata": {},
   "source": [
    "Turn on alignment and locking for these three images using the `align_images` method, with `lock_match=True`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1196f246-af9f-4300-9c4a-42a8cd2c94bb",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Fill in code to align images"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "966174eb",
   "metadata": {},
   "source": [
    "Note that zooming or panning one image will do the same in the other images as well."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cc4be978-d51f-421f-8484-f5e1ff40ff7e",
   "metadata": {},
   "source": [
    "**Exercise**: For the J band image, set the zoom factor to 4.0, and pan the image to coordinates `x=target.ra.deg`, `y=target.dec.deg`, using `coord=\"j2000\"`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1c916fa6-4c0b-47ed-a71b-61fd7c6cd254",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Fill in zoom and pan commands"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8d349a1d",
   "metadata": {},
   "source": [
    "***"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9d76e3c2",
   "metadata": {},
   "source": [
    "## About This Notebook"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "66dabe59",
   "metadata": {},
   "source": [
    "**Author:** David Shupe (IRSA Scientist) and Joyce Kim (IRSA Scientist) in conjunction with the IRSA Science Team\n",
    "\n",
    "**Updated:** 2024-12-19; modified 2025-01-09 to make exercises\n",
    "\n",
    "**Contact:** [the IRSA Helpdesk](https://irsa.ipac.caltech.edu/docs/help_desk.html) with questions or reporting problems."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}