source: CCCC/trunk/ceda_cc/checkT.py @ 407

Subversion URL: http://proj.badc.rl.ac.uk/svn/exarch/CCCC/trunk/ceda_cc/checkT.py
Revision 407, 5.7 KB checked in by mjuckes, 4 years ago (diff)

relaxed check on missing_values attribute for SPECS .. check went beyond requirement

Line 
1
2import string
3
4
5t1 = '19900101'
6
7t2 = ['days since 1980-01-01', '360-day', 3600]
8
9lm = [31,28,31,30,31,30,31,31,30,31,30,31]
10
11## gregorian or standard
12## mixed gregorian/julian
13## leap if after 1582: gregorian leap-year rule, otherwise julian.
14
15## proleptic_gregorian
16
17  ##  A Gregorian calendar extended to dates before 1582-10-15. That is, a year is a leap year if either (i) it is divisible by 4 but not by 100 or (ii) it is divisible by 400.
18## noleap or 365_day; leap = False
19
20## all_leap or 366_day; leap = True
21
22## 360_day; lm=12, no leap
23
24## julian: every 4th year a leap year
25
26class cfCalSupport:
27
28  def __init__(self,calendar):
29    self.clist = ['gregorian','proleptic_gregorian','all_leap','360_day','365_day','julian']
30    self.alist = ['standard','366_day','noleap']
31    self.calias = {'standard':'gregorian', '366_day':'all_leap', 'noleap':'365_day' }
32    assert calendar in self.clist + self.alist, 'Calendar %s not recognised' % calendar
33    self.calendar = calendar
34    self.cal = self.calias.get( calendar, calendar )
35    self.lm = lm
36    if self.cal == 'all_leap':
37      self.lfixed = True
38      self.leap = True
39      self.ylen = 366
40    elif self.cal == '360_day':
41      self.lfixed = True
42      self.leap = False
43      self.ylen = 360
44      self.lm = [30,30,30,30,30,30,30,30,30,30,30,30,]
45    elif self.cal == '365_day':
46      self.lfixed = True
47      self.leap = False
48      self.ylen = 365
49    else:
50      self.lfixed = False
51      self.ylen0 = sum( self.lm) 
52      if self.cal == 'julian':
53        self.ylen = 365.25
54      elif self.cal == 'proleptic_gregorian':
55        self.ylen = 365 + 97./400.
56      elif self.cal == 'gregorian':
57## compromise between gregorian and julian.
58        self.ylen = 365 + 98.5/400.
59        self.ylenx = [365.25,365 + 97./400.]
60     
61    self.alm = [0,]
62    for k in range(11):
63      self.alm.append( self.alm[k] + self.lm[k] )
64
65  def isleap( self, year ):
66    if self.lfixed:
67      return self.leap
68    elif self.cal == 'julian':
69      return self.isJulianLeap(year)
70    elif self.cal == 'proleptic_gregorian':
71      return self.isProlepticGregorianLeap(year)
72    elif year > 1582:
73      return self.isProlepticGregorianLeap(year)
74    else:
75      return self.isJulianLeap(year)
76
77  def isJulianLeap(self,year):
78    iy = int(year)
79    return (iy/4)*4 == iy
80
81  def isProlepticGregorianLeap(self,year):
82    iy = int(year)
83    return ( (iy/4)*4 == iy and (iy/100)*100 != iy ) or ( (iy/400)*400 == iy )
84
85  def dayOff(self,year):
86    """Returns the number of days since year zero for given year"""
87    if self.cal in ['all_leap','360_day','365_day']:
88      return year*self.ylen
89    else:
90      if self.cal == 'julian':
91        return  int(1 + round( 365.25*year - 0.51 )  )
92      elif self.cal == 'proleptic_gregorian':
93## see http://quasar.as.utexas.edu/BillInfo/JulianDatesG.html
94        y = year - 1
95        a = round(y/100 - 0.5)
96        return  int(1 - a + round(a/4 - 0.5) + round( 365.25*year - 0.51 )  )
97      elif self.cal == 'gregorian':
98        if year <= 1583: 
99          return  int(1 + round( 365.25*year - 0.51 )  )
100        else:
101          y = year - 1
102          a = round(y/100 - 0.5)
103          return  int(1 - a + round(a/4 - 0.5) + round( 365.25*year - 0.51 )  ) + 12
104## proleptic gregorian has 12 fewer leap years than julian.
105      else:
106        assert False, "not suported yet"
107 
108  def d2y(self,ybase, d):
109    dy = d/float(self.ylen)
110    if not self.lfixed:
111      #####
112      pass
113      ############ need some work here to get correct day .....
114      ############ need to use known periodicities (4,400 and mixed) as a starting point ###########
115   
116
117  def dayDiff(self,year1,year0):
118    return self.dayOff(year1) - self.dayOff(year0)
119
120  def dayInYear(self,year,month,day):
121    if self.cal == '360_day':
122      return (month-1)*30 + day
123    x = 0
124    if self.isleap(year) and month > 2:
125      x = 1
126    return x + self.alm[month-1] + day
127
128  def yearLen(self,year):
129    if self.lfixed:
130      return self.ylen
131    x=0
132    if self.isleap(year):
133      x = 1
134    return self.ylen0 + x
135
136
137  def setBase( self, base ):
138    bits = string.split( base )
139    assert bits[0] == 'days'
140## only support days here -- for more general support coudl use cf-python, introducing dependency on numpy and python netCDF4
141    assert bits[1] == 'since'
142    assert len(bits) > 2,'Not enough elements in time base specification string %s' % base
143    s = re.compile( '^([0-9]{1,4}-[0-9]{1,2}-[0-9]{1,2})($|T(.*))' )
144    t = re.compile( '^([0-9]{1,2}):([0-9]{1,2}):([0-9.]*)($|Z$)' )
145    bb = s.findall( bits[2] )
146    assert len(bb) > 0, 'Failed to parse time base reference year/month/day: %s' % base
147    y,m,d = map( int, string.split( bb[0], '-' ) )
148    if len(bits) == 3 and len(bb[2]) ==  0:
149      h,mn,s = 0,0,0.
150    else:
151      if len(bits) == 4:
152         hms = bits[3]
153      else:
154         hms = bb[2]
155      cc = t.findall( hms )
156      h = int(cc[0])
157      mn = int(cc[1])
158      s = float(cc[2])
159    self.base = (y,m,d,h,mn,s)
160
161
162## method to obtain extended index -- (time-base time)/dt or (time-base time)*freq
163
164def timeUnitsScan( c1 ):
165    bits = string.split( c1 )
166    assert bits[0] == 'days'
167## only support days here -- for more general support coudl use cf-python, introducing dependency on numpy and python netCDF4
168    assert bits[1] == 'since'
169
170def jabs( y,m,d,cal ):
171  if cal == '360-day':
172    return y*360 + m*30 + d
173
174def c1(t1,t2):
175  y,m,d = map(int,[t1[:4],t1[4:6],t1[6:]])
176  bits = string.split( t2[0] )
177  assert bits[0] == 'days'
178  assert bits[1] == 'since'
179
180  calendar = t2[1]
181  if calendar == '360-day':
182    ly = 360
183    lm = 30
184
185  bt = bits[2]
186  by,bm,bd = map(int,[bt[:4],bt[5:7],bt[8:]])
187  j0 = jabs(y,m,d, calendar )
188  j1 = jabs(by,bm,bd, calendar ) + t2[2]
189  if j0 != j1:
190    print 'dates do not match: %s, %s:: %s' % (t1,str(t2),j1-j0)
191  else:
192    print 'OK'
193
194c1(t1,t2)
Note: See TracBrowser for help on using the repository browser.