SunPy: Python for Solar Physics

Stuart J. Mumford

sunpy_banner

The Sun!

  • The Sun is our local star

  • Lots of data!

  • Process data

    – Compensate for instrumental effects.

    – Calibrate the data with some physical coords.

    – etc.

  • Analyse data

    – Do science!

    – Manipulate and analyse.

Data Flow

dataflow

Data Flow

dataflow-server

Data Flow

dataflow-processing

Data Flow

dataflow-analysis

SolarSoft

SolarSoft is a primarily IDL based set of programs and routines for the download, processing and analysis of solar data.

IDL - Interactive Data Language

  • IT'S HUGE! (9.5GB installed)

  • The source is freely available

  • No clear way to contribute

  • Badly Documented

  • No version control

  • Complex to update

  • Binary plugins written in C/C++/FORTRAN etc.

SunPy

SunPy was created out of frustration with IDL / SolarSoft.

It aims to provide a simple free and open-source alternative to SolarSoft.

Takes full advantage of the scientific python workflow:

  • Git / GitHub for version control / collaboration

  • Minimal data in the repo

  • (Reasonably) well documented

  • Makes use of Python packaging and PyPI

  • Some C-API code.

  • Python 2.6 & 2.7

  • Python 3.x support coming

  • Depends upon:

    • NumPy, SciPy & matplotlib

    • pandas

    • PyFITS / AstroPy.io.fits

    • suds & beautifulsoup4 (net module)

SunPy's Features

  • Integraded data retrevial and download: Virtual Solar Observatory(VSO) / Heliophysics Events Knowledgebase (HEK)

  • Core data types: Map, Lightcurve, Spectra

  • Physical coordinate transforms: World Coordinate System (WCS)

  • Transparent file handling and JPEG2000 support: I/O

  • Visualisation via matplotlib wrappers

SunPy Maps

SunPy Maps are the image data type in SunPy, maps are created using a factory called

In []:
import sunpy
sunpy.Map()

Maps can be created from many forms of input:

  • Preloaded tuples of (data, header) pairs

mymap = sunpy.Map((data, header))

  • data, header pairs, not in tuples

mymap = sunpy.Map(data, header)

  • File names

mymap = sunpy.Map('file1.fits')

  • All fits files in a directory by giving a directory

mymap = sunpy.Map('local_dir/sub_dir')

  • Some regex globs

mymap = sunpy.Map('eit_*.fits')

  • URLs

mymap = sunpy.Map(url_str)

  • Lists of any of the above

mymap = sunpy.Map(['file1.fits', 'file2.fits', 'file3.fits', 'directory1/'])

  • Any mixture of the above not in a list

mymap = sunpy.Map((data, header), data2, header2, 'file1.fits', url_str, 'eit_*.fits')

Instrument Specialisation

There are specialisations of Map for different instruments to compensate for variations in the meta data and other homogeneties between instruments. The following instruments currently have a SunPy Map subclass:

  • Hinode / XRT
  • Proba2 / SWAP
  • SDO / AIA
  • SDO / HMI
  • SoHO / EIT
  • SoHO / LASCO
  • SoHO / MDI
  • STEREO / EUVI
  • STEREO / COR
  • Yohkoh / SXT

The Map factory class automatically detects if there is a subclass of Map applicable to the data file it is processing.

Visualisation

As of SunPy 0.2 we made most of the plotting routines function in a similar manner to the matplotlib API:

In [5]:
%pylab inline
Welcome to pylab, a matplotlib-based Python environment [backend: module://IPython.kernel.zmq.pylab.backend_inline].
For more information, type 'help(pylab)'.
In [27]:
import matplotlib.pyplot as plt
import sunpy

aia = sunpy.Map(sunpy.AIA_171_IMAGE)

plt.figure(figsize=(8,8))

im = aia.plot()

plt.title("Hello!")
cbar = plt.colorbar(im)
plt.show()

LightCurves

Lightcurves are the 1D timeseries data type, currently lightcurves have specialisations for GOES, LYRA and EVE. Currently the lightcurve objects implement their own data grabbers using the .create() method:

In [14]:
import sunpy.lightcurve
goes = sunpy.lightcurve.GOESLightCurve.create('2012/06/01', '2012/06/05')
goes.plot()
Out[14]:
<matplotlib.axes.AxesSubplot at 0x2111eb90>

Spectra

The final core data type is the two spectral datatypes Spectrogram and Spectra. SunPy has an implementation for the e-Callisto network of instrumets:

In [26]:
import sunpy
from sunpy.spectra.sources.callisto import CallistoSpectrogram
spec = CallistoSpectrogram.read(sunpy.CALLISTO_IMAGE)
spec.plot()
Out[26]:

Data Download and Animation

SunPy can download data and also produce animations using matplotlib's animation routines:

In []:
import matplotlib.pyplot as plt 
import sunpy 
from sunpy.net import vso
	
client = vso.VSOClient()
	
result = client.query(
    vso.attrs.Time((2013, 6, 1, 8, 30, 0),
                   (2013, 6, 1, 8, 35, 0)),
    vso.attrs.Instrument('aia'),
    vso.attrs.Wave(304,304))

print "Number of records found: %d " % result.num_records() 
res = client.get(result, path="../AIA/{file}").wait()

cube = sunpy.Map(res, cube=True)
ani = cube.plot(controls=False)
In []:
Writer = animation.writers['ffmpeg']
writer = Writer(fps=5, metadata=dict(artist='Me'), bitrate=1800)

ani.save('my_animation.mp4', writer=writer)
plt.show()

Helioviewer.org

"Helioviewer.org is an open-source project for the visualization of solar and heliospheric data. The project is funded by ESA and NASA." It allows the download of a lot of pre-processed images as well as the exploration of a large archive of solar data on the web.

SunPy provides a interface to the Helioviewer API

In [2]:
from sunpy.net.helioviewer import HelioviewerClient
from IPython.core.display import Image 

hv = HelioviewerClient()
datasources = hv.get_data_sources()

res = hv.download_png('2099/01/01', 4.8, "[SDO,AIA,AIA,304,1,100]", x0=0, y0=0, width=512, height=512)
Image(res)
Out[2]:
In [6]:
res = hv.download_png('2099/01/01', 40,
                      "[SDO,AIA,AIA,304,1,100],[SDO,AIA,AIA,193,1,50],[SOHO,LASCO,C2,white-light,1,100],[SOHO,LASCO,C3,white-light,1,100]",
                      x0=0, y0=0, width=512, height=512)
Image(res)
Out[6]: