Changeset 2457 in lmdz_wrf

Timestamp:

Apr 22, 2019, 8:53:32 PM (6 years ago)

Author:

lfita

Message:

Adding `fix_CFdates'

File:

• trunk/tools/generic_tools.py (modified) (3 diffs)

trunk/tools/generic_tools.py

@@ -727 +727 @@
         if not searchInlist(availcalendar , calendar):
             print errormsg
-            print '  ' + fname + ": calendar '" + + "' not ready !!"
+            print '  ' + fname + ": calendar '" + calendar + "' not ready !!"
             print '    avaialable ones:', availcalendar
             quit(-1)
@@ -754 +754 @@
               '_00:00:00')
 
+        removeleap = False
         if calendar == 'noleap' or calendar == '365d':
-            removeleap = False
             newdate = advance_cfDate(refdate, timeval, tunits)
             if int(yrref) <= newdate.year:
@@ -15803 +15803 @@
 #        print '  ' , '(', inds[i,j,0], ',', inds[i,j,1], ')'
 
+def fix_CFdates(timevals, origcftimeu, origcal,                                      \
+  newcftimeu= 'hours since 1949-12-01 00:00:00', newcal='gregorian'):
+    """ Fixing CF-time values with wrong setting
+       NOTE: CF-times does not authorize timeu = 'month'
+             python datetime package does not work for dates < 1900
+      timevals: actual CF time-values
+      origcftimeu: current original CF-time values [timeu] since [Refdate]
+      newcftimeu: new correct CF and python time units
+    >>> times = [15.5, 45, 74.5, 105, 135.5, 166, 196.5, 227.5, 258, 288.5, 319]
+    >>> fix_CFdates(times, 'days since 1850-01-01 00:00:00', 'noleap', 
+      'days since 1949-12-01 00:00:00')
+    [-36477.5 -36448.  -36418.5 -36388.  -36357.5 -36327.  -36296.5 -36265.5
+     -36235.  -36204.5 -36174. ]
+    """
+    import datetime as dt
+    import subprocess as sub
+
+    fname = 'fix_CFdates'
+
+    availtu = ['years', 'weeks', 'days', 'hours', 'minutes', 'seconds']
+
+    # original
+    origtv = origcftimeu.split(' ')
+    origtunits = origtv[0]
+    origSrefdate = origtv[2]
+
+    yrref=origSrefdate[0:4]
+    monref=origSrefdate[5:7]
+    dayref=origSrefdate[8:10]
+
+    if int(yrref) < 1970: origsing = -1.
+    else: origsing = 1.
+
+# Does original reference date contain a time value [YYYY]-[MM]-[DD] [HH]:[MI]:[SS]
+##
+    trefT = origcftimeu.find(':')
+    if not trefT == -1:
+        if len(origtv) == 3:
+            horref=int(origSrefdate[11:13])
+            minref=int(origSrefdate[14:16])
+            secref=int(origSrefdate[17:19])
+        else:
+            origSreftime = origtv[3]
+            horref=int(origSreftime[0:2])
+            minref=int(origSreftime[3:5])
+            secref=int(origSreftime[6:8])
+    else:
+        horref = 0
+        minref = 0
+        secref = 0
+
+    addsecs = horref*3600 + minref*60 + secref
+    addsS = str(addsecs)
+
+    # Using shell date `coreutil'
+    seconds = sub.Popen('/bin/date +%s -u -d"' + origSrefdate[0:10] + ' ' + addsS +  \
+      ' seconds"', shell=True, stdout=sub.PIPE)
+    origsecondsdiff = np.float64(seconds.communicate()[0].replace('\n',''))
+
+    # new
+    newtv = newcftimeu.split(' ')
+    newtunits = newtv[0]
+    newSrefdate = newtv[2]
+
+    yrref=newSrefdate[0:4]
+    monref=newSrefdate[5:7]
+    dayref=newSrefdate[8:10]  
+
+    if int(yrref) < 1970: newsing = -1.
+    else: newsing = 1.
+
+# Does reference date contain a time value [YYYY]-[MM]-[DD] [HH]:[MI]:[SS]
+##
+    trefT = newcftimeu.find(':')
+    if not trefT == -1:
+        if len(newtv) == 3:
+            horref=int(newSrefdate[11:13])
+            minref=int(newSrefdate[14:16])
+            secref=int(newSrefdate[17:19])
+        else:
+            newSreftime = newtv[3]
+            horref=int(newSreftime[0:2])
+            minref=int(newSreftime[3:5])
+            secref=int(newSreftime[6:8])
+    else:
+        horref = 0
+        minref = 0
+        secref = 0
+
+    addsecs = horref*3600 + minref*60 + secref
+    addsS = str(addsecs)
+
+    # Using shell date `coreutil'
+    seconds = sub.Popen('/bin/date +%s -u -d"' + newSrefdate[0:10] + ' ' + addsS +   \
+      ' seconds"', shell=True, stdout=sub.PIPE)
+    newsecondsdiff = np.float64(seconds.communicate()[0].replace('\n',''))
+
+    # Seconds difference respect new date
+    secsdiff = origsecondsdiff - newsecondsdiff
+
+    if origtunits == 'years':
+        secsdtr = 365 * 24 * 3600.
+    elif origtunits == 'weeks':
+        secsdtr = 7 * 24 * 3600.
+    elif origtunits == 'days':
+        secsdtr = 24 * 3600.
+    elif origtunits == 'hours':
+        secsdtr = 3600.
+    elif origtunits == 'minutes':
+        secsdtr = 60.
+    elif origtunits == 'seconds':
+        secsdtr = 1.
+    else:
+        print errormsg
+        print '  ' + fname + ": original time-units '" + origtunits + "' not ready !!"
+        print '    available ones: ', availtu
+        quit(-1)
+
+    if type(timevals) != type(np.ones((2), dtype=int)):
+        newtimevals = np.array(timevals, dtype=np.float64)
+        newtimevals = newtimevals*secsdtr
+    else:
+        newtimevals = timevals*secsdtr
+    newtimevals = newtimevals + origsecondsdiff
+    newtimevals = newtimevals - newsecondsdiff
+
+    # Re-scaling to the new time units
+    if newtunits == 'years':
+        secsdtr = 365 * 24 * 3600.
+    elif newtunits == 'weeks':
+        secsdtr = 7 * 24 * 3600.
+    elif newtunits == 'days':
+        secsdtr = 24 * 3600.
+    elif newtunits == 'hours':
+        secsdtr = 3600.
+    elif newtunits == 'minutes':
+        secsdtr = 60.
+    elif newtunits == 'seconds':
+        secsdtr = 1.
+    else:
+        print errormsg
+        print '  ' + fname + ": new time-units '" + newtunits + "' not ready !!"
+        print '    available ones: ', availtu
+        quit(-1)
+    newtimevals = newtimevals / secsdtr
+
+    if origcal == 'noleap3' or origcal == '365d':
+        print '  ' + fname + ": found non-standard calendar '" + origcal + "' !!"
+        # Imposing to gregorian calendar
+        newtimevals = impose_gregorian(newtimevals, newcftimeu, origcal)
+
+    return newtimevals
+
+times = [15.5, 45, 74.5, 105, 135.5, 166, 196.5, 227.5, 258, 288.5, 319]
+times = [56649.5, 56680, 56710.5, 56741, 56771.5, 56802.5, 56833, 56863.5, 56894, 56924.5]
+print fix_CFdates(times, 'days since 1850-01-01 00:00:00', 'noleap', \
+  'days since 1949-12-01 00:00:00')
 
 #quit()