brandon-rhodes
commented
11 months ago It turns out that this is a very big and interesting question for Skyfield, because its initial design was around having lots of times — a Time vector of every day in a whole decade, for example — and trying to make those operations NumPy-efficient. Your problem is a challenge because it's orthogonal to that: given a single time, what are n objects doing?
Instead of trying to write up a huge answer that tackles everything at once, I'm going to try making a little daily comment that moves one step along in our understanding of how we're going to do this efficiently in Skyfield! Otherwise it might be weeks before you hear anything from me. 🙂 But one little step per day should be doable for my schedule, and let you see the progress we can make towards making your problem faster.
So, today, a first observation is that Skyfield's frame methods are a little unfortunate in your case, because .frame_latlon_and_rates() calls .frame_xyz_and_velocity() (with the plan of converting its x,y,z coordinates to angles) which then thinks, okay, where was the frame pointing at the time .t that this vector was measured, so it calls the frame's .rotation_at() method to find out.
So your opos object is getting hammered with n requests, where n is the number of satellites, for the answer to the question “how exactly were you oriented at time t?” And each time it performs exactly the same expensive computation and returns the same answer. Which is wasteful: because Skyfield doesn't see the big picture that you're going to ask this over and over for n satellites, it duplicates the same work over and over.
Skyfield's design actually tries to side-step this problem in one particular case: when .altaz() is called, Skyfield accesses a cached private property named _altaz_rotation on the observer's geographic location, instead of asking it for its orientation directly. This means that it's computed only on the first try, instead of every time.
But that doesn't work for (a) the general case where the frame doesn't happen to be alt-az, or (b) even for the specific case you yourself have fallen into where the frame is in fact an alt-az horizon-based frame but you don't just want the angles altitude and azimuth and distance, but where you want their rates of change as well, and so have to call the general-purpose .frame_latlon_and_rates() method instead of the specific special .altaz() method.
Maybe Skyfield needs to switch to a more general caching mechanism that works for any frame? We'll have to see, there's lots of investigation to do here! For today, I'll just suggest you try adding a cache yourself to the opos, and see if you get a speedup. Try this:
from functools import lru_cache
opos.rotation_at = lru_cache(maxsize=None)(opos.rotation_at)
opos._dRdt_times_RT_at = lru_cache(maxsize=None)(opos._dRdt_times_RT_at)That teaches the opos object to cache repeated calls to its frame-rotation methods if the same time gets asked about multiple times. Add it at the top of your code, and let me know if you see a difference — that would validate for us our initial idea here, that a big expense is that the frame is being asked for its orientation repeatedly!
Sec42
Hi,
I have a list of reception results (timestamp, frequency) that I am trying to match a specific satellite within a TLE file by matching the doppler shift of the signal to the speed of the satellite.
I was able to produce a working example using
frame_latlon_and_rates, however the result is quite slow (in the order of 30 matches / second on my test server).My input data contains about 14M measurements from a 24h period, so I am at least one order of magnitude away from it being usable.
Am I making some obvious mistake here? Most of the satellites are below the horizon, but as elevation is also a result of
frame_latlon_and_rateschecking it doesn't really give an improvement.I am hoping you can help here.
The snippet requires an TLE file (e.g. iridium-NEXT.txt) and the input file looks like this: