cdfepoch.breakdown_tt2000 returns different results for scalar vs array

[jameswilburlewis]

Consider the following test code, which makes some calls to cdfepoch.breakdown_tt2000 for TT2000 timestamps just before (t1), during (t2), and after (t3) the leap second that was added at the end of 2016:

    def test_cdfepoch_leap(self):
        from cdflib import cdfepoch
        t1=np.int64(536500867184000000)
        t2=t1+1000000000
        t3=t2+1000000000
        b3=cdfepoch.breakdown_tt2000(t3)
        print("t3 ",b3)
        b4=cdfepoch.breakdown_tt2000([t2,t3])
        print("t2 t3",b4)
        b5=cdfepoch.breakdown_tt2000([t1,t3])
        print("t1 t3", b5)

Output:

t3  [2017    1    1    0    0    0    0    0    0]
t2 t3 [[2016   12   31   23   60    0    0    0    0]
 [2016   12   31   23   60    0    0    0    0]]
t1 t3 [[2016   12   31   23   59   59    0    0    0]
 [2016   12   31   23   60    0    0    0    0]]

Note that if t3 is provided as a scalar, it converts as expected (2017-01-01/00:00:00), but if it's an element of an array that spans the leap second boundary, it's broken down incorrectly as 2016-12-31/23:60:00. And Python's datetime.datetime will choke on it.

Here's an excerpt from breakdown_tt200:

        datxzero = datxs[:, 1] == 0.0
        epochs[post72 & ~datxzero] -= 1
        xdates = CDFepoch._EPOCHbreakdownTT2000(epochs)

        # If 1 second was subtracted, add 1 second back in
        # Be careful not to go 60 or above
        xdates[5, post72 & ~datxzero] += 1
        xdates[4, post72 & ~datxzero] += np.floor(xdates[5, post72 & ~datxzero] / 60.0)
        xdates[5, post72 & ~datxzero] = xdates[5, post72 & ~datxzero] % 60

I think the intention here might to prevent the seconds field from getting a value of 60 during a leap second. But it just pushes the problem from the seconds field to the minutes field. And it was weird to see the different results for the same value, depending on whether it was a scalar or an element of an array. To explain that, we need to look at _LeapSecondsfromJ2000(), which sets datxs. The datxs[:,1] component is an overflow flag that I think is supposed to mark leap seconds, and it's getting set differently for the scalar and array cases:

        if CDFepoch.NST is None:
            CDFepoch._LoadLeapNanoSecondsTable()
        for i, _ in reversed(list(enumerate(CDFepoch.NST))):
            idxs = (j == -1) & (nanosecs >= CDFepoch.NST[i])
            j[idxs] = i
            if i < (CDFepoch.NDAT - 1):
                overflow = nanosecs + 1000000000 >= CDFepoch.NST[i + 1]
                da[overflow, 1] = 1.0
            if np.all(j > 0):
                break

For this particular test case, If the input is a scalar, this branch will not be executed:

            if i < (CDFepoch.NDAT - 1):
                overflow = nanosecs + 1000000000 >= CDFepoch.NST[i + 1]
                da[overflow, 1] = 1.0

However, if the input is an array that spans the most recent leap second boundary, that branch does get executed, and it looks to me like the overflow flag will be set for ALL values during or after the final leap second, not just the values close to the boundary.

So, I think several changes are needed:

1. Add a clause to the overflow test to prevent dates beyond the final leap second from being marked
2. If the seconds field reaches 60, don't just roll it over to the minutes field; that change might need to propagate all the way up to the year (but not always, because leap seconds can also occur in June)
3. (maybe) Alternatively, let the seconds field remain at 60, but test for that before making a datetime out of it. If necessary, change the seconds field to 59, then make the datetime, then use the datetime methods to add back the extra second. This might need to be checked in several places, but it might be easier and clearer than replicating the rollover logic that's already built into datetime.
4. Add some regression tests to catch these and similar scenarios

For the very near term, I will probably just work around it with some exception handling in the PySPEDAS TT2000 code, but in a few weeks I might have time to look at this a bit more, and hopefully put a PR together to address the issues.

[jameswilburlewis]

Closely related issue: #296