1# -*- coding: utf-8 -*- 2# Licensed under a 3-clause BSD style license - see LICENSE.rst 3""" 4The astropy.time package provides functionality for manipulating times and 5dates. Specific emphasis is placed on supporting time scales (e.g. UTC, TAI, 6UT1) and time representations (e.g. JD, MJD, ISO 8601) that are used in 7astronomy. 8""" 9 10import os 11import copy 12import enum 13import operator 14import threading 15from datetime import datetime, date, timedelta 16from time import strftime 17from warnings import warn 18 19import numpy as np 20import erfa 21 22from astropy import units as u, constants as const 23from astropy.units import UnitConversionError 24from astropy.utils import ShapedLikeNDArray 25from astropy.utils.compat.misc import override__dir__ 26from astropy.utils.data_info import MixinInfo, data_info_factory 27from astropy.utils.exceptions import AstropyWarning 28from .utils import day_frac 29from .formats import (TIME_FORMATS, TIME_DELTA_FORMATS, 30 TimeJD, TimeUnique, TimeAstropyTime, TimeDatetime) 31# Import TimeFromEpoch to avoid breaking code that followed the old example of 32# making a custom timescale in the documentation. 33from .formats import TimeFromEpoch # noqa 34 35from astropy.extern import _strptime 36 37__all__ = ['TimeBase', 'Time', 'TimeDelta', 'TimeInfo', 'update_leap_seconds', 38 'TIME_SCALES', 'STANDARD_TIME_SCALES', 'TIME_DELTA_SCALES', 39 'ScaleValueError', 'OperandTypeError'] 40 41 42STANDARD_TIME_SCALES = ('tai', 'tcb', 'tcg', 'tdb', 'tt', 'ut1', 'utc') 43LOCAL_SCALES = ('local',) 44TIME_TYPES = dict((scale, scales) for scales in (STANDARD_TIME_SCALES, LOCAL_SCALES) 45 for scale in scales) 46TIME_SCALES = STANDARD_TIME_SCALES + LOCAL_SCALES 47MULTI_HOPS = {('tai', 'tcb'): ('tt', 'tdb'), 48 ('tai', 'tcg'): ('tt',), 49 ('tai', 'ut1'): ('utc',), 50 ('tai', 'tdb'): ('tt',), 51 ('tcb', 'tcg'): ('tdb', 'tt'), 52 ('tcb', 'tt'): ('tdb',), 53 ('tcb', 'ut1'): ('tdb', 'tt', 'tai', 'utc'), 54 ('tcb', 'utc'): ('tdb', 'tt', 'tai'), 55 ('tcg', 'tdb'): ('tt',), 56 ('tcg', 'ut1'): ('tt', 'tai', 'utc'), 57 ('tcg', 'utc'): ('tt', 'tai'), 58 ('tdb', 'ut1'): ('tt', 'tai', 'utc'), 59 ('tdb', 'utc'): ('tt', 'tai'), 60 ('tt', 'ut1'): ('tai', 'utc'), 61 ('tt', 'utc'): ('tai',), 62 } 63GEOCENTRIC_SCALES = ('tai', 'tt', 'tcg') 64BARYCENTRIC_SCALES = ('tcb', 'tdb') 65ROTATIONAL_SCALES = ('ut1',) 66TIME_DELTA_TYPES = dict((scale, scales) 67 for scales in (GEOCENTRIC_SCALES, BARYCENTRIC_SCALES, 68 ROTATIONAL_SCALES, LOCAL_SCALES) for scale in scales) 69TIME_DELTA_SCALES = GEOCENTRIC_SCALES + BARYCENTRIC_SCALES + ROTATIONAL_SCALES + LOCAL_SCALES 70# For time scale changes, we need L_G and L_B, which are stored in erfam.h as 71# /* L_G = 1 - d(TT)/d(TCG) */ 72# define ERFA_ELG (6.969290134e-10) 73# /* L_B = 1 - d(TDB)/d(TCB), and TDB (s) at TAI 1977/1/1.0 */ 74# define ERFA_ELB (1.550519768e-8) 75# These are exposed in erfa as erfa.ELG and erfa.ELB. 76# Implied: d(TT)/d(TCG) = 1-L_G 77# and d(TCG)/d(TT) = 1/(1-L_G) = 1 + (1-(1-L_G))/(1-L_G) = 1 + L_G/(1-L_G) 78# scale offsets as second = first + first * scale_offset[(first,second)] 79SCALE_OFFSETS = {('tt', 'tai'): None, 80 ('tai', 'tt'): None, 81 ('tcg', 'tt'): -erfa.ELG, 82 ('tt', 'tcg'): erfa.ELG / (1. - erfa.ELG), 83 ('tcg', 'tai'): -erfa.ELG, 84 ('tai', 'tcg'): erfa.ELG / (1. - erfa.ELG), 85 ('tcb', 'tdb'): -erfa.ELB, 86 ('tdb', 'tcb'): erfa.ELB / (1. - erfa.ELB)} 87 88# triple-level dictionary, yay! 89SIDEREAL_TIME_MODELS = { 90 'mean': { 91 'IAU2006': {'function': erfa.gmst06, 'scales': ('ut1', 'tt')}, 92 'IAU2000': {'function': erfa.gmst00, 'scales': ('ut1', 'tt')}, 93 'IAU1982': {'function': erfa.gmst82, 'scales': ('ut1',), 'include_tio': False} 94 }, 95 'apparent': { 96 'IAU2006A': {'function': erfa.gst06a, 'scales': ('ut1', 'tt')}, 97 'IAU2000A': {'function': erfa.gst00a, 'scales': ('ut1', 'tt')}, 98 'IAU2000B': {'function': erfa.gst00b, 'scales': ('ut1',)}, 99 'IAU1994': {'function': erfa.gst94, 'scales': ('ut1',), 'include_tio': False} 100 }} 101 102 103class _LeapSecondsCheck(enum.Enum): 104 NOT_STARTED = 0 # No thread has reached the check 105 RUNNING = 1 # A thread is running update_leap_seconds (_LEAP_SECONDS_LOCK is held) 106 DONE = 2 # update_leap_seconds has completed 107 108 109_LEAP_SECONDS_CHECK = _LeapSecondsCheck.NOT_STARTED 110_LEAP_SECONDS_LOCK = threading.RLock() 111 112 113class TimeInfo(MixinInfo): 114 """ 115 Container for meta information like name, description, format. This is 116 required when the object is used as a mixin column within a table, but can 117 be used as a general way to store meta information. 118 """ 119 attr_names = MixinInfo.attr_names | {'serialize_method'} 120 _supports_indexing = True 121 122 # The usual tuple of attributes needed for serialization is replaced 123 # by a property, since Time can be serialized different ways. 124 _represent_as_dict_extra_attrs = ('format', 'scale', 'precision', 125 'in_subfmt', 'out_subfmt', 'location', 126 '_delta_ut1_utc', '_delta_tdb_tt') 127 128 # When serializing, write out the `value` attribute using the column name. 129 _represent_as_dict_primary_data = 'value' 130 131 mask_val = np.ma.masked 132 133 @property 134 def _represent_as_dict_attrs(self): 135 method = self.serialize_method[self._serialize_context] 136 if method == 'formatted_value': 137 out = ('value',) 138 elif method == 'jd1_jd2': 139 out = ('jd1', 'jd2') 140 else: 141 raise ValueError("serialize method must be 'formatted_value' or 'jd1_jd2'") 142 143 return out + self._represent_as_dict_extra_attrs 144 145 def __init__(self, bound=False): 146 super().__init__(bound) 147 148 # If bound to a data object instance then create the dict of attributes 149 # which stores the info attribute values. 150 if bound: 151 # Specify how to serialize this object depending on context. 152 # If ``True`` for a context, then use formatted ``value`` attribute 153 # (e.g. the ISO time string). If ``False`` then use float jd1 and jd2. 154 self.serialize_method = {'fits': 'jd1_jd2', 155 'ecsv': 'formatted_value', 156 'hdf5': 'jd1_jd2', 157 'yaml': 'jd1_jd2', 158 'parquet': 'jd1_jd2', 159 None: 'jd1_jd2'} 160 161 def get_sortable_arrays(self): 162 """ 163 Return a list of arrays which can be lexically sorted to represent 164 the order of the parent column. 165 166 Returns 167 ------- 168 arrays : list of ndarray 169 """ 170 parent = self._parent 171 jd_approx = parent.jd 172 jd_remainder = (parent - parent.__class__(jd_approx, format='jd')).jd 173 return [jd_approx, jd_remainder] 174 175 @property 176 def unit(self): 177 return None 178 179 info_summary_stats = staticmethod( 180 data_info_factory(names=MixinInfo._stats, 181 funcs=[getattr(np, stat) for stat in MixinInfo._stats])) 182 # When Time has mean, std, min, max methods: 183 # funcs = [lambda x: getattr(x, stat)() for stat_name in MixinInfo._stats]) 184 185 def _construct_from_dict_base(self, map): 186 if 'jd1' in map and 'jd2' in map: 187 # Initialize as JD but revert to desired format and out_subfmt (if needed) 188 format = map.pop('format') 189 out_subfmt = map.pop('out_subfmt', None) 190 map['format'] = 'jd' 191 map['val'] = map.pop('jd1') 192 map['val2'] = map.pop('jd2') 193 out = self._parent_cls(**map) 194 out.format = format 195 if out_subfmt is not None: 196 out.out_subfmt = out_subfmt 197 198 else: 199 map['val'] = map.pop('value') 200 out = self._parent_cls(**map) 201 202 return out 203 204 def _construct_from_dict(self, map): 205 delta_ut1_utc = map.pop('_delta_ut1_utc', None) 206 delta_tdb_tt = map.pop('_delta_tdb_tt', None) 207 208 out = self._construct_from_dict_base(map) 209 210 if delta_ut1_utc is not None: 211 out._delta_ut1_utc = delta_ut1_utc 212 if delta_tdb_tt is not None: 213 out._delta_tdb_tt = delta_tdb_tt 214 215 return out 216 217 def new_like(self, cols, length, metadata_conflicts='warn', name=None): 218 """ 219 Return a new Time instance which is consistent with the input Time objects 220 ``cols`` and has ``length`` rows. 221 222 This is intended for creating an empty Time instance whose elements can 223 be set in-place for table operations like join or vstack. It checks 224 that the input locations and attributes are consistent. This is used 225 when a Time object is used as a mixin column in an astropy Table. 226 227 Parameters 228 ---------- 229 cols : list 230 List of input columns (Time objects) 231 length : int 232 Length of the output column object 233 metadata_conflicts : str ('warn'|'error'|'silent') 234 How to handle metadata conflicts 235 name : str 236 Output column name 237 238 Returns 239 ------- 240 col : Time (or subclass) 241 Empty instance of this class consistent with ``cols`` 242 243 """ 244 # Get merged info attributes like shape, dtype, format, description, etc. 245 attrs = self.merge_cols_attributes(cols, metadata_conflicts, name, 246 ('meta', 'description')) 247 attrs.pop('dtype') # Not relevant for Time 248 col0 = cols[0] 249 250 # Check that location is consistent for all Time objects 251 for col in cols[1:]: 252 # This is the method used by __setitem__ to ensure that the right side 253 # has a consistent location (and coerce data if necessary, but that does 254 # not happen in this case since `col` is already a Time object). If this 255 # passes then any subsequent table operations via setitem will work. 256 try: 257 col0._make_value_equivalent(slice(None), col) 258 except ValueError: 259 raise ValueError('input columns have inconsistent locations') 260 261 # Make a new Time object with the desired shape and attributes 262 shape = (length,) + attrs.pop('shape') 263 jd2000 = 2451544.5 # Arbitrary JD value J2000.0 that will work with ERFA 264 jd1 = np.full(shape, jd2000, dtype='f8') 265 jd2 = np.zeros(shape, dtype='f8') 266 tm_attrs = {attr: getattr(col0, attr) 267 for attr in ('scale', 'location', 268 'precision', 'in_subfmt', 'out_subfmt')} 269 out = self._parent_cls(jd1, jd2, format='jd', **tm_attrs) 270 out.format = col0.format 271 272 # Set remaining info attributes 273 for attr, value in attrs.items(): 274 setattr(out.info, attr, value) 275 276 return out 277 278 279class TimeDeltaInfo(TimeInfo): 280 _represent_as_dict_extra_attrs = ('format', 'scale') 281 282 def _construct_from_dict(self, map): 283 return self._construct_from_dict_base(map) 284 285 def new_like(self, cols, length, metadata_conflicts='warn', name=None): 286 """ 287 Return a new TimeDelta instance which is consistent with the input Time objects 288 ``cols`` and has ``length`` rows. 289 290 This is intended for creating an empty Time instance whose elements can 291 be set in-place for table operations like join or vstack. It checks 292 that the input locations and attributes are consistent. This is used 293 when a Time object is used as a mixin column in an astropy Table. 294 295 Parameters 296 ---------- 297 cols : list 298 List of input columns (Time objects) 299 length : int 300 Length of the output column object 301 metadata_conflicts : str ('warn'|'error'|'silent') 302 How to handle metadata conflicts 303 name : str 304 Output column name 305 306 Returns 307 ------- 308 col : Time (or subclass) 309 Empty instance of this class consistent with ``cols`` 310 311 """ 312 # Get merged info attributes like shape, dtype, format, description, etc. 313 attrs = self.merge_cols_attributes(cols, metadata_conflicts, name, 314 ('meta', 'description')) 315 attrs.pop('dtype') # Not relevant for Time 316 col0 = cols[0] 317 318 # Make a new Time object with the desired shape and attributes 319 shape = (length,) + attrs.pop('shape') 320 jd1 = np.zeros(shape, dtype='f8') 321 jd2 = np.zeros(shape, dtype='f8') 322 out = self._parent_cls(jd1, jd2, format='jd', scale=col0.scale) 323 out.format = col0.format 324 325 # Set remaining info attributes 326 for attr, value in attrs.items(): 327 setattr(out.info, attr, value) 328 329 return out 330 331 332class TimeBase(ShapedLikeNDArray): 333 """Base time class from which Time and TimeDelta inherit.""" 334 335 # Make sure that reverse arithmetic (e.g., TimeDelta.__rmul__) 336 # gets called over the __mul__ of Numpy arrays. 337 __array_priority__ = 20000 338 339 # Declare that Time can be used as a Table column by defining the 340 # attribute where column attributes will be stored. 341 _astropy_column_attrs = None 342 343 def __getnewargs__(self): 344 return (self._time,) 345 346 def _init_from_vals(self, val, val2, format, scale, copy, 347 precision=None, in_subfmt=None, out_subfmt=None): 348 """ 349 Set the internal _format, scale, and _time attrs from user 350 inputs. This handles coercion into the correct shapes and 351 some basic input validation. 352 """ 353 if precision is None: 354 precision = 3 355 if in_subfmt is None: 356 in_subfmt = '*' 357 if out_subfmt is None: 358 out_subfmt = '*' 359 360 # Coerce val into an array 361 val = _make_array(val, copy) 362 363 # If val2 is not None, ensure consistency 364 if val2 is not None: 365 val2 = _make_array(val2, copy) 366 try: 367 np.broadcast(val, val2) 368 except ValueError: 369 raise ValueError('Input val and val2 have inconsistent shape; ' 370 'they cannot be broadcast together.') 371 372 if scale is not None: 373 if not (isinstance(scale, str) 374 and scale.lower() in self.SCALES): 375 raise ScaleValueError("Scale {!r} is not in the allowed scales " 376 "{}".format(scale, 377 sorted(self.SCALES))) 378 379 # If either of the input val, val2 are masked arrays then 380 # find the masked elements and fill them. 381 mask, val, val2 = _check_for_masked_and_fill(val, val2) 382 383 # Parse / convert input values into internal jd1, jd2 based on format 384 self._time = self._get_time_fmt(val, val2, format, scale, 385 precision, in_subfmt, out_subfmt) 386 self._format = self._time.name 387 388 # Hack from #9969 to allow passing the location value that has been 389 # collected by the TimeAstropyTime format class up to the Time level. 390 # TODO: find a nicer way. 391 if hasattr(self._time, '_location'): 392 self.location = self._time._location 393 del self._time._location 394 395 # If any inputs were masked then masked jd2 accordingly. From above 396 # routine ``mask`` must be either Python bool False or an bool ndarray 397 # with shape broadcastable to jd2. 398 if mask is not False: 399 mask = np.broadcast_to(mask, self._time.jd2.shape) 400 self._time.jd1[mask] = 2451544.5 # Set to JD for 2000-01-01 401 self._time.jd2[mask] = np.nan 402 403 def _get_time_fmt(self, val, val2, format, scale, 404 precision, in_subfmt, out_subfmt): 405 """ 406 Given the supplied val, val2, format and scale try to instantiate 407 the corresponding TimeFormat class to convert the input values into 408 the internal jd1 and jd2. 409 410 If format is `None` and the input is a string-type or object array then 411 guess available formats and stop when one matches. 412 """ 413 414 if (format is None 415 and (val.dtype.kind in ('S', 'U', 'O', 'M') or val.dtype.names)): 416 # Input is a string, object, datetime, or a table-like ndarray 417 # (structured array, recarray). These input types can be 418 # uniquely identified by the format classes. 419 formats = [(name, cls) for name, cls in self.FORMATS.items() 420 if issubclass(cls, TimeUnique)] 421 422 # AstropyTime is a pseudo-format that isn't in the TIME_FORMATS registry, 423 # but try to guess it at the end. 424 formats.append(('astropy_time', TimeAstropyTime)) 425 426 elif not (isinstance(format, str) 427 and format.lower() in self.FORMATS): 428 if format is None: 429 raise ValueError("No time format was given, and the input is " 430 "not unique") 431 else: 432 raise ValueError("Format {!r} is not one of the allowed " 433 "formats {}".format(format, 434 sorted(self.FORMATS))) 435 else: 436 formats = [(format, self.FORMATS[format])] 437 438 assert formats 439 problems = {} 440 for name, cls in formats: 441 try: 442 return cls(val, val2, scale, precision, in_subfmt, out_subfmt) 443 except UnitConversionError: 444 raise 445 except (ValueError, TypeError) as err: 446 # If ``format`` specified then there is only one possibility, so raise 447 # immediately and include the upstream exception message to make it 448 # easier for user to see what is wrong. 449 if len(formats) == 1: 450 raise ValueError( 451 f'Input values did not match the format class {format}:' 452 + os.linesep 453 + f'{err.__class__.__name__}: {err}' 454 ) from err 455 else: 456 problems[name] = err 457 else: 458 raise ValueError(f'Input values did not match any of the formats ' 459 f'where the format keyword is optional: ' 460 f'{problems}') from problems[formats[0][0]] 461 462 @property 463 def writeable(self): 464 return self._time.jd1.flags.writeable & self._time.jd2.flags.writeable 465 466 @writeable.setter 467 def writeable(self, value): 468 self._time.jd1.flags.writeable = value 469 self._time.jd2.flags.writeable = value 470 471 @property 472 def format(self): 473 """ 474 Get or set time format. 475 476 The format defines the way times are represented when accessed via the 477 ``.value`` attribute. By default it is the same as the format used for 478 initializing the `Time` instance, but it can be set to any other value 479 that could be used for initialization. These can be listed with:: 480 481 >>> list(Time.FORMATS) 482 ['jd', 'mjd', 'decimalyear', 'unix', 'unix_tai', 'cxcsec', 'gps', 'plot_date', 483 'stardate', 'datetime', 'ymdhms', 'iso', 'isot', 'yday', 'datetime64', 484 'fits', 'byear', 'jyear', 'byear_str', 'jyear_str'] 485 """ 486 return self._format 487 488 @format.setter 489 def format(self, format): 490 """Set time format""" 491 if format not in self.FORMATS: 492 raise ValueError(f'format must be one of {list(self.FORMATS)}') 493 format_cls = self.FORMATS[format] 494 495 # Get the new TimeFormat object to contain time in new format. Possibly 496 # coerce in/out_subfmt to '*' (default) if existing subfmt values are 497 # not valid in the new format. 498 self._time = format_cls( 499 self._time.jd1, self._time.jd2, 500 self._time._scale, self.precision, 501 in_subfmt=format_cls._get_allowed_subfmt(self.in_subfmt), 502 out_subfmt=format_cls._get_allowed_subfmt(self.out_subfmt), 503 from_jd=True) 504 505 self._format = format 506 507 def __repr__(self): 508 return ("<{} object: scale='{}' format='{}' value={}>" 509 .format(self.__class__.__name__, self.scale, self.format, 510 getattr(self, self.format))) 511 512 def __str__(self): 513 return str(getattr(self, self.format)) 514 515 def __hash__(self): 516 517 try: 518 loc = getattr(self, 'location', None) 519 if loc is not None: 520 loc = loc.x.to_value(u.m), loc.y.to_value(u.m), loc.z.to_value(u.m) 521 522 return hash((self.jd1, self.jd2, self.scale, loc)) 523 524 except TypeError: 525 if self.ndim != 0: 526 reason = '(must be scalar)' 527 elif self.masked: 528 reason = '(value is masked)' 529 else: 530 raise 531 532 raise TypeError(f"unhashable type: '{self.__class__.__name__}' {reason}") 533 534 @property 535 def scale(self): 536 """Time scale""" 537 return self._time.scale 538 539 def _set_scale(self, scale): 540 """ 541 This is the key routine that actually does time scale conversions. 542 This is not public and not connected to the read-only scale property. 543 """ 544 545 if scale == self.scale: 546 return 547 if scale not in self.SCALES: 548 raise ValueError("Scale {!r} is not in the allowed scales {}" 549 .format(scale, sorted(self.SCALES))) 550 551 if scale == 'utc' or self.scale == 'utc': 552 # If doing a transform involving UTC then check that the leap 553 # seconds table is up to date. 554 _check_leapsec() 555 556 # Determine the chain of scale transformations to get from the current 557 # scale to the new scale. MULTI_HOPS contains a dict of all 558 # transformations (xforms) that require intermediate xforms. 559 # The MULTI_HOPS dict is keyed by (sys1, sys2) in alphabetical order. 560 xform = (self.scale, scale) 561 xform_sort = tuple(sorted(xform)) 562 multi = MULTI_HOPS.get(xform_sort, ()) 563 xforms = xform_sort[:1] + multi + xform_sort[-1:] 564 # If we made the reverse xform then reverse it now. 565 if xform_sort != xform: 566 xforms = tuple(reversed(xforms)) 567 568 # Transform the jd1,2 pairs through the chain of scale xforms. 569 jd1, jd2 = self._time.jd1, self._time.jd2_filled 570 for sys1, sys2 in zip(xforms[:-1], xforms[1:]): 571 # Some xforms require an additional delta_ argument that is 572 # provided through Time methods. These values may be supplied by 573 # the user or computed based on available approximations. The 574 # get_delta_ methods are available for only one combination of 575 # sys1, sys2 though the property applies for both xform directions. 576 args = [jd1, jd2] 577 for sys12 in ((sys1, sys2), (sys2, sys1)): 578 dt_method = '_get_delta_{}_{}'.format(*sys12) 579 try: 580 get_dt = getattr(self, dt_method) 581 except AttributeError: 582 pass 583 else: 584 args.append(get_dt(jd1, jd2)) 585 break 586 587 conv_func = getattr(erfa, sys1 + sys2) 588 jd1, jd2 = conv_func(*args) 589 590 jd1, jd2 = day_frac(jd1, jd2) 591 if self.masked: 592 jd2[self.mask] = np.nan 593 594 self._time = self.FORMATS[self.format](jd1, jd2, scale, self.precision, 595 self.in_subfmt, self.out_subfmt, 596 from_jd=True) 597 598 @property 599 def precision(self): 600 """ 601 Decimal precision when outputting seconds as floating point (int 602 value between 0 and 9 inclusive). 603 """ 604 return self._time.precision 605 606 @precision.setter 607 def precision(self, val): 608 del self.cache 609 if not isinstance(val, int) or val < 0 or val > 9: 610 raise ValueError('precision attribute must be an int between ' 611 '0 and 9') 612 self._time.precision = val 613 614 @property 615 def in_subfmt(self): 616 """ 617 Unix wildcard pattern to select subformats for parsing string input 618 times. 619 """ 620 return self._time.in_subfmt 621 622 @in_subfmt.setter 623 def in_subfmt(self, val): 624 self._time.in_subfmt = val 625 del self.cache 626 627 @property 628 def out_subfmt(self): 629 """ 630 Unix wildcard pattern to select subformats for outputting times. 631 """ 632 return self._time.out_subfmt 633 634 @out_subfmt.setter 635 def out_subfmt(self, val): 636 # Setting the out_subfmt property here does validation of ``val`` 637 self._time.out_subfmt = val 638 del self.cache 639 640 @property 641 def shape(self): 642 """The shape of the time instances. 643 644 Like `~numpy.ndarray.shape`, can be set to a new shape by assigning a 645 tuple. Note that if different instances share some but not all 646 underlying data, setting the shape of one instance can make the other 647 instance unusable. Hence, it is strongly recommended to get new, 648 reshaped instances with the ``reshape`` method. 649 650 Raises 651 ------ 652 ValueError 653 If the new shape has the wrong total number of elements. 654 AttributeError 655 If the shape of the ``jd1``, ``jd2``, ``location``, 656 ``delta_ut1_utc``, or ``delta_tdb_tt`` attributes cannot be changed 657 without the arrays being copied. For these cases, use the 658 `Time.reshape` method (which copies any arrays that cannot be 659 reshaped in-place). 660 """ 661 return self._time.jd1.shape 662 663 @shape.setter 664 def shape(self, shape): 665 del self.cache 666 667 # We have to keep track of arrays that were already reshaped, 668 # since we may have to return those to their original shape if a later 669 # shape-setting fails. 670 reshaped = [] 671 oldshape = self.shape 672 673 # In-place reshape of data/attributes. Need to access _time.jd1/2 not 674 # self.jd1/2 because the latter are not guaranteed to be the actual 675 # data, and in fact should not be directly changeable from the public 676 # API. 677 for obj, attr in ((self._time, 'jd1'), 678 (self._time, 'jd2'), 679 (self, '_delta_ut1_utc'), 680 (self, '_delta_tdb_tt'), 681 (self, 'location')): 682 val = getattr(obj, attr, None) 683 if val is not None and val.size > 1: 684 try: 685 val.shape = shape 686 except Exception: 687 for val2 in reshaped: 688 val2.shape = oldshape 689 raise 690 else: 691 reshaped.append(val) 692 693 def _shaped_like_input(self, value): 694 if self._time.jd1.shape: 695 if isinstance(value, np.ndarray): 696 return value 697 else: 698 raise TypeError( 699 f"JD is an array ({self._time.jd1!r}) but value " 700 f"is not ({value!r})") 701 else: 702 # zero-dimensional array, is it safe to unbox? 703 if (isinstance(value, np.ndarray) 704 and not value.shape 705 and not np.ma.is_masked(value)): 706 if value.dtype.kind == 'M': 707 # existing test doesn't want datetime64 converted 708 return value[()] 709 elif value.dtype.fields: 710 # Unpack but keep field names; .item() doesn't 711 # Still don't get python types in the fields 712 return value[()] 713 else: 714 return value.item() 715 else: 716 return value 717 718 @property 719 def jd1(self): 720 """ 721 First of the two doubles that internally store time value(s) in JD. 722 """ 723 jd1 = self._time.mask_if_needed(self._time.jd1) 724 return self._shaped_like_input(jd1) 725 726 @property 727 def jd2(self): 728 """ 729 Second of the two doubles that internally store time value(s) in JD. 730 """ 731 jd2 = self._time.mask_if_needed(self._time.jd2) 732 return self._shaped_like_input(jd2) 733 734 def to_value(self, format, subfmt='*'): 735 """Get time values expressed in specified output format. 736 737 This method allows representing the ``Time`` object in the desired 738 output ``format`` and optional sub-format ``subfmt``. Available 739 built-in formats include ``jd``, ``mjd``, ``iso``, and so forth. Each 740 format can have its own sub-formats 741 742 For built-in numerical formats like ``jd`` or ``unix``, ``subfmt`` can 743 be one of 'float', 'long', 'decimal', 'str', or 'bytes'. Here, 'long' 744 uses ``numpy.longdouble`` for somewhat enhanced precision (with 745 the enhancement depending on platform), and 'decimal' 746 :class:`decimal.Decimal` for full precision. For 'str' and 'bytes', the 747 number of digits is also chosen such that time values are represented 748 accurately. 749 750 For built-in date-like string formats, one of 'date_hms', 'date_hm', or 751 'date' (or 'longdate_hms', etc., for 5-digit years in 752 `~astropy.time.TimeFITS`). For sub-formats including seconds, the 753 number of digits used for the fractional seconds is as set by 754 `~astropy.time.Time.precision`. 755 756 Parameters 757 ---------- 758 format : str 759 The format in which one wants the time values. Default: the current 760 format. 761 subfmt : str or None, optional 762 Value or wildcard pattern to select the sub-format in which the 763 values should be given. The default of '*' picks the first 764 available for a given format, i.e., 'float' or 'date_hms'. 765 If `None`, use the instance's ``out_subfmt``. 766 767 """ 768 # TODO: add a precision argument (but ensure it is keyword argument 769 # only, to make life easier for TimeDelta.to_value()). 770 if format not in self.FORMATS: 771 raise ValueError(f'format must be one of {list(self.FORMATS)}') 772 773 cache = self.cache['format'] 774 # Try to keep cache behaviour like it was in astropy < 4.0. 775 key = format if subfmt is None else (format, subfmt) 776 if key not in cache: 777 if format == self.format: 778 tm = self 779 else: 780 tm = self.replicate(format=format) 781 782 # Some TimeFormat subclasses may not be able to handle being passes 783 # on a out_subfmt. This includes some core classes like 784 # TimeBesselianEpochString that do not have any allowed subfmts. But 785 # those do deal with `self.out_subfmt` internally, so if subfmt is 786 # the same, we do not pass it on. 787 kwargs = {} 788 if subfmt is not None and subfmt != tm.out_subfmt: 789 kwargs['out_subfmt'] = subfmt 790 try: 791 value = tm._time.to_value(parent=tm, **kwargs) 792 except TypeError as exc: 793 # Try validating subfmt, e.g. for formats like 'jyear_str' that 794 # do not implement out_subfmt in to_value() (because there are 795 # no allowed subformats). If subfmt is not valid this gives the 796 # same exception as would have occurred if the call to 797 # `to_value()` had succeeded. 798 tm._time._select_subfmts(subfmt) 799 800 # Subfmt was valid, so fall back to the original exception to see 801 # if it was lack of support for out_subfmt as a call arg. 802 if "unexpected keyword argument 'out_subfmt'" in str(exc): 803 raise ValueError( 804 f"to_value() method for format {format!r} does not " 805 f"support passing a 'subfmt' argument") from None 806 else: 807 # Some unforeseen exception so raise. 808 raise 809 810 value = tm._shaped_like_input(value) 811 cache[key] = value 812 return cache[key] 813 814 @property 815 def value(self): 816 """Time value(s) in current format""" 817 return self.to_value(self.format, None) 818 819 @property 820 def masked(self): 821 return self._time.masked 822 823 @property 824 def mask(self): 825 return self._time.mask 826 827 def insert(self, obj, values, axis=0): 828 """ 829 Insert values before the given indices in the column and return 830 a new `~astropy.time.Time` or `~astropy.time.TimeDelta` object. 831 832 The values to be inserted must conform to the rules for in-place setting 833 of ``Time`` objects (see ``Get and set values`` in the ``Time`` 834 documentation). 835 836 The API signature matches the ``np.insert`` API, but is more limited. 837 The specification of insert index ``obj`` must be a single integer, 838 and the ``axis`` must be ``0`` for simple row insertion before the 839 index. 840 841 Parameters 842 ---------- 843 obj : int 844 Integer index before which ``values`` is inserted. 845 values : array-like 846 Value(s) to insert. If the type of ``values`` is different 847 from that of quantity, ``values`` is converted to the matching type. 848 axis : int, optional 849 Axis along which to insert ``values``. Default is 0, which is the 850 only allowed value and will insert a row. 851 852 Returns 853 ------- 854 out : `~astropy.time.Time` subclass 855 New time object with inserted value(s) 856 857 """ 858 # Validate inputs: obj arg is integer, axis=0, self is not a scalar, and 859 # input index is in bounds. 860 try: 861 idx0 = operator.index(obj) 862 except TypeError: 863 raise TypeError('obj arg must be an integer') 864 865 if axis != 0: 866 raise ValueError('axis must be 0') 867 868 if not self.shape: 869 raise TypeError('cannot insert into scalar {} object' 870 .format(self.__class__.__name__)) 871 872 if abs(idx0) > len(self): 873 raise IndexError('index {} is out of bounds for axis 0 with size {}' 874 .format(idx0, len(self))) 875 876 # Turn negative index into positive 877 if idx0 < 0: 878 idx0 = len(self) + idx0 879 880 # For non-Time object, use numpy to help figure out the length. (Note annoying 881 # case of a string input that has a length which is not the length we want). 882 if not isinstance(values, self.__class__): 883 values = np.asarray(values) 884 n_values = len(values) if values.shape else 1 885 886 # Finally make the new object with the correct length and set values for the 887 # three sections, before insert, the insert, and after the insert. 888 out = self.__class__.info.new_like([self], len(self) + n_values, name=self.info.name) 889 890 out._time.jd1[:idx0] = self._time.jd1[:idx0] 891 out._time.jd2[:idx0] = self._time.jd2[:idx0] 892 893 # This uses the Time setting machinery to coerce and validate as necessary. 894 out[idx0:idx0 + n_values] = values 895 896 out._time.jd1[idx0 + n_values:] = self._time.jd1[idx0:] 897 out._time.jd2[idx0 + n_values:] = self._time.jd2[idx0:] 898 899 return out 900 901 def __setitem__(self, item, value): 902 if not self.writeable: 903 if self.shape: 904 raise ValueError('{} object is read-only. Make a ' 905 'copy() or set "writeable" attribute to True.' 906 .format(self.__class__.__name__)) 907 else: 908 raise ValueError('scalar {} object is read-only.' 909 .format(self.__class__.__name__)) 910 911 # Any use of setitem results in immediate cache invalidation 912 del self.cache 913 914 # Setting invalidates transform deltas 915 for attr in ('_delta_tdb_tt', '_delta_ut1_utc'): 916 if hasattr(self, attr): 917 delattr(self, attr) 918 919 if value is np.ma.masked or value is np.nan: 920 self._time.jd2[item] = np.nan 921 return 922 923 value = self._make_value_equivalent(item, value) 924 925 # Finally directly set the jd1/2 values. Locations are known to match. 926 if self.scale is not None: 927 value = getattr(value, self.scale) 928 self._time.jd1[item] = value._time.jd1 929 self._time.jd2[item] = value._time.jd2 930 931 def isclose(self, other, atol=None): 932 """Returns a boolean or boolean array where two Time objects are 933 element-wise equal within a time tolerance. 934 935 This evaluates the expression below:: 936 937 abs(self - other) <= atol 938 939 Parameters 940 ---------- 941 other : `~astropy.time.Time` 942 Time object for comparison. 943 atol : `~astropy.units.Quantity` or `~astropy.time.TimeDelta` 944 Absolute tolerance for equality with units of time (e.g. ``u.s`` or 945 ``u.day``). Default is two bits in the 128-bit JD time representation, 946 equivalent to about 40 picosecs. 947 """ 948 if atol is None: 949 # Note: use 2 bits instead of 1 bit based on experience in precision 950 # tests, since taking the difference with a UTC time means one has 951 # to do a scale change. 952 atol = 2 * np.finfo(float).eps * u.day 953 954 if not isinstance(atol, (u.Quantity, TimeDelta)): 955 raise TypeError("'atol' argument must be a Quantity or TimeDelta instance, got " 956 f'{atol.__class__.__name__} instead') 957 958 try: 959 # Separate these out so user sees where the problem is 960 dt = self - other 961 dt = abs(dt) 962 out = dt <= atol 963 except Exception as err: 964 raise TypeError("'other' argument must support subtraction with Time " 965 f"and return a value that supports comparison with " 966 f"{atol.__class__.__name__}: {err}") 967 968 return out 969 970 def copy(self, format=None): 971 """ 972 Return a fully independent copy the Time object, optionally changing 973 the format. 974 975 If ``format`` is supplied then the time format of the returned Time 976 object will be set accordingly, otherwise it will be unchanged from the 977 original. 978 979 In this method a full copy of the internal time arrays will be made. 980 The internal time arrays are normally not changeable by the user so in 981 most cases the ``replicate()`` method should be used. 982 983 Parameters 984 ---------- 985 format : str, optional 986 Time format of the copy. 987 988 Returns 989 ------- 990 tm : Time object 991 Copy of this object 992 """ 993 return self._apply('copy', format=format) 994 995 def replicate(self, format=None, copy=False, cls=None): 996 """ 997 Return a replica of the Time object, optionally changing the format. 998 999 If ``format`` is supplied then the time format of the returned Time 1000 object will be set accordingly, otherwise it will be unchanged from the 1001 original. 1002 1003 If ``copy`` is set to `True` then a full copy of the internal time arrays 1004 will be made. By default the replica will use a reference to the 1005 original arrays when possible to save memory. The internal time arrays 1006 are normally not changeable by the user so in most cases it should not 1007 be necessary to set ``copy`` to `True`. 1008 1009 The convenience method copy() is available in which ``copy`` is `True` 1010 by default. 1011 1012 Parameters 1013 ---------- 1014 format : str, optional 1015 Time format of the replica. 1016 copy : bool, optional 1017 Return a true copy instead of using references where possible. 1018 1019 Returns 1020 ------- 1021 tm : Time object 1022 Replica of this object 1023 """ 1024 return self._apply('copy' if copy else 'replicate', format=format, cls=cls) 1025 1026 def _apply(self, method, *args, format=None, cls=None, **kwargs): 1027 """Create a new time object, possibly applying a method to the arrays. 1028 1029 Parameters 1030 ---------- 1031 method : str or callable 1032 If string, can be 'replicate' or the name of a relevant 1033 `~numpy.ndarray` method. In the former case, a new time instance 1034 with unchanged internal data is created, while in the latter the 1035 method is applied to the internal ``jd1`` and ``jd2`` arrays, as 1036 well as to possible ``location``, ``_delta_ut1_utc``, and 1037 ``_delta_tdb_tt`` arrays. 1038 If a callable, it is directly applied to the above arrays. 1039 Examples: 'copy', '__getitem__', 'reshape', `~numpy.broadcast_to`. 1040 args : tuple 1041 Any positional arguments for ``method``. 1042 kwargs : dict 1043 Any keyword arguments for ``method``. If the ``format`` keyword 1044 argument is present, this will be used as the Time format of the 1045 replica. 1046 1047 Examples 1048 -------- 1049 Some ways this is used internally:: 1050 1051 copy : ``_apply('copy')`` 1052 replicate : ``_apply('replicate')`` 1053 reshape : ``_apply('reshape', new_shape)`` 1054 index or slice : ``_apply('__getitem__', item)`` 1055 broadcast : ``_apply(np.broadcast, shape=new_shape)`` 1056 """ 1057 new_format = self.format if format is None else format 1058 1059 if callable(method): 1060 apply_method = lambda array: method(array, *args, **kwargs) 1061 1062 else: 1063 if method == 'replicate': 1064 apply_method = None 1065 else: 1066 apply_method = operator.methodcaller(method, *args, **kwargs) 1067 1068 jd1, jd2 = self._time.jd1, self._time.jd2 1069 if apply_method: 1070 jd1 = apply_method(jd1) 1071 jd2 = apply_method(jd2) 1072 1073 # Get a new instance of our class and set its attributes directly. 1074 tm = super().__new__(cls or self.__class__) 1075 tm._time = TimeJD(jd1, jd2, self.scale, precision=0, 1076 in_subfmt='*', out_subfmt='*', from_jd=True) 1077 1078 # Optional ndarray attributes. 1079 for attr in ('_delta_ut1_utc', '_delta_tdb_tt', 'location'): 1080 try: 1081 val = getattr(self, attr) 1082 except AttributeError: 1083 continue 1084 1085 if apply_method: 1086 # Apply the method to any value arrays (though skip if there is 1087 # only an array scalar and the method would return a view, 1088 # since in that case nothing would change). 1089 if getattr(val, 'shape', ()): 1090 val = apply_method(val) 1091 elif method == 'copy' or method == 'flatten': 1092 # flatten should copy also for a single element array, but 1093 # we cannot use it directly for array scalars, since it 1094 # always returns a one-dimensional array. So, just copy. 1095 val = copy.copy(val) 1096 1097 setattr(tm, attr, val) 1098 1099 # Copy other 'info' attr only if it has actually been defined and the 1100 # time object is not a scalar (issue #10688). 1101 # See PR #3898 for further explanation and justification, along 1102 # with Quantity.__array_finalize__ 1103 if 'info' in self.__dict__: 1104 tm.info = self.info 1105 1106 # Make the new internal _time object corresponding to the format 1107 # in the copy. If the format is unchanged this process is lightweight 1108 # and does not create any new arrays. 1109 if new_format not in tm.FORMATS: 1110 raise ValueError(f'format must be one of {list(tm.FORMATS)}') 1111 1112 NewFormat = tm.FORMATS[new_format] 1113 1114 tm._time = NewFormat( 1115 tm._time.jd1, tm._time.jd2, 1116 tm._time._scale, 1117 precision=self.precision, 1118 in_subfmt=NewFormat._get_allowed_subfmt(self.in_subfmt), 1119 out_subfmt=NewFormat._get_allowed_subfmt(self.out_subfmt), 1120 from_jd=True) 1121 tm._format = new_format 1122 tm.SCALES = self.SCALES 1123 1124 return tm 1125 1126 def __copy__(self): 1127 """ 1128 Overrides the default behavior of the `copy.copy` function in 1129 the python stdlib to behave like `Time.copy`. Does *not* make a 1130 copy of the JD arrays - only copies by reference. 1131 """ 1132 return self.replicate() 1133 1134 def __deepcopy__(self, memo): 1135 """ 1136 Overrides the default behavior of the `copy.deepcopy` function 1137 in the python stdlib to behave like `Time.copy`. Does make a 1138 copy of the JD arrays. 1139 """ 1140 return self.copy() 1141 1142 def _advanced_index(self, indices, axis=None, keepdims=False): 1143 """Turn argmin, argmax output into an advanced index. 1144 1145 Argmin, argmax output contains indices along a given axis in an array 1146 shaped like the other dimensions. To use this to get values at the 1147 correct location, a list is constructed in which the other axes are 1148 indexed sequentially. For ``keepdims`` is ``True``, the net result is 1149 the same as constructing an index grid with ``np.ogrid`` and then 1150 replacing the ``axis`` item with ``indices`` with its shaped expanded 1151 at ``axis``. For ``keepdims`` is ``False``, the result is the same but 1152 with the ``axis`` dimension removed from all list entries. 1153 1154 For ``axis`` is ``None``, this calls :func:`~numpy.unravel_index`. 1155 1156 Parameters 1157 ---------- 1158 indices : array 1159 Output of argmin or argmax. 1160 axis : int or None 1161 axis along which argmin or argmax was used. 1162 keepdims : bool 1163 Whether to construct indices that keep or remove the axis along 1164 which argmin or argmax was used. Default: ``False``. 1165 1166 Returns 1167 ------- 1168 advanced_index : list of arrays 1169 Suitable for use as an advanced index. 1170 """ 1171 if axis is None: 1172 return np.unravel_index(indices, self.shape) 1173 1174 ndim = self.ndim 1175 if axis < 0: 1176 axis = axis + ndim 1177 1178 if keepdims and indices.ndim < self.ndim: 1179 indices = np.expand_dims(indices, axis) 1180 1181 index = [indices 1182 if i == axis 1183 else np.arange(s).reshape( 1184 (1,) * (i if keepdims or i < axis else i - 1) 1185 + (s,) 1186 + (1,) * (ndim - i - (1 if keepdims or i > axis else 2)) 1187 ) 1188 for i, s in enumerate(self.shape)] 1189 1190 return tuple(index) 1191 1192 def argmin(self, axis=None, out=None): 1193 """Return indices of the minimum values along the given axis. 1194 1195 This is similar to :meth:`~numpy.ndarray.argmin`, but adapted to ensure 1196 that the full precision given by the two doubles ``jd1`` and ``jd2`` 1197 is used. See :func:`~numpy.argmin` for detailed documentation. 1198 """ 1199 # First get the minimum at normal precision. 1200 jd1, jd2 = self.jd1, self.jd2 1201 approx = np.min(jd1 + jd2, axis, keepdims=True) 1202 1203 # Approx is very close to the true minimum, and by subtracting it at 1204 # full precision, all numbers near 0 can be represented correctly, 1205 # so we can be sure we get the true minimum. 1206 # The below is effectively what would be done for 1207 # dt = (self - self.__class__(approx, format='jd')).jd 1208 # which translates to: 1209 # approx_jd1, approx_jd2 = day_frac(approx, 0.) 1210 # dt = (self.jd1 - approx_jd1) + (self.jd2 - approx_jd2) 1211 dt = (jd1 - approx) + jd2 1212 1213 return dt.argmin(axis, out) 1214 1215 def argmax(self, axis=None, out=None): 1216 """Return indices of the maximum values along the given axis. 1217 1218 This is similar to :meth:`~numpy.ndarray.argmax`, but adapted to ensure 1219 that the full precision given by the two doubles ``jd1`` and ``jd2`` 1220 is used. See :func:`~numpy.argmax` for detailed documentation. 1221 """ 1222 # For procedure, see comment on argmin. 1223 jd1, jd2 = self.jd1, self.jd2 1224 approx = np.max(jd1 + jd2, axis, keepdims=True) 1225 1226 dt = (jd1 - approx) + jd2 1227 1228 return dt.argmax(axis, out) 1229 1230 def argsort(self, axis=-1): 1231 """Returns the indices that would sort the time array. 1232 1233 This is similar to :meth:`~numpy.ndarray.argsort`, but adapted to ensure 1234 that the full precision given by the two doubles ``jd1`` and ``jd2`` 1235 is used, and that corresponding attributes are copied. Internally, 1236 it uses :func:`~numpy.lexsort`, and hence no sort method can be chosen. 1237 """ 1238 # For procedure, see comment on argmin. 1239 jd1, jd2 = self.jd1, self.jd2 1240 approx = jd1 + jd2 1241 remainder = (jd1 - approx) + jd2 1242 1243 if axis is None: 1244 return np.lexsort((remainder.ravel(), approx.ravel())) 1245 else: 1246 return np.lexsort(keys=(remainder, approx), axis=axis) 1247 1248 def min(self, axis=None, out=None, keepdims=False): 1249 """Minimum along a given axis. 1250 1251 This is similar to :meth:`~numpy.ndarray.min`, but adapted to ensure 1252 that the full precision given by the two doubles ``jd1`` and ``jd2`` 1253 is used, and that corresponding attributes are copied. 1254 1255 Note that the ``out`` argument is present only for compatibility with 1256 ``np.min``; since `Time` instances are immutable, it is not possible 1257 to have an actual ``out`` to store the result in. 1258 """ 1259 if out is not None: 1260 raise ValueError("Since `Time` instances are immutable, ``out`` " 1261 "cannot be set to anything but ``None``.") 1262 return self[self._advanced_index(self.argmin(axis), axis, keepdims)] 1263 1264 def max(self, axis=None, out=None, keepdims=False): 1265 """Maximum along a given axis. 1266 1267 This is similar to :meth:`~numpy.ndarray.max`, but adapted to ensure 1268 that the full precision given by the two doubles ``jd1`` and ``jd2`` 1269 is used, and that corresponding attributes are copied. 1270 1271 Note that the ``out`` argument is present only for compatibility with 1272 ``np.max``; since `Time` instances are immutable, it is not possible 1273 to have an actual ``out`` to store the result in. 1274 """ 1275 if out is not None: 1276 raise ValueError("Since `Time` instances are immutable, ``out`` " 1277 "cannot be set to anything but ``None``.") 1278 return self[self._advanced_index(self.argmax(axis), axis, keepdims)] 1279 1280 def ptp(self, axis=None, out=None, keepdims=False): 1281 """Peak to peak (maximum - minimum) along a given axis. 1282 1283 This is similar to :meth:`~numpy.ndarray.ptp`, but adapted to ensure 1284 that the full precision given by the two doubles ``jd1`` and ``jd2`` 1285 is used. 1286 1287 Note that the ``out`` argument is present only for compatibility with 1288 `~numpy.ptp`; since `Time` instances are immutable, it is not possible 1289 to have an actual ``out`` to store the result in. 1290 """ 1291 if out is not None: 1292 raise ValueError("Since `Time` instances are immutable, ``out`` " 1293 "cannot be set to anything but ``None``.") 1294 return (self.max(axis, keepdims=keepdims) 1295 - self.min(axis, keepdims=keepdims)) 1296 1297 def sort(self, axis=-1): 1298 """Return a copy sorted along the specified axis. 1299 1300 This is similar to :meth:`~numpy.ndarray.sort`, but internally uses 1301 indexing with :func:`~numpy.lexsort` to ensure that the full precision 1302 given by the two doubles ``jd1`` and ``jd2`` is kept, and that 1303 corresponding attributes are properly sorted and copied as well. 1304 1305 Parameters 1306 ---------- 1307 axis : int or None 1308 Axis to be sorted. If ``None``, the flattened array is sorted. 1309 By default, sort over the last axis. 1310 """ 1311 return self[self._advanced_index(self.argsort(axis), axis, 1312 keepdims=True)] 1313 1314 @property 1315 def cache(self): 1316 """ 1317 Return the cache associated with this instance. 1318 """ 1319 return self._time.cache 1320 1321 @cache.deleter 1322 def cache(self): 1323 del self._time.cache 1324 1325 def __getattr__(self, attr): 1326 """ 1327 Get dynamic attributes to output format or do timescale conversion. 1328 """ 1329 if attr in self.SCALES and self.scale is not None: 1330 cache = self.cache['scale'] 1331 if attr not in cache: 1332 if attr == self.scale: 1333 tm = self 1334 else: 1335 tm = self.replicate() 1336 tm._set_scale(attr) 1337 if tm.shape: 1338 # Prevent future modification of cached array-like object 1339 tm.writeable = False 1340 cache[attr] = tm 1341 return cache[attr] 1342 1343 elif attr in self.FORMATS: 1344 return self.to_value(attr, subfmt=None) 1345 1346 elif attr in TIME_SCALES: # allowed ones done above (self.SCALES) 1347 if self.scale is None: 1348 raise ScaleValueError("Cannot convert TimeDelta with " 1349 "undefined scale to any defined scale.") 1350 else: 1351 raise ScaleValueError("Cannot convert {} with scale " 1352 "'{}' to scale '{}'" 1353 .format(self.__class__.__name__, 1354 self.scale, attr)) 1355 1356 else: 1357 # Should raise AttributeError 1358 return self.__getattribute__(attr) 1359 1360 @override__dir__ 1361 def __dir__(self): 1362 result = set(self.SCALES) 1363 result.update(self.FORMATS) 1364 return result 1365 1366 def _match_shape(self, val): 1367 """ 1368 Ensure that `val` is matched to length of self. If val has length 1 1369 then broadcast, otherwise cast to double and make sure shape matches. 1370 """ 1371 val = _make_array(val, copy=True) # be conservative and copy 1372 if val.size > 1 and val.shape != self.shape: 1373 try: 1374 # check the value can be broadcast to the shape of self. 1375 val = np.broadcast_to(val, self.shape, subok=True) 1376 except Exception: 1377 raise ValueError('Attribute shape must match or be ' 1378 'broadcastable to that of Time object. ' 1379 'Typically, give either a single value or ' 1380 'one for each time.') 1381 1382 return val 1383 1384 def _time_comparison(self, other, op): 1385 """If other is of same class as self, compare difference in self.scale. 1386 Otherwise, return NotImplemented 1387 """ 1388 if other.__class__ is not self.__class__: 1389 try: 1390 other = self.__class__(other, scale=self.scale) 1391 except Exception: 1392 # Let other have a go. 1393 return NotImplemented 1394 1395 if(self.scale is not None and self.scale not in other.SCALES 1396 or other.scale is not None and other.scale not in self.SCALES): 1397 # Other will also not be able to do it, so raise a TypeError 1398 # immediately, allowing us to explain why it doesn't work. 1399 raise TypeError("Cannot compare {} instances with scales " 1400 "'{}' and '{}'".format(self.__class__.__name__, 1401 self.scale, other.scale)) 1402 1403 if self.scale is not None and other.scale is not None: 1404 other = getattr(other, self.scale) 1405 1406 return op((self.jd1 - other.jd1) + (self.jd2 - other.jd2), 0.) 1407 1408 def __lt__(self, other): 1409 return self._time_comparison(other, operator.lt) 1410 1411 def __le__(self, other): 1412 return self._time_comparison(other, operator.le) 1413 1414 def __eq__(self, other): 1415 """ 1416 If other is an incompatible object for comparison, return `False`. 1417 Otherwise, return `True` if the time difference between self and 1418 other is zero. 1419 """ 1420 return self._time_comparison(other, operator.eq) 1421 1422 def __ne__(self, other): 1423 """ 1424 If other is an incompatible object for comparison, return `True`. 1425 Otherwise, return `False` if the time difference between self and 1426 other is zero. 1427 """ 1428 return self._time_comparison(other, operator.ne) 1429 1430 def __gt__(self, other): 1431 return self._time_comparison(other, operator.gt) 1432 1433 def __ge__(self, other): 1434 return self._time_comparison(other, operator.ge) 1435 1436 1437class Time(TimeBase): 1438 """ 1439 Represent and manipulate times and dates for astronomy. 1440 1441 A `Time` object is initialized with one or more times in the ``val`` 1442 argument. The input times in ``val`` must conform to the specified 1443 ``format`` and must correspond to the specified time ``scale``. The 1444 optional ``val2`` time input should be supplied only for numeric input 1445 formats (e.g. JD) where very high precision (better than 64-bit precision) 1446 is required. 1447 1448 The allowed values for ``format`` can be listed with:: 1449 1450 >>> list(Time.FORMATS) 1451 ['jd', 'mjd', 'decimalyear', 'unix', 'unix_tai', 'cxcsec', 'gps', 'plot_date', 1452 'stardate', 'datetime', 'ymdhms', 'iso', 'isot', 'yday', 'datetime64', 1453 'fits', 'byear', 'jyear', 'byear_str', 'jyear_str'] 1454 1455 See also: http://docs.astropy.org/en/stable/time/ 1456 1457 Parameters 1458 ---------- 1459 val : sequence, ndarray, number, str, bytes, or `~astropy.time.Time` object 1460 Value(s) to initialize the time or times. Bytes are decoded as ascii. 1461 val2 : sequence, ndarray, or number; optional 1462 Value(s) to initialize the time or times. Only used for numerical 1463 input, to help preserve precision. 1464 format : str, optional 1465 Format of input value(s) 1466 scale : str, optional 1467 Time scale of input value(s), must be one of the following: 1468 ('tai', 'tcb', 'tcg', 'tdb', 'tt', 'ut1', 'utc') 1469 precision : int, optional 1470 Digits of precision in string representation of time 1471 in_subfmt : str, optional 1472 Unix glob to select subformats for parsing input times 1473 out_subfmt : str, optional 1474 Unix glob to select subformat for outputting times 1475 location : `~astropy.coordinates.EarthLocation` or tuple, optional 1476 If given as an tuple, it should be able to initialize an 1477 an EarthLocation instance, i.e., either contain 3 items with units of 1478 length for geocentric coordinates, or contain a longitude, latitude, 1479 and an optional height for geodetic coordinates. 1480 Can be a single location, or one for each input time. 1481 If not given, assumed to be the center of the Earth for time scale 1482 transformations to and from the solar-system barycenter. 1483 copy : bool, optional 1484 Make a copy of the input values 1485 """ 1486 SCALES = TIME_SCALES 1487 """List of time scales""" 1488 1489 FORMATS = TIME_FORMATS 1490 """Dict of time formats""" 1491 1492 def __new__(cls, val, val2=None, format=None, scale=None, 1493 precision=None, in_subfmt=None, out_subfmt=None, 1494 location=None, copy=False): 1495 1496 if isinstance(val, Time): 1497 self = val.replicate(format=format, copy=copy, cls=cls) 1498 else: 1499 self = super().__new__(cls) 1500 1501 return self 1502 1503 def __init__(self, val, val2=None, format=None, scale=None, 1504 precision=None, in_subfmt=None, out_subfmt=None, 1505 location=None, copy=False): 1506 1507 if location is not None: 1508 from astropy.coordinates import EarthLocation 1509 if isinstance(location, EarthLocation): 1510 self.location = location 1511 else: 1512 self.location = EarthLocation(*location) 1513 if self.location.size == 1: 1514 self.location = self.location.squeeze() 1515 else: 1516 if not hasattr(self, 'location'): 1517 self.location = None 1518 1519 if isinstance(val, Time): 1520 # Update _time formatting parameters if explicitly specified 1521 if precision is not None: 1522 self._time.precision = precision 1523 if in_subfmt is not None: 1524 self._time.in_subfmt = in_subfmt 1525 if out_subfmt is not None: 1526 self._time.out_subfmt = out_subfmt 1527 self.SCALES = TIME_TYPES[self.scale] 1528 if scale is not None: 1529 self._set_scale(scale) 1530 else: 1531 self._init_from_vals(val, val2, format, scale, copy, 1532 precision, in_subfmt, out_subfmt) 1533 self.SCALES = TIME_TYPES[self.scale] 1534 1535 if self.location is not None and (self.location.size > 1 1536 and self.location.shape != self.shape): 1537 try: 1538 # check the location can be broadcast to self's shape. 1539 self.location = np.broadcast_to(self.location, self.shape, 1540 subok=True) 1541 except Exception as err: 1542 raise ValueError('The location with shape {} cannot be ' 1543 'broadcast against time with shape {}. ' 1544 'Typically, either give a single location or ' 1545 'one for each time.' 1546 .format(self.location.shape, self.shape)) from err 1547 1548 def _make_value_equivalent(self, item, value): 1549 """Coerce setitem value into an equivalent Time object""" 1550 1551 # If there is a vector location then broadcast to the Time shape 1552 # and then select with ``item`` 1553 if self.location is not None and self.location.shape: 1554 self_location = np.broadcast_to(self.location, self.shape, subok=True)[item] 1555 else: 1556 self_location = self.location 1557 1558 if isinstance(value, Time): 1559 # Make sure locations are compatible. Location can be either None or 1560 # a Location object. 1561 if self_location is None and value.location is None: 1562 match = True 1563 elif ((self_location is None and value.location is not None) 1564 or (self_location is not None and value.location is None)): 1565 match = False 1566 else: 1567 match = np.all(self_location == value.location) 1568 if not match: 1569 raise ValueError('cannot set to Time with different location: ' 1570 'expected location={} and ' 1571 'got location={}' 1572 .format(self_location, value.location)) 1573 else: 1574 try: 1575 value = self.__class__(value, scale=self.scale, location=self_location) 1576 except Exception: 1577 try: 1578 value = self.__class__(value, scale=self.scale, format=self.format, 1579 location=self_location) 1580 except Exception as err: 1581 raise ValueError('cannot convert value to a compatible Time object: {}' 1582 .format(err)) 1583 return value 1584 1585 @classmethod 1586 def now(cls): 1587 """ 1588 Creates a new object corresponding to the instant in time this 1589 method is called. 1590 1591 .. note:: 1592 "Now" is determined using the `~datetime.datetime.utcnow` 1593 function, so its accuracy and precision is determined by that 1594 function. Generally that means it is set by the accuracy of 1595 your system clock. 1596 1597 Returns 1598 ------- 1599 nowtime : :class:`~astropy.time.Time` 1600 A new `Time` object (or a subclass of `Time` if this is called from 1601 such a subclass) at the current time. 1602 """ 1603 # call `utcnow` immediately to be sure it's ASAP 1604 dtnow = datetime.utcnow() 1605 return cls(val=dtnow, format='datetime', scale='utc') 1606 1607 info = TimeInfo() 1608 1609 @classmethod 1610 def strptime(cls, time_string, format_string, **kwargs): 1611 """ 1612 Parse a string to a Time according to a format specification. 1613 See `time.strptime` documentation for format specification. 1614 1615 >>> Time.strptime('2012-Jun-30 23:59:60', '%Y-%b-%d %H:%M:%S') 1616 <Time object: scale='utc' format='isot' value=2012-06-30T23:59:60.000> 1617 1618 Parameters 1619 ---------- 1620 time_string : str, sequence, or ndarray 1621 Objects containing time data of type string 1622 format_string : str 1623 String specifying format of time_string. 1624 kwargs : dict 1625 Any keyword arguments for ``Time``. If the ``format`` keyword 1626 argument is present, this will be used as the Time format. 1627 1628 Returns 1629 ------- 1630 time_obj : `~astropy.time.Time` 1631 A new `~astropy.time.Time` object corresponding to the input 1632 ``time_string``. 1633 1634 """ 1635 time_array = np.asarray(time_string) 1636 1637 if time_array.dtype.kind not in ('U', 'S'): 1638 err = "Expected type is string, a bytes-like object or a sequence"\ 1639 " of these. Got dtype '{}'".format(time_array.dtype.kind) 1640 raise TypeError(err) 1641 1642 to_string = (str if time_array.dtype.kind == 'U' else 1643 lambda x: str(x.item(), encoding='ascii')) 1644 iterator = np.nditer([time_array, None], 1645 op_dtypes=[time_array.dtype, 'U30']) 1646 1647 for time, formatted in iterator: 1648 tt, fraction = _strptime._strptime(to_string(time), format_string) 1649 time_tuple = tt[:6] + (fraction,) 1650 formatted[...] = '{:04}-{:02}-{:02}T{:02}:{:02}:{:02}.{:06}'\ 1651 .format(*time_tuple) 1652 1653 format = kwargs.pop('format', None) 1654 out = cls(*iterator.operands[1:], format='isot', **kwargs) 1655 if format is not None: 1656 out.format = format 1657 1658 return out 1659 1660 def strftime(self, format_spec): 1661 """ 1662 Convert Time to a string or a numpy.array of strings according to a 1663 format specification. 1664 See `time.strftime` documentation for format specification. 1665 1666 Parameters 1667 ---------- 1668 format_spec : str 1669 Format definition of return string. 1670 1671 Returns 1672 ------- 1673 formatted : str or numpy.array 1674 String or numpy.array of strings formatted according to the given 1675 format string. 1676 1677 """ 1678 formatted_strings = [] 1679 for sk in self.replicate('iso')._time.str_kwargs(): 1680 date_tuple = date(sk['year'], sk['mon'], sk['day']).timetuple() 1681 datetime_tuple = (sk['year'], sk['mon'], sk['day'], 1682 sk['hour'], sk['min'], sk['sec'], 1683 date_tuple[6], date_tuple[7], -1) 1684 fmtd_str = format_spec 1685 if '%f' in fmtd_str: 1686 fmtd_str = fmtd_str.replace('%f', '{frac:0{precision}}'.format( 1687 frac=sk['fracsec'], precision=self.precision)) 1688 fmtd_str = strftime(fmtd_str, datetime_tuple) 1689 formatted_strings.append(fmtd_str) 1690 1691 if self.isscalar: 1692 return formatted_strings[0] 1693 else: 1694 return np.array(formatted_strings).reshape(self.shape) 1695 1696 def light_travel_time(self, skycoord, kind='barycentric', location=None, ephemeris=None): 1697 """Light travel time correction to the barycentre or heliocentre. 1698 1699 The frame transformations used to calculate the location of the solar 1700 system barycentre and the heliocentre rely on the erfa routine epv00, 1701 which is consistent with the JPL DE405 ephemeris to an accuracy of 1702 11.2 km, corresponding to a light travel time of 4 microseconds. 1703 1704 The routine assumes the source(s) are at large distance, i.e., neglects 1705 finite-distance effects. 1706 1707 Parameters 1708 ---------- 1709 skycoord : `~astropy.coordinates.SkyCoord` 1710 The sky location to calculate the correction for. 1711 kind : str, optional 1712 ``'barycentric'`` (default) or ``'heliocentric'`` 1713 location : `~astropy.coordinates.EarthLocation`, optional 1714 The location of the observatory to calculate the correction for. 1715 If no location is given, the ``location`` attribute of the Time 1716 object is used 1717 ephemeris : str, optional 1718 Solar system ephemeris to use (e.g., 'builtin', 'jpl'). By default, 1719 use the one set with ``astropy.coordinates.solar_system_ephemeris.set``. 1720 For more information, see `~astropy.coordinates.solar_system_ephemeris`. 1721 1722 Returns 1723 ------- 1724 time_offset : `~astropy.time.TimeDelta` 1725 The time offset between the barycentre or Heliocentre and Earth, 1726 in TDB seconds. Should be added to the original time to get the 1727 time in the Solar system barycentre or the Heliocentre. 1728 Also, the time conversion to BJD will then include the relativistic correction as well. 1729 """ 1730 1731 if kind.lower() not in ('barycentric', 'heliocentric'): 1732 raise ValueError("'kind' parameter must be one of 'heliocentric' " 1733 "or 'barycentric'") 1734 1735 if location is None: 1736 if self.location is None: 1737 raise ValueError('An EarthLocation needs to be set or passed ' 1738 'in to calculate bary- or heliocentric ' 1739 'corrections') 1740 location = self.location 1741 1742 from astropy.coordinates import (UnitSphericalRepresentation, CartesianRepresentation, 1743 HCRS, ICRS, GCRS, solar_system_ephemeris) 1744 1745 # ensure sky location is ICRS compatible 1746 if not skycoord.is_transformable_to(ICRS()): 1747 raise ValueError("Given skycoord is not transformable to the ICRS") 1748 1749 # get location of observatory in ITRS coordinates at this Time 1750 try: 1751 itrs = location.get_itrs(obstime=self) 1752 except Exception: 1753 raise ValueError("Supplied location does not have a valid `get_itrs` method") 1754 1755 with solar_system_ephemeris.set(ephemeris): 1756 if kind.lower() == 'heliocentric': 1757 # convert to heliocentric coordinates, aligned with ICRS 1758 cpos = itrs.transform_to(HCRS(obstime=self)).cartesian.xyz 1759 else: 1760 # first we need to convert to GCRS coordinates with the correct 1761 # obstime, since ICRS coordinates have no frame time 1762 gcrs_coo = itrs.transform_to(GCRS(obstime=self)) 1763 # convert to barycentric (BCRS) coordinates, aligned with ICRS 1764 cpos = gcrs_coo.transform_to(ICRS()).cartesian.xyz 1765 1766 # get unit ICRS vector to star 1767 spos = (skycoord.icrs.represent_as(UnitSphericalRepresentation). 1768 represent_as(CartesianRepresentation).xyz) 1769 1770 # Move X,Y,Z to last dimension, to enable possible broadcasting below. 1771 cpos = np.rollaxis(cpos, 0, cpos.ndim) 1772 spos = np.rollaxis(spos, 0, spos.ndim) 1773 1774 # calculate light travel time correction 1775 tcor_val = (spos * cpos).sum(axis=-1) / const.c 1776 return TimeDelta(tcor_val, scale='tdb') 1777 1778 def earth_rotation_angle(self, longitude=None): 1779 """Calculate local Earth rotation angle. 1780 1781 Parameters 1782 --------------- 1783 longitude : `~astropy.units.Quantity`, `~astropy.coordinates.EarthLocation`, str, or None; optional 1784 The longitude on the Earth at which to compute the Earth rotation 1785 angle (taken from a location as needed). If `None` (default), taken 1786 from the ``location`` attribute of the Time instance. If the special 1787 string 'tio', the result will be relative to the Terrestrial 1788 Intermediate Origin (TIO) (i.e., the output of `~erfa.era00`). 1789 1790 Returns 1791 ------- 1792 `~astropy.coordinates.Longitude` 1793 Local Earth rotation angle with units of hourangle. 1794 1795 See Also 1796 -------- 1797 astropy.time.Time.sidereal_time 1798 1799 References 1800 ---------- 1801 IAU 2006 NFA Glossary 1802 (currently located at: https://syrte.obspm.fr/iauWGnfa/NFA_Glossary.html) 1803 1804 Notes 1805 ----- 1806 The difference between apparent sidereal time and Earth rotation angle 1807 is the equation of the origins, which is the angle between the Celestial 1808 Intermediate Origin (CIO) and the equinox. Applying apparent sidereal 1809 time to the hour angle yields the true apparent Right Ascension with 1810 respect to the equinox, while applying the Earth rotation angle yields 1811 the intermediate (CIRS) Right Ascension with respect to the CIO. 1812 1813 The result includes the TIO locator (s'), which positions the Terrestrial 1814 Intermediate Origin on the equator of the Celestial Intermediate Pole (CIP) 1815 and is rigorously corrected for polar motion. 1816 (except when ``longitude='tio'``). 1817 1818 """ 1819 if isinstance(longitude, str) and longitude == 'tio': 1820 longitude = 0 1821 include_tio = False 1822 else: 1823 include_tio = True 1824 1825 return self._sid_time_or_earth_rot_ang(longitude=longitude, 1826 function=erfa.era00, scales=('ut1',), 1827 include_tio=include_tio) 1828 1829 def sidereal_time(self, kind, longitude=None, model=None): 1830 """Calculate sidereal time. 1831 1832 Parameters 1833 --------------- 1834 kind : str 1835 ``'mean'`` or ``'apparent'``, i.e., accounting for precession 1836 only, or also for nutation. 1837 longitude : `~astropy.units.Quantity`, `~astropy.coordinates.EarthLocation`, str, or None; optional 1838 The longitude on the Earth at which to compute the Earth rotation 1839 angle (taken from a location as needed). If `None` (default), taken 1840 from the ``location`` attribute of the Time instance. If the special 1841 string 'greenwich' or 'tio', the result will be relative to longitude 1842 0 for models before 2000, and relative to the Terrestrial Intermediate 1843 Origin (TIO) for later ones (i.e., the output of the relevant ERFA 1844 function that calculates greenwich sidereal time). 1845 model : str or None; optional 1846 Precession (and nutation) model to use. The available ones are: 1847 - {0}: {1} 1848 - {2}: {3} 1849 If `None` (default), the last (most recent) one from the appropriate 1850 list above is used. 1851 1852 Returns 1853 ------- 1854 `~astropy.coordinates.Longitude` 1855 Local sidereal time, with units of hourangle. 1856 1857 See Also 1858 -------- 1859 astropy.time.Time.earth_rotation_angle 1860 1861 References 1862 ---------- 1863 IAU 2006 NFA Glossary 1864 (currently located at: https://syrte.obspm.fr/iauWGnfa/NFA_Glossary.html) 1865 1866 Notes 1867 ----- 1868 The difference between apparent sidereal time and Earth rotation angle 1869 is the equation of the origins, which is the angle between the Celestial 1870 Intermediate Origin (CIO) and the equinox. Applying apparent sidereal 1871 time to the hour angle yields the true apparent Right Ascension with 1872 respect to the equinox, while applying the Earth rotation angle yields 1873 the intermediate (CIRS) Right Ascension with respect to the CIO. 1874 1875 For the IAU precession models from 2000 onwards, the result includes the 1876 TIO locator (s'), which positions the Terrestrial Intermediate Origin on 1877 the equator of the Celestial Intermediate Pole (CIP) and is rigorously 1878 corrected for polar motion (except when ``longitude='tio'`` or ``'greenwich'``). 1879 1880 """ # docstring is formatted below 1881 1882 if kind.lower() not in SIDEREAL_TIME_MODELS.keys(): 1883 raise ValueError('The kind of sidereal time has to be {}'.format( 1884 ' or '.join(sorted(SIDEREAL_TIME_MODELS.keys())))) 1885 1886 available_models = SIDEREAL_TIME_MODELS[kind.lower()] 1887 1888 if model is None: 1889 model = sorted(available_models.keys())[-1] 1890 elif model.upper() not in available_models: 1891 raise ValueError( 1892 'Model {} not implemented for {} sidereal time; ' 1893 'available models are {}' 1894 .format(model, kind, sorted(available_models.keys()))) 1895 1896 model_kwargs = available_models[model.upper()] 1897 1898 if isinstance(longitude, str) and longitude in ('tio', 'greenwich'): 1899 longitude = 0 1900 model_kwargs = model_kwargs.copy() 1901 model_kwargs['include_tio'] = False 1902 1903 return self._sid_time_or_earth_rot_ang(longitude=longitude, **model_kwargs) 1904 1905 if isinstance(sidereal_time.__doc__, str): 1906 sidereal_time.__doc__ = sidereal_time.__doc__.format( 1907 'apparent', sorted(SIDEREAL_TIME_MODELS['apparent'].keys()), 1908 'mean', sorted(SIDEREAL_TIME_MODELS['mean'].keys())) 1909 1910 def _sid_time_or_earth_rot_ang(self, longitude, function, scales, include_tio=True): 1911 """Calculate a local sidereal time or Earth rotation angle. 1912 1913 Parameters 1914 ---------- 1915 longitude : `~astropy.units.Quantity`, `~astropy.coordinates.EarthLocation`, str, or None; optional 1916 The longitude on the Earth at which to compute the Earth rotation 1917 angle (taken from a location as needed). If `None` (default), taken 1918 from the ``location`` attribute of the Time instance. 1919 function : callable 1920 The ERFA function to use. 1921 scales : tuple of str 1922 The time scales that the function requires on input. 1923 include_tio : bool, optional 1924 Whether to includes the TIO locator corrected for polar motion. 1925 Should be `False` for pre-2000 IAU models. Default: `True`. 1926 1927 Returns 1928 ------- 1929 `~astropy.coordinates.Longitude` 1930 Local sidereal time or Earth rotation angle, with units of hourangle. 1931 1932 """ 1933 from astropy.coordinates import Longitude, EarthLocation 1934 from astropy.coordinates.builtin_frames.utils import get_polar_motion 1935 from astropy.coordinates.matrix_utilities import rotation_matrix 1936 1937 if longitude is None: 1938 if self.location is None: 1939 raise ValueError('No longitude is given but the location for ' 1940 'the Time object is not set.') 1941 longitude = self.location.lon 1942 elif isinstance(longitude, EarthLocation): 1943 longitude = longitude.lon 1944 else: 1945 # Sanity check on input; default unit is degree. 1946 longitude = Longitude(longitude, u.degree, copy=False) 1947 1948 theta = self._call_erfa(function, scales) 1949 1950 if include_tio: 1951 # TODO: this duplicates part of coordinates.erfa_astrom.ErfaAstrom.apio; 1952 # maybe posisble to factor out to one or the other. 1953 sp = self._call_erfa(erfa.sp00, ('tt',)) 1954 xp, yp = get_polar_motion(self) 1955 # Form the rotation matrix, CIRS to apparent [HA,Dec]. 1956 r = (rotation_matrix(longitude, 'z') 1957 @ rotation_matrix(-yp, 'x', unit=u.radian) 1958 @ rotation_matrix(-xp, 'y', unit=u.radian) 1959 @ rotation_matrix(theta+sp, 'z', unit=u.radian)) 1960 # Solve for angle. 1961 angle = np.arctan2(r[..., 0, 1], r[..., 0, 0]) << u.radian 1962 1963 else: 1964 angle = longitude + (theta << u.radian) 1965 1966 return Longitude(angle, u.hourangle) 1967 1968 def _call_erfa(self, function, scales): 1969 # TODO: allow erfa functions to be used on Time with __array_ufunc__. 1970 erfa_parameters = [getattr(getattr(self, scale)._time, jd_part) 1971 for scale in scales 1972 for jd_part in ('jd1', 'jd2_filled')] 1973 1974 result = function(*erfa_parameters) 1975 1976 if self.masked: 1977 result[self.mask] = np.nan 1978 1979 return result 1980 1981 def get_delta_ut1_utc(self, iers_table=None, return_status=False): 1982 """Find UT1 - UTC differences by interpolating in IERS Table. 1983 1984 Parameters 1985 ---------- 1986 iers_table : `~astropy.utils.iers.IERS`, optional 1987 Table containing UT1-UTC differences from IERS Bulletins A 1988 and/or B. Default: `~astropy.utils.iers.earth_orientation_table` 1989 (which in turn defaults to the combined version provided by 1990 `~astropy.utils.iers.IERS_Auto`). 1991 return_status : bool 1992 Whether to return status values. If `False` (default), iers 1993 raises `IndexError` if any time is out of the range 1994 covered by the IERS table. 1995 1996 Returns 1997 ------- 1998 ut1_utc : float or float array 1999 UT1-UTC, interpolated in IERS Table 2000 status : int or int array 2001 Status values (if ``return_status=`True```):: 2002 ``astropy.utils.iers.FROM_IERS_B`` 2003 ``astropy.utils.iers.FROM_IERS_A`` 2004 ``astropy.utils.iers.FROM_IERS_A_PREDICTION`` 2005 ``astropy.utils.iers.TIME_BEFORE_IERS_RANGE`` 2006 ``astropy.utils.iers.TIME_BEYOND_IERS_RANGE`` 2007 2008 Notes 2009 ----- 2010 In normal usage, UT1-UTC differences are calculated automatically 2011 on the first instance ut1 is needed. 2012 2013 Examples 2014 -------- 2015 To check in code whether any times are before the IERS table range:: 2016 2017 >>> from astropy.utils.iers import TIME_BEFORE_IERS_RANGE 2018 >>> t = Time(['1961-01-01', '2000-01-01'], scale='utc') 2019 >>> delta, status = t.get_delta_ut1_utc(return_status=True) # doctest: +REMOTE_DATA 2020 >>> status == TIME_BEFORE_IERS_RANGE # doctest: +REMOTE_DATA 2021 array([ True, False]...) 2022 """ 2023 if iers_table is None: 2024 from astropy.utils.iers import earth_orientation_table 2025 iers_table = earth_orientation_table.get() 2026 2027 return iers_table.ut1_utc(self.utc, return_status=return_status) 2028 2029 # Property for ERFA DUT arg = UT1 - UTC 2030 def _get_delta_ut1_utc(self, jd1=None, jd2=None): 2031 """ 2032 Get ERFA DUT arg = UT1 - UTC. This getter takes optional jd1 and 2033 jd2 args because it gets called that way when converting time scales. 2034 If delta_ut1_utc is not yet set, this will interpolate them from the 2035 the IERS table. 2036 """ 2037 # Sec. 4.3.1: the arg DUT is the quantity delta_UT1 = UT1 - UTC in 2038 # seconds. It is obtained from tables published by the IERS. 2039 if not hasattr(self, '_delta_ut1_utc'): 2040 from astropy.utils.iers import earth_orientation_table 2041 iers_table = earth_orientation_table.get() 2042 # jd1, jd2 are normally set (see above), except if delta_ut1_utc 2043 # is access directly; ensure we behave as expected for that case 2044 if jd1 is None: 2045 self_utc = self.utc 2046 jd1, jd2 = self_utc._time.jd1, self_utc._time.jd2_filled 2047 scale = 'utc' 2048 else: 2049 scale = self.scale 2050 # interpolate UT1-UTC in IERS table 2051 delta = iers_table.ut1_utc(jd1, jd2) 2052 # if we interpolated using UT1 jds, we may be off by one 2053 # second near leap seconds (and very slightly off elsewhere) 2054 if scale == 'ut1': 2055 # calculate UTC using the offset we got; the ERFA routine 2056 # is tolerant of leap seconds, so will do this right 2057 jd1_utc, jd2_utc = erfa.ut1utc(jd1, jd2, delta.to_value(u.s)) 2058 # calculate a better estimate using the nearly correct UTC 2059 delta = iers_table.ut1_utc(jd1_utc, jd2_utc) 2060 2061 self._set_delta_ut1_utc(delta) 2062 2063 return self._delta_ut1_utc 2064 2065 def _set_delta_ut1_utc(self, val): 2066 del self.cache 2067 if hasattr(val, 'to'): # Matches Quantity but also TimeDelta. 2068 val = val.to(u.second).value 2069 val = self._match_shape(val) 2070 self._delta_ut1_utc = val 2071 2072 # Note can't use @property because _get_delta_tdb_tt is explicitly 2073 # called with the optional jd1 and jd2 args. 2074 delta_ut1_utc = property(_get_delta_ut1_utc, _set_delta_ut1_utc) 2075 """UT1 - UTC time scale offset""" 2076 2077 # Property for ERFA DTR arg = TDB - TT 2078 def _get_delta_tdb_tt(self, jd1=None, jd2=None): 2079 if not hasattr(self, '_delta_tdb_tt'): 2080 # If jd1 and jd2 are not provided (which is the case for property 2081 # attribute access) then require that the time scale is TT or TDB. 2082 # Otherwise the computations here are not correct. 2083 if jd1 is None or jd2 is None: 2084 if self.scale not in ('tt', 'tdb'): 2085 raise ValueError('Accessing the delta_tdb_tt attribute ' 2086 'is only possible for TT or TDB time ' 2087 'scales') 2088 else: 2089 jd1 = self._time.jd1 2090 jd2 = self._time.jd2_filled 2091 2092 # First go from the current input time (which is either 2093 # TDB or TT) to an approximate UT1. Since TT and TDB are 2094 # pretty close (few msec?), assume TT. Similarly, since the 2095 # UT1 terms are very small, use UTC instead of UT1. 2096 njd1, njd2 = erfa.tttai(jd1, jd2) 2097 njd1, njd2 = erfa.taiutc(njd1, njd2) 2098 # subtract 0.5, so UT is fraction of the day from midnight 2099 ut = day_frac(njd1 - 0.5, njd2)[1] 2100 2101 if self.location is None: 2102 # Assume geocentric. 2103 self._delta_tdb_tt = erfa.dtdb(jd1, jd2, ut, 0., 0., 0.) 2104 else: 2105 location = self.location 2106 # Geodetic params needed for d_tdb_tt() 2107 lon = location.lon 2108 rxy = np.hypot(location.x, location.y) 2109 z = location.z 2110 self._delta_tdb_tt = erfa.dtdb( 2111 jd1, jd2, ut, lon.to_value(u.radian), 2112 rxy.to_value(u.km), z.to_value(u.km)) 2113 2114 return self._delta_tdb_tt 2115 2116 def _set_delta_tdb_tt(self, val): 2117 del self.cache 2118 if hasattr(val, 'to'): # Matches Quantity but also TimeDelta. 2119 val = val.to(u.second).value 2120 val = self._match_shape(val) 2121 self._delta_tdb_tt = val 2122 2123 # Note can't use @property because _get_delta_tdb_tt is explicitly 2124 # called with the optional jd1 and jd2 args. 2125 delta_tdb_tt = property(_get_delta_tdb_tt, _set_delta_tdb_tt) 2126 """TDB - TT time scale offset""" 2127 2128 def __sub__(self, other): 2129 # T - Tdelta = T 2130 # T - T = Tdelta 2131 other_is_delta = not isinstance(other, Time) 2132 if other_is_delta: # T - Tdelta 2133 # Check other is really a TimeDelta or something that can initialize. 2134 if not isinstance(other, TimeDelta): 2135 try: 2136 other = TimeDelta(other) 2137 except Exception: 2138 return NotImplemented 2139 2140 # we need a constant scale to calculate, which is guaranteed for 2141 # TimeDelta, but not for Time (which can be UTC) 2142 out = self.replicate() 2143 if self.scale in other.SCALES: 2144 if other.scale not in (out.scale, None): 2145 other = getattr(other, out.scale) 2146 else: 2147 if other.scale is None: 2148 out._set_scale('tai') 2149 else: 2150 if self.scale not in TIME_TYPES[other.scale]: 2151 raise TypeError("Cannot subtract Time and TimeDelta instances " 2152 "with scales '{}' and '{}'" 2153 .format(self.scale, other.scale)) 2154 out._set_scale(other.scale) 2155 # remove attributes that are invalidated by changing time 2156 for attr in ('_delta_ut1_utc', '_delta_tdb_tt'): 2157 if hasattr(out, attr): 2158 delattr(out, attr) 2159 2160 else: # T - T 2161 # the scales should be compatible (e.g., cannot convert TDB to LOCAL) 2162 if other.scale not in self.SCALES: 2163 raise TypeError("Cannot subtract Time instances " 2164 "with scales '{}' and '{}'" 2165 .format(self.scale, other.scale)) 2166 self_time = (self._time if self.scale in TIME_DELTA_SCALES 2167 else self.tai._time) 2168 # set up TimeDelta, subtraction to be done shortly 2169 out = TimeDelta(self_time.jd1, self_time.jd2, format='jd', 2170 scale=self_time.scale) 2171 2172 if other.scale != out.scale: 2173 other = getattr(other, out.scale) 2174 2175 jd1 = out._time.jd1 - other._time.jd1 2176 jd2 = out._time.jd2 - other._time.jd2 2177 2178 out._time.jd1, out._time.jd2 = day_frac(jd1, jd2) 2179 2180 if other_is_delta: 2181 # Go back to left-side scale if needed 2182 out._set_scale(self.scale) 2183 2184 return out 2185 2186 def __add__(self, other): 2187 # T + Tdelta = T 2188 # T + T = error 2189 if isinstance(other, Time): 2190 raise OperandTypeError(self, other, '+') 2191 2192 # Check other is really a TimeDelta or something that can initialize. 2193 if not isinstance(other, TimeDelta): 2194 try: 2195 other = TimeDelta(other) 2196 except Exception: 2197 return NotImplemented 2198 2199 # ideally, we calculate in the scale of the Time item, since that is 2200 # what we want the output in, but this may not be possible, since 2201 # TimeDelta cannot be converted arbitrarily 2202 out = self.replicate() 2203 if self.scale in other.SCALES: 2204 if other.scale not in (out.scale, None): 2205 other = getattr(other, out.scale) 2206 else: 2207 if other.scale is None: 2208 out._set_scale('tai') 2209 else: 2210 if self.scale not in TIME_TYPES[other.scale]: 2211 raise TypeError("Cannot add Time and TimeDelta instances " 2212 "with scales '{}' and '{}'" 2213 .format(self.scale, other.scale)) 2214 out._set_scale(other.scale) 2215 # remove attributes that are invalidated by changing time 2216 for attr in ('_delta_ut1_utc', '_delta_tdb_tt'): 2217 if hasattr(out, attr): 2218 delattr(out, attr) 2219 2220 jd1 = out._time.jd1 + other._time.jd1 2221 jd2 = out._time.jd2 + other._time.jd2 2222 2223 out._time.jd1, out._time.jd2 = day_frac(jd1, jd2) 2224 2225 # Go back to left-side scale if needed 2226 out._set_scale(self.scale) 2227 2228 return out 2229 2230 # Reverse addition is possible: <something-Tdelta-ish> + T 2231 # but there is no case of <something> - T, so no __rsub__. 2232 def __radd__(self, other): 2233 return self.__add__(other) 2234 2235 def to_datetime(self, timezone=None): 2236 # TODO: this could likely go through to_value, as long as that 2237 # had an **kwargs part that was just passed on to _time. 2238 tm = self.replicate(format='datetime') 2239 return tm._shaped_like_input(tm._time.to_value(timezone)) 2240 2241 to_datetime.__doc__ = TimeDatetime.to_value.__doc__ 2242 2243 2244class TimeDelta(TimeBase): 2245 """ 2246 Represent the time difference between two times. 2247 2248 A TimeDelta object is initialized with one or more times in the ``val`` 2249 argument. The input times in ``val`` must conform to the specified 2250 ``format``. The optional ``val2`` time input should be supplied only for 2251 numeric input formats (e.g. JD) where very high precision (better than 2252 64-bit precision) is required. 2253 2254 The allowed values for ``format`` can be listed with:: 2255 2256 >>> list(TimeDelta.FORMATS) 2257 ['sec', 'jd', 'datetime'] 2258 2259 Note that for time differences, the scale can be among three groups: 2260 geocentric ('tai', 'tt', 'tcg'), barycentric ('tcb', 'tdb'), and rotational 2261 ('ut1'). Within each of these, the scales for time differences are the 2262 same. Conversion between geocentric and barycentric is possible, as there 2263 is only a scale factor change, but one cannot convert to or from 'ut1', as 2264 this requires knowledge of the actual times, not just their difference. For 2265 a similar reason, 'utc' is not a valid scale for a time difference: a UTC 2266 day is not always 86400 seconds. 2267 2268 See also: 2269 2270 - https://docs.astropy.org/en/stable/time/ 2271 - https://docs.astropy.org/en/stable/time/index.html#time-deltas 2272 2273 Parameters 2274 ---------- 2275 val : sequence, ndarray, number, `~astropy.units.Quantity` or `~astropy.time.TimeDelta` object 2276 Value(s) to initialize the time difference(s). Any quantities will 2277 be converted appropriately (with care taken to avoid rounding 2278 errors for regular time units). 2279 val2 : sequence, ndarray, number, or `~astropy.units.Quantity`; optional 2280 Additional values, as needed to preserve precision. 2281 format : str, optional 2282 Format of input value(s) 2283 scale : str, optional 2284 Time scale of input value(s), must be one of the following values: 2285 ('tdb', 'tt', 'ut1', 'tcg', 'tcb', 'tai'). If not given (or 2286 ``None``), the scale is arbitrary; when added or subtracted from a 2287 ``Time`` instance, it will be used without conversion. 2288 copy : bool, optional 2289 Make a copy of the input values 2290 """ 2291 SCALES = TIME_DELTA_SCALES 2292 """List of time delta scales.""" 2293 2294 FORMATS = TIME_DELTA_FORMATS 2295 """Dict of time delta formats.""" 2296 2297 info = TimeDeltaInfo() 2298 2299 def __new__(cls, val, val2=None, format=None, scale=None, 2300 precision=None, in_subfmt=None, out_subfmt=None, 2301 location=None, copy=False): 2302 2303 if isinstance(val, TimeDelta): 2304 self = val.replicate(format=format, copy=copy, cls=cls) 2305 else: 2306 self = super().__new__(cls) 2307 2308 return self 2309 2310 def __init__(self, val, val2=None, format=None, scale=None, copy=False): 2311 if isinstance(val, TimeDelta): 2312 if scale is not None: 2313 self._set_scale(scale) 2314 else: 2315 if format is None: 2316 format = 'datetime' if isinstance(val, timedelta) else 'jd' 2317 2318 self._init_from_vals(val, val2, format, scale, copy) 2319 2320 if scale is not None: 2321 self.SCALES = TIME_DELTA_TYPES[scale] 2322 2323 def replicate(self, *args, **kwargs): 2324 out = super().replicate(*args, **kwargs) 2325 out.SCALES = self.SCALES 2326 return out 2327 2328 def to_datetime(self): 2329 """ 2330 Convert to ``datetime.timedelta`` object. 2331 """ 2332 tm = self.replicate(format='datetime') 2333 return tm._shaped_like_input(tm._time.value) 2334 2335 def _set_scale(self, scale): 2336 """ 2337 This is the key routine that actually does time scale conversions. 2338 This is not public and not connected to the read-only scale property. 2339 """ 2340 2341 if scale == self.scale: 2342 return 2343 if scale not in self.SCALES: 2344 raise ValueError("Scale {!r} is not in the allowed scales {}" 2345 .format(scale, sorted(self.SCALES))) 2346 2347 # For TimeDelta, there can only be a change in scale factor, 2348 # which is written as time2 - time1 = scale_offset * time1 2349 scale_offset = SCALE_OFFSETS[(self.scale, scale)] 2350 if scale_offset is None: 2351 self._time.scale = scale 2352 else: 2353 jd1, jd2 = self._time.jd1, self._time.jd2 2354 offset1, offset2 = day_frac(jd1, jd2, factor=scale_offset) 2355 self._time = self.FORMATS[self.format]( 2356 jd1 + offset1, jd2 + offset2, scale, 2357 self.precision, self.in_subfmt, 2358 self.out_subfmt, from_jd=True) 2359 2360 def _add_sub(self, other, op): 2361 """Perform common elements of addition / subtraction for two delta times""" 2362 # If not a TimeDelta then see if it can be turned into a TimeDelta. 2363 if not isinstance(other, TimeDelta): 2364 try: 2365 other = TimeDelta(other) 2366 except Exception: 2367 return NotImplemented 2368 2369 # the scales should be compatible (e.g., cannot convert TDB to TAI) 2370 if(self.scale is not None and self.scale not in other.SCALES 2371 or other.scale is not None and other.scale not in self.SCALES): 2372 raise TypeError("Cannot add TimeDelta instances with scales " 2373 "'{}' and '{}'".format(self.scale, other.scale)) 2374 2375 # adjust the scale of other if the scale of self is set (or no scales) 2376 if self.scale is not None or other.scale is None: 2377 out = self.replicate() 2378 if other.scale is not None: 2379 other = getattr(other, self.scale) 2380 else: 2381 out = other.replicate() 2382 2383 jd1 = op(self._time.jd1, other._time.jd1) 2384 jd2 = op(self._time.jd2, other._time.jd2) 2385 2386 out._time.jd1, out._time.jd2 = day_frac(jd1, jd2) 2387 2388 return out 2389 2390 def __add__(self, other): 2391 # If other is a Time then use Time.__add__ to do the calculation. 2392 if isinstance(other, Time): 2393 return other.__add__(self) 2394 2395 return self._add_sub(other, operator.add) 2396 2397 def __sub__(self, other): 2398 # TimeDelta - Time is an error 2399 if isinstance(other, Time): 2400 raise OperandTypeError(self, other, '-') 2401 2402 return self._add_sub(other, operator.sub) 2403 2404 def __radd__(self, other): 2405 return self.__add__(other) 2406 2407 def __rsub__(self, other): 2408 out = self.__sub__(other) 2409 return -out 2410 2411 def __neg__(self): 2412 """Negation of a `TimeDelta` object.""" 2413 new = self.copy() 2414 new._time.jd1 = -self._time.jd1 2415 new._time.jd2 = -self._time.jd2 2416 return new 2417 2418 def __abs__(self): 2419 """Absolute value of a `TimeDelta` object.""" 2420 jd1, jd2 = self._time.jd1, self._time.jd2 2421 negative = jd1 + jd2 < 0 2422 new = self.copy() 2423 new._time.jd1 = np.where(negative, -jd1, jd1) 2424 new._time.jd2 = np.where(negative, -jd2, jd2) 2425 return new 2426 2427 def __mul__(self, other): 2428 """Multiplication of `TimeDelta` objects by numbers/arrays.""" 2429 # Check needed since otherwise the self.jd1 * other multiplication 2430 # would enter here again (via __rmul__) 2431 if isinstance(other, Time): 2432 raise OperandTypeError(self, other, '*') 2433 elif ((isinstance(other, u.UnitBase) 2434 and other == u.dimensionless_unscaled) 2435 or (isinstance(other, str) and other == '')): 2436 return self.copy() 2437 2438 # If other is something consistent with a dimensionless quantity 2439 # (could just be a float or an array), then we can just multiple in. 2440 try: 2441 other = u.Quantity(other, u.dimensionless_unscaled, copy=False) 2442 except Exception: 2443 # If not consistent with a dimensionless quantity, try downgrading 2444 # self to a quantity and see if things work. 2445 try: 2446 return self.to(u.day) * other 2447 except Exception: 2448 # The various ways we could multiply all failed; 2449 # returning NotImplemented to give other a final chance. 2450 return NotImplemented 2451 2452 jd1, jd2 = day_frac(self.jd1, self.jd2, factor=other.value) 2453 out = TimeDelta(jd1, jd2, format='jd', scale=self.scale) 2454 2455 if self.format != 'jd': 2456 out = out.replicate(format=self.format) 2457 return out 2458 2459 def __rmul__(self, other): 2460 """Multiplication of numbers/arrays with `TimeDelta` objects.""" 2461 return self.__mul__(other) 2462 2463 def __truediv__(self, other): 2464 """Division of `TimeDelta` objects by numbers/arrays.""" 2465 # Cannot do __mul__(1./other) as that looses precision 2466 if ((isinstance(other, u.UnitBase) 2467 and other == u.dimensionless_unscaled) 2468 or (isinstance(other, str) and other == '')): 2469 return self.copy() 2470 2471 # If other is something consistent with a dimensionless quantity 2472 # (could just be a float or an array), then we can just divide in. 2473 try: 2474 other = u.Quantity(other, u.dimensionless_unscaled, copy=False) 2475 except Exception: 2476 # If not consistent with a dimensionless quantity, try downgrading 2477 # self to a quantity and see if things work. 2478 try: 2479 return self.to(u.day) / other 2480 except Exception: 2481 # The various ways we could divide all failed; 2482 # returning NotImplemented to give other a final chance. 2483 return NotImplemented 2484 2485 jd1, jd2 = day_frac(self.jd1, self.jd2, divisor=other.value) 2486 out = TimeDelta(jd1, jd2, format='jd', scale=self.scale) 2487 2488 if self.format != 'jd': 2489 out = out.replicate(format=self.format) 2490 return out 2491 2492 def __rtruediv__(self, other): 2493 """Division by `TimeDelta` objects of numbers/arrays.""" 2494 # Here, we do not have to worry about returning NotImplemented, 2495 # since other has already had a chance to look at us. 2496 return other / self.to(u.day) 2497 2498 def to(self, unit, equivalencies=[]): 2499 """ 2500 Convert to a quantity in the specified unit. 2501 2502 Parameters 2503 ---------- 2504 unit : unit-like 2505 The unit to convert to. 2506 equivalencies : list of tuple 2507 A list of equivalence pairs to try if the units are not directly 2508 convertible (see :ref:`astropy:unit_equivalencies`). If `None`, no 2509 equivalencies will be applied at all, not even any set globallyq 2510 or within a context. 2511 2512 Returns 2513 ------- 2514 quantity : `~astropy.units.Quantity` 2515 The quantity in the units specified. 2516 2517 See also 2518 -------- 2519 to_value : get the numerical value in a given unit. 2520 """ 2521 return u.Quantity(self._time.jd1 + self._time.jd2, 2522 u.day).to(unit, equivalencies=equivalencies) 2523 2524 def to_value(self, *args, **kwargs): 2525 """Get time delta values expressed in specified output format or unit. 2526 2527 This method is flexible and handles both conversion to a specified 2528 ``TimeDelta`` format / sub-format AND conversion to a specified unit. 2529 If positional argument(s) are provided then the first one is checked 2530 to see if it is a valid ``TimeDelta`` format, and next it is checked 2531 to see if it is a valid unit or unit string. 2532 2533 To convert to a ``TimeDelta`` format and optional sub-format the options 2534 are:: 2535 2536 tm = TimeDelta(1.0 * u.s) 2537 tm.to_value('jd') # equivalent of tm.jd 2538 tm.to_value('jd', 'decimal') # convert to 'jd' as a Decimal object 2539 tm.to_value('jd', subfmt='decimal') 2540 tm.to_value(format='jd', subfmt='decimal') 2541 2542 To convert to a unit with optional equivalencies, the options are:: 2543 2544 tm.to_value('hr') # convert to u.hr (hours) 2545 tm.to_value('hr', []) # specify equivalencies as a positional arg 2546 tm.to_value('hr', equivalencies=[]) 2547 tm.to_value(unit='hr', equivalencies=[]) 2548 2549 The built-in `~astropy.time.TimeDelta` options for ``format`` are: 2550 {'jd', 'sec', 'datetime'}. 2551 2552 For the two numerical formats 'jd' and 'sec', the available ``subfmt`` 2553 options are: {'float', 'long', 'decimal', 'str', 'bytes'}. Here, 'long' 2554 uses ``numpy.longdouble`` for somewhat enhanced precision (with the 2555 enhancement depending on platform), and 'decimal' instances of 2556 :class:`decimal.Decimal` for full precision. For the 'str' and 'bytes' 2557 sub-formats, the number of digits is also chosen such that time values 2558 are represented accurately. Default: as set by ``out_subfmt`` (which by 2559 default picks the first available for a given format, i.e., 'float'). 2560 2561 Parameters 2562 ---------- 2563 format : str, optional 2564 The format in which one wants the `~astropy.time.TimeDelta` values. 2565 Default: the current format. 2566 subfmt : str, optional 2567 Possible sub-format in which the values should be given. Default: as 2568 set by ``out_subfmt`` (which by default picks the first available 2569 for a given format, i.e., 'float' or 'date_hms'). 2570 unit : `~astropy.units.UnitBase` instance or str, optional 2571 The unit in which the value should be given. 2572 equivalencies : list of tuple 2573 A list of equivalence pairs to try if the units are not directly 2574 convertible (see :ref:`astropy:unit_equivalencies`). If `None`, no 2575 equivalencies will be applied at all, not even any set globally or 2576 within a context. 2577 2578 Returns 2579 ------- 2580 value : ndarray or scalar 2581 The value in the format or units specified. 2582 2583 See also 2584 -------- 2585 to : Convert to a `~astropy.units.Quantity` instance in a given unit. 2586 value : The time value in the current format. 2587 2588 """ 2589 if not (args or kwargs): 2590 raise TypeError('to_value() missing required format or unit argument') 2591 2592 # TODO: maybe allow 'subfmt' also for units, keeping full precision 2593 # (effectively, by doing the reverse of quantity_day_frac)? 2594 # This way, only equivalencies could lead to possible precision loss. 2595 if ('format' in kwargs 2596 or (args != () and (args[0] is None or args[0] in self.FORMATS))): 2597 # Super-class will error with duplicate arguments, etc. 2598 return super().to_value(*args, **kwargs) 2599 2600 # With positional arguments, we try parsing the first one as a unit, 2601 # so that on failure we can give a more informative exception. 2602 if args: 2603 try: 2604 unit = u.Unit(args[0]) 2605 except ValueError as exc: 2606 raise ValueError("first argument is not one of the known " 2607 "formats ({}) and failed to parse as a unit." 2608 .format(list(self.FORMATS))) from exc 2609 args = (unit,) + args[1:] 2610 2611 return u.Quantity(self._time.jd1 + self._time.jd2, 2612 u.day).to_value(*args, **kwargs) 2613 2614 def _make_value_equivalent(self, item, value): 2615 """Coerce setitem value into an equivalent TimeDelta object""" 2616 if not isinstance(value, TimeDelta): 2617 try: 2618 value = self.__class__(value, scale=self.scale, format=self.format) 2619 except Exception as err: 2620 raise ValueError('cannot convert value to a compatible TimeDelta ' 2621 'object: {}'.format(err)) 2622 return value 2623 2624 def isclose(self, other, atol=None, rtol=0.0): 2625 """Returns a boolean or boolean array where two TimeDelta objects are 2626 element-wise equal within a time tolerance. 2627 2628 This effectively evaluates the expression below:: 2629 2630 abs(self - other) <= atol + rtol * abs(other) 2631 2632 Parameters 2633 ---------- 2634 other : `~astropy.units.Quantity` or `~astropy.time.TimeDelta` 2635 Quantity or TimeDelta object for comparison. 2636 atol : `~astropy.units.Quantity` or `~astropy.time.TimeDelta` 2637 Absolute tolerance for equality with units of time (e.g. ``u.s`` or 2638 ``u.day``). Default is one bit in the 128-bit JD time representation, 2639 equivalent to about 20 picosecs. 2640 rtol : float 2641 Relative tolerance for equality 2642 """ 2643 try: 2644 other_day = other.to_value(u.day) 2645 except Exception as err: 2646 raise TypeError(f"'other' argument must support conversion to days: {err}") 2647 2648 if atol is None: 2649 atol = np.finfo(float).eps * u.day 2650 2651 if not isinstance(atol, (u.Quantity, TimeDelta)): 2652 raise TypeError("'atol' argument must be a Quantity or TimeDelta instance, got " 2653 f'{atol.__class__.__name__} instead') 2654 2655 return np.isclose(self.to_value(u.day), other_day, 2656 rtol=rtol, atol=atol.to_value(u.day)) 2657 2658 2659class ScaleValueError(Exception): 2660 pass 2661 2662 2663def _make_array(val, copy=False): 2664 """ 2665 Take ``val`` and convert/reshape to an array. If ``copy`` is `True` 2666 then copy input values. 2667 2668 Returns 2669 ------- 2670 val : ndarray 2671 Array version of ``val``. 2672 """ 2673 if isinstance(val, (tuple, list)) and len(val) > 0 and isinstance(val[0], Time): 2674 dtype = object 2675 else: 2676 dtype = None 2677 2678 val = np.array(val, copy=copy, subok=True, dtype=dtype) 2679 2680 # Allow only float64, string or object arrays as input 2681 # (object is for datetime, maybe add more specific test later?) 2682 # This also ensures the right byteorder for float64 (closes #2942). 2683 if val.dtype.kind == "f" and val.dtype.itemsize >= np.dtype(np.float64).itemsize: 2684 pass 2685 elif val.dtype.kind in 'OSUMaV': 2686 pass 2687 else: 2688 val = np.asanyarray(val, dtype=np.float64) 2689 2690 return val 2691 2692 2693def _check_for_masked_and_fill(val, val2): 2694 """ 2695 If ``val`` or ``val2`` are masked arrays then fill them and cast 2696 to ndarray. 2697 2698 Returns a mask corresponding to the logical-or of masked elements 2699 in ``val`` and ``val2``. If neither is masked then the return ``mask`` 2700 is ``None``. 2701 2702 If either ``val`` or ``val2`` are masked then they are replaced 2703 with filled versions of themselves. 2704 2705 Parameters 2706 ---------- 2707 val : ndarray or MaskedArray 2708 Input val 2709 val2 : ndarray or MaskedArray 2710 Input val2 2711 2712 Returns 2713 ------- 2714 mask, val, val2: ndarray or None 2715 Mask: (None or bool ndarray), val, val2: ndarray 2716 """ 2717 def get_as_filled_ndarray(mask, val): 2718 """ 2719 Fill the given MaskedArray ``val`` from the first non-masked 2720 element in the array. This ensures that upstream Time initialization 2721 will succeed. 2722 2723 Note that nothing happens if there are no masked elements. 2724 """ 2725 fill_value = None 2726 2727 if np.any(val.mask): 2728 # Final mask is the logical-or of inputs 2729 mask = mask | val.mask 2730 2731 # First unmasked element. If all elements are masked then 2732 # use fill_value=None from above which will use val.fill_value. 2733 # As long as the user has set this appropriately then all will 2734 # be fine. 2735 val_unmasked = val.compressed() # 1-d ndarray of unmasked values 2736 if len(val_unmasked) > 0: 2737 fill_value = val_unmasked[0] 2738 2739 # Fill the input ``val``. If fill_value is None then this just returns 2740 # an ndarray view of val (no copy). 2741 val = val.filled(fill_value) 2742 2743 return mask, val 2744 2745 mask = False 2746 if isinstance(val, np.ma.MaskedArray): 2747 mask, val = get_as_filled_ndarray(mask, val) 2748 if isinstance(val2, np.ma.MaskedArray): 2749 mask, val2 = get_as_filled_ndarray(mask, val2) 2750 2751 return mask, val, val2 2752 2753 2754class OperandTypeError(TypeError): 2755 def __init__(self, left, right, op=None): 2756 op_string = '' if op is None else f' for {op}' 2757 super().__init__( 2758 "Unsupported operand type(s){}: " 2759 "'{}' and '{}'".format(op_string, 2760 left.__class__.__name__, 2761 right.__class__.__name__)) 2762 2763 2764def _check_leapsec(): 2765 global _LEAP_SECONDS_CHECK 2766 if _LEAP_SECONDS_CHECK != _LeapSecondsCheck.DONE: 2767 from astropy.utils import iers 2768 with _LEAP_SECONDS_LOCK: 2769 # There are three ways we can get here: 2770 # 1. First call (NOT_STARTED). 2771 # 2. Re-entrant call (RUNNING). We skip the initialisation 2772 # and don't worry about leap second errors. 2773 # 3. Another thread which raced with the first call 2774 # (RUNNING). The first thread has relinquished the 2775 # lock to us, so initialization is complete. 2776 if _LEAP_SECONDS_CHECK == _LeapSecondsCheck.NOT_STARTED: 2777 _LEAP_SECONDS_CHECK = _LeapSecondsCheck.RUNNING 2778 update_leap_seconds() 2779 _LEAP_SECONDS_CHECK = _LeapSecondsCheck.DONE 2780 2781 2782def update_leap_seconds(files=None): 2783 """If the current ERFA leap second table is out of date, try to update it. 2784 2785 Uses `astropy.utils.iers.LeapSeconds.auto_open` to try to find an 2786 up-to-date table. See that routine for the definition of "out of date". 2787 2788 In order to make it safe to call this any time, all exceptions are turned 2789 into warnings, 2790 2791 Parameters 2792 ---------- 2793 files : list of path-like, optional 2794 List of files/URLs to attempt to open. By default, uses defined by 2795 `astropy.utils.iers.LeapSeconds.auto_open`, which includes the table 2796 used by ERFA itself, so if that is up to date, nothing will happen. 2797 2798 Returns 2799 ------- 2800 n_update : int 2801 Number of items updated. 2802 2803 """ 2804 try: 2805 from astropy.utils import iers 2806 2807 table = iers.LeapSeconds.auto_open(files) 2808 return erfa.leap_seconds.update(table) 2809 2810 except Exception as exc: 2811 warn("leap-second auto-update failed due to the following " 2812 f"exception: {exc!r}", AstropyWarning) 2813 return 0 2814