This notebook presents the merge of the various pristine catalogues to produce HELP mater catalogue on Lockman SWIRE.
from herschelhelp_internal import git_version
print("This notebook was run with herschelhelp_internal version: \n{}".format(git_version()))
%matplotlib inline
#%config InlineBackend.figure_format = 'svg'
import matplotlib.pyplot as plt
plt.rc('figure', figsize=(10, 6))
import os
import time
from astropy import units as u
from astropy.coordinates import SkyCoord
from astropy.table import Column, Table
import numpy as np
from pymoc import MOC
from herschelhelp_internal.masterlist import merge_catalogues, nb_merge_dist_plot, specz_merge
from herschelhelp_internal.utils import coords_to_hpidx, ebv, gen_help_id, inMoc
TMP_DIR = os.environ.get('TMP_DIR', "./data_tmp")
OUT_DIR = os.environ.get('OUT_DIR', "./data")
SUFFIX = os.environ.get('SUFFIX', time.strftime("_%Y%m%d"))
try:
os.makedirs(OUT_DIR)
except FileExistsError:
pass
wfc = Table.read("{}/INT-WFC.fits".format(TMP_DIR))
rcs = Table.read("{}/RCSLenS.fits".format(TMP_DIR))
ps1 = Table.read("{}/PS1.fits".format(TMP_DIR))
sparcs = Table.read("{}/SpARCS.fits".format(TMP_DIR))
dxs = Table.read("{}/UKIDSS-DXS.fits".format(TMP_DIR))
swire= Table.read("{}/SWIRE.fits".format(TMP_DIR))
servs = Table.read("{}/SERVS.fits".format(TMP_DIR))
uhs = Table.read("{}/UHS.fits".format(TMP_DIR))
We first merge the optical catalogues and then add the infrared ones: WFC, DXS, SpARCS, HSC, PS1, SERVS, SWIRE.
At every step, we look at the distribution of the distances separating the sources from one catalogue to the other (within a maximum radius) to determine the best cross-matching radius.
master_catalogue = wfc
master_catalogue['wfc_ra'].name = 'ra'
master_catalogue['wfc_dec'].name = 'dec'
nb_merge_dist_plot(
SkyCoord(master_catalogue['ra'], master_catalogue['dec']),
SkyCoord(rcs['rcs_ra'], rcs['rcs_dec'])
)
# Given the graph above, we use 0.8 arc-second radius
master_catalogue = merge_catalogues(master_catalogue, rcs, "rcs_ra", "rcs_dec", radius=0.8*u.arcsec)
nb_merge_dist_plot(
SkyCoord(master_catalogue['ra'], master_catalogue['dec']),
SkyCoord(ps1['ps1_ra'], ps1['ps1_dec'])
)
# Given the graph above, we use 0.8 arc-second radius
master_catalogue = merge_catalogues(master_catalogue, ps1, "ps1_ra", "ps1_dec", radius=0.8*u.arcsec)
nb_merge_dist_plot(
SkyCoord(master_catalogue['ra'], master_catalogue['dec']),
SkyCoord(sparcs['sparcs_ra'], sparcs['sparcs_dec'])
)
# Given the graph above, we use 0.8 arc-second radius
master_catalogue = merge_catalogues(master_catalogue, sparcs, "sparcs_ra", "sparcs_dec", radius=0.8*u.arcsec)
nb_merge_dist_plot(
SkyCoord(master_catalogue['ra'], master_catalogue['dec']),
SkyCoord(dxs['dxs_ra'], dxs['dxs_dec'])
)
# Given the graph above, we use 0.8 arc-second radius
master_catalogue = merge_catalogues(master_catalogue, dxs, "dxs_ra", "dxs_dec", radius=0.8*u.arcsec)
nb_merge_dist_plot(
SkyCoord(master_catalogue['ra'], master_catalogue['dec']),
SkyCoord(servs['servs_ra'], servs['servs_dec'])
)
# Given the graph above, we use 1 arc-second radius
master_catalogue = merge_catalogues(master_catalogue, servs, "servs_ra", "servs_dec", radius=1.*u.arcsec)
nb_merge_dist_plot(
SkyCoord(master_catalogue['ra'], master_catalogue['dec']),
SkyCoord(swire['swire_ra'], swire['swire_dec'])
)
# Given the graph above, we use 1 arc-second radius
master_catalogue = merge_catalogues(master_catalogue, swire, "swire_ra", "swire_dec", radius=1.*u.arcsec)
nb_merge_dist_plot(
SkyCoord(master_catalogue['ra'], master_catalogue['dec']),
SkyCoord(uhs['uhs_ra'], uhs['uhs_dec'])
)
# Given the graph above, we use 1 arc-second radius
master_catalogue = merge_catalogues(master_catalogue, uhs, "uhs_ra", "uhs_dec", radius=1.*u.arcsec)
When we merge the catalogues, astropy masks the non-existent values (e.g. when a row comes only from a catalogue and has no counterparts in the other, the columns from the latest are masked for that row). We indicate to use NaN for masked values for floats columns, False for flag columns and -1 for ID columns.
for col in master_catalogue.colnames:
if "m_" in col or "merr_" in col or "f_" in col or "ferr_" in col or "stellarity" in col:
master_catalogue[col].fill_value = np.nan
elif "flag" in col:
master_catalogue[col].fill_value = 0
elif "id" in col:
master_catalogue[col].fill_value = -1
master_catalogue = master_catalogue.filled()
master_catalogue[:10].show_in_notebook()
Each pristine catalogue contains a flag indicating if the source was associated to a another nearby source that was removed during the cleaning process. We merge these flags in a single one.
flag_cleaned_columns = [column for column in master_catalogue.colnames
if 'flag_cleaned' in column]
flag_column = np.zeros(len(master_catalogue), dtype=bool)
for column in flag_cleaned_columns:
flag_column |= master_catalogue[column]
master_catalogue.add_column(Column(data=flag_column, name="flag_cleaned"))
master_catalogue.remove_columns(flag_cleaned_columns)
Each pristine catalogue contains a flag indicating the probability of a source being a Gaia object (0: not a Gaia object, 1: possibly, 2: probably, 3: definitely). We merge these flags taking the highest value.
flag_gaia_columns = [column for column in master_catalogue.colnames
if 'flag_gaia' in column]
master_catalogue.add_column(Column(
data=np.max([master_catalogue[column] for column in flag_gaia_columns], axis=0),
name="flag_gaia"
))
master_catalogue.remove_columns(flag_gaia_columns)
Each prisitine catalogue may contain one or several stellarity columns indicating the probability (0 to 1) of each source being a star. We merge these columns taking the highest value.
stellarity_columns = [column for column in master_catalogue.colnames
if 'stellarity' in column]
master_catalogue.add_column(Column(
data=np.nanmax([master_catalogue[column] for column in stellarity_columns], axis=0),
name="stellarity"
))
master_catalogue.remove_columns(stellarity_columns)
master_catalogue.add_column(
ebv(master_catalogue['ra'], master_catalogue['dec'])
)
master_catalogue.add_column(Column(gen_help_id(master_catalogue['ra'], master_catalogue['dec']),
name="help_id"))
master_catalogue.add_column(Column(np.full(len(master_catalogue), "Lockman SWIRE", dtype='<U18'),
name="field"))
# Check that the HELP Ids are unique
if len(master_catalogue) != len(np.unique(master_catalogue['help_id'])):
print("The HELP IDs are not unique!!!")
else:
print("OK!")
specz = Table.read("../../dmu23/dmu23_Lockman-SWIRE/data/Lockman_SWIRE-specz-v2.1.fits")
nb_merge_dist_plot(
SkyCoord(master_catalogue['ra'], master_catalogue['dec']),
SkyCoord(specz['ra'] * u.deg, specz['dec'] * u.deg)
)
master_catalogue = specz_merge(master_catalogue, specz, radius=1. * u.arcsec)
Both SERVS and SWIRE provide IRAC1 and IRAC2 fluxes. SERVS is deeper but tends to under-estimate flux of bright sources (Mattia said over 2000 µJy) as illustrated by this comparison of SWIRE, SERVS, and Spitzer-EIP fluxes.
seip = Table.read("../../dmu0/dmu0_SEIP/data/SEIP_Lockman-SWIRE.fits")
seip_coords = SkyCoord(seip['ra'], seip['dec'])
idx, d2d, _ = seip_coords.match_to_catalog_sky(SkyCoord(master_catalogue['ra'], master_catalogue['dec']))
mask = d2d <= 2 * u.arcsec
fig, ax = plt.subplots()
ax.scatter(seip['i1_f_ap1'][mask], master_catalogue[idx[mask]]['f_ap_servs_irac1'], label="SERVS", s=2.)
ax.scatter(seip['i1_f_ap1'][mask], master_catalogue[idx[mask]]['f_ap_swire_irac1'], label="SWIRE", s=2.)
ax.set_xscale('log')
ax.set_yscale('log')
ax.set_xlabel("SEIP flux [μJy]")
ax.set_ylabel("SERVS/SWIRE flux [μJy]")
ax.set_title("IRAC 1")
ax.legend()
ax.axvline(2000, color="black", linestyle="--", linewidth=1.)
ax.plot(seip['i1_f_ap1'][mask], seip['i1_f_ap1'][mask], linewidth=.1, color="black", alpha=.5);
fig, ax = plt.subplots()
ax.scatter(seip['i2_f_ap1'][mask], master_catalogue[idx[mask]]['f_ap_servs_irac2'], label="SERVS", s=2.)
ax.scatter(seip['i2_f_ap1'][mask], master_catalogue[idx[mask]]['f_ap_swire_irac2'], label="SWIRE", s=2.)
ax.set_xscale('log')
ax.set_yscale('log')
ax.set_xlabel("SEIP flux [μJy]")
ax.set_ylabel("SERVS/SWIRE flux [μJy]")
ax.set_title("IRAC 2")
ax.legend()
ax.axvline(2000, color="black", linestyle="--", linewidth=1.)
ax.plot(seip['i1_f_ap2'][mask], seip['i1_f_ap2'][mask], linewidth=.1, color="black", alpha=.5);
When both SWIRE and SERVS fluxes are provided, we use the SERVS flux below 2000 μJy and the SWIRE flux over.
We create a table indicating for each source the origin on the IRAC1 and IRAC2 fluxes that will be saved separately.
irac_origin = Table()
irac_origin.add_column(master_catalogue['help_id'])
# IRAC1 aperture flux and magnitudes
has_servs = ~np.isnan(master_catalogue['f_ap_servs_irac1'])
has_swire = ~np.isnan(master_catalogue['f_ap_swire_irac1'])
has_both = has_servs & has_swire
print("{} sources with SERVS flux".format(np.sum(has_servs)))
print("{} sources with SWIRE flux".format(np.sum(has_swire)))
print("{} sources with SERVS and SWIRE flux".format(np.sum(has_both)))
has_servs_above_limit = has_servs.copy()
has_servs_above_limit[has_servs] = master_catalogue['f_ap_servs_irac1'][has_servs] > 2000
use_swire = (has_swire & ~has_servs) | (has_both & has_servs_above_limit)
use_servs = (has_servs & ~(has_both & has_servs_above_limit))
print("{} sources for which we use SERVS".format(np.sum(use_servs)))
print("{} sources for which we use SWIRE".format(np.sum(use_swire)))
f_ap_irac = np.full(len(master_catalogue), np.nan)
f_ap_irac[use_servs] = master_catalogue['f_ap_servs_irac1'][use_servs]
f_ap_irac[use_swire] = master_catalogue['f_ap_swire_irac1'][use_swire]
ferr_ap_irac = np.full(len(master_catalogue), np.nan)
ferr_ap_irac[use_servs] = master_catalogue['ferr_ap_servs_irac1'][use_servs]
ferr_ap_irac[use_swire] = master_catalogue['ferr_ap_swire_irac1'][use_swire]
m_ap_irac = np.full(len(master_catalogue), np.nan)
m_ap_irac[use_servs] = master_catalogue['m_ap_servs_irac1'][use_servs]
m_ap_irac[use_swire] = master_catalogue['m_ap_swire_irac1'][use_swire]
merr_ap_irac = np.full(len(master_catalogue), np.nan)
merr_ap_irac[use_servs] = master_catalogue['merr_ap_servs_irac1'][use_servs]
merr_ap_irac[use_swire] = master_catalogue['merr_ap_swire_irac1'][use_swire]
master_catalogue.add_column(Column(data=f_ap_irac, name="f_ap_irac_i1"))
master_catalogue.add_column(Column(data=ferr_ap_irac, name="ferr_ap_irac_i1"))
master_catalogue.add_column(Column(data=m_ap_irac, name="m_ap_irac_i1"))
master_catalogue.add_column(Column(data=merr_ap_irac, name="merr_ap_irac_i1"))
master_catalogue.remove_columns(['f_ap_servs_irac1', 'f_ap_swire_irac1', 'ferr_ap_servs_irac1',
'ferr_ap_swire_irac1', 'm_ap_servs_irac1', 'm_ap_swire_irac1',
'merr_ap_servs_irac1', 'merr_ap_swire_irac1'])
origin = np.full(len(master_catalogue), ' ', dtype='<U5')
origin[use_servs] = "SERVS"
origin[use_swire] = "SWIRE"
irac_origin.add_column(Column(data=origin, name="IRAC1_app"))
# IRAC1 total flux and magnitudes
has_servs = ~np.isnan(master_catalogue['f_servs_irac1'])
has_swire = ~np.isnan(master_catalogue['f_swire_irac1'])
has_both = has_servs & has_swire
print("{} sources with SERVS flux".format(np.sum(has_servs)))
print("{} sources with SWIRE flux".format(np.sum(has_swire)))
print("{} sources with SERVS and SWIRE flux".format(np.sum(has_both)))
has_servs_above_limit = has_servs.copy()
has_servs_above_limit[has_servs] = master_catalogue['f_servs_irac1'][has_servs] > 2000
use_swire = (has_swire & ~has_servs) | (has_both & has_servs_above_limit)
use_servs = (has_servs & ~(has_both & has_servs_above_limit))
print("{} sources for which we use SERVS".format(np.sum(use_servs)))
print("{} sources for which we use SWIRE".format(np.sum(use_swire)))
f_ap_irac = np.full(len(master_catalogue), np.nan)
f_ap_irac[use_servs] = master_catalogue['f_servs_irac1'][use_servs]
f_ap_irac[use_swire] = master_catalogue['f_swire_irac1'][use_swire]
ferr_ap_irac = np.full(len(master_catalogue), np.nan)
ferr_ap_irac[use_servs] = master_catalogue['ferr_servs_irac1'][use_servs]
ferr_ap_irac[use_swire] = master_catalogue['ferr_swire_irac1'][use_swire]
flag_irac = np.full(len(master_catalogue), False, dtype=bool)
flag_irac[use_servs] = master_catalogue['flag_servs_irac1'][use_servs]
flag_irac[use_swire] = master_catalogue['flag_swire_irac1'][use_swire]
m_ap_irac = np.full(len(master_catalogue), np.nan)
m_ap_irac[use_servs] = master_catalogue['m_servs_irac1'][use_servs]
m_ap_irac[use_swire] = master_catalogue['m_swire_irac1'][use_swire]
merr_ap_irac = np.full(len(master_catalogue), np.nan)
merr_ap_irac[use_servs] = master_catalogue['merr_servs_irac1'][use_servs]
merr_ap_irac[use_swire] = master_catalogue['merr_swire_irac1'][use_swire]
master_catalogue.add_column(Column(data=f_ap_irac, name="f_irac_i1"))
master_catalogue.add_column(Column(data=ferr_ap_irac, name="ferr_irac_i1"))
master_catalogue.add_column(Column(data=m_ap_irac, name="m_irac_i1"))
master_catalogue.add_column(Column(data=merr_ap_irac, name="merr_irac_i1"))
master_catalogue.add_column(Column(data=flag_irac, name="flag_irac_i1"))
master_catalogue.remove_columns(['f_servs_irac1', 'f_swire_irac1', 'ferr_servs_irac1',
'ferr_swire_irac1', 'm_servs_irac1', 'flag_servs_irac1', 'm_swire_irac1',
'merr_servs_irac1', 'merr_swire_irac1', 'flag_swire_irac1'])
origin = np.full(len(master_catalogue), ' ', dtype='<U5')
origin[use_servs] = "SERVS"
origin[use_swire] = "SWIRE"
irac_origin.add_column(Column(data=origin, name="IRAC1_total"))
# IRAC2 aperture flux and magnitudes
has_servs = ~np.isnan(master_catalogue['f_ap_servs_irac2'])
has_swire = ~np.isnan(master_catalogue['f_ap_swire_irac2'])
has_both = has_servs & has_swire
print("{} sources with SERVS flux".format(np.sum(has_servs)))
print("{} sources with SWIRE flux".format(np.sum(has_swire)))
print("{} sources with SERVS and SWIRE flux".format(np.sum(has_both)))
has_servs_above_limit = has_servs.copy()
has_servs_above_limit[has_servs] = master_catalogue['f_ap_servs_irac2'][has_servs] > 2000
use_swire = (has_swire & ~has_servs) | (has_both & has_servs_above_limit)
use_servs = (has_servs & ~(has_both & has_servs_above_limit))
print("{} sources for which we use SERVS".format(np.sum(use_servs)))
print("{} sources for which we use SWIRE".format(np.sum(use_swire)))
f_ap_irac = np.full(len(master_catalogue), np.nan)
f_ap_irac[use_servs] = master_catalogue['f_ap_servs_irac2'][use_servs]
f_ap_irac[use_swire] = master_catalogue['f_ap_swire_irac2'][use_swire]
ferr_ap_irac = np.full(len(master_catalogue), np.nan)
ferr_ap_irac[use_servs] = master_catalogue['ferr_ap_servs_irac2'][use_servs]
ferr_ap_irac[use_swire] = master_catalogue['ferr_ap_swire_irac2'][use_swire]
m_ap_irac = np.full(len(master_catalogue), np.nan)
m_ap_irac[use_servs] = master_catalogue['m_ap_servs_irac2'][use_servs]
m_ap_irac[use_swire] = master_catalogue['m_ap_swire_irac2'][use_swire]
merr_ap_irac = np.full(len(master_catalogue), np.nan)
merr_ap_irac[use_servs] = master_catalogue['merr_ap_servs_irac2'][use_servs]
merr_ap_irac[use_swire] = master_catalogue['merr_ap_swire_irac2'][use_swire]
master_catalogue.add_column(Column(data=f_ap_irac, name="f_ap_irac_i2"))
master_catalogue.add_column(Column(data=ferr_ap_irac, name="ferr_ap_irac_i2"))
master_catalogue.add_column(Column(data=m_ap_irac, name="m_ap_irac_i2"))
master_catalogue.add_column(Column(data=merr_ap_irac, name="merr_ap_irac_i2"))
master_catalogue.remove_columns(['f_ap_servs_irac2', 'f_ap_swire_irac2', 'ferr_ap_servs_irac2',
'ferr_ap_swire_irac2', 'm_ap_servs_irac2', 'm_ap_swire_irac2',
'merr_ap_servs_irac2', 'merr_ap_swire_irac2'])
origin = np.full(len(master_catalogue), ' ', dtype='<U5')
origin[use_servs] = "SERVS"
origin[use_swire] = "SWIRE"
irac_origin.add_column(Column(data=origin, name="IRAC2_app"))
# IRAC2 total flux and magnitudes
has_servs = ~np.isnan(master_catalogue['f_servs_irac2'])
has_swire = ~np.isnan(master_catalogue['f_swire_irac2'])
has_both = has_servs & has_swire
print("{} sources with SERVS flux".format(np.sum(has_servs)))
print("{} sources with SWIRE flux".format(np.sum(has_swire)))
print("{} sources with SERVS and SWIRE flux".format(np.sum(has_both)))
has_servs_above_limit = has_servs.copy()
has_servs_above_limit[has_servs] = master_catalogue['f_servs_irac2'][has_servs] > 2000
use_swire = (has_swire & ~has_servs) | (has_both & has_servs_above_limit)
use_servs = (has_servs & ~(has_both & has_servs_above_limit))
print("{} sources for which we use SERVS".format(np.sum(use_servs)))
print("{} sources for which we use SWIRE".format(np.sum(use_swire)))
f_irac = np.full(len(master_catalogue), np.nan)
f_irac[use_servs] = master_catalogue['f_servs_irac2'][use_servs]
f_irac[use_swire] = master_catalogue['f_swire_irac2'][use_swire]
ferr_irac = np.full(len(master_catalogue), np.nan)
ferr_irac[use_servs] = master_catalogue['ferr_servs_irac2'][use_servs]
ferr_irac[use_swire] = master_catalogue['ferr_swire_irac2'][use_swire]
flag_irac = np.full(len(master_catalogue), False, dtype=bool)
flag_irac[use_servs] = master_catalogue['flag_servs_irac2'][use_servs]
flag_irac[use_swire] = master_catalogue['flag_swire_irac2'][use_swire]
m_irac = np.full(len(master_catalogue), np.nan)
m_irac[use_servs] = master_catalogue['m_servs_irac2'][use_servs]
m_irac[use_swire] = master_catalogue['m_swire_irac2'][use_swire]
merr_irac = np.full(len(master_catalogue), np.nan)
merr_irac[use_servs] = master_catalogue['merr_servs_irac2'][use_servs]
merr_irac[use_swire] = master_catalogue['merr_swire_irac2'][use_swire]
master_catalogue.add_column(Column(data=f_irac, name="f_irac_i2"))
master_catalogue.add_column(Column(data=ferr_irac, name="ferr_irac_i2"))
master_catalogue.add_column(Column(data=m_irac, name="m_irac_i2"))
master_catalogue.add_column(Column(data=merr_irac, name="merr_irac_i2"))
master_catalogue.add_column(Column(data=flag_irac, name="flag_irac_i2"))
master_catalogue.remove_columns(['f_servs_irac2', 'f_swire_irac2', 'ferr_servs_irac2',
'ferr_swire_irac2', 'm_servs_irac2', 'flag_servs_irac2', 'm_swire_irac2',
'merr_servs_irac2', 'merr_swire_irac2', 'flag_swire_irac2'])
origin = np.full(len(master_catalogue), ' ', dtype='<U5')
origin[use_servs] = "SERVS"
origin[use_swire] = "SWIRE"
irac_origin.add_column(Column(data=origin, name="IRAC2_total"))
irac_origin.write("{}/lockman-swire_irac_fluxes_origins{}.fits".format(OUT_DIR, SUFFIX), overwrite=True)
These both come from wfcam. If UHS is to be included (it currently is not) we will have to decide between the two. DXS is deeper so it would probably be a case of preferentially taking DXS J fluxes.
wfcam_origin = Table()
wfcam_origin.add_column(master_catalogue['help_id'])
#Total J flux
use_dxs = ~np.isnan(master_catalogue['f_ukidss_j'])
use_uhs = np.isnan(master_catalogue['f_ukidss_j']) & ~np.isnan(master_catalogue['f_uhs_j'])
f_wfcam = np.full(len(master_catalogue), np.nan)
f_wfcam[use_uhs] = master_catalogue['f_uhs_j'][use_uhs]
f_wfcam[use_dxs] = master_catalogue['f_ukidss_j'][use_dxs]
ferr_wfcam = np.full(len(master_catalogue), np.nan)
ferr_wfcam[use_uhs] = master_catalogue['ferr_uhs_j'][use_uhs]
ferr_wfcam[use_dxs] = master_catalogue['ferr_ukidss_j'][use_dxs]
m_wfcam = np.full(len(master_catalogue), np.nan)
m_wfcam[use_uhs] = master_catalogue['m_uhs_j'][use_uhs]
m_wfcam[use_dxs] = master_catalogue['m_ukidss_j'][use_dxs]
merr_wfcam = np.full(len(master_catalogue), np.nan)
merr_wfcam[use_uhs] = master_catalogue['merr_uhs_j'][use_uhs]
merr_wfcam[use_dxs] = master_catalogue['merr_ukidss_j'][use_dxs]
flag_wfcam = np.full(len(master_catalogue), np.nan)
flag_wfcam[use_uhs] = master_catalogue['flag_uhs_j'][use_uhs]
flag_wfcam[use_dxs] = master_catalogue['flag_ukidss_j'][use_dxs]
master_catalogue.add_column(Column(data=f_wfcam, name="f_wfcam_j"))
master_catalogue.add_column(Column(data=ferr_wfcam, name="ferr_wfcam_j"))
master_catalogue.add_column(Column(data=m_wfcam, name="m_wfcam_j"))
master_catalogue.add_column(Column(data=merr_wfcam, name="merr_wfcam_j"))
master_catalogue.add_column(Column(data=flag_wfcam, name="flag_wfcam_j"))
master_catalogue.remove_columns(['f_uhs_j', 'f_ukidss_j',
'ferr_uhs_j', 'ferr_ukidss_j',
'm_uhs_j', 'm_ukidss_j',
'merr_uhs_j', 'merr_ukidss_j',
'flag_uhs_j', 'flag_ukidss_j'])
origin = np.full(len(master_catalogue), ' ', dtype='<U5')
origin[use_uhs] = "UHS"
origin[use_dxs] = "DXS"
wfcam_origin.add_column(Column(data=origin, name="wfcam_j"))
#Aperture J flux
use_ap_dxs = ~np.isnan(master_catalogue['f_ap_ukidss_j'])
use_ap_uhs = np.isnan(master_catalogue['f_ap_ukidss_j']) & ~np.isnan(master_catalogue['f_ap_uhs_j'])
f_ap_wfcam = np.full(len(master_catalogue), np.nan)
f_ap_wfcam[use_ap_uhs] = master_catalogue['f_ap_uhs_j'][use_ap_uhs]
f_ap_wfcam[use_ap_dxs] = master_catalogue['f_ap_ukidss_j'][use_ap_dxs]
ferr_ap_wfcam = np.full(len(master_catalogue), np.nan)
ferr_ap_wfcam[use_ap_uhs] = master_catalogue['ferr_ap_uhs_j'][use_ap_uhs]
ferr_ap_wfcam[use_ap_dxs] = master_catalogue['ferr_ap_ukidss_j'][use_ap_dxs]
m_ap_wfcam = np.full(len(master_catalogue), np.nan)
m_ap_wfcam[use_ap_uhs] = master_catalogue['m_ap_uhs_j'][use_ap_uhs]
m_ap_wfcam[use_ap_dxs] = master_catalogue['m_ap_ukidss_j'][use_ap_dxs]
merr_ap_wfcam = np.full(len(master_catalogue), np.nan)
merr_ap_wfcam[use_ap_uhs] = master_catalogue['merr_ap_uhs_j'][use_ap_uhs]
merr_ap_wfcam[use_ap_dxs] = master_catalogue['merr_ap_ukidss_j'][use_ap_dxs]
master_catalogue.add_column(Column(data=f_ap_wfcam, name="f_ap_wfcam_j"))
master_catalogue.add_column(Column(data=ferr_ap_wfcam, name="ferr_ap_wfcam_j"))
master_catalogue.add_column(Column(data=m_ap_wfcam, name="m_ap_wfcam_j"))
master_catalogue.add_column(Column(data=merr_ap_wfcam, name="merr_ap_wfcam_j"))
master_catalogue.remove_columns(['f_ap_uhs_j', 'f_ap_ukidss_j',
'ferr_ap_uhs_j', 'ferr_ap_ukidss_j',
'm_ap_uhs_j', 'm_ap_ukidss_j',
'merr_ap_uhs_j', 'merr_ap_ukidss_j'])
origin = np.full(len(master_catalogue), ' ', dtype='<U5')
origin[use_ap_uhs] = "UHS"
origin[use_ap_dxs] = "DXS"
wfcam_origin.add_column(Column(data=origin, name="wfcam_ap_j"))
wfcam_origin.write("{}/lockman-swire_wfcam_fluxes_origins{}.fits".format(OUT_DIR, SUFFIX), overwrite=True)
master_catalogue['m_ukidss_k'].name = 'm_wfcam_k'
master_catalogue['merr_ukidss_k'].name = 'merr_wfcam_k'
master_catalogue['f_ukidss_k'].name = 'f_wfcam_k'
master_catalogue['ferr_ukidss_k'].name = 'ferr_wfcam_k'
master_catalogue['flag_ukidss_k'].name = 'flag_wfcam_k'
master_catalogue['m_ap_ukidss_k'].name = 'm_ap_wfcam_k'
master_catalogue['merr_ap_ukidss_k'].name = 'merr_ap_wfcam_k'
master_catalogue['f_ap_ukidss_k'].name = 'f_ap_wfcam_k'
master_catalogue['ferr_ap_ukidss_k'].name = 'ferr_ap_wfcam_k'
SpARCS apears to be significantly deeper and contains both total and aperture magnitudes so we take SpARCS over RCSLenS if both are available
Survey | Bands |
---|---|
SpARCS | ugrz |
RCSLenS | grizy |
megacam_origin = Table()
megacam_origin.add_column(master_catalogue['help_id'])
megacam_stats = Table()
megacam_stats.add_column(Column(data=['u','g','r','i','z','y'], name="Band"))
for col in [ "SpARCS", "RCSLenS"]:
megacam_stats.add_column(Column(data=np.full(6, 0), name="{}".format(col)))
megacam_stats.add_column(Column(data=np.full(6, 0), name="use {}".format(col)))
megacam_bands = ['u','g','r','i','z','y'] # Lowercase naming convention (k is Ks)
for band in ['g', 'r', 'z']: #the bands in both catalogues
# Megacam total flux
has_sparcs = ~np.isnan(master_catalogue['f_cfht_megacam_' + band])
has_rcs = ~np.isnan(master_catalogue['f_rcs_' + band])
use_sparcs = has_sparcs
use_rcs = has_rcs & ~has_sparcs
f_megacam = np.full(len(master_catalogue), np.nan)
f_megacam[use_sparcs] = master_catalogue['f_cfht_megacam_' + band][use_sparcs]
f_megacam[use_rcs] = master_catalogue['f_rcs_' + band][use_rcs]
ferr_megacam = np.full(len(master_catalogue), np.nan)
ferr_megacam[use_sparcs] = master_catalogue['ferr_cfht_megacam_' + band][use_sparcs]
ferr_megacam[use_rcs] = master_catalogue['ferr_rcs_' + band][use_rcs]
m_megacam = np.full(len(master_catalogue), np.nan)
m_megacam[use_sparcs] = master_catalogue['m_cfht_megacam_' + band][use_sparcs]
m_megacam[use_rcs] = master_catalogue['m_rcs_' + band][use_rcs]
merr_megacam = np.full(len(master_catalogue), np.nan)
merr_megacam[use_sparcs] = master_catalogue['merr_cfht_megacam_' + band][use_sparcs]
merr_megacam[use_rcs] = master_catalogue['merr_rcs_' + band][use_rcs]
flag_megacam = np.full(len(master_catalogue), np.nan)
flag_megacam[use_sparcs] = master_catalogue['flag_cfht_megacam_' + band][use_sparcs]
flag_megacam[use_rcs] = master_catalogue['flag_rcs_' + band][use_rcs]
master_catalogue.add_column(Column(data=f_megacam, name="f_megacam_" + band))
master_catalogue.add_column(Column(data=ferr_megacam, name="ferr_megacam_" + band))
master_catalogue.add_column(Column(data=m_megacam, name="m_megacam_" + band))
master_catalogue.add_column(Column(data=merr_megacam, name="merr_megacam_" + band))
master_catalogue.add_column(Column(data=flag_megacam, name="flag_megacam_" + band))
old_columns = ['f_cfht_megacam_' + band,
'ferr_cfht_megacam_' + band,
'm_cfht_megacam_' + band,
'merr_cfht_megacam_' + band,
'flag_cfht_megacam_' + band,
'f_rcs_' + band,
'ferr_rcs_' + band,
'm_rcs_' + band,
'merr_rcs_' + band,
'flag_rcs_' + band]
master_catalogue.remove_columns(old_columns)
origin = np.full(len(master_catalogue), ' ', dtype='<U5')
origin[use_sparcs] = "SpARCS"
origin[use_rcs] = "RCSLenS"
megacam_origin.add_column(Column(data=origin, name= 'f_megacam_' + band ))
megacam_stats['RCSLenS'][megacam_stats['Band'] == band] = np.sum(has_rcs)
megacam_stats['SpARCS'][megacam_stats['Band'] == band] = np.sum(has_sparcs)
megacam_stats['use RCSLenS'][megacam_stats['Band'] == band] = np.sum(use_rcs)
megacam_stats['use SpARCS'][megacam_stats['Band'] == band] = np.sum(use_sparcs)
#No aperture fluxes in rcs so just rename all the sparcs columns
for col in master_catalogue.colnames:
if '_cfht_megacam_' in col:
master_catalogue[col].name = col.replace('_cfht_megacam_','_megacam_')
if ('_rcs_' in col) and ('_ap_rcs_' not in col):
master_catalogue[col].name = col.replace('_rcs_','_megacam_')
if '_ap_rcs_' in col:
master_catalogue.remove_columns(col)
megacam_stats.show_in_notebook()
megacam_origin.write("{}/lockman-swire_megacam_fluxes_origins{}.fits".format(OUT_DIR, SUFFIX), overwrite=True)
master_catalogue.colnames
We add a binary flag_optnir_obs
indicating that a source was observed in a given wavelenght domain:
It's an integer binary flag, so a source observed both in optical and near-infrared by not in mid-infrared would have this flag at 1 + 2 = 3.
Note 1: The observation flag is based on the creation of multi-order coverage maps from the catalogues, this may not be accurate, especially on the edges of the coverage.
Note 2: Being on the observation coverage does not mean having fluxes in that wavelength domain. For sources observed in one domain but having no flux in it, one must take into consideration de different depths in the catalogue we are using.
wfc_moc = MOC(filename="../../dmu0/dmu0_INTWFC/data/lh_intwfc_v2.1_HELP_coverage_MOC.fits")
sparcs_moc = MOC(filename="../../dmu0/dmu0_SpARCS/data/SpARCS_HELP_Lockman-SWIRE_MOC.fits")
rcs_moc = MOC(filename="../../dmu0/dmu0_RCSLenS/data/RCSLenS_Lockman-SWIRE_MOC.fits")
ps1_moc = MOC(filename="../../dmu0/dmu0_PanSTARRS1-3SS/data/PanSTARRS1-3SS_Lockman-SWIRE_MOC.fits")
dxs_moc = MOC(filename="../../dmu0/dmu0_UKIDSS-DXS_DR10plus/data/UKIDSS-DR10plus_Lockman-SWIRE_MOC.fits")
servs_moc = MOC(filename="../../dmu0/dmu0_DataFusion-Spitzer/data/DF-SERVS_Lockman-SWIRE_MOC.fits")
swire_moc = MOC(filename="../../dmu0/dmu0_DataFusion-Spitzer/data/DF-SWIRE_Lockman-SWIRE_MOC.fits")
uhs_moc = MOC(filename="../../dmu0/dmu0_UHS/data/UHS-DR1_Lockman-SWIRE_MOC.fits")
was_observed_optical = inMoc(
master_catalogue['ra'], master_catalogue['dec'],
wfc_moc + sparcs_moc + rcs_moc + ps1_moc)
was_observed_nir = inMoc(
master_catalogue['ra'], master_catalogue['dec'],
dxs_moc + uhs_moc
)
was_observed_mir = inMoc(
master_catalogue['ra'], master_catalogue['dec'],
servs_moc + swire_moc
)
master_catalogue.add_column(
Column(
1 * was_observed_optical + 2 * was_observed_nir + 4 * was_observed_mir,
name="flag_optnir_obs")
)
We add a binary flag_optnir_det
indicating that a source was detected in a given wavelenght domain:
It's an integer binary flag, so a source detected both in optical and near-infrared by not in mid-infrared would have this flag at 1 + 2 = 3.
Note 1: We use the total flux columns to know if the source has flux, in some catalogues, we may have aperture flux and no total flux.
To get rid of artefacts (chip edges, star flares, etc.) we consider that a source is detected in one wavelength domain when it has a flux value in at least two bands. That means that good sources will be excluded from this flag when they are on the coverage of only one band.
# SpARCS is a catalogue of sources detected in r (with fluxes measured at
# this prior position in the other bands). Thus, we are only using the r
# CFHT band.
# Check to use catalogue flags from HSC and PanSTARRS.
nb_optical_flux = (
1 * ~np.isnan(master_catalogue['f_wfc_u']) +
1 * ~np.isnan(master_catalogue['f_wfc_g']) +
1 * ~np.isnan(master_catalogue['f_wfc_r']) +
1 * ~np.isnan(master_catalogue['f_wfc_i']) +
1 * ~np.isnan(master_catalogue['f_wfc_z']) +
1 * ~np.isnan(master_catalogue['f_gpc1_g']) +
1 * ~np.isnan(master_catalogue['f_gpc1_r']) +
1 * ~np.isnan(master_catalogue['f_gpc1_i']) +
1 * ~np.isnan(master_catalogue['f_gpc1_z']) +
1 * ~np.isnan(master_catalogue['f_gpc1_y']) +
1 * ~np.isnan(master_catalogue['f_megacam_u']) +
1 * ~np.isnan(master_catalogue['f_megacam_g']) +
1 * ~np.isnan(master_catalogue['f_megacam_r']) +
1 * ~np.isnan(master_catalogue['f_megacam_i']) +
1 * ~np.isnan(master_catalogue['f_megacam_z']) +
1 * ~np.isnan(master_catalogue['f_megacam_y'])
)
nb_nir_flux = (
1 * ~np.isnan(master_catalogue['f_wfcam_j']) +
1 * ~np.isnan(master_catalogue['f_wfcam_k'])
)
nb_mir_flux = (
1 * ~np.isnan(master_catalogue['f_irac_i1']) +
1 * ~np.isnan(master_catalogue['f_irac_i2']) +
1 * ~np.isnan(master_catalogue['f_irac_i3']) +
1 * ~np.isnan(master_catalogue['f_irac_i4'])
)
has_optical_flux = nb_optical_flux >= 2
has_nir_flux = nb_nir_flux >= 2
has_mir_flux = nb_mir_flux >= 2
master_catalogue.add_column(
Column(
1 * has_optical_flux + 2 * has_nir_flux + 4 * has_mir_flux,
name="flag_optnir_det")
)
We are producing a table associating to each HELP identifier, the identifiers of the sources in the pristine catalogue. This can be used to easily get additional information from them.
id_names = []
for col in master_catalogue.colnames:
if '_id' in col:
id_names += [col]
if '_intid' in col:
id_names += [col]
print(id_names)
master_catalogue[id_names].write(
"{}/master_list_cross_ident_lockman-swire{}.fits".format(OUT_DIR, SUFFIX), overwrite=True)
id_names.remove('help_id')
master_catalogue.remove_columns(id_names)
We are adding a column with a HEALPix index at order 13 associated with each source.
master_catalogue.add_column(Column(
data=coords_to_hpidx(master_catalogue['ra'], master_catalogue['dec'], order=13),
name="hp_idx"
))
columns = ["help_id", "field", "ra", "dec", "hp_idx"]
bands = [column[5:] for column in master_catalogue.colnames if 'f_ap' in column]
for band in bands:
columns += ["f_ap_{}".format(band), "ferr_ap_{}".format(band),
"m_ap_{}".format(band), "merr_ap_{}".format(band),
"f_{}".format(band), "ferr_{}".format(band),
"m_{}".format(band), "merr_{}".format(band)] #No flag_{}
columns += [ 'm_megacam_i', 'ferr_megacam_i', 'merr_megacam_i', 'flag_megacam_i', 'f_megacam_i',
'm_megacam_y', 'ferr_megacam_y', 'merr_megacam_y', 'flag_megacam_y', 'f_megacam_y']
columns += ["stellarity", "flag_cleaned", "flag_merged", "flag_gaia", "flag_optnir_obs", "flag_optnir_det",
"zspec", "zspec_qual", "zspec_association_flag", "ebv"]
# We check for columns in the master catalogue that we will not save to disk.
print("Missing columns: {}".format(set(master_catalogue.colnames) - set(columns)))
master_catalogue[columns].write("{}/master_catalogue_lockman-swire{}.fits".format(OUT_DIR, SUFFIX), overwrite=True)