1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
|
"""
============================================
Coordinates computations using SPICE kernels
============================================
How to use SPICE kernels provided by space missions to perform coordinates
computations.
The `SPICE <https://naif.jpl.nasa.gov/naif/>`__ observation geometry information
system is being increasingly used by space missions to describe the locations of
spacecraft and the time-varying orientations of reference frames. Here are two
examples of mission SPICE kernels:
* `Solar Orbiter <https://www.cosmos.esa.int/web/spice/solar_orbiter>`__
* `Parker Solar Probe <https://psp-gateway.jhuapl.edu/website/Ancillary/SclkFiles>`__
The `sunpy.coordinates.spice` module enables the use of the
`~astropy.coordinates.SkyCoord` API to perform SPICE computations such as the
location of bodies or the transformation of a vector from one coordinate frame
to another coordinate frame. Although SPICE kernels can define coordinate
frames that are very similar to the frames that `sunpy.coordinates` already
provides, there will very likely be slight differences. Using
`sunpy.coordinates.spice` will ensure that the definitions are exactly what the
mission specifies and that the results are identical to other implementations
of SPICE (e.g., CSPICE or Icy).
.. note::
This example requires the optional dependency `~spiceypy.spiceypy` to be
installed.
"""
import matplotlib.pyplot as plt
import numpy as np
from matplotlib.dates import DateFormatter
import astropy.units as u
from astropy.coordinates import SkyCoord
from sunpy.coordinates import spice
from sunpy.data import cache
from sunpy.time import parse_time
###############################################################################
# Download a small subset (~30 MB) of the Solar Orbiter SPICE kernel set that
# corresponds to about a day of the mission.
obstime = parse_time('2022-10-12') + np.arange(720) * u.min
kernel_urls = [
"ck/solo_ANC_soc-sc-fof-ck_20180930-21000101_V03.bc",
"ck/solo_ANC_soc-stix-ck_20180930-21000101_V03.bc",
"ck/solo_ANC_soc-flown-att_20221011T142135-20221012T141817_V01.bc",
"fk/solo_ANC_soc-sc-fk_V09.tf",
"fk/solo_ANC_soc-sci-fk_V08.tf",
"ik/solo_ANC_soc-stix-ik_V02.ti",
"lsk/naif0012.tls",
"pck/pck00010.tpc",
"sclk/solo_ANC_soc-sclk_20231015_V01.tsc",
"spk/de421.bsp",
"spk/solo_ANC_soc-orbit-stp_20200210-20301120_280_V1_00288_V01.bsp",
]
kernel_urls = [f"http://spiftp.esac.esa.int/data/SPICE/SOLAR-ORBITER/kernels/{url}"
for url in kernel_urls]
kernel_files = [cache.download(url) for url in kernel_urls]
###############################################################################
# Initialize `sunpy.coordinates.spice` with these kernels, which will create
# classes for `~astropy.coordinates.SkyCoord` to use.
spice.initialize(kernel_files)
###############################################################################
# The above call automatically installs all SPICE frames defined in the
# kernels, but you may also want to use one of the built-in SPICE frames (e.g.,
# inertial frames or body-fixed frames). Here, we manually install the
# 'IAU_SUN' built-in SPICE frame for potential later use.
spice.install_frame('IAU_SUN')
###############################################################################
# We can request the location of the spacecraft in any SPICE frame. Here, we
# request it in 'SOLO_HEEQ', which is Stonyhurst heliographic coordinates.
spacecraft = spice.get_body('Solar Orbiter', obstime, spice_frame='SOLO_HEEQ')
print(spacecraft[:4])
###############################################################################
# Plot the radial distance from the Sun over the time range.
fig = plt.figure()
ax = fig.add_subplot()
ax.plot(obstime.datetime64, spacecraft.distance.to('AU'))
ax.xaxis.set_major_formatter(DateFormatter('%H:%M'))
ax.set_xlabel('2022 October 12 (UTC)')
ax.set_ylabel('Radial distance (AU)')
ax.set_title('Solar Orbiter distance from Sun center')
###############################################################################
# We can then transform the coordinate to a different SPICE frame. When
# specifying the frame for `~astropy.coordinates.SkyCoord`, SPICE frame names
# should be prepended with ``'spice_'``. Here, we transform it to 'SOLO_GAE'.
spacecraft_gae = spacecraft.transform_to("spice_SOLO_GAE")
print(spacecraft_gae[:4])
###############################################################################
# Plot the radial distance from the Earth over the time range.
fig = plt.figure()
ax = fig.add_subplot()
ax.plot(obstime.datetime64, spacecraft_gae.distance.to('AU'))
ax.xaxis.set_major_formatter(DateFormatter('%H:%M'))
ax.set_xlabel('2022 October 12 (UTC)')
ax.set_ylabel('Radial distance (AU)')
ax.set_title('Solar Orbiter distance from Earth center')
###############################################################################
# We can also leverage the Solar Orbiter SPICE kernels to look at instrument
# pointing. Let's define a coordinate that points directly along the line of
# sight of the STIX instrument. For the 'SOLO_STIX_ILS' frame, 0 degrees
# longitude is in the anti-Sun direction, while 180 degrees longitude is in the
# Sun direction.
stix_ils = SkyCoord(np.repeat(0*u.deg, len(obstime)),
np.repeat(0*u.deg, len(obstime)),
frame='spice_SOLO_STIX_ILS', obstime=obstime)
print(stix_ils[:4])
###############################################################################
# We can transform that line of sight to the SPICE frame 'SOLO_SUN_RTN', which
# is similar to helioprojective coordinates as observed from Solar Orbiter,
# except that the disk center is at 180 degrees longitude instead of 0 degrees
# longitude. Given how the line-of-sight coordinate is defined above, the
# latitude and longitude values of the resulting coordinate are the pitch and
# yaw offsets from disk center, respectively.
stix_ils_rtn = stix_ils.transform_to("spice_SOLO_SUN_RTN")
print(stix_ils_rtn[:4])
###############################################################################
# Plot the pitch/yaw offsets over the time range.
fig, ax = plt.subplots(2, 1)
ax[0].plot(obstime.datetime64, stix_ils_rtn.lat.to('arcsec'))
ax[0].xaxis.set_major_formatter(DateFormatter('%H:%M'))
ax[0].set_xlabel('2022 October 12 (UTC)')
ax[0].set_ylabel('Pitch offset (arcsec)')
ax[1].plot(obstime.datetime64, stix_ils_rtn.lon.to('arcsec'))
ax[1].xaxis.set_major_formatter(DateFormatter('%H:%M'))
ax[1].set_xlabel('2022 October 12 (UTC)')
ax[1].set_ylabel('Yaw offset (arcsec)')
ax[0].set_title('Pointing offset of STIX from disk center')
plt.show()
###############################################################################
# Finally, we can query the instrument field of view (FOV) via SPICE, which
# will be in the 'SOLO_STIX_ILS' frame. This call returns the corners of the
# rectangular FOV of the STIX instrument, and you can see they are centered
# around 180 degrees longitude, which is the direction of the Sun in this
# frame.
stix_fov = spice.get_fov('SOLO_STIX', obstime[0])
print(stix_fov)
###############################################################################
# More usefully, every coordinate in a SPICE frame has a
# :meth:`~sunpy.coordinates.spice.SpiceBaseCoordinateFrame.to_helioprojective`
# method that converts the coordinate to `~sunpy.coordinates.Helioprojective`
# with the ``observer`` at the center of te SPICE frame. For the
# 'SOLO_STIX_ILS' frame, the center is Solar Orbiter, which is exactly what we
# want.
print(stix_fov.to_helioprojective())
# sphinx_gallery_thumbnail_number = 3
|
