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

Minor Planet Center interface for retrieving observations #110

Merged
merged 20 commits into from
Oct 30, 2023

Conversation

tristandijkstra
Copy link
Contributor

@tristandijkstra tristandijkstra commented Sep 6, 2023

This pull request adds an interface for retrieving observations from the Minor Planet Center (MPC) for estimation. Most notably, this PR includes a class called BatchMPC() which enables automatic retrieval of observations for multiple bodies, the filtering of said observations and automatic conversion to a Tudat ObservationCollection and link settings. Estimation for one or more objects observed by any number of terestrial satellites has been tested. Support for space telescope observations is included but has not been tested yet.

This PR also adds a new submodule data (from tudatpy.data.mpc import BatchMPC). This new submodule is python exclusive and will include more interfaces in the future.

Two examples for BatchMPC will be included here: tudat-team/tudatpy-examples#28

@DominicDirkx
Copy link
Member

Additional note: astropy should now be added as a dependency in the environment.yaml files

@DominicDirkx
Copy link
Member

Small note, this addition:
90bc0b2
is required to make the new module accessible in the 'regular' setup


"""

def __init__(self) -> None:
Copy link
Member

Choose a reason for hiding this comment

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

I would recommend using a SystemOfBodies object (bodies variable in most applications) as input here. This would then be an input to the _add_observatory_positions function, which would make sure that the station position definition is consistent with the Earth shape model that is provided (see TODO entry in that function)

Copy link
Contributor Author

Choose a reason for hiding this comment

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

I have tried this solution but I found it detracting from the user experience. Most notably, it felt unintuitive to set up systemofbodies and its dependencies before having at least a cursory look at your observations. As an alternative, I have moved the _add_observatory_positions() method call to the to_tudat(). The only downside is that users would have to run to_tudat() first to see carthesian positions of the observatories (has been documented). However it does mean that the use of BatchMPC remains the same as before.

-------
estimation.ObservationCollection
ObservationCollection with the observations found in the batch
Dict[str, observation.LinkDefinition]
Copy link
Member

Choose a reason for hiding this comment

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

Why do we need the dict[str, LinkDefinition] here? This information should also be in the ObservationCollection (?)

Copy link
Contributor Author

Choose a reason for hiding this comment

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

This has been fixed, users must now retrieve links from the ObservationCollection

@tristandijkstra
Copy link
Contributor Author

In addition to the above comments 2 features have been suggested in the companion examples PR (tudat-team/tudatpy-examples#28):

  1. A plot with right ascension as function of time, and declination as function of time.
  2. included_satellites is now a mandatory parameter of to_tudat() to avoid confusion of where satellite data is filtered.

With that this PR should be complete pending the completion and testing of the estimation example of tudat-team/tudatpy-examples#28

@tristandijkstra
Copy link
Contributor Author

tristandijkstra commented Oct 18, 2023

One more thing!

I've added a quick implementation of the JPL small body database (https://ssd.jpl.nasa.gov/tools/sbdb_lookup.html#/) under data.sbdb, this has been done to make it easier to retrieve spk codes and names of the asteroids, which MPC does not provide.

Please have a look here. In the future this could also be used to add masses to the bodies.

I will use it in the estimation example of this PR to reduce the amount of hardcoding required: tudat-team/tudatpy-examples#28

@DominicDirkx
Copy link
Member

Hi Tristan, thanks for the updates. Great to also have the link to the SBDB :) Two remaining question:

  • Concerning where to provide the bodies to the BatchMPC class, it's fine to add this later than the constructor. Could you verify that, also in this new solution, the geographic coordinates of the stations in the MPC database will always be consistent with the Earth model that is used?
  • I think we discussed this (but I may have been mistaken, and it may have been in the context of the Horizons database), that it would be good to have a small unit test for this code. For instance: retrieve data from some object, compare positions to Spice, and verify that the RMS difference is small.

@tristandijkstra
Copy link
Contributor Author

tristandijkstra commented Oct 28, 2023

Hi Dominic,

  • Yes, the radius is consistent with the Earth body provided to to_tudat(), thus updating anytime the user retrieves new observations.
  • We have discussed this, I think we decided it was a bad idea to make a unit test dependant on estimation, as we wouldnt know if the error is due to a change in the estimation module or MPC, similarly we might change the spice kernel. The alternative was to compare to RA+DEC to JPL horizons directly (in a future PR). For now, some basic tests checking the results of a retrieval, filter and to_tudat() could also be good.

In general, I have a list of future additions to the data submodule (both MPC and Horizons) that will be implemented in one or more future PRs. I will later formulate this list into an issue.

def codes_300_spkid(self):
"""Returns spice kernel number for the codes_300ast_20100725.bsp spice kernel"""
spkid = self.spkid[0] + self.spkid[2:]
if spkid == "2000001":
Copy link
Member

Choose a reason for hiding this comment

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

Small question: why are some of the asteroids returned by name, and other by ID?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

It's this really annoying thing with the codes_300ast_20100725.bsp spice kernel. For some asteroids, specifically those with an if statement in this fucntion, only identify by their name E.g. Eros. All other objects in this kernel identify with a number. For some reason, there is one zero extra in the id JPL gives, compared to the codes_300ast_20100725.bsp kernel, thats why I do line 0. To allow for for an automatic process where you only have to change the list of MPC codes in my examples, I have made this which retrieves the correct code for any body in codes_300ast_20100725.bsp. The function above (spkid) gives the general ID, this function is for codes_300ast_20100725.bsp specifically. Otherwise I would have to do some if statement in the examples or a try except.

See also: https://naif.jpl.nasa.gov/pub/naif/generic_kernels/spk/asteroids/aa_summaries.txt
And: tudat-team/tudatpy-examples#28 (comment), point 2.5

Copy link
Member

Choose a reason for hiding this comment

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

Right, thanks! Could you add a note to the comments of the code on why these specific cases are checked. Someone looking at these some months or years from now will have the same question

@DominicDirkx DominicDirkx merged commit 0bb610e into tudat-team:develop Oct 30, 2023
1 of 2 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

2 participants