In [2]:
import os
import numpy as np
import numpy.ma as ma
from astropy.io import ascii
import matplotlib.pyplot as plt
from astroML.stats import binned_statistic_2d
import seaborn as sns
#from scipy.stats import gaussian_kde
%matplotlib inline

sns.set(style='white', font_scale=1.6, palette='deep')


:0: FutureWarning: IPython widgets are experimental and may change in the future.

In [4]:
# Read the full catalog
dat = ascii.read(os.path.join(os.getenv('HOME'),'dr2-grz.txt'), 
                 format='fixed_width', delimiter='|', data_start=2)#, data_end=10000)
dat


Out[4]:
<Table masked=True length=4471514>
idradecgfluxrfluxzfluxgflux_ivarrflux_ivarzflux_ivar
str11str17float64float64float64float64float64float64float64
328187035.563262479-5.0759236563253.7093146.618280.66848.306242.950616.8985
328202635.5448541701-5.066507045231.791727.669827.6593185.748154.07154.2688
328203335.5574103117-5.065634432372.776724.340685.83656181.597159.67256.8928
328204035.5344433281-5.064822051035.0967922.5321110.0169.803110.53446.4268
328205435.577066248-5.0633803081837.461880.7769143.24462.76760.740421.8262
328218135.5557897738-5.057487781172.728557.0215110.25695.43386.717627.0506
328221835.5891039988-5.0513569478929.857678.8151171.00659.50253.458319.385
328222735.5921434121-5.054230690161.684026.8567815.8476173.971147.29450.8141
328227135.5878112882-5.052667110553.453829.9903819.8505165.662142.18450.3153
328227735.5595275713-5.0509455659725.183940.914465.73595.41397.48937.8232
...........................
262426803318.94188087910.12261746641.394274.849949.69338502.873407.10474.3937
262426806319.0555371910.123102280511.172925.655940.9123369.168245.99981.3191
262426814319.09743518210.12366668751.32757.1375217.3142545.588338.34161.3776
262426817319.10306277410.12408178993.047518.1576819.963662.328638.7335.48922
262426821319.07464008910.124764516437.0958.387675.3267174.318140.46362.8838
262426822319.07635544610.12486900551.57634.349827.27934699.338553.738115.718
262426827319.04200278210.12361328485.6302821.9844120.753503.369275.15149.7688
262426838319.13508023310.12443738791.956626.6577414.3401707.312484.282103.135
262426847318.96611555810.12472185943.3522413.352438.5873590.178360.4883.4146
(4471513 ros)--------------

In [3]:
# Do some juggling to deal with masked values and then define the g-r and r-z colors.
rcut = 20.0
dat['gflux'].fill_value = -99
dat['rflux'].fill_value = -99
dat['zflux'].fill_value = -99
good = np.where((dat['gflux']>0)*(dat['rflux']>0)*(dat['zflux']>0)*(dat['rflux']>10**(0.4*(22.5-rcut)))*1)[0]
gr = -2.5*np.log10(np.array(dat['gflux'][good]/dat['rflux'][good]))
rz = -2.5*np.log10(np.array(dat['rflux'][good]/dat['zflux'][good]))
print(len(gr))


---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
<ipython-input-3-c95ef5a2e9a3> in <module>()
      1 # Do some juggling to deal with masked values and then define the g-r and r-z colors.
      2 rcut = 20.0
----> 3 dat['gflux'].fill_value = -99
      4 dat['rflux'].fill_value = -99
      5 dat['zflux'].fill_value = -99

NameError: name 'dat' is not defined

In [1]:
# Make a fancy histogram - http://www.astroml.org/examples/datasets/plot_SDSS_SSPP.html
xlim = (-0.7, 2.0)
ylim = (-0.7, 2.0)
N, xedges, yedges = binned_statistic_2d(rz, gr, gr, 'count', 
                                        bins=100, range=[xlim, ylim])
N[N<10] = 0
print(np.log10(N).max())
# Define custom colormaps: Set pixels with no sources to white
cmap = plt.cm.jet
cmap.set_bad('w', 1.)
cmap_multicolor = plt.cm.jet
cmap_multicolor.set_bad('w', 1.)


---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
<ipython-input-1-8e89e1dac3ee> in <module>()
      2 xlim = (-0.7, 2.0)
      3 ylim = (-0.7, 2.0)
----> 4 N, xedges, yedges = binned_statistic_2d(rz, gr, gr, 'count', 
      5                                         bins=100, range=[xlim, ylim])
      6 N[N<10] = 0

NameError: name 'binned_statistic_2d' is not defined

In [35]:
# Create figure and subplots
fig = plt.figure(figsize=(10, 6))
fig.subplots_adjust(wspace=0.25, left=0.1, right=0.95,
                    bottom=0.07, top=0.95)

#ax.scatter(rz, gr, c=nn, s=50, edgecolor='', cmap=cm.Blues)
#sns.kdeplot(rz, gr, clip=(xlim, ylim), ax=ax, gridsize=40,
#           cmap="Reds", shade=True, cut=0, shade_lowest=False)
#ax.scatter(gr[:100], rz[:100], marker='s')
plt.subplot(111, xticks=[4000, 5000, 6000, 7000, 8000])
plt.imshow(np.log10(N.T), origin='lower', 
          extent=[xedges[0], xedges[-1], yedges[0], yedges[-1]],
          aspect='auto', interpolation='nearest', cmap=cmap)
levels = np.linspace(0, np.log10(N.max()), 7)[2:]
plt.contour(np.log10(N.T), levels, colors='k', linewidths=1,
            extent=[xedges[0], xedges[-1], yedges[0], yedges[-1]])
plt.xlim(xedges[0], xedges[-1])
plt.ylim(yedges[0], yedges[-1])
plt.xlabel('r - z')
plt.ylabel('g - r')
plt.colorbar()
#cb = plt.colorbar(ticks=[0, 1, 2, 3], format=r'$10^{%i}$')
#cb.set_label(r'$\mathrm{Number\ of\ objects}$')
#plt.clim(0, 3)
#plt.tight_layout()
plt.show()



In [ ]: