Doing Astronomy with Python

One of the things that makes Python so powerful is that you can find a module for almostanything. In this article, I cover Astropy, which was originallydeveloped by the Space Telescope Science Institute for doing astronomy calculations likeimage processing and observatory calculations. Because this is a Pythonprogram, you can install it with eitherpip oreasy_install. Your Linuxdistribution already should have a package included. For example, in Debian-baseddistributions, you can install it with this:

sudo apt-get install python-astropy

There is also a separate package, python-astropy-doc, that containsextra documentation for Astropy. Because Astropy is a fairly large system,it is broken into smaller sub-packages. This should be familiar toanyone who has worked with packages like SciPy or NumPy before. So,simply using the following actually isn't very useful:

import astropy

You more likely will need to importindividual sub-packages with commands like this:

from astropy.io import fits

There are sub-packages to handle file IO, cosmological calculations andcoordinate systems, among many other topics. Let me provide a basicintroduction to some of the available functionality so you have an ideaof what you can do.

First, let's look at how to deal with data files. Thecommon file format used in astrophysics and astronomy is the FITS fileformat. The PyFITS Python package was written to read and write FITSfiles. This code is actually the same as the code in the sub-packageastropy.io.fits, so you can use it in the same way. You actuallycan even just drop Astropy in as a plugin with the following:

from astropy.io import fits as pyfits

This way, you can use existing file management code without having to makeany changes.

The first thing you need to do is open your data file with:

from astropy.io import fitshdulist = fits.open("My_File.fit")

This returns an object that behaves like a list. Each element of thereturned object maps to a Header-Data Unit (HDU) in the file. You canget more information on the file with this command:

hdulist.info()

Each of the individual elements has a header and data portion. Youcan access them to see details about the data you are about toprocess.

Along with all of the library functions, Astropy includes aseries of command-line utilities to work with FITS files. You can checkthe headers of a FITS file with thefitsheader utility. You can checkyour FITS file withfitscheck, and you even can find the differences betweentwo files withfitsdiff.

A common computational process in astronomy is imageprocessing. The convolution sub-package provides two categories ofconvolution operations: direct and FFT. You can do one-, two- and three-dimensionalconvolutions. The visualization sub-package handles more basic imageprocessing like normalization and stretching. You can combine multipletransformations very easily. The + operator is overloaded to applytransformations that are "added" together in series. So, acommand like this:

transform = SqrtStretch() + PercentileInterval(90.)

gives you a new function,transform, that combines the two separatetransformations in a single step. This sub-package also includes ascript,fits2bitmap, that can do conversions between different fileformats.

A second common computational task in astronomy is doing statistics based onobservations, and Astropy provides a sub-package called stats.Although the scipy.stats sub-package provides alot of functionality, some astronomy-focused functions are missing, so the astropy.stats sub-package fills in those missingfunctions.

Once you have your data loaded, you can use the modelingsub-package. You can do 1D and 2D modeling withastropy.modeling. Thisincludes curve-fitting functionality, where you can do linear andnonlinear function fitting. There are built-in functions to fit Gaussiancurves and polynomials. This fitting is handled with a least-squaresmethod. With version 1.0, you can build compound models by combiningexisting models with arithmetic operators.

When you are ready to start doing calculations, you will need to useconstants. Instead of remembering them or usingthem with potential typos, Astropy includes a complete list of allthe standard scientific constants that you will need when doing numericalwork. You can import the entire list with this:

from astropy import constants

If you need only a few of the constants, like maybe the speed of light,you can import them individually with this:

from astropy.constants import c

The really cool thing about these constants is that they are actually"Quantity" objects. This means you can do something like change theunits being used with a command like the following:

c.to('km/s')

Because it is so prevalent, you can use CGS units withc.cgs.

Thereare also two sub-packages to handle coordinate systems. Astronomicalcoordinate systems are handled by the coordinates sub-package, andworld coordinate systems are handled by the wcs sub-package. Inthe coordinates sub-package, the core object is theSkyCoordobject. There are methods of this object to handle conversions betweencoordinate systems or distances from one point to the origin within agiven coordinate system. The wcs sub-package allows for mapping datawithin a FITS file onto a "real-world" coordinate system in order to analyzethem correctly. This includes functionality to deal with complications,like projections onto the sphere of the sky.

You even can do cosmological calculations with Astropy. The cosmologysub-package actually includes functionality to model the evolutionof the entire cosmos based on a set of initial conditions thatyou set. Although you can set your own cosmology, severalbuilt-in cosmologies are available. These are based on the WMAP andPlanck satellite data.

Most functionality is built off a core FLRWobject. This object represents a homogeneous, isotropic cosmologydefined by the Friedmann-Lemaitre-Robertson-Walker metric from GeneralRelativity. However, this class can't be used directly. You need to subclass it first. There are several included in the cosmology sub-package,such as the FlatLambdaCDM object that includes dark energy. You can lookat various values, like the Hubble constant, at various points duringthe evolution of the cosmology. You also can include contributions to theenergy density from matter, dark energy and even photons and neutrinos.

Now that you've seen a bit of what you can do with astropy, ifastronomical calculations are on your radar, there is much moreavailable. Additionally, there is the concept of affiliated packages. Theseare packages that basically are built on top of the core functionalityprovided by Astropy. Although they aren't part of Astropy, they are beingbuilt up into a community-driven environment for doing astrophysics. Itdefinitely will be worth your while to take a look at the extended world ofavailable packages.