1"""
2Performs conversions of netCDF time coordinate data to/from datetime objects.
3"""
4
5from cpython.object cimport (PyObject_RichCompare, Py_LT, Py_LE, Py_EQ,
6                             Py_NE, Py_GT, Py_GE)
7import cython
8import numpy as np
9import re
10import sys
11import time
12from datetime import datetime as datetime_python
13from datetime import timedelta, MINYEAR
14import time                     # strftime
15try:
16    from itertools import izip as zip
17except ImportError:  # python 3.x
18    pass
19
20microsec_units = ['microseconds','microsecond', 'microsec', 'microsecs']
21millisec_units = ['milliseconds', 'millisecond', 'millisec', 'millisecs']
22sec_units =      ['second', 'seconds', 'sec', 'secs', 's']
23min_units =      ['minute', 'minutes', 'min', 'mins']
24hr_units =       ['hour', 'hours', 'hr', 'hrs', 'h']
25day_units =      ['day', 'days', 'd']
26month_units =    ['month', 'months'] # only allowed for 360_day calendar
27_units = microsec_units+millisec_units+sec_units+min_units+hr_units+day_units
28# supported calendars. Includes synonyms ('standard'=='gregorian',
29# '366_day'=='all_leap','365_day'=='noleap')
30# see http://cfconventions.org/cf-conventions/cf-conventions.html#calendar
31# for definitions.
32_calendars = ['standard', 'gregorian', 'proleptic_gregorian',
33              'noleap', 'julian', 'all_leap', '365_day', '366_day', '360_day']
34# Following are number of Days Per Month
35cdef int[12] _dpm      = [31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31]
36cdef int[12] _dpm_leap = [31, 29, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31]
37cdef int[12] _dpm_360  = [30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30]
38# Same as above, but SUM of previous months (no leap years).
39cdef int[13] _spm_365day = [0, 31, 59, 90, 120, 151, 181, 212, 243, 273, 304, 334, 365]
40cdef int[13] _spm_366day = [0, 31, 60, 91, 121, 152, 182, 213, 244, 274, 305, 335, 366]
41# Reverse operator lookup for datetime.__richcmp__
42_rop_lookup = {Py_LT: '__gt__', Py_LE: '__ge__', Py_EQ: '__eq__',
43               Py_GT: '__lt__', Py_GE: '__le__', Py_NE: '__ne__'}
44
45__version__ = '1.0.3.4'
46
47# Adapted from http://delete.me.uk/2005/03/iso8601.html
48# Note: This regex ensures that all ISO8601 timezone formats are accepted - but, due to legacy support for other timestrings, not all incorrect formats can be rejected.
49#       For example, the TZ spec "+01:0" will still work even though the minutes value is only one character long.
50ISO8601_REGEX = re.compile(r"(?P<year>[+-]?[0-9]{1,4})(-(?P<month>[0-9]{1,2})(-(?P<day>[0-9]{1,2})"
51                           r"(((?P<separator1>.)(?P<hour>[0-9]{1,2}):(?P<minute>[0-9]{1,2})(:(?P<second>[0-9]{1,2})(\.(?P<fraction>[0-9]+))?)?)?"
52                           r"((?P<separator2>.?)(?P<timezone>Z|(([-+])([0-9]{2})((:([0-9]{2}))|([0-9]{2}))?)))?)?)?)?"
53                           )
54# Note: The re module apparently does not support branch reset groups that allow redifinition of the same group name in alternative branches as PCRE does.
55#       Using two different group names is also somewhat ugly, but other solutions might hugely inflate the expression. feel free to contribute a better solution.
56TIMEZONE_REGEX = re.compile(
57    "(?P<prefix>[+-])(?P<hours>[0-9]{2})(?:(?::(?P<minutes1>[0-9]{2}))|(?P<minutes2>[0-9]{2}))?")
58
59class real_datetime(datetime_python):
60    """add dayofwk and dayofyr attributes to python datetime instance"""
61    @property
62    def dayofwk(self):
63        # 0=Monday, 6=Sunday
64        return self.weekday()
65    @property
66    def dayofyr(self):
67        return self.timetuple().tm_yday
68    nanosecond = 0 # workaround for pandas bug (cftime issue #77)
69
70# start of the gregorian calendar
71gregorian = real_datetime(1582,10,15)
72
73def _datesplit(timestr):
74    """split a time string into two components, units and the remainder
75    after 'since'
76    """
77    try:
78        (units, sincestring, remainder) = timestr.split(None,2)
79    except ValueError as e:
80        raise ValueError('Incorrectly formatted CF date-time unit_string')
81
82    if sincestring.lower() != 'since':
83        raise ValueError("no 'since' in unit_string")
84
85    return units.lower(), remainder
86
87def _dateparse(timestr):
88    """parse a string of the form time-units since yyyy-mm-dd hh:mm:ss,
89    return a datetime instance"""
90    # same as version in cftime, but returns a timezone naive
91    # python datetime instance with the utc_offset included.
92
93    (units, isostring) = _datesplit(timestr)
94
95    # parse the date string.
96    year, month, day, hour, minute, second, utc_offset =\
97        _parse_date( isostring.strip() )
98    if year >= MINYEAR:
99        basedate = real_datetime(year, month, day, hour, minute, second)
100        # subtract utc_offset from basedate time instance (which is timezone naive)
101        basedate -= timedelta(days=utc_offset/1440.)
102    else:
103        if not utc_offset:
104            basedate = datetime(year, month, day, hour, minute, second)
105        else:
106            raise ValueError('cannot use utc_offset for reference years <= 0')
107    return basedate
108
109def _round_half_up(x):
110    # 'round half up' so 0.5 rounded to 1 (instead of 0 as in numpy.round)
111    return np.ceil(np.floor(2.*x)/2.)
112
113cdef _parse_date_and_units(timestr,calendar='standard'):
114    """parse a string of the form time-units since yyyy-mm-dd hh:mm:ss
115    return a tuple (units,utc_offset, datetimeinstance)"""
116    (units, isostring) = _datesplit(timestr)
117    if not ((units in month_units and calendar=='360_day') or units in _units):
118        if units in month_units and calendar != '360_day':
119            raise ValueError("'months since' units only allowed for '360_day' calendar")
120        else:
121            raise ValueError(
122            "units must be one of 'seconds', 'minutes', 'hours' or 'days' (or singular version of these), got '%s'" % units)
123    # parse the date string.
124    year, month, day, hour, minute, second, utc_offset = _parse_date(
125        isostring.strip())
126    return units, utc_offset, datetime(year, month, day, hour, minute, second)
127
128
129def date2num(dates,units,calendar='standard'):
130        """date2num(dates,units,calendar='standard')
131
132    Return numeric time values given datetime objects. The units
133    of the numeric time values are described by the `units` argument
134    and the `calendar` keyword. The datetime objects must
135    be in UTC with no time-zone offset.  If there is a
136    time-zone offset in `units`, it will be applied to the
137    returned numeric values.
138
139    **`dates`**: A datetime object or a sequence of datetime objects.
140    The datetime objects should not include a time-zone offset.
141
142    **`units`**: a string of the form `<time units> since <reference time>`
143    describing the time units. `<time units>` can be days, hours, minutes,
144    seconds, milliseconds or microseconds. `<reference time>` is the time
145    origin. `months_since` is allowed *only* for the `360_day` calendar.
146
147    **`calendar`**: describes the calendar used in the time calculations.
148    All the values currently defined in the
149    [CF metadata convention](http://cfconventions.org)
150    Valid calendars `'standard', 'gregorian', 'proleptic_gregorian'
151    'noleap', '365_day', '360_day', 'julian', 'all_leap', '366_day'`.
152    Default is `'standard'`, which is a mixed Julian/Gregorian calendar.
153
154    returns a numeric time value, or an array of numeric time values
155    with approximately 100 microsecond accuracy.
156        """
157        calendar = calendar.lower()
158        basedate = _dateparse(units)
159        (unit, ignore) = _datesplit(units)
160        # real-world calendars limited to positive reference years.
161        if calendar in ['julian', 'standard', 'gregorian', 'proleptic_gregorian']:
162            if basedate.year == 0:
163                msg='zero not allowed as a reference year, does not exist in Julian or Gregorian calendars'
164                raise ValueError(msg)
165
166        if (calendar == 'proleptic_gregorian' and basedate.year >= MINYEAR) or \
167           (calendar in ['gregorian','standard'] and basedate > gregorian):
168            # use python datetime module,
169            isscalar = False
170            try:
171                dates[0]
172            except:
173                isscalar = True
174            if isscalar:
175                dates = np.array([dates])
176            else:
177                dates = np.array(dates)
178                shape = dates.shape
179            ismasked = False
180            if np.ma.isMA(dates) and np.ma.is_masked(dates):
181                mask = dates.mask
182                ismasked = True
183            times = []
184            for date in dates.flat:
185                if ismasked and not date:
186                    times.append(None)
187                else:
188                    td = date - basedate
189                    # total time in microseconds.
190                    totaltime = td.microseconds + (td.seconds + td.days * 24 * 3600) * 1.e6
191                    if unit in microsec_units:
192                        times.append(totaltime)
193                    elif unit in millisec_units:
194                        times.append(totaltime/1.e3)
195                    elif unit in sec_units:
196                        times.append(totaltime/1.e6)
197                    elif unit in min_units:
198                        times.append(totaltime/1.e6/60)
199                    elif unit in hr_units:
200                        times.append(totaltime/1.e6/3600)
201                    elif unit in day_units:
202                        times.append(totaltime/1.e6/3600./24.)
203                    else:
204                        raise ValueError('unsupported time units')
205            if isscalar:
206                return times[0]
207            else:
208                return np.reshape(np.array(times), shape)
209        else: # use cftime module for other calendars
210            cdftime = utime(units,calendar=calendar)
211            return cdftime.date2num(dates)
212
213
214def num2date(times,units,calendar='standard',only_use_cftime_datetimes=False):
215    """num2date(times,units,calendar='standard')
216
217    Return datetime objects given numeric time values. The units
218    of the numeric time values are described by the `units` argument
219    and the `calendar` keyword. The returned datetime objects represent
220    UTC with no time-zone offset, even if the specified
221    `units` contain a time-zone offset.
222
223    **`times`**: numeric time values.
224
225    **`units`**: a string of the form `<time units> since <reference time>`
226    describing the time units. `<time units>` can be days, hours, minutes,
227    seconds, milliseconds or microseconds. `<reference time>` is the time
228    origin. `months_since` is allowed *only* for the `360_day` calendar.
229
230    **`calendar`**: describes the calendar used in the time calculations.
231    All the values currently defined in the
232    [CF metadata convention](http://cfconventions.org)
233    Valid calendars `'standard', 'gregorian', 'proleptic_gregorian'
234    'noleap', '365_day', '360_day', 'julian', 'all_leap', '366_day'`.
235    Default is `'standard'`, which is a mixed Julian/Gregorian calendar.
236
237    **`only_use_cftime_datetimes`**: if False (default), datetime.datetime
238    objects are returned from num2date where possible; if True dates which
239    subclass cftime.datetime are returned for all calendars.
240
241    returns a datetime instance, or an array of datetime instances with
242    approximately 100 microsecond accuracy.
243
244    ***Note***: The datetime instances returned are 'real' python datetime
245    objects if `calendar='proleptic_gregorian'`, or
246    `calendar='standard'` or `'gregorian'`
247    and the date is after the breakpoint between the Julian and
248    Gregorian calendars (1582-10-15). Otherwise, they are 'phony' datetime
249    objects which support some but not all the methods of 'real' python
250    datetime objects. The datetime instances
251    do not contain a time-zone offset, even if the specified `units`
252    contains one.
253    """
254    calendar = calendar.lower()
255    basedate = _dateparse(units)
256    (unit, ignore) = _datesplit(units)
257    # real-world calendars limited to positive reference years.
258    if calendar in ['julian', 'standard', 'gregorian', 'proleptic_gregorian']:
259        if basedate.year == 0:
260            msg='zero not allowed as a reference year, does not exist in Julian or Gregorian calendars'
261            raise ValueError(msg)
262
263    postimes =  (np.asarray(times) > 0).all() # netcdf4-python issue #659
264    if only_use_cftime_datetimes:
265        cdftime = utime(units, calendar=calendar,
266                        only_use_cftime_datetimes=only_use_cftime_datetimes)
267        return cdftime.num2date(times)
268    elif postimes and ((calendar == 'proleptic_gregorian' and basedate.year >= MINYEAR) or \
269       (calendar in ['gregorian','standard'] and basedate > gregorian)):
270        # use python datetime module,
271        isscalar = False
272        try:
273            times[0]
274        except:
275            isscalar = True
276        if isscalar:
277            times = np.array([times],dtype='d')
278        else:
279            times = np.array(times, dtype='d')
280            shape = times.shape
281        ismasked = False
282        if np.ma.isMA(times) and np.ma.is_masked(times):
283            mask = times.mask
284            ismasked = True
285        dates = []
286        for time in times.flat:
287            if ismasked and not time:
288                dates.append(None)
289            else:
290                # convert to total seconds
291                if unit in microsec_units:
292                    tsecs = time/1.e6
293                elif unit in millisec_units:
294                    tsecs = time/1.e3
295                elif unit in sec_units:
296                    tsecs = time
297                elif unit in min_units:
298                    tsecs = time*60.
299                elif unit in hr_units:
300                    tsecs = time*3600.
301                elif unit in day_units:
302                    tsecs = time*86400.
303                else:
304                    raise ValueError('unsupported time units')
305                # compute time delta.
306                days = tsecs // 86400.
307                msecsd = tsecs*1.e6 - days*86400.*1.e6
308                secs = msecsd // 1.e6
309                msecs = np.round(msecsd - secs*1.e6)
310                td = timedelta(days=days,seconds=secs,microseconds=msecs)
311                # add time delta to base date.
312                date = basedate + td
313                dates.append(date)
314        if isscalar:
315            return dates[0]
316        else:
317            return np.reshape(np.array(dates), shape)
318    else: # use cftime for other calendars
319        cdftime = utime(units,calendar=calendar)
320        return cdftime.num2date(times)
321
322
323def date2index(dates, nctime, calendar=None, select='exact'):
324    """date2index(dates, nctime, calendar=None, select='exact')
325
326    Return indices of a netCDF time variable corresponding to the given dates.
327
328    **`dates`**: A datetime object or a sequence of datetime objects.
329    The datetime objects should not include a time-zone offset.
330
331    **`nctime`**: A netCDF time variable object. The nctime object must have a
332    `units` attribute.
333
334    **`calendar`**: describes the calendar used in the time calculations.
335    All the values currently defined in the
336    [CF metadata convention](http://cfconventions.org)
337    Valid calendars `'standard', 'gregorian', 'proleptic_gregorian'
338    'noleap', '365_day', '360_day', 'julian', 'all_leap', '366_day'`.
339    Default is `'standard'`, which is a mixed Julian/Gregorian calendar.
340    If `calendar` is None, its value is given by `nctime.calendar` or
341    `standard` if no such attribute exists.
342
343    **`select`**: `'exact', 'before', 'after', 'nearest'`
344    The index selection method. `exact` will return the indices perfectly
345    matching the dates given. `before` and `after` will return the indices
346    corresponding to the dates just before or just after the given dates if
347    an exact match cannot be found. `nearest` will return the indices that
348    correspond to the closest dates.
349
350    returns an index (indices) of the netCDF time variable corresponding
351    to the given datetime object(s).
352    """
353    try:
354        nctime.units
355    except AttributeError:
356        raise AttributeError("netcdf time variable is missing a 'units' attribute")
357    if calendar == None:
358        calendar = getattr(nctime, 'calendar', 'standard')
359    calendar = calendar.lower()
360    basedate = _dateparse(nctime.units)
361    # real-world calendars limited to positive reference years.
362    if calendar in ['julian', 'standard', 'gregorian', 'proleptic_gregorian']:
363        if basedate.year == 0:
364            msg='zero not allowed as a reference year, does not exist in Julian or Gregorian calendars'
365            raise ValueError(msg)
366
367    if (calendar == 'proleptic_gregorian' and basedate.year >= MINYEAR) or \
368       (calendar in ['gregorian','standard'] and basedate > gregorian):
369        # use python datetime
370        times = date2num(dates,nctime.units,calendar=calendar)
371        return time2index(times, nctime, calendar, select)
372    else: # use cftime module for other cases
373        return _date2index(dates, nctime, calendar, select)
374
375
376def JulianDayFromDate(date, calendar='standard'):
377    """JulianDayFromDate(date, calendar='standard')
378
379    creates a Julian Day from a 'datetime-like' object.  Returns the fractional
380    Julian Day (approximately 100 microsecond accuracy).
381
382    if calendar='standard' or 'gregorian' (default), Julian day follows Julian
383    Calendar on and before 1582-10-5, Gregorian calendar after 1582-10-15.
384
385    if calendar='proleptic_gregorian', Julian Day follows gregorian calendar.
386
387    if calendar='julian', Julian Day follows julian calendar.
388    """
389
390    # check if input was scalar and change return accordingly
391    isscalar = False
392    try:
393        date[0]
394    except:
395        isscalar = True
396
397    date = np.atleast_1d(np.array(date))
398    year = np.empty(len(date), dtype=np.int32)
399    month = year.copy()
400    day = year.copy()
401    hour = year.copy()
402    minute = year.copy()
403    second = year.copy()
404    microsecond = year.copy()
405    jd = np.empty(year.shape, np.longdouble)
406    cdef long double[:] jd_view = jd
407    cdef Py_ssize_t i_max = len(date)
408    cdef Py_ssize_t i
409    for i in range(i_max):
410        d = date[i]
411        year[i] = d.year
412        month[i] = d.month
413        day[i] = d.day
414        hour[i] = d.hour
415        minute[i] = d.minute
416        second[i] = d.second
417        microsecond[i] = d.microsecond
418        jd_view[i] = <double>_IntJulianDayFromDate(<int>year[i],<int>month[i],<int>day[i],calendar)
419
420    # at this point jd is an integer representing noon UTC on the given
421    # year,month,day.
422    # compute fractional day from hour,minute,second,microsecond
423    fracday = hour / 24.0 + minute / 1440.0 + (second + microsecond/1.e6) / 86400.0
424    jd = jd - 0.5 + fracday
425
426    if isscalar:
427        return jd[0]
428    else:
429        return jd
430
431def DateFromJulianDay(JD, calendar='standard', only_use_cftime_datetimes=False,
432                      return_tuple=False):
433    """
434
435    returns a 'datetime-like' object given Julian Day. Julian Day is a
436    fractional day with approximately 100 microsecond accuracy.
437
438    if calendar='standard' or 'gregorian' (default), Julian day follows Julian
439    Calendar on and before 1582-10-5, Gregorian calendar after  1582-10-15.
440
441    if calendar='proleptic_gregorian', Julian Day follows gregorian calendar.
442
443    if calendar='julian', Julian Day follows julian calendar.
444
445    If only_use_cftime_datetimes is set to True, then cftime.datetime
446    objects are returned for all calendars.  Otherwise the datetime object is a
447    'real' datetime object if the date falls in the Gregorian calendar
448    (i.e. calendar='proleptic_gregorian', or  calendar = 'standard'/'gregorian'
449    and the date is after 1582-10-15).  In all other cases a 'phony' datetime
450    objects are used, which are actually instances of cftime.datetime.
451    """
452
453    julian = np.atleast_1d(np.array(JD, dtype=np.longdouble))
454
455    def getdateinfo(julian):
456        # get the day (Z) and the fraction of the day (F)
457        # use 'round half up' rounding instead of numpy's even rounding
458        # so that 0.5 is rounded to 1.0, not 0 (cftime issue #49)
459        Z = np.atleast_1d(np.int32(_round_half_up(julian)))
460        F = (julian + 0.5 - Z).astype(np.longdouble)
461
462        cdef Py_ssize_t i_max = len(Z)
463        year = np.empty(i_max, dtype=np.int32)
464        month = np.empty(i_max, dtype=np.int32)
465        day = np.empty(i_max, dtype=np.int32)
466        dayofyr = np.zeros(i_max,dtype=np.int32)
467        dayofwk = np.zeros(i_max,dtype=np.int32)
468        cdef int ijd
469        cdef Py_ssize_t i
470        for i in range(i_max):
471            ijd = Z[i]
472            year[i],month[i],day[i],dayofwk[i],dayofyr[i] = _IntJulianDayToDate(ijd,calendar)
473
474        if calendar in ['standard', 'gregorian']:
475            ind_before = np.where(julian < 2299160.5)[0]
476        else:
477            ind_before = None
478
479        # compute hour, minute, second, microsecond, convert to int32
480        hour = np.clip((F * 24.).astype(np.int64), 0, 23)
481        F   -= hour / 24.
482        minute = np.clip((F * 1440.).astype(np.int64), 0, 59)
483        second = np.clip((F - minute / 1440.) * 86400., 0, None)
484        microsecond = (second % 1)*1.e6
485        hour = hour.astype(np.int32)
486        minute = minute.astype(np.int32)
487        second = second.astype(np.int32)
488        microsecond = microsecond.astype(np.int32)
489
490        return year,month,day,hour,minute,second,microsecond,dayofyr,dayofwk,ind_before
491
492    year,month,day,hour,minute,second,microsecond,dayofyr,dayofwk,ind_before =\
493    getdateinfo(julian)
494    # round to nearest second if within ms_eps microseconds
495    # (to avoid ugly errors in datetime formatting - alternative
496    # to adding small offset all the time as was done previously)
497    # see netcdf4-python issue #433 and cftime issue #78
498    # this is done by rounding microsends up or down, then
499    # recomputing year,month,day etc
500    # ms_eps is proportional to julian day,
501    # about 47 microseconds in 2000 for Julian base date in -4713
502    ms_eps = np.atleast_1d(np.array(np.finfo(np.float64).eps,np.longdouble))
503    ms_eps = 86400000000.*np.maximum(ms_eps*julian, ms_eps)
504    microsecond = np.where(microsecond < ms_eps, 0, microsecond)
505    indxms = microsecond > 1000000-ms_eps
506    if indxms.any():
507        julian[indxms] = julian[indxms] + 2*ms_eps[indxms]/86400000000.
508        year,month,day,hour,minute,second,microsecond,dayofyr,dayofwk,ind_before =\
509        getdateinfo(julian)
510        microsecond[indxms] = 0
511
512    # check if input was scalar and change return accordingly
513    isscalar = False
514    try:
515        JD[0]
516    except:
517        isscalar = True
518
519    is_real_dateime = False
520    if calendar == 'proleptic_gregorian':
521        # datetime.datetime does not support years < 1
522        #if year < 0:
523        if only_use_cftime_datetimes:
524            datetime_type = DatetimeProlepticGregorian
525        else:
526            if (year < 0).any(): # netcdftime issue #28
527               datetime_type = DatetimeProlepticGregorian
528            else:
529               is_real_datetime = True
530               datetime_type = real_datetime
531    elif calendar in ('standard', 'gregorian'):
532        # return a 'real' datetime instance if calendar is proleptic
533        # Gregorian or Gregorian and all dates are after the
534        # Julian/Gregorian transition
535        if len(ind_before) == 0 and not only_use_cftime_datetimes:
536            is_real_datetime = True
537            datetime_type = real_datetime
538        else:
539            is_real_datetime = False
540            datetime_type = DatetimeGregorian
541    elif calendar == "julian":
542        datetime_type = DatetimeJulian
543    elif calendar in ["noleap","365_day"]:
544        datetime_type = DatetimeNoLeap
545    elif calendar in ["all_leap","366_day"]:
546        datetime_type = DatetimeAllLeap
547    elif calendar == "360_day":
548        datetime_type = Datetime360Day
549    else:
550        raise ValueError("unsupported calendar: {0}".format(calendar))
551
552    if not isscalar:
553        if return_tuple:
554            return np.array([args for args in
555                            zip(year, month, day, hour, minute, second,
556                                microsecond,dayofwk,dayofyr)])
557        else:
558            if is_real_datetime:
559                return np.array([datetime_type(*args)
560                                 for args in
561                                 zip(year, month, day, hour, minute, second,
562                                     microsecond)])
563            else:
564                return np.array([datetime_type(*args)
565                                 for args in
566                                 zip(year, month, day, hour, minute, second,
567                                     microsecond,dayofwk,dayofyr)])
568
569    else:
570        if return_tuple:
571            return (year[0], month[0], day[0], hour[0],
572                    minute[0], second[0], microsecond[0],
573                    dayofwk[0], dayofyr[0])
574        else:
575            if is_real_datetime:
576                return datetime_type(year[0], month[0], day[0], hour[0],
577                                     minute[0], second[0], microsecond[0])
578            else:
579                return datetime_type(year[0], month[0], day[0], hour[0],
580                                     minute[0], second[0], microsecond[0],
581                                     dayofwk[0], dayofyr[0])
582
583class utime:
584
585    """
586Performs conversions of netCDF time coordinate
587data to/from datetime objects.
588
589To initialize: C{t = utime(unit_string,calendar='standard')}
590
591where
592
593B{C{unit_string}} is a string of the form
594C{'time-units since <time-origin>'} defining the time units.
595
596Valid time-units are days, hours, minutes and seconds (the singular forms
597are also accepted). An example unit_string would be C{'hours
598since 0001-01-01 00:00:00'}. months is allowed as a time unit
599*only* for the 360_day calendar.
600
601The B{C{calendar}} keyword describes the calendar used in the time calculations.
602All the values currently defined in the U{CF metadata convention
603<http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.1/cf-conventions.html#time-coordinate>}
604are accepted. The default is C{'standard'}, which corresponds to the mixed
605Gregorian/Julian calendar used by the C{udunits library}. Valid calendars
606are:
607
608C{'gregorian'} or C{'standard'} (default):
609
610Mixed Gregorian/Julian calendar as defined by udunits.
611
612C{'proleptic_gregorian'}:
613
614A Gregorian calendar extended to dates before 1582-10-15. That is, a year
615is a leap year if either (i) it is divisible by 4 but not by 100 or (ii)
616it is divisible by 400.
617
618C{'noleap'} or C{'365_day'}:
619
620Gregorian calendar without leap years, i.e., all years are 365 days long.
621all_leap or 366_day Gregorian calendar with every year being a leap year,
622i.e., all years are 366 days long.
623
624C{'360_day'}:
625
626All years are 360 days divided into 30 day months.
627
628C{'julian'}:
629
630Proleptic Julian calendar, extended to dates after 1582-10-5. A year is a
631leap year if it is divisible by 4.
632
633The C{L{num2date}} and C{L{date2num}} class methods can used to convert datetime
634instances to/from the specified time units using the specified calendar.
635
636The datetime instances returned by C{num2date} are 'real' python datetime
637objects if the date falls in the Gregorian calendar (i.e.
638C{calendar='proleptic_gregorian', 'standard'} or C{'gregorian'} and
639the date is after 1582-10-15). Otherwise, they are 'phony' datetime
640objects which are actually instances of C{L{cftime.datetime}}.  This is
641because the python datetime module cannot handle the weird dates in some
642calendars (such as C{'360_day'} and C{'all_leap'}) which don't exist in any real
643world calendar.
644
645
646Example usage:
647
648>>> from cftime import utime
649>>> from datetime import  datetime
650>>> cdftime = utime('hours since 0001-01-01 00:00:00')
651>>> date = datetime.now()
652>>> print date
6532016-10-05 08:46:27.245015
654>>>
655>>> t = cdftime.date2num(date)
656>>> print t
65717669840.7742
658>>>
659>>> date = cdftime.num2date(t)
660>>> print date
6612016-10-05 08:46:27.244996
662>>>
663
664The resolution of the transformation operation is approximately a millisecond.
665
666Warning:  Dates between 1582-10-5 and 1582-10-15 do not exist in the
667C{'standard'} or C{'gregorian'} calendars.  An exception will be raised if you pass
668a 'datetime-like' object in that range to the C{L{date2num}} class method.
669
670Words of Wisdom from the British MetOffice concerning reference dates:
671
672"udunits implements the mixed Gregorian/Julian calendar system, as
673followed in England, in which dates prior to 1582-10-15 are assumed to use
674the Julian calendar. Other software cannot be relied upon to handle the
675change of calendar in the same way, so for robustness it is recommended
676that the reference date be later than 1582. If earlier dates must be used,
677it should be noted that udunits treats 0 AD as identical to 1 AD."
678
679@ivar origin: datetime instance defining the origin of the netCDF time variable.
680@ivar calendar:  the calendar used (as specified by the C{calendar} keyword).
681@ivar unit_string:  a string defining the the netCDF time variable.
682@ivar units:  the units part of C{unit_string} (i.e. 'days', 'hours', 'seconds').
683    """
684
685    def __init__(self, unit_string, calendar='standard',
686                 only_use_cftime_datetimes=False):
687        """
688@param unit_string: a string of the form
689C{'time-units since <time-origin>'} defining the time units.
690
691Valid time-units are days, hours, minutes and seconds (the singular forms
692are also accepted). An example unit_string would be C{'hours
693since 0001-01-01 00:00:00'}. months is allowed as a time unit
694*only* for the 360_day calendar.
695
696@keyword calendar: describes the calendar used in the time calculations.
697All the values currently defined in the U{CF metadata convention
698<http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.1/cf-conventions.html#time-coordinate>}
699are accepted. The default is C{'standard'}, which corresponds to the mixed
700Gregorian/Julian calendar used by the C{udunits library}. Valid calendars
701are:
702 - C{'gregorian'} or C{'standard'} (default):
703 Mixed Gregorian/Julian calendar as defined by udunits.
704 - C{'proleptic_gregorian'}:
705 A Gregorian calendar extended to dates before 1582-10-15. That is, a year
706 is a leap year if either (i) it is divisible by 4 but not by 100 or (ii)
707 it is divisible by 400.
708 - C{'noleap'} or C{'365_day'}:
709 Gregorian calendar without leap years, i.e., all years are 365 days long.
710 - C{'all_leap'} or C{'366_day'}:
711 Gregorian calendar with every year being a leap year, i.e.,
712 all years are 366 days long.
713 -C{'360_day'}:
714 All years are 360 days divided into 30 day months.
715 -C{'julian'}:
716 Proleptic Julian calendar, extended to dates after 1582-10-5. A year is a
717 leap year if it is divisible by 4.
718
719@keyword only_use_cftime_datetimes: if False (default), datetime.datetime
720objects are returned from num2date where possible; if True dates which subclass
721cftime.datetime are returned for all calendars.
722
723@returns: A class instance which may be used for converting times from netCDF
724units to datetime objects.
725        """
726        calendar = calendar.lower()
727        if calendar in _calendars:
728            self.calendar = calendar
729        else:
730            raise ValueError(
731                "calendar must be one of %s, got '%s'" % (str(_calendars), calendar))
732        units, tzoffset, self.origin =\
733        _parse_date_and_units(unit_string,calendar)
734        # real-world calendars limited to positive reference years.
735        if self.calendar in ['julian', 'standard', 'gregorian', 'proleptic_gregorian']:
736            if self.origin.year == 0:
737                msg='zero not allowed as a reference year, does not exist in Julian or Gregorian calendars'
738                raise ValueError(msg)
739        self.tzoffset = np.array(tzoffset,dtype=np.longdouble)  # time zone offset in minutes
740        self.units = units
741        self.unit_string = unit_string
742        if self.calendar in ['noleap', '365_day'] and self.origin.month == 2 and self.origin.day == 29:
743            raise ValueError(
744                'cannot specify a leap day as the reference time with the noleap calendar')
745        if self.calendar == '360_day' and self.origin.day > 30:
746            raise ValueError(
747                'there are only 30 days in every month with the 360_day calendar')
748        self._jd0 = JulianDayFromDate(self.origin, calendar=self.calendar)
749        self.only_use_cftime_datetimes = only_use_cftime_datetimes
750
751    def date2num(self, date):
752        """
753        Returns C{time_value} in units described by L{unit_string}, using
754        the specified L{calendar}, given a 'datetime-like' object.
755
756        The datetime object must represent UTC with no time-zone offset.
757        If there is a time-zone offset implied by L{unit_string}, it will
758        be applied to the returned numeric values.
759
760        Resolution is approximately a millisecond.
761
762        If C{calendar = 'standard'} or C{'gregorian'} (indicating
763        that the mixed Julian/Gregorian calendar is to be used), an
764        exception will be raised if the 'datetime-like' object describes
765        a date between 1582-10-5 and 1582-10-15.
766
767        Works for scalars, sequences and numpy arrays.
768        Returns a scalar if input is a scalar, else returns a numpy array.
769        """
770        isscalar = False
771        try:
772            date[0]
773        except:
774            isscalar = True
775        if not isscalar:
776            date = np.array(date)
777            shape = date.shape
778        if isscalar:
779            jdelta = JulianDayFromDate(date, self.calendar)-self._jd0
780        else:
781            jdelta = JulianDayFromDate(date.flat, self.calendar)-self._jd0
782        if not isscalar:
783            jdelta = np.array(jdelta)
784        # convert to desired units, add time zone offset.
785        if self.units in microsec_units:
786            jdelta = jdelta * 86400. * 1.e6  + self.tzoffset * 60. * 1.e6
787        elif self.units in millisec_units:
788            jdelta = jdelta * 86400. * 1.e3  + self.tzoffset * 60. * 1.e3
789        elif self.units in sec_units:
790            jdelta = jdelta * 86400. + self.tzoffset * 60.
791        elif self.units in min_units:
792            jdelta = jdelta * 1440. + self.tzoffset
793        elif self.units in hr_units:
794            jdelta = jdelta * 24. + self.tzoffset / 60.
795        elif self.units in day_units:
796            jdelta = jdelta + self.tzoffset / 1440.
797        elif self.units in month_units and self.calendar == '360_day':
798            jdelta = jdelta/30. + self.tzoffset / (30. * 1440.)
799        else:
800            raise ValueError('unsupported time units')
801        if isscalar:
802            return jdelta.astype(np.float64)
803        else:
804            return np.reshape(jdelta.astype(np.float64), shape)
805
806    def num2date(self, time_value):
807        """
808        Return a 'datetime-like' object given a C{time_value} in units
809        described by L{unit_string}, using L{calendar}.
810
811        dates are in UTC with no offset, even if L{unit_string} contains
812        a time zone offset from UTC.
813
814        Resolution is approximately a millisecond.
815
816        Works for scalars, sequences and numpy arrays.
817        Returns a scalar if input is a scalar, else returns a numpy array.
818
819        The datetime instances returned by C{num2date} are 'real' python datetime
820        objects if the date falls in the Gregorian calendar (i.e.
821        C{calendar='proleptic_gregorian'}, or C{calendar = 'standard'/'gregorian'} and
822        the date is after 1582-10-15). Otherwise, they are 'phony' datetime
823        objects which are actually instances of cftime.datetime.  This is
824        because the python datetime module cannot handle the weird dates in some
825        calendars (such as C{'360_day'} and C{'all_leap'}) which
826        do not exist in any real world calendar.
827        """
828        isscalar = False
829        try:
830            time_value[0]
831        except:
832            isscalar = True
833        ismasked = False
834        if np.ma.isMA(time_value) and np.ma.is_masked(time_value):
835            mask = time_value.mask
836            ismasked = True
837        if not isscalar:
838            time_value = np.array(time_value, dtype='d')
839            shape = time_value.shape
840        # convert to desired units, subtract time zone offset.
841        if self.units in microsec_units:
842            jdelta = time_value / 86400000000. - self.tzoffset / 1440.
843        elif self.units in millisec_units:
844            jdelta = time_value / 86400000. - self.tzoffset / 1440.
845        elif self.units in sec_units:
846            jdelta = time_value / 86400. - self.tzoffset / 1440.
847        elif self.units in min_units:
848            jdelta = time_value / 1440. - self.tzoffset / 1440.
849        elif self.units in hr_units:
850            jdelta = time_value / 24. - self.tzoffset / 1440.
851        elif self.units in day_units:
852            jdelta = time_value - self.tzoffset / 1440.
853        elif self.units in month_units and self.calendar == '360_day':
854            # only allowed for 360_day calendar
855            jdelta = time_value * 30. - self.tzoffset / 1440.
856        else:
857            raise ValueError('unsupported time units')
858        jd = self._jd0 + jdelta
859        if not isscalar:
860            if ismasked:
861                date = []
862                for j, m in zip(jd.flat, mask.flat):
863                    if not m:
864                        date.append(DateFromJulianDay(j, self.calendar,
865                                                      self.only_use_cftime_datetimes))
866                    else:
867                        date.append(None)
868            else:
869                date = DateFromJulianDay(jd.flat, self.calendar,
870                                         self.only_use_cftime_datetimes)
871        else:
872            if ismasked and mask.item():
873                date = None
874            else:
875                date = DateFromJulianDay(jd, self.calendar,
876                                         self.only_use_cftime_datetimes)
877        if isscalar:
878            return date
879        else:
880            return np.reshape(np.array(date), shape)
881
882
883cdef _parse_timezone(tzstring):
884    """Parses ISO 8601 time zone specs into tzinfo offsets
885
886    Adapted from pyiso8601 (http://code.google.com/p/pyiso8601/)
887    """
888    if tzstring == "Z":
889        return 0
890    # This isn't strictly correct, but it's common to encounter dates without
891    # time zones so I'll assume the default (which defaults to UTC).
892    if tzstring is None:
893        return 0
894    m = TIMEZONE_REGEX.match(tzstring)
895    prefix, hours, minutes1, minutes2 = m.groups()
896    hours = int(hours)
897    # Note: Minutes don't have to be specified in tzstring, so if the group is not found it means minutes is 0.
898    #       Also, due to the timezone regex definition, there are two mutually exclusive groups that might hold the minutes value, so check both.
899    minutes = int(minutes1) if minutes1 is not None else int(minutes2) if minutes2 is not None else 0
900    if prefix == "-":
901        hours = -hours
902        minutes = -minutes
903    return minutes + hours * 60.
904
905
906cpdef _parse_date(datestring):
907    """Parses ISO 8601 dates into datetime objects
908
909    The timezone is parsed from the date string, assuming UTC
910    by default.
911
912    Adapted from pyiso8601 (http://code.google.com/p/pyiso8601/)
913    """
914    if not isinstance(datestring, str) and not isinstance(datestring, unicode):
915        raise ValueError("Expecting a string %r" % datestring)
916    m = ISO8601_REGEX.match(datestring.strip())
917    if not m:
918        raise ValueError("Unable to parse date string %r" % datestring)
919    groups = m.groupdict()
920    tzoffset_mins = _parse_timezone(groups["timezone"])
921    if groups["hour"] is None:
922        groups["hour"] = 0
923    if groups["minute"] is None:
924        groups["minute"] = 0
925    if groups["second"] is None:
926        groups["second"] = 0
927    # if groups["fraction"] is None:
928    #    groups["fraction"] = 0
929    # else:
930    #    groups["fraction"] = int(float("0.%s" % groups["fraction"]) * 1e6)
931    iyear = int(groups["year"])
932    return iyear, int(groups["month"]), int(groups["day"]),\
933        int(groups["hour"]), int(groups["minute"]), int(groups["second"]),\
934        tzoffset_mins
935
936cdef _check_index(indices, times, nctime, calendar, select):
937    """Return True if the time indices given correspond to the given times,
938    False otherwise.
939
940    Parameters:
941
942    indices : sequence of integers
943    Positive integers indexing the time variable.
944
945    times : sequence of times.
946    Reference times.
947
948    nctime : netCDF Variable object
949    NetCDF time object.
950
951    calendar : string
952    Calendar of nctime.
953
954    select : string
955    Index selection method.
956    """
957    N = nctime.shape[0]
958    if (indices < 0).any():
959        return False
960
961    if (indices >= N).any():
962        return False
963
964    try:
965        t = nctime[indices]
966        nctime = nctime
967    # WORKAROUND TO CHANGES IN SLICING BEHAVIOUR in 1.1.2
968    # this may be unacceptably slow...
969    # if indices are unsorted, or there are duplicate
970    # values in indices, read entire time variable into numpy
971    # array so numpy slicing rules can be used.
972    except IndexError:
973        nctime = nctime[:]
974        t = nctime[indices]
975# if fancy indexing not available, fall back on this.
976#   t=[]
977#   for ind in indices:
978#       t.append(nctime[ind])
979
980    if select == 'exact':
981        return np.all(t == times)
982
983    elif select == 'before':
984        ta = nctime[np.clip(indices + 1, 0, N - 1)]
985        return np.all(t <= times) and np.all(ta > times)
986
987    elif select == 'after':
988        tb = nctime[np.clip(indices - 1, 0, N - 1)]
989        return np.all(t >= times) and np.all(tb < times)
990
991    elif select == 'nearest':
992        ta = nctime[np.clip(indices + 1, 0, N - 1)]
993        tb = nctime[np.clip(indices - 1, 0, N - 1)]
994        delta_after = ta - t
995        delta_before = t - tb
996        delta_check = np.abs(times - t)
997        return np.all(delta_check <= delta_after) and np.all(delta_check <= delta_before)
998
999
1000def _date2index(dates, nctime, calendar=None, select='exact'):
1001    """
1002    _date2index(dates, nctime, calendar=None, select='exact')
1003
1004    Return indices of a netCDF time variable corresponding to the given dates.
1005
1006    @param dates: A datetime object or a sequence of datetime objects.
1007    The datetime objects should not include a time-zone offset.
1008
1009    @param nctime: A netCDF time variable object. The nctime object must have a
1010    C{units} attribute. The entries are assumed to be stored in increasing
1011    order.
1012
1013    @param calendar: Describes the calendar used in the time calculation.
1014    Valid calendars C{'standard', 'gregorian', 'proleptic_gregorian'
1015    'noleap', '365_day', '360_day', 'julian', 'all_leap', '366_day'}.
1016    Default is C{'standard'}, which is a mixed Julian/Gregorian calendar
1017    If C{calendar} is None, its value is given by C{nctime.calendar} or
1018    C{standard} if no such attribute exists.
1019
1020    @param select: C{'exact', 'before', 'after', 'nearest'}
1021    The index selection method. C{exact} will return the indices perfectly
1022    matching the dates given. C{before} and C{after} will return the indices
1023    corresponding to the dates just before or just after the given dates if
1024    an exact match cannot be found. C{nearest} will return the indices that
1025    correspond to the closest dates.
1026    """
1027    try:
1028        nctime.units
1029    except AttributeError:
1030        raise AttributeError("netcdf time variable is missing a 'units' attribute")
1031    # Setting the calendar.
1032    if calendar == None:
1033        calendar = getattr(nctime, 'calendar', 'standard')
1034    cdftime = utime(nctime.units,calendar=calendar)
1035    times = cdftime.date2num(dates)
1036    return time2index(times, nctime, calendar=calendar, select=select)
1037
1038
1039def time2index(times, nctime, calendar=None, select='exact'):
1040    """
1041    time2index(times, nctime, calendar=None, select='exact')
1042
1043    Return indices of a netCDF time variable corresponding to the given times.
1044
1045    @param times: A numeric time or a sequence of numeric times.
1046
1047    @param nctime: A netCDF time variable object. The nctime object must have a
1048    C{units} attribute. The entries are assumed to be stored in increasing
1049    order.
1050
1051    @param calendar: Describes the calendar used in the time calculation.
1052    Valid calendars C{'standard', 'gregorian', 'proleptic_gregorian'
1053    'noleap', '365_day', '360_day', 'julian', 'all_leap', '366_day'}.
1054    Default is C{'standard'}, which is a mixed Julian/Gregorian calendar
1055    If C{calendar} is None, its value is given by C{nctime.calendar} or
1056    C{standard} if no such attribute exists.
1057
1058    @param select: C{'exact', 'before', 'after', 'nearest'}
1059    The index selection method. C{exact} will return the indices perfectly
1060    matching the times given. C{before} and C{after} will return the indices
1061    corresponding to the times just before or just after the given times if
1062    an exact match cannot be found. C{nearest} will return the indices that
1063    correspond to the closest times.
1064    """
1065    try:
1066        nctime.units
1067    except AttributeError:
1068        raise AttributeError("netcdf time variable is missing a 'units' attribute")
1069    # Setting the calendar.
1070    if calendar == None:
1071        calendar = getattr(nctime, 'calendar', 'standard')
1072
1073    num = np.atleast_1d(times)
1074    N = len(nctime)
1075
1076    # Trying to infer the correct index from the starting time and the stride.
1077    # This assumes that the times are increasing uniformly.
1078    if len(nctime) >= 2:
1079        t0, t1 = nctime[:2]
1080        dt = t1 - t0
1081    else:
1082        t0 = nctime[0]
1083        dt = 1.
1084    if select in ['exact', 'before']:
1085        index = np.array((num - t0) / dt, int)
1086    elif select == 'after':
1087        index = np.array(np.ceil((num - t0) / dt), int)
1088    else:
1089        index = np.array(np.around((num - t0) / dt), int)
1090
1091    # Checking that the index really corresponds to the given time.
1092    # If the times do not correspond, then it means that the times
1093    # are not increasing uniformly and we try the bisection method.
1094    if not _check_index(index, times, nctime, calendar, select):
1095
1096        # Use the bisection method. Assumes nctime is ordered.
1097        import bisect
1098        index = np.array([bisect.bisect_right(nctime, n) for n in num], int)
1099        before = index == 0
1100
1101        index = np.array([bisect.bisect_left(nctime, n) for n in num], int)
1102        after = index == N
1103
1104        if select in ['before', 'exact'] and np.any(before):
1105            raise ValueError(
1106                'Some of the times given are before the first time in `nctime`.')
1107
1108        if select in ['after', 'exact'] and np.any(after):
1109            raise ValueError(
1110                'Some of the times given are after the last time in `nctime`.')
1111
1112        # Find the times for which the match is not perfect.
1113        # Use list comprehension instead of the simpler `nctime[index]` since
1114        # not all time objects support numpy integer indexing (eg dap).
1115        index[after] = N - 1
1116        ncnum = np.squeeze([nctime[i] for i in index])
1117        mismatch = np.nonzero(ncnum != num)[0]
1118
1119        if select == 'exact':
1120            if len(mismatch) > 0:
1121                raise ValueError(
1122                    'Some of the times specified were not found in the `nctime` variable.')
1123
1124        elif select == 'before':
1125            index[after] = N
1126            index[mismatch] -= 1
1127
1128        elif select == 'after':
1129            pass
1130
1131        elif select == 'nearest':
1132            nearest_to_left = num[mismatch] < np.array(
1133                [float(nctime[i - 1]) + float(nctime[i]) for i in index[mismatch]]) / 2.
1134            index[mismatch] = index[mismatch] - 1 * nearest_to_left
1135
1136        else:
1137            raise ValueError(
1138                "%s is not an option for the `select` argument." % select)
1139
1140        # Correct for indices equal to -1
1141        index[before] = 0
1142
1143    # convert numpy scalars or single element arrays to python ints.
1144    return _toscalar(index)
1145
1146
1147cdef _toscalar(a):
1148    if a.shape in [(), (1,)]:
1149        return a.item()
1150    else:
1151        return a
1152
1153cdef to_tuple(dt):
1154    """Turn a datetime.datetime instance into a tuple of integers. Elements go
1155    in the order of decreasing significance, making it easy to compare
1156    datetime instances. Parts of the state that don't affect ordering
1157    are omitted. Compare to datetime.timetuple()."""
1158    return (dt.year, dt.month, dt.day, dt.hour, dt.minute,
1159            dt.second, dt.microsecond)
1160
1161# a cache of converters (utime instances) for different calendars
1162cdef dict _converters
1163_converters = {}
1164for calendar in _calendars:
1165    _converters[calendar] = utime("seconds since 1-1-1", calendar)
1166
1167@cython.embedsignature(True)
1168cdef class datetime(object):
1169    """
1170The base class implementing most methods of datetime classes that
1171mimic datetime.datetime but support calendars other than the proleptic
1172Gregorial calendar.
1173    """
1174    cdef readonly int year, month, day, hour, minute, dayofwk, dayofyr
1175    cdef readonly int second, microsecond
1176    cdef readonly str calendar
1177
1178    # Python's datetime.datetime uses the proleptic Gregorian
1179    # calendar. This boolean is used to decide whether a
1180    # cftime.datetime instance can be converted to
1181    # datetime.datetime.
1182    cdef readonly bint datetime_compatible
1183
1184    def __init__(self, int year, int month, int day, int hour=0, int minute=0, int second=0,
1185                 int microsecond=0, int dayofwk=-1, int dayofyr=1):
1186        """dayofyr set to 1 by default - otherwise time.strftime will complain"""
1187
1188        self.year = year
1189        self.month = month
1190        self.day = day
1191        self.hour = hour
1192        self.minute = minute
1193        self.dayofwk = dayofwk # 0 is Monday, 6 is Sunday
1194        self.dayofyr = dayofyr
1195        self.second = second
1196        self.microsecond = microsecond
1197        self.calendar = ""
1198
1199        self.datetime_compatible = True
1200
1201    @property
1202    def format(self):
1203        return '%Y-%m-%d %H:%M:%S'
1204
1205    def strftime(self, format=None):
1206        """
1207        Return a string representing the date, controlled by an explicit format
1208        string. For a complete list of formatting directives, see section
1209        'strftime() and strptime() Behavior' in the base Python documentation.
1210        """
1211        if format is None:
1212            format = self.format
1213        return _strftime(self, format)
1214
1215    def replace(self, **kwargs):
1216        """Return datetime with new specified fields."""
1217        args = {"year": self.year,
1218                "month": self.month,
1219                "day": self.day,
1220                "hour": self.hour,
1221                "minute": self.minute,
1222                "second": self.second,
1223                "microsecond": self.microsecond,
1224                "dayofwk": self.dayofwk,
1225                "dayofyr": self.dayofyr}
1226
1227        for name, value in kwargs.items():
1228            args[name] = value
1229
1230        return self.__class__(**args)
1231
1232    def timetuple(self):
1233        """
1234        Return a time.struct_time such as returned by time.localtime().
1235        The DST flag is -1. d.timetuple() is equivalent to
1236        time.struct_time((d.year, d.month, d.day, d.hour, d.minute,
1237        d.second, d.weekday(), yday, dst)), where yday is the
1238        day number within the current year starting with 1 for January 1st.
1239        """
1240        return time.struct_time((self.year, self.month, self.day, self.hour,
1241                self.minute, self.second, self.dayofwk, self.dayofyr, -1))
1242
1243    cpdef _to_real_datetime(self):
1244        return real_datetime(self.year, self.month, self.day,
1245                             self.hour, self.minute, self.second,
1246                             self.microsecond)
1247
1248    def __repr__(self):
1249        return "{0}.{1}{2}".format('cftime',
1250                                   self.__class__.__name__,
1251                                   self._getstate())
1252
1253    def __str__(self):
1254        return "{:04d}-{:02d}-{:02d} {:02d}:{:02d}:{:02d}".format(
1255            self.year, self.month, self.day, self.hour, self.minute, self.second)
1256
1257    def __hash__(self):
1258        try:
1259            d = self._to_real_datetime()
1260        except ValueError:
1261            return hash(self.timetuple())
1262        return hash(d)
1263
1264    cdef to_tuple(self):
1265        return (self.year, self.month, self.day, self.hour, self.minute,
1266                self.second, self.microsecond)
1267
1268    def __richcmp__(self, other, int op):
1269        cdef datetime dt, dt_other
1270        dt = self
1271        if isinstance(other, datetime):
1272            dt_other = other
1273            # comparing two datetime instances
1274            if dt.calendar == dt_other.calendar:
1275                return PyObject_RichCompare(dt.to_tuple(), dt_other.to_tuple(), op)
1276            else:
1277                # Note: it *is* possible to compare datetime
1278                # instances that use difference calendars by using
1279                # utime.date2num(), but this implementation does
1280                # not attempt it.
1281                raise TypeError("cannot compare {0!r} and {1!r} (different calendars)".format(dt, dt_other))
1282        elif isinstance(other, datetime_python):
1283            # comparing datetime and real_datetime
1284            if not dt.datetime_compatible:
1285                raise TypeError("cannot compare {0!r} and {1!r} (different calendars)".format(self, other))
1286            return PyObject_RichCompare(dt.to_tuple(), to_tuple(other), op)
1287        else:
1288            # With Python3 we can simply return NotImplemented. If the other
1289            # object does not support rich comparison for cftime then a
1290            # TypeError will be automatically raised. However, Python2 is not
1291            # consistent with this Python3 behaviour. In Python2, we only
1292            # delegate the comparison operation to the other object iff it has
1293            # suitable rich comparison support available. This is deduced by
1294            # introspection of the other object. Otherwise, we explicitly raise
1295            # a TypeError to avoid Python2 defaulting to using either __cmp__
1296            # comparision on the other object, or worst still, object ID
1297            # comparison. Either way, at this point the comparision is deemed
1298            # not valid from our perspective.
1299            if sys.version_info.major == 2:
1300                rop = _rop_lookup[op]
1301                if (hasattr(other, '__richcmp__') or hasattr(other, rop)):
1302                    # The other object potentially has the smarts to handle
1303                    # the comparision, so allow the Python machinery to hand
1304                    # the operation off to the other object.
1305                    return NotImplemented
1306                # Otherwise, the comparison is not valid.
1307                emsg = "cannot compare {0!r} and {1!r}"
1308                raise TypeError(emsg.format(self, other))
1309            else:
1310                # Delegate responsibility of comparison to the other object.
1311                return NotImplemented
1312
1313    cdef _getstate(self):
1314        return (self.year, self.month, self.day, self.hour,
1315                self.minute, self.second, self.microsecond,
1316                self.dayofwk, self.dayofyr)
1317
1318    def __reduce__(self):
1319        """special method that allows instance to be pickled"""
1320        return (self.__class__, self._getstate())
1321
1322    cdef _add_timedelta(self, other):
1323        return NotImplemented
1324
1325    def __add__(self, other):
1326        cdef datetime dt
1327        if isinstance(self, datetime) and isinstance(other, timedelta):
1328            dt = self
1329            delta = other
1330        elif isinstance(self, timedelta) and isinstance(other, datetime):
1331            dt = other
1332            delta = self
1333        else:
1334            return NotImplemented
1335        return dt._add_timedelta(delta)
1336
1337    def __sub__(self, other):
1338        cdef datetime dt
1339        if isinstance(self, datetime): # left arg is a datetime instance
1340            dt = self
1341            if isinstance(other, datetime):
1342                # datetime - datetime
1343                if dt.calendar != other.calendar:
1344                    raise ValueError("cannot compute the time difference between dates with different calendars")
1345                if dt.calendar == "":
1346                    raise ValueError("cannot compute the time difference between dates that are not calendar-aware")
1347                converter = _converters[dt.calendar]
1348                return timedelta(seconds=converter.date2num(dt) - converter.date2num(other))
1349            elif isinstance(other, datetime_python):
1350                # datetime - real_datetime
1351                if not dt.datetime_compatible:
1352                    raise ValueError("cannot compute the time difference between dates with different calendars")
1353                return dt._to_real_datetime() - other
1354            elif isinstance(other, timedelta):
1355                # datetime - timedelta
1356                return dt._add_timedelta(-other)
1357            else:
1358                return NotImplemented
1359        else:
1360            if isinstance(self, datetime_python):
1361                # real_datetime - datetime
1362                if not other.datetime_compatible:
1363                    raise ValueError("cannot compute the time difference between dates with different calendars")
1364                return self - other._to_real_datetime()
1365            else:
1366                return NotImplemented
1367
1368@cython.embedsignature(True)
1369cdef class DatetimeNoLeap(datetime):
1370    """
1371Phony datetime object which mimics the python datetime object,
1372but uses the "noleap" ("365_day") calendar.
1373    """
1374    def __init__(self, *args, **kwargs):
1375        datetime.__init__(self, *args, **kwargs)
1376        self.calendar = "noleap"
1377        self.datetime_compatible = False
1378        assert_valid_date(self, no_leap, False, has_year_zero=True)
1379        # if dayofwk, dayofyr not set, calculate them.
1380        if self.dayofwk < 0:
1381            jd = JulianDayFromDate(self,calendar='365_day')
1382            year,month,day,hour,mn,sec,ms,dayofwk,dayofyr =\
1383            DateFromJulianDay(jd,return_tuple=True,calendar='365_day')
1384            self.dayofwk = dayofwk
1385            self.dayofyr = dayofyr
1386
1387    cdef _add_timedelta(self, delta):
1388        return DatetimeNoLeap(*add_timedelta(self, delta, no_leap, False))
1389
1390@cython.embedsignature(True)
1391cdef class DatetimeAllLeap(datetime):
1392    """
1393Phony datetime object which mimics the python datetime object,
1394but uses the "all_leap" ("366_day") calendar.
1395    """
1396    def __init__(self, *args, **kwargs):
1397        datetime.__init__(self, *args, **kwargs)
1398        self.calendar = "all_leap"
1399        self.datetime_compatible = False
1400        assert_valid_date(self, all_leap, False, has_year_zero=True)
1401        # if dayofwk, dayofyr not set, calculate them.
1402        if self.dayofwk < 0:
1403            jd = JulianDayFromDate(self,calendar='366_day')
1404            year,month,day,hour,mn,sec,ms,dayofwk,dayofyr =\
1405            DateFromJulianDay(jd,return_tuple=True,calendar='366_day')
1406            self.dayofwk = dayofwk
1407            self.dayofyr = dayofyr
1408
1409    cdef _add_timedelta(self, delta):
1410        return DatetimeAllLeap(*add_timedelta(self, delta, all_leap, False))
1411
1412@cython.embedsignature(True)
1413cdef class Datetime360Day(datetime):
1414    """
1415Phony datetime object which mimics the python datetime object,
1416but uses the "360_day" calendar.
1417    """
1418    def __init__(self, *args, **kwargs):
1419        datetime.__init__(self, *args, **kwargs)
1420        self.calendar = "360_day"
1421        self.datetime_compatible = False
1422        assert_valid_date(self, no_leap, False, has_year_zero=True, is_360_day=True)
1423        # if dayofwk, dayofyr not set, calculate them.
1424        if self.dayofwk < 0:
1425            jd = JulianDayFromDate(self,calendar='360_day')
1426            year,month,day,hour,mn,sec,ms,dayofwk,dayofyr =\
1427            DateFromJulianDay(jd,return_tuple=True,calendar='360_day')
1428            self.dayofwk = dayofwk
1429            self.dayofyr = dayofyr
1430
1431    cdef _add_timedelta(self, delta):
1432        return Datetime360Day(*add_timedelta_360_day(self, delta))
1433
1434@cython.embedsignature(True)
1435cdef class DatetimeJulian(datetime):
1436    """
1437Phony datetime object which mimics the python datetime object,
1438but uses the "julian" calendar.
1439    """
1440    def __init__(self, *args, **kwargs):
1441        datetime.__init__(self, *args, **kwargs)
1442        self.calendar = "julian"
1443        self.datetime_compatible = False
1444        assert_valid_date(self, is_leap_julian, False)
1445        # if dayofwk, dayofyr not set, calculate them.
1446        if self.dayofwk < 0:
1447            jd = JulianDayFromDate(self,calendar='julian')
1448            year,month,day,hour,mn,sec,ms,dayofwk,dayofyr =\
1449            DateFromJulianDay(jd,return_tuple=True,calendar='julian')
1450            self.dayofwk = dayofwk
1451            self.dayofyr = dayofyr
1452
1453    cdef _add_timedelta(self, delta):
1454        return DatetimeJulian(*add_timedelta(self, delta, is_leap_julian, False))
1455
1456@cython.embedsignature(True)
1457cdef class DatetimeGregorian(datetime):
1458    """
1459Phony datetime object which mimics the python datetime object,
1460but uses the mixed Julian-Gregorian ("standard", "gregorian") calendar.
1461
1462The last date of the Julian calendar is 1582-10-4, which is followed
1463by 1582-10-15, using the Gregorian calendar.
1464
1465Instances using the date after 1582-10-15 can be compared to
1466datetime.datetime instances and used to compute time differences
1467(datetime.timedelta) by subtracting a DatetimeGregorian instance from
1468a datetime.datetime instance or vice versa.
1469    """
1470    def __init__(self, *args, **kwargs):
1471        datetime.__init__(self, *args, **kwargs)
1472        self.calendar = "gregorian"
1473
1474        # dates after 1582-10-15 can be converted to and compared to
1475        # proleptic Gregorian dates
1476        if self.to_tuple() >= (1582, 10, 15, 0, 0, 0, 0):
1477            self.datetime_compatible = True
1478        else:
1479            self.datetime_compatible = False
1480        assert_valid_date(self, is_leap_gregorian, True)
1481        # if dayofwk, dayofyr not set, calculate them.
1482        if self.dayofwk < 0:
1483            jd = JulianDayFromDate(self,calendar='gregorian')
1484            year,month,day,hour,mn,sec,ms,dayofwk,dayofyr =\
1485            DateFromJulianDay(jd,return_tuple=True,calendar='gregorian')
1486            self.dayofwk = dayofwk
1487            self.dayofyr = dayofyr
1488
1489    cdef _add_timedelta(self, delta):
1490        return DatetimeGregorian(*add_timedelta(self, delta, is_leap_gregorian, True))
1491
1492@cython.embedsignature(True)
1493cdef class DatetimeProlepticGregorian(datetime):
1494    """
1495Phony datetime object which mimics the python datetime object,
1496but allows for dates that don't exist in the proleptic gregorian calendar.
1497
1498Supports timedelta operations by overloading + and -.
1499
1500Has strftime, timetuple, replace, __repr__, and __str__ methods. The
1501format of the string produced by __str__ is controlled by self.format
1502(default %Y-%m-%d %H:%M:%S). Supports comparisons with other phony
1503datetime instances using the same calendar; comparison with
1504datetime.datetime instances is possible for cftime.datetime
1505instances using 'gregorian' and 'proleptic_gregorian' calendars.
1506
1507Instance variables are year,month,day,hour,minute,second,microsecond,dayofwk,dayofyr,
1508format, and calendar.
1509    """
1510    def __init__(self, *args, **kwargs):
1511        datetime.__init__(self, *args, **kwargs)
1512        self.calendar = "proleptic_gregorian"
1513        self.datetime_compatible = True
1514        assert_valid_date(self, is_leap_proleptic_gregorian, False)
1515        # if dayofwk, dayofyr not set, calculate them.
1516        if self.dayofwk < 0:
1517            jd = JulianDayFromDate(self,calendar='proleptic_gregorian')
1518            year,month,day,hour,mn,sec,ms,dayofwk,dayofyr =\
1519            DateFromJulianDay(jd,return_tuple=True,calendar='proleptic_gregorian')
1520            self.dayofwk = dayofwk
1521            self.dayofyr = dayofyr
1522
1523    cdef _add_timedelta(self, delta):
1524        return DatetimeProlepticGregorian(*add_timedelta(self, delta,
1525                                                         is_leap_proleptic_gregorian, False))
1526
1527_illegal_s = re.compile(r"((^|[^%])(%%)*%s)")
1528
1529
1530cdef _findall(text, substr):
1531    # Also finds overlaps
1532    sites = []
1533    i = 0
1534    while 1:
1535        j = text.find(substr, i)
1536        if j == -1:
1537            break
1538        sites.append(j)
1539        i = j + 1
1540    return sites
1541
1542# Every 28 years the calendar repeats, except through century leap
1543# years where it's 6 years.  But only if you're using the Gregorian
1544# calendar.  ;)
1545
1546
1547cdef _strftime(datetime dt, fmt):
1548    if _illegal_s.search(fmt):
1549        raise TypeError("This strftime implementation does not handle %s")
1550    # don't use strftime method at all.
1551    # if dt.year > 1900:
1552    #    return dt.strftime(fmt)
1553
1554    year = dt.year
1555    # For every non-leap year century, advance by
1556    # 6 years to get into the 28-year repeat cycle
1557    delta = 2000 - year
1558    off = 6 * (delta // 100 + delta // 400)
1559    year = year + off
1560
1561    # Move to around the year 2000
1562    year = year + ((2000 - year) // 28) * 28
1563    timetuple = dt.timetuple()
1564    s1 = time.strftime(fmt, (year,) + timetuple[1:])
1565    sites1 = _findall(s1, str(year))
1566
1567    s2 = time.strftime(fmt, (year + 28,) + timetuple[1:])
1568    sites2 = _findall(s2, str(year + 28))
1569
1570    sites = []
1571    for site in sites1:
1572        if site in sites2:
1573            sites.append(site)
1574
1575    s = s1
1576    syear = "%4d" % (dt.year,)
1577    for site in sites:
1578        s = s[:site] + syear + s[site + 4:]
1579    return s
1580
1581cdef bint is_leap_julian(int year):
1582    "Return 1 if year is a leap year in the Julian calendar, 0 otherwise."
1583    return _is_leap(year, calendar='julian')
1584
1585cdef bint is_leap_proleptic_gregorian(int year):
1586    "Return 1 if year is a leap year in the Proleptic Gregorian calendar, 0 otherwise."
1587    return _is_leap(year, calendar='proleptic_gregorian')
1588
1589cdef bint is_leap_gregorian(int year):
1590    "Return 1 if year is a leap year in the Gregorian calendar, 0 otherwise."
1591    return _is_leap(year, calendar='standard')
1592
1593cdef bint all_leap(int year):
1594    "Return True for all years."
1595    return True
1596
1597cdef bint no_leap(int year):
1598    "Return False for all years."
1599    return False
1600
1601cdef int * month_lengths(bint (*is_leap)(int), int year):
1602    if is_leap(year):
1603        return _dpm_leap
1604    else:
1605        return _dpm
1606
1607cdef void assert_valid_date(datetime dt, bint (*is_leap)(int),
1608                            bint julian_gregorian_mixed,
1609                            bint has_year_zero=False,
1610                            bint is_360_day=False) except *:
1611    cdef int[12] month_length
1612
1613    if not has_year_zero:
1614        if dt.year == 0:
1615            raise ValueError("invalid year provided in {0!r}".format(dt))
1616    if is_360_day:
1617        month_length = _dpm_360
1618    else:
1619        month_length = month_lengths(is_leap, dt.year)
1620
1621    if dt.month < 1 or dt.month > 12:
1622        raise ValueError("invalid month provided in {0!r}".format(dt))
1623
1624    if dt.day < 1 or dt.day > month_length[dt.month-1]:
1625        raise ValueError("invalid day number provided in {0!r}".format(dt))
1626
1627    if julian_gregorian_mixed and dt.year == 1582 and dt.month == 10 and dt.day > 4 and dt.day < 15:
1628        raise ValueError("{0!r} is not present in the mixed Julian/Gregorian calendar".format(dt))
1629
1630    if dt.hour < 0 or dt.hour > 23:
1631        raise ValueError("invalid hour provided in {0!r}".format(dt))
1632
1633    if dt.minute < 0 or dt.minute > 59:
1634        raise ValueError("invalid minute provided in {0!r}".format(dt))
1635
1636    if dt.second < 0 or dt.second > 59:
1637        raise ValueError("invalid second provided in {0!r}".format(dt))
1638
1639    if dt.microsecond < 0 or dt.microsecond > 999999:
1640        raise ValueError("invalid microsecond provided in {0!r}".format(dt))
1641
1642# Add a datetime.timedelta to a cftime.datetime instance. Uses
1643# integer arithmetic to avoid rounding errors and preserve
1644# microsecond accuracy.
1645#
1646# The argument is_leap is the pointer to a function returning 1 for leap years and 0 otherwise.
1647#
1648# This implementation supports 365_day (no_leap), 366_day (all_leap),
1649# julian, proleptic_gregorian, and the mixed Julian/Gregorian
1650# (standard, gregorian) calendars by using different is_leap and
1651# julian_gregorian_mixed arguments.
1652#
1653# The date of the transition from the Julian to Gregorian calendar and
1654# the number of invalid dates are hard-wired (1582-10-4 is the last day
1655# of the Julian calendar, after which follows 1582-10-15).
1656cdef tuple add_timedelta(datetime dt, delta, bint (*is_leap)(int), bint julian_gregorian_mixed):
1657    cdef int microsecond, second, minute, hour, day, month, year
1658    cdef int delta_microseconds, delta_seconds, delta_days
1659    cdef int* month_length
1660    cdef int extra_days, n_invalid_dates
1661
1662    # extract these inputs here to avoid type conversion in the code below
1663    delta_microseconds = delta.microseconds
1664    delta_seconds = delta.seconds
1665    delta_days = delta.days
1666
1667    # shift microseconds, seconds, days
1668    microsecond = dt.microsecond + delta_microseconds
1669    second = dt.second + delta_seconds
1670    minute = dt.minute
1671    hour = dt.hour
1672    day = dt.day
1673    month = dt.month
1674    year = dt.year
1675
1676    month_length = month_lengths(is_leap, year)
1677
1678    n_invalid_dates = 10 if julian_gregorian_mixed else 0
1679
1680    # Normalize microseconds, seconds, minutes, hours.
1681    second += microsecond // 1000000
1682    microsecond = microsecond % 1000000
1683    minute += second // 60
1684    second = second % 60
1685    hour += minute // 60
1686    minute = minute % 60
1687    extra_days = hour // 24
1688    hour = hour % 24
1689
1690    delta_days += extra_days
1691
1692    while delta_days < 0:
1693        if year == 1582 and month == 10 and day > 14 and day + delta_days < 15:
1694            delta_days -= n_invalid_dates    # skip over invalid dates
1695        if day + delta_days < 1:
1696            delta_days += day
1697            # decrement month
1698            month -= 1
1699            if month < 1:
1700                month = 12
1701                year -= 1
1702                if year == 0:
1703                    year = -1
1704                month_length = month_lengths(is_leap, year)
1705            day = month_length[month-1]
1706        else:
1707            day += delta_days
1708            delta_days = 0
1709
1710    while delta_days > 0:
1711        if year == 1582 and month == 10 and day < 5 and day + delta_days > 4:
1712            delta_days += n_invalid_dates    # skip over invalid dates
1713        if day + delta_days > month_length[month-1]:
1714            delta_days -= month_length[month-1] - (day - 1)
1715            # increment month
1716            month += 1
1717            if month > 12:
1718                month = 1
1719                year += 1
1720                if year == 0:
1721                    year = 1
1722                month_length = month_lengths(is_leap, year)
1723            day = 1
1724        else:
1725            day += delta_days
1726            delta_days = 0
1727
1728    return (year, month, day, hour, minute, second, microsecond, -1, 1)
1729
1730# Add a datetime.timedelta to a cftime.datetime instance with the 360_day calendar.
1731#
1732# Assumes that the 360_day calendar (unlike the rest of supported
1733# calendars) has the year 0. Also, there are no leap years and all
1734# months are 30 days long, so we can compute month and year by using
1735# "//" and "%".
1736cdef tuple add_timedelta_360_day(datetime dt, delta):
1737    cdef int microsecond, second, minute, hour, day, month, year
1738    cdef int delta_microseconds, delta_seconds, delta_days
1739
1740    # extract these inputs here to avoid type conversion in the code below
1741    delta_microseconds = delta.microseconds
1742    delta_seconds = delta.seconds
1743    delta_days = delta.days
1744
1745    # shift microseconds, seconds, days
1746    microsecond = dt.microsecond + delta_microseconds
1747    second = dt.second + delta_seconds
1748    minute = dt.minute
1749    hour = dt.hour
1750    day = dt.day + delta_days
1751    month = dt.month
1752    year = dt.year
1753
1754    # Normalize microseconds, seconds, minutes, hours, days, and months.
1755    second += microsecond // 1000000
1756    microsecond = microsecond % 1000000
1757    minute += second // 60
1758    second = second % 60
1759    hour += minute // 60
1760    minute = minute % 60
1761    day += hour // 24
1762    hour = hour % 24
1763    # day and month are counted from 1; all months have 30 days
1764    month += (day - 1) // 30
1765    day = (day - 1) % 30 + 1
1766    # all years have 12 months
1767    year += (month - 1) // 12
1768    month = (month - 1) % 12 + 1
1769
1770    return (year, month, day, hour, minute, second, microsecond, -1, 1)
1771
1772# Calendar calculations base on calcals.c by David W. Pierce
1773# http://meteora.ucsd.edu/~pierce/calcalcs
1774
1775cdef _is_leap(int year, calendar):
1776    cdef int tyear
1777    cdef bint leap
1778    calendar = _check_calendar(calendar)
1779    if year == 0:
1780        raise ValueError('year zero does not exist in the %s calendar' %\
1781                calendar)
1782    # Because there is no year 0 in the Julian calendar, years -1, -5, -9, etc
1783    # are leap years.
1784    if year < 0:
1785        tyear = year + 1
1786    else:
1787        tyear = year
1788    if calendar == 'proleptic_gregorian' or (calendar == 'standard' and year > 1581):
1789        if tyear % 4: # not divisible by 4
1790            leap = False
1791        elif year % 100: # not divisible by 100
1792            leap = True
1793        elif year % 400: # not divisible by 400
1794            leap = False
1795        else:
1796            leap = True
1797    elif calendar == 'julian' or (calendar == 'standard' and year < 1582):
1798        leap = tyear % 4 == 0
1799    elif calendar == '366_day':
1800        leap = True
1801    else:
1802        leap = False
1803    return leap
1804
1805cdef _IntJulianDayFromDate(int year,int month,int day,calendar,skip_transition=False):
1806    """Compute integer Julian Day from year,month,day and calendar.
1807
1808    Allowed calendars are 'standard', 'gregorian', 'julian',
1809    'proleptic_gregorian','360_day', '365_day', '366_day', 'noleap',
1810    'all_leap'.
1811
1812    'noleap' is a synonym for '365_day'
1813    'all_leap' is a synonym for '366_day'
1814    'gregorian' is a synonym for 'standard'
1815
1816    Negative years allowed back to -4714
1817    (proleptic_gregorian) or -4713 (standard or gregorian calendar).
1818
1819    Negative year values are allowed in 360_day,365_day,366_day calendars.
1820
1821    Integer julian day is number of days since noon UTC -4713-1-1
1822    in the julian or mixed julian/gregorian calendar, or noon UTC
1823    -4714-11-24 in the proleptic_gregorian calendar. Reference
1824    date is noon UTC 0-1-1 for other calendars.
1825
1826    There is no year zero in standard (mixed), julian, or proleptic_gregorian
1827    calendars.
1828
1829    Subtract 0.5 to get 00 UTC on that day.
1830
1831    optional kwarg 'skip_transition':  When True, leave a 10-day
1832    gap in Julian day numbers between Oct 4 and Oct 15 1582 (the transition
1833    from Julian to Gregorian calendars).  Default False, ignored
1834    unless calendar = 'standard'."""
1835    cdef int jday, jday_jul, jday_greg
1836    cdef bint leap
1837    cdef int[12] dpm2use
1838
1839    # validate inputs.
1840    calendar = _check_calendar(calendar)
1841    if month < 1 or month > 12 or day < 1 or day > 31:
1842        msg = "date %04d-%02d-%02d does not exist in the %s calendar" %\
1843        (year,month,day,calendar)
1844        raise ValueError(msg)
1845
1846    # handle all calendars except standard, julian, proleptic_gregorian.
1847    if calendar == '360_day':
1848        return _IntJulianDayFromDate_360day(year,month,day)
1849    elif calendar == '365_day':
1850        return _IntJulianDayFromDate_365day(year,month,day)
1851    elif calendar == '366_day':
1852        return _IntJulianDayFromDate_366day(year,month,day)
1853
1854    # handle standard, julian, proleptic_gregorian calendars.
1855    if year == 0:
1856        raise ValueError('year zero does not exist in the %s calendar' %\
1857                calendar)
1858    if (calendar == 'proleptic_gregorian'         and year < -4714) or\
1859       (calendar in ['julian','standard']  and year < -4713):
1860        raise ValueError('year out of range for %s calendar' % calendar)
1861    leap = _is_leap(year,calendar)
1862    if not leap and month == 2 and day == 29:
1863        raise ValueError('%s is not a leap year' % year)
1864
1865    # add year offset
1866    if year < 0:
1867        year += 4801
1868    else:
1869        year += 4800
1870
1871    if leap:
1872        dpm2use = _dpm_leap
1873    else:
1874        dpm2use = _dpm
1875
1876    jday = day
1877    for m in range(month-1,0,-1):
1878        jday += dpm2use[m-1]
1879
1880    jday_greg = jday + 365*(year-1) + (year-1)//4 - (year-1)//100 + (year-1)//400
1881    jday_greg -= 31739 # fix year offset
1882    jday_jul = jday + 365*(year-1) + (year-1)//4
1883    jday_jul -= 31777 # fix year offset
1884    if calendar == 'julian':
1885        return jday_jul
1886    elif calendar == 'proleptic_gregorian':
1887        return jday_greg
1888    elif calendar == 'standard':
1889        # check for invalid days in mixed calendar (there are 10 missing)
1890        if jday_jul >= 2299161 and jday_jul < 2299171:
1891            raise ValueError('invalid date in mixed calendar')
1892        if jday_jul < 2299161: # 1582 October 15
1893            return jday_jul
1894        else:
1895            if skip_transition:
1896                return jday_greg+10
1897            else:
1898                return jday_greg
1899
1900    return jday
1901
1902cdef _IntJulianDayToDate(int jday,calendar,skip_transition=False):
1903    """Compute the year,month,day,dow,doy given the integer Julian day.
1904    and calendar. (dow = day of week with 0=Mon,6=Sun and doy is day of year).
1905
1906    Allowed calendars are 'standard', 'gregorian', 'julian',
1907    'proleptic_gregorian','360_day', '365_day', '366_day', 'noleap',
1908    'all_leap'.
1909
1910    'noleap' is a synonym for '365_day'
1911    'all_leap' is a synonym for '366_day'
1912    'gregorian' is a synonym for 'standard'
1913
1914    optional kwarg 'skip_transition':  When True, assume a 10-day
1915    gap in Julian day numbers between Oct 4 and Oct 15 1582 (the transition
1916    from Julian to Gregorian calendars).  Default False, ignored
1917    unless calendar = 'standard'."""
1918    cdef int year,month,day,dow,doy,yp1,tjday
1919    cdef int[12] dpm2use
1920    cdef int[13] spm2use
1921
1922    # validate inputs.
1923    calendar = _check_calendar(calendar)
1924
1925    # handle all calendars except standard, julian, proleptic_gregorian.
1926    if calendar == '360_day':
1927        return _IntJulianDayToDate_360day(jday)
1928    elif calendar == '365_day':
1929        return _IntJulianDayToDate_365day(jday)
1930    elif calendar == '366_day':
1931        return _IntJulianDayToDate_366day(jday)
1932
1933    # handle standard, julian, proleptic_gregorian calendars.
1934    if jday < 0:
1935        raise ValueError('julian day must be a positive integer')
1936    # Make first estimate for year. subtract 4714 or 4713 because Julian Day number
1937    # 0 occurs in year 4714 BC in the Gregorian calendar and 4713 BC in the
1938    # Julian calendar.
1939    if calendar == 'proleptic_gregorian':
1940        year = jday//366 - 4714
1941    elif calendar in ['standard','julian']:
1942        year = jday//366 - 4713
1943
1944    # compute day of week.
1945    dow = _get_dow(jday)
1946
1947    if not skip_transition and calendar == 'standard' and jday > 2299160: jday += 10
1948
1949    # Advance years until we find the right one
1950    yp1 = year + 1
1951    if yp1 == 0:
1952       yp1 = 1 # no year 0
1953    tjday = _IntJulianDayFromDate(yp1,1,1,calendar,skip_transition=True)
1954    while jday >= tjday:
1955        year += 1
1956        if year == 0:
1957            year = 1
1958        yp1 = year + 1
1959        if yp1 == 0:
1960            yp1 = 1
1961        tjday = _IntJulianDayFromDate(yp1,1,1,calendar,skip_transition=True)
1962    if _is_leap(year, calendar):
1963        dpm2use = _dpm_leap
1964        spm2use = _spm_366day
1965    else:
1966        dpm2use = _dpm
1967        spm2use = _spm_365day
1968    month = 1
1969    tjday =\
1970    _IntJulianDayFromDate(year,month,dpm2use[month-1],calendar,skip_transition=True)
1971    while jday > tjday:
1972        month += 1
1973        tjday =\
1974        _IntJulianDayFromDate(year,month,dpm2use[month-1],calendar,skip_transition=True)
1975    tjday = _IntJulianDayFromDate(year,month,1,calendar,skip_transition=True)
1976    day = jday - tjday + 1
1977    if month == 1:
1978        doy = day
1979    else:
1980        doy = spm2use[month-1]+day
1981    return year,month,day,dow,doy
1982
1983cdef _get_dow(int jday):
1984    """compute day of week.
1985    0 = Sunday, 6 = Sat, valid after noon UTC"""
1986    cdef int dow
1987    dow = (jday + 1) % 7
1988    # convert to ISO 8601 (0 = Monday, 6 = Sunday), like python datetime
1989    dow -= 1
1990    if dow == -1: dow = 6
1991    return dow
1992
1993cdef _check_calendar(calendar):
1994    """validate calendars, convert to subset of names to get rid of synonyms"""
1995    if calendar not in _calendars:
1996        raise ValueError('unsupported calendar')
1997    calout = calendar
1998    # remove 'gregorian','noleap','all_leap'
1999    if calendar in ['gregorian','standard']:
2000        calout = 'standard'
2001    if calendar == 'noleap':
2002        calout = '365_day'
2003    if calendar == 'all_leap':
2004        calout = '366_day'
2005    return calout
2006
2007cdef _IntJulianDayFromDate_360day(int year,int month,int day):
2008    """Compute integer Julian Day from year,month,day in
2009    360_day calendar"""
2010    return year*360 + (month-1)*30 + day - 1
2011
2012cdef _IntJulianDayFromDate_365day(int year,int month,int day):
2013    """Compute integer Julian Day from year,month,day in
2014    365_day calendar"""
2015    if month == 2 and day == 29:
2016        raise ValueError('no leap days in 365_day calendar')
2017    return year*365 + _spm_365day[month-1] + day - 1
2018
2019cdef _IntJulianDayFromDate_366day(int year,int month,int day):
2020    """Compute integer Julian Day from year,month,day in
2021    366_day calendar"""
2022    return year*366 + _spm_366day[month-1] + day - 1
2023
2024cdef _IntJulianDayToDate_365day(int jday):
2025    """Compute the year,month,day given the integer Julian day
2026    for 365_day calendar."""
2027    cdef int year,month,day,nextra,dow
2028
2029    year = jday//365
2030    nextra = jday - year*365
2031    doy    = nextra + 1 # Julday numbering starts at 0, doy starts at 1
2032    month = 1
2033    while doy > _spm_365day[month]:
2034        month += 1
2035    day = doy - _spm_365day[month-1]
2036
2037    # compute day of week.
2038    dow = _get_dow(jday)
2039
2040    return year,month,day,dow,doy
2041
2042cdef _IntJulianDayToDate_366day(int jday):
2043    """Compute the year,month,day given the integer Julian day
2044    for 366_day calendar."""
2045    cdef int year,month,day,nextra,dow
2046
2047    year = jday//366
2048    nextra = jday - year*366
2049    doy    = nextra + 1 # Julday numbering starts at 0, doy starts at 1
2050    month = 1
2051    while doy > _spm_366day[month]:
2052        month += 1
2053    day = doy - _spm_366day[month-1]
2054
2055    # compute day of week.
2056    dow = _get_dow(jday)
2057
2058    return year,month,day,dow,doy
2059
2060cdef _IntJulianDayToDate_360day(int jday):
2061    """Compute the year,month,day given the integer Julian day
2062    for 360_day calendar."""
2063    cdef int year,month,day,nextra,dow
2064
2065    year = jday//360
2066    nextra = jday - year*360
2067    doy    = nextra + 1 # Julday numbering starts at 0, doy starts at 1
2068    month = nextra//30 + 1
2069    day   = doy - (month-1)*30
2070
2071    # compute day of week.
2072    dow = _get_dow(jday)
2073
2074    return year,month,day,dow,doy
2075