Changeset 654 in lmdz_wrf

Timestamp:

Oct 15, 2015, 11:00:04 AM (10 years ago)

Author:

lfita

Message:

Adding 'WRFzhgt', WRF heights retrieved from the geopotential

Location:

trunk/tools

Files:

• diagnostics.py (modified) (15 diffs)
• variables_values.dat (modified) (1 diff)

trunk/tools/diagnostics.py ¶

@@ -12 +12 @@
 import re
 import nc_var_tools as ncvar
-import datetime as dt
+import datetime as dtime
 
 main = 'diagnostics.py'
 errormsg = 'ERROR -- error -- ERROR -- error'
 warnmsg = 'WARNING -- warning -- WARNING -- warning'
+
+# Constants
+grav = 9.81
 
 # Gneral information
@@ -797 +800 @@
 availdiags = ['ACRAINTOT', 'accum', 'clt', 'cllmh', 'deaccum', 'LMDZrh', 'mslp',     \
   'OMEGAw', 'RAINTOT',   \
-  'rvors', 'td', 'turbulence', 'WRFrvors', 'wds', 'wss']
+  'rvors', 'td', 'turbulence', 'WRFgeop', 'WRFp', 'WRFrvors', 'wds', 'wss',          \
+  'WRFheight']
 
 methods = ['accum', 'deaccum']
@@ -806 +810 @@
   'WRFp', 'WRFtd',                                                                   \
   'WRFpos', 'WRFprc', 'WRFprls', 'WRFrh', 'LMDZrh', 'LMDZrhs', 'WRFrhs', 'WRFrvors', \
-  'WRFt', 'WRFtime', 'WRFwds', 'WRFwss']
+  'WRFt', 'WRFtime', 'WRFwds', 'WRFwss', 'WRFheight']
 
 ofile = 'diagnostics.nc'
@@ -846 +850 @@
 
 # Looking for specific variables that might be use in more than one diagnostic
+WRFgeop_compute = False
 WRFp_compute = False
 WRFt_compute = False
@@ -852 +857 @@
 WRFdens_compute = False
 WRFpos_compute = False
+WRFtime_compute = False
 
 for idiag in range(Ndiags):
     if diags[idiag].split('|')[1].find('@') == -1:
         depvars = diags[idiag].split('|')[1]
+        if depvars == 'WRFgeop':WRFgeop_compute = True
         if depvars == 'WRFp': WRFp_compute = True
         if depvars == 'WRFt': WRFt_compute = True
@@ -862 +869 @@
         if depvars == 'WRFdens': WRFdens_compute = True
         if depvars == 'WRFpos': WRFpos_compute = True
-
+        if depvars == 'WRFtime': WRFtime_compute = True
     else:
         depvars = diags[idiag].split('|')[1].split('@')
+        if ncvar.searchInlist(depvars, 'WRFgeop'): WRFgeop_compute = True
         if ncvar.searchInlist(depvars, 'WRFp'): WRFp_compute = True
         if ncvar.searchInlist(depvars, 'WRFt'): WRFt_compute = True
@@ -871 +879 @@
         if ncvar.searchInlist(depvars, 'WRFdens'): WRFdens_compute = True
         if ncvar.searchInlist(depvars, 'WRFpos'): WRFpos_compute = True
+        if ncvar.searchInlist(depvars, 'WRFtime'): WRFtime_compute = True
+
+if WRFgeop_compute:
+    print '  ' + main + ': Retrieving geopotential value from WRF as PH + PHB'
+    dimv = ncobj.variables['PH'].shape
+    WRFgeop = ncobj.variables['PH'][:] + ncobj.variables['PHB'][:]
 
 if WRFp_compute:
@@ -904 +918 @@
     print '    ' + main + ': computing air density from WRF as ((MU + MUB) * ' +     \
       'DNW)/g ...'
-    grav = 9.81
 
 # Just we need in in absolute values: Size of the central grid cell
@@ -955 +968 @@
     for it in range(1,dt):
         WRFpos[it,:,:] = WRFpos[0,:,:]
+
+if WRFtime_compute:
+    print '    ' + main + ': computing time from WRF as CFtime(Times) ...'
+
+    refdate='19491201000000'
+    tunitsval='minutes'
+
+    timeobj = ncobj.variables['Times']
+    timewrfv = timeobj[:]
+
+    yrref=refdate[0:4]
+    monref=refdate[4:6]
+    dayref=refdate[6:8]
+    horref=refdate[8:10]
+    minref=refdate[10:12]
+    secref=refdate[12:14]
+
+    refdateS = yrref + '-' + monref + '-' + dayref + ' ' + horref + ':' + minref +   \
+      ':' + secref
+
+    dt = timeobj.shape[0]
+    WRFtime = np.zeros((dt), dtype=np.float)
+
+    for it in range(dt):
+        wrfdates = ncvar.datetimeStr_conversion(timewrfv[it,:],'WRFdatetime', 'matYmdHMS')
+        WRFtime[it] = ncvar.realdatetime1_CFcompilant(wrfdates, refdate, tunitsval)
+
+    tunits = tunitsval + ' since ' + refdateS
+
 
 ### ## #
@@ -1143 +1185 @@
         if depvars[2] != 'WRFtime':
             var2 = ncobj.variables[depvars[2]]
+        else:
+            var2 = np.arange(var0.shape[0], dtype=int)
 
         var = var0[:] + var1[:]
@@ -1152 +1196 @@
 
 # Transforming to a flux
-        if var2.shape[0] > 0:
+        if var2.shape[0] > 1:
             if depvars[2] != 'WRFtime':
                 dtimeunits = var2.getncattr('units')
@@ -1167 +1211 @@
                     tmf1 = tmf1 + time1[ic]
                     tmf2 = tmf2 + time2[ic]
-                dtdate1 = dt.datetime.strptime(tmf1,"%Y-%m-%d_%H:%M:%S")
-                dtdate2 = dt.datetime.strptime(tmf2,"%Y-%m-%d_%H:%M:%S")
+                dtdate1 = dtime.datetime.strptime(tmf1,"%Y-%m-%d_%H:%M:%S")
+                dtdate2 = dtime.datetime.strptime(tmf2,"%Y-%m-%d_%H:%M:%S")
                 diffdate12 = dtdate2 - dtdate1
                 dtime = diffdate12.total_seconds()
@@ -1280 +1324 @@
         ncvar.insert_variable(ncobj, 'bils', diagout, dnames, dvnames, newnc)
 
+# WRFgeop geopotential from WRF as PH + PHB
+    elif diag == 'WRFgeop':
+            
+        diagout = WRFgeop
+
+        ncvar.insert_variable(ncobj, 'zg', diagout, dnames, dvnames, newnc)
+
 # WRFp pressure from WRF as P + PB
     elif diag == 'WRFp':
@@ -1341 +1392 @@
 
         ncvar.insert_variable(ncobj, 'rvors', diagout, dnamesvar, dvnamesvar, newnc)
+
+
+# WRFtime (u10, v10, WRFpos)
+    elif diag == 'WRFtime':
+            
+        diagout = WRFtime
+
+        dnamesvar = ['Time']
+        dvnamesvar = ['Times']
+
+        ncvar.insert_variable(ncobj, 'time', diagout, dnamesvar, dvnamesvar, newnc)
 
 # wss (u10, v10)
@@ -1359 +1421 @@
 
         ncvar.insert_variable(ncobj, 'wss', diagout, dnamesvar, dvnamesvar, newnc)
+
+# WRFheight height from WRF geopotential as WRFGeop/g
+    elif diag == 'WRFheight':
+            
+        diagout = WRFgeop/grav
+
+        ncvar.insert_variable(ncobj, 'zhgt', diagout, dnames, dvnames, newnc)
 
     else:

trunk/tools/variables_values.dat ¶

@@ -484 +484 @@
 z, z, z,  0., 100., z, -, Greens
 Z, z, z,  0., 100., z, -, Greens
+zhgt, zhgt, height,  0., 50000., altitude|from|surface, m, rainbow
 zg, zg, geopotential_height, 0., 80000., geopotential|height, m2s-2, rainbow
 WRFght, zg, geopotential_height, 0., 80000., geopotential|height, m2s-2, rainbow