Skip to content

Commit

Permalink
Merge pull request #7526 from Deus1704/solar_rotate_fix
Browse files Browse the repository at this point in the history
fixes the `solar_rotate_coordinate()`
  • Loading branch information
nabobalis committed Jun 12, 2024
2 parents 767cd79 + 0aca03f commit 558e55e
Show file tree
Hide file tree
Showing 3 changed files with 18 additions and 2 deletions.
2 changes: 2 additions & 0 deletions changelog/7526.bugfix.rst
Original file line number Diff line number Diff line change
@@ -0,0 +1,2 @@
Added a check in `sunpy.physics.differential_rotation.solar_rotate_coordinate` to ensure the input frame has an "observer" attribute before replicating frame
attributes, preventing potential issues with frames lacking this attribute.
6 changes: 5 additions & 1 deletion sunpy/physics/differential_rotation.py
Original file line number Diff line number Diff line change
Expand Up @@ -255,8 +255,12 @@ def solar_rotate_coordinate(coordinate, observer=None, time=None, **diff_rot_kwa
# Calculate where the rotated coordinate appears as seen by new observer
# for the coordinate system of the input coordinate. The translational
# motion of the Sun will be ignored for the transformation.
frame_newobs = coordinate.frame.replicate_without_data(observer=new_observer,

if "observer" in coordinate.frame.frame_attributes.keys():
frame_newobs = coordinate.frame.replicate_without_data(observer=new_observer,
obstime=new_observer.obstime)
else:
frame_newobs = coordinate.frame.replicate_without_data(obstime=new_observer.obstime)
with transform_with_sun_center():
return heliographic_rotated.transform_to(frame_newobs)

Expand Down
12 changes: 11 additions & 1 deletion sunpy/physics/tests/test_differential_rotation.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@
import astropy.units as u
from astropy.coordinates import SkyCoord
from astropy.tests.helper import assert_quantity_allclose
from astropy.time import TimeDelta
from astropy.time import Time, TimeDelta

import sunpy.map
from sunpy.coordinates import frames, transform_with_sun_center
Expand Down Expand Up @@ -132,6 +132,16 @@ def test_solar_rotate_coordinate():
# Test that the SkyCoordinate is Helioprojective
assert isinstance(d.frame, frames.Helioprojective)

# Test that the function works correctly with a HGS coordinate.
earth_coord = get_earth(Time("2022-03-30"))
coord_hpc = SkyCoord(100*u.arcsec, 100*u.arcsec, frame=frames.Helioprojective(observer=earth_coord))

coord_hgs = coord_hpc.transform_to(frames.HeliographicStonyhurst)
with pytest.warns(UserWarning, match="Using 'time' assumes an Earth-based observer"):
rotated_coord_hgs = solar_rotate_coordinate(coord_hgs, time=Time("2022-03-31"))

assert isinstance(rotated_coord_hgs.frame, frames.HeliographicStonyhurst)


def test_consistency_with_rotatedsunframe():
old_observer = frames.HeliographicStonyhurst(10*u.deg, 20*u.deg, 1*u.AU, obstime='2001-01-01')
Expand Down

0 comments on commit 558e55e

Please sign in to comment.