Import some packages for data handling and plotting. Also import datalab code.
# Python imports
import numpy as N
import pandas as pd
from cStringIO import StringIO
import pylab as p
import matplotlib
%matplotlib inline
# Datalab imports
from dl import authClient, queryClient
# authentication token; use it for all your queries below
token = authClient.login('anonymous')
# query DB for a list of all SMASH fields
query = "SELECT ra,dec,nexp,nobj FROM smash_dr1.field WHERE nexp IS NOT NULL AND nobj IS NOT NULL" # define sql query as a string
result = queryClient.query(token,adql=query,fmt='csv') # submit the query, using your authentication token
print "Inspect the returned result: it's a csv-formatted string:"
print result[:105]
Inspect the returned result: it's a csv-formatted string: ra,dec,nexp,nobj 103.81954,-36.918889,29,1751546 145.212,-26.9897,32,1538828 149.082,-27.9336,46,2313654
For easier handling, load the columns into arrays.
fi = StringIO(result) # creates an in-RAM file (not on disk, just in RAM)
ra,dec,nexp,nobj = N.loadtxt(fi,unpack=True,skiprows=1,delimiter=',') # reads columns from that 'file'
import astropy.units as u # units module
from astropy.coordinates import SkyCoord
# use a helper, don't re-invent the wheel
c = SkyCoord(ra=ra*u.degree, dec=dec*u.degree, frame='icrs')
ra = c.ra.wrap_at(180 * u.deg).radian # ra and ...
dec = c.dec.radian # ... dec now in radians
p.figure(figsize=(8,4.2))
p.subplot(111, projection="aitoff")
im = p.scatter(ra, dec, marker='o', s=nexp, c=N.log10(nobj), edgecolors='none',alpha=1)
cb = p.colorbar(im)
cb.set_label('log10 number of sources in field')
p.title(r"SMASH field positions. Blob size $\propto$ nexp, color $\propto$ log10(nobj)",y=1.08)
p.grid(True)
p.xlabel('RA')
p.ylabel('DEC')
<matplotlib.text.Text at 0x7fa1dc9d6310>
Select from the object
table all g and r mags of all sources in a single SMASH field. Apply some quality criteria.
fieldID = 169 # a specific SMASH field
# note: the DB encodes 'no measurement' values as 99.99
query = "SELECT gmag,rmag FROM smash_dr1.object WHERE ABS(sharp)<0.5 AND fieldid=%d AND gmag<99 AND rmag<99" % fieldID
result = queryClient.query(token,adql=query,fmt='csv') # query the DB
The returned result is a long, csv-formatted string. Let's place it into a Pandas dataframe for easier access.
df = pd.read_csv(StringIO(result))
print df.columns # lists the names of the columns in dataframe
Index([u'gmag', u'rmag'], dtype='object')
You can easily save the table to a csv file, for later use.
filename = 'smash_avg_photo_field%d.csv' % fieldID
df.to_csv(filename,index=False) # don't write the index column
# get the file size
import os
filesize = os.path.getsize(filename) / 1024.**2 # in MB
# display a download link
from IPython.display import HTML
HTML('<a href="%s">Click on this link to download file: %s (%.1f MB)</a>' % (filename,filename,filesize))
# make a figure
fig = p.figure(figsize=(11,5))
ax1 = fig.add_subplot(121)
ax2 = fig.add_subplot(122)
# mag histograms
ax1.hist(df['gmag'],bins=100,histtype='step',color='b',label='g',lw=2,alpha=0.6)
ax1.hist(df['rmag'],bins=100,histtype='step',color='r',label='r',lw=2,alpha=0.6)
ax1.set_yscale('log')
ax1.set_xlabel('mag')
ax1.set_ylabel('number of sources')
ax1.set_title('Magnitude histograms')
ax1.legend(loc='upper left',frameon=False)
# color-magnitude diagram (Hess diagram)
im = ax2.hexbin(df['gmag']-df['rmag'],df['gmag'],gridsize=200,cmap=matplotlib.cm.viridis,norm=matplotlib.colors.LogNorm())
cb = p.colorbar(im,aspect=40,pad=0.01)
cb.set_label('number of sources per 2D bin')
ax2.set_xlabel('g - r')
ax2.set_ylabel('g')
ax2.set_title('Color-magnitude (Hess) diagram')
ax2.set_xlim(-3,3)
ax2.set_ylim(28,13)
(28, 13)
We can retrieve all calibrated magnitudes of a specific object, for all bands, and all epochs, from the source
table.
objID = '169.429960' # this happens to be a RR Lyrae star in Hydra II
# select columns Julian date, calibrated mag, filter
# note: the DB encodes 'no measurement' values as 99.99
query = "SELECT mjd,cmag,filter FROM smash_dr1.source WHERE id='%s' AND cmag<99 ORDER BY mjd ASC" % objID
result = queryClient.query(token,adql=query,fmt='csv') # run the query
df = pd.read_csv(StringIO(result)) # load to Pandas dataframe
Plot the light-curves of this source separately for every band. We also zoom in on a range of dates with interesting measurements.
colors = ('b','g','y','m','r') # define plotting colors
mjd0 = 57111 # zero-point Julian date
for j,filt in enumerate(list('ugriz')):
sel = (df['filter'] == filt) # selection: only one filter
p.plot(df['mjd'][sel]-mjd0,df['cmag'][sel],marker='.',ms=4,ls='-',lw=1,color=colors[j],alpha=0.4,label=filt)
p.xlabel('Julian date - %d' % mjd0)
p.ylabel('cmag')
p.xlim(0.9,2.4)
p.ylim(22.2,20.9)
p.legend(loc='lower center',frameon=True)
<matplotlib.legend.Legend at 0x7fa1df02cf10>
In this range of dates there were no u and z-band measurements.