In [1]:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from astropy.io import fits
import seaborn as sns
%matplotlib notebook
The following query was run to select sources between declination -10 and -9.66666666666, in order to match these sources to those that were classified by the PS1 RF model:
select ps1.*,
sqrt(power(pmra,2) + power(pmdec,2)) as pm,
sqrt( power(pmra,2)*power(pmra_error, 2)/(power(pmra,2) + power(pmdec,2))
+ power(pmdec,2)*power(pmdec_error, 2)/(power(pmra,2) + power(pmdec,2))
+ 2*pmra*pmdec/(power(pmra,2) + power(pmdec,2))*pmra_pmdec_corr*pmra_error*pmdec_error) as pm_unc,
gaia.parallax_over_error, gaia.phot_g_mean_mag, gaia.b
from gaiadr2.panstarrs1_best_neighbour as ps1
inner join gaiadr2.gaia_source as gaia
on ps1.source_id = gaia.source_id
where gaia.astrometric_params_solved > 3
and gaia.dec between -10 and -9.6666666666666667
and phot_bp_rp_excess_factor < 1.3+0.06*power(phot_bp_mean_mag-phot_rp_mean_mag,2)
AND phot_bp_rp_excess_factor > 1.0+0.015*power(phot_bp_mean_mag-phot_rp_mean_mag,2)
AND astrometric_chi2_al/(astrometric_n_good_obs_al-5) < 1.44*greatest(1,exp(-0.4*(phot_g_mean_mag-19.5)))
order by ps1.original_ext_source_id
The results are stored in the file neg10_match-result.fits
In [16]:
gaia = fits.getdata("neg10_match-result.fits")
In [17]:
ps1 = pd.read_hdf("/Users/adamamiller/Desktop/PS1_fits/hdf5/dec_neg10_0_classifications.h5")
In [18]:
gaia_objid = np.array(gaia['original_ext_source_id'], dtype="int64") #.byteswap().newbyteorder()
gaia_parallax = np.array(gaia['parallax_over_error'], dtype="float") #.byteswap().newbyteorder()
gaia_pm = np.array(gaia['pm']/gaia['pm_unc'], dtype="float") #.byteswap().newbyteorder()
gaia_g = np.array(gaia["phot_g_mean_mag"], dtype="float")
gaia_b = np.array(gaia["b"], dtype="float")
In [19]:
gaia_df = pd.DataFrame(gaia_objid, columns=["objid"])
gaia_df["parallax_signif"] = gaia_parallax
gaia_df["pm_signif"] = gaia_pm
gaia_df["G"] = gaia_g
gaia_df["b"] = gaia_b
In [20]:
ps1_objid = np.array(ps1["objid"], dtype="int64") #.byteswap().newbyteorder()
ps1_sg_score = np.array(ps1["rf_score"], dtype="float") #.byteswap().newbyteorder()
In [21]:
ps1_df = pd.DataFrame(ps1_objid, columns=["objid"])
ps1_df["sg_score"] = ps1_sg_score
Merge the two DFs in order to directly compare the sg_score to the parallax measurements
In [22]:
merge = pd.merge(ps1_df, gaia_df, on="objid")
In [23]:
merge
Out[23]:
In [52]:
fig, ax = plt.subplots()
hax = ax.hexbin(merge["sg_score"], merge["parallax_signif"],
bins="log", mincnt=1, gridsize=1000,
extent=(0,1,-20,120))
ax.plot([-1,2], [6,6], 'DarkOrange', lw=2)
ax.set_xlim(0,1.02)
ax.set_ylim(-10,40)
ax.set_xlabel("sg_score")
ax.set_ylabel("parallax SNR")
fig.colorbar(hax)
fig.tight_layout()
In [22]:
# hax = sns.jointplot(merge["sg_score"], merge["parallax_signif"],
# kind="hex", stat_func=None,
# bins="log", mincnt=1, gridsize=1000,
# xlim=(0,1), ylim=(-20,120),
# cmap="viridis",
# marginal_kws={'bins':1000})
# hax.ax_joint.set_ylim(-10,40)
In [53]:
fig, ax = plt.subplots()
hax = ax.hexbin(merge["G"], merge["parallax_signif"],
bins="log", mincnt=1, gridsize=1000,
extent=(10,21,-20,120))
ax.plot([-1,22], [6,6], 'DarkOrange', lw=2)
ax.set_xlim(12,21)
ax.set_ylim(-10,40)
ax.set_xlabel("G (mag)")
ax.set_ylabel("parallax SNR")
fig.colorbar(hax)
fig.tight_layout()
In [54]:
fig, ax = plt.subplots()
hax = ax.hexbin(merge["b"], merge["parallax_signif"],
bins="log", mincnt=1, gridsize=1000,
extent=(-75,54,-20,120))
ax.plot([-75,54], [6,6], 'DarkOrange', lw=2)
ax.set_xlim(-75,54)
ax.set_ylim(-10,40)
ax.set_xlabel("Galactic latitude")
ax.set_ylabel("parallax SNR")
fig.colorbar(hax)
fig.tight_layout()
In [44]:
plt.hexbin?
Create the same 3 plots, but this time showing the proper motion SNR as a function of sg_score, G mag, and galactic latitude in the existing PS1 catalog.
In [55]:
fig, ax = plt.subplots()
hax = ax.hexbin(merge["sg_score"], merge["pm_signif"],
bins="log", mincnt=1, gridsize=1000,
extent=(0,1,-1,150))
ax.plot([-1,2], [54,54], 'DarkOrange', lw=2)
ax.set_xlim(0,1)
ax.set_ylim(-0.5,150)
ax.set_xlabel("sg_score")
ax.set_ylabel("pm SNR")
fig.colorbar(hax)
fig.tight_layout()
In [56]:
fig, ax = plt.subplots()
hax = ax.hexbin(merge["G"], merge["pm_signif"],
bins="log", mincnt=1, gridsize=1000,
extent=(10,21,-1,150))
# ax.plot([-1,22], [80,80], 'DarkOrange', lw=2)
ax.plot([-1,22], [54,54], 'DarkOrange', lw=2)
ax.set_xlim(12,21)
ax.set_ylim(-0.5,150)
ax.set_xlabel("G (mag)")
ax.set_ylabel("pm SNR")
fig.colorbar(hax)
fig.tight_layout()