<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Astropy Smorgasbord
Astropy is a the python package designed to contain the basic utilities needed for astronomy.  It has been many years in the making and replicated some of the functionaily of other astronomy packages from IRAF and IDL.  However, it still is missing many packages and utilities needed for say reduction, handling spectra, etc. As it is a community driven effort, its unclear how it will take for additional functionality to be included.
Tables
Data tables are convenient way of dealing with large data sets, particularly if you need to create, read, modifiy and/or write files (i.e. ASCII, FITS, HDF5, etc...).  Both Pandas and Astropy have their own unique support of data tables and it is possible to convert back and forth betweeen them.  The question arises why are there two different types of data table formats and what are the benifits of using one over the other. 
Pandas strengths are within its statistical machinary available for data analysis. Whereas the Astropy table has a few astronomy specific advantages; 1) being able to work with astronomy formats (i.e. FITS, sextractor, CDS, VOTables, etc) and 2) being able to use units within the table.  It may becomes necessary to switch back and forth between tables depending on the task at hand.  
Here we will give you a very brief introduction to Astropy tables.
First you need to load the following module inorder to utilize an Astropy table.
Step1:  Here we will create a QTable which is a quantity table.  The first column is the redshift of the source, the second is the star formation rate (SFR), the third is the stellar mass and the fourth is the rotational velocity. Since it is a QTable, we can define the units for each of the entries, if we like.
Step2:  When we create the table, we can assign the names to each column
Step3:  We can display the table by typing the variable into the notebook
Step4:  We can display the information about the table (i.e. data type, units, name)
Step5:  We can view a specific column by calling is name.
Step6:  We can also group subsets of columns.
Step7:  We can call a specific row, in this case the first row, which shows all the entries for the first object.  Remember python indexing starts with zero!
Step8:  We can write to a standard format if we like, binary FITS tables are often used for this. If could be an ASCII file (either sapace, tab or comma delimited) or whatever your favorite sharable file format is.  There is a whole list of supported formats that can be read in, many astronomy specific.
Step9:  You can read it in just as easily too. As you can see, the units are stored as well, as long as you are using a QTable.
Step10:  We can add new entries to the table.  In this case specific star formation rate (sSFR).
Step11:  Maybe we are unhappy using years and we want to convert to Gyr. There are a whole list of available units here.
Step12:  If we want to estimate the H-alpha line lunminosity from the SFR using an emperical relation (such as one from Kennicutt et al. 1998).
Step13:  Say we want to go from erg to photons, which are both "energy".  You can see it is not quite as straight forward by the follwoing error.
Step14:  In order to convert to photons/s properly, we need to check the units. Our beginning and ending units match (i.e. energy/time), we just need to acount for the energy of the photon which is depends on the wavelength.
$$ E_{phot} = \frac{h c}{\lambda} $$
(1) As an exercise for the reader. Modify the above table to estimate the number of photons/s from the H-alpha line luminosity the energy of the photon instead.
Photometry with Photutils
What is photometry?  It is the measurement of light from a source, either represented as flux or magnitude. There are several ways to compute the flux from a source. In this tutorial we will go over aperture photometry, where you draw a circular region around the source and measure the flux inside.
Photoutils is a useful astropy package for performing photometry.  This includes but not limited to aperture photometry, PSF photometry, star finding, source isophotes, background subtraction, etc...
Large parts of this tutorial are borrowed from Getting Started with Photutils documentation.
Here we load in a star field image provided with Photutils.  We are interested in a small subsection of the image in which we want to do photometry.
Step15:  The type function is useful for determining what data type you are dealing with, especially if it is an unfamilar one.
Step16:  Assuming your background is uniform (and it most certianly is not), here is a rough background subtraction. 
This may not work for all astronomical fields (e.g. clusters, low surface brightness features, etc). There are more sophisticated ways to deal with the background subtraction, see the Background Estimation documentation.
Step17:  We use DAOStarFinder to find all of the objects in the image. 
The FWHM (Full Width Half Maximum) is a way of measureing the size of the source. Often for unresolved sources, like stars, the FWHM should be matched the resolution of the telescope (or the seeing).  
The threshold is the value above the background in which the star finding algorithm will find sources. Often this should be done in a statistical sense in which a real detection is n-sigma above the background (or n standard deviations). We can measure the standard deviation of background. 
Statistics aside
Step18:  Astroquery
Astroquery is a Astropy package for querying astronomical databases.  There many, many, many  catalogs and image services.  Astroquery can be used to find catalogs and images.  In this tutorial, we will show a few different databases and and ways to search for data.
NED query
NED is the NASA/IPAC Extragalacitc Database.  It can be useful for finding measurements, properties and papers about
your favorite galaxy (or galaxies).  Like many databases you can search for galaxies via their name or coordinates.
First we need to import the modules that we will use for this tutorial.
Step19:  Print the columns returned from the search.
Step20:  Print the redshift of the source
Step21:  Lets now search for objects around M82.
Step22:  Gaia query
Gaia is a space observatory designed to observe the precise positions of stars in order to determine their distances and motions within the galaxy.  Gaia star coordinates are the new golden standard for astrometry.  Being able to querry their database to correct your catalogs and images to their coordinate system would be highly beneficial.
Step23:  IRSA query of 2MASS
IRSA is the NASA/IPAC Infrared Science Archive.  IRSA typically has a infrared mission focus, such as 2MASS, Herschel, Planck, Spitzer, WISE, etc...
Step24:  SDSS query
SDSS or the Sloan Digital Sky Survey is a massive survey of the sky which originally was optical imaging and spectroscopy.  It has been operating for over 20 years and the well beyond the original planned survey as it has expanded to near-infrared observations of stars and IFUs (integral field units) for observing galaxies and eventually the whole sky LVM (Local Volume Mapper). It is an immense data set to draw from.
Step25:  Coordinates of your favorite object, in this case M82 in degrees.
Step26:  Skyview
Skyview is a virtual observatory that is connected to several surveys and imaging data sets spanning multiple wavelengths. It can be a handy way to get quick access to multi-wavelength imaging data for your favorite object.
Step27:  Plot image
Step28:  Plot the in RA and Declination coordinates space, instead of pixels. The WCS (or World Coordinates System) in a FITS file has a set of keywords that define the coordinate system.  
Typical keywords are the following
Step29:  Use the FITS header to define the WCS of the image.
<ASSISTANT_TASK:>
Python Code:
from astropy.table import QTable
import astropy.units as u
import numpy as np
Explanation: Astropy Smorgasbord
Astropy is a the python package designed to contain the basic utilities needed for astronomy.  It has been many years in the making and replicated some of the functionaily of other astronomy packages from IRAF and IDL.  However, it still is missing many packages and utilities needed for say reduction, handling spectra, etc. As it is a community driven effort, its unclear how it will take for additional functionality to be included.
Tables
Data tables are convenient way of dealing with large data sets, particularly if you need to create, read, modifiy and/or write files (i.e. ASCII, FITS, HDF5, etc...).  Both Pandas and Astropy have their own unique support of data tables and it is possible to convert back and forth betweeen them.  The question arises why are there two different types of data table formats and what are the benifits of using one over the other. 
Pandas strengths are within its statistical machinary available for data analysis. Whereas the Astropy table has a few astronomy specific advantages; 1) being able to work with astronomy formats (i.e. FITS, sextractor, CDS, VOTables, etc) and 2) being able to use units within the table.  It may becomes necessary to switch back and forth between tables depending on the task at hand.  
Here we will give you a very brief introduction to Astropy tables.
First you need to load the following module inorder to utilize an Astropy table.
End of explanation
a = [0.10, 0.15, 0.2]
b = [10.0, 2.0,  100.0] * u.M_sun / u.yr
c = [1e10, 1e9,  1e11]  * u.M_sun
d = [150., 100., 2000.] * u.km / u.s
Explanation: Here we will create a QTable which is a quantity table.  The first column is the redshift of the source, the second is the star formation rate (SFR), the third is the stellar mass and the fourth is the rotational velocity. Since it is a QTable, we can define the units for each of the entries, if we like.
End of explanation
t = QTable([a, b, c, d],
           names=('redshift', 'sfr', 'stellar_mass', 'velocity'))
Explanation: When we create the table, we can assign the names to each column
End of explanation
t
Explanation: We can display the table by typing the variable into the notebook
End of explanation
t.info
Explanation: We can display the information about the table (i.e. data type, units, name)
End of explanation
t['sfr']
Explanation: We can view a specific column by calling is name.
End of explanation
t['sfr','velocity']
Explanation: We can also group subsets of columns.
End of explanation
t[0]
Explanation: We can call a specific row, in this case the first row, which shows all the entries for the first object.  Remember python indexing starts with zero!
End of explanation
t.write('my_table.fits', overwrite=True) 
Explanation: We can write to a standard format if we like, binary FITS tables are often used for this. If could be an ASCII file (either sapace, tab or comma delimited) or whatever your favorite sharable file format is.  There is a whole list of supported formats that can be read in, many astronomy specific.
End of explanation
all_targets = QTable.read("my_table.fits")
all_targets
Explanation: You can read it in just as easily too. As you can see, the units are stored as well, as long as you are using a QTable.
End of explanation
t['ssfr'] = t['sfr'] / t['stellar_mass']
t
Explanation: We can add new entries to the table.  In this case specific star formation rate (sSFR).
End of explanation
t['ssfr'].to(1./u.gigayear)
Explanation: Maybe we are unhappy using years and we want to convert to Gyr. There are a whole list of available units here.
End of explanation
t['lum'] = (t['sfr'] * u.erg / u.s )/(7.9e-42 * u.Msun / u.yr)
t
Explanation: If we want to estimate the H-alpha line lunminosity from the SFR using an emperical relation (such as one from Kennicutt et al. 1998).
End of explanation
(t['lum']).to(u.ph / u.s)
Explanation: Say we want to go from erg to photons, which are both "energy".  You can see it is not quite as straight forward by the follwoing error.
End of explanation
import numpy as np
from photutils import datasets
hdu = datasets.load_star_image()  
image = hdu.data[500:700, 500:700].astype(float)
head = hdu.header
Explanation: In order to convert to photons/s properly, we need to check the units. Our beginning and ending units match (i.e. energy/time), we just need to acount for the energy of the photon which is depends on the wavelength.
$$ E_{phot} = \frac{h c}{\lambda} $$
(1) As an exercise for the reader. Modify the above table to estimate the number of photons/s from the H-alpha line luminosity the energy of the photon instead.
Photometry with Photutils
What is photometry?  It is the measurement of light from a source, either represented as flux or magnitude. There are several ways to compute the flux from a source. In this tutorial we will go over aperture photometry, where you draw a circular region around the source and measure the flux inside.
Photoutils is a useful astropy package for performing photometry.  This includes but not limited to aperture photometry, PSF photometry, star finding, source isophotes, background subtraction, etc...
Large parts of this tutorial are borrowed from Getting Started with Photutils documentation.
Here we load in a star field image provided with Photutils.  We are interested in a small subsection of the image in which we want to do photometry.
End of explanation
print(type(hdu))
print(type(image))
print(type(head))
Explanation: The type function is useful for determining what data type you are dealing with, especially if it is an unfamilar one.
End of explanation
image -= np.median(image) 
Explanation: Assuming your background is uniform (and it most certianly is not), here is a rough background subtraction. 
This may not work for all astronomical fields (e.g. clusters, low surface brightness features, etc). There are more sophisticated ways to deal with the background subtraction, see the Background Estimation documentation.
End of explanation
from photutils import DAOStarFinder
from astropy.stats import mad_std
bkg_sigma = mad_std(image)
daofind = DAOStarFinder(fwhm=4., threshold=3.*bkg_sigma)  
sources = daofind(image)
for col in sources.colnames: sources[col].info.format = '%.8g'  # for consistent table output
print(sources[:10])
print(sources.info)
from photutils import aperture_photometry, CircularAperture
positions = np.transpose((sources['xcentroid'], sources['ycentroid']))  
apertures = CircularAperture(positions, r=4.)  
phot_table = aperture_photometry(image, apertures)
for col in phot_table.colnames: phot_table[col].info.format = '%.8g'  # for consistent table output
print(phot_table[:10])
print(phot_table.info)
import matplotlib.pyplot as plt
plt.imshow(image, cmap='gray_r', origin='lower')
apertures.plot(color='blue', lw=1.5, alpha=0.5)
print(head.keys())
Explanation: We use DAOStarFinder to find all of the objects in the image. 
The FWHM (Full Width Half Maximum) is a way of measureing the size of the source. Often for unresolved sources, like stars, the FWHM should be matched the resolution of the telescope (or the seeing).  
The threshold is the value above the background in which the star finding algorithm will find sources. Often this should be done in a statistical sense in which a real detection is n-sigma above the background (or n standard deviations). We can measure the standard deviation of background. 
Statistics aside:
Assuming a Gaussian distribution for the noise, 1-sigma contains 68% of the distribution, 2-sigma contains 95% of the distribution and 3-sigma contains 99.7% of the distribution. Often 5-sigma is the golden standard in physics and astronomy for a "real" detection as it is contains 99.99994% background distribution.  
Here we use 3-sigma (3 standard devations) above the background for a "real" detection.
End of explanation
from astroquery.ned import Ned
from astropy import coordinates
import astropy.units as u
result_table = Ned.query_object("m82")
print(result_table)
Explanation: Astroquery
Astroquery is a Astropy package for querying astronomical databases.  There many, many, many  catalogs and image services.  Astroquery can be used to find catalogs and images.  In this tutorial, we will show a few different databases and and ways to search for data.
NED query
NED is the NASA/IPAC Extragalacitc Database.  It can be useful for finding measurements, properties and papers about
your favorite galaxy (or galaxies).  Like many databases you can search for galaxies via their name or coordinates.
First we need to import the modules that we will use for this tutorial.
End of explanation
print(result_table.keys())
Explanation: Print the columns returned from the search.
End of explanation
print(result_table['Redshift'])
Explanation: Print the redshift of the source
End of explanation
result_table = Ned.query_region("m82", radius=1 * u.arcmin)
print(result_table)
Explanation: Lets now search for objects around M82.
End of explanation
from astroquery.gaia import Gaia
coord = coordinates.SkyCoord(ra=280, dec=-60, unit=(u.degree, u.degree), frame='icrs')
radius = u.Quantity(1.0, u.arcmin)
j = Gaia.cone_search_async(coord, radius)
r = j.get_results()
r.pprint()
print(r.keys())
Explanation: Gaia query
Gaia is a space observatory designed to observe the precise positions of stars in order to determine their distances and motions within the galaxy.  Gaia star coordinates are the new golden standard for astrometry.  Being able to querry their database to correct your catalogs and images to their coordinate system would be highly beneficial.
End of explanation
from astroquery.irsa import Irsa
Irsa.list_catalogs()
table = Irsa.query_region("m82", catalog="fp_psc", spatial="Cone",radius=10 * u.arcmin)
table
table.keys()
Explanation: IRSA query of 2MASS
IRSA is the NASA/IPAC Infrared Science Archive.  IRSA typically has a infrared mission focus, such as 2MASS, Herschel, Planck, Spitzer, WISE, etc...
End of explanation
from astroquery.sdss import SDSS
Explanation: SDSS query
SDSS or the Sloan Digital Sky Survey is a massive survey of the sky which originally was optical imaging and spectroscopy.  It has been operating for over 20 years and the well beyond the original planned survey as it has expanded to near-infrared observations of stars and IFUs (integral field units) for observing galaxies and eventually the whole sky LVM (Local Volume Mapper). It is an immense data set to draw from.
End of explanation
ra, dec = 148.969687, 69.679383
co = coordinates.SkyCoord(ra=ra, dec=dec,unit=(u.deg, u.deg), frame='fk5')
xid = SDSS.query_region(co, radius=10 * u.arcmin)
# print the first 10 entries
print(xid[:10])
print(xid.keys())
print(xid['ra','dec'][:10]) # print the first 10 entries
Explanation: Coordinates of your favorite object, in this case M82 in degrees.
End of explanation
from astroquery.skyview import SkyView
SkyView.list_surveys()
pflist = SkyView.get_images(position='M82', survey=['SDSSr'],radius=10 * u.arcmin)
ext = 0
pf = pflist[0] # first element of the list, might need a loop if multiple images
m82_image = pf[ext].data
Explanation: Skyview
Skyview is a virtual observatory that is connected to several surveys and imaging data sets spanning multiple wavelengths. It can be a handy way to get quick access to multi-wavelength imaging data for your favorite object.
End of explanation
ax = plt.subplot()
ax.imshow(m82_image, cmap='gray_r', origin='lower', vmin=-10, vmax=20)
ax.set_xlabel('X (pixels)')
ax.set_ylabel('Y (pixels)')
Explanation: Plot image
End of explanation
from astropy.wcs import WCS
Explanation: Plot the in RA and Declination coordinates space, instead of pixels. The WCS (or World Coordinates System) in a FITS file has a set of keywords that define the coordinate system.  
Typical keywords are the following:
CRVAL1, CRVAL2 - RA and Declination (in degrees) coordinates of CRPIX1 and CRPIX2.
CRPIX1, CRPIX2 - Pixel positions of CRVAL1 and CRVAL2.
CDELT1, CDELT2 - pixel scale (or plate scale) in degrees/pixel 
(or CD matrix: CD1_1, CD1_2, CD2_1, CD2_2 which represent the pixel scale and rotation with respect to the position angle)
End of explanation
head = pf[ext].header
wcs = WCS(head)
ax = plt.subplot(projection=wcs)
ax.imshow(m82_image, cmap='gray_r', origin='lower', vmin=-10, vmax=20)
#ax.grid(color='white', ls='solid')
ax.set_xlabel('Right Ascension (J2000)')
ax.set_ylabel('Declination (J2000)')
#ax.scatter(xid['ra'],xid['dec'],marker="o",s=50,transform=ax.get_transform('fk5'),edgecolor='b', facecolor='none')
Explanation: Use the FITS header to define the WCS of the image.
End of explanation
<END_TASK>