{ "cells": [ { "cell_type": "markdown", "id": "51b115be", "metadata": {}, "source": [ "# Automated Multiband Forced Photometry on Large Datasets\n", "***\n", "\n", "## Learning Goals:\n", "By the end of this tutorial, you will be able to:\n", "- get catalogs and images from NASA archives in the cloud where possible\n", "- measure fluxes at any location by running forced photometry using \"The Tractor\" \n", "- employ parallel processing to make this as fast as possible\n", "- cross match large catalogs\n", "- plot results\n", "\n", "## Introduction:\n", "This code performs photometry in an automated fashion at all locations in an input catalog on 4 bands of IRAC data from IRSA and 2 bands of Galex data from MAST. The resulting catalog is then cross-matched with a Chandra catalog from HEASARC to generate a multiband catalog to facilitate galaxy evolution studies.\n", "\n", "The code will run on 2 different science platforms and makes full use of multiple processors to optimize run time on large datasets.\n", "\n", "## Input:\n", "- RA and DEC within COSMOS catalog\n", "- desired catalog radius in arcminutes\n", "- mosaics of that region for IRAC and Galex\n", "\n", "## Output:\n", "- merged, multiband, science ready pandas dataframe\n", "- IRAC color color plot for identifying interesting populations\n", "\n", "## Non-standard Imports\n", "- `tractor` code which does the forced photometry from Lang et al., 2016\n", "- `astroquery` to interface with archives APIs\n", "- `astropy` to work with coordinates/units and data structures\n", "- `skimage` to work with the images\n", "\n", "\n", "## Authors:\n", "Jessica Krick, David Shupe, Marziye JafariYazani, Brigitta Sipőcz, Vandana Desai, Steve Groom, Troy Raen\n", "\n", "\n", "## Acknowledgements:\n", "Kristina Nyland for the workflow of the tractor wrapper.\\\n", "MAST, HEASARC, & IRSA Fornax teams" ] }, { "cell_type": "code", "execution_count": null, "id": "fdaec1e6", "metadata": {}, "outputs": [], "source": [ "#ensure all dependencies are installed\n", "!pip install -r requirements.txt" ] }, { "cell_type": "code", "execution_count": null, "id": "982b8779", "metadata": {}, "outputs": [], "source": [ "# standard lib imports\n", "\n", "import math\n", "import time\n", "import warnings\n", "import concurrent.futures\n", "import sys\n", "import os\n", "import re\n", "from typing import NamedTuple\n", "\n", "# Third party imports\n", "\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "from skimage.transform import rotate\n", "import pandas as pd\n", "import seaborn as sns\n", "import statsmodels\n", "import mpld3\n", "from tqdm import tqdm\n", "\n", "from astropy.nddata import Cutout2D\n", "from astropy.io import fits\n", "from astropy.table import Table\n", "from astropy.coordinates import SkyCoord\n", "import astropy.units as u\n", "\n", "from astroquery.ipac.irsa import Irsa\n", "from astroquery.heasarc import Heasarc\n", "from astroquery.mast import Observations\n", "import pyvo\n", "\n", "# Local code imports\n", "sys.path.append('code_src/')\n", "\n", "from display_images import display_images\n", "import cutout\n", "from exceptions import TractorError\n", "import photometry\n", "from plot_SED import plot_SED\n", "from nway_write_header import nway_write_header\n", "from photometry import Band\n", "from photometry import lookup_img_pair\n", "#from prepare_prf import prepare_prf\n", "\n", "# temporarily let the notebook start without tractor as dependency\n", "try:\n", " from find_nconfsources import find_nconfsources\n", "except ImportError:\n", " print(\"tractor is missing\")\n", " pass\n", "\n", "\n", "%matplotlib inline" ] }, { "cell_type": "markdown", "id": "475affc8", "metadata": {}, "source": [ "## 1. Retrieve Initial Catalog from IRSA\n", "- Automatically set up a catalog with ra, dec, photometric redshifts, fiducial band fluxes, & probability that it is a star \n", "- Catalog we are using is COSMOS2015 (Laigle et al. 2016) \n", "- Data exploration" ] }, { "cell_type": "code", "execution_count": null, "id": "d4f57a27", "metadata": {}, "outputs": [], "source": [ "#pull a COSMOS catalog from IRSA using pyVO\n", "\n", "#what is the central RA and DEC of the desired catalog\n", "coords = SkyCoord('150.01d 2.2d', frame='icrs') #COSMOS center acording to Simbad\n", "\n", "#how large is the search radius, in arcmin\n", "radius = 48.0 * u.arcmin \n", "\n", "#specify only select columns to limit the size of the catalog\n", "cols = [\n", " 'ra', 'dec', 'id', 'Ks_FLUX_APER2', 'Ks_FLUXERR_APER2', 'PHOTOZ', 'SPLASH_1_MAG',\n", " 'SPLASH_1_MAGERR', 'SPLASH_1_FLUX', 'SPLASH_1_FLUX_ERR', 'SPLASH_2_FLUX',\n", " 'SPLASH_2_FLUX_ERR', 'SPLASH_3_FLUX', 'SPLASH_3_FLUX_ERR', 'SPLASH_4_FLUX',\n", " 'SPLASH_4_FLUX_ERR', 'FLUX_GALEX_NUV', 'FLUX_GALEX_FUV', 'FLUX_CHANDRA_05_2',\n", " 'FLUX_CHANDRA_2_10', 'FLUX_CHANDRA_05_10', 'ID_CHANDRA09 ', 'type', 'r_MAG_AUTO',\n", " 'r_MAGERR_AUTO', 'FLUX_24', 'FLUXERR_24', 'MAG_GALEX_NUV', 'MAGERR_GALEX_NUV',\n", " 'MAG_GALEX_FUV', 'MAGERR_GALEX_FUV'\n", "]\n", "\n", "tap = pyvo.dal.TAPService('https://irsa.ipac.caltech.edu/TAP')\n", "result = tap.run_sync(\"\"\"\n", " SELECT {}\n", " FROM cosmos2015\n", " WHERE CONTAINS(POINT('ICRS',ra, dec), CIRCLE('ICRS',{}, {}, {}))=1\n", " \"\"\".format(','.join(cols), coords.ra.value, coords.dec.value, radius.to(u.deg).value))\n", "cosmos_table = result.to_table()\n", "\n", "\n", "print(\"Number of objects: \", len(cosmos_table))" ] }, { "cell_type": "markdown", "id": "082f3bce", "metadata": {}, "source": [ "### 1a. Filter Catalog\n", "- if desired could filter the initial catalog to only include desired sources" ] }, { "cell_type": "code", "execution_count": null, "id": "1e15d6a1", "metadata": {}, "outputs": [], "source": [ "#an example of how to filter the catalog to \n", "#select those rows with either chandra fluxes or Galex NUV fluxes\n", "\n", "#example_table = cosmos_table[(cosmos_table['flux_chandra_05_10']> 0) | (cosmos_table['flux_galex_fuv'] > 0)]\n" ] }, { "cell_type": "markdown", "id": "a8e55ad6", "metadata": {}, "source": [ "## 2. Retrieve Image Datasets from the Cloud" ] }, { "cell_type": "markdown", "id": "b02092e4", "metadata": {}, "source": [ "#### Use the fornax cloud access API to obtain the IRAC data from the IRSA S3 bucket. \n", "\n", "Details here may change as the prototype code is being added to the appropriate libraries, as well as the data holding to the appropriate NGAP storage as opposed to IRSA resources." ] }, { "cell_type": "code", "execution_count": null, "id": "eed28ee8", "metadata": {}, "outputs": [], "source": [ "# Temporary solution, remove when the fornax API is added to the image\n", "# This relies on the assumption that https://github.com/fornax-navo/fornax-cloud-access-API is being cloned to this environment. \n", "# If it's not, then run a ``git clone https://github.com/fornax-navo/fornax-cloud-access-API --depth=1`` from a terminal at the highest directory root.\n", "# You may need to update the fork if you forked it in the past\n", "\n", "import os\n", "if not os.path.exists('../../fornax-cloud-access-API'):\n", " ! git clone https://github.com/fornax-navo/fornax-cloud-access-API --depth=1 ../../fornax-cloud-access-API" ] }, { "cell_type": "code", "execution_count": null, "id": "cddaf0c3", "metadata": {}, "outputs": [], "source": [ "sys.path.append('../../fornax-cloud-access-API')\n", "\n", "import pyvo\n", "import fornax" ] }, { "cell_type": "code", "execution_count": null, "id": "9d26add7", "metadata": {}, "outputs": [], "source": [ "# Getting the COSMOS address from the registry to follow PyVO user case approach. We could hardwire it.\n", "image_services = pyvo.regsearch(servicetype='image')\n", "irsa_cosmos = [s for s in image_services if 'irsa' in s.ivoid and 'cosmos' in s.ivoid][0]\n", "\n", "# The search returns 11191 entries, but unfortunately we cannot really filter efficiently in the query\n", "# itself (https://irsa.ipac.caltech.edu/applications/Atlas/AtlasProgramInterface.html#inputparam)\n", "# to get only the Spitzer IRAC results from COSMOS as a mission. We will do the filtering in a next step before download.\n", "cosmos_results = irsa_cosmos.search(coords).to_table()\n", "\n", "spitzer = cosmos_results[cosmos_results['dataset'] == 'IRAC']" ] }, { "cell_type": "code", "execution_count": null, "id": "cfb4de4f", "metadata": {}, "outputs": [], "source": [ "# Temporarily add the cloud_access metadata to the Atlas response. \n", "# This dataset has limited acces, thus 'region' should be used instead of 'open'.\n", "# S3 access should be available from the daskhub and those who has their IRSA token set up.\n", "\n", "fname = spitzer['fname']\n", "spitzer['cloud_access'] = [(f'{{\"aws\": {{ \"bucket_name\": \"irsa-mast-tike-spitzer-data\",'\n", " f' \"region\": \"us-east-1\",'\n", " f' \"access\": \"restricted\",'\n", " f' \"key\": \"data/COSMOS/{fn}\" }} }}') for fn in fname]" ] }, { "cell_type": "code", "execution_count": null, "id": "52857198", "metadata": {}, "outputs": [], "source": [ "# Adding function to download multiple files using the fornax API. \n", "# Requires https://github.com/fornax-navo/fornax-cloud-access-API/pull/4\n", "def fornax_download(data_table, data_subdirectory, access_url_column='access_url',\n", " fname_filter=None, verbose=False):\n", " working_dir = os.getcwd()\n", " \n", " data_directory = f\"data/{data_subdirectory}\"\n", " os.chdir(data_directory)\n", " for row in data_table:\n", " if fname_filter is not None and fname_filter not in row['fname']:\n", " continue\n", " handler = fornax.get_data_product(row, 'aws', access_url_column=access_url_column, verbose=verbose)\n", " temp_file = handler.download()\n", " # on-prem download returns temp file path, s3 download downloads file\n", " if temp_file:\n", " os.rename(temp_file, os.path.basename(row['fname']))\n", " \n", " os.chdir(working_dir)" ] }, { "cell_type": "code", "execution_count": null, "id": "beae3a47", "metadata": {}, "outputs": [], "source": [ "fornax_download(spitzer, access_url_column='sia_url', fname_filter='go2_sci', \n", " data_subdirectory='IRAC', verbose=False)" ] }, { "cell_type": "markdown", "id": "2af22c0b", "metadata": {}, "source": [ "#### Use IVOA image search and Fornax download to obtain Galex from the MAST archive" ] }, { "cell_type": "code", "execution_count": null, "id": "26069f47", "metadata": {}, "outputs": [], "source": [ "#the Galex mosaic of COSMOS is broken into 4 seperate images\n", "#need to know which Galex image the targets are nearest to.\n", "#make a new column in dataframe which figures this out\n", "\n", "#four centers for 1, 2, 3, 4 are\n", "ra_center=[150.369,150.369,149.869,149.869]\n", "dec_center=[2.45583,1.95583,2.45583,1.95583]\n", "\n", "galex = SkyCoord(ra = ra_center*u.degree, dec = dec_center*u.degree)\n", "catalog = SkyCoord(ra = cosmos_table['ra'], dec = cosmos_table['dec'])\n", "#idx, d2d, d3d = match_coordinates_sky(galex, catalog) #only finds the nearest one\n", "#idx, d2d, d3d = galex.match_to_catalog_sky(catalog) #only finds the nearest one\n", "\n", "cosmos_table['COSMOS_01'] = galex[0].separation(catalog)\n", "cosmos_table['COSMOS_02'] = galex[1].separation(catalog)\n", "cosmos_table['COSMOS_03'] = galex[2].separation(catalog)\n", "cosmos_table['COSMOS_04'] = galex[3].separation(catalog)\n", "\n", "#convert to pandas\n", "df = cosmos_table.to_pandas()\n", "\n", "#which row has the minimum value of distance to the galex images\n", "df['galex_image'] = df[['COSMOS_01','COSMOS_02','COSMOS_03','COSMOS_04']].idxmin(axis = 1)" ] }, { "cell_type": "code", "execution_count": null, "id": "ac10566d", "metadata": {}, "outputs": [], "source": [ "# 76k with 15arcmin diameter IRAC images\n", "df.describe()" ] }, { "cell_type": "code", "execution_count": null, "id": "9961ddec", "metadata": {}, "outputs": [], "source": [ "# Assign the endpoint for the MAST Galex Simple Image Access service.\n", "galex_sia_url = 'https://mast.stsci.edu/portal_vo/Mashup/VoQuery.asmx/SiaV1?MISSION=GALEX&'\n", "\n", "# Query the Galex SIA service using the COSMOS coordinates defined above.\n", "query_result = pyvo.dal.sia.search(galex_sia_url, pos=coords, size=0.0)\n", "galex_image_products = query_result.to_table()\n", "\n", "# Filter the products so we only download SCIENCE products with exposure time > 40000\n", "filtered_products = galex_image_products[(galex_image_products['timExposure'] > 40000.0)]\n", "filtered_products = filtered_products[(filtered_products['productType'] == \"SCIENCE\")]\n", "\n", "# The column containing the on-prem access to fall back on if cloud is unavailable\n", "access_url_column = query_result.fieldname_with_ucd('VOX:Image_AccessReference')\n", "\n", "# Download filtered products from AWS\n", "download_subdir = 'Galex'\n", "fornax_download(filtered_products, access_url_column=access_url_column, data_subdirectory=download_subdir)" ] }, { "cell_type": "code", "execution_count": null, "id": "13aca8e3", "metadata": {}, "outputs": [], "source": [ "# Get the GALEX sky background fits files in addition to the mosaics\n", "skybg_products = []\n", "skybkg_pattern = re.compile(r\"COSMOS_0[1-4]-[fn]d-skybg\")\n", "\n", "for row in galex_image_products:\n", " if skybkg_pattern.search(row['name']):\n", " skybg_products.append(row)\n", "\n", "# Download skybg products from AWS\n", "download_subdir = 'Galex'\n", "fornax_download(skybg_products, access_url_column=access_url_column, data_subdirectory=download_subdir)" ] }, { "cell_type": "code", "execution_count": null, "id": "faa3e0df", "metadata": {}, "outputs": [], "source": [ "#make sure there aren't any troublesome rows in the catalog\n", "#are there missing values in any rows?\n", "df.isna().sum()\n", "\n", "#don't mind that there are missing values for some of the fluxes\n", "#The rest of the rows are complete" ] }, { "cell_type": "code", "execution_count": null, "id": "cb00427d", "metadata": {}, "outputs": [], "source": [ "#out of curiosity how many of each type of source are in this catalog\n", "#Type: 0 = galaxy, 1 = star, 2 = X-ray source, -9 is failure to fit\n", "df.type.value_counts()" ] }, { "cell_type": "markdown", "id": "f74c5ee7", "metadata": {}, "source": [ "## 3. Run Forced Photometry\n", "- Calculate a flux at a given position in 2 IRAC and 2 Galex bands" ] }, { "cell_type": "markdown", "id": "e22c2c61", "metadata": {}, "source": [ "### 3a. Setup\n", "- initialize data frame columns to hold the results\n", "- collect the parameters for each band/channel\n", "- collect the input images" ] }, { "cell_type": "code", "execution_count": null, "id": "931c55fc", "metadata": {}, "outputs": [], "source": [ "# initialize columns in data frame for photometry results\n", "cols = [\"ch1flux\", \"ch1flux_unc\", \"ch2flux\", \"ch2flux_unc\", \"ch3flux\", \"ch3flux_unc\",\n", " \"ch4flux\", \"ch4flux_unc\", \"ch5flux\", \"ch5flux_unc\", \"ch6flux\", \"ch6flux_unc\"]\n", "df[cols] = 0.0\n", "\n", "# list to collect all the bands\n", "all_bands = []" ] }, { "cell_type": "code", "execution_count": null, "id": "b6a3712e", "metadata": {}, "outputs": [], "source": [ "# IRAC channels\n", "\n", "irac_band_indexes = [\n", " 0, # ch1\n", " 1, # ch2\n", " 2, # ch3\n", " 3, # ch4\n", "]\n", "\n", "irac_fluxconversion = (1E12) / (4.254517E10) * (0.6) *(0.6)\n", "\n", "irac_mosaic_pix_scale = 0.6\n", "\n", "irac_cutout_width = 10 # in arcseconds, taken from Nyland et al. 2017\n", "\n", "irac_prfs = [\n", " fits.open('data/IRAC/PRF_IRAC_ch1.fits')[0].data,\n", " fits.open('data/IRAC/PRF_IRAC_ch2.fits')[0].data,\n", " fits.open('data/IRAC/PRF_IRAC_ch3.fits')[0].data,\n", " fits.open('data/IRAC/PRF_IRAC_ch4.fits')[0].data,\n", "]\n", "\n", "# zip parameters for each band into a container and append to the master list\n", "irac_bands = [\n", " Band(\n", " idx, prf, irac_cutout_width, irac_fluxconversion, irac_mosaic_pix_scale\n", " )\n", " for idx, prf in zip(irac_band_indexes, irac_prfs)\n", "]\n", "all_bands += irac_bands" ] }, { "cell_type": "code", "execution_count": null, "id": "85ead361", "metadata": {}, "outputs": [], "source": [ "# GALEX bands\n", "\n", "galex_band_indexes = [\n", " 4, # nuv\n", " 5, # fuv\n", "]\n", "\n", "galex_cutout_width = 40\n", "\n", "galex_fluxconversions = [\n", " 3.373E1, # uJy. fudging this to make the numbers bigger for plotting later\n", " 1.076E2, # uJy. fudging this to make the numbers bigger for plotting later\n", "]\n", "\n", "galex_mosaic_pix_scale = 1.5\n", "\n", "prf_nuv = fits.open(\"data/Galex/PSFnuv_faint.fits\")[0].data\n", "prf_fuv = fits.open(\"data/Galex/PSFfuv.fits\")[0].data\n", "prf_nuv = prf_nuv[0:119, 0:119]\n", "prf_fuv = prf_fuv[0:119, 0:119]\n", "\n", "#these are much larger than the cutouts we are using, so only keep the central region which is the size of our cutouts\n", "ngalex_pix = galex_cutout_width / galex_mosaic_pix_scale\n", "prf_cen = int(60)\n", "prf_nuv = prf_nuv[(prf_cen - int(ngalex_pix / 2) - 1) : (prf_cen + int(ngalex_pix / 2)),\n", " (prf_cen - int(ngalex_pix / 2) - 1) : (prf_cen + int(ngalex_pix / 2))]\n", "prf_fuv = prf_fuv[(prf_cen - int(ngalex_pix / 2) - 1) : (prf_cen + int(ngalex_pix / 2)),\n", " (prf_cen - int(ngalex_pix / 2) - 1) : (prf_cen + int(ngalex_pix / 2))]\n", "galex_prfs = [prf_nuv, prf_fuv]\n", "\n", "# zip parameters for each band into a container and append to the master list\n", "galex_bands = [\n", " Band(\n", " idx, prf, galex_cutout_width, flux_conv, galex_mosaic_pix_scale\n", " )\n", " for idx, prf, flux_conv in zip(galex_band_indexes, galex_prfs, galex_fluxconversions)\n", "]\n", "all_bands += galex_bands" ] }, { "cell_type": "code", "execution_count": null, "id": "0642f2ac", "metadata": {}, "outputs": [], "source": [ "#Collect input images\n", "# collect the files in pairs: (science image, sky-background image)\n", "# if the same file should be used for both, just send it once\n", "sci_bkg_pairs = [\n", " # IRAC. use the science image to calculate the background\n", " ('data/IRAC/irac_ch1_go2_sci_10.fits', ),\n", " ('data/IRAC/irac_ch2_go2_sci_10.fits', ),\n", " ('data/IRAC/irac_ch3_go2_sci_10.fits', ),\n", " ('data/IRAC/irac_ch4_go2_sci_10.fits', ),\n", " # GALEX. calculate the background from a dedicated file\n", " ('data/Galex/COSMOS_01-nd-int.fits.gz', 'data/Galex/COSMOS_01-nd-skybg.fits.gz'),\n", " ('data/Galex/COSMOS_01-fd-int.fits.gz', 'data/Galex/COSMOS_01-fd-skybg.fits.gz'),\n", " ('data/Galex/COSMOS_02-nd-int.fits.gz', 'data/Galex/COSMOS_02-nd-skybg.fits.gz'),\n", " ('data/Galex/COSMOS_02-fd-int.fits.gz', 'data/Galex/COSMOS_02-fd-skybg.fits.gz'),\n", " ('data/Galex/COSMOS_03-nd-int.fits.gz', 'data/Galex/COSMOS_03-nd-skybg.fits.gz'),\n", " ('data/Galex/COSMOS_03-fd-int.fits.gz', 'data/Galex/COSMOS_03-fd-skybg.fits.gz'),\n", " ('data/Galex/COSMOS_04-nd-int.fits.gz', 'data/Galex/COSMOS_04-nd-skybg.fits.gz'),\n", " ('data/Galex/COSMOS_04-fd-int.fits.gz', 'data/Galex/COSMOS_04-fd-skybg.fits.gz'),\n", "]\n" ] }, { "cell_type": "markdown", "id": "b5db63b2", "metadata": {}, "source": [ "### 3b. Main Function to do the Forced Photometry" ] }, { "cell_type": "code", "execution_count": null, "id": "44a79795", "metadata": {}, "outputs": [], "source": [ "def calc_instrflux(band, ra, dec, stype, ks_flux_aper2, img_pair, df):\n", " \"\"\"\n", " Calculate single-band instrumental fluxes and uncertainties at the given ra, dec \n", " using tractor.\n", "\n", " Parameters:\n", " -----------\n", " band : `Band`\n", " Collection of parameters for a single band. \n", " A `Band` is a named tuple with the following attributes:\n", " idx : int\n", " Identifier for the band/channel.\n", " (integer in [0, 1, 2, 3, 4, 5] for the four IRAC bands and two Galex bands)\n", " prf : np.ndarray\n", " Point spread function for the band/channel.\n", " cutout_width : int\n", " width of desired cutout in arcseconds\n", " flux_conv : float\n", " factor used to convert tractor result to microjanskies\n", " mosaic_pix_scale : float\n", " Pixel scale of the image\n", " ra, dec : float\n", " celestial coordinates for measuring photometry\n", " stype : int\n", " 0, 1, 2, -9 for star, galaxy, x-ray source\n", " ks_flux_aper_2 : float\n", " flux in aperture 2\n", " img_pair : tuple\n", " Pair of images for science and background respectively.\n", " If the tuple only contains one element it will serve double duty.\n", " A tuple element can be a `fits.ImageHDU` or the path to a FITS file as a `str`.\n", " df : pd.DataFrame\n", " Source catalog.\n", " Previous arguments (ra, dec, stype, ks_flux_aper_2) come from a single row of this df.\n", " However, we must also pass the entire dataframe in order to find nearby sources which are possible contaminates.\n", "\n", " Returns:\n", " --------\n", " outband : int\n", " reflects the input band index for identification purposes\n", " flux : float\n", " Measured flux in microJansky.\n", " NaN if the forced photometery failed.\n", " flux_unc : float\n", " Flux uncertainty in microJansky, calculated from the tractor results.\n", " NaN if the forced photometery failed or if tractor didn't report a flux variance.\n", " \"\"\"\n", " # cutout a small region around the object of interest\n", " subimage, bkgsubimage, x1, y1, subimage_wcs = cutout.extract_pair(\n", " ra, dec, img_pair=img_pair, cutout_width=band.cutout_width, mosaic_pix_scale=band.mosaic_pix_scale\n", " )\n", " \n", " # find nearby sources that are possible contaminants\n", " objsrc, nconfsrcs = find_nconfsources(\n", " ra, dec, stype, ks_flux_aper2, x1, y1, band.cutout_width, subimage_wcs, df\n", " )\n", "\n", " # estimate the background\n", " skymean, skynoise = photometry.calc_background(bkgsubimage=bkgsubimage)\n", "\n", " # do the forced photometry\n", " # if tractor fails to converge, just return NaNs\n", " try:\n", " flux_var = photometry.run_tractor(\n", " subimage=subimage, prf=band.prf, objsrc=objsrc, skymean=skymean, skynoise=skynoise\n", " )\n", " except TractorError:\n", " return (band.idx, np.nan, np.nan)\n", "\n", " # convert the results\n", " microJy_flux, microJy_unc = photometry.interpret_tractor_results(\n", " flux_var=flux_var, flux_conv=band.flux_conv, objsrc=objsrc, nconfsrcs=nconfsrcs\n", " )\n", "\n", " return (band.idx, microJy_flux, microJy_unc)" ] }, { "cell_type": "markdown", "id": "7aaddfba", "metadata": {}, "source": [ "### 3c. Calculate Forced Photometry with Straightforward but Slow Method\n", "- no longer in use but keeping to demonstrate this capability\\\n", "as well as the increase in speed when going to parallelization" ] }, { "cell_type": "raw", "id": "6618b40a", "metadata": {}, "source": [ "%%time\n", "#do the calculation without multiprocessing for benchmarking\n", "\n", "#make a copy for parallel computation\n", "pl_df = df.copy(deep=True)\n", "\n", "t0 = time.time()\n", "#for each object\n", "for row in df.itertuples():\n", " #for each band\n", " for band in range(6):\n", " #measure the flux with tractor\n", " outband, flux, unc = calc_instrflux(band, row.ra, row.dec, row.type, row.ks_flux_aper2)\n", " #put the results back into the dataframe\n", " df.loc[row.Index, 'ch{:d}flux'.format(outband+1)] = flux\n", " df.loc[row.Index, 'ch{:d}flux_unc'.format(outband+1)] = unc\n", " #print(row.ra, row.dec, row.type, row.ks_flux_aper2, band+1,\n", " # outband, flux, unc)\n", "t1 = time.time()\n", "\n", "\n", "#10,000 sources took 1.5 hours with this code on the IPAC SP" ] }, { "cell_type": "markdown", "id": "60533433", "metadata": {}, "source": [ "### 3d. Calculate Forced Photometry - Parallelization\n", "- Parallelization: we can either interate over the rows of the dataframe and run the four bands in parallel; or we could zip together the row index, band, ra, dec," ] }, { "cell_type": "code", "execution_count": null, "id": "de32f88d", "metadata": {}, "outputs": [], "source": [ "paramlist = []\n", "for row in df.itertuples():\n", " for band in all_bands:\n", " img_pair = lookup_img_pair(sci_bkg_pairs, band.idx, row.galex_image) # file paths only\n", " paramlist.append(\n", " [row.Index, band, row.ra, row.dec, row.type, row.ks_flux_aper2, img_pair, df]\n", " )\n", "print('paramlist: ', len(paramlist))" ] }, { "cell_type": "code", "execution_count": null, "id": "e82a0353", "metadata": {}, "outputs": [], "source": [ "#proove we can do this for one object\n", "calc_instrflux(*paramlist[0][1:])\n", "\n", "#same thing, different syntax\n", "# calc_instrflux(paramlist[0][1], paramlist[0][2], paramlist[0][3], paramlist[0][4], paramlist[0][5], paramlist[0][6])" ] }, { "cell_type": "code", "execution_count": null, "id": "6a3d0411", "metadata": {}, "outputs": [], "source": [ "#wrapper to measure the photometry on a single object, single band\n", "def calculate_flux(args):\n", " \"\"\"Calculate flux.\"\"\"\n", " val = calc_instrflux(*args[1:])\n", " # add simple logging\n", " if (args[0] % 100) == 0 and val[0] == 0:\n", " with open('output/output.log', 'a') as fp: fp.write(f'{args[0]}\\n')\n", " return(args[0], val)" ] }, { "cell_type": "code", "execution_count": null, "id": "09723553", "metadata": {}, "outputs": [], "source": [ "# if results were previously saved to this location, load them\n", "# else start a pool of workers to calculate results in parallel, and save them here\n", "fname = f'output/results_{radius.value}.npz'\n", "\n", "if os.path.exists(fname):\n", " results = np.load(fname, allow_pickle=True)['results']\n", "\n", "else:\n", " from multiprocessing import Pool\n", " t0 = time.time()\n", " with open('output/output.log', 'w') as fp: fp.write('')\n", " with Pool() as pool:\n", " results = pool.map(calculate_flux, paramlist)\n", " dtime = time.time() - t0\n", " np.savez(fname, results=results)\n", " print(f'Parallel calculation took {dtime} seconds')" ] }, { "cell_type": "code", "execution_count": null, "id": "0fd87078", "metadata": {}, "outputs": [], "source": [ "## put the results into the main daraframe\n", "for res in results:\n", " idx,(ich, val, err) = res\n", " df.loc[idx, f'ch{ich+1}flux'] = val\n", " df.loc[idx, f'ch{ich+1}flux_unc'] = err" ] }, { "cell_type": "code", "execution_count": null, "id": "0107b549", "metadata": {}, "outputs": [], "source": [ "#Count the number of non-zero ch1 fluxes\n", "print('Parallel calculation: number of ch1 fluxes filled in =',\n", " np.sum(df.ch1flux > 0))" ] }, { "cell_type": "markdown", "id": "fd9036a2", "metadata": {}, "source": [ "### 3e. Cleanup" ] }, { "cell_type": "code", "execution_count": null, "id": "e7075f48", "metadata": {}, "outputs": [], "source": [ "#had to call the galex flux columns ch5 and ch6\n", "#fix that by renaming them now\n", "cols = {'ch5flux':'nuvflux', 'ch5flux_unc':'nuvflux_unc','ch6flux':'fuvflux', 'ch6flux_unc':'fuvflux_unc'}\n", "df.rename(columns=cols, inplace = True)" ] }, { "cell_type": "code", "execution_count": null, "id": "678acb95", "metadata": {}, "outputs": [], "source": [ "#When doing a large run of a large area, save the dataframe with the forced photometry \n", "# so we don't have to do the forced photometry every time\n", "\n", "df.to_pickle(f'output/COSMOS_{radius.value}arcmin.pkl')" ] }, { "cell_type": "code", "execution_count": null, "id": "91a24a05", "metadata": {}, "outputs": [], "source": [ "#If you are not runnig the forced photometry, then read in the catalog from a previous run\n", "\n", "#df = pd.read_pickle('output/COSMOS_48.0arcmin.pkl')" ] }, { "cell_type": "markdown", "id": "29d7b40c", "metadata": {}, "source": [ "### 3f. Plot to Confirm our Photometry Results \n", "- compare to published COSMOS 2015 catalog" ] }, { "cell_type": "code", "execution_count": null, "id": "65c183c4", "metadata": {}, "outputs": [], "source": [ "%%time\n", "#plot tractor fluxes vs. catalog splash fluxes\n", "#should see a straightline with a slope of 1\n", "\n", "#setup to plot\n", "fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2)\n", "fluxmax = 200\n", "ymax = 100\n", "xmax = 100\n", "#ch1 \n", "#first shrink the dataframe to only those rows where I have tractor photometry \n", "df_tractor = df[(df.splash_1_flux> 0) & (df.splash_1_flux < fluxmax)] #200\n", "#sns.regplot(data = df_tractor, x = \"splash_1_flux\", y = \"ch1flux\", ax = ax1, robust = True)\n", "sns.scatterplot(data = df_tractor, x = \"splash_1_flux\", y = \"ch1flux\", ax = ax1)\n", "\n", "#add a diagonal line with y = x\n", "lims = [\n", " np.min([ax1.get_xlim(), ax1.get_ylim()]), # min of both axes\n", " np.max([ax1.get_xlim(), ax1.get_ylim()]), # max of both axes\n", "]\n", "\n", "# now plot both limits against eachother\n", "ax1.plot(lims, lims, 'k-', alpha=0.75, zorder=0)\n", "ax1.set(xlabel = 'COSMOS 2015 flux ($\\mu$Jy)', ylabel = 'tractor flux ($\\mu$Jy)', title = 'IRAC 3.6')\n", "ax1.set_ylim([0, ymax])\n", "ax1.set_xlim([0, xmax])\n", "\n", "\n", "#ch2 \n", "#first shrink the dataframe to only those rows where I have tractor photometry \n", "df_tractor = df[(df.splash_2_flux> 0) & (df.splash_2_flux < fluxmax)]\n", "#sns.regplot(data = df_tractor, x = \"splash_2_flux\", y = \"ch2flux\", ax = ax2, robust = True)\n", "sns.scatterplot(data = df_tractor, x = \"splash_2_flux\", y = \"ch2flux\", ax = ax2)\n", "\n", "#add a diagonal line with y = x\n", "lims = [\n", " np.min([ax2.get_xlim(), ax2.get_ylim()]), # min of both axes\n", " np.max([ax2.get_xlim(), ax2.get_ylim()]), # max of both axes\n", "]\n", "\n", "# now plot both limits against eachother\n", "ax2.plot(lims, lims, 'k-', alpha=0.75, zorder=0)\n", "ax2.set(xlabel = 'COSMOS 2015 flux ($\\mu$Jy)', ylabel = 'tractor flux ($\\mu$Jy)', title = 'IRAC 4.5')\n", "ax2.set_ylim([0, ymax])\n", "ax2.set_xlim([0, xmax])\n", "\n", "\n", "#ch3 \n", "#first shrink the dataframe to only those rows where I have tractor photometry\n", "df_tractor = df[(df.splash_3_flux> 0) & (df.splash_3_flux < fluxmax)]\n", "\n", "#sns.regplot(data = df_tractor, x = \"splash_3_flux\", y = \"ch3flux\", ax = ax3, robust = True)\n", "sns.scatterplot(data = df_tractor, x = \"splash_3_flux\", y = \"ch3flux\", ax = ax3)\n", "\n", "#add a diagonal line with y = x\n", "lims = [\n", " np.min([ax3.get_xlim(), ax3.get_ylim()]), # min of both axes\n", " np.max([ax3.get_xlim(), ax3.get_ylim()]), # max of both axes\n", "]\n", "\n", "# now plot both limits against eachother\n", "ax3.plot(lims, lims, 'k-', alpha=0.75, zorder=0)\n", "ax3.set(xlabel = 'COSMOS 2015 flux ($\\mu$Jy)', ylabel = 'tractor flux ($\\mu$Jy)', title = 'IRAC 5.8')\n", "ax3.set_ylim([0, ymax])\n", "ax3.set_xlim([0, xmax])\n", "\n", "\n", "#ch4 \n", "#first shrink the dataframe to only those rows where I have tractor photometry \n", "df_tractor = df[(df.splash_4_flux> 0) & (df.splash_4_flux < fluxmax)]\n", "\n", "#sns.regplot(data = df_tractor, x = \"splash_4_flux\", y = \"ch4flux\", ax = ax4, robust = True)\n", "sns.scatterplot(data = df_tractor, x = \"splash_4_flux\", y = \"ch4flux\", ax = ax4)\n", "\n", "#add a diagonal line with y = x\n", "lims = [\n", " np.min([ax4.get_xlim(), ax4.get_ylim()]), # min of both axes\n", " np.max([ax4.get_xlim(), ax4.get_ylim()]), # max of both axes\n", "]\n", "\n", "# now plot both limits against eachother\n", "ax4.plot(lims, lims, 'k-', alpha=0.75, zorder=0)\n", "ax4.set(xlabel = 'COSMOS 2015 flux ($\\mu$Jy)', ylabel = 'tractor flux ($\\mu$Jy)', title = 'IRAC 8.0')\n", "ax4.set_ylim([0, ymax])\n", "ax4.set_xlim([0, xmax])\n", "\n", "\n", "plt.tight_layout()\n", "\n", "fig.subplots_adjust( hspace=0.5)\n", "fig.set_size_inches(8, 12)\n", "\n", "#plt.savefig('output/flux_comparison.png')" ] }, { "cell_type": "markdown", "id": "2ca3d900", "metadata": {}, "source": [ "Tractor is working for IRAC; Comparison of tractor derived fluxes with COSMOS 2015 fluxes for all four Spitzer IRAC channels. Blue points represent each object from the subset of the COSMOS 2015 catalog. The blue line is a linear regression robust fit to the data with uncertainties shown as the light blue wedge. The black line is a y = x line plotted to guide the eye." ] }, { "cell_type": "markdown", "id": "7cafae11", "metadata": {}, "source": [ "## 4. Cross Match our New Photometry Catalog with an X-ray archival Catalog\n", "We are using `nway` as the tool to do the cross match (Salvato et al. 2017).\n", "`nway` expects input as two fits table files and outputs a third table file with all the possible matches and their probabilities of being the correct match. We then sort that catalog and take only the best matches to be the true matches." ] }, { "cell_type": "markdown", "id": "0bd009f3", "metadata": {}, "source": [ "### 4a. Retrieve the HEASARC Catalog" ] }, { "cell_type": "code", "execution_count": null, "id": "25306dad", "metadata": {}, "outputs": [], "source": [ "#first get an X-ray catalog from Heasarc\n", "heasarc = Heasarc()\n", "table = heasarc.query_mission_list()\n", "mask = (table['Mission'] ==\"CHANDRA\")\n", "chandratable = table[mask] \n", "\n", "#find out which tables exist on Heasarc\n", "#chandratable.pprint(max_lines = 200, max_width = 130)\n", "\n", "#want ccosmoscat\n", "mission = 'ccosmoscat'\n", "#coords already defined above where I pull the original COSMOS catalog\n", "ccosmoscat_rad = 1 #radius of chandra cosmos catalog\n", "ccosmoscat = heasarc.query_region(coords, mission=mission, radius='1 degree', resultmax = 5000, fields = \"ALL\")" ] }, { "cell_type": "markdown", "id": "a3a9241a", "metadata": {}, "source": [ "### 4b. Run `nway` to do the Cross-Match" ] }, { "cell_type": "code", "execution_count": null, "id": "d9f70a3a", "metadata": {}, "outputs": [], "source": [ "#setup:\n", "\n", "#astropy doesn't recognize capitalized units\n", "#so there might be some warnings here on writing out the file, but we can safely ignore those\n", "\n", "#need to make the chandra catalog into a fits table \n", "#and needs to include area of the survey.\n", "nway_write_header('data/Chandra/COSMOS_chandra.fits', 'CHANDRA', float(ccosmoscat_rad**2) )\n", "\n", "\n", "#also need to transform the main pandas dataframe into fits table for nway\n", "#make an index column for tracking later\n", "df['ID'] = range(1, len(df) + 1)\n", "\n", "#need this to be a fits table and needs to include area of the survey.\n", "rad_in_arcmin = radius.value #units attached to this are confusing nway down the line\n", "nway_write_header('data/multiband_phot.fits', 'OPT', float((2*rad_in_arcmin/60)**2) )\n" ] }, { "cell_type": "code", "execution_count": null, "id": "6a9a77f4", "metadata": {}, "outputs": [], "source": [ "# call nway\n", "!/home/jovyan/.local/bin/nway.py 'data/Chandra/COSMOS_chandra.fits' :ERROR_RADIUS 'data/multiband_phot.fits' 0.1 --out=data/Chandra/chandra_multiband.fits --radius 15 --prior-completeness 0.9\n", " \n", "#!/opt/conda/bin/nway.py 'data/Chandra/COSMOS_chandra.fits' :ERROR_RADIUS 'data/multiband_phot.fits' 0.1 --out=data/Chandra/chandra_multiband.fits --radius 15 --prior-completeness 0.9" ] }, { "cell_type": "code", "execution_count": null, "id": "e13906ce", "metadata": {}, "outputs": [], "source": [ "#Clean up the cross match results and merge them back into main pandas dataframe\n", "\n", "#read in the nway matched catalog\n", "xmatch = Table.read('data/Chandra/chandra_multiband.fits', hdu = 1)\n", "df_xmatch = xmatch.to_pandas()\n", "\n", "#manual suggests that p_i should be greater than 0.1 for a pure catalog.\n", "#The matched catalog has multiple optical associations for some of the XMM detections.\n", "#simplest thing to do is only keep match_flag = 1\n", "matched = df_xmatch.loc[(df_xmatch['p_i']>=0.1) & df_xmatch['match_flag']==1]\n", "\n", "#merge this info back into the df_optical dataframe.\n", "merged = pd.merge(df, matched, 'outer',left_on='ID', right_on = 'OPT_ID')\n", "\n", "#remove all the rows which start with \"OPT\" because they are duplications of the original catalog\n", "merged = merged.loc[:, ~merged.columns.str.startswith('OPT')]\n", "\n", "#somehow the matching is giving negative fluxes in the band where there is no detection \n", "#if there is a detection in the other band\n", "#clean that up to make those negative fluxes = 0\n", "\n", "merged.loc[merged['flux_chandra_2_10'] < 0, 'flux_chandra_2_10'] = 0\n", "merged.loc[merged['flux_chandra_05_2'] < 0, 'flux_chandra_05_2'] = 0" ] }, { "cell_type": "code", "execution_count": null, "id": "989f9911", "metadata": {}, "outputs": [], "source": [ "#How many Chandra sources are there?\n", "#How many Galex sources are there?\n", "\n", "#make a new column which is a bool of existing chandra measurements\n", "merged['cosmos_chandra_detect'] = 0\n", "merged.loc[merged.flux_chandra_2_10 > 0,'cosmos_chandra_detect']=1\n", "\n", "#make one for Galex too\n", "merged['galex_detect'] = 0\n", "merged.loc[merged.flux_galex_nuv > 0,'galex_detect']=1\n", "\n", "#make chandra hardness ratio column:\n", "#hard = 'flux_chandra_2_10', soft = flux_chandra_05_2, HR = (H-S)/(H+S)\n", "merged['chandra_HR'] = (merged['flux_chandra_2_10'] - merged['flux_chandra_05_2']) / (merged['flux_chandra_2_10'] + merged['flux_chandra_05_2'])\n", "\n", "\n", "print('number of Chandra detections =',np.sum(merged.cosmos_chandra_detect > 0))\n", "print('number of Galex detections =',np.sum(merged.galex_detect > 0))" ] }, { "cell_type": "markdown", "id": "ab18392c", "metadata": {}, "source": [ "## 5. Plot Final Results\n", "- We want to understand something about populations based on their colors" ] }, { "cell_type": "code", "execution_count": null, "id": "63c295b2", "metadata": {}, "outputs": [], "source": [ "#IRAC color color plots akin to Lacy et al. 2004\n", "#overplot galex sources\n", "#overplot xray sources\n", "\n", "#first select on 24 micron \n", "merged_24 = merged[(merged.flux_24 >= 0) ] \n", "\n", "#negative Galex fluxes are causing problems, so set those to zero\n", "merged_24.loc[merged_24.fuvflux < 0, 'fuvflux'] = 0\n", "merged_24.loc[merged_24.nuvflux < 0, 'nuvflux'] = 0\n", "\n", "#make color columns\n", "merged_24['F5.8divF3.6'] = merged_24.ch3flux / merged_24.ch1flux\n", "merged_24['F8.0divF4.5'] = merged_24.ch4flux / merged_24.ch2flux\n", "\n", "#detected in all IRAC bands\n", "merged_allirac = merged_24[(merged_24['F8.0divF4.5'] > 0) & (merged_24['F5.8divF3.6'] > 0)]\n", "\n", "#plot all the points\n", "fig, ax = plt.subplots()\n", "sns.scatterplot(data = merged_allirac, x = 'F5.8divF3.6', y = 'F8.0divF4.5',\n", " ax = ax, alpha = 0.5, label = 'all')\n", "\n", "#plot only those points with Galex detections\n", "galex_detect = merged_allirac[merged_allirac.galex_detect > 0]\n", "sns.scatterplot(data = galex_detect, x = 'F5.8divF3.6', y = 'F8.0divF4.5',\n", " ax = ax, alpha = 0.5, label = 'Galex detect')\n", "\n", "#plot only those points with chandra detections\n", "chandra_detect = merged_allirac[merged_allirac.cosmos_chandra_detect > 0]\n", "sns.scatterplot(data = chandra_detect, x = 'F5.8divF3.6', y = 'F8.0divF4.5',\n", " ax = ax, label = 'Chandra detect')\n", "\n", "\n", "\n", "ax.set(xscale=\"log\", yscale=\"log\")\n", "ax.set_ylim([0.1, 10])\n", "ax.set_xlim([0.1, 10])\n", "\n", "ax.set(xlabel = 'log F5.8/F3.6', ylabel = 'log F8.0/F4.5')\n", "ax.legend(loc='lower right')\n", "plt.title('IRAC Color Color Plot')" ] }, { "cell_type": "markdown", "id": "8f25846e", "metadata": {}, "source": [ "This figure shows an IRAC color color plot akin to the seminal work by Lacy et al. 2004. Points are color coded for those with Galex UV detections and those with Chandra x-ray detections. Note that the different populations are seperating out in this color color space." ] }, { "cell_type": "code", "execution_count": null, "id": "daffe78c", "metadata": {}, "outputs": [], "source": [ "#UV IR color color plot akin to Bouquin et al. 2015\n", "fig, ax = plt.subplots()\n", "merged['FUV-NUV'] = merged.mag_galex_fuv - merged.mag_galex_nuv\n", "merged['NUV-3.6'] = merged.mag_galex_nuv - merged.splash_1_mag\n", "\n", "\n", "#plot all the points\n", "#sns.scatterplot(data = merged, x = 'NUV-3.6', y = 'FUV-NUV',\n", "# ax = ax, alpha = 0.5)\n", "\n", "#plot only those points with Galex detections\n", "galex_detect = merged[merged.galex_detect > 0]\n", "sns.kdeplot(data = galex_detect, x = 'NUV-3.6', y = 'FUV-NUV',\n", " ax = ax, fill = True, levels = 15)#scatterplot , alpha = 0.5\n", "\n", "#plot only those points with chandra detections\n", "#now with color coding Chandra sources by hardness ratio a la Moutard et al. 2020\n", "chandra_detect = merged[merged.cosmos_chandra_detect > 0]\n", "sns.scatterplot(data = chandra_detect, x = 'NUV-3.6', y = 'FUV-NUV',\n", " ax = ax, hue= 'chandra_HR',palette=\"flare\")\n", "\n", "#whew that legend for the hue is terrible\n", "#try making it into a colorbar outside the plot instead\n", "norm = plt.Normalize(merged['chandra_HR'].min(), merged['chandra_HR'].max())\n", "sm = plt.cm.ScalarMappable(cmap=\"flare\", norm=norm)\n", "sm.set_array([])\n", "\n", "# Remove the legend and add a colorbar\n", "ax.get_legend().remove()\n", "ax.figure.colorbar(sm)\n", "\n", "#ax.set(xscale=\"log\", yscale=\"log\")\n", "ax.set_ylim([-0.5, 3.5])\n", "ax.set_xlim([-1, 7])\n", "\n", "ax.set(xlabel = 'NUV - [3.6]', ylabel = 'FUV - NUV')\n", "#plt.legend([],[], frameon=False)\n", "\n", "#fig.savefig(\"output/color_color.png\")\n", "#mpld3.display(fig) " ] }, { "cell_type": "markdown", "id": "944ab7d8", "metadata": {}, "source": [ "We extend the works of Bouquin et al. 2015 and Moutard et al. 2020 by showing a GALEX - Spitzer color color diagram over plotted with Chandra detections. Blue galaxies in these colors are generated by O and B stars and so must currently be forming stars. We find a tight blue cloud in this color space identifying those star forming galaxies. Galaxies off of the blue cloud have had their star formation quenched, quite possibly by the existence of an AGN through removal of the gas reservoir required for star formation. Chandra detected galaxies host AGN, and while those are more limited in number, can be shown here to be a hosted by all kinds of galaxies, including quiescent galaxies which would be in the upper right of this plot. This likely implies that AGN are indeed involved in quenching star formation. Additionally, we show the Chandra hardness ratio (HR) color coded according to the vertical color bar on the right side of the plot. Those AGN with higher hardness ratios have their soft x-ray bands heavily obscured and appear to reside preferentially toward the quiescent galaxies." ] }, { "cell_type": "markdown", "id": "1a8c12c9", "metadata": {}, "source": [ "## References\n", "\n", "This work made use of:\n", "\n", "- Astroquery; Ginsburg et al., 2019, 2019AJ....157...98G\n", "\n", "- Astropy; Astropy Collaboration 2022, Astropy Collaboration 2018, Astropy Collaboration 2013, 2022ApJ...935..167A, 2018AJ....156..123A, 2013A&A...558A..33A\n", "\n", "- The Tractor; Lang et al. 2016, 2016AJ....151...36L\n", "\n", "- Nyland et al. 2017 , 2017ApJS..230....9N\n", "\n", "- Salvato et al. 2018, 2018MNRAS.473.4937S\n", "\n", "- Laigle et al. 2016, 2016ApJS..224...24L" ] }, { "cell_type": "code", "execution_count": null, "id": "f6857fa3", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "jupytext": { "text_representation": { "extension": ".md", "format_name": "myst", "format_version": 0.13, "jupytext_version": "1.15.2" } }, "kernelspec": { "display_name": "Python [conda env:tractor]", "language": "python", "name": "conda-env-tractor-py" }, "source_map": [ 12, 56, 61, 117, 124, 154, 159, 165, 169, 175, 186, 193, 206, 218, 239, 242, 246, 272, 277, 297, 311, 320, 324, 329, 336, 346, 379, 420, 441, 445, 521, 527, 550, 555, 566, 574, 585, 604, 612, 616, 620, 627, 634, 638, 643, 737, 741, 747, 751, 766, 770, 791, 798, 824, 843, 848, 891, 895, 936, 940, 958 ] }, "nbformat": 4, "nbformat_minor": 5 }