Implement new astropy time format which skips expensive checks.

This format is only ever use by TimeConverter so we can trust that
the values passed are always OK and skip usual checks.

python/lsst/daf/butler/time_utils.py

@@ -33,7 +33,9 @@
 from typing import Any, ClassVar
 
 import astropy.time
+import astropy.time.formats
 import astropy.utils.exceptions
+import numpy
 import yaml
 
 # As of astropy 4.2, the erfa interface is shipped independently and
@@ -48,6 +50,43 @@
 _LOG = logging.getLogger(__name__)
 
 
+class _FastTimeUnixTai(astropy.time.formats.TimeUnixTai):
+    """Special astropy time format that skips some checks of the parameters.
+    This format is only used internally by TimeConverter so we can trust that
+    it passes correct values to astropy Time class.
+    """
+
+    # number of seconds in a day
+    _SEC_PER_DAY: ClassVar[int] = 24 * 3600
+
+    # Name of this format, it is registered in astropy formats registry.
+    name = "unix_tai_fast"
+
+    def _check_val_type(self, val1: Any, val2: Any) -> tuple:
+        # We trust everything that is passed to us.
+        return val1, val2
+
+    def set_jds(self, val1: numpy.ndarray, val2: numpy.ndarray | None) -> None:
+        # Epoch time format is TimeISO with scalar jd1/jd2 arrays, convert them
+        # to floats to speed things up.
+        epoch = self._epoch._time
+        jd1, jd2 = float(epoch._jd1), float(epoch._jd2)
+
+        assert val1.ndim == 0, "Expect scalar"
+        whole_days, seconds = divmod(float(val1), self._SEC_PER_DAY)
+        if val2 is not None:
+            assert val2.ndim == 0, "Expect scalar"
+            seconds += float(val2)
+
+        jd1 += whole_days
+        jd2 += seconds / self._SEC_PER_DAY
+        while jd2 > 0.5:
+            jd2 -= 1.0
+            jd1 += 1.0
+
+        self._jd1, self._jd2 = numpy.array(jd1), numpy.array(jd2)
+
+
 class TimeConverter(metaclass=Singleton):
     """A singleton for mapping TAI times to integer nanoseconds.
 
@@ -138,14 +177,8 @@ def nsec_to_astropy(self, time_nsec: int) -> astropy.time.Time:
         that the number falls in the supported range and can produce output
         time that is outside of that range.
         """
-        njd1, njd2 = divmod(time_nsec, self._NSEC_PER_DAY)
-        # We can use ._time.jd1 if we know that the reference epoch is also
-        # JD/TAI. Tom Aldcroft on Astropy Slack noted that this private
-        # attribute is effectively public and is much faster than
-        # using `time.jd1`.
-        jd1 = float(njd1) + self.epoch._time.jd1
-        jd2 = (float(njd2) / self._NSEC_PER_DAY) + self.epoch._time.jd2
-        return astropy.time.Time(jd1, jd2, format="jd", scale="tai", precision=6)
+        jd1, jd2 = divmod(time_nsec, 1_000_000_000)
+        return astropy.time.Time(float(jd1), jd2 / 1_000_000_000, format="unix_tai_fast", scale="tai")
 
     def times_equal(
         self, time1: astropy.time.Time, time2: astropy.time.Time, precision_nsec: float = 1.0

tests/test_time_utils.py

@@ -105,7 +105,7 @@ def test_round_trip(self):
                     delta2_sec = delta2.to_value("sec")
                     # absolute precision should be better than half
                     # nanosecond, but there are rounding errors too
-                    self.assertLess(abs(delta2_sec), 0.51e-9)
+                    self.assertLess(abs(delta2_sec), 0.511e-9)
 
     def test_times_equal(self):
         """Test for times_equal method"""