Changeset 773 in lmdz_wrf

Timestamp:

May 24, 2016, 11:05:49 AM (9 years ago)

Author:

lfita

Message:

Adding draw_vertical_levels' and plot_vertical_lev', plotting vertical levels distribution

Location:

trunk/tools

Files:

• drawing.py (modified) (9 diffs)
• drawing_tools.py (modified) (1 diff)

trunk/tools/drawing.py

@@ -3 +3 @@
 from netCDF4 import Dataset as NetCDFFile
 import drawing_tools as drw
+import generic_tools as gen
 from optparse import OptionParser
 import sys
@@ -44 +45 @@
   'draw_topo_geogrid',                                                               \
   'draw_topo_geogrid_boxes', 'draw_trajectories', 'draw_vals_trajectories',          \
-  'draw_vectors',  'list_graphics', 'variable_values']
+  'draw_vectors',  'draw_vertical_levels', 'list_graphics', 'variable_values']
 
 def draw_2D_shad(ncfile, values, varn):
@@ -1345 +1346 @@
         [min/max]Topo: minimum and maximum values of topography to draw
         lonlatL: limits of longitudes and latitudes [lonmin, latmin, lonmax, latmax] or None
-        title: title of the graph
+        title: title of the graph ('!' for spaces)
         graphic_kind: kind of figure (jpg, pdf, png)
         mapvalues: map characteristics [proj],[res]
@@ -1366 +1367 @@
         quit()
 
-    expectargs = '[minTopo]:[maxTopo]:[lonlatL]:[title]:[graphic_kind]:[mapvalues]'
+    expectargs = '[minTopo],[maxTopo]:[lonlatL]:[title]:[graphic_kind]:[mapvalues]'
  
     drw.check_arguments(fname,values,expectargs,':')
@@ -1385 +1386 @@
         lonlatL[3] = np.float(lonlatLS.split(',')[3])
 
-    grtit = values.split(':')[2]
+    grtit = values.split(':')[2].replace('!',' ')
     kindfig = values.split(':')[3]
     mapvalues = values.split(':')[4]
@@ -1416 +1417 @@
         [min/max]Topo: minimum and maximum values of topography to draw
         lonlatL: limits of longitudes and latitudes [lonmin, latmin, lonmax, latmax] or None
-        title: title of the graph
+        title: title of the graph ('!' for spaces)
         graphic_kind: kind of figure (jpg, pdf, png)
         mapvalues: map characteristics [proj],[res]
@@ -1461 +1462 @@
         lonlatL[3] = np.float(lonlatLS.split(',')[3])
 
-    grtit = values.split(':')[2]
+    grtit = values.split(':')[2].replace('!', ' ')
     kindfig = values.split(':')[3]
     mapvalues = values.split(':')[4]
@@ -3662 +3663 @@
     of.close()
 
+def draw_vertical_levels(ncfile, values, varn):
+    """ plotting vertical levels distribution
+    draw_topo_geogrid(ncfile, values, varn)
+      ncfile= file to use
+      values= [zlog]:[dzlog]:[title]:[graphic_kind]:[legloc]
+        zlog: to use logarithmic scale on the height axis ('true/false')
+        dzlog: to use logarithmic scale on the difference of height between levels axis ('true/false')
+        title: title of the graph ('!' for spaces)
+        graphic_kind: kind of figure (jpg, pdf, png)
+        legloc= location of the legend (0, automatic)
+          1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right',
+          5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
+          9: 'upper center', 10: 'center'      kfig= kind of figure
+      varn= name of the variable with the vertical levels
+        'WRFz': for WRF z-levels (computed as (PH + PHB)/g, from a PHB(0,i,j) = 0)
+    """
+    fname = 'draw_vertical_levels'
+
+    if values == 'h':
+        print fname + '_____________________________________________________________'
+        print draw_vertical_levels.__doc__
+        quit()
+
+    expectargs = '[zlog]:[dzlog]:[title]:[graphic_kind]:[legloc]'
+ 
+    drw.check_arguments(fname,values,expectargs,':')
+
+    zlog = values.split(':')[0]
+    dzlog = values.split(':')[1]
+    title = values.split(':')[2].replace('!',' ')
+    kindfig = values.split(':')[3]
+    legloc = int(values.split(':')[4])
+
+    if not os.path.isfile(ncfile):
+        print errormsg
+        print '  ' + fname + ': file "' + ncfile + '" does not exist !!'
+        quit(-1)    
+
+    objf = NetCDFFile(ncfile, 'r')
+
+    if varn == 'WRFz': 
+        if not gen.searchInlist(objf.variables,'PH'):
+            print errormsg
+            print '  ' + fname + ": WRF file '" + ncfile + "' does not have " +      \
+              "variable 'PH' !!"
+            quit(-1)
+        if not gen.searchInlist(objf.variables,'PHB'):
+            print errormsg
+            print '  ' + fname + ": WRF file '" + ncfile + "' does not have " +      \
+              "variable 'PHB' !!"
+            quit(-1)
+
+        objph = objf.variables['PH']
+        objphb = objf.variables['PHB']
+        geop = objph[:] + objphb[:]
+
+        ijz0 = gen.index_mat(geop[0,], 0.)
+        zvals = geop[0, :, ijz0[0], ijz0[1]] / 9.8
+    else:
+        if not objf.has_variable(varn):
+            print errormsg
+            print '  ' + fname + ": file '" + ncfile + "' does not have variable '" +\
+              varn + "' !!"
+            quit(-1)
+        objvar = objf.variables[varn]
+        if len(objvar.shape) == 4:
+            print warnmsg
+            print '  ' + fname + ": assuming that variable '" + varn + "' with " +   \
+              "shape: dt,dz,dy,dx. Tacking first time-step"
+
+            ijz0 = gen.index_mat(objvar[0,0,], 0.)
+            zvals = objvar[0, :, ijz0[0], ijz0[1]]
+        elif len(objvar.shape) == 3:
+            print warnmsg
+            print '  ' + fname + ": assuming that variable '" + varn + "' with " +   \
+              "shape: dz,dy,dx"
+
+            ijz0 = gen.index_mat(objvar[0,], 0.)
+            zvals = objvar[:, ijz0[0], ijz0[1]]
+        
+        elif len(objvar.shape) == 2:
+            print warnmsg
+            print '  ' + fname + ": assuming that variable '" + varn + "' with " +   \
+              "shape: dz,dyx"
+
+            ijz0 = gen.index_mat(objvar[0,], 0.)
+            zvals = objvar[:, ijz0[0]]
+        else:
+            zvals = objvar[:]
+
+    if zlog == 'true':
+        zlogv = True
+    elif zlog == 'false':
+        zlogv = False
+    else:
+        print errormsg
+        print '  ' + fname + ": wrong value for zlog: '" + zlog + "' !!"
+        print "    must be either: 'true' or 'false'"
+        quit(-1)
+    if dzlog == 'true':
+        dzlogv = True
+    elif dzlog == 'false':
+        dzlogv = False
+    else:
+        print errormsg
+        print '  ' + fname + ": wrong value for dzlog: '" + dzlog + "' !!"
+        print "    must be either: 'true' or 'false'"
+        quit(-1)
+
+    drw.plot_vertical_lev(zvals, zlogv, dzlogv, title, kindfig, legloc)
+
+    objf.close()
+
+    return
+
 #quit()
 
@@ -3760 +3876 @@
 elif oper == 'draw_vectors':
     draw_vectors(opts.ncfile, opts.values, opts.varname)
+elif oper == 'draw_vertical_levels':
+    draw_vertical_levels(opts.ncfile, opts.values, opts.varname)
 elif oper == 'list_graphics':
 # From: http://www.diveintopython.net/power_of_introspection/all_together.html

trunk/tools/drawing_tools.py

@@ -6561 +6561 @@
     return
 
+
+def plot_vertical_lev(vertz, zlog, dzlog, gtit, kfig, lloc):
+    """ plotting vertical levels distribution
+    plot_vertical_lev(vertz, gtit, kfig, lloc)
+      vertz= distribution of vertical heights [z]
+      zlog: to use logarithmic scale on the height axis ('true/false')
+      dzlog: to use logarithmic scale on the difference of height between levels axis ('true/false')
+      gtit= title of the graph ('!' for spaces)
+      kfig= kind of figure (jpg, pdf, png)
+      lloc= location of the legend (0, automatic)
+        1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right',
+        5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
+        9: 'upper center', 10: 'center'      kfig= kind of figure
+    """
+    fname = 'plot_vertical_lev'
+
+    figname = 'vertical_lev'
+
+    if vertz == 'h':
+        print fname + '_____________________________________________________________'
+        print plot_vertical_lev.__doc__
+        quit()
+
+    Nlev = len(vertz)
+    dvertz = vertz[1:Nlev] - vertz[0:Nlev-1]
+
+# From: http://stackoverflow.com/questions/14762181/adding-a-y-axis-label-to-secondary-y-axis-in-matplotlib
+
+    plt.rc('text', usetex=True)
+
+    fig, ax1 = plt.subplots()
+    ax2 = ax1.twinx()
+
+    plt.xlim(0,Nlev)
+
+    l1 = ax1.plot(range(Nlev), vertz, 'r-x', label='height')
+    l2 = ax2.plot(range(1,Nlev), dvertz, 'b-x', label='dheight')
+
+    ax1.set_xlabel('level (\#)')
+    ax1.set_ylabel('height (m)', color='r')
+    ax1.set_ylim(1,np.max(vertz))
+    ax1.grid()
+    ax2.set_ylabel('difference between levels (m)', color='b')
+
+    plt.title(gtit)
+    if zlog: ax1.set_yscale('log')
+    if dzlog: ax2.set_yscale('log')
+
+    output_kind(kfig, figname, True)
+
+    return
+