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

Quicker versions of target_rise_time & target_set_time #423

Closed
iwellaway opened this issue Oct 16, 2019 · 9 comments
Closed

Quicker versions of target_rise_time & target_set_time #423

iwellaway opened this issue Oct 16, 2019 · 9 comments

Comments

@iwellaway
Copy link

@iwellaway iwellaway commented Oct 16, 2019

We are using the Astroplan Observer object extensively in the development of our robotic scheduler for the INT. We need to run a 10 year simulation so any way we can speed up the calculations would be a great help.

Is it possible to produce a faster, albeit probably slightly less accurate, version of the Observer object methods: target_rise_time, target_set_time and target_meridian_transit_time?

@bmorris3

This comment has been minimized.

Copy link
Contributor

@bmorris3 bmorris3 commented Oct 16, 2019

Thanks for the note @iwellaway. I've implemented a faster "trigonometric solver" in pull request #327, which is currently under review. Would you please leave comments on that pull request if you have any feedback on the implementation? I'd like to make sure that it will work well for folks like you!

@iwellaway

This comment has been minimized.

Copy link
Author

@iwellaway iwellaway commented Oct 16, 2019

@eteq

This comment has been minimized.

Copy link
Member

@eteq eteq commented Oct 16, 2019

@iwellaway - is there a specific test case you can provide that might let us do some profiling on your use case? Ideally a compact snippet that we can just run with astroplan, but if it's going to be hard to do that, a pointer to how to run a specific scheduler example with your code would also work.

I have a suspicion that there are areas of these functions that could be much better optimized at least for some use cases without going the "lower-accuracy" direction. As I said in #327 it's not really clear to me that the coordinates machinery itself is the "hot spot" in these functions, especially given that there have been several major performance enhancements since #327 was first created. So I'd like to investigate that a bit first before going in direction like #327.

@bmorris3

This comment has been minimized.

Copy link
Contributor

@bmorris3 bmorris3 commented Oct 17, 2019

To begin this discussion, I've timed the Observer.target_rise_time function with and without the trig approximation, see gist here, which shows a 2.4x speed-up (for a single target). The only machinery change between these functions is the coordinate engine swap.

@bmorris3

This comment has been minimized.

Copy link
Contributor

@bmorris3 bmorris3 commented Oct 17, 2019

For the sake of transparency in this process, I'm going to briefly outline what astroplan is currently doing to compute rise/set times, mostly in the Observer._calc_riseset function:

  1. Generate a grid of times with 5 min resolution over the 24 hours before or after the specified reference time [code] (fast!)
  2. Compute the alt/az of the target at each time on the 24 hour grid [code] (slower)
  3. Identify two grid-points where the target rises/sets above/below the user-defined horizon altitude, interpolate between the two to find the rise/set time [code] (fast!)

There are several imperfect ways to speed up step (2). They might include:

  • reducing the number of grid-points on which the altitude is computed in step (1-2)
  • identifying existing options within coordinate machinery that can be toggled to speed things up (i.e. if we set the pressure to 0*u.bar, does turning atmospheric refraction off actually speed things up?) (hacky)
  • speeding up the existing astropy coordinates submodule (difficult)
  • swapping the engine computation machinery (inconsistent)

After thinking about it a bit, perhaps it'd be worth exploring the first option in my list, since it's the most consistent approach. I can do a small study to see how precision in rise/set times varies as a function of the number of grid points – if it turns out we can reduce the number of grid points by a factor of a few for a rise/set time accuracy tradeoff of a few minutes, we may be able to similarly decrease the runtime.

@bmorris3

This comment has been minimized.

Copy link
Contributor

@bmorris3 bmorris3 commented Oct 17, 2019

I've created a branch of astroplan which has an extra keyword argument to the rise/set functions called N which exposes to the user the number of grid-points on which the altitudes are computed in steps (1-2). That way you can easily say: "compute rise/set times on a grid of 10 points per 24 hours or 1000 points per 24 hours," for coarse and fine precisions, respectively. By default, we run N=150 within astroplan right now, but we left that keyword argument hanging around rather than hard-coding N=150 because we originally intended for users to be able to change it!

Here's some code that probes how much precision you lose when you go from N=1000 to 10:

%matplotlib inline
import matplotlib.pyplot as plt
from astropy.time import Time
import astropy.units as u
import numpy as np

from astroplan import Observer, FixedTarget

target = FixedTarget.from_name('Sirius')
apo = Observer.at_site('APO')

time_grid = Time.now() + np.arange(0, 50) * u.day

accurate_riseset = apo.target_rise_time(time_grid, target, N=1000, grid_times_targets=True)
medium_riseset = apo.target_rise_time(time_grid, target, N=100, grid_times_targets=True)
coarse_riseset = apo.target_rise_time(time_grid, target, N=10, grid_times_targets=True)

plt.plot((medium_riseset - accurate_riseset).to(u.min).value)
plt.plot((coarse_riseset - accurate_riseset).to(u.min).value)
plt.ylim([-1, 4])
plt.ylabel('Computed - Expected [min]')
plt.xlabel('time [days]')

astroplan_precision

The blue line shows a comparison between N=1000 and N=100, showing that our typical uncertainties on rise/set times are much better than one minute, and the orange curve shows the rise/set time precisions for N=10, demonstrating that even if there are only ten grid points, we're still predicting rise/set times with precisions better than 5 minutes – which is the rough specified tolerance of folks like @iwellaway!

This is great news, because the runtime scales roughly linearly with the number of grid points, demonstrated below,

%timeit apo.target_rise_time(time_grid, target, N=150, grid_times_targets=True)
%timeit apo.target_rise_time(time_grid, target, N=10, grid_times_targets=True)
1.43 s ± 15.5 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
111 ms ± 597 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

so we can get order of magnitude speedups with precisions still good to a few minutes simply by exposing a keyword argument to the user. @eteq – this looks like a much more palatable enhancement to me, because:

  • the N-hancement requires no meaningful changes to astropy or astroplan! It's a three-line fix.
  • N is intuitive to explain, it's simply "the number of grid-points on which to search for the horizon crossings of the target over a 24 hour period"
  • N is an intuitive parameter to adjust – small N yields low precision at high speeds, high N yields high precision at slow speeds

Would you agree?

@iwellaway

This comment has been minimized.

Copy link
Author

@iwellaway iwellaway commented Oct 17, 2019

@eteq We pre-load our database with all of the rise and set times for each target for each day of the project. This allows us to quickly construct the actual observable period when the telescope requests a target from the scheduler program.

Essentially, we just iterate over each target, running the rise, set and meridian functions for each day and adding a new record to a database table with the results. There are over 1000 targets and the we are currently running a one year simulation so this process takes quite a long time.

Hope this helps.

@iwellaway

This comment has been minimized.

Copy link
Author

@iwellaway iwellaway commented Oct 17, 2019

@bmorris3 This sounds like a really good solution. Thanks for your work on this!

@bmorris3

This comment has been minimized.

Copy link
Contributor

@bmorris3 bmorris3 commented Dec 8, 2019

With #424 merged, I'm closing this issue for now, until further great speed improvements are necessary. Thanks @iwellaway!

@bmorris3 bmorris3 closed this Dec 8, 2019
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Linked pull requests

Successfully merging a pull request may close this issue.

None yet
3 participants
You can’t perform that action at this time.