Section 1.5: Surveys and Data-sets Used in Examples
Adrian M. Price-Whelan
To see what datasets are available, first import the datasets module:
import astroML.datasets as datasets
Remove the comment below, then go to the end of the line (after the period) and hit tab
:
# datasets.
SDSS
Photometry
from IPython.display import Image
Image("http://www.sdss3.org/images/dr8_spectro_coverage.png", width=800)
sdss_data = datasets.fetch_imaging_sample()
figure(figsize=(4,4))
plot(sdss_data['ra'], sdss_data['dec'],
marker='.', linestyle='none', alpha=0.1)
xlim(360, 0)
ylim(-90, 90)
(-90, 90)
figure(figsize=(12,6))
plot(sdss_data['ra'], sdss_data['dec'],
marker='.', linestyle='none', alpha=0.05)
title("Source density")
xlim(10.5, -0.5)
ylim(-2, 2)
(-2, 2)
Now, as an example, we separate the stars and galaxies and plot a color-magnitude diagram:
stars = sdss_data[sdss_data['type'] == 6]
galaxies = sdss_data[sdss_data['type'] == 3]
stars.size, galaxies.size
(82003, 248750)
fig,axes = subplots(1,2,figsize=(12,4),
sharex=True,sharey=True)
axes[0].plot(stars['gRaw']-stars['rRaw'],
stars['rRaw']-stars['iRaw'],
alpha=0.02, marker='.', linestyle='none')
axes[0].set_xlabel(r"$g-r$")
axes[0].set_ylabel(r"$r-i$")
axes[0].set_title("Stars")
axes[1].plot(galaxies['gRaw']-galaxies['rRaw'],
galaxies['rRaw']-galaxies['iRaw'],
alpha=0.02, marker='.', linestyle='none')
axes[1].set_xlabel(r"$g-r$")
axes[1].set_title("Galaxies")
axes[0].set_xlim(-0.5, 3)
axes[0].set_ylim(-1, 2)
(-1, 2)
Example: Transforming to Galactic coordinates with Astropy
from astropy.coordinates import ICRS
import astropy.units as u
eq = ICRS(stars['ra']*u.deg, stars['dec']*u.deg)
l,b = eq.galactic.l, eq.galactic.b
plot(l.degree, b.degree,
marker='.',linestyle='none',alpha=0.05)
xlim(95,117)
ylim(-64,-59)
xlabel('$l$ [deg]')
ylabel('$b$ [deg]')
<matplotlib.text.Text at 0x1103847d0>
Spectroscopy
spectrum = datasets.fetch_sdss_spectrum(2872, 54533, 71)
from scipy.ndimage import gaussian_filter
figure(figsize=(12,6))
plot(spectrum.wavelength(), gaussian_filter(spectrum.spectrum,2),
drawstyle='steps', alpha=0.75, color='k')
xlim(min(spectrum.wavelength()), max(spectrum.wavelength()))
xlabel(r"Wavelength [$\AA$]")
ylabel(r"Flux")
<matplotlib.text.Text at 0x1105ae3d0>
DR7 quasar catalog example
qso = datasets.fetch_dr7_quasar()
First we plot the redshit
figure(figsize=(8,5))
n,bins,patches = hist(qso['redshift'], bins=100,
alpha=0.75)
xlabel('Redshift')
<matplotlib.text.Text at 0x1106c5c90>
SEGUE
segue = datasets.fetch_sdss_sspp()
eq = ICRS(segue['ra']*u.deg, segue['dec']*u.deg)
l,b = eq.galactic.l, eq.galactic.b
By default, galactic longitude wraps at 360 -- here, we switch to (-180,180)
l = l.wrap_at(180*u.deg)
figure(figsize=(12,8))
subplot(projection='hammer')
plot(l.radian, b.radian,
marker='.',linestyle='none',alpha=0.05)
xlabel('$l$ [deg]')
ylabel('$b$ [deg]')
<matplotlib.text.Text at 0x110772e90>
/Users/adrian/projects/matplotlib/build/lib.macosx-10.9-intel-2.7/matplotlib/cbook.py:1730: DeprecationWarning: using a non-integer number instead of an integer will result in an error in the future result = np.zeros(new_shape, a.dtype)
Example: scatter plot of RV vs. Metallicty
repr(segue.dtype.names)
"('ra', 'dec', 'Ar', 'upsf', 'uErr', 'gpsf', 'gErr', 'rpsf', 'rErr', 'ipsf', 'iErr', 'zpsf', 'zErr', 'pmL', 'pmB', 'pmErr', 'mjd', 'plate', 'fiber', 'SNR', 'FeH', 'FeHErr', 'Teff', 'TeffErr', 'logg', 'loggErr', 'alphFe', 'alphFeErr', 'radVel', 'radVelErr')"
good = ((segue['radVelErr']/segue['radVel']) < 0.1) & \
((segue['FeHErr']/segue['FeH']) < 0.1) & \
((segue['alphFeErr']/segue['alphFe']) < 0.1) & \
(segue['rpsf'] > 15) & (segue['rpsf'] < 19)
good_segue = segue[good]
good_segue.size
203204
Next I'm going to use some of my own code to transform these heliocentric radial velocities into Galactocentric:
from streams.coordinates import vhel_to_vgsr
eq = ICRS(good_segue['ra']*u.deg, good_segue['dec']*u.deg)
l,b = eq.galactic.l, eq.galactic.b
vgsr = vhel_to_vgsr(l, b, good_segue['radVel']*u.km/u.s)
plot(good_segue['FeH'], vgsr.value,
marker='.',alpha=0.02,linestyle='none')
xlim(-3, 0.1)
ylim(-400, 400)
xlabel(r'[Fe/H]')
ylabel(r'$V_{\rm gsr}$')
<matplotlib.text.Text at 0x1108981d0>
Here I plot the alpha abundance versus metallcity:
plot(good_segue['FeH'], good_segue['alphFe'],
marker='.',alpha=0.03,linestyle='none')
xlabel(r'[Fe/H]')
ylabel(r'[$\alpha$/Fe]')
xlim(-3,0)
ylim(0.05,0.6)
(0.05, 0.6)
LINEAR
linear = datasets.fetch_LINEAR_sample()
trg = linear.targets
sesar2013 = ((trg['raLIN'] > 309) | (trg['raLIN'] < 60)) & \
(np.abs(trg['decLIN']) < 1.23) & \
(trg['ug'] > 0.75) & (trg['ug'] < 1.45) & \
(trg['gr'] > -0.25) & (trg['gr'] < 0.4) & \
(trg['ri'] > -0.2) & (trg['ri'] < 0.2) & \
(trg['iz'] > -0.3) & (trg['iz'] < 0.3)
sesar2013 = trg[sesar2013]
sesar2013.size
154
rrlyr = sesar2013[0]
P = 10**rrlyr['LP1']
t,mag,err = linear[rrlyr['objectID']].T
phase = (t / P + 0.1) % 1
figure(figsize=(12,6))
errorbar(phase, mag, err, linestyle='none',
marker='o', ecolor='#888888')
xlabel("Phase")
ylabel("mag")
<matplotlib.text.Text at 0x10f68ab10>