differential_rotation#
- sunpy.sun.models.differential_rotation(duration: Unit('s'), latitude: Unit('deg'), *, model='howard', frame_time='sidereal')[source]#
Computes the change in longitude over a duration for a given latitude.
Since the Sun is not a rigid body, different heliographic latitudes rotate with different periods. This is known as solar differential rotation.
- Parameters:
duration (
Quantity
) – Amount of time to rotate over.latitude (
Quantity
) – Heliographic latitude.model (
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 fromsidereal_rotation_rate
.frame_time (
str
) – If'sidereal'
, returns the change in longitude as referenced to distant stars. If'synodic'
, returns the apparent change in longitude as observed by the average orbital motion of Earth, which results in a slower rotation rate. Defaults to'sidereal'
.
- Returns:
Quantity
– The change in longitude
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
Solar surface velocity fields determined from small magnetic features (Howard et al. 1990)
A comparison of differential rotation measurements (Beck 2000, includes Snodgrass values)
Examples
Simple Differential RotationDefault rotation calculation over two days at 30 degrees latitude:
>>> import numpy as np >>> import astropy.units as u >>> from sunpy.sun.models import differential_rotation >>> differential_rotation(2 * u.day, 30 * u.deg) <Longitude 27.36432679 deg>
Default rotation over two days for a number of latitudes:
>>> differential_rotation(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 model ‘allen’:
>>> differential_rotation(2 * u.day, np.linspace(-70, 70, 20) * u.deg, model='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>