In [1]:
import splat
import wisps
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from os.path import expanduser
homedir = expanduser("~")
from astropy.io import ascii
%matplotlib inline
from scipy import stats
import ast

In [2]:
df=ascii.read('/users/caganze/aegis_3dhst.v4.1.cats/Catalog/aegis_3dhst.v4.1.cat').to_pandas()

In [3]:
df['id']


Out[3]:
0            1
1            2
2            3
3            4
4            5
         ...  
41195    41196
41196    41197
41197    41198
41198    41199
41199    41200
Name: id, Length: 41200, dtype: int64

In [4]:
def strip_diction(s):
        #from string to dictionary
        return pd.Series(ast.literal_eval(((s.strip().replace('nan', "'nan'").replace('-inf', "'nan'")).replace('inf', "'nan'"))))
    
def reformat_index_table(df):
    #assign dictionary keys
    ids_list=[]
    snr_list=[]
    for k in df.columns:
        first=df[k].iloc[0]
        print (k)
        if isinstance(first, str):
            if first.strip().startswith("{"):
                if not k=='indices':
                    snr_list.append(df[k].apply(lambda x: strip_diction(x)))
                else:
                    ids_list.append(df[k].apply(lambda x: strip_diction(x)).applymap(tuple))
      
    
    return snr_list, ids_list

def replace(x):
    #print (x)
    if (isinstance(x, float)):
        if (np.isnan(x)):
            return tuple([np.nan, np.nan])
    else:
        return x
    
def comp_f_test(tupl):
    num=tupl[0]
    den=tupl[1]
    if den ==0: return np.nan
    else:
        x=num/den
        return stats.f.cdf(x, 2, 1, 0, scale=1)

In [5]:
#get photometry catalogs
w_phot=pd.read_csv(wisps.OUTPUT_FILES+'/wisp_photometry.csv')
ht_phot=pd.read_csv(wisps.OUTPUT_FILES+'/hst3d_photometry_all.csv')

In [6]:
hst_phot=ht_phot[~ ht_phot.grism_id.isin([ht_phot.grism_id.iloc[-1]])]

In [7]:
#use the indices
hst_3d_indices=pd.read_pickle('/users/caganze/new_indices_measurements_october.pkl')
wisp_indices=pd.read_pickle('/users/caganze/new_wisps_indices_measurements_october.pkl')

In [8]:
#drop nan columns
hst_3d_indices=(hst_3d_indices[~ hst_3d_indices.isnull().all(1)])

In [9]:
#combine both
new_indices=pd.concat([wisp_indices,hst_3d_indices ])

In [10]:
new_indices=new_indices[~ new_indices.isnull().all(1)]

In [11]:
new_indices=new_indices.reset_index(drop=True)

In [12]:
new_indices.columns


Out[12]:
Index(['CH_4/H-Cont', 'CH_4/H_2O-1', 'CH_4/H_2O-2', 'CH_4/J-Cont',
       'H-cont/H_2O-1', 'H-cont/H_2O-2', 'H-cont/J-Cont', 'H_2O-1/J-Cont',
       'H_2O-2/H_2O-1', 'H_2O-2/J-Cont', 'cdf_snr', 'f_test', 'line_chi',
       'name', 'snr1', 'snr2', 'snr3', 'snr4', 'spex_chi', 'spt'],
      dtype='object')

In [13]:
#cols=list(wisps.INDEX_NAMES)

In [14]:
def combined_wisp_hst_catalogs(hst3d_phot,wisp_phot, indices):
    """
    combine both hst-3d and wisps into one big file with all the information
    """
    #hst_3d does not have 110 photometry
    hst3d_phot['F110_mag']=np.nan
    hst3d_phot['F110_mag_er']=np.nan

    
    #combine flags into one flag
    flgs=hst3d_phot[['use_phot_x', 'f_cover', 'f_flagged', 'f_negative']].values
    hst3d_phot['phot_flags']= pd.Series([i for i in flgs])
    
    hst3d_phot['survey']='HST3D'
    wisp_phot['survey']='WISP'
    wisp_phot=wisp_phot.rename(columns={'EXTRACTION_FLAG':'phot_flags'})
   
    #rename some columns
    indices=indices.rename(columns={'name':'grism_id'})
    
    ##drop .ascii from hst_phot
    indices['grism_id']=indices['grism_id']
    
    #combined_photometry (the order matters: HST3D+WISPP
    comb_phot=pd.DataFrame()
    grism_ids=hst3d_phot['grism_id'].append(wisp_phot['grism_id'])
    star_flags=hst3d_phot['star_flag'].append(wisp_phot['star_flag'])
    comb_phot['grism_id']=grism_ids
    comb_phot['star_flag']=star_flags
    comb_phot['class_star']=hst3d_phot['class_star'].append(wisp_phot['star_flag'])
    
    print(comb_phot.columns)
    for flt in ['110', '140', '160']:
        mag_tuple1=hst3d_phot[['F'+flt+'_mag', 'F'+flt+'_mag_er']].apply(tuple, axis=1)
        mag_tuple2=wisp_phot[['F'+flt+'W', 'F'+flt+'W_ER']].apply(tuple, axis=1)
        mags=mag_tuple1.append(mag_tuple2)
        comb_phot['F'+flt]=mags
        
    ras=hst3d_phot['ra_x'].append(wisp_phot['RA'])
    decs=hst3d_phot['dec_x'].append(wisp_phot['DEC'])
        
    comb_phot['RA']=ras
    comb_phot['DEC']=decs
    comb_phot['survey']=hst3d_phot['survey'].append(wisp_phot['survey'])
    #comb_phot['flags']=hst3d_phot['flags'].append(wisp_phot['flags'])
    
    #strip white spaces from grism_ids #the combination might pose problems
    comb_phot['grism_id']=comb_phot['grism_id']
    indices['grism_id']=indices['grism_id']
    
    indices=indices.drop_duplicates(subset='grism_id')
    comb_phot=comb_phot.drop_duplicates(subset='grism_id')
    
    comb_phot.grism_id=comb_phot.grism_id.apply(lambda x: str(x).replace('.1D.ascii', ' ').strip().lower())
    indices.grism_id=indices.grism_id.apply(lambda x: str(x).replace('.1D.ascii', ' ').strip().lower())
    
    print (comb_phot.grism_id.iloc[0], indices.grism_id.iloc[0])
    master_table=pd.merge(indices, comb_phot, on='grism_id', validate='one_to_one')
    
    # I probably lost tons of objects with grism id ='0000'
    print (master_table.shape, comb_phot.shape, indices.shape)

    
    master_table[['spex_chi', 'line_chi']]=master_table[['spex_chi', 'line_chi']].applymap(np.float)
    master_table['x']=master_table.spex_chi/master_table.line_chi
    master_table['f_test']=master_table.f_test

    #save the result
    #master_table=master_table.join(df)

    #drop the spectrum column because it makes the file heavier
    #master_table=master_table.drop(columns='spectra')

    #make the cut 
    

    return master_table

In [15]:
w_phot.columns


Out[15]:
Index(['Unnamed: 0', 'F140W', 'F140W_ER', 'EXTRACTION_FLAG', 'star_flag',
       'F160W', 'F160W_ER', 'F110W', 'F110W_ER', 'index', 'RA_DEC_NAME',
       'NUMBER', '2', '3', '4', '5', '6', '7', '8', '9', '10', '11', '12',
       '13', '14', 'FLAGS', 'RA', 'DEC', 'FIELD', 'grism_id'],
      dtype='object')

In [16]:
mdf=combined_wisp_hst_catalogs(hst_phot, w_phot, new_indices)


Index(['grism_id', 'star_flag', 'class_star'], dtype='object')
aegis-26-g141_00002 par321-00514
(271915, 28) (279596, 9) (430098, 20)

In [17]:
mdf.shape


Out[17]:
(271915, 29)

In [18]:
#get all the observation info for each field
obs=pd.read_csv(wisps.OUTPUT_FILES+'/observation_log.csv')
obs=obs.drop(columns=['Unnamed: 0']).drop_duplicates(subset='POINTING').reindex()

In [19]:
obs.columns=[x.lower() for x in obs.columns]

In [20]:
obs.columns


Out[20]:
Index(['ra (deg)', 'dec(deg)', 'l (deg)', 'b (deg)', 'exposure (s)',
       'observation date (ut)', 'pointing'],
      dtype='object')

In [21]:
(obs[obs.pointing=='par1'])['exposure (s)'].values[0]


Out[21]:
2609.0

In [22]:
def get_proper_pointing(grism_id):
    grism_id=grism_id.lower()
    if grism_id.startswith('par'):
        return grism_id.split('-')[0]
    else:
        return grism_id.split('-g141')[0]
    
def get_time(pntg):
    try:
        g=(obs[obs.pointing==pntg])[['exposure (s)', 'observation date (ut)']]
        return g.values[0]
    except:
        return (np.nan, np.nan)

def add_pointing_information(row):
    
    pntg=get_proper_pointing(row.grism_id)
    obst,   obsdate = get_time(pntg)
    s3 = pd.Series({'pointing':pntg, 'exposure_time': obst, 'observation_date':obsdate})
    row=row.append(s3)
    return row

In [23]:
obs.columns=obs.columns.str.lower()

In [24]:
mt=mdf.reset_index(drop=True).apply(add_pointing_information, axis=1)

In [25]:
mt.exposure_time.values, mt.shape


Out[25]:
(array([4915., 4915., 4915., ..., 3812., 3812., 3812.]), (271915, 32))

In [26]:
mt.to_hdf(wisps.COMBINED_PHOTO_SPECTRO_FILE, key='all_phot_spec_data')

In [40]:
cands_good=pd.read_pickle(wisps.LIBRARIES+'/candidates.pkl')

In [54]:
#df_ids=wisps.UCD_SPECTRA.grism_id
df_ids=pd.read_pickle(wisps.LIBRARIES+'/candidates_ids.pkl')[0]

In [56]:
cand_df=df=pd.DataFrame()
cand_df['grism_id']=df_ids.apply(lambda x: x.replace('g141', 'G141')).values

In [60]:
objects_ofinter=((cand_df[ cand_df.isin(mt.grism_id.apply(lambda x: x.replace('g141', 'G141')).values)])[['grism_id']])

In [63]:
objects_ofinter.shape, df_ids.shape


Out[63]:
((367, 1), (367,))

In [64]:
objects_ofinter.to_hdf(wisps.LIBRARIES+'/objects_of_interest.hdf', key='good')

In [28]:
stars=mt[mt.class_star !=0.]

In [29]:
stars.shape


Out[29]:
(114477, 32)

In [30]:
stars.to_hdf(wisps.COMBINED_PHOTO_SPECTRO_FILE, key='stars')

In [31]:
st=wisps.Annotator.reformat_table(mt)

In [32]:
snrs=(st[st.F140.between(15.0, 26.)]).snr1.values
mags_F140=(st[st.F140.between(15.0, 26.)]).F140.values
mags_F110=(st[st.F140.between(15.0, 26.)]).F110.values
mags_F160=(st[st.F140.between(15.0, 26.)]).F160.values
xs=(st[st.F140.between(15.0, 26.)]).class_star

In [33]:
idx=np.logical_and(snrs<200., snrs>0.1)

In [34]:
import seaborn as sns

In [35]:
cmap=sns.diverging_palette(150, 275, s=80, l=55, n=9, as_cmap=True)

In [36]:
fig, ax=plt.subplots(figsize=(8, 6))
c=plt.scatter(np.log10(snrs[idx]), mags_F140[idx], marker='+', alpha=0.3, c=xs[idx], cmap=cmap)
plt.xlabel('Log SNR', fontsize=18)
plt.ylabel('F140W', fontsize=18)
#plt.axvline(np.log10(3.), c='k')
plt.minorticks_on()

cbar=plt.colorbar()
plt.tight_layout()
cbar.ax.set_title('CLASS STAR FLAG', fontsize=18)


Out[36]:
Text(0.5, 1, 'CLASS STAR FLAG')

In [37]:
fig, ax=plt.subplots(figsize=(8, 6))
c=plt.scatter(np.log10(snrs[idx]), mags_F110[idx], marker='+', alpha=0.3, c=xs[idx], cmap=cmap)
plt.xlabel('Log SNR', fontsize=18)
plt.ylabel('F110W', fontsize=18)
#plt.axvline(np.log10(3.), c='k')
plt.minorticks_on()

cbar=plt.colorbar()
plt.tight_layout()
cbar.ax.set_title('CLASS STAR FLAG', fontsize=18)


Out[37]:
Text(0.5, 1, 'CLASS STAR FLAG')

In [38]:
fig, ax=plt.subplots(figsize=(8, 6))
c=plt.scatter(np.log10(snrs[idx]), mags_F160[idx], marker='+', alpha=0.3, c=xs[idx], cmap=cmap)
plt.xlabel('Log SNR', fontsize=18)
plt.ylabel('F110W', fontsize=18)
#plt.axvline(np.log10(3.), c='k')
plt.minorticks_on()

cbar=plt.colorbar()
plt.tight_layout()
cbar.ax.set_title('CLASS STAR FLAG', fontsize=18)


Out[38]:
Text(0.5, 1, 'CLASS STAR FLAG')

In [39]:
mt.shape


Out[39]:
(271915, 32)

In [ ]: