In [1]:
# This changes the current directory to the base saga directory - make sure to run this first!
# This is necessary to be able to import the py files and use the right directories,
# while keeping all the notebooks in their own directory.
import os
import sys
if 'saga_base_dir' not in locals():
saga_base_dir = os.path.abspath('..')
if saga_base_dir not in sys.path:
os.chdir(saga_base_dir)
In [2]:
%matplotlib inline
from matplotlib import pyplot as plt
from matplotlib import rcParams
rcParams['figure.figsize'] = (16, 10)
rcParams['image.interpolation'] = 'none'
rcParams['image.origin'] = 'lower'
In [3]:
for module in ['hosts', 'targeting', 'magellan']:
if module in globals():
reload(globals()[module])
else:
globals()[module] = __import__(module)
#g = targeting.get_gama() #re-caches the gama catalog
In [4]:
hostlst = hosts.get_saga_hosts_from_google('etollerud', )
In [190]:
from astropy import table
from astropy.table import Table, Column
from astropy.coordinates import SkyCoord
In [6]:
aen = [h for h in hostlst if h.name.lower() == 'aeneid'][0]
aen.shortname = 'Ae15B'
In [292]:
aencatall = aen.load_and_reprocess_sdss_catalog('catalogs/base_sql_nsa148734.fits.gz')
aencat = aencatall[aencatall['REMOVE']==-1]
aencat
Out[292]:
In [293]:
#this is the OLD one, without WISE or AAT inputs. Below we replace this with the new one
mlcat = Table.read('catalogs/Nsaid_148734_Prediction.txt.gz', format='ascii.commented_header')
mlcat
Out[293]:
In [383]:
#there are two ML catalogs, one with WISE and the other without
mlcatwowise = Table.read('catalogs/SAGAobjidPrediction_SDSS_148724.fit')
mlcatwwise = Table.read('catalogs/SAGAobjidPrediction_SDSSwise_148724.fit')
#their OBJID columns are strings and are weirdly broken - need to fix:
for cat in (mlcatwowise, mlcatwwise):
newoid = [int(i) for i in cat['OBJID']]
del cat['OBJID']
cat['OBJID'] = newoid
mlcatwowisetokeep = np.array([i not in mlcatwwise['OBJID'] for i in mlcatwowise['OBJID']])
mlcat = table.vstack((mlcatwwise, mlcatwowise[mlcatwowisetokeep]))
mlcat
Out[383]:
In [399]:
#oh wait, now there's a catalog that's pre-combined them in a better way. Use that.
mlcat = Table.read('catalogs/SAGAobjidMaxPrediction_148724.fit')
mlcat
Out[399]:
In [400]:
#now xmatch in the sense of ML things with matches in aencat
catsc = SkyCoord(aencat['ra'], aencat['dec'], unit=u.deg)
mlsc = SkyCoord(mlcat['RA'], mlcat['DEC'], unit=u.deg)
idxsc, d2d, _ = mlsc.match_to_catalog_sky(catsc)
In [401]:
plt.hist(d2d.arcsec, bins=100, histtype='step',range=(0,60))
plt.tight_layout()
np.sum(d2d<1*u.arcsec),np.sum(d2d>1*u.arcsec), mlcat[d2d>1*u.arcsec]['PROBABILITY_CLASS_1']
Out[401]:
In [402]:
#add probabilities to catalog
probs = -np.ones(len(aencat), dtype=float) #-1 for those w/o
probs[idxsc[d2d<1*u.arcsec]] = mlcat['PROBABILITY_CLASS_1'][d2d<1*u.arcsec]
if 'ML_prob' in aencat.colnames:
aencat['ML_prob'] = probs
else:
aencat.add_column(Column(name='ML_prob', data=probs))
In [403]:
gricolorcuts = {'g-r': (None, 0.8, 2),
'r-i': (None, 0.5, 2)}
sagacolorcuts = gricolorcuts.copy()
sagacolorcuts['r-K'] = (None, 2.0, 2)
sagacolorcuts['r-w1'] = (None, 2.6, 2)
def uggrline_cut(cat, sl=1.5, inter=-0.2):
gt = cat['u']-cat['g']+2*(cat['u_err']+cat['g_err'])
lt = sl*(cat['g']-cat['r']-2*(cat['g_err']+cat['r_err'])) + inter
return gt > lt
sagacolorcuts['funcs'] = [uggrline_cut]
In [404]:
aen.environskpc
Out[404]:
In [405]:
magcut = 22
aentargs = targeting.select_targets(aen, colorcuts=gricolorcuts, outercutrad=aen.environskpc*u.kpc,
galvsallcutoff=magcut, faintlimit=magcut, verbose=True,
removespecstars=False, removegalsathighz =False,
catalog=aencat) #do these because the fits catalogs don't have spec_class
len(aentargs)
Out[405]:
In [406]:
goodprob = 0.05
closegoodprob_aencatmsk = (aencat['ML_prob']>goodprob)&(aencat['rhost_kpc']<aen.environskpc)&(aencat['r']<magcut)
np.sum(aentargs['ML_prob']>goodprob), np.sum(closegoodprob_aencatmsk)
Out[406]:
In [407]:
#so need to add in the missing ones...
print(len(aentargs))
torem = aentargs['ML_prob']>goodprob # first remove them
aentargs = aentargs[~torem]
print(len(aentargs))
#then add them back in
aentargs = table.vstack((aentargs, aencat[closegoodprob_aencatmsk]))
print(len(aentargs))
In [408]:
pris = np.zeros(len(aentargs), dtype=int)
pris[aentargs['r']<22] = 4
pris[aentargs['r']<21.5] = 3
pris[aentargs['r']<21] = 2
pris[aentargs['ML_prob']>goodprob] = 1
if 'imacs_pri' in aentargs.colnames:
aentargs['imacs_pri'] = pris
else:
aentargs.add_column(Column(name='imacs_pri', data=pris))
np.bincount(pris)
Out[408]:
In [409]:
allspectaken = Table.read('SAGADropbox/data/allspectaken_v3.fits.gz')
allspectaken = allspectaken[allspectaken['HOST_NSAID'] == aen.nsaid]
allspectaken
Out[409]:
In [410]:
# now ID those with a good enough zquality
okzqual = 3
goodspectaken = allspectaken[allspectaken['ZQUALITY'] >= okzqual]
In [411]:
#x-match with the target list
goodspecsc = SkyCoord(goodspectaken['ra'], goodspectaken['dec'], unit=u.deg)
aentargssc = SkyCoord(aentargs['ra'], aentargs['dec'], unit=u.deg)
idx, d2d, _ = aentargssc.match_to_catalog_sky(goodspecsc)
plt.hist(d2d.arcsec, bins=100,range=(0,10),log=True,histtype='step')
alreadyobserved = d2d<1*u.arcsec
np.sum(alreadyobserved), np.sum(~alreadyobserved)
Out[411]:
In [412]:
!rm imacs_targets/Ae15B.cat imacs_targets/Ae15B_ini.obs
magellan.build_imacs_targeting_files(aen, 'Marchi/Munoz', targs=aentargs[~alreadyobserved], date='2015-9-5',
refmagrange={'r':(17, 17.5)}, inclhost=False, onlygals=False)
In [414]:
plt.figure(figsize=(10,8))
magellan.plot_imacs_masks(aen, eastleft=True, plotpris =True)
plt.tight_layout()