In [ ]:
from os import path
# Third-party
from astropy.io import fits
import astropy.coordinates as coord
import astropy.units as u
import matplotlib as mpl
import matplotlib.pyplot as plt
import numpy as np
%matplotlib inline
from sqlalchemy import func
from twoface import TWOFACE_CACHE_PATH
from twoface.db import (db_connect, AllStar, AllVisit, AllVisitToAllStar,
StarResult, Status, JokerRun, initialize_db)
from twoface.data import APOGEERVData
from twoface.plot import plot_data_orbits
from twoface.mass import m2_func
from config import FIGURES_PATH1
In [ ]:
allvisit = fits.getdata('/Users/adrian/data/APOGEE_DR14/allVisit-l31c.2.fits')
allstar = fits.getdata('/Users/adrian/data/APOGEE_DR14/allStar-l31c.2.fits')
print("{0} unique stars in APOGEE DR14".format(len(np.unique(allstar['APOGEE_ID']))))
print("{0} visits in APOGEE DR14".format(len(allvisit)))
In [ ]:
(allstar['NVISITS'] == 1).sum(), (allstar['NVISITS'] > 1).sum()
In [ ]:
Session, _ = db_connect(path.join(TWOFACE_CACHE_PATH, 'apogee.sqlite'))
session = Session()
In [ ]:
print("{0} unique stars in TwoFace DB".format(session.query(AllStar.apogee_id).distinct().count()))
print("{0} visits in TwoFace DB".format(session.query(AllVisit.id).distinct().count()))
In [ ]:
nvisits = np.array(session.query(func.count(AllVisit.id)).join(AllVisitToAllStar, AllStar)\
.group_by(AllStar.apogee_id)\
.having(func.count(AllVisit.id) >= 3).all())
nvisits = nvisits[:,0]
In [ ]:
fig, ax = plt.subplots(1, 1, figsize=(6, 4))
_ = ax.hist(nvisits, bins=np.logspace(0.5, 7, 22, base=2.), rasterized=True)
ax.set_xscale('log', basex=2)
ax.set_yscale('log')
ax.set_xlim(2.5, 140)
ax.set_ylim(0.8, 10**5.1)
ax.xaxis.set_ticks(2**np.arange(2, 7+1, 1))
ax.xaxis.set_ticklabels([str(x) for x in ax.get_xticks()])
ax.yaxis.set_ticks(10**np.arange(0, 5+1, 1))
ax.set_xlabel('$N$ visits')
ax.set_ylabel('$N$ stars')
fig.tight_layout()
fig.savefig(str(FIGURES_PATH1 / 'nvisits.pdf'), rasterized=True, dpi=250)
In [ ]:
(nvisits < 8).sum(), (nvisits >= 8).sum()
In [ ]:
(nvisits < 8).sum() / len(nvisits)
logg - Teff and metallicity of sample:
In [ ]:
res = np.array(session.query(AllStar.logg, AllStar.teff, AllStar.fe_h).all())
logg,teff,feh = res.T
In [ ]:
fig, axes = plt.subplots(1, 2, figsize=(10, 4.5))
axes[0].hist(feh, bins=np.linspace(-3, 1, 64), rasterized=True);
axes[0].set_yscale('log')
axes[0].set_xlabel(r'[Fe/H]')
axes[0].set_ylabel(r'$N$')
axes[0].xaxis.set_ticks([-3, -2, -1, 0, 1])
c = axes[1].scatter(teff, logg, c=feh, marker='.', alpha=0.2,
linewidth=0, vmin=-1, vmax=0.75,
cmap='RdYlBu_r', s=8, rasterized=True)
axes[1].set_xlim(5600, 3600)
axes[1].set_ylim(4.5, -0.5)
axes[1].set_xlabel(r'$T_{\rm eff}$ [K]')
axes[1].set_ylabel(r'$\log g$')
axes[1].xaxis.set_ticks(np.arange(4000, 5500+1, 500))
fig.tight_layout()
fig.subplots_adjust(right=0.82)
cax = fig.add_axes([0.85, 0.185, 0.02, 0.78])
cb = fig.colorbar(c, cax=cax, drawedges=False)
cb.set_label('[Fe/H]')
cb.solids.set_edgecolor("face")
cb.solids.set_rasterized(True)
fig.savefig(str(FIGURES_PATH1 / 'logg_teff_feh.pdf'), rasterized=True, dpi=250)
In [ ]:
print('{0} total visits'.format(len(allvisit)))
mask1 = (np.isfinite(allvisit['VHELIO']) &
np.isfinite(allvisit['VRELERR']) &
(allvisit['VRELERR'] < 100.))
print('mask 1: removes {0}'.format(np.logical_not(mask1).sum()))
skip_mask = np.sum(2 ** np.array([9, 12, 13])) # PERSIST_HIGH, PERSIST_JUMP_POS, PERSIST_JUMP_NEG
mask2 = (allvisit['STARFLAG'] & skip_mask) == 0
print('mask 2: removes {0} more'.format((np.logical_not(mask2) & mask1).sum()))
skip_mask += np.sum(2 ** np.array([3, 4])) # VERY_BRIGHT_NEIGHBOR, LOW_SNR
mask3 = (allvisit['STARFLAG'] & skip_mask) == 0
print('mask 3: removes {0} more'.format((np.logical_not(mask3) & mask2 & mask1).sum()))
# Remove STAR_BAD stars and restrict logg range:
apogee_ids = allstar['APOGEE_ID'][((allstar['ASPCAPFLAG'] & np.sum(2 ** np.array([23]))) == 0) &
(allstar['LOGG'] > 0) & (allstar['LOGG'] < 4)]
print('(star cut removes {0} stars)'.format(len(allstar)-len(apogee_ids)))
mask4 = np.isin(allvisit['APOGEE_ID'], apogee_ids)
print('mask 4: removes {0} more'.format((np.logical_not(mask4) & mask3 & mask2 & mask1).sum()))
print('\n total removed: {0}'.format(np.logical_not(mask4 & mask3 & mask2 & mask1).sum()))
In [ ]:
tmp_visits = allvisit[mask1 & mask2 & mask3 & mask4]
v_apogee_ids, counts = np.unique(tmp_visits['APOGEE_ID'], return_counts=True)
stars = allstar[np.isin(allstar['APOGEE_ID'], v_apogee_ids[counts >= 3])]
print(len(np.unique(stars['APOGEE_ID'])))
visits = tmp_visits[np.isin(tmp_visits['APOGEE_ID'], stars['APOGEE_ID'])]
print(len(visits))
In [ ]: