# 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.map
import sunpy.coordinates
import sunpy.coordinates.wcs_utils
from sunpy.net 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 = Fido.search(wave, aia | stereo)
```

```files = Fido.fetch(res)
print(files)
```

Out:

```Files Downloaded:   0%|          | 0/2 [00:00<?, ?file/s]
['/home/docs/sunpy/data/aia_lev1_304a_2011_01_01t00_00_08_12z_image_lev1.fits', '/home/docs/sunpy/data/20110101_000615_n4eub.fts']
```

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,
frame=m.coordinate_frame))
for m in sunpy.map.Map(files)}
```

Out:

```/home/docs/checkouts/readthedocs.org/user_builds/sunpy/conda/stable/lib/python3.7/site-packages/astropy/io/fits/card.py:1003: VerifyWarning: Card is too long, comment will be truncated.
VerifyWarning)
keyword looks very much like CROTAn but isn't.
colsel=colsel)
keyword looks very much like CROTAn but isn't.
colsel=colsel)
keyword looks very much like CROTAn but isn't.
colsel=colsel)
keyword looks very much like CROTAn but isn't.
colsel=colsel)
```

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)
m.plot(axes=ax)
``` Out:

```/home/docs/checkouts/readthedocs.org/user_builds/sunpy/conda/stable/lib/python3.7/site-packages/astropy/wcs/wcs.py:466: FITSFixedWarning: CROTA = -4.13136390199
keyword looks very much like CROTAn but isn't.
colsel=colsel)
```

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,
frame=maps['AIA'].coordinate_frame)
aia_top_right = SkyCoord(aia_bottom_left.Tx + aia_width,
aia_bottom_left.Ty + aia_height,
frame=maps['AIA'].coordinate_frame)
```

Plot a rectangle around the region we want to crop

```m = maps['AIA']
fig = plt.figure()
m.plot(axes=ax)
m.draw_rectangle(aia_bottom_left, aia_width, aia_height)
``` Out:

```[<matplotlib.patches.Rectangle object at 0x7f142a02ea10>]
```

Create a submap of this area

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

```<matplotlib.image.AxesImage object at 0x7f142a038190>
```

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)

print(hpc_aia)
```

Out:

```<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)
(-0.01598098, -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)
print(hpc_B)
```

Out:

```/home/docs/checkouts/readthedocs.org/user_builds/sunpy/conda/stable/lib/python3.7/site-packages/astropy/wcs/wcs.py:466: FITSFixedWarning: CROTA = -4.13136390199
keyword looks very much like CROTAn but isn't.
colsel=colsel)
<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, m)
[[(422.14140721, -387.56099958, 1.57805355e+11)],
[(642.54613611, -376.44172931, 1.57947290e+11)],
[(513.08620817, -161.10016561, 1.57784422e+11)],
[(709.38834261, -149.2088207 , 1.57926317e+11)]]>
```

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)
m.plot(axes=ax)

# coord is the top-right corner coord is the bottom-left corner.
w = (coord.Tx - coord.Tx)
h = (coord.Ty - coord.Ty)
m.draw_rectangle(coord, w, h,
transform=ax.get_transform('world'))
``` Out:

```/home/docs/checkouts/readthedocs.org/user_builds/sunpy/conda/stable/lib/python3.7/site-packages/astropy/wcs/wcs.py:466: FITSFixedWarning: CROTA = -4.13136390199
keyword looks very much like CROTAn but isn't.
colsel=colsel)
keyword looks very much like CROTAn but isn't.
colsel=colsel)
```

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

```subeuvi = maps['EUVI'].submap(hpc_B, hpc_B)
fig = plt.figure()
plt.subplot(projection=subeuvi)
subeuvi.plot()
``` Out:

```/home/docs/checkouts/readthedocs.org/user_builds/sunpy/conda/stable/lib/python3.7/site-packages/astropy/wcs/wcs.py:466: FITSFixedWarning: CROTA = -4.13136390199
keyword looks very much like CROTAn but isn't.
colsel=colsel)
keyword looks very much like CROTAn but isn't.
colsel=colsel)
keyword looks very much like CROTAn but isn't.
colsel=colsel)
keyword looks very much like CROTAn but isn't.
colsel=colsel)
keyword looks very much like CROTAn but isn't.
colsel=colsel)

<matplotlib.image.AxesImage object at 0x7f1403ac0610>
```

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)
m.plot(axes=ax)
plt.show()
``` Out:

```/home/docs/checkouts/readthedocs.org/user_builds/sunpy/conda/stable/lib/python3.7/site-packages/astropy/wcs/wcs.py:466: FITSFixedWarning: CROTA = -4.13136390199
keyword looks very much like CROTAn but isn't.
colsel=colsel)
```

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

Gallery generated by Sphinx-Gallery