{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This notebook prepare the catalogues that will be analysed by CIGALE for SED fitting and physical parameter estimation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import os\n",
    "os.environ['LOG_LEVEL'] = 'INFO'\n",
    "\n",
    "from astropy.table import Table\n",
    "\n",
    "from herschelhelp.filters import correct_galactic_extinction\n",
    "from herschelhelp.external import convert_table_for_cigale"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "SUFFIX=\"20180613\"\n",
    "master_catalogue = Table.read(\"../../dmu32/dmu32_CDFS-SWIRE/data/CDFS-SWIRE_{}_cigale.fits\".format(SUFFIX))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2171051"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "len(master_catalogue)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Best sources\n",
    "\n",
    "Define a good far-IR measurement as:\n",
    "- an existing flux in the band;\n",
    "- the flag from XID+ must not be set;\n",
    "- the signal to noise ratio must be over 2."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "good = {}\n",
    "for band in ['pacs_green', 'pacs_red', 'spire_250', 'spire_350', 'spire_500']:\n",
    "    good[band] = (~np.isnan(master_catalogue['f_{}'.format(band)]) & \n",
    "                                   ~master_catalogue['flag_{}'.format(band)])\n",
    "    good[band][good[band]] &= (master_catalogue[good[band]]['f_{}'.format(band)] /\n",
    "                              master_catalogue[good[band]]['ferr_{}'.format(band)] >= 2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We will keep only sources with at leat 2 good far-IR measurements (we may actually use less sources are not all may have a redshift)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "combined_good = np.sum(list(good.values()), axis=0) >= 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of good sources: 11086\n"
     ]
    }
   ],
   "source": [
    "print(\"Number of good sources: {}\".format(np.sum(combined_good)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Only sources with at least two optical and at least two near infrared detections\n",
    "optnir = ((master_catalogue['flag_optnir_det'] == 3) \n",
    "          | (master_catalogue['flag_optnir_det'] == 7))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Main catalogue for CIGALE"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "best_catalogue = master_catalogue[combined_good].copy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Correction for galactic extinction\n",
    "best_catalogue = correct_galactic_extinction(best_catalogue, inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": true,
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:herschelhelp.external:For 3 sources, the band omegacam_u should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 3 sources, the band wfi_u should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 1 sources, the band acs_f435w should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n"
     ]
    }
   ],
   "source": [
    "# Convertion to CIGALE format\n",
    "best_catalogue = convert_table_for_cigale(best_catalogue, inplace=True, remove_zerofluxes=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Band selection\n",
    "\n",
    "We want to use only one filter for similar bands. We define an order of preference and set to NaN the flux in the lower prefered bands when a prefered band is available. Some band may have a 0 flux, we set there values to NaN."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "u_bands = [                      \"omegacam_u\"]\n",
    "g_bands = [           \"decam_g\", \"omegacam_g\", \"gpc1_g\"]\n",
    "r_bands = [           \"decam_r\", \"omegacam_r\", \"gpc1_r\", \"wfi_r\"]\n",
    "i_bands = [           \"decam_i\", \"omegacam_i\", \"gpc1_i\"]\n",
    "z_bands = [\"vista_z\", \"decam_z\", \"omegacam_z\", \"gpc1_z\",]\n",
    "y_bands = [\"vista_y\", \"decam_y\",               \"gpc1_y\" ]\n",
    "\n",
    "def remove_unneeded_fluxes(list_of_bands):\n",
    "    for band_idx, band in enumerate(list_of_bands[:-1]):\n",
    "        mask = ~np.isnan(best_catalogue[band])\n",
    "        for lower_band in list_of_bands[band_idx+1:]:\n",
    "            best_catalogue[lower_band][mask] = np.nan\n",
    "            best_catalogue[\"{}_err\".format(lower_band)][mask] = np.nan"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "remove_unneeded_fluxes(g_bands)\n",
    "remove_unneeded_fluxes(u_bands)\n",
    "remove_unneeded_fluxes(r_bands)\n",
    "remove_unneeded_fluxes(i_bands)\n",
    "remove_unneeded_fluxes(z_bands)\n",
    "remove_unneeded_fluxes(y_bands)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "best_catalogue.write(\"data_tmp/CDFS-SWIRE_cigale_best_extcor_{}.fits\".format(SUFFIX), overwrite=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# All optnir objects with photo-z"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of optnir sources with z-phot: 110939\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:herschelhelp.external:For 123 sources, the band omegacam_u should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 117 sources, the band wfi_u should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 68 sources, the band wfi_416nm should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 82 sources, the band acs_f435w should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 68 sources, the band wfi_b should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 60 sources, the band wfi_461nm should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 68 sources, the band omegacam_g should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 68 sources, the band decam_g should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 59 sources, the band wfi_485nm should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 68 sources, the band gpc1_g should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 54 sources, the band wfi_518nm should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 58 sources, the band wfi_v should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 39 sources, the band wfi_571nm should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 59 sources, the band acs_f606w should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 37 sources, the band wfi_604nm should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 52 sources, the band gpc1_r should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 43 sources, the band omegacam_r should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 40 sources, the band decam_r should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 29 sources, the band wfi_646nm should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 39 sources, the band wfi_r should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 15 sources, the band wfi_696nm should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 15 sources, the band gpc1_i should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 22 sources, the band omegacam_i should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 15 sources, the band acs_f775w should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 13 sources, the band decam_i should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 14 sources, the band acs_f814w should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n"
     ]
    }
   ],
   "source": [
    "best_catalogue = master_catalogue[optnir].copy()\n",
    "best_catalogue = best_catalogue[~np.isnan(best_catalogue[\"redshift\"])]\n",
    "print(\"Number of optnir sources with z-phot: {}\".format(len(best_catalogue)))\n",
    "# Correction for galactic extinction\n",
    "best_catalogue = correct_galactic_extinction(best_catalogue, inplace=True)\n",
    "# Convertion to CIGALE format\n",
    "os.environ['LOG_LEVEL'] = 'INFO'\n",
    "best_catalogue = convert_table_for_cigale(best_catalogue, inplace=True, remove_zerofluxes=True)\n",
    "remove_unneeded_fluxes(g_bands)\n",
    "remove_unneeded_fluxes(u_bands)\n",
    "remove_unneeded_fluxes(r_bands)\n",
    "remove_unneeded_fluxes(i_bands)\n",
    "remove_unneeded_fluxes(z_bands)\n",
    "remove_unneeded_fluxes(y_bands)\n",
    "best_catalogue.write(\"data_tmp/CDFS-SWIRE_cigale_optnir_extcor_{}.fits\".format(SUFFIX), overwrite=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Catalogue using spectroscopic redshift"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "best_catalogue = master_catalogue[optnir].copy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "best_catalogue.remove_column(\"redshift\")\n",
    "best_catalogue[\"zspec\"].name = \"redshift\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "best_catalogue = best_catalogue[~np.isnan(best_catalogue[\"redshift\"])]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of sources with z-spec: 37625\n"
     ]
    }
   ],
   "source": [
    "print(\"Number of sources with z-spec: {}\".format(len(best_catalogue)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Correction for galactic extinction\n",
    "best_catalogue = correct_galactic_extinction(best_catalogue, inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "collapsed": true
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "INFO:herschelhelp.external:For 240 sources, the band omegacam_u should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 219 sources, the band wfi_u should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 57 sources, the band wfi_416nm should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 115 sources, the band acs_f435w should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 59 sources, the band wfi_b should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 13 sources, the band wfi_461nm should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 59 sources, the band omegacam_g should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 59 sources, the band decam_g should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 10 sources, the band wfi_485nm should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 59 sources, the band gpc1_g should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 3 sources, the band wfi_518nm should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 7 sources, the band wfi_v should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n",
      "INFO:herschelhelp.external:For 10 sources, the band acs_f606w should not be used because it overlaps or is below the Lyman limit at the redshift of these sources. These fluxes were set to NaN.\n"
     ]
    }
   ],
   "source": [
    "# Convertion to CIGALE format\n",
    "os.environ['LOG_LEVEL'] = 'INFO'\n",
    "best_catalogue = convert_table_for_cigale(best_catalogue, inplace=True, remove_zerofluxes=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "remove_unneeded_fluxes(g_bands)\n",
    "remove_unneeded_fluxes(u_bands)\n",
    "remove_unneeded_fluxes(r_bands)\n",
    "remove_unneeded_fluxes(i_bands)\n",
    "remove_unneeded_fluxes(z_bands)\n",
    "remove_unneeded_fluxes(y_bands)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "best_catalogue.write(\"data_tmp/CDFS-SWIRE_cigale_optnir_extcor_zspec_{}.fits\".format(SUFFIX), overwrite=True)"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python (herschelhelp_internal)",
   "language": "python",
   "name": "helpint"
  },
  "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.7.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}