diff_rot#

sunpy.physics.differential_rotation.diff_rot(duration: Unit('s'), latitude: Unit('deg'), rot_type='howard', frame_time='sidereal')[source]#

This function computes the change in longitude over days in degrees.

Parameters:
  • duration (Quantity) – Number of seconds to rotate over.

  • latitude (Quantity) – heliographic coordinate latitude in Degrees.

  • rot_type (str) – The differential rotation model to use.

    One of:

    howard : Use values from Howard et al. (1990)
    snodgrass : Use values from Snodgrass et. al. (1983)
    allen : Use values from Allen’s Astrophysical Quantities, and simpler equation.
    rigid : Use values from sidereal_rotation_rate.
  • frame_time (str) – One of : 'sidereal' or 'synodic'. Choose ‘type of day’ time reference frame.

Returns:

longitude_delta (Quantity) – The change in longitude over days (units=degrees)

Notes

The rotation rate at a heliographic latitude \(\theta\) is given by

\[A + B \sin^{2} \left (\theta \right ) + C \sin^{4} \left ( \theta \right )\]

where \(A, B, C\) are constants that depend on the model:

Model

A

B

C

Unit

howard

2.894

-0.428

-0.370

microrad/s

snodgrass

2.851

-0.343

-0.474

microrad/s

allen

14.44

-3.0

0

deg/day

rigid

14.1844

0

0

deg/day

1 microrad/s is approximately 4.95 deg/day.

References

Examples

Simple Differential Rotation

Simple Differential Rotation

Comparing differential-rotation models

Comparing differential-rotation models

Default rotation calculation over two days at 30 degrees latitude:

>>> import numpy as np
>>> import astropy.units as u
>>> from sunpy.physics.differential_rotation import diff_rot
>>> diff_rot(2 * u.day, 30 * u.deg)
<Longitude 27.36432679 deg>

Default rotation over two days for a number of latitudes:

>>> diff_rot(2 * u.day, np.linspace(-70, 70, 20) * u.deg)
<Longitude [22.05449682, 23.03214991, 24.12033958, 25.210281  ,
            26.21032832, 27.05716463, 27.71932645, 28.19299667,
            28.49196765, 28.63509765, 28.63509765, 28.49196765,
            28.19299667, 27.71932645, 27.05716463, 26.21032832,
            25.210281  , 24.12033958, 23.03214991, 22.05449682] deg>

With rotation type ‘allen’:

>>> diff_rot(2 * u.day, np.linspace(-70, 70, 20) * u.deg, 'allen')
<Longitude [23.58186667, 24.14800185, 24.82808733, 25.57737945,
        26.34658134, 27.08508627, 27.74430709, 28.28087284,
        28.6594822 , 28.85522599, 28.85522599, 28.6594822 ,
        28.28087284, 27.74430709, 27.08508627, 26.34658134,
        25.57737945, 24.82808733, 24.14800185, 23.58186667] deg>