AIA to STEREO coordinate conversion

How to convert a point of a source on an AIA image to a position on a STEREO image.

import matplotlib.pyplot as plt

import astropy.units as u
from astropy.coordinates import SkyCoord

import sunpy.coordinates
import sunpy.coordinates.wcs_utils
from import Fido, attrs as a

The first step is to download some data, we are going to get an image from early 2011 when the STEREO spacecraft were roughly 90 deg seperated from the Earth.

stereo = (a.vso.Source('STEREO_B') &
          a.Instrument('EUVI') &
          a.Time('2011-01-01', '2011-01-01T00:10:00'))

aia = (a.Instrument('AIA') &
       a.vso.Sample(24 * u.hour) &
       a.Time('2011-01-01', '2011-01-02'))

wave = a.Wavelength(30 * u.nm, 31 * u.nm)

res =, aia | stereo)

Download the files:

files = Fido.fetch(res)


['/home/docs/sunpy/data/secchi_l0_b_img_euvi_20110101_20110101_000615_n4eub.fts', '/home/docs/sunpy/data/aia_lev1_304a_2011_01_01t00_00_08_12z_image_lev1.fits']

Create a dictionary with the two maps, cropped down to full disk.

maps = {m.detector: m.submap(SkyCoord([-1100, 1100]*u.arcsec,
                                      [-1100, 1100]*u.arcsec,
        for m in}

Plot both maps

fig = plt.figure(figsize=(10, 4))
for i, m in enumerate(maps.values()):
    ax = fig.add_subplot(1, 2, i+1, projection=m)

We are now going to pick out a region around the south west corner:

aia_width = 200 * u.arcsec
aia_height = 250 * u.arcsec
aia_bottom_left = SkyCoord([[-800, -300]] * u.arcsec,
aia_top_right = SkyCoord(aia_bottom_left.Tx + aia_width,
                         aia_bottom_left.Ty + aia_height,

Plot a rectangle around the region we want to crop

m = maps['AIA']
fig = plt.figure()
ax = fig.add_subplot(111, projection=m)
m.draw_rectangle(aia_bottom_left, aia_width, aia_height)

Create a submap of this area

subaia = maps['AIA'].submap(aia_bottom_left, aia_top_right)
fig = plt.figure()

We now want to crop out this same area on the STEREO EUVI image. First, we create a SkyCoord object with the four corners of the box. When we create this object, we use Map.coordinate_frame so that the location parameters of SDO are correctly set.

corners = ([aia_bottom_left.Tx, aia_bottom_left.Ty],
           [aia_bottom_left.Tx + aia_width, aia_bottom_left.Ty],
           [aia_bottom_left.Tx, aia_bottom_left.Ty + aia_height],
           [aia_top_right.Tx, aia_top_right.Ty])

hpc_aia = SkyCoord(corners, frame=maps['AIA'].coordinate_frame)



<SkyCoord (Helioprojective: obstime=2011-01-01T00:00:08.120, rsun=696000000.0 m, observer=<HeliographicStonyhurst Coordinate (obstime=2011-01-01T00:00:08.120): (lon, lat, radius) in (deg, deg, m)
    (359.98401902, -2.97460829, 1.47100674e+11)>): (Tx, Ty) in arcsec
    [[(-800., -300.)],
     [(-600., -300.)],
     [(-800.,  -50.)],
     [(-600.,  -50.)]]>

We can now transform to from the AIA frame to the EUVI frame. This transformation first transforms to Heliographic Stonyhurst coordinates and then into the EUVI frame.

hpc_B = hpc_aia.transform_to(maps['EUVI'].coordinate_frame)


<SkyCoord (Helioprojective: obstime=2011-01-01T00:07:01.811, rsun=695700000.0 m, observer=<HeliographicStonyhurst Coordinate (obstime=2011-01-01T00:07:01.811): (lon, lat, radius) in (deg, deg, m)
    (-89.58098377, 6.76626853, 1.58345292e+11)>): (Tx, Ty, distance) in (arcsec, arcsec, km)
    [[(422.21047096, -387.55701629, 1.57805382e+08)],
     [(642.59938241, -376.43551334, 1.57947332e+08)],
     [(513.1552792 , -161.09530724, 1.57784455e+08)],
     [(709.44158877, -149.2019928 , 1.57926363e+08)]]>

Now we can plot this box on both the AIA and EUVI images:

fig = plt.figure(figsize=(10, 4))
for i, (m, coord) in enumerate(zip([maps['EUVI'], maps['AIA']],
                                   [hpc_B, hpc_aia])):
    ax = fig.add_subplot(1, 2, i+1, projection=m)

    # coord[3] is the top-right corner coord[0] is the bottom-left corner.
    w = (coord[3].Tx - coord[0].Tx)
    h = (coord[3].Ty - coord[0].Ty)
    m.draw_rectangle(coord[0], w, h,

We can now zoom in on the region in the EUVI image:

subeuvi = maps['EUVI'].submap(hpc_B[0], hpc_B[3])
fig = plt.figure()

Putting them together:

fig = plt.figure(figsize=(15, 5))
for i, m in enumerate((subeuvi, subaia)):
    ax = fig.add_subplot(1, 2, i+1, projection=m)

Total running time of the script: ( 0 minutes 12.684 seconds)

Gallery generated by Sphinx-Gallery