Changeset 1872 in lmdz_wrf for trunk/tools

Timestamp:

Mar 27, 2018, 5:50:52 PM (7 years ago)

Author:

lfita

Message:

Adding:

• `draw_multi_2D_shad': plotting multiple 2D fields with same projection with shading and sharing colorbar

Location:

trunk/tools

Files:

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

trunk/tools/drawing.py

@@ -78 +78 @@
 # draw_lines: Function to draw different lines at the same time from different files
 # draw_lines_time: Function to draw different lines at the same time from different files with times
+# draw_multi_2D_shad: plotting multiple 2D fields with same projection with shading and sharing colorbar
 # draw_Neighbourghood_evol: Function to draw the temporal evolution of a neighbourghood around a point
 # draw_points: Function to plot a series of points
@@ -108 +109 @@
   'draw_2D_shad_line_time', 'draw_bar', 'draw_bar_line', 'draw_bar_line_time',       \
   'draw_bar_time', 'draw_barbs', 'draw_basins', 'draw_cycle', 'draw_ensembles_time', \
-  'draw_2lines', 'draw_2lines_time', 'draw_lines', 'draw_lines_time',                \
+  'draw_2lines', 'draw_2lines_time', 'draw_multi_2D_shad', 'draw_lines',             \
+  'draw_lines_time',                                                                 \
   'draw_Neighbourghood_evol',                                                        \
   'draw_points', 'draw_points_lonlat',                                               \
@@ -8393 +8395 @@
     return
 
+
+def draw_multi_2D_shad(ncfiles, values, axfig=None, fig=None):
+    """ plotting multiple 2D fields with same projection with shading and sharing colorbar
+    draw_multi_2D_shad(ncfile, values, varn)
+      ncfiles= ';' list of files to use as [filen]@[varn]@[dimvals]
+        [filen]: name of the file
+        [varn]: name of the variable in the file to be shaded
+        [dimvals]: ',' list of [dimname]|[value] telling at which dimension of the 
+          variable a given value is required:
+            * [integer]: which value of the dimension
+            * -1: all along the dimension
+            * -9: last value of the dimension
+            * [beg]@[end]@[inc] slice from [beg] to [end] every [inc]
+            * NOTE, no dim name all the dimension size
+      values=[vnamefs]:[dimxvn]:[dimyvn]:[dimxyfmt]:[colorbarvals]:[sminv],[smaxv]:
+       [kindfig]:[reverse]:[mapv]:[labels]:[Nrow]:[Ncol]:[globaltitle]:[close]
+        [vnamefs]: Name in the figure of the variable to be shaded
+        [dimx/yvn]: name of the variables with the values of the final dimensions (x,y)
+        [dimxyfmt]=[dxs],[dxf],[Ndx],[ordx],[dys],[dyf],[Ndy],[ordx]: format of the values at each axis (or 'auto')
+          [dxs]: style of x-axis ('auto' for 'pretty')
+            'Nfix', values computed at even 'Ndx'
+            'Vfix', values computed at even 'Ndx' increments
+            'pretty', values computed following aprox. 'Ndx' at 'pretty' intervals (2.,2.5,4,5,10)
+          [dxf]: format of the labels at the x-axis ('auto' for '%5g')
+          [Ndx]: Number of ticks at the x-axis ('auto' for 5)
+          [ordx]: angle of orientation of ticks at the x-axis ('auto' for horizontal)
+          [dys]: style of y-axis ('auto' for 'pretty')
+          [dyf]: format of the labels at the y-axis ('auto' for '%5g')
+          [Ndy]: Number of ticks at the y-axis ('auto' for 5)
+          [ordy]: angle of orientation of ticks at the y-axis ('auto' for horizontal)
+        [colorbarvals]=[colbarn],[fmtcolorbar],[orientation]
+          [colorbarn]: name of the color bar
+          [fmtcolorbar]: format of the numbers in the color bar 'C'-like ('auto' for %6g)
+          [orientation]: orientation of the colorbar ('vertical' (default, by 'auto'), 'horizontal')
+        [smin/axv]: minimum and maximum value for the shading or:
+          'Srange': for full range
+          'Saroundmean@val': for mean-xtrm,mean+xtrm where xtrm = np.min(mean-min@val,max@val-mean)
+          'Saroundminmax@val': for min*val,max*val
+          'Saroundpercentile@val': for median-xtrm,median+xtrm where xtrm = np.min(median-percentile_(val),
+            percentile_(100-val)-median)
+          'Smean@val': for -xtrm,xtrm where xtrm = np.min(mean-min*@val,max*@val-mean)
+          'Smedian@val': for -xtrm,xtrm where xtrm = np.min(median-min@val,max@val-median)
+          'Spercentile@val': for -xtrm,xtrm where xtrm = np.min(median-percentile_(val),
+             percentile_(100-val)-median)
+        [kindfig]: kind of figure output file: png, pdf, ...
+        [reverse]: Transformation of the values
+          * 'transpose': reverse the axes (x-->y, y-->x)
+          * 'flip'@[x/y]: flip the axis x or y
+        [mapv]: map characteristics: [proj],[res]
+          see full documentation: http://matplotlib.org/basemap/
+          [proj]: projection
+            * 'cyl', cilindric
+            * 'lcc', lambert conformal
+          [res]: resolution:
+            * 'c', crude
+            * 'l', low
+            * 'i', intermediate
+            * 'h', high
+            * 'f', full
+        [labels]: ',' list of consecutive labels for each individual panel ('!' for spaces)
+        [Nrow]: number of rows in the figure
+        [Ncol]: number of columns in the figure
+        [globaltitle]: global title of the figure ('!' for spaces)
+        [close]: Whether figure should be finished or not
+      valules= 'rh:z|-1,x|-1:z|-1,x|-1:lat:pressure:BuPu:0.,100.:rh:pdf:flip@y:None'
+      varn= 
+    """
+
+    fname = 'draw_multi_2D_shad'
+    if values == 'h':
+        print fname + '_____________________________________________________________'
+        print draw_2D_shad.__doc__
+        quit()
+
+    expectargs = '[vnamefs]:[dimxvn]:[dimyvn]:[dimxyf]:[colbarvals]:' +             \
+      '[sminv],[smaxv]:[kindfig]:[reverse]:[mapv]:[labels]:[Nrow]:[Ncol]:' +        \
+      '[globaltitle]:[close]'
+ 
+    drw.check_arguments(fname,values,expectargs,':')
+
+    vnamesfig = values.split(':')[0]
+    vdimxn = values.split(':')[1]
+    vdimyn = values.split(':')[2]
+    dimxyf = values.split(':')[3]
+    colorbarvals = values.split(':')[4]
+    shadminmax = values.split(':')[5]
+    figkind = values.split(':')[6]
+    revals = values.split(':')[7]
+    mapvalue = values.split(':')[8]
+    labels = values.split(':')[9].split(',')
+    Nrow = int(values.split(':')[10])
+    Ncol = int(values.split(':')[11])
+    globaltitle = values.split(':')[12].replace('!', ' ')
+    close = gen.Str_Bool(values.split(':')[13])
+
+    fvds = ncfiles.split(';')
+    lvalshad = []
+
+    for fvd in fvds:
+        ncfile = fvd.split('@')[0]
+        varn = fvd.split('@')[1]
+        dimvals = fvd.split('@')[2]
+
+        if not os.path.isfile(ncfile):
+            print errormsg
+            print '  ' + fname + ': shading file "' + ncfile + '" does not exist !!'
+            quit(-1)    
+
+        objsf = NetCDFFile(ncfile, 'r')
+    
+        varns = varn.split(',')[0]
+
+        if  not objsf.variables.has_key(varns):
+            print errormsg
+            print '  ' + fname + ': shading file "' + ncfiles +                      \
+              '" does not have variable "' +  varns + '" !!'
+            varns = sorted(objsf.variables.keys())
+            print '    available ones:', varns
+            quit(-1)
+
+# Variables' values
+        objvars = objsf.variables[varns]
+
+        valshad, dimsshad = drw.slice_variable(objvars,                              \
+          dimvals.replace('|',':').replace(',','|'))
+        lvalshad.append(valshad)
+
+        # Getting variable information from first file
+        if fvd == fvds[0]:
+            dimnamesv = [vdimxn, vdimyn]
+    
+            if drw.searchInlist(objvars.ncattrs(),'units'):
+                varunits = objvars.getncattr('units')
+            else:
+                print warnmsg
+                print '  ' + fname + ": variable '" + varn + "' without units!!"
+                varunits = '-'
+
+            if  not objsf.variables.has_key(vdimxn):
+                print errormsg
+                print '  ' + fname + ': shading file "' + ncfiles +                  \
+                  '" does not have dimension variable "' +  vdimxn + '" !!'
+                quit(-1)
+            if  not objsf.variables.has_key(vdimyn):
+                print errormsg
+                print '  ' + fname + ': shading file "' + ncfiles +                  \
+                  '" does not have dimension variable "' +  vdimyn + '" !!'
+                quit(-1)
+
+            objdimx = objsf.variables[vdimxn]
+            objdimy = objsf.variables[vdimyn]
+            if drw.searchInlist(objdimx.ncattrs(),'units'):
+                odimxu = objdimx.getncattr('units')
+            else:
+                print warnmsg
+                print '  ' + fname + ": variable dimension '" + vdimxn +             \
+                  "' without units!!"
+                odimxu = '-'
+
+            if drw.searchInlist(objdimy.ncattrs(),'units'):
+                odimyu = objdimy.getncattr('units')
+            else:
+                print warnmsg
+                print '  ' + fname + ": variable dimension '" + vdimyn +             \
+                  "' without units!!"
+                odimyu = '-'
+    
+            odimxv, odimyv = drw.dxdy_lonlatDIMS(objdimx[:], objdimy[:],             \
+              objdimx.dimensions, objdimy.dimensions,                                \
+              dimvals.replace(':','|').split(','))
+
+    shading_nx = []
+    if shadminmax.split(',')[0][0:1] != 'S':
+            shading_nx.append(np.float(shadminmax.split(',')[0]))
+    else:
+        shading_nx.append(shadminmax.split(',')[0])
+
+    if shadminmax.split(',')[1][0:1] != 'S':
+        shading_nx.append(np.float(shadminmax.split(',')[1]))
+    else:
+        shading_nx.append(shadminmax.split(',')[1])
+
+    if mapvalue == 'None': mapvalue = None
+
+    colbarn, fmtcolbar, colbaror = drw.colorbar_vals(colorbarvals,',')
+    colormapv = [colbarn, fmtcolbar, colbaror]
+
+    xstyl, xaxf, Nxax, xaxor, ystyl, yaxf, Nyax, yaxor = drw.format_axes(dimxyf,',')
+    xaxis = [xstyl, xaxf, Nxax, xaxor]
+    yaxis = [ystyl, yaxf, Nyax, yaxor]
+
+    if revals == 'None':
+        revals = None
+
+    drw.multi_plot_2D_shadow(lvalshad, odimxv, odimyv, odimxu, odimyu, xaxis, yaxis, \
+      dimnamesv, colormapv, shading_nx, varunits, figkind, revals, mapvalue, labels, \
+      Nrow, Ncol, globaltitle, close)
+
+    return
+#python drawing.py -o draw_multi_2D_shad -f '../PY/wrfout_d01_1995-01-01_00:00:00@T2@west_east|-1,south_north|-1,Time|0;../PY/wrfout_d01_1995-01-01_00:00:00@T2@west_east|-1,south_north|-1,Time|1;../PY/wrfout_d01_1995-01-01_00:00:00@T2@west_east|-1,south_north|-1,Time|2;../PY/wrfout_d01_1995-01-01_00:00:00@T2@west_east|-1,south_north|-1,Time|3' -S 'tas:XLONG:XLAT:auto:rainbow,auto,auto:Srange,Srange:png:None:cyl,l:0!UTC,1!UTC,2!UTC,3!UTC:2:2:tas!at!2001-11-11:True'
+
 #quit()
 
@@ -8416 +8619 @@
   'draw_2D_shad_line', 'draw_2D_shad_line_time', 'draw_2lines', 'draw_2lines_time',  \
   'draw_bar', 'draw_bar_line', 'draw_bar_line_time', 'draw_bar_time',                \
-  'draw_cycle', 'draw_ensembles_time', 'draw_lines',                                 \
+  'draw_cycle', 'draw_ensembles_time', 'draw_multi_2D_shad', 'draw_lines',           \
   'draw_lines_time', 'draw_points', 'draw_Taylor', 'draw_time_lag',                  \
   'draw_topo_geogrid_boxes',                                                         \
@@ -8504 +8707 @@
     elif oper == 'draw_lines_time':
         draw_lines_time(opts.ncfile, opts.values, opts.varname)
+    elif oper == 'draw_multi_2D_shad':
+        draw_multi_2D_shad(opts.ncfile, opts.values)
     elif oper == 'draw_points':
         draw_points(opts.ncfile, opts.values)

trunk/tools/drawing_tools.py

@@ -10843 +10843 @@
 
     return
+
+def multi_plot_2D_shadow(lvarsv, dimxv, dimyv, dimxu, dimyu, xaxv, yaxv, dimn,       \
+  cbarv, vs, uts, kfig, reva, mapv, labs, Ncol, Nrow, ftit, ifclose):
+    """ plot of multiple 2D shadowing variables sharing colorbar
+      lvarsv= list of 2D values to plot with shading
+      dim[x/y]v = list of values at the axes of x and y
+      dim[x/y]u = units at the axes of x and y
+      xaxv= list with the x-axis paramteres [style, format, number and orientation]
+      yaxv= list with the y-axis paramteres [style, format, number and orientation]
+      dimn= dimension names to plot
+      cbarv= list with the parameters of the color bar [colorbar, cbarfmt, cbaror]
+        colorbar: name of the color bar to use
+        cbarfmt: format of the numbers in the colorbar
+        cbaror: orientation of the colorbar
+      vs= minmum and maximum values to plot in shadow or:
+        'Srange': for full range
+        'Saroundmean@val': for mean-xtrm,mean+xtrm where xtrm = np.min(mean-min@val,max@val-mean)
+        'Saroundminmax@val': for min*val,max*val
+        'Saroundpercentile@val': for median-xtrm,median+xtrm where xtrm = np.min(median-percentile_(val),
+          percentile_(100-val)-median)
+        'Smean@val': for -xtrm,xtrm where xtrm = np.min(mean-min*@val,max*@val-mean)
+        'Smedian@val': for -xtrm,xtrm where xtrm = np.min(median-min@val,max@val-median)
+        'Spercentile@val': for -xtrm,xtrm where xtrm = np.min(median-percentile_(val),
+           percentile_(100-val)-median)
+      uts= units of the variable to shadow
+      vtit= title of the variable
+      kfig= kind of figure (jpg, pdf, png)
+      reva= ('|' for combination)
+        * 'transpose': reverse the axes (x-->y, y-->x)
+        * 'flip'@[x/y]: flip the axis x or y
+      mapv= map characteristics: [proj],[res]
+        see full documentation: http://matplotlib.org/basemap/
+        [proj]: projection
+          * 'cyl', cilindric
+          * 'lcc', lambert conformal
+        [res]: resolution:
+          * 'c', crude
+          * 'l', low
+          * 'i', intermediate
+          * 'h', high
+          * 'f', full
+      labs: list of consecutive labels for each figure
+      ftit: global title of the figure
+      Ncol: number of column panels
+      Nrow: number of rows of panels
+      ifclose= boolean value whether figure should be close (finish) or not
+    """
+##    import matplotlib as mpl
+##    mpl.use('Agg')
+##    import matplotlib.pyplot as plt
+    fname = 'multi_plot_2D_shadow'
+
+#    print dimyv[73,21]
+#    dimyv[73,21] = -dimyv[73,21]
+#    print 'Lluis dimsv: ',np.min(dimxv), np.max(dimxv), ':', np.min(dimyv), np.max(dimyv)
+
+    if type(lvarsv) == type('S') and varsv == 'h':
+        print fname + '_____________________________________________________________'
+        print multi_plot_2D_shadow.__doc__
+        quit()
+ 
+    # Getting the right lon values for plotting
+    if mapv is not None:
+        dimxv = np.where(dimxv > 180., dimxv-360., dimxv)
+
+    # Axis ticks
+    # Usually axis > x must be the lon, thus...
+    dimxv0 = dimxv.copy()
+    dimyv0 = dimyv.copy()
+
+    dxn = dimxv.min()
+    dxx = dimxv.max()
+    dyn = dimyv.min()
+    dyx = dimyv.max()
+
+    if xaxv[0] == 'pretty':
+        dimxt0 = np.array(gen.pretty_int(dxn,dxx,xaxv[2]))
+    elif xaxv[0] == 'Nfix':
+        dimxt0 = np.arange(dxn,dxx,(dxx-dxn)/(1.*xaxv[2]))
+    elif xaxv[0] == 'Vfix':
+        dimxt0 = np.arange(0,dxx,xaxv[2])
+    if yaxv[0] == 'pretty':
+        dimyt0 = np.array(gen.pretty_int(dyn,dyx,yaxv[2]))
+    elif yaxv[0] == 'Nfix':
+        dimyt0 = np.arange(dyn,dyx,(dyx-dyn)/(1.*yaxv[2]))
+    elif yaxv[0] == 'Vfix':
+        dimyt0 = np.arange(0,dyx,yaxv[2])
+
+    dimxl0 = []
+    for i in range(len(dimxt0)): dimxl0.append('{:{style}}'.format(dimxt0[i], style=xaxv[1]))
+    dimyl0 = []
+    for i in range(len(dimyt0)): dimyl0.append('{:{style}}'.format(dimyt0[i], style=yaxv[1]))
+
+    dimxT0 = variables_values(dimn[0])[0] + ' (' + units_lunits(dimxu) + ')'
+    dimyT0 = variables_values(dimn[1])[0] + ' (' + units_lunits(dimyu) + ')'
+
+    if mapv is not None:
+        pixkind = 'data'
+    else:
+        # No following data values
+        dimxt0 = np.arange(len(dimxt0),dtype=np.float)/(len(dimxt0))
+        dimyt0 = np.arange(len(dimyt0),dtype=np.float)/(len(dimyt0))
+        pixkind = 'fixpixel'
+
+    if reva is not None:
+        varsv, dimxv, dimyv, dimxt, dimyt, dimxl, dimyl, dimxT, dimyT =              \
+          transform(varsv, reva, dxv=dimxv0, dyv=dimyv0, dxt=dimxt0, dyt=dimyt0,     \
+          dxl=dimxl0, dyl=dimyl0, dxtit=dimxT0, dytit=dimyT0)
+    else:
+        dimxv = dimxv0
+        dimyv = dimyv0
+        dimxt = dimxt0
+        dimyt = dimyt0
+        dimxl = dimxl0
+        dimyl = dimyl0
+        dimxT = dimxT0
+        dimyT = dimyT0
+
+    if len(dimxv[:].shape) == 3:
+        xdims = '1,2'
+    elif len(dimxv[:].shape) == 2:
+        xdims = '0,1'
+    elif len(dimxv[:].shape) == 1:
+        xdims = '0'
+
+    if len(dimyv[:].shape) == 3:
+        ydims = '1,2'
+    elif len(dimyv[:].shape) == 2:
+        ydims = '0,1'
+    elif len(dimyv[:].shape) == 1:
+        ydims = '0'
+
+    lon0, lat0 = dxdy_lonlat(dimxv,dimyv, xdims, ydims)
+
+    if not mapv is None:
+        map_proj=mapv.split(',')[0]
+        map_res=mapv.split(',')[1]
+
+        dx = lon0.shape[1]
+        dy = lat0.shape[0]
+
+        nlon = np.min(lon0)
+        xlon = np.max(lon0)
+        nlat = np.min(lat0)
+        xlat = np.max(lat0)
+
+        lon2 = lon0[dy/2,dx/2]
+        lat2 = lat0[dy/2,dx/2]
+
+        print 'lon2:', lon2, 'lat2:', lat2, 'SW pt:', nlon, ',', nlat, 'NE pt:',     \
+          xlon, ',', xlat
+
+        if map_proj == 'cyl':
+            m = Basemap(projection=map_proj, llcrnrlon=nlon, llcrnrlat=nlat,         \
+              urcrnrlon=xlon, urcrnrlat= xlat, resolution=map_res)
+        elif map_proj == 'lcc':
+            m = Basemap(projection=map_proj, lat_0=lat2, lon_0=lon2, llcrnrlon=nlon, \
+              llcrnrlat=nlat, urcrnrlon=xlon, urcrnrlat= xlat, resolution=map_res)
+        else:
+            print errormsg
+            print '  ' + fname + ": map projection '" + map_proj + "' not defined!!!"
+            print '    available: cyl, lcc'
+            quit(-1)
+ 
+        x,y = m(lon0,lat0)
+        if cbarv[0] == 'gist_gray':
+            m.drawcoastlines(color="red")
+        else:
+            m.drawcoastlines()
+        meridians = gen.pretty_int(nlon,xlon,xaxv[2])
+        m.drawmeridians(meridians,labels=[True,False,False,True])
+        parallels = gen.pretty_int(nlat,xlat,yaxv[2])
+        m.drawparallels(parallels,labels=[False,True,True,False])
+
+    else:
+        # No following data values
+        x = (dimxv-np.min(dimxv))/(np.max(dimxv) - np.min(dimxv))
+        y = (dimyv-np.min(dimyv))/(np.max(dimyv) - np.min(dimyv))
+
+# Changing limits of the colors
+    vsend = graphic_range(vs,lvarsv)
+
+    # Number of panels
+    Npanels = len(lvarsv)
+    if Ncol*Nrow != Npanels:
+        print errormsg
+        print '  ' + fname + ': no coincident number of panels:', Npanels,           \
+          'with provied distribution of them:', Ncol, Nrow
+
+    plt.rc('text', usetex=True)
+    fig, axes = plt.subplots(Nrow, Ncol)
+
+    irow=0
+    icol=0
+    for iv in range(Npanels):
+        varsv = lvarsv[iv]
+
+        if len(varsv.shape) != 2:
+            print errormsg
+            print '  ' + fname + ': wrong variable shadow rank:', varsv.shape,       \
+              'is has to be 2D!!'
+            quit(-1)
+        axes[irow,icol] = plt.subplot(Nrow ,Ncol, iv+1)
+        im = plt.pcolormesh(x, y, varsv, cmap=plt.get_cmap(cbarv[0]), vmin=vsend[0], \
+          vmax=vsend[1])
+    
+        if not mapv is None:
+            m.drawcoastlines()
+            m.drawmeridians(meridians,labels=[True,False,False,True])
+            m.drawparallels(parallels,labels=[False,True,True,False])
+            if irow == Nrow -1: plt.xlabel('W-E')
+            if icol == 0: plt.ylabel('S-N')
+        else:
+            plt.xlabel(dimxT)
+            plt.ylabel(dimyT)
+            plt.xticks(dimxt, dimxl, rotation=xaxv[3])
+            plt.yticks(dimyt, dimyl, rotation=yaxv[3])
+
+        plt.axis([x.min(), x.max(), y.min(), y.max()])
+
+        graphtit = gen.latex_text(labs[iv].replace('!', ' '))
+        plt.title(gen.latex_text(graphtit))
+        icol = icol + 1
+        if icol > Ncol-1:
+            irow = irow + 1
+            icol = 0
+
+    fig.subplots_adjust(right=0.8, hspace=0.4)
+    cbar_ax = fig.add_axes([0.86, 0.25, 0.03, 0.50])
+
+    if cbarv[2] == 'horizontal':
+        cbar = fig.colorbar(im, cax=cbar_ax, format=cbarv[1],orientation=cbarv[2])
+        # From: http://stackoverflow.com/questions/32050030/rotation-of-colorbar-tick-labels-in-matplotlib
+        ticklabels= cbar.ax.get_xticklabels()
+        Nticks = len(ticklabels)
+        ticklabs = []
+        for itick in range(Nticks): ticklabs.append(ticklabels[itick].get_text())
+        cbar.ax.set_xticklabels(ticklabs,rotation=90)
+    else:
+        cbar = fig.colorbar(im, cax=cbar_ax, format=cbarv[1],orientation=cbarv[2])
+
+    #fig.colorbar(im, cax=cbar_ax)
+# units labels
+    ulab = units_lunits(uts)
+    cbar_ax.set_label(ulab)
+    plt.annotate(ulab, xy=(0.96,0.5), xycoords='figure fraction', rotation=90)
+    cbar_ax.tick_params(labelsize=8)
+
+    plt.suptitle(ftit)
+    figname = 'multi_2Dfields_shadow'
+
+    output_kind(kfig, figname, ifclose)
+
+    return
+