[Spice_discussion] Re: how accurate is conics?

[Nat Bachman]

Hi Trent,

Thanks for your professional-quality messages.
Getting test code and a test meta-kernel makes life a
lot easier.

Your code is fine; it's just telling you something about
reality you may not have expected.

When you extract osculating elements from an ephemeris,
you get, at best, just that: an "osculating" solution: the
state of the osculating solution matches that of the input
ephemeris at the epoch of the elements, and the rest of the
osculating solution is a kinematic extrapolation.
Any noise in the data will be reflected in the elements,
so long term extrapolation of osculating elements should
be expected to give poor to horrible results unless you're starting
from ideal, artificial data that describe true two-body motion
and don't describe a near-singular case.

Elements that have been smoothed can give much better results
if the orbit is well approximated by a two-body model. But in reality,
two-body models yield at best crude approximations to actual
orbits of objects within the solar system.

Although I don't see a problem in this case, orbits whose
elements are near the singularities of the mapping from
states to elements are cases where OSCELT can't be expected
to give good results. Examples include near-parabolic orbits,
and very nearly circular or equatorial orbits.

Finally, there *is* a numerical problem with CONICS, but I
don't believe it's an issue here. The Stumpff function routine
upon which CONICS relies produces too much round-off error.
I've linked your program against a corrected routine (some variant
of which will be released in the N0065 SPICE Toolkit), but I
saw no significant change in the results.

Here are some tests that should give you a feel for what OSCELT
does and how well extrapolation of elements can work for your
test case:

    1) Run your program with a short time step, say 1 hour.
       If the elements make any sense at all, you should get
       position agreement at the few-meter level.

    2) Get an extrapolated state from CONICS, use OSCELT on that
       state, and see whether you can recover the time-independent
       elements obtained from your initial OSCELT call. Also
       extrapolate backwards to the epoch of your initial OSCELT
       call to try to recover the initial state.

    3) Run your program with elements provided by the JPL Horizons
       program instead of those from OSCELT. These elements
       likely have been selected to fit all of the observational data.

You also can use the SPICE utility MKSPK to produce your own SPK files
from conic elements; your test should work well with such files unless
the elements represent near-singular cases. I don't think there's anything
new to learn about osculating elements from this; it's just another kind
of consistency check.

A few remarks on your test code:

    1) You can get the number of seconds in a Julian year
       from the SPICE routine jyear_c. (This is the number of
       seconds in 365.25 Julian days.)

    2) You can get NAIF ID codes from body names by calling bodn2c_c
       or bods2c_c.

       You can create new name-ID code associations at run time
       via assignments in your meta-kernel (or other text kernel)
       See the NAIF IDS tutorial for details.

    3) You can get solar and planetary GM values from the kernel

          de-403-masses.tpc

       which is available in the "PCK" data area on the NAIF
       web site.

       You can fetch these values into your code using

          gdpool_c


Regarding two-body approximations for the lunar orbit: try
comparing the magnitudes of the gravitational forces exerted
on the moon by the earth and by the sun. Then answer the
question "is this two-body motion"? :)

Best regards,

   -Nat

Nathaniel.Bachman at jpl.nasa.gov

















Trent Waddington wrote:
> On Thu, Sep 2, 2010 at 1:00 PM, Trent Waddington
> <trent.waddington at gmail.com> wrote:
>> I'm trying to get a feel for how much ephemeris generated by the
>> HORIZONS system differs from the mathematical ideal.  I've been
>> looking at asteroids but, as a bit of a sanity check, I decided to
>> write some code to compare the output of conics with the precalculated
>> ephemeris of the Moon.
> 
> It occurs to me that this is anything but sane :P
> 
> Here's the same code but using heliocentric osculating elements, for
> an arbitrary SPK.  I've been using ones generated for asteroids..
> 
>     http://quantumg.net/pconics.c
>     http://quantumg.net/2006_RH120.bsp
> 
> with a step of 1 year the position is 43902km wrong.. that seems like a lot.
> 
> Thanks,
> 
> Trent
> 
> _______________________________________________
> Spice_discussion mailing list
> Spice_discussion at naif.jpl.nasa.gov
> http://naif.jpl.nasa.gov/mailman/listinfo/spice_discussion