In [1]:
# GENERAL PURPOSE PACKAGES
import os
import glob
import tarfile
from urllib.request import urlretrieve
from datetime import date
import timeit
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
In [2]:
# Import ZI tools
%load_ext autoreload
%autoreload 2
# importing ZI tools:
# import ZItools as zit
# importing ZI tools with KT changes:
import KTtools as ktt
In [3]:
# NOT NEEDED ANY MORE, NO MATCHING OF CATALOGS DONE
# import pyspherematch
In [4]:
# MATCH RAD FOR PYSPHEREMATCH
tol_asec = 3. # matching radius in arc.sec
tol_deg = tol_asec/3600.
# MAX PERMITTED MAG ERROR
max_MAGerr = 0.05
In [5]:
homedir = '/home/user/SDSS_SSC/'
datadir = homedir+'DataDir/'
# matched SDSS + PS1 catalogs
nSSC2PS1 = datadir+'nSSC2PS1_matchedv1.csv'
# SDSS 2020: USE THIS VERSION ONLY - FROM ZI ON JUNE, 13, 2020
newSSC2020_cat = datadir+'stripe82calibStars_v3.2.dat'
# DO NOT USE THese earlier VERSIONs
# SDSS 2020: USE THIS VERSION ONLY - FROM ZI ON JUNE, 9, 2020
# newSSC2020_cat = 'stripe82calibStars_v3.1.dat'
# newSSC2020_cat = 'N2020_stripe82calibStars.dat'
# PS1
PS1_cat = datadir+'PS1_STARS_Stripe82_area.csv'
In [6]:
colnames = ['ra','dec','nEpochs','g_Nobs','g_mMed','g_mErr',
'r_Nobs','r_mMed','r_mErr','i_Nobs','i_mMed','i_mErr',
'z_Nobs','z_mMed','z_mErr','raMean','decMean',
'gMeanPSFMag','gMeanPSFMagErr','gMeanPSFMagNpt',
'rMeanPSFMag','rMeanPSFMagErr','rMeanPSFMagNpt',
'iMeanPSFMag','iMeanPSFMagErr','iMeanPSFMagNpt',
'zMeanPSFMag','zMeanPSFMagErr','zMeanPSFMagNpt']
Ndtype = {'nEpochs':'int64','g_Nobs':'int64',\
'r_Nobs':'int64','i_Nobs':'int64','z_Nobs':'int64',\
'gMeanPSFMagNpt':'int64','rMeanPSFMagNpt':'int64',\
'iMeanPSFMagNpt':'int64','zMeanPSFMagNpt':'int64'}
In [7]:
%%time
n2020PS1 = pd.read_csv(nSSC2PS1,header=0,dtype=Ndtype,comment='#')
nrows,ncols = n2020PS1.shape
print('N2020+PS1, as read: num rows, cols: ',nrows,ncols)
In [8]:
n2020PS1.info()
In [9]:
n2020PS1.head()
Out[9]:
In [10]:
### code for generating new quantities, such as dra, ddec, colors, differences in mags, etc
### NOTE: matches IS A DATAFRAME WITH COL NAMES AS GIVEN IN THIS FUNCTIONS
def derivedColumns(matches):
matches['dra'] = (matches['ra']-matches['raMean'])*3600
matches['ddec'] = (matches['dec']-matches['decMean'])*3600
ra = matches['ra']
matches['raW'] = np.where(ra > 180, ra-360, ra)
matches['dg'] = matches['g_mMed'] - matches['gMeanPSFMag']
matches['dr'] = matches['r_mMed'] - matches['rMeanPSFMag']
matches['di'] = matches['i_mMed'] - matches['iMeanPSFMag']
matches['dz'] = matches['z_mMed'] - matches['zMeanPSFMag']
matches['gr'] = matches['g_mMed'] - matches['r_mMed']
matches['ri'] = matches['r_mMed'] - matches['i_mMed']
matches['gi'] = matches['g_mMed'] - matches['i_mMed']
matches['iz'] = matches['i_mMed'] - matches['z_mMed']
matches['dgr'] = matches['dg'] - matches['dr']
matches['dri'] = matches['dr'] - matches['di']
matches['diz'] = matches['di'] - matches['dz']
matches['drz'] = matches['dr'] - matches['dz']
matches['dgi'] = matches['dg'] - matches['di']
return
In [11]:
derivedColumns(n2020PS1)
n2020PS1.head()
Out[11]:
In [12]:
def doOneColor(d, kw):
print('=========== WORKING ON:', kw['Ystr'], '===================')
xVec = d[kw['Xstr']]
yVec = d[kw['Ystr']]
#
# this is where the useable range of color specified by Xstr is defined
# it's really a hack - these limits should be passed via kw...
xMin = 0.5
xMax = 3.1
nBinX = 52
xBin, nPts, medianBin, sigGbin = ktt.fitMedians(xVec, yVec, xMin, xMax,nBinX, 0)
fig,ax = plt.subplots(1,1,figsize=(8,6))
# This is for a scatter plot
# ax.scatter(xVec, yVec, s=0.01, c='blue')
# Replace with hist2D
# get x,y-binning
yMax,yMin = np.quantile(yVec, [0.99,0.01])
yMax,yMin,nBinY = round(yMax,2),round(yMin,2),25
xedges = np.linspace(xMin,xMax,nBinX+1,endpoint=False,dtype=float)
yedges = np.linspace(yMin,yMax,nBinY+1,endpoint=False,dtype=float)
Histo2D, xedges, yedges = np.histogram2d(xVec,yVec, bins=(xedges, yedges))
Histo2D = Histo2D.T
X, Y = np.meshgrid(xedges, yedges)
# cs = ax.pcolormesh(X,Y, Histo2D, cmap='Greys')
# cs = ax.pcolormesh(X,Y, Histo2D, cmap='plasma')
cs = ax.pcolormesh(X,Y, Histo2D, cmap=kw['cmap'])
#ax.scatter(xBin, medianBin, s=5.2, c='yellow')
ax.scatter(xBin, medianBin, s=5.2, c='black', marker='+')
# ax.set_xlim(0.4,3.2)
# ax.set_ylim(-0.5,0.5)
ax.set_xlim(xMin,xMax)
ax.set_ylim(yMin,yMax)
ax.set_xlabel(kw['Xstr'])
ax.set_ylabel(kw['Ystr'])
fig.colorbar(cs, ax=ax, shrink=0.9)
# THERE IS NO ANALYTIC COLOR TERM: linear interpolation of the binned medians!
d['colorfit'] = np.interp(xVec, xBin, medianBin)
# the following line corrects the trend given by the binned medians
d['colorresid'] = d[kw['Ystr']] - d['colorfit']
# note that we only use the restricted range of color for ZP analysis
# 0.3 mag limit is to reject gross outliers
goodC = d[(np.abs(d['colorresid'])<0.3)&(xVec>xMin)&(xVec<xMax)]
### plots
# RA
print(' stats for RA binning medians:')
plotNameRoot = kw['plotNameRoot'] + kw['Ystr']
plotName = plotNameRoot + '_RA.png'
Ylabel =kw['Ystr'] + ' residuals'
kwOC = {"Xstr":'raW', "Xmin":-45, "Xmax":47, "Xlabel":'R.A. (deg)', \
"Ystr":'colorresid', "Ymin":-0.07, "Ymax":0.07,"nBinY":25, "Ylabel":Ylabel, \
"XminBin":-43, "XmaxBin":45, "nBinX":56, "Nsigma":3, "offset":0.01, \
"plotName":plotName, "symbSize":kw['symbSize'], "cmap":kw['cmap']}
ktt.plotdelMag_KT(goodC, kwOC)
print('made plot', plotName)
# Dec
print('-----------')
print(' stats for Dec binning medians:')
plotName = plotNameRoot + '_Dec.png'
kwOC = {"Xstr":'dec', "Xmin":-1.3, "Xmax":1.3, "Xlabel":'Declination (deg)', \
"Ystr":'colorresid', "Ymin":-0.07, "Ymax":0.07,"nBinY":25, "Ylabel":Ylabel, \
"XminBin":-1.26, "XmaxBin":1.26, "nBinX":52, "Nsigma":3, "offset":0.01, \
"plotName":plotName, "symbSize":kw['symbSize'], "cmap":kw['cmap']}
ktt.plotdelMag_KT(goodC, kwOC)
print('made plot', plotName)
# r SDSS
print('-----------')
print(' stats for SDSS r binning medians:')
plotName = plotNameRoot + '_rmag.png'
kwOC = {"Xstr":'r_mMed', "Xmin":14.3, "Xmax":22.2, "Xlabel":'SDSS r (mag)', \
"Ystr":'colorresid', "Ymin":-0.07, "Ymax":0.07,"nBinY":25, "Ylabel":Ylabel, \
"XminBin":14.5, "XmaxBin":21.5, "nBinX":55, "Nsigma":3, "offset":0.01, \
"plotName":plotName, "symbSize":kw['symbSize'], "cmap":kw['cmap']}
ktt.plotdelMag_KT(goodC, kwOC)
print('made plot', plotName)
print('------------------------------------------------------------------')
In [14]:
# These keywords what you want to plot on x
keywords = {"Xstr":'gi', "plotNameRoot":'colorResidPS12_'}
# These keywords are for plot style
keywords["symbSize"] = 0.05
keywords["cmap"] = 'plasma' # gives a good dark background for the overplots
# keywords["cmap"] = 'PuBuGn' # the overplots look washed out
# create a series of plots
# for color in ('dg', 'dr', 'di', 'dz', 'dgr', 'dri', 'diz'):
for color in ('dr', 'dgr'): # use this for testing
keywords["Ystr"] = color
doOneColor(n2020PS1, keywords)
In [ ]: