SunPy–The Community-developed, Free, and Open-source Solar Data Analysis Environment for Python

VAST Seminar Series / Virtual / 15 May 2024

Will Barnes

American University/NASA GSFC

on behalf of The SunPy Project

What is solar physics?

• Astronomy/astrophysics of the Sun
• Can study our star in greater detail than any other
• Impact on space weather on the earth
• A laboratory for astrophysics, fundamental plasma physics, and more!
• Primary focus is the solar atmosphere

Solar data are dyanmic

• Magnetic field holds the key to why the corona is heated to million degree temperatures
• Plasma trapped along field lines
• Twisted magnetic field releases energy into the plasma, thereby heating it
• We can start to understand this system by observing these areas of intense magnetic activity at different wavelengths as a function of time

…high resolution

• Collect many different kinds of data
• Spectra allow us to see more details: more precise temperatures, how things are moving, where things are forming
• We observe different structures at different wavelengths

…observed from many viewpoints

Challenges

Data…

• …are spread across many different providers
• …are often large and split across many files
• …can be very heterogeneous (may or may not follow FITS standard)
• …must be held together with the metadata
• …are expressed in complex coordinate systems

Software…

• …is typically developed by instrument teams on a best effort basis
• …development is not coordinated across instrument teams
• …has traditionally been written in IDL, though with a centralized “core” codebase (SSW)

sunpy: Solar Data Analysis in Python

• Began in March of 2011 at NASA GSFC
• Frustration with licensing fees, fragility of SSW
• Early attempts using GDL\(\to\)Python
• v0.1 released 9/2011
• v1.0 in 6/2019-paper currently has >200 citations
• Now at v5.1 with v6.0 imminent
• Some early skepticism due to heritage of IDL/SSW
• Now the default choice, especially for ECRs
• Open-source and openly developed–by the community, for the community
• Built on SciPy ecosystem, especially astropy

sunpy: Solar Data Analysis in Python

Searching for Data with sunpy

from sunpy.net import Fido, attrs as a
import astropy.units as u

1query = Fido.search(
  a.Time('2018-05-29 18:00', '2018-05-29 18:00:10'),
  a.Wavelength(171*u.angstrom),
  a.Instrument.aia
)
print(query)

1

Search remote repository for AIA 171 Å images between 18:00 and 18:10 on May 29 2018.

Results from 1 Provider:

1 Results from the VSOClient:
Source: http://vso.stanford.edu/cgi-bin/search
Total estimated size: 67.789 Mbyte

       Start Time               End Time        Source ... Extent Type   Size  
                                                       ...              Mibyte 
----------------------- ----------------------- ------ ... ----------- --------
2018-05-29 18:00:09.000 2018-05-29 18:00:10.000    SDO ...    FULLDISK 64.64844

2files = Fido.fetch(query, path='data/{instrument}')
print(files)

2

Download all FITS files corresponding to this search result to data/AIA.

['data/AIA/aia_lev1_171a_2018_05_29t18_00_09_35z_image_lev1.fits']

Loading Data with sunpy

import sunpy.map

m = sunpy.map.Map(files)
m.peek()

Inspecting Data with sunpy

print(m.detector)
print(m.wavelength)
print(m.date)
print(m.observer_coordinate)

AIA
171.0 Angstrom
2018-05-29T18:00:09.350
<SkyCoord (HeliographicStonyhurst: obstime=2018-05-29T18:00:09.350, rsun=696000.0 km): (lon, lat, radius) in (deg, deg, m)
    (0.00562552, -0.97774891, 1.51600359e+11)>

from astropy.coordinates import SkyCoord

corner = SkyCoord(Tx=-375*u.arcsec, Ty=0*u.arcsec,
                  frame=m.coordinate_frame)
print(m.world_to_pixel(corner))
m_cutout = m.submap(corner,
                    width=500*u.arcsec,
                    height=500*u.arcsec)
m_cutout.peek()

PixelPair(x=<Quantity 1425.73662106 pix>, y=<Quantity 2049.01122848 pix>)

Transforming Data with sunpy

from sunpy.coordinates import propagate_with_solar_surface
m_seq = sunpy.map.Map('data/sequence/*.fits', sequence=True)
fig = plt.figure(figsize=(16, 4), layout='constrained')
for i, m in enumerate(m_seq):
  ax = fig.add_subplot(1, len(m_seq), i+1, projection=m)
  m.plot(axes=ax)
  with propagate_with_solar_surface():
    blc = m_cutout.bottom_left_coord.transform_to(m.coordinate_frame)
    trc = m_cutout.top_right_coord.transform_to(m.coordinate_frame)
  m.draw_quadrangle(blc, top_right=trc)

Transforming Data with sunpy

with propagate_with_solar_surface():
  m_seq_aligned = sunpy.map.Map([m.reproject_to(m_cutout.wcs) for m in m_seq], sequence=True)
fig = plt.figure(figsize=(16, 4), layout='constrained')
for i, m in enumerate(m_seq_aligned):
  ax = fig.add_subplot(1,len(m_seq_aligned), i+1, projection=m)
  m.plot(axes=ax, cmap='sdoaia171', title=m_seq[i].date)

The SunPy Project

The SunPy Ecosystem

Functionality

Integration

Documentation

Testing

Duplication

Community

Development Status

Summary

• Solar data are large, complex, and heterogenous, but hold great value in combination
• sunpy is a community-developed Python package for solar data analysis
• The sunpy package provides the ability to search, load, and transform solar data
• The SunPy Project is composed of not just the code, but also the people and the community
• The SunPy Project maintains an interoperable ecosystem of tools for working with solar data

Resources

• Website: sunpy.org
• GitHub: github.com/sunpy
• Documentation: docs.sunpy.org
• Matrix Chat