Format question: difference between segment.start_jd and initial_epoch

[ChristopherRabotin]

This might be a file format question, so I understand if you prefer me to offset this to the NAIF team.

In spk.generate, you use the initial epoch to compute the position of the object. In spk.init, you compute the start JD from the start second, which is ephemeris time I believe.

How come there is a difference of 4 days between both values in the de438.bsp file (2415016.5 != 2415020.5)? Similar question for the 3 day difference between the end_jde and the end_second.

Thanks

[brandon-rhodes]

Which number is the smaller date? Is it possible that they've delivered a set of coefficients whose parameters, for the sake of better fitting or something, extend to a few days beyond the end of the time that the ephemeris is actually valid for? I'd be happy if the promised dates in the header were inside of the dates given along with the coefficients; but unhappy if they weren't.

[ChristopherRabotin]

The initial epoch is the the smaller of the two dates. This is the latest planetary ephemeris data from JPL, de438. I can check the dates in the comments of that file to check if there is an issue in the file itself.

[brandon-rhodes]

Good, so it is in the direction that I had hoped! They have supplied data for a curve between two dates, but only promised in the segment header that it's valid for a narrower range of dates. I'm not sure of their technical reason for doing so (it's really not valid over the whole span? or so that it doesn't conflict with the next file they release?) but I don't think any software using it would have any problem if the segment-header dates are within the range of the polynomial end dates.

[ChristopherRabotin]

Okay. The confusion I'm having now is that I have extracted the start and end modified Julian dates as they are in the file, but I don't know how to perform the interpolation if I've got the wrong dates. Should I use the "initial epoch" of the segment instead of the start_jd value then?

On Sat, Jan 5, 2019, 17:27 Brandon Rhodes notifications@github.com wrote:

Good, so it is in the direction that I had hoped! They have supplied data for a curve between two dates, but only promised in the segment header that it's valid for a narrower range of dates. I'm not sure of their technical reason for doing so (it's really not valid over the whole span? or so that it doesn't conflict with the next file they release?) but I don't think any software using it would have any problem if the segment-header dates are within the range of the polynomial end dates.

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/brandon-rhodes/python-jplephem/issues/29#issuecomment-451704533, or mute the thread https://github.com/notifications/unsubscribe-auth/AEma6EOVizOJtjFEY3SZmeqAIaPudMJ6ks5vAUMHgaJpZM4ZpTMT .

[brandon-rhodes]

My guess is that the initial epoch is the correct date from which to form the parameter that you'll insert into the interpolation formulae. But — could you maybe try it both ways, and see if one answer makes sense and the other looks obviously wrong?

[ChristopherRabotin]

So it turned out that the problem I was having is due to, I think, a difference in how the Chebyshev polynominals are computed between the NAIF BSP way and the traditional way. I've just now fixed the issue (I had to work on other stuff until now).

Here is the Python implementation for BSPs:

            T = [0] * coefficient_count
            T[0] = 1.0
            T[1] = t1 = 2.0 * offset / interval_length - 1.0
            twot1 = t1 + t1
            for i in range(2, coefficient_count):
                T[i] = twot1 * T[i - 1] - T[i - 2]

And here is the chebval implementation of numpy:

        x2 = 2*x
        c0 = c[-2]
        c1 = c[-1]
        for i in range(3, len(c) + 1):
            tmp = c0
            c0 = c[-i] - c1
        c1 = tmp + c1*x2

[brandon-rhodes]

I'm glad you were able to resolve the problem! If you're ever interested in Skyfield having native support for the kind of file you were trying to read, I'd be happy to look at a pull request that includes your code. Either way, I'm glad you're no longer blocked.