<jupyter_start><jupyter_text>Comparisons with APW et al<jupyter_code>import numpy as np
import matplotlib
import matplotlib.pyplot as plt
import pickle
from astropy.table import Table

import sys
sys.path.insert(0, "../")

import velociraptor

%matplotlib inline
# TODO: cross-match this
#Get data from APOGEE
apw_lnk_percentiles = Table.read("../data/lnK-percentiles-xm-rv-all.fits")
x = apw_lnk_percentiles["lnK_per_1"]
y = apw_lnk_percentiles["rv_excess_variance"]**0.5

fig, ax = plt.subplots(figsize=(12, 12))
ax.scatter(x, y, s=5, alpha=0.5, c="#000000")
ax.set_xlabel(r"$\log{(K\,\,/\,\,\mathrm{km\,s}^{-1})}$ $(\mathrm{APW})$")
ax.set_ylabel(r"$\sigma_{rv,excess}\,\,/\,\,\mathrm{km\,s}^{-1}$")

x = np.exp(apw_lnk_percentiles["lnK_per_1"])
y = apw_lnk_percentiles["rv_excess_variance"]**0.5
z = apw_lnk_percentiles["p_sb_50"]

fig, ax = plt.subplots(figsize=(12, 12))
scat = ax.scatter(x, y, c=z, s=5, alpha=1, cmap="viridis")
cbar = plt.colorbar(scat)
cbar.set_label(r"\textrm{binary probability}")


ax.set_xlabel(r"$K_{1\%}\,\,/\,\,\mathrm{km\,s}^{-1}$ $(\mathrm{APW})$")
ax.set_ylabel(r"$\sigma_{rv,excess}\,\,/\,\,\mathrm{km\,s}^{-1}$")
ax.semilogx()
ax.semilogy()
limits = (1e-3, 1e2)
ax.plot(limits, limits, c="#666666", linestyle=":", linewidth=0.5)
ax.set_xlim(limits)
ax.set_ylim(limits)
x = np.exp(apw_lnk_percentiles["lnK_per_1"])
y = apw_lnk_percentiles["rv_excess_variance"]**0.5
z = apw_lnk_percentiles["phot_g_mean_mag"]

fig, ax = plt.subplots(figsize=(12, 12))
scat = ax.scatter(x, y, c=z, s=5, alpha=1, cmap="viridis", vmax=13, vmin=10)
cbar = plt.colorbar(scat)


ax.set_xlabel(r"$K_{1\%}\,\,/\,\,\mathrm{km\,s}^{-1}$ $(\mathrm{APW})$")
ax.set_ylabel(r"$\sigma_{rv,excess}\,\,/\,\,\mathrm{km\,s}^{-1}$")
ax.semilogx()
ax.semilogy()
limits = (1e-3, 1e2)
ax.plot(limits, limits, c="#666666", linestyle=":", linewidth=0.5)
ax.set_xlim(limits)
ax.set_ylim(limits)
print(apw_lnk_percentiles.dtype.names)
# Comparisons with unimodal catalog.
unimodal = Table.read("../data/highK-unimodal-xm-rv-all.fits")

print(unimodal.dtype.names)
denominator = (unimodal["P"] * np.sqrt(1 - unimodal["e"]**2))
x = unimodal["K"]/denominator
y = np.sqrt(unimodal["rv_excess_variance"])/denominator

fig, ax = plt.subplots(1, 1, figsize=(12, 12))
ax.scatter(x, y)
ax.set_xlabel(r"${K}/{P\sqrt{1-e^2}}\,\,/\,\mathrm{km\,s}^{-1}\,\mathrm{day}^{-1}$")
ax.set_ylabel(r"${\sigma_\mathrm{rv,excess}}/{P\sqrt{1-e^2}}\,\,/\,\mathrm{km\,s}^{-1}\,\mathrm{day}^{-1}$")
ax.semilogx()
ax.semilogy()

_ = np.array([x, y])
limits = (1e-3, 1.5 * np.nanmax(_))
ax.plot(limits, limits, c="#666666", zorder=-1, lw=0.5, linestyle=":")
ax.set_xlim(limits)
ax.set_ylim(limits)

# What *kinds* of sources did we detect or not detect?
consistent_with_single = (unimodal["rv_excess_variance"] == 0)
def what_did_we_find(data, consistent_with_single,
                     latex_label_name=None, ax=None, **kwargs):
    if ax is None:
        fig, ax = plt.subplots()
    else:
        fig = ax.figure
    
    kwds = dict(histtype="bar", stacked="True", alpha=0.5)
    kwds.update(kwargs)

    finite = np.isfinite(data)
    ax.hist([
        data[finite * consistent_with_single],
        data[finite * (~consistent_with_single)],
    ], label=("single", "binary"), **kwargs)
    
    legend = plt.legend()
    
    if latex_label_name is not None:
        ax.set_xlabel(latex_label_name)
    
    ax.set_ylabel(r"\textrm{count}")
    
    fig.tight_layout()
    
    return fig

fig = what_did_we_find(unimodal["K"], consistent_with_single,
                       latex_label_name=r"$K\,\,/\,\,{\rm km\,s}^{-1}$",
                       bins=30)
fig = what_did_we_find(unimodal["P"], consistent_with_single,
                       latex_label_name=r"$P\,\,/\,\,{\rm days}^{-1}$",
                       bins=30)
fig = what_did_we_find(unimodal["e"], consistent_with_single,
                       latex_label_name=r"$e$",
                       bins=30)
label_names = ("P", "M0", "M1", "K", "v0", "omega", "e")
K = len(label_names)
fig, axes = plt.subplots(1, K, figsize=(4*K, 4))
for k, (ax, label_name) in enumerate(zip(axes, label_names)):
    _ =  what_did_we_find(unimodal[label_name], consistent_with_single,
                         latex_label_name=label_name, bins=30, ax=ax)
    
denominator = (unimodal["P"] * np.sqrt(1 - unimodal["e"]**2))
x = unimodal["K"]/denominator
y = np.sqrt(unimodal["rv_excess_variance"])/denominator
z = unimodal["p_sb_50"]

fig, ax = plt.subplots(1, 1, figsize=(12, 12))
scat = ax.scatter(x, y, c=z, s=50, cmap="viridis", vmin=0, vmax=1)
ax.set_xlabel(r"${K}/{P\sqrt{1-e^2}}\,\,/\,\mathrm{km\,s}^{-1}\,\mathrm{day}^{-1}$")
ax.set_ylabel(r"${\sigma_\mathrm{rv,excess}}/{P\sqrt{1-e^2}}\,\,/\,\mathrm{km\,s}^{-1}\,\mathrm{day}^{-1}$")
ax.semilogx()
ax.semilogy()
cbar = plt.colorbar(scat)
cbar.set_label(r"\textrm{binary probability}")

_ = np.array([x, y])
limits = (1e-3, 1.5 * np.nanmax(_))
ax.plot(limits, limits, c="#666666", zorder=-1, lw=0.5, linestyle=":")
ax.set_xlim(limits)
ax.set_ylim(limits)

denominator = (unimodal["P"] * np.sqrt(1 - unimodal["e"]**2))
x = unimodal["K"]/denominator
y = np.sqrt(unimodal["rv_excess_variance"])/denominator
z = np.log10(unimodal["P"])

fig, ax = plt.subplots(1, 1, figsize=(12, 12))
scat = ax.scatter(x, y, c=z, s=50, cmap="viridis")
ax.set_xlabel(r"${K}/{P\sqrt{1-e^2}}\,\,/\,\mathrm{km\,s}^{-1}\,\mathrm{day}^{-1}$")
ax.set_ylabel(r"${\sigma_\mathrm{rv,excess}}/{P\sqrt{1-e^2}}\,\,/\,\mathrm{km\,s}^{-1}\,\mathrm{day}^{-1}$")
ax.semilogx()
ax.semilogy()
cbar = plt.colorbar(scat)
cbar.set_label(r"$\log_{10}(P\,\,/\,\,\mathrm{days}^{-1})$")

_ = np.array([x, y])
limits = (1e-3, 1.5 * np.nanmax(_))
ax.plot(limits, limits, c="#666666", zorder=-1, lw=0.5, linestyle=":")
ax.set_xlim(limits)
ax.set_ylim(limits)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.23355149844987944, 'fraction_numerical': 0.033069238718567, 'mean_word_length': 4.38682745825603, 'pattern_counts': {'":': 5, '<': 5, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '1', 'dir': 'jupyter-structured-clean-dedup', 'id': '140777', 'license': 'MIT', 'n_tokens_mistral': 2598, 'n_tokens_neox': 2466, 'n_words': 356, 'path': 'notebooks/comparison-apw.ipynb', 'repo_name': 'andycasey/velociraptor', 'size': '5806'}