We need about 10,000 galaxies from SDSS, in a narrow redshift range around z = 0.1, to match with our halos.
We can grab this data from the SDSS skyserver SQL server from python, using the mechanize
web browser, and then manipulate the catalog with pandas
.
In [1]:
# Querying skyserver, using code from Eduardo Martin's blog post at
# http://balbuceosastropy.blogspot.com/2013/10/an-easy-way-to-make-sql-queries-from.html
import pandas as pd
def SDSS_select(sql):
from StringIO import StringIO # To read a string like a file
import mechanize
url = "http://skyserver.sdss3.org/dr10/en/tools/search/sql.aspx"
br = mechanize.Browser()
br.open(url)
br.select_form(name="sql")
br['cmd'] = sql
br['format']=['csv']
response = br.submit()
file_like = StringIO(response.get_data())
return pd.read_csv(file_like, skiprows=1)
We want a sample of galaxies and their 5-band photometry, plus their redshifts and possibly positions.
In [2]:
galaxies = "SELECT top 10000 \
objid, ra, dec, z, \
dered_u AS mag_u, \
dered_g AS mag_g, \
dered_r AS mag_r, \
dered_i AS mag_i, \
dered_z AS mag_z \
FROM SpecPhoto \
WHERE \
(class = 'Galaxy')"
print galaxies
This SQL checks out at the SkyServer SQL box, which has a syntax checking option.
In [3]:
data = SDSS_select(galaxies)
data
Out[3]:
Let's see if we can see the cosmic web. Initial plotting code from Josh Bloom's AstroHackWeek 2014 lecture.
In [4]:
# For pretty plotting
!pip install --upgrade seaborn
In [5]:
%pylab inline
import seaborn
seaborn.set()
import copy
First let's see the redshift distribution.
In [6]:
bins = hist(data['z'].values,bins=100) ; xlabel("redshift") ; ylabel("N")
Out[6]:
Now let's look at everything else, colored by redshift.
In [7]:
import matplotlib as mpl
import matplotlib.cm as cm
data = data[(data['mag_r'] > -9999)]
## truncate the color at z=0.6 just to keep some contrast.
norm = mpl.colors.Normalize(vmin=min(data['z'].values), vmax=0.6)
cmap = cm.jet_r
m = cm.ScalarMappable(norm=norm, cmap=cmap)
rez = pd.scatter_matrix(data[0:2000], alpha=0.2,figsize=[15,15],color=m.to_rgba(data['z'].values))
Save this for future use.
In [8]:
data.to_csv('galaxies.clean.csv')
In [ ]: