Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Implemented the vector transformation from one frame to another #7452

Open
wants to merge 15 commits into
base: main
Choose a base branch
from
Open
1 change: 1 addition & 0 deletions changelog/7452.feature.rst
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
Implemented the `~sunpy.coordinates.spice.transform_vector_field` which takes a vector field in one frame and transforms it to another frame.
77 changes: 76 additions & 1 deletion sunpy/coordinates/spice.py
Original file line number Diff line number Diff line change
Expand Up @@ -77,7 +77,7 @@
from sunpy.time.time import _variables_for_parse_time_docstring
from sunpy.util.decorators import add_common_docstring

__all__ = ['SpiceBaseCoordinateFrame', 'get_body', 'get_fov', 'initialize', 'install_frame']
__all__ = ['SpiceBaseCoordinateFrame', 'get_body', 'get_fov', 'initialize', 'install_frame', 'transform_vector_field']


# Note that this epoch is very slightly different from the typical definition of J2000.0 (in TT)
Expand Down Expand Up @@ -467,3 +467,78 @@ def get_fov(instrument, time, *, resolution=100):
obstime=obstime,
representation_type='unitspherical')
return fov

def transform_vector_field(components, source_frame, target_frame, from_time, to_time=None):
"""
Transforms a vector field (``components``) from ``source_frame`` to ``target_frame`` at a specified
time using the SPICE toolkit, which provides accurate transformations between
different astronomical frames.

Parameters
----------
components : `~astropy.units.Quantity`
A 3-element Quantity representing the vector field in the source frame.

source_frame : `str`, `int`, `~astropy.coordinates.BaseCoordinateFrame`
The frame of the input vector field, specified as a frame name, SPICE ID, or Astropy BaseCoordinateFrame.

target_frame : `str`, `int`, `~astropy.coordinates.BaseCoordinateFrame`
The frame to which the vector field will be transformed, specified in the same formats as ``source_frame``.

from_time : `str`, `~astropy.time.Time`, `~datetime.datetime`, `~datetime.date`, `~numpy.datetime64`, `~pandas.Series`, `~pandas.DatetimeIndex`, `~pandas.DataFrame`
The time at which the vector field is defined in the source frame, in any format accepted by :func:`sunpy.time.parse_time`.

to_time : `str`, `~astropy.time.Time`, `~datetime.datetime`, `~datetime.date`, `~numpy.datetime64`, `~pandas.Series`, `~pandas.DatetimeIndex`, `~pandas.DataFrame`, optional
The time at which the vector field should be defined in the target frame. Defaults to ``from_time``, resulting in a spatial-only transformation.

Returns
-------
`~astropy.units.Quantity`
The components of the vector field in the target frame. This is a 3-element Quantity with the same units
as the input ``components``.

Examples
--------
>>> from astropy import units as u
>>> from sunpy.coordinates.spice import transform_vector_field

>>> vec_components = [1, 0, 0] * u.T

>>> source_frame = "J2000"
>>> target_frame = "Galactic"
>>> from_time = '2001-01-01T00:00:00'

>>> transformed_vector = transform_vector_field(vec_components, source_frame, target_frame, from_time)

* The transform_vector_field function will first determine the transformation matrix that converts vectors from the source frame to the target frame at the specified time.
* This transformation matrix is a 3x3 matrix that represents the rotation between the two frames.

* The function then multiplies this transformation matrix with the input vector.
* This is a matrix-vector multiplication, which results in a new vector in the target frame.

>>> transformed_vector
<Quantity [-0.05487554, 0.49410945, -0.86766614] T>

* transformed_vector is now the input vector expressed in the target frame.

"""
# Convert sunpy frames to SPICE strings
source_frame_spice = source_frame.name.upper() if hasattr(source_frame, 'name') else str(source_frame).upper()
target_frame_spice = target_frame.name.upper() if hasattr(target_frame, 'name') else str(target_frame).upper()

from_time = parse_time(from_time)
to_time = parse_time(to_time) if to_time else from_time

from_time_et = _convert_to_et(from_time)
to_time_et = _convert_to_et(to_time)

# Compute state transformation matrix using spiceypy
state_matrix = spiceypy.sxform(source_frame_spice, target_frame_spice, to_time_et or from_time_et)

# The state transformation matrix is a 6x6 matrix. The upper left 3x3 block is the rotation matrix.
rotation_matrix = state_matrix[:3, :3]

# Apply rotation matrix to vector field components
transformed_components = rotation_matrix @ components

return transformed_components
14 changes: 14 additions & 0 deletions sunpy/coordinates/tests/test_spice.py
Original file line number Diff line number Diff line change
Expand Up @@ -216,3 +216,17 @@
sunpy_coord = hpc.geocentricsolarecliptic
assert_quantity_allclose(sunpy_coord.lon, spice_coord.lon, rtol=1e-6)
assert_quantity_allclose(sunpy_coord.lat, spice_coord.lat, rtol=1e-6)

def test_transform_vector_field(spice_test):
components = np.array([1, 0, 0]) * u.T
source_frame = "J2000"
target_frame = "Galactic"
from_time = '2001-01-01T00:00:00'

Check warning on line 224 in sunpy/coordinates/tests/test_spice.py

View check run for this annotation

Codecov / codecov/patch

sunpy/coordinates/tests/test_spice.py#L221-L224

Added lines #L221 - L224 were not covered by tests

result = spice.transform_vector_field(components, source_frame, target_frame, from_time)

Check warning on line 226 in sunpy/coordinates/tests/test_spice.py

View check run for this annotation

Codecov / codecov/patch

sunpy/coordinates/tests/test_spice.py#L226

Added line #L226 was not covered by tests

assert isinstance(result, u.Quantity)

Check warning on line 228 in sunpy/coordinates/tests/test_spice.py

View check run for this annotation

Codecov / codecov/patch

sunpy/coordinates/tests/test_spice.py#L228

Added line #L228 was not covered by tests

assert result.shape == components.shape
Deus1704 marked this conversation as resolved.
Show resolved Hide resolved
expected_result = np.array([-0.054875539, 0.49410945, -0.86766614])* u.T
assert_quantity_allclose(result.value, expected_result.value)

Check warning on line 232 in sunpy/coordinates/tests/test_spice.py

View check run for this annotation

Codecov / codecov/patch

sunpy/coordinates/tests/test_spice.py#L230-L232

Added lines #L230 - L232 were not covered by tests