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

Add a function to calculate the solar eclipse amount for an observer #7142

Merged
merged 7 commits into from Nov 1, 2023

Conversation

ayshih
Copy link
Member

@ayshih ayshih commented Aug 11, 2023

Go go solar eclipse!

Checking with other resources indicates that we can use our coordinates framework to calculate the start/end of total solar eclipses to within a few seconds.

To do:

  • Verify accuracy
  • Unit tests
  • Gallery example
  • Changelog entry
  • Add a docstring note recommending the use of a JPL ephemeris
  • Add an option to get_body_heliographic_stonyhurst() to silence the info message for the light-travel-time correction?

@ayshih ayshih added coordinates Affects the coordinates submodule No Backport A PR that isn't to be backported to any release branch. (To be used as a flag to other maintainers) labels Aug 11, 2023
@ayshih
Copy link
Member Author

ayshih commented Aug 11, 2023

import matplotlib.pyplot as plt
import numpy as np
from matplotlib.dates import DateFormatter

import astropy.units as u
from astropy.coordinates import EarthLocation
from astropy.time import Time

from sunpy.coordinates import sun

obstime = Time('2024-04-08 18:43') + np.arange(-100, 100) * u.min
observer = EarthLocation.from_geodetic(-96.808891*u.deg, 32.779167*u.deg).get_itrs(obstime)
amount = sun.eclipse_amount(observer)

plt.figure()
plt.plot(obstime.datetime64, amount)
plt.gca().xaxis.set_major_formatter(DateFormatter('%H:%M'))
plt.title('From Dallas on 2024 April 8')
plt.ylabel('Eclipse percentage')
plt.xlabel('UTC')
plt.grid()
plt.show()

Figure_1

@wtbarnes
Copy link
Member

wtbarnes commented Aug 11, 2023

Want to make this an example gallery entry? 😁

@ayshih
Copy link
Member Author

ayshih commented Aug 11, 2023

It's already on the to-do list. =). Appropriately spiffed up, of course

@ayshih
Copy link
Member Author

ayshih commented Aug 11, 2023

Accuracy looks excellent once I remember to use the JPL ephemeris rather than Astropy's built-in ephemeris

@wtbarnes wtbarnes added this to the 5.1.0 milestone Oct 4, 2023
# the array of observation times.

location = EarthLocation.from_geodetic(-96.808891*u.deg, 32.779167*u.deg)
observer = location.get_itrs(obstime)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Somewhat off on a tangent, but can we pass EarthLocation as an observer? i.e. skip this line. If not why not? 😆

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Because EarthLocation doesn't include obstime?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Drives the PR car off down another street

Couldn't we still do the conversion though? Wont most transforms pull the obstime from the frame if it's missing on the observer?

Copy link
Member Author

@ayshih ayshih Oct 17, 2023

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think you're not looking at how the observer variable is later used. Despite its name, it's not being used as the value for an observer frame attribute. The function eclipse_amount() takes a single argument, which is the coordinate of the observer including the obstime. If the function were to accept EarthLocation for the observer location, the function would need to additionally require the obstime to be supplied as a second argument.

.. minigallery:: sunpy.coordinates.sun.eclipse_amount
"""
# The radius of the Moon to use
k = 0.272281 if minimum else 0.2725076
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should we not put this in constants somewhere rather than in here?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

And also define it better. It's not a radius but a fraction of Earth's radius.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It's of course a dimensionless fraction, but it's also the value of the lunar radius in units of Earth radii.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't know where to stick such constants, since we don't really have a home for non-Sun constants. I've tacked on a TODO comment, and someone in the future can figure it out.

@ehsteve
Copy link
Member

ehsteve commented Oct 13, 2023

I tested this code and it worked great. This is a cool example and this capability could be used for an update to NumFocus. To make it more relevant, how about choosing a point where both 2023 and 2024 eclipses pass and show both the total eclipse and annular eclipse curves on the same plot?

@ayshih ayshih force-pushed the eclipse branch 3 times, most recently from 831543e to acdc703 Compare October 17, 2023 12:08
@ayshih
Copy link
Member Author

ayshih commented Oct 17, 2023

To make it more relevant, how about choosing a point where both 2023 and 2024 eclipses pass and show both the total eclipse and annular eclipse curves on the same plot?

Okay, I've picked a location on the outskirts of San Antonio that is on both eclipse paths.

@ayshih ayshih marked this pull request as ready for review October 17, 2023 13:13
@ayshih ayshih requested review from a team as code owners October 17, 2023 13:13
@ayshih ayshih added the Needs Review Needs reviews before merge label Oct 17, 2023
Copy link
Member

@dstansby dstansby left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Looks good, I've left a few comments 👍

sunpy/coordinates/sun.py Outdated Show resolved Hide resolved
sunpy/coordinates/sun.py Outdated Show resolved Hide resolved
sunpy/coordinates/sun.py Outdated Show resolved Hide resolved
sunpy/coordinates/sun.py Outdated Show resolved Hide resolved
sunpy/coordinates/sun.py Outdated Show resolved Hide resolved
@dstansby dstansby added the Whats New? Needs a section added to the current Whats New? page. label Oct 18, 2023
@ayshih ayshih force-pushed the eclipse branch 2 times, most recently from dac22ab to 57ce27b Compare October 19, 2023 00:39
sunpy/coordinates/ephemeris.py Show resolved Hide resolved
sunpy/coordinates/sun.py Outdated Show resolved Hide resolved
Copy link
Member

@dstansby dstansby left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should fix the test...

sunpy/coordinates/tests/test_sun.py Outdated Show resolved Hide resolved
@ayshih ayshih force-pushed the eclipse branch 2 times, most recently from c4efd83 to 5ac2d71 Compare October 25, 2023 20:17
Copy link
Member

@Cadair Cadair left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

A couple of very minor tweaks.

sunpy/coordinates/ephemeris.py Show resolved Hide resolved
sunpy/coordinates/sun.py Outdated Show resolved Hide resolved
@ayshih ayshih force-pushed the eclipse branch 2 times, most recently from 846b886 to ac9e625 Compare October 27, 2023 06:05
@nabobalis nabobalis removed the Needs Review Needs reviews before merge label Nov 1, 2023
@nabobalis nabobalis merged commit f57452c into sunpy:main Nov 1, 2023
25 of 26 checks passed
@ayshih ayshih deleted the eclipse branch November 3, 2023 01:26
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
coordinates Affects the coordinates submodule No Backport A PR that isn't to be backported to any release branch. (To be used as a flag to other maintainers) Whats New? Needs a section added to the current Whats New? page.
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

6 participants