diff --git a/src/LIB/Python/Panel_Plot.py b/src/LIB/Python/Panel_Plot.py index 2072d1cf9225c27f14b79849e8b11e5e49885ce8..e4588ddde7852b734d0972de097edcbf41bf894f 100644 --- a/src/LIB/Python/Panel_Plot.py +++ b/src/LIB/Python/Panel_Plot.py @@ -9,7 +9,7 @@ MNH_LIC for details. version 1. @author: 07/2021 Quentin Rodier """ import matplotlib as mpl -#mpl.use('Agg') +# mpl.use('Agg') import matplotlib.pyplot as plt from matplotlib import cm from matplotlib.colors import ListedColormap @@ -17,9 +17,10 @@ import numpy as np import cartopy import cartopy.feature as cfeature + class PanelPlot(): - - def __init__(self, nb_l, nb_c, Lfigsize, bigtitle, titlepad=40, minmaxpad=1.03, timepad=-0.06, lateralminmaxpad=0.86, + + def __init__(self, nb_l, nb_c, Lfigsize, bigtitle, bigtitlepad=0.95, titlepad=40, minmaxpad=1.03, timepad=-0.06, lateralminmaxpad=0.86, labelcolorbarpad=6.0, colorbaraspect=20, colorbarpad=0.04, tickspad=0.8, minmaxTextSize=10, bigtitleSize=13, titleSize=12, legendSize=10, xlabelSize=11, ylabelSize=11, timeSize=11, cbTicksLabelSize=11, cbTitleSize=11, xyTicksLabelSize=10, figBoxLinewidth=1, @@ -31,6 +32,7 @@ class PanelPlot(): self.nb_c = nb_c # Panel number of rows self.nb_graph = 0 # New independent graph within the subplot + self.bigtitlepad = bigtitlepad # Panel title vertical position self.titlepad = titlepad # Title pad (vertical shift) from graph self.tickspad = tickspad # Ticks pad (between ticks and axis label) self.minmaxpad = minmaxpad # Min/Max print pad (vertical shift) @@ -56,873 +58,1007 @@ class PanelPlot(): self.xyTicksLength = xyTicksLength # Ticks length # Initialization of the panel plots - self.fig = plt.figure(figsize=(self.Lfigsize[0],self.Lfigsize[1])) + self.fig = plt.figure(figsize=(self.Lfigsize[0], self.Lfigsize[1])) self.fig.set_dpi(125) - self.fig.suptitle(self.bigtitle,fontsize=bigtitleSize) + self.fig.suptitle(self.bigtitle, fontsize=self.bigtitleSize, y=self.bigtitlepad) def save_graph(self, iplt, fig, fig_name='tempgraph'): - """ - Create a temporary png file of the panel plot which can be converted to PDF - """ - self.iplt = iplt - self.fig = fig - self.fig_name=fig_name # .png figure prefix name - - self.fig.savefig(self.fig_name+str(self.iplt)) - self.iplt+=1 - return self.iplt - - def draw_Backmap(self,drawCoastLines, ax, projo): - """ - Handle drawing of the background plot (coastlines, departements, grid lines and labels) - """ - self.drawCoastLines = drawCoastLines - self.projo = projo - - # Grid lines and labels - if 'PlateCarree' in str(projo): + """ + Create a temporary png file of the panel plot which can be converted to PDF + """ + self.iplt = iplt + self.fig = fig + self.fig_name = fig_name # .png figure prefix name + + self.fig.savefig(self.fig_name + str(self.iplt)) + self.iplt += 1 + return self.iplt + + def draw_Backmap(self, drawCoastLines, ax, projo): + """ + Handle drawing of the background plot (coastlines, departements, grid lines and labels) + """ + self.drawCoastLines = drawCoastLines + self.projo = projo + + # Grid lines and labels gl = ax.gridlines(crs=self.projo, draw_labels=True, linewidth=1, color='gray') if float(cartopy.__version__[:4]) >= 0.18: - gl.top_labels = False - gl.right_labels = False + from cartopy.mpl.ticker import (LatitudeLocator, LongitudeLocator,LongitudeFormatter, LatitudeFormatter) + gl.top_labels = False + gl.right_labels = False + gl.xlines = True + gl.ylines = True + gl.xlocator = LongitudeLocator() + gl.ylocator = LatitudeLocator() + gl.xformatter = LongitudeFormatter() + gl.yformatter = LatitudeFormatter() else: - gl.xlabels_top = False - gl.ylabels_right = False - - # Coastlines - if self.drawCoastLines and 'GeoAxes' in str(type(ax)): - ax.coastlines(resolution='10m') - - # Countries border - ax.add_feature(cfeature.BORDERS) - ax.add_feature(cfeature.LAKES, alpha=0.7) - - def addWhitecm(self, colormap_in, nb_level,whiteTop=False): + gl.xlabels_top = False + gl.ylabels_right = False + gl.xlabel_style = {'size': self.xyTicksLabelSize, 'color': 'black'} + gl.ylabel_style = {'size': self.xyTicksLabelSize, 'color': 'black'} + + # Coastlines + if self.drawCoastLines and 'GeoAxes' in str(type(ax)): + ax.coastlines(resolution='10m', color='black', linewidth=1) + + # Countries border + ax.add_feature(cfeature.BORDERS) + ax.add_feature(cfeature.LAKES, alpha=0.7) + + def addWhitecm(self, colormap_in, nb_level, whiteTop=False): """ Add a white color at the top (whiteTop=True) or bottom of the colormap w.r.t. the number of independent colors used """ color_map = cm.get_cmap(colormap_in, 256) newcolor_map = color_map(np.linspace(0, 1, 256)) - whites = np.array([1, 1, 1, 1]) #RBG code + opacity + whites = np.array([1, 1, 1, 1]) # RBG code + opacity if whiteTop: - for i in range(int(256/nb_level)): newcolor_map[-i, :] = whites + for i in range(int(256 / nb_level)): + newcolor_map[-i, :] = whites else: - for i in range(int(256/nb_level)): newcolor_map[:i, :] = whites + for i in range(int(256 / nb_level)): + newcolor_map[:i, :] = whites newcmp = ListedColormap(newcolor_map) return newcmp - def set_Title(self, ax, i, title, Lid_overlap, xlab, ylab): - """ - Handle top title of each graph - Parameters : - - titlepad : global attribute for vertical shift placement w.r.t the graph - """ - self.ax = ax - self.title = title - self.Lid_overlap = Lid_overlap - self.i = i - #self.ax[self.i].set_xlabel("test", fontweight='bold') - if not Lid_overlap: - self.ax[self.i].set_title(title, pad=self.titlepad, fontsize=self.titleSize) - else: # If graph overlap, title is concatenated - new_title = self.ax[self.i].get_title() + ' and ' + title - self.ax[self.i].set_title(new_title, pad=self.titlepad, fontsize=self.titleSize) + """ + Handle top title of each graph + Parameters : + - titlepad : global attribute for vertical shift placement w.r.t the graph + """ + self.ax = ax + self.title = title + self.Lid_overlap = Lid_overlap + self.i = i + #self.ax[self.i].set_xlabel("test", fontweight='bold') + if not Lid_overlap: + self.ax[self.i].set_title(title, pad=self.titlepad, fontsize=self.titleSize) + else: # If graph overlap, title is concatenated + new_title = self.ax[self.i].get_title() + ' and ' + title + self.ax[self.i].set_title(new_title, pad=self.titlepad, fontsize=self.titleSize) def set_xlim(self, ax, i, xlim): - """ - Handle x limits plotted if necessary - """ - self.ax = ax - self.xlim = xlim - self.i = i - - self.ax[self.i].set_xlim(xlim[0],xlim[1])#, fontsize=self.xlabelSize) - + """ + Handle x limits plotted if necessary + """ + self.ax = ax + self.xlim = xlim + self.i = i + + self.ax[self.i].set_xlim(xlim[0], xlim[1]) # , fontsize=self.xlabelSize) + def set_ylim(self, ax, i, ylim): - """ - Handle x limits plotted if necessary - """ - self.ax = ax - self.ylim = ylim - self.i = i - - self.ax[self.i].set_ylim(ylim[0],ylim[1])#, fontsize=self.ylabelSize) + """ + Handle x limits plotted if necessary + """ + self.ax = ax + self.ylim = ylim + self.i = i + + self.ax[self.i].set_ylim(ylim[0], ylim[1]) # , fontsize=self.ylabelSize) def set_XYaxislab(self, ax, i, xlab, ylab): - """ - Handle x and y axis labels - """ - self.ax = ax - self.xlab = xlab - self.ylab = ylab - self.i = i - - # This handing label is a known issue with GeoAxes of cartopy - # https://stackoverflow.com/questions/35479508/cartopy-set-xlabel-set-ylabel-not-ticklabels - # https://github.com/SciTools/cartopy/issues/1332 - if 'GeoAxes' in str(type(self.ax[self.i])): - self.ax[self.i].text(-0.11, 0.45, ylab, verticalalignment='top', horizontalalignment='left', - rotation='vertical', rotation_mode='anchor', transform=self.ax[self.i].transAxes, color='black', fontsize=self.ylabelSize) - self.ax[self.i].text(0.45, -0.06, xlab, verticalalignment='top', horizontalalignment='left', - rotation='horizontal', rotation_mode='anchor', transform=self.ax[self.i].transAxes, color='black', fontsize=self.xlabelSize) - else: - self.ax[self.i].set_xlabel(xlab, fontsize=self.xlabelSize, labelpad=0.1) - self.ax[self.i].set_ylabel(ylab, fontsize=self.ylabelSize, labelpad=0.1) - + """ + Handle x and y axis labels + """ + self.ax = ax + self.xlab = xlab + self.ylab = ylab + self.i = i + + # This handing label is a known issue with GeoAxes of cartopy + # https://stackoverflow.com/questions/35479508/cartopy-set-xlabel-set-ylabel-not-ticklabels + # https://github.com/SciTools/cartopy/issues/1332 + if 'GeoAxes' in str(type(self.ax[self.i])): + self.ax[self.i].text(-0.11, 0.45, ylab, verticalalignment='top', horizontalalignment='left', + rotation='vertical', rotation_mode='anchor', transform=self.ax[self.i].transAxes, color='black', fontsize=self.ylabelSize) + self.ax[self.i].text(0.45, -0.06, xlab, verticalalignment='top', horizontalalignment='left', + rotation='horizontal', rotation_mode='anchor', transform=self.ax[self.i].transAxes, color='black', fontsize=self.xlabelSize) + else: + self.ax[self.i].set_xlabel(xlab, fontsize=self.xlabelSize, labelpad=0.1) + self.ax[self.i].set_ylabel(ylab, fontsize=self.ylabelSize, labelpad=0.1) + def addLine(self, ax, beg_coord, end_coord, color='black', linewidth=0.2): - self.ax = ax - self.beg_coord = beg_coord - self.end_coord = end_coord - self.color = color - self.linewidth = linewidth - - x1, y1 = [self.beg_coord[0],self.end_coord[0]], [self.beg_coord[1], self.end_coord[1]] - ax.plot(x1, y1,color=self.color,linewidth=self.linewidth) + self.ax = ax + self.beg_coord = beg_coord + self.end_coord = end_coord + self.color = color + self.linewidth = linewidth + + x1, y1 = [self.beg_coord[0], self.end_coord[0]], [self.beg_coord[1], self.end_coord[1]] + ax.plot(x1, y1, color=self.color, linewidth=self.linewidth) def print_minmax(self, var, title): print(str(title) + " min = " + str(np.nanmin(var)) + " max = " + str(np.nanmax(var))) def set_minmaxText(self, ax, i, var, title, Lid_overlap, facconv): - """ - Show min and max value Text in the plot - If overlap variable, text is align to the right - TODO : handle more than 2 overlap variables - """ - self.ax = ax - self.var = var - self.i = i - self.title = title - self.facconv = facconv - - strtext = " min = " + "{:.3e}".format(np.nanmin(var*facconv)) + " max = " + "{:.3e}".format(np.nanmax(var*facconv)) - if not Lid_overlap: - self.ax[self.i].text(0.01, self.minmaxpad, strtext, verticalalignment='top', horizontalalignment='left', - transform=self.ax[self.i].transAxes, color='black', fontsize=self.minmaxTextSize) - else: - self.ax[self.i].text(self.lateralminmaxpad, self.minmaxpad, strtext, verticalalignment='top', horizontalalignment='right', - transform=self.ax[self.i].transAxes, color='black', fontsize=self.minmaxTextSize) - # Print to help choose min/max value for ticks - self.print_minmax(var*facconv, title) - + """ + Show min and max value Text in the plot + If overlap variable, text is align to the right + TODO : handle more than 2 overlap variables + """ + self.ax = ax + self.var = var + self.i = i + self.title = title + self.facconv = facconv + + strtext = " min = " + "{:.3e}".format(np.nanmin(var * facconv)) + " max = " + "{:.3e}".format(np.nanmax(var * facconv)) + if not Lid_overlap: + self.ax[self.i].text(0.01, self.minmaxpad, strtext, verticalalignment='top', horizontalalignment='left', + transform=self.ax[self.i].transAxes, color='black', fontsize=self.minmaxTextSize) + else: + self.ax[self.i].text(self.lateralminmaxpad, self.minmaxpad, strtext, verticalalignment='top', horizontalalignment='right', + transform=self.ax[self.i].transAxes, color='black', fontsize=self.minmaxTextSize) + # Print to help choose min/max value for ticks + self.print_minmax(var * facconv, title) + def showTimeText(self, ax, i, timetxt): - """ - Show time validity - """ - self.ax = ax - self.i = i - self.timetxt = timetxt - - strtext = "Time = " + timetxt - self.ax[self.i].text(0.0, self.timepad, strtext, verticalalignment='top', horizontalalignment='left', - transform=self.ax[self.i].transAxes, color='black', fontsize=self.timeSize) - - - def psectionV(self, Lxx=[], Lzz=[], Lvar=[], Lxlab=[], Lylab=[], Ltitle=[], Lminval=[], Lmaxval=[], - Lstep=[], Lstepticks=[], Lcolormap=[], Lcbarlabel=[], LcolorLine=[], - Lfacconv=[], ax=[], Lid_overlap=[], colorbar=True, orog=[], Lxlim=[], Lylim=[], Ltime=[], Lpltype=[], LaddWhite_cm=[], LwhiteTop=[]): - """ - Vertical cross section plot - Parameters : - - Lxx : List of x or y coordinate variable or time axis - - Lzz : List of z coordinates variable - - Lvar : List of variables to plot - - Lxlab : List of x-axis label - - Lylab : List of y-axis label - - Lxlim : List of x (min, max) value plotted - - Lylim : List of y (min, max) value plotted - - Ltime : List of time (validity) - - Ltitle : List of sub-title - - Lminval: List of minimum value for each colorbar - - Lmaxval: List of maximum value for each colorbar - - Lstep : List of color-steps for each colorbar - - Lstepticks : List of value of labels for each colorbar - - Lcolormap : List of colormap - - LcolorLine : List of colors for colors arg of contour (color line only) - - Lcbarlabel : List of colorbar label legend (units) - - Lfacconv : List of factors for unit conversion of each variables - - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel - - Lid_overlap: List of number index of plot to overlap current variables - - Lpltype : List of types of plot 'cf' or 'c'. cf=contourf, c=contour (lines only) - - colorbar : show colorbar or not - - LaddWhite_cm : List of boolean to add white color to a colormap at the last bottom (low value) tick colorbar - - LwhiteTop : List of boolean to add the white color at the first top (high value). If false, the white is added at the bottom if Laddwhite_cm=T - - orog : Orography variable - """ - self.ax = ax - firstCall = (len(self.ax) == 0) - - # Initialize default value w.r.t to the number of plots -# D={'Lxlab':Lxlab, 'Lylab':Lylab, 'Ltitle':Ltitle,'Lminval':Lminval, 'Lmaxval':Lmaxval, + """ + Show time validity + """ + self.ax = ax + self.i = i + self.timetxt = timetxt + + strtext = "Time = " + timetxt + self.ax[self.i].text(0.0, self.timepad, strtext, verticalalignment='top', horizontalalignment='left', + transform=self.ax[self.i].transAxes, color='black', fontsize=self.timeSize) + + def psectionV(self, Lxx=[], Lzz=[], Lvar=[], Lxlab=[], Lylab=[], Ltitle=[], Lminval=[], Lmaxval=[], + Lstep=[], Lstepticks=[], Lcolormap=[], Lcbarlabel=[], LcolorLine=[], Lcbformatlabel=[], Llinewidth=[], + Lfacconv=[], ax=[], Lid_overlap=[], colorbar=True, orog=[], Lxlim=[], Lylim=[], Ltime=[], Lpltype=[], LaddWhite_cm=[], LwhiteTop=[]): + """ + Vertical cross section plot + Parameters : + - Lxx : List of x or y coordinate variable or time axis + - Lzz : List of z coordinates variable + - Lvar : List of variables to plot + - Lxlab : List of x-axis label + - Lylab : List of y-axis label + - Lxlim : List of x (min, max) value plotted + - Lylim : List of y (min, max) value plotted + - Ltime : List of time (validity) + - Ltitle : List of sub-title + - Lminval: List of minimum value for each colorbar + - Lmaxval: List of maximum value for each colorbar + - Lstep : List of color-steps for each colorbar + - Lstepticks : List of value of labels for each colorbar + - Lcolormap : List of colormap + - Llinewidth : List of lines thickness of contour + - LcolorLine : List of colors for colors arg of contour (color line only) + - Lcbarlabel : List of colorbar label legend (units) + - Lfacconv : List of factors for unit conversion of each variables + - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel + - Lid_overlap: List of number index of plot to overlap current variables + - Lpltype : List of types of plot 'cf' or 'c'. cf=contourf, c=contour (lines only) + - colorbar : show colorbar or not + - LaddWhite_cm : List of boolean to add white color to a colormap at the last bottom (low value) tick colorbar + - LwhiteTop : List of boolean to add the white color at the first top (high value). If false, the white is added at the bottom if Laddwhite_cm=T + - Lcbformatlabel: List of boolean to reduce the format to exponential 1.1E+02 format colorbar label + - orog : Orography variable + """ + self.ax = ax + firstCall = (len(self.ax) == 0) + + # Initialize default value w.r.t to the number of plots +# D={'Lxlab':Lxlab, 'Lylab':Lylab, 'Ltitle':Ltitle,'Lminval':Lminval, 'Lmaxval':Lmaxval, # 'Lstep':Lstep, 'Lstepticks':Lstepticks, 'Lcolormap':Lcolormap, 'Lcbarlabel':Lcbarlabel, 'Lfacconv':Lfacconv, 'Ltime':Ltime, # 'LaddWhite_cm':LaddWhite_cm, 'Lpltype':Lpltype} # D = initialize_default_val(Lvar, D) - # Default values - if not Lfacconv: Lfacconv = [1.0]*len(Lvar) - if not Lcolormap and not LcolorLine: Lcolormap=['gist_rainbow_r']*len(Lvar) # If no color given, a cmap is given - if not Lcolormap: LcolorLine=['black']*len(Lvar) - if not Lpltype: Lpltype=['cf']*len(Lvar) - if not LaddWhite_cm: LaddWhite_cm=[False]*len(Lvar) - if not LwhiteTop: LwhiteTop=[False]*len(Lvar) - if not Lylab: Lylab = ['']*len(Lvar) - if not Lxlab: Lxlab = ['']*len(Lvar) - # Add an extra percentage of the top max value for forcing the colorbar show the true user maximum value (correct a bug) - Lmaxval = list(map(lambda x, y: x + 1E-6*y, Lmaxval, Lstep) ) #The extra value is 1E-6 times the step ticks of the colorbar - - # On all variables to plot - for i,var in enumerate(Lvar): - if firstCall: #1st call - iax = i - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1)) - self.nb_graph+=1 - elif Lid_overlap != []: #overlapping plot - iax = Lid_overlap[i] - else: #existing ax with no overlapping (graphd appended to existing panel) - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1)) - self.nb_graph+=1 - iax = len(self.ax)-1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter - - # Colors normalization - norm = mpl.colors.Normalize(vmax=Lmaxval[i], vmin=Lminval[i]) - - # Print min/max (stout and on plot) - self.set_minmaxText(self.ax, iax, var, Ltitle[i], Lid_overlap, Lfacconv[i]) - - # Print time validity - if Ltime: self.showTimeText(self.ax, iax, str(Ltime[i])) - - # Number of contours level - if not Lstep[i]: # Default value of number of steps is 20 - Lstep[i] = (Lmaxval[i] - Lminval[i])/20 - Lstepticks[i] = Lstep[i] - - levels_contour = np.arange(Lminval[i],Lmaxval[i],step=Lstep[i]) - - # Add White to colormap - if LaddWhite_cm[i] and Lcolormap: Lcolormap[i]=self.addWhitecm(Lcolormap[i], len(levels_contour), LwhiteTop[i]) - - # Plot - if Lpltype[i]=='c': # Contour - if LcolorLine: - cf = self.ax[iax].contour(Lxx[i], Lzz[i], var*Lfacconv[i], levels=levels_contour, - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], colors=LcolorLine[i], linewidths=0.1) - else: - cf = self.ax[iax].contour(Lxx[i], Lzz[i], var*Lfacconv[i], levels=levels_contour, - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i], linewidths=0.1) - else: # Contourf - cf = self.ax[iax].contourf(Lxx[i], Lzz[i], var*Lfacconv[i], levels=levels_contour, - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) - - # Title - self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap,Lxlab[i], Lylab[i]) - - # X/Y Axis label - self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) - - # Ticks label - self.ax[iax].tick_params(axis='both', labelsize=self.xyTicksLabelSize, width=self.xyTicksWidth, length=self.xyTicksLength, pad=self.tickspad) - - # Bounding box of the plot line width - for axis in ['top','bottom','left','right']: - self.ax[iax].spines[axis].set_linewidth(self.figBoxLinewidth) - - # X/Y Axis limits value - if Lxlim: - try: - self.set_xlim(self.ax, iax, Lxlim[i]) - except: - pass - if Lylim: - try: - self.set_ylim(self.ax, iax, Lylim[i]) - except: - pass - - # Color label on contour-line - if Lpltype[i]=='c': # Contour - self.ax[iax].clabel(cf, fontsize=self.cbTicksLabelSize) - #self.ax[iax].clabel(cf, levels=np.arange(Lminval[i],Lmaxval[i],step=Lstep[i]), fontsize=self.cbTicksLabelSize) #TODO bug, levels not recognized - - #Filling area under topography - if not orog==[]: - self.ax[iax].fill_between(Lxx[i][0,:], orog, color='black', linewidth=0.2) - - # Colorbar - if colorbar: - cb=plt.colorbar(cf, ax=self.ax[iax], fraction=0.031, pad=self.colorbarpad, ticks=np.arange(Lminval[i],Lmaxval[i],Lstepticks[i]), aspect=self.colorbaraspect) - cb.ax.set_title(Lcbarlabel[i], pad=self.labelcolorbarpad, loc='left', fontsize=self.cbTitleSize) # This creates a new AxesSubplot only for the colorbar y=0 ==> location at the bottom - - return self.fig + # Default values + if not Lfacconv: + Lfacconv = [1.0] * len(Lvar) + if not Lcolormap and not LcolorLine: + Lcolormap = ['gist_rainbow_r'] * len(Lvar) # If no color given, a cmap is given + if not Lcolormap: + LcolorLine = ['black'] * len(Lvar) + if not Lpltype: + Lpltype = ['cf'] * len(Lvar) + if not LaddWhite_cm: + LaddWhite_cm = [False] * len(Lvar) + if not LwhiteTop: + LwhiteTop = [False] * len(Lvar) + if not Lylab: + Lylab = [''] * len(Lvar) + if not Lxlab: + Lxlab = [''] * len(Lvar) + if not Lcbformatlabel: + Lcbformatlabel = [False] * len(Lvar) + if not Llinewidth: + Llinewidth = [1.0] * len(Lvar) + + # Add an extra percentage of the top max value for forcing the colorbar show the true user maximum value (correct a bug) + Lmaxval = list(map(lambda x, y: x + 1E-6 * y, Lmaxval, Lstep)) # The extra value is 1E-6 times the step ticks of the colorbar + + # On all variables to plot + for i, var in enumerate(Lvar): + if firstCall: # 1st call + iax = i + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1)) + self.nb_graph += 1 + elif Lid_overlap != []: # overlapping plot + iax = Lid_overlap[i] + else: # existing ax with no overlapping (graphd appended to existing panel) + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1)) + self.nb_graph += 1 + iax = len(self.ax) - 1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter + + # Colors normalization + norm = mpl.colors.Normalize(vmax=Lmaxval[i], vmin=Lminval[i]) + + # Print min/max (stout and on plot) + self.set_minmaxText(self.ax, iax, var, Ltitle[i], Lid_overlap, Lfacconv[i]) + + # Print time validity + if Ltime: + self.showTimeText(self.ax, iax, str(Ltime[i])) + + # Number of contours level + if not Lstep[i]: # Default value of number of steps is 20 + Lstep[i] = (Lmaxval[i] - Lminval[i]) / 20 + Lstepticks[i] = Lstep[i] + + levels_contour = np.arange(Lminval[i], Lmaxval[i], step=Lstep[i]) + + # Add White to colormap + if LaddWhite_cm[i] and Lcolormap: + Lcolormap[i] = self.addWhitecm(Lcolormap[i], len(levels_contour), LwhiteTop[i]) + + # Plot + if Lpltype[i] == 'c': # Contour + if LcolorLine: + cf = self.ax[iax].contour(Lxx[i], Lzz[i], var * Lfacconv[i], levels=levels_contour, + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], colors=LcolorLine[i], linewidths=Llinewidth[i]) + else: + cf = self.ax[iax].contour(Lxx[i], Lzz[i], var * Lfacconv[i], levels=levels_contour, + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i], linewidths=Llinewidth[i]) + else: # Contourf + cf = self.ax[iax].contourf(Lxx[i], Lzz[i], var * Lfacconv[i], levels=levels_contour, + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) + + # Title + self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap, Lxlab[i], Lylab[i]) + + # X/Y Axis label + self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) + + # Ticks label + self.ax[iax].tick_params( + axis='both', + labelsize=self.xyTicksLabelSize, + width=self.xyTicksWidth, + length=self.xyTicksLength, + pad=self.tickspad) + + # Bounding box of the plot line width + for axis in ['top', 'bottom', 'left', 'right']: + self.ax[iax].spines[axis].set_linewidth(self.figBoxLinewidth) + + # X/Y Axis limits value + if Lxlim: + try: + self.set_xlim(self.ax, iax, Lxlim[i]) + except BaseException: + pass + if Lylim: + try: + self.set_ylim(self.ax, iax, Lylim[i]) + except BaseException: + pass + + # Color label on contour-line + if Lpltype[i] == 'c': # Contour + self.ax[iax].clabel(cf, fontsize=self.cbTicksLabelSize) + # self.ax[iax].clabel(cf, + # levels=np.arange(Lminval[i],Lmaxval[i],step=Lstep[i]), + # fontsize=self.cbTicksLabelSize) #TODO bug, levels not recognized + + # Filling area under topography + if not orog == []: + if Lxx[i].ndim == 1: + self.ax[iax].fill_between(Lxx[i], orog, color='black', linewidth=0.2) + else: + self.ax[iax].fill_between(Lxx[i][0, :], orog, color='black', linewidth=0.2) + + # Colorbar + if colorbar: + cb = plt.colorbar( + cf, + ax=self.ax[iax], + fraction=0.031, + pad=self.colorbarpad, + ticks=np.arange( + Lminval[i], + Lmaxval[i], + Lstepticks[i]), + aspect=self.colorbaraspect) + cb.ax.tick_params(labelsize=self.cbTicksLabelSize) + # This creates a new AxesSubplot only for the colorbar y=0 ==> location at the bottom + cb.ax.set_title(Lcbarlabel[i], pad=self.labelcolorbarpad, loc='left', fontsize=self.cbTitleSize) + if Lcbformatlabel[i]: + cb.ax.set_yticklabels(["{:.1E}".format(i) for i in cb.get_ticks()]) + + return self.fig def pXY_lines(self, Lxx=[], Lyy=[], Lxlab=[], Lylab=[], Ltitle=[], Llinetype=[], Llinewidth=[], - Llinecolor=[], Llinelabel=[], LfacconvX=[], LfacconvY=[], ax=[], id_overlap=None, Lxlim=[], - Lylim=[], Ltime=[], LaxisColor=[]): - """ - XY (multiple)-lines plot - Parameters : - - Lxx : List of variables to plot or coordinates along the X axis - - Lyy : List of variables to plot or coordinates along the Y axis - - Lxlab : List of x-axis label - - Lylab : List of y-axis label - - Lxlim : List of x (min, max) value plotted - - Lylim : List of y (min, max) value plotted - - Ltime : List of time (validity) - - Ltitle : List of sub-title - - Llinewidth : List of lines thickness - - Llinetype : List of line types - - Lcolorlines: List of color lines - - Llinelabel : List of legend label lines - - LfacconvX/Y: List of factors for unit conversion of the variables/coordinates to plot on X and Y axis - - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel - - Lid_overlap: List of number index of plot to overlap current variables - - LaxisColor : List of colors for multiple x-axis overlap - """ - self.ax = ax - firstCall = (len(self.ax) == 0) - # Defaults value convert to x number of variables list - if not LfacconvX: LfacconvX = [1.0]*len(Lxx) - if not LfacconvY: LfacconvY = [1.0]*len(Lxx) - if not Llinewidth: Llinewidth = [1.0]*len(Lxx) - if not Llinecolor: Llinecolor = ['blue']*len(Lxx) - if not Llinetype: Llinetype = ['-']*len(Lxx) - if not Llinelabel: Llinelabel = ['']*len(Lxx) - if not LaxisColor: LaxisColor = ['black']*len(Lxx) - if not Lylab: Lylab = ['']*len(Lxx) - if not Lxlab: Lxlab = ['']*len(Lxx) - - if firstCall: # 1st call - iax = 0 - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,1, label='graph axe x down')) - self.nb_graph+=1 - elif id_overlap: # overlapping plot with a different x-axis - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph, label='graph axe x top',frame_on=False)) - iax = len(self.ax)-1 - else: # existing ax with no overlapping (graph appended to existing panel) - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1)) - self.nb_graph+=1 - iax = len(self.ax)-1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter - - # On all variables to plot - for i,var in enumerate(Lxx): - # Print time validity - if Ltime: self.showTimeText(self.ax, iax, str(Ltime[i])) - - # Plot - cf = self.ax[iax].plot(Lxx[i]*LfacconvX[i], Lyy[i]*LfacconvY[i], color=Llinecolor[i], ls=Llinetype[i], - label=Llinelabel[i], linewidth=Llinewidth[i]) - # Legend - #TODO : Handling legend with overlap two axis lines in the same box. For now, placement is by hand - if not id_overlap: - self.ax[iax].legend(loc='upper right', bbox_to_anchor=(1, 0.95),fontsize=self.legendSize) - else: - self.ax[iax].legend(loc='upper right', bbox_to_anchor=(1, 0.90),fontsize=self.legendSize) + Llinecolor=[], Llinelabel=[], LfacconvX=[], LfacconvY=[], ax=[], id_overlap=None, Lxlim=[], + Lylim=[], Ltime=[], LaxisColor=[], LlocLegend=[]): + """ + XY (multiple)-lines plot + Parameters : + - Lxx : List of variables to plot or coordinates along the X axis + - Lyy : List of variables to plot or coordinates along the Y axis + - Lxlab : List of x-axis label + - Lylab : List of y-axis label + - Lxlim : List of x (min, max) value plotted + - Lylim : List of y (min, max) value plotted + - Ltime : List of time (validity) + - Ltitle : List of sub-title + - Llinewidth : List of lines thickness + - Llinetype : List of line types + - Lcolorlines: List of color lines + - Llinelabel : List of legend label lines + - LfacconvX/Y: List of factors for unit conversion of the variables/coordinates to plot on X and Y axis + - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel + - Lid_overlap: List of number index of plot to overlap current variables + - LaxisColor : List of colors for multiple x-axis overlap + - LlocLegend : List of localisation of the legend : 'best', 'upper left', 'upper right', 'lower left', 'lower right', + 'upper center', 'lower center', 'center left', 'center right', 'center' + """ + self.ax = ax + firstCall = (len(self.ax) == 0) + # Defaults value convert to x number of variables list + if not LfacconvX: + LfacconvX = [1.0] * len(Lxx) + if not LfacconvY: + LfacconvY = [1.0] * len(Lxx) + if not Llinewidth: + Llinewidth = [1.0] * len(Lxx) + if not Llinecolor: + Llinecolor = ['blue'] * len(Lxx) + if not Llinetype: + Llinetype = ['-'] * len(Lxx) + if not Llinelabel: + Llinelabel = [''] * len(Lxx) + if not LaxisColor: + LaxisColor = ['black'] * len(Lxx) + if not Lylab: + Lylab = [''] * len(Lxx) + if not Lxlab: + Lxlab = [''] * len(Lxx) + if not LlocLegend: + LlocLegend = ['upper right'] * len(Lxx) - # Title - if Ltitle: self.set_Title(self.ax, iax, Ltitle[i], id_overlap,Lxlab[i], Lylab[i]) - - # Ticks label - self.ax[iax].tick_params(axis='both', labelsize=self.xyTicksLabelSize, width=self.xyTicksWidth, length=self.xyTicksLength, pad=self.tickspad) - - # X/Y Axis label - if id_overlap: - self.ax[iax].xaxis.tick_top() - self.ax[iax].xaxis.set_label_position('top') - self.ax[iax].set_xlabel(Lxlab[i]) - self.ax[iax].set_ylabel(Lylab[i]) - else: - self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) - - self.ax[iax].tick_params(axis='x', colors=LaxisColor[i]) - - # X/Y Axis limits value - if Lxlim: - try: - self.set_xlim(self.ax, iax, Lxlim[i]) - except: - pass - if Lylim: - try: - self.set_ylim(self.ax, iax, Lylim[i]) - except: - pass - return self.fig - - def psectionH(self, lon=[],lat=[], Lvar=[], Lcarte=[], Llevel=[], Lxlab=[], Lylab=[], Ltitle=[], Lminval=[], Lmaxval=[], - Lstep=[], Lstepticks=[], Lcolormap=[], LcolorLine=[], Lcbarlabel=[], Lproj=[], Lfacconv=[], coastLines=True, ax=[], - Lid_overlap=[], colorbar=True, Ltime=[], LaddWhite_cm=[], LwhiteTop=[], Lpltype=[], Lcbformatlabel=[]): - """ - Horizontal cross section plot - Parameters : - - lon : longitude 2D array - - lat : latitude 2D array - - Lvar : List of variables to plot - - Lcarte : Zooming [lonmin, lonmax, latmin, latmax] - - Llevel : List of k-level value for the section plot (ignored if variable is already 2D) - - Lxlab : List of x-axis label - - Lylab : List of y-axis label - - Ltitle : List of sub-title - - Ltime : List of time (validity) - - Lminval: List of minimum value for each colorbar - - Lmaxval: List of maximum value for each colorbar - - Lstep : List of color-steps for each colorbar - - Lstepticks : List of value of labels for each colorbar - - Lcolormap : List of colormap - - LcolorLine : List of colors for colors arg of contour (color line only) - - Lcbarlabel : List of colorbar label legend (units) - - Lproj : List of ccrs cartopy projection ([] for cartesian coordinates) - - Lfacconv : List of factors for unit conversion of each variables - - coastLines : Boolean to plot coast lines and grid lines - - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel - - Lid_overlap : List of number index of plot to overlap current variables - - colorbar : show colorbar or not - - Lpltype : List of types of plot 'cf' or 'c'. cf=contourf, c=contour (lines only) - - LaddWhite_cm : List of boolean to add white color to a colormap at the first (low value) tick colorbar - - LwhiteTop : List of boolean to add the white color at the first top (high value). If false, the white is added at the bottom if Laddwhite_cm=T - - Lcbformatlabel: List of boolean to reduce the format to exponential 1.1E+02 format colorbar label - """ - self.ax = ax - firstCall = (len(self.ax) == 0) - - # Initialize default value w.r.t to the number of plots -# D={'lon':lon, 'lat':lat, 'Lcarte':Lcarte, 'Llevel':Llevel, 'Lxlab':Lxlab, 'Lylab':Lylab, 'Ltitle':Ltitle,'Lminval':Lminval, 'Lmaxval':Lmaxval, + if firstCall: # 1st call + iax = 0 + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, 1, label='graph axe x down')) + self.nb_graph += 1 + elif id_overlap: # overlapping plot with a different x-axis + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph, label='graph axe x top', frame_on=False)) + iax = len(self.ax) - 1 + else: # existing ax with no overlapping (graph appended to existing panel) + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1)) + self.nb_graph += 1 + iax = len(self.ax) - 1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter + + # On all variables to plot + for i, var in enumerate(Lxx): + # Print time validity + if Ltime: + self.showTimeText(self.ax, iax, str(Ltime[i])) + + # Plot + cf = self.ax[iax].plot(Lxx[i] * LfacconvX[i], Lyy[i] * LfacconvY[i], color=Llinecolor[i], ls=Llinetype[i], + label=Llinelabel[i], linewidth=Llinewidth[i]) + # Legend + # TODO : Handling legend with overlap two axis lines in the same box. For now, placement is by hand + if not id_overlap: + self.ax[iax].legend(loc=LlocLegend[i], bbox_to_anchor=(1, 0.95), fontsize=self.legendSize) + else: + self.ax[iax].legend(loc=LlocLegend[i], bbox_to_anchor=(1, 0.90), fontsize=self.legendSize) + + # Title + if Ltitle: + self.set_Title(self.ax, iax, Ltitle[i], id_overlap, Lxlab[i], Lylab[i]) + + # Ticks label + self.ax[iax].tick_params( + axis='both', + labelsize=self.xyTicksLabelSize, + width=self.xyTicksWidth, + length=self.xyTicksLength, + pad=self.tickspad) + + # X/Y Axis label + if id_overlap: + self.ax[iax].xaxis.tick_top() + self.ax[iax].xaxis.set_label_position('top') + self.ax[iax].set_xlabel(Lxlab[i]) + self.ax[iax].set_ylabel(Lylab[i]) + else: + self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) + + self.ax[iax].tick_params(axis='x', colors=LaxisColor[i]) + + # X/Y Axis limits value + if Lxlim: + try: + self.set_xlim(self.ax, iax, Lxlim[i]) + except BaseException: + pass + if Lylim: + try: + self.set_ylim(self.ax, iax, Lylim[i]) + except BaseException: + pass + return self.fig + + def psectionH(self, lon=[], lat=[], Lvar=[], Lcarte=[], Llevel=[], Lxlab=[], Lylab=[], Ltitle=[], Lminval=[], Lmaxval=[], + Lstep=[], Lstepticks=[], Lcolormap=[], LcolorLine=[], Lcbarlabel=[], Lproj=[], Lfacconv=[], coastLines=True, ax=[], + Lid_overlap=[], colorbar=True, Ltime=[], LaddWhite_cm=[], LwhiteTop=[], Lpltype=[], Lcbformatlabel=[], Llinewidth=[]): + """ + Horizontal cross section plot + Parameters : + - lon : longitude 2D array + - lat : latitude 2D array + - Lvar : List of variables to plot + - Lcarte : Zooming [lonmin, lonmax, latmin, latmax] + - Llevel : List of k-level value for the section plot (ignored if variable is already 2D) + - Lxlab : List of x-axis label + - Lylab : List of y-axis label + - Ltitle : List of sub-title + - Ltime : List of time (validity) + - Lminval: List of minimum value for each colorbar + - Lmaxval: List of maximum value for each colorbar + - Lstep : List of color-steps for each colorbar + - Lstepticks : List of value of labels for each colorbar + - Llinewidth : List of lines thickness of contour + - Lcolormap : List of colormap + - LcolorLine : List of colors for colors arg of contour (color line only) + - Lcbarlabel : List of colorbar label legend (units) + - Lproj : List of ccrs cartopy projection ([] for cartesian coordinates) + - Lfacconv : List of factors for unit conversion of each variables + - coastLines : Boolean to plot coast lines and grid lines + - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel + - Lid_overlap : List of number index of plot to overlap current variables + - colorbar : show colorbar or not + - Lpltype : List of types of plot 'cf' or 'c'. cf=contourf, c=contour (lines only) + - LaddWhite_cm : List of boolean to add white color to a colormap at the first (low value) tick colorbar + - LwhiteTop : List of boolean to add the white color at the first top (high value). If false, the white is added at the bottom if Laddwhite_cm=T + - Lcbformatlabel: List of boolean to reduce the format to exponential 1.1E+02 format colorbar label + """ + self.ax = ax + firstCall = (len(self.ax) == 0) + + # Initialize default value w.r.t to the number of plots +# D={'lon':lon, 'lat':lat, 'Lcarte':Lcarte, 'Llevel':Llevel, 'Lxlab':Lxlab, 'Lylab':Lylab, 'Ltitle':Ltitle,'Lminval':Lminval, 'Lmaxval':Lmaxval, # 'Lstep':Lstep, 'Lstepticks':Lstepticks, 'Lcolormap':Lcolormap, 'Lcbarlabel':Lcbarlabel, 'Lproj':Lproj, 'Lfacconv':Lfacconv, 'Ltime':Ltime, # 'LaddWhite_cm':LaddWhite_cm, 'Lpltype':Lpltype} # D = initialize_default_val(Lvar, D) - - # If all plots are not using conversion factor, convert it to List - if not Lfacconv: Lfacconv = [1.0]*len(Lvar) - if not Lylab: Lylab = ['']*len(Lvar) - if not Lxlab: Lxlab = ['']*len(Lvar) - if not Lcolormap and not LcolorLine: Lcolormap=['gist_rainbow_r']*len(Lvar) # If no color given, a cmap is given - if not Lcolormap: LcolorLine=['black']*len(Lvar) - if not Lpltype: Lpltype=['cf']*len(Lvar) - if not LaddWhite_cm: LaddWhite_cm=[False]*len(Lvar) - if not LwhiteTop: LwhiteTop=[False]*len(Lvar) - if not Lcbformatlabel: Lcbformatlabel=[False]*len(Lvar) - # Add an extra percentage of the top max value for forcing the colorbar show the true user maximum value (correct a bug) - if Lstep: Lmaxval = list(map(lambda x, y: x + 1E-6*y, Lmaxval, Lstep) ) #The extra value is 1E-6 times the step ticks of the colorbar - - # This following must be after the extra percentage top value addition - if not Lstep: Lstep=[None]*len(Lvar) - if not Lstepticks: Lstepticks=[None]*len(Lvar) - - # On all variables to plot - for i,var in enumerate(Lvar): - if firstCall: #1st call - iax = i - if Lproj: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1, projection = Lproj[i])) - else: # Cartesian coordinates - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1)) - self.nb_graph+=1 - elif Lid_overlap != []: #overlapping plot - iax = Lid_overlap[i] - else: #existing ax with no overlapping (graphd appended to existing panel) - if Lproj: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1, projection = Lproj[i])) - else: # Cartesian coordinates - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1)) - self.nb_graph+=1 - iax = len(self.ax)-1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter - - #Colors normalization - norm = mpl.colors.Normalize(vmax=Lmaxval[i], vmin=Lminval[i]) - - # Zooming - if len(Lcarte) == 4: #zoom - self.ax[iax].set_xlim(Lcarte[0], Lcarte[1]) - self.ax[iax].set_ylim(Lcarte[2], Lcarte[3]) - # Variable to plot w.r.t dimensions - if var.ndim==2: - vartoPlot = var[:,:] - else: - vartoPlot = var[Llevel[i],:,:] - - # Print min/max (stout and on plot) - self.set_minmaxText(self.ax, iax, vartoPlot, Ltitle[i], Lid_overlap, Lfacconv[i]) - - # Print time validity - if Ltime: self.showTimeText(self.ax, iax, str(Ltime[i])) - - # Number of contours level - if not Lstep[i]: # Default value of number of steps is 20 - Lstep[i] = (Lmaxval[i] - Lminval[i])/20 - Lstepticks[i] = Lstep[i] - - levels_contour = np.arange(Lminval[i],Lmaxval[i],step=Lstep[i]) - - # Add White to colormap - if LaddWhite_cm[i] and Lcolormap: Lcolormap[i]=self.addWhitecm(Lcolormap[i], len(levels_contour), LwhiteTop[i]) - - # Plot - if Lproj: - if Lpltype[i]=='c': # Contour - if LcolorLine: - cf = self.ax[iax].contour(lon[i], lat[i], vartoPlot*Lfacconv[i], levels=levels_contour,transform=Lproj[i], - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], colors=LcolorLine[i]) - else: - cf = self.ax[iax].contour(lon[i], lat[i], vartoPlot*Lfacconv[i], levels=levels_contour,transform=Lproj[i], - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) - else: - cf = self.ax[iax].contourf(lon[i], lat[i], vartoPlot*Lfacconv[i], levels=levels_contour,transform=Lproj[i], - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) - else: # Cartesian coordinates - if Lpltype[i]=='c': # Contour - cf = self.ax[iax].contour(lon[i], lat[i], vartoPlot*Lfacconv[i], levels=levels_contour, - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], colors=LcolorLine[i]) - else: - cf = self.ax[iax].contourf(lon[i], lat[i], vartoPlot*Lfacconv[i], levels=levels_contour, - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) - # Title - self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap,Lxlab[i], Lylab[i]) - - # Coastlines / Grid lines and labels - if Lproj: self.draw_Backmap(coastLines, self.ax[iax], Lproj[i]) - - # X/Y Axis - self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) - # Ticks label - self.ax[iax].tick_params(axis='both', labelsize=self.xyTicksLabelSize, width=self.xyTicksWidth, length=self.xyTicksLength, pad=self.tickspad) - - # Color label on contour-line - if Lpltype[i]=='c': # Contour - if 'GeoAxes' in str(type(self.ax[self.i])): # cartopy does not like the levels arguments in clabel, known issue - self.ax[iax].clabel(cf, fontsize=self.cbTicksLabelSize) - else: - self.ax[iax].clabel(cf, levels=np.arange(Lminval[i],Lmaxval[i],step=Lstep[i]), fontsize=self.cbTicksLabelSize) - - # Colorbar - if colorbar: - cb=plt.colorbar(cf, ax=self.ax[iax], fraction=0.031, pad=self.colorbarpad, ticks=np.arange(Lminval[i],Lmaxval[i],Lstepticks[i]), aspect=self.colorbaraspect) - cb.ax.set_title(Lcbarlabel[i], pad = self.labelcolorbarpad, loc='left', fontsize=self.cbTitleSize) #This creates a new AxesSubplot only for the colorbar y=0 ==> location at the bottom - if Lcbformatlabel[i]: cb.ax.set_yticklabels(["{:.1E}".format(i) for i in cb.get_ticks()]) - - return self.fig - - def pvector(self, Lxx=[], Lyy=[], Lvar1=[], Lvar2=[], Lcarte=[], Llevel=[], Lxlab=[], Lylab=[], - Ltitle=[], Lwidth=[], Larrowstep=[], Lcolor=[], Llegendval=[], Llegendlabel=[], + + # If all plots are not using conversion factor, convert it to List + if not Lfacconv: + Lfacconv = [1.0] * len(Lvar) + if not Lylab: + Lylab = [''] * len(Lvar) + if not Lxlab: + Lxlab = [''] * len(Lvar) + if not Lcolormap and not LcolorLine: + Lcolormap = ['gist_rainbow_r'] * len(Lvar) # If no color given, a cmap is given + if not Lcolormap: + LcolorLine = ['black'] * len(Lvar) + if not Lpltype: + Lpltype = ['cf'] * len(Lvar) + if not LaddWhite_cm: + LaddWhite_cm = [False] * len(Lvar) + if not LwhiteTop: + LwhiteTop = [False] * len(Lvar) + if not Lcbformatlabel: + Lcbformatlabel = [False] * len(Lvar) + if not Llinewidth: + Llinewidth = [1.0] * len(Lvar) + + # Add an extra percentage of the top max value for forcing the colorbar show the true user maximum value (correct a bug) + if Lstep: + Lmaxval = list(map(lambda x, y: x + 1E-6 * y, Lmaxval, Lstep)) # The extra value is 1E-6 times the step ticks of the colorbar + + # This following must be after the extra percentage top value addition + if not Lstep: + Lstep = [None] * len(Lvar) + if not Lstepticks: + Lstepticks = [None] * len(Lvar) + + # On all variables to plot + for i, var in enumerate(Lvar): + if firstCall: # 1st call + iax = i + if Lproj: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1, projection=Lproj[i])) + else: # Cartesian coordinates + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1)) + self.nb_graph += 1 + elif Lid_overlap != []: # overlapping plot + iax = Lid_overlap[i] + else: # existing ax with no overlapping (graphd appended to existing panel) + if Lproj: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1, projection=Lproj[i])) + else: # Cartesian coordinates + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1)) + self.nb_graph += 1 + iax = len(self.ax) - 1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter + + # Colors normalization + norm = mpl.colors.Normalize(vmax=Lmaxval[i], vmin=Lminval[i]) + + # Zooming + if len(Lcarte) == 4: # zoom + self.ax[iax].set_xlim(Lcarte[0], Lcarte[1]) + self.ax[iax].set_ylim(Lcarte[2], Lcarte[3]) + # Variable to plot w.r.t dimensions + if var.ndim == 2: + vartoPlot = var[:, :] + else: + vartoPlot = var[Llevel[i], :, :] + + # Print min/max (stout and on plot) + self.set_minmaxText(self.ax, iax, vartoPlot, Ltitle[i], Lid_overlap, Lfacconv[i]) + + # Print time validity + if Ltime: + self.showTimeText(self.ax, iax, str(Ltime[i])) + + # Number of contours level + if not Lstep[i]: # Default value of number of steps is 20 + Lstep[i] = (Lmaxval[i] - Lminval[i]) / 20 + Lstepticks[i] = Lstep[i] + + levels_contour = np.arange(Lminval[i], Lmaxval[i], step=Lstep[i]) + + # Add White to colormap + if LaddWhite_cm[i] and Lcolormap: + Lcolormap[i] = self.addWhitecm(Lcolormap[i], len(levels_contour), LwhiteTop[i]) + + # Plot + if Lproj: + if Lpltype[i] == 'c': # Contour + if LcolorLine: + cf = self.ax[iax].contour(lon[i], lat[i], vartoPlot * Lfacconv[i], levels=levels_contour, transform=Lproj[i], + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], colors=LcolorLine[i],linewidths=Llinewidth[i]) + else: + cf = self.ax[iax].contour(lon[i], lat[i], vartoPlot * Lfacconv[i], levels=levels_contour, transform=Lproj[i], + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i],linewidths=Llinewidth[i]) + else: + cf = self.ax[iax].contourf(lon[i], lat[i], vartoPlot * Lfacconv[i], levels=levels_contour, transform=Lproj[i], + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) + else: # Cartesian coordinates + if Lpltype[i] == 'c': # Contour + cf = self.ax[iax].contour(lon[i], lat[i], vartoPlot * Lfacconv[i], levels=levels_contour, + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], colors=LcolorLine[i],linewidths=Llinewidth[i]) + else: + cf = self.ax[iax].contourf(lon[i], lat[i], vartoPlot * Lfacconv[i], levels=levels_contour, + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) + # Title + self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap, Lxlab[i], Lylab[i]) + + # Coastlines / Grid lines and labels + if Lproj: + self.draw_Backmap(coastLines, self.ax[iax], Lproj[i]) + + # X/Y Axis + self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) + # Ticks label + self.ax[iax].tick_params( + axis='both', + labelsize=self.xyTicksLabelSize, + width=self.xyTicksWidth, + length=self.xyTicksLength, + pad=self.tickspad) + + # Color label on contour-line + if Lpltype[i] == 'c': # Contour + if 'GeoAxes' in str(type(self.ax[self.i])): # cartopy does not like the levels arguments in clabel, known issue + self.ax[iax].clabel(cf, fontsize=self.cbTicksLabelSize) + else: + self.ax[iax].clabel(cf, levels=np.arange(Lminval[i], Lmaxval[i], step=Lstep[i]), fontsize=self.cbTicksLabelSize) + + # Colorbar + if colorbar: + cb = plt.colorbar( + cf, + ax=self.ax[iax], + fraction=0.031, + pad=self.colorbarpad, + ticks=np.arange( + Lminval[i], + Lmaxval[i], + Lstepticks[i]), + aspect=self.colorbaraspect) + cb.ax.tick_params(labelsize=self.cbTicksLabelSize) + # This creates a new AxesSubplot only for the colorbar y=0 ==> location at the bottom + cb.ax.set_title(Lcbarlabel[i], pad=self.labelcolorbarpad, loc='left', fontsize=self.cbTitleSize) + if Lcbformatlabel[i]: + cb.ax.set_yticklabels(["{:.1E}".format(i) for i in cb.get_ticks()]) + + return self.fig + + def pvector(self, Lxx=[], Lyy=[], Lvar1=[], Lvar2=[], Lcarte=[], Llevel=[], Lxlab=[], Lylab=[], + Ltitle=[], Lwidth=[], Larrowstep=[], Lcolor=[], Llegendval=[], Llegendlabel=[], Lproj=[], Lfacconv=[], ax=[], coastLines=True, Lid_overlap=[], Ltime=[], Lscale=[], Lylim=[], Lxlim=[]): - """ - Vectors - Parameters : - - Lxx : List of x or y coordinate variable (lat or ni or nm) - - Lyy : List of y coordinates variable (lon or level) - - Lvar1 : List of wind-component along x/y or oblic axis (3D for hor. section, 2D for vertical section) - - Lvar2 : List of wind-component along y-axis : v-component for horizontal section / w-component for vertical section - - Lcarte : Zooming [lonmin, lonmax, latmin, latmax] - - Llevel : List of k-level value for the horizontal section plot (ignored if variable is already 2D) - - Lxlab : List of x-axis label - - Lylab : List of y-axis label - - Lxlim : List of x (min, max) value plotted - - Lylim : List of y (min, max) value plotted - - Ltitle : List of sub-titles - - Ltime : List of time (validity) - - Lwidth : List of thickness of the arrows - - Lscale : List of scale for the length of the arrows (high value <=> small length) - - Larrowstep : List of sub-sample (frequency) if too much arrows - - Lcolor : List of colors for the arrows (default: black) - - Llegendval : List of value for the legend of the default arrow - - Llegendlabel : List of labels for the legend of the default arrow - - Lproj : List of ccrs cartopy projection - - Lfacconv : List of factors for unit conversion of each variables - - coastLines : Boolean to plot coast lines and grid lines - - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel - - Lid_overlap : List of number index of plot to overlap current variables - """ - self.ax = ax - firstCall = (len(self.ax) == 0) - # If all plots are not using conversion factor, convert it to List - if not Lfacconv: Lfacconv = [1.0]*len(Lvar1) - if not Lcolor: Lcolor = ['black']*len(Lvar1) - if not Lscale : Lscale = [None]*len(Lvar1) - if not Lylab: Lylab = ['']*len(Lvar1) - if not Lxlab: Lxlab = ['']*len(Lvar1) - # On all variables to plot - for i,var1 in enumerate(Lvar1): - if firstCall: #1st call - iax = i - if Lproj: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1, projection = Lproj[i])) - else: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1)) - self.nb_graph+=1 - elif Lid_overlap != []: #overlapping plot - iax = Lid_overlap[i] - else: #existing ax with no overlapping (graphd appended to existing panel) - if Lproj: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1, projection = Lproj[i])) - else: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1)) - self.nb_graph+=1 - iax = len(self.ax)-1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter - - # Zooming - if len(Lcarte) == 4: #zoom - self.ax[iax].set_xlim(Lcarte[0], Lcarte[1]) - self.ax[iax].set_ylim(Lcarte[2], Lcarte[3]) - - # Variable to plot w.r.t dimensions - if var1.ndim==2: - vartoPlot1 = var1[:,:] - vartoPlot2 = Lvar2[i][:,:] - else: # Variable is 3D : only for horizontal section - vartoPlot1 = var1[Llevel[i],:,:] - vartoPlot2 = Lvar2[i][Llevel[i],:,:] - - # Print min/max val to help choose colorbar steps - self.set_minmaxText(self.ax, iax, np.sqrt(vartoPlot1**2 + vartoPlot2**2), Ltitle[i], Lid_overlap, Lfacconv[i]) - - # Print time validity - if Ltime: self.showTimeText(self.ax, iax, str(Ltime[i])) - - # Plot - axeX = Lxx[i] - axeY = Lyy[i] - if Lxx[i].ndim == 2: - cf = self.ax[iax].quiver(axeX[::Larrowstep[i],::Larrowstep[i]], axeY[::Larrowstep[i],::Larrowstep[i]], - vartoPlot1[::Larrowstep[i],::Larrowstep[i]], vartoPlot2[::Larrowstep[i],::Larrowstep[i]], - width=Lwidth[i] ,angles='uv', color=Lcolor[i],scale=Lscale[i]) - else: - cf = self.ax[iax].quiver(axeX[::Larrowstep[i]], axeY[::Larrowstep[i]], - vartoPlot1[::Larrowstep[i],::Larrowstep[i]], vartoPlot2[::Larrowstep[i],::Larrowstep[i]], - width=Lwidth[i] ,angles='uv', color=Lcolor[i], scale=Lscale[i]) - # Title - self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap,Lxlab[i], Lylab[i]) - - # X/Y Axis Label - self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) - - # X/Y Axis limits value - if Lxlim: - try: - self.set_xlim(self.ax, iax, Lxlim[i]) - except: - pass - if Lylim: - try: - self.set_ylim(self.ax, iax, Lylim[i]) - except: - pass - - # Coastlines: - if Lproj: self.draw_Backmap(coastLines, self.ax[iax], Lproj[i]) - - # Arrow legend key - qk = self.ax[iax].quiverkey(cf, 1.0, -0.05, Llegendval[i], str(Llegendval[i]) + Llegendlabel[i], labelpos='E', color='black') - - return self.fig + """ + Vectors + Parameters : + - Lxx : List of x or y coordinate variable (lat or ni or nm) + - Lyy : List of y coordinates variable (lon or level) + - Lvar1 : List of wind-component along x/y or oblic axis (3D for hor. section, 2D for vertical section) + - Lvar2 : List of wind-component along y-axis : v-component for horizontal section / w-component for vertical section + - Lcarte : Zooming [lonmin, lonmax, latmin, latmax] + - Llevel : List of k-level value for the horizontal section plot (ignored if variable is already 2D) + - Lxlab : List of x-axis label + - Lylab : List of y-axis label + - Lxlim : List of x (min, max) value plotted + - Lylim : List of y (min, max) value plotted + - Ltitle : List of sub-titles + - Ltime : List of time (validity) + - Lwidth : List of thickness of the arrows + - Lscale : List of scale for the length of the arrows (high value <=> small length) + - Larrowstep : List of sub-sample (frequency) if too much arrows + - Lcolor : List of colors for the arrows (default: black) + - Llegendval : List of value for the legend of the default arrow + - Llegendlabel : List of labels for the legend of the default arrow + - Lproj : List of ccrs cartopy projection + - Lfacconv : List of factors for unit conversion of each variables + - coastLines : Boolean to plot coast lines and grid lines + - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel + - Lid_overlap : List of number index of plot to overlap current variables + """ + self.ax = ax + firstCall = (len(self.ax) == 0) + # If all plots are not using conversion factor, convert it to List + if not Lfacconv: + Lfacconv = [1.0] * len(Lvar1) + if not Lcolor: + Lcolor = ['black'] * len(Lvar1) + if not Lscale: + Lscale = [None] * len(Lvar1) + if not Lylab: + Lylab = [''] * len(Lvar1) + if not Lxlab: + Lxlab = [''] * len(Lvar1) + # On all variables to plot + for i, var1 in enumerate(Lvar1): + if firstCall: # 1st call + iax = i + if Lproj: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1, projection=Lproj[i])) + else: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1)) + self.nb_graph += 1 + elif Lid_overlap != []: # overlapping plot + iax = Lid_overlap[i] + else: # existing ax with no overlapping (graphd appended to existing panel) + if Lproj: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1, projection=Lproj[i])) + else: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1)) + self.nb_graph += 1 + iax = len(self.ax) - 1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter + + # Zooming + if len(Lcarte) == 4: # zoom + self.ax[iax].set_xlim(Lcarte[0], Lcarte[1]) + self.ax[iax].set_ylim(Lcarte[2], Lcarte[3]) + + # Variable to plot w.r.t dimensions + if var1.ndim == 2: + vartoPlot1 = var1[:, :] + vartoPlot2 = Lvar2[i][:, :] + else: # Variable is 3D : only for horizontal section + vartoPlot1 = var1[Llevel[i], :, :] + vartoPlot2 = Lvar2[i][Llevel[i], :, :] + + # Print min/max val to help choose colorbar steps + self.set_minmaxText(self.ax, iax, np.sqrt(vartoPlot1**2 + vartoPlot2**2), Ltitle[i], Lid_overlap, Lfacconv[i]) + + # Print time validity + if Ltime: + self.showTimeText(self.ax, iax, str(Ltime[i])) + + # Plot + axeX = Lxx[i] + axeY = Lyy[i] + if Lxx[i].ndim == 2: + cf = self.ax[iax].quiver(axeX[::Larrowstep[i], ::Larrowstep[i]], axeY[::Larrowstep[i], ::Larrowstep[i]], + vartoPlot1[::Larrowstep[i], ::Larrowstep[i]], vartoPlot2[::Larrowstep[i], ::Larrowstep[i]], + width=Lwidth[i], angles='uv', color=Lcolor[i], scale=Lscale[i]) + else: + cf = self.ax[iax].quiver(axeX[::Larrowstep[i]], axeY[::Larrowstep[i]], + vartoPlot1[::Larrowstep[i], ::Larrowstep[i]], vartoPlot2[::Larrowstep[i], ::Larrowstep[i]], + width=Lwidth[i], angles='uv', color=Lcolor[i], scale=Lscale[i]) + # Title + self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap, Lxlab[i], Lylab[i]) + + # X/Y Axis Label + self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) + + # X/Y Axis limits value + if Lxlim: + try: + self.set_xlim(self.ax, iax, Lxlim[i]) + except BaseException: + pass + if Lylim: + try: + self.set_ylim(self.ax, iax, Lylim[i]) + except BaseException: + pass + + # Coastlines: + if Lproj: + self.draw_Backmap(coastLines, self.ax[iax], Lproj[i]) + + # Arrow legend key + qk = self.ax[iax].quiverkey(cf, 1.0, -0.05, Llegendval[i], str(Llegendval[i]) + Llegendlabel[i], + labelpos='E', color='black', fontproperties={'size': self.legendSize}) + + return self.fig def pstreamline(self, Lxx=[], Lyy=[], Lvar1=[], Lvar2=[], Lcarte=[], Llevel=[], Lxlab=[], Lylab=[], Llinewidth=[], Ldensity=[], - Ltitle=[], Lcolor=[], Lproj=[], Lfacconv=[], ax=[], coastLines=True, Lid_overlap=[], Ltime=[], - Lylim=[], Lxlim=[]): - """ - Wind stream lines - Parameters : - - Lxx : List of x or y coordinate variable (lat or ni or nm) - - Lyy : List of y coordinates variable (lon or level) - - Lvar1 : List of wind-component along x/y or oblic axis (3D for hor. section, 2D for vertical section) - - Lvar2 : List of wind-component along y-axis : v-component for horizontal section / w-component for vertical section - - Lcarte : Zooming [lonmin, lonmax, latmin, latmax] - - Llevel : List of k-level value for the horizontal section plot (ignored if variable is already 2D) - - Lxlab : List of x-axis label - - Lylab : List of y-axis label - - Lxlim : List of x (min, max) value plotted - - Lylim : List of y (min, max) value plotted - - Ltitle : List of sub-titles - - Ltime : List of time (validity) - - Llinewidth : List of lines thickness - - Ldensity : List of density that control the closeness of streamlines - - Lcolor : List of colors for the streamline (default: black) - - Lproj : List of ccrs cartopy projection - - Lfacconv : List of factors for unit conversion of each variables - - coastLines : Boolean to plot coast lines and grid lines - - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel - - Lid_overlap : List of number index of plot to overlap current variables - """ - self.ax = ax - firstCall = (len(self.ax) == 0) - # If all plots are not using conversion factor, convert it to List - if not Lfacconv: Lfacconv = [1.0]*len(Lvar1) - if not Lcolor: Lcolor = ['black']*len(Lvar1) - if not Lylab: Lylab = ['']*len(Lvar1) - if not Lxlab: Lxlab = ['']*len(Lvar1) - if not Llinewidth: Llinewidth = [1.0]*len(Lvar1) - if not Ldensity: Ldensity = [1.0]*len(Lvar1) - - # On all variables to plot - for i,var1 in enumerate(Lvar1): - if firstCall: #1st call - iax = i - if Lproj: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1, projection = Lproj[i])) - else: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1)) - self.nb_graph+=1 - elif Lid_overlap != []: #overlapping plot - iax = Lid_overlap[i] - else: #existing ax with no overlapping (graphd appended to existing panel) - if Lproj: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1, projection = Lproj[i])) - else: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1)) - self.nb_graph+=1 - iax = len(self.ax)-1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter - - # Zooming - if len(Lcarte) == 4: #zoom - self.ax[iax].set_xlim(Lcarte[0], Lcarte[1]) - self.ax[iax].set_ylim(Lcarte[2], Lcarte[3]) - - # Variable to plot w.r.t dimensions - if var1.ndim==2: - vartoPlot1 = var1[:,:] - vartoPlot2 = Lvar2[i][:,:] - else: # Variable is 3D : only for horizontal section - vartoPlot1 = var1[Llevel[i],:,:] - vartoPlot2 = Lvar2[i][Llevel[i],:,:] - - # Print min/max val to help choose steps - self.set_minmaxText(self.ax, iax, np.sqrt(vartoPlot1**2 + vartoPlot2**2), Ltitle[i], Lid_overlap, Lfacconv[i]) - - # Print time validity - if Ltime: self.showTimeText(self.ax, iax, str(Ltime[i])) + Ltitle=[], Lcolor=[], Lproj=[], Lfacconv=[], ax=[], coastLines=True, Lid_overlap=[], Ltime=[], + Lylim=[], Lxlim=[]): + """ + Wind stream lines + Parameters : + - Lxx : List of x or y coordinate variable (lat or ni or nm) + - Lyy : List of y coordinates variable (lon or level) + - Lvar1 : List of wind-component along x/y or oblic axis (3D for hor. section, 2D for vertical section) + - Lvar2 : List of wind-component along y-axis : v-component for horizontal section / w-component for vertical section + - Lcarte : Zooming [lonmin, lonmax, latmin, latmax] + - Llevel : List of k-level value for the horizontal section plot (ignored if variable is already 2D) + - Lxlab : List of x-axis label + - Lylab : List of y-axis label + - Lxlim : List of x (min, max) value plotted + - Lylim : List of y (min, max) value plotted + - Ltitle : List of sub-titles + - Ltime : List of time (validity) + - Llinewidth : List of lines thickness + - Ldensity : List of density that control the closeness of streamlines + - Lcolor : List of colors for the streamline (default: black) + - Lproj : List of ccrs cartopy projection + - Lfacconv : List of factors for unit conversion of each variables + - coastLines : Boolean to plot coast lines and grid lines + - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel + - Lid_overlap : List of number index of plot to overlap current variables + """ + self.ax = ax + firstCall = (len(self.ax) == 0) + # If all plots are not using conversion factor, convert it to List + if not Lfacconv: + Lfacconv = [1.0] * len(Lvar1) + if not Lcolor: + Lcolor = ['black'] * len(Lvar1) + if not Lylab: + Lylab = [''] * len(Lvar1) + if not Lxlab: + Lxlab = [''] * len(Lvar1) + if not Llinewidth: + Llinewidth = [1.0] * len(Lvar1) + if not Ldensity: + Ldensity = [1.0] * len(Lvar1) - # Plot - cf = self.ax[iax].streamplot(Lxx[i], Lyy[i], vartoPlot1, vartoPlot2, density=Ldensity[i], linewidth=Llinewidth[i], color=Lcolor[i]) - - # Title - self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap,Lxlab[i], Lylab[i]) + # On all variables to plot + for i, var1 in enumerate(Lvar1): + if firstCall: # 1st call + iax = i + if Lproj: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1, projection=Lproj[i])) + else: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1)) + self.nb_graph += 1 + elif Lid_overlap != []: # overlapping plot + iax = Lid_overlap[i] + else: # existing ax with no overlapping (graphd appended to existing panel) + if Lproj: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1, projection=Lproj[i])) + else: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1)) + self.nb_graph += 1 + iax = len(self.ax) - 1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter + + # Zooming + if len(Lcarte) == 4: # zoom + self.ax[iax].set_xlim(Lcarte[0], Lcarte[1]) + self.ax[iax].set_ylim(Lcarte[2], Lcarte[3]) + + # Variable to plot w.r.t dimensions + if var1.ndim == 2: + vartoPlot1 = var1[:, :] + vartoPlot2 = Lvar2[i][:, :] + else: # Variable is 3D : only for horizontal section + vartoPlot1 = var1[Llevel[i], :, :] + vartoPlot2 = Lvar2[i][Llevel[i], :, :] + + # Print min/max val to help choose steps + self.set_minmaxText(self.ax, iax, np.sqrt(vartoPlot1**2 + vartoPlot2**2), Ltitle[i], Lid_overlap, Lfacconv[i]) + + # Print time validity + if Ltime: + self.showTimeText(self.ax, iax, str(Ltime[i])) + + # Plot + cf = self.ax[iax].streamplot(Lxx[i], Lyy[i], vartoPlot1, vartoPlot2, density=Ldensity[i], linewidth=Llinewidth[i], color=Lcolor[i]) + + # Title + self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap, Lxlab[i], Lylab[i]) + + # X/Y Axis Label + self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) + + # X/Y Axis limits value + if Lxlim: + try: + self.set_xlim(self.ax, iax, Lxlim[i]) + except BaseException: + pass + if Lylim: + try: + self.set_ylim(self.ax, iax, Lylim[i]) + except BaseException: + pass + + # Coastlines: + if Lproj: + self.draw_Backmap(coastLines, self.ax[iax], Lproj[i]) + + return self.fig - # X/Y Axis Label - self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) - - # X/Y Axis limits value - if Lxlim: - try: - self.set_xlim(self.ax, iax, Lxlim[i]) - except: - pass - if Lylim: - try: - self.set_ylim(self.ax, iax, Lylim[i]) - except: - pass - - # Coastlines: - if Lproj: self.draw_Backmap(coastLines, self.ax[iax], Lproj[i]) - - return self.fig - def pXY_bar(self, Lbins=[], Lvar=[], Lxlab=[], Lylab=[], Ltitle=[], Lcolor=[], Lwidth=[], - Llinecolor=[], Llinewidth=[], Lfacconv=[], ax=[], id_overlap=None, Lxlim=[], - Lylim=[], Ltime=[], LaxisColor=[]): - """ - XY Histogram - Parameters : - - Lbins : List of bins - - Lvar : List of the value for each bin - - Lxlab : List of x-axis label - - Lylab : List of y-axis label - - Ltitle : List of sub-title - - Lcolor : List of color (or sequence of colors for each value) of the bars - - Lwidth : List of width of the bars - - Llinecolor : List of line color of the bar edges - - Llinewidth : List of lines thickness of the bar edges - - Lfacconv : List of factors for unit conversion of each variables - - Lxlim : List of x (min, max) value plotted - - Lylim : List of y (min, max) value plotted - - Ltime : List of time (validity) - - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel - - Lid_overlap: List of number index of plot to overlap current variables - - LaxisColor : List of colors for multiple x-axis overlap - """ - self.ax = ax - firstCall = (len(self.ax) == 0) - # Defaults value convert to x number of variables list - if not Lfacconv: Lfacconv = [1.0]*len(Lvar) - if not LaxisColor: LaxisColor = ['black']*len(Lvar) - if not Lylab: Lylab = ['']*len(Lvar) - if not Lxlab: Lxlab = ['']*len(Lvar) - if not Lcolor: Lcolor = ['black']*len(Lvar) - if not Lwidth: Lwidth = [1]*len(Lvar) - if not Llinecolor: Llinecolor = ['black']*len(Lvar) - if not Llinewidth: Llinewidth = [0]*len(Lvar) - - # On all variables to plot - for i,var in enumerate(Lvar): - if firstCall: # 1st call - iax = i - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1, label='graph axe x down')) - self.nb_graph+=1 - elif id_overlap: # overlapping plot with a different x-axis - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph, label='graph axe x top',frame_on=False)) - iax = len(self.ax)-1 - else: # existing ax with no overlapping (graph appended to existing panel) - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1)) - self.nb_graph+=1 - print(len(self.ax)) - iax = len(self.ax)-1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter - - # Print time validity - if Ltime: self.showTimeText(self.ax, iax, str(Ltime[i])) - - # Bins by labels - labels = np.array([str(L) for L in Lbins[i]]) - - # Plot - cf = self.ax[iax].bar(labels, var*Lfacconv[i], width=Lwidth[i], color=Lcolor[i], linewidth=Llinewidth[i], edgecolor=Llinecolor[i]) - - # Legend - #TODO : Handling legend with overlap two axis lines in the same box. For now, placement is by hand - if not id_overlap: - self.ax[iax].legend(loc='upper right', bbox_to_anchor=(1, 0.95)) - else: - self.ax[iax].legend(loc='upper right', bbox_to_anchor=(1, 0.90)) - - # Title - if Ltitle: self.set_Title(self.ax, iax, Ltitle[i], id_overlap,Lxlab[i], Lylab[i]) - - # X/Y Axis label - if id_overlap: - self.ax[iax].xaxis.tick_top() - self.ax[iax].xaxis.set_label_position('top') - self.ax[iax].set_xlabel(Lxlab[i]) - self.ax[iax].set_ylabel(Lylab[i]) - else: - self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) - - self.ax[iax].tick_params(axis='x', colors=LaxisColor[i]) - - # X/Y Axis limits value - if Lxlim: - try: - self.set_xlim(self.ax, iax, Lxlim[i]) - except: - pass - if Lylim: - try: - self.set_ylim(self.ax, iax, Lylim[i]) - except: - pass - - return self.fig -#def initialize_default_val(Lvar, Dparam): + Llinecolor=[], Llinewidth=[], Lfacconv=[], ax=[], id_overlap=None, Lxlim=[], + Lylim=[], Ltime=[], LaxisColor=[], LlocLegend=[]): + """ + XY Histogram + Parameters : + - Lbins : List of bins + - Lvar : List of the value for each bin + - Lxlab : List of x-axis label + - Lylab : List of y-axis label + - Ltitle : List of sub-title + - Lcolor : List of color (or sequence of colors for each value) of the bars + - Lwidth : List of width of the bars + - Llinecolor : List of line color of the bar edges + - Llinewidth : List of lines thickness of the bar edges + - Lfacconv : List of factors for unit conversion of each variables + - Lxlim : List of x (min, max) value plotted + - Lylim : List of y (min, max) value plotted + - Ltime : List of time (validity) + - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel + - Lid_overlap: List of number index of plot to overlap current variables + - LaxisColor : List of colors for multiple x-axis overlap + - LlocLegend : List of localisation of the legend : 'best', 'upper left', 'upper right', 'lower left', 'lower right', + 'upper center', 'lower center', 'center left', 'center right', 'center' + """ + self.ax = ax + firstCall = (len(self.ax) == 0) + # Defaults value convert to x number of variables list + if not Lfacconv: + Lfacconv = [1.0] * len(Lvar) + if not LaxisColor: + LaxisColor = ['black'] * len(Lvar) + if not Lylab: + Lylab = [''] * len(Lvar) + if not Lxlab: + Lxlab = [''] * len(Lvar) + if not Lcolor: + Lcolor = ['black'] * len(Lvar) + if not Lwidth: + Lwidth = [1] * len(Lvar) + if not Llinecolor: + Llinecolor = ['black'] * len(Lvar) + if not Llinewidth: + Llinewidth = [0] * len(Lvar) + if not LlocLegend: + LlocLegend = ['upper right'] * len(Lvar) + + # On all variables to plot + for i, var in enumerate(Lvar): + if firstCall: # 1st call + iax = i + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1, label='graph axe x down')) + self.nb_graph += 1 + elif id_overlap: # overlapping plot with a different x-axis + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph, label='graph axe x top', frame_on=False)) + iax = len(self.ax) - 1 + else: # existing ax with no overlapping (graph appended to existing panel) + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1)) + self.nb_graph += 1 + print(len(self.ax)) + iax = len(self.ax) - 1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter + + # Print time validity + if Ltime: + self.showTimeText(self.ax, iax, str(Ltime[i])) + + # Bins by labels + labels = np.array([str(L) for L in Lbins[i]]) + + # Plot + cf = self.ax[iax].bar(labels, var * Lfacconv[i], width=Lwidth[i], color=Lcolor[i], linewidth=Llinewidth[i], edgecolor=Llinecolor[i]) + + # Legend + # TODO : Handling legend with overlap two axis lines in the same box. For now, placement is by hand + if not id_overlap: + self.ax[iax].legend(loc=LlocLegend[i], bbox_to_anchor=(1, 0.95)) + else: + self.ax[iax].legend(loc=LlocLegend[i], bbox_to_anchor=(1, 0.90)) + + # Title + if Ltitle: + self.set_Title(self.ax, iax, Ltitle[i], id_overlap, Lxlab[i], Lylab[i]) + + # X/Y Axis label + if id_overlap: + self.ax[iax].xaxis.tick_top() + self.ax[iax].xaxis.set_label_position('top') + self.ax[iax].set_xlabel(Lxlab[i]) + self.ax[iax].set_ylabel(Lylab[i]) + else: + self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) + + self.ax[iax].tick_params(axis='x', colors=LaxisColor[i]) + + # X/Y Axis limits value + if Lxlim: + try: + self.set_xlim(self.ax, iax, Lxlim[i]) + except BaseException: + pass + if Lylim: + try: + self.set_ylim(self.ax, iax, Lylim[i]) + except BaseException: + pass + + return self.fig +# def initialize_default_val(Lvar, Dparam): # #TODO : initialize default value for all parameter of all type of graphs # #Returns All the parameters given in Dparam where : # "- If no value is found (empty list []) : return the default value (if exist) * number of graph # #- If ONE value only is found : return the value copied x times the number of graph # #- If the list is complete, nothing is done # #CURRENT PROBLEM -# #The returned value do not change the true referenced variable given as argument +# #The returned value do not change the true referenced variable given as argument # # Number of graphs # l = len(Lvar) -# +# # # Add an extra percentage of the top max value for forcing the colorbar show the true user maximum value (correct a bug) # #if Dparam['Lstep'] and Dparam['Lmaxval']: Dparam['Lmaxval'] = list(map(lambda x, y: x + 1E-6*y, Dparam['Lmaxval'], Dparam['Lstep']) ) #The extra value is 1E-6 times the step ticks of the colorbar # print(Dparam.items()) -# +# # # Read on which parameters initialize the default values # for args_t in list(Dparam.items()): # Test on all arguments present, if they are empty list, default values apply for each graph # args = list(args_t) @@ -936,16 +1072,16 @@ class PanelPlot(): # elif args[0] == 'Lstepticks' and not args[1]: args[1] = [None]*l # Dparam[args[0]] = args[1] # print(args) -# +# # # Check if there is no value for a parameter -## for args_t in list(Dparam.items()): +# for args_t in list(Dparam.items()): ## args = list(args_t) -## if args[1] +# if args[1] # # Number of contours level -## for i in range(l): -## if 'Lstepticks' in arguments.args and 'Lmaxval' in arguments.args and 'Lminval' in arguments.args: -## Lstep[i] = (Lmaxval[i] - Lminval[i])/20 # Default value of number of steps is 20 -## Lstepticks[i] = Lstep[i] # Default value is stepticks are the same as steps values +# for i in range(l): +# if 'Lstepticks' in arguments.args and 'Lmaxval' in arguments.args and 'Lminval' in arguments.args: +# Lstep[i] = (Lmaxval[i] - Lminval[i])/20 # Default value of number of steps is 20 +# Lstepticks[i] = Lstep[i] # Default value is stepticks are the same as steps values +# # -# # return Dparam diff --git a/src/LIB/Python/misc_functions.py b/src/LIB/Python/misc_functions.py index efc0916e04a60a7ee3cd2a276d1bb1eba3de0182..53577de6ee71fbc12524c19af6daff7657bd5fe9 100644 --- a/src/LIB/Python/misc_functions.py +++ b/src/LIB/Python/misc_functions.py @@ -13,207 +13,214 @@ from scipy.interpolate import RectBivariateSpline import numpy as np import math + def convert_date(datesince, time_in_sec): - return str(time_in_sec) + datesince[:33] + return str(time_in_sec) + datesince[:33] + class mean_operator(): - def MYM(self,var): + def MYM(self, var): ny = var.shape[1] out = copy.deepcopy(var) - for j in range(ny-1): - out[:,j,:] = (var[:,j,:] + var[:,j+1,:])*0.5 + for j in range(ny - 1): + out[:, j, :] = (var[:, j, :] + var[:, j + 1, :]) * 0.5 return out - - def MXM(self,var): + + def MXM(self, var): nx = var.shape[2] out = copy.deepcopy(var) - for i in range(nx-1): - out[:,:,i] = (var[:,:,i] + var[:,:,i+1])*0.5 + for i in range(nx - 1): + out[:, :, i] = (var[:, :, i] + var[:, :, i + 1]) * 0.5 return out - - def MZM(self,var): + + def MZM(self, var): nz = var.shape[0] out = copy.deepcopy(var) - for k in range(nz-1): - out[k,:,:] = (var[k,:,:] + var[k+1,:,:])*0.5 + for k in range(nz - 1): + out[k, :, :] = (var[k, :, :] + var[k + 1, :, :]) * 0.5 return out + def windvec_verti_proj(u, v, level, angle): """Compute the projected horizontal wind vector on an axis with a given angle w.r.t. the x/ni axes (West-East) - + Parameters ---------- u : array 3D U-wind component - + v : array 3D V-wind component level : array 1D level dimension array - + angle : float angle (radian) of the new axe w.r.t the x/ni axes (West-East). angle = 0 for (z,x) sections, angle=pi/2 for (z,y) sections - + Returns ------- - + projected_wind : array 3D a 3D wind component projected on the axe to be used with Panel_Plot.pvector as Lvar1 - """ + """ projected_wind = copy.deepcopy(u) for k in range(len(level)): - projected_wind[k,:,:] = u[k,:,:]*math.cos(angle) + v[k,:,:]*math.sin(angle) + projected_wind[k, :, :] = u[k, :, :] * math.cos(angle) + v[k, :, :] * math.sin(angle) return projected_wind + def oblique_proj(var, ni, nj, lvl, i_beg, j_beg, i_end, j_end): """Compute an oblique projection of a variable w.r.t. its axes - + Parameters ---------- var : array 3D or 2D the variable to project (e.g. THT, ZS) - + ni : array 1D 1D x-axis of the 3D dimension - + nj : array 1D 1D y-axis of the 3D dimension - + level : array 1D 1D z-axe of the 3D dimension - + i_beg, j_beg : int coordinate of the begin point of the new axe - + i_end, j_end : int - coordinate of the end point of the new axe - + coordinate of the end point of the new axe + Returns ------- angle_proj : float the angle (radian) of the new axe w.r.t the x/ni axes (West-East) - + out_var : array 2D or 1D a 2D (z,m) or 1D (m) variable projected on the oblique axe - + axe_m : array 1D a 1D m new axe (distance from the beggining point) - + """ - dist_seg=np.sqrt((i_end-i_beg)**2.0 + (j_end-j_beg)**2.0) # Distance de la section oblique m - if var.ndim ==3: - out_var = np.zeros((len(lvl),int(dist_seg)+1)) # Initialisation du nouveau champs projeté dans la coupe (z,m) + dist_seg = np.sqrt((i_end - i_beg)**2.0 + (j_end - j_beg)**2.0) # Distance de la section oblique m + if var.ndim == 3: + out_var = np.zeros((len(lvl), int(dist_seg) + 1)) # Initialisation du nouveau champs projeté dans la coupe (z,m) else: # 2D - out_var = np.zeros(int(dist_seg)+1) # Initialisation du nouveau champs projeté dans la coupe (m) + out_var = np.zeros(int(dist_seg) + 1) # Initialisation du nouveau champs projeté dans la coupe (m) - axe_m = np.zeros(int(dist_seg)+1) # Axe des abscisses qui sera tracé selon la coupe - axe_m_coord = [] # Coordonnées x,y des points qui composent l'axe - axe_m_coord.append( (ni[i_beg],nj[j_beg]) ) # Le premier point est celui donné par l'utilisateur + axe_m = np.zeros(int(dist_seg) + 1) # Axe des abscisses qui sera tracé selon la coupe + axe_m_coord = [] # Coordonnées x,y des points qui composent l'axe + axe_m_coord.append((ni[i_beg], nj[j_beg])) # Le premier point est celui donné par l'utilisateur - for m in range(int(dist_seg)): # Discrétisation selon distance de la coupe / int(distance_de_la_coupe) - axe_m_coord.append( (axe_m_coord[0][0] + (ni[i_end]-ni[i_beg])/(int(dist_seg))*(m+1), - axe_m_coord[0][1] + (nj[j_end]-nj[j_beg])/(int(dist_seg))*(m+1) )) - axe_m[m+1] = np.sqrt((ni[i_beg]-axe_m_coord[m+1][0])**2 + (nj[j_beg]-axe_m_coord[m+1][1])**2) + for m in range(int(dist_seg)): # Discrétisation selon distance de la coupe / int(distance_de_la_coupe) + axe_m_coord.append((axe_m_coord[0][0] + (ni[i_end] - ni[i_beg]) / (int(dist_seg)) * (m + 1), + axe_m_coord[0][1] + (nj[j_end] - nj[j_beg]) / (int(dist_seg)) * (m + 1))) + axe_m[m + 1] = np.sqrt((ni[i_beg] - axe_m_coord[m + 1][0])**2 + (nj[j_beg] - axe_m_coord[m + 1][1])**2) - if var.ndim ==3: # 3D variables to project + if var.ndim == 3: # 3D variables to project for k in range(len(lvl)): - a=RectBivariateSpline(ni, nj,var[k,:,:],kx=1,ky=1) # Interpolation par niveau à l'ordre 1 pour éviter des valeurs négatives de champs strictement > 0 - for m in range(int(dist_seg)+1): - out_var[k,m] = a.ev(axe_m_coord[m][0],axe_m_coord[m][1]) # La fonction ev de RectBivariate retourne la valeur la plus proche du point considéré + a = RectBivariateSpline(nj, ni, var[k, :, :], kx=1, ky=1) + for m in range(int(dist_seg) + 1): + # La fonction ev de RectBivariate retourne la valeur la plus proche du point considéré + out_var[k, m] = a.ev(axe_m_coord[m][1], axe_m_coord[m][0]) else: # 2D variables to project - a=RectBivariateSpline(ni, nj,var[:,:],kx=1,ky=1) - for m in range(int(dist_seg)+1): - out_var[m] = a.ev(axe_m_coord[m][0],axe_m_coord[m][1]) - - angle_proj = math.acos((ni[i_end]-ni[i_beg])/axe_m[-1]) + a = RectBivariateSpline(nj, ni, var[:, :], kx=1, ky=1) + for m in range(int(dist_seg) + 1): + out_var[m] = a.ev(axe_m_coord[m][1], axe_m_coord[m][0]) + + angle_proj = math.acos((ni[i_end] - ni[i_beg]) / axe_m[-1]) return angle_proj, out_var, axe_m + def comp_altitude1DVar(oneVar2D, orography, ztop, level, n_xory): """Compute and returns an altitude and x or y grid mesh variable in 2D following the topography in 1D To be used with 2D simulations - + Parameters ---------- oneVar2D : array 2D a 2D array (e.g. UT, THT) - + orography : array 1D 1D orography (ZS) - + ztop : real scalar of the top height of the model (ZTOP) - + level : array 1D 1D level variable (level or level_w) - - n_xory : array 1D + + n_xory : array 1D 1D directionnal grid variable (ni_u, nj_u, ni_v or nj_v) - + Returns ------- altitude a 2D altitude variable with topography taken into account - + n_xory_2D a 2D directionnal variable duplicated from n_xory """ n_xory_2D = copy.deepcopy(oneVar2D) - altitude = copy.deepcopy(oneVar2D) - + altitude = copy.deepcopy(oneVar2D) + for k in range(len(level)): - n_xory_2D[k,:] = n_xory + n_xory_2D[k, :] = n_xory for j in range(len(n_xory)): for k in range(len(level)): - altitude[k,j] = orography[j] + level[k]*((ztop-orography[j])/ztop) + altitude[k, j] = orography[j] + level[k] * ((ztop - orography[j]) / ztop) return altitude, n_xory_2D + def comp_altitude2DVar(oneVar3D, orography, ztop, level, n_y, n_x): """Compute and returns an altitude and x or y grid mesh variable in 3D following the topography in 2D To be used with 3D simulations - + Parameters ---------- oneVar3D : array 3D a 3D array (e.g. UT, THT) - + orography : array 2D 2D orography (ZS) - + ztop : real scalar of the top height of the model (ZTOP) - + level : array 1D 1D level variable (level or level_w) - - n_x : array 1D + + n_x : array 1D 1D directionnal grid variable along i (ni_u, or ni_v) - - n_y : array 1D + + n_y : array 1D 1D directionnal grid variable along j (nj_u, or nj_v) - + Returns ------- altitude a 3D altitude variable with topography taken into account - + n_x3D a 3D directionnal variable duplicated from n_x - + n_y3D a 3D directionnal variable duplicated from n_y """ n_x3D = copy.deepcopy(oneVar3D) n_y3D = copy.deepcopy(oneVar3D) - altitude = copy.deepcopy(oneVar3D) + altitude = copy.deepcopy(oneVar3D) for i in range(len(level)): - n_y3D[i,:] = n_y - n_x3D[i,:] = n_x + n_y3D[i, :] = n_y + n_x3D[i, :] = n_x for i in range(oneVar3D.shape[2]): for j in range(oneVar3D.shape[1]): - if ztop==0: - altitude[:,i,j] = level[:] + if ztop == 0: + altitude[:, i, j] = level[:] else: for k in range(len(level)): - altitude[k,j,i] = orography[j,i] + level[k]*((ztop-orography[j,i])/ztop) + altitude[k, j, i] = orography[j, i] + level[k] * ((ztop - orography[j, i]) / ztop) return altitude, n_x3D, n_y3D diff --git a/src/LIB/Python/read_MNHfile.py b/src/LIB/Python/read_MNHfile.py index 38072becffe3b1d8c066c6f370683b4988cf7ccf..7e240f65c91f1abbb80576f99655cdad016f080d 100644 --- a/src/LIB/Python/read_MNHfile.py +++ b/src/LIB/Python/read_MNHfile.py @@ -11,33 +11,67 @@ MNH_LIC for details. version 1. import netCDF4 as nc import numpy as np -def read_withEPY(LnameFiles,Dvar_input, Dvar_output={}, path='.'): + +def read_withEPY(LnameFiles, Dvar_input, Dvar_output={}, path='.'): + """Read a netCDF4 Meso-NH file, LFI, FA or GRIB 1/2 file with the Epygram Meteo-France library + Parameters + ---------- + LnameFiles : list of str + list of Meso-NH netCDF4 file (diachronic or synchronous) + + Dvar_input : Dict{'keyFile' : [field_identifier]} + where + 'keyFile' is a str corresponding to a key for the file number in LnameFiles (by order) + [field_identifier] is + - GRIB1 [indicatorOfParameter, paramId, indicatorOfTypeOfLevel, level, casual name] for GRIB1 + - GRIB2 [discipline, parameterCategory, typeOfFirstFixedSurface, parameterNumber, level, casual name] for GRIB2 + - MNH netcdf ('variable string', level) for netcdf MNH; set level=0 for 2D variables + - LFI ('variable string', level) for LFI; set level=0 for 2D variables + - FA 'variable string' + + path : str or list + path(s) of the files to read + + Returns + ------- + Dvar : Dict + Dvar[ifile]['var_name'] an epygram list of resources + + """ import epygram epygram.init_env() - for i,keyFiles in enumerate(Dvar_input.keys()): + for i, keyFiles in enumerate(Dvar_input.keys()): + if isinstance(path, list): + ipath = path[i] + else: + ipath = path print('Reading file ' + keyFiles) - theFile = epygram.formats.resource(LnameFiles[i],'r') - Dvar_output[keyFiles] = {} #initialize dic for each files - for var in Dvar_input[keyFiles]: #For each files + theFile = epygram.formats.resource(ipath + LnameFiles[i], 'r') + Dvar_output[keyFiles] = {} # initialize dic for each files + for var in Dvar_input[keyFiles]: # For each files # Read variables if(theFile.format == 'FA'): Dvar_output[keyFiles][var] = theFile.readfield(var) elif(theFile.format == 'LFI'): - if(var[1]==None or var[1]==0): # 2D Field + if(var[1] is None or var[1] == 0): # 2D Field Dvar_output[keyFiles][var[0]] = theFile.readfield(var) - else: # 3D Field - Dvar_output[keyFiles][var[0]+str(var[1])] = theFile.readfield(var).getlevel(k=var[1]) + else: # 3D Field + Dvar_output[keyFiles][var[0] + str(var[1])] = theFile.readfield(var).getlevel(k=var[1]) elif(theFile.format == 'netCDFMNH'): - if(var[1]==None or var[1]==0): # 2D Field + if(var[1] is None or var[1] == 0): # 2D Field Dvar_output[keyFiles][var[0]] = theFile.readfield(var[0]) else: - Dvar_output[keyFiles][var[0]+str(var[1])] = theFile.readfield(var[0]).getlevel(k=var[1]) - elif(theFile.format == 'GRIB'): - if len(var)==6: # GRIB2 - Dvar_output[keyFiles][var[5]] = theFile.readfield({'discipline': var[0], 'parameterCategory': var[1], 'typeOfFirstFixedSurface': var[2],'parameterNumber': var[3], 'level': var[4]}) - elif len(var)==5: # GRIB1 - Dvar_output[keyFiles][var[4]] = theFile.readfield({'indicatorOfParameter': var[0], 'paramId': var[1], 'indicatorOfTypeOfLevel': var[2], 'level': var[3]}) - else: epygramError("GRIB format error. GRIB1 expects 4 values : [indicatorOfParameter, paramId, indicatorOfTypeOfLevel, level, 'casual name'], GRIB2 expects 5 values [discipline, parameterCategory, typeOfFirstFixedSurface, parameterNumber, level, casual name]") + Dvar_output[keyFiles][var[0] + str(var[1])] = theFile.readfield(var[0]).getlevel(k=var[1]) + elif(theFile.format == 'GRIB'): + if len(var) == 6: # GRIB2 + Dvar_output[keyFiles][var[5]] = theFile.readfield( + {'discipline': var[0], 'parameterCategory': var[1], 'typeOfFirstFixedSurface': var[2], 'parameterNumber': var[3], 'level': var[4]}) + elif len(var) == 5: # GRIB1 + Dvar_output[keyFiles][var[4]] = theFile.readfield( + {'indicatorOfParameter': var[0], 'paramId': var[1], 'indicatorOfTypeOfLevel': var[2], 'level': var[3]}) + else: + epygramError( + "GRIB format error. GRIB1 expects 4 values : [indicatorOfParameter, paramId, indicatorOfTypeOfLevel, level, 'casual name'], GRIB2 expects 5 values [discipline, parameterCategory, typeOfFirstFixedSurface, parameterNumber, level, casual name]") else: raise epygramError("Unknown format file, please use FA, LFI, GRIB or MNH NetCDF") theFile.close() @@ -45,22 +79,23 @@ def read_withEPY(LnameFiles,Dvar_input, Dvar_output={}, path='.'): # Transform spectral data to physics space (for AROME and ARPEGE) for f in Dvar_output: for var in Dvar_output[f]: - try: + try: if(Dvar_output[f][var].spectral): Dvar_output[f][var].sp2gp() - except: + except BaseException: break return Dvar_output + def read_netcdf(LnameFiles, Dvar_input, path='.', get_data_only=True, del_empty_dim=True, removeHALO=True): """Read a netCDF4 Meso-NH file For each file, call functions to read diachronic or synchronous file - + Parameters ---------- LnameFiles : list of str list of Meso-NH netCDF4 file (diachronic or synchronous) - + Dvar_input : Dict{'keyFile' : 'var_name',('group_name','var_name')} where 'keyFile' is a str corresponding to a key for the file number in LnameFiles (by order) @@ -69,211 +104,225 @@ def read_netcdf(LnameFiles, Dvar_input, path='.', get_data_only=True, del_empty_ e.g. : {'f1':['ZS', 'WT','ni', 'level'], 'f2':[('/LES_budgets/Cartesian/Not_time_averaged/Not_normalized/cart/',MEAN_TH'),('/Budgets/RI','AVEF')] } - - path : str - unique path of the files - + + path : str or list + path(s) of the files to read + get_data_only : bool, default: True if True, the function returns Dvar as masked_array (only data) if False, the function returns Dvar as netCDF4._netCDF4.Variable - + del_empty_dim : bool, default: True if get_data_only=True and del_empty_dim=True, returns Dvar as an array without dimensions with size 1 and 0 e.g. : an array of dimensions (time_budget, cart_level, cart_nj, cart_ni) with shape (180,1,50,1) is returned (180,50) - + removeHALO : bool, default: True - if True, remove first and last (NHALO=1) point [1:-1] if get_data_only=True on each + if True, remove first and last (NHALO=1) point [1:-1] if get_data_only=True on each level, level_w, ni, ni_u, ni_v, nj, nj_u, nj_v dimensions - + Returns ------- - Dvar : Dict + Dvar : Dict Dvar[ifile]['var_name'] if the group contains only one variable Dvar[ifile][('group_name','var_name')] if the group contains more than one variable """ Dvar = {} - for i,keyFiles in enumerate(Dvar_input.keys()): + for i, keyFiles in enumerate(Dvar_input.keys()): + if isinstance(path, list): + ipath = path[i] + else: + ipath = path print('Reading file ' + keyFiles) - print(path + LnameFiles[i]) - theFile = nc.Dataset(path + LnameFiles[i],'r') + print(ipath + LnameFiles[i]) + theFile = nc.Dataset(ipath + LnameFiles[i], 'r') Dvar[keyFiles] = {} - if '000' in LnameFiles[i][-6:-3]: + if '000' in LnameFiles[i][-6:-3]: if theFile['MASDEV'][0] <= 54: raise TypeError('The python lib is available for MNH >= 5.5') else: Dvar[keyFiles] = read_TIMESfiles_55(theFile, Dvar_input[keyFiles], Dvar[keyFiles], get_data_only, del_empty_dim, removeHALO) else: - Dvar[keyFiles]= read_BACKUPfile(theFile, Dvar_input[keyFiles], Dvar[keyFiles], get_data_only, del_empty_dim, removeHALO) - #theFile.close() + Dvar[keyFiles] = read_BACKUPfile(theFile, Dvar_input[keyFiles], Dvar[keyFiles], get_data_only, del_empty_dim, removeHALO) + # theFile.close() return Dvar + def read_var(theFile, Dvar, var_name, get_data_only, del_empty_dim=True, removeHALO=True): """Read a netCDF4 variable - + Parameters ---------- theFile : netCDF4._netCDF4.Dataset a Meso-NH diachronic netCDF4 file - + var_name : str a Meso-NH netCDF4 variable name - + get_data_only : bool, default: True if True, the function returns Dvar as masked_array (only data) if False, the function returns Dvar as netCDF4._netCDF4.Variable - + del_empty_dim : bool, default: True if get_data_only=True and del_empty_dim=True, returns Dvar as an array without dimensions with size 1 and 0 e.g. : an array of dimensions (time_budget, cart_level, cart_nj, cart_ni) with shape (180,1,50,1) is returned (180,50) - + removeHALO : bool, default: True - if True, remove first and last (NHALO=1) point [1:-1] if get_data_only=True on each + if True, remove first and last (NHALO=1) point [1:-1] if get_data_only=True on each level, level_w, ni, ni_u, ni_v, nj, nj_u, nj_v dimensions - + Returns ------- - Dvar : Dict + Dvar : Dict Dvar['var_name'] if the group contains only one variable Dvar[('group_name','var_name')] if the group contains more than one variable """ try: var_dim = theFile.variables[var_name].ndim var_dim_name = theFile.variables[var_name].dimensions - except: - raise KeyError("Group and variable name not found in the file, check the group/variable name with ncdump -h YourMNHFile.000.nc. You asked for variable : " + var_name) + except BaseException: + raise KeyError( + "Group and variable name not found in the file, check the group/variable name with ncdump -h YourMNHFile.000.nc. You asked for variable : " + + var_name) if not get_data_only: Dvar[var_name] = theFile.variables[var_name] - else: + else: if var_dim == 0: Dvar[var_name] = theFile.variables[var_name][0].data elif var_dim == 1: Dvar[var_name] = theFile.variables[var_name][:] elif var_dim == 2: - Dvar[var_name] = theFile.variables[var_name][:,:] + Dvar[var_name] = theFile.variables[var_name][:, :] elif var_dim == 3: - Dvar[var_name] = theFile.variables[var_name][:,:,:] + Dvar[var_name] = theFile.variables[var_name][:, :, :] elif var_dim == 4: - Dvar[var_name] = theFile.variables[var_name][:,:,:,:] + Dvar[var_name] = theFile.variables[var_name][:, :, :, :] elif var_dim == 5: - Dvar[var_name] = theFile.variables[var_name][:,:,:,:,:] + Dvar[var_name] = theFile.variables[var_name][:, :, :, :, :] elif var_dim == 6: - Dvar[var_name] = theFile.variables[var_name][:,:,:,:,:,:] + Dvar[var_name] = theFile.variables[var_name][:, :, :, :, :, :] elif var_dim == 7: - Dvar[var_name] = theFile.variables[var_name][:,:,:,:,:,:,:] + Dvar[var_name] = theFile.variables[var_name][:, :, :, :, :, :, :] if removeHALO: for i in range(8): try: - if var_dim_name[i]=='level' or var_dim_name[i]=='level_w' or \ - var_dim_name[i]=='ni' or var_dim_name[i]=='ni_u' or var_dim_name[i]=='ni_v' or \ - var_dim_name[i]=='nj' or var_dim_name[i]=='nj_u' or var_dim_name[i]=='nj_v': + if var_dim_name[i] == 'level' or var_dim_name[i] == 'level_w' or \ + var_dim_name[i] == 'ni' or var_dim_name[i] == 'ni_u' or var_dim_name[i] == 'ni_v' or \ + var_dim_name[i] == 'nj' or var_dim_name[i] == 'nj_u' or var_dim_name[i] == 'nj_v': if var_dim != 0: - Dvar[var_name] = removetheHALO(i+1, Dvar[var_name]) - except: + Dvar[var_name] = removetheHALO(i + 1, Dvar[var_name]) + except BaseException: break if del_empty_dim: - Ldimtosqueeze=[] + Ldimtosqueeze = [] var_shape = theFile.variables[var_name].shape for i in range(8): try: - if var_shape[i]==1: Ldimtosqueeze.append(i) + if var_shape[i] == 1: + Ldimtosqueeze.append(i) except IndexError: break - Ldimtosqueeze=tuple(Ldimtosqueeze) + Ldimtosqueeze = tuple(Ldimtosqueeze) Dvar[var_name] = np.squeeze(Dvar[var_name], axis=Ldimtosqueeze) return Dvar -def read_from_group(theFile, Dvar, group_name, var_name, get_data_only, del_empty_dim=True,removeHALO=True): - """Read a variable from a netCDF4 group - + +def read_from_group(theFile, Dvar, group_name, var_name, get_data_only, del_empty_dim=True, removeHALO=True): + """Read a variable from a netCDF4 group + Parameters ---------- theFile : netCDF4._netCDF4.Dataset a Meso-NH diachronic netCDF4 file - + group_name : str a Meso-NH netCDF4 group name - + var_name : str a Meso-NH netCDF4 variable name - + get_data_only : bool, default: True if True, the function returns Dvar as masked_array (only data) if False, the function returns Dvar as netCDF4._netCDF4.Variable - + del_empty_dim : bool, default: True if get_data_only=True and del_empty_dim=True, returns Dvar as an array without dimensions with size 1 and 0 e.g. : an array of dimensions (time_budget, cart_level, cart_nj, cart_ni) with shape (180,1,50,1) is returned (180,50) - + removeHALO : bool, default: True - if True, remove first and last (NHALO=1) point [1:-1] if get_data_only=True on each + if True, remove first and last (NHALO=1) point [1:-1] if get_data_only=True on each level, level_w, ni, ni_u, ni_v, nj, nj_u, nj_v dimensions - + Returns ------- - Dvar : Dict + Dvar : Dict Dvar['var_name'] if the group contains only one variable Dvar[('group_name','var_name')] if the group contains more than one variable """ try: var_dim = theFile[group_name].variables[var_name].ndim var_dim_name = theFile[group_name].variables[var_name].dimensions - except: - raise KeyError("Group and variable name not found in the file, check the group/variable name with ncdump -h YourMNHFile.000.nc. You asked for group/variable : " + group_name + var_name) - + except BaseException: + raise KeyError( + "Group and variable name not found in the file, check the group/variable name with ncdump -h YourMNHFile.000.nc. You asked for group/variable : " + + group_name + + var_name) + if not get_data_only: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name] else: if var_dim == 0: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][0].data + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][0].data if var_dim == 1: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:] elif var_dim == 2: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:,:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:, :] elif var_dim == 3: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:,:,:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:, :, :] elif var_dim == 4: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:,:,:,:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:, :, :, :] elif var_dim == 5: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:,:,:,:,:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:, :, :, :, :] elif var_dim == 6: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:,:,:,:,:,:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:, :, :, :, :, :] elif var_dim == 7: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:,:,:,:,:,:,:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:, :, :, :, :, :, :] if removeHALO: for i in range(8): try: - if var_dim_name[i]=='level' or var_dim_name[i]=='level_w' or \ - var_dim_name[i]=='ni' or var_dim_name[i]=='ni_u' or var_dim_name[i]=='ni_v' or \ - var_dim_name[i]=='nj' or var_dim_name[i]=='nj_u' or var_dim_name[i]=='nj_v': + if var_dim_name[i] == 'level' or var_dim_name[i] == 'level_w' or \ + var_dim_name[i] == 'ni' or var_dim_name[i] == 'ni_u' or var_dim_name[i] == 'ni_v' or \ + var_dim_name[i] == 'nj' or var_dim_name[i] == 'nj_u' or var_dim_name[i] == 'nj_v': if var_dim != 0: - Dvar[(group_name,var_name)] = removetheHALO(i+1, Dvar[(group_name,var_name)]) - except: + Dvar[(group_name, var_name)] = removetheHALO(i + 1, Dvar[(group_name, var_name)]) + except BaseException: break if del_empty_dim: - Ldimtosqueeze=[] - var_shape = Dvar[(group_name,var_name)].shape + Ldimtosqueeze = [] + var_shape = Dvar[(group_name, var_name)].shape for i in range(8): try: - if var_shape[i]==1: Ldimtosqueeze.append(i) + if var_shape[i] == 1: + Ldimtosqueeze.append(i) except IndexError: break - Ldimtosqueeze=tuple(Ldimtosqueeze) - Dvar[(group_name,var_name)] = np.squeeze(Dvar[(group_name,var_name)], axis=Ldimtosqueeze) - + Ldimtosqueeze = tuple(Ldimtosqueeze) + Dvar[(group_name, var_name)] = np.squeeze(Dvar[(group_name, var_name)], axis=Ldimtosqueeze) + # LES budget, ZTSERIES needs to be transposed to use psection functions without specifying .T each time if 'LES_budget' in group_name or 'ZTSERIES' in group_name or 'XTSERIES' in group_name: - Dvar[(group_name,var_name)] = Dvar[(group_name,var_name)].T + Dvar[(group_name, var_name)] = Dvar[(group_name, var_name)].T return Dvar + def read_BACKUPfile(theFile, Dvar_input, Dvar, get_data_only, del_empty_dim=True, removeHALO=True): """Read variables from Meso-NH MASDEV >= 5.5.0 synchronous file For all variables in Dvar_input of one file, call functions to read the variable of the group+variable - + Parameters ---------- theFile : netCDF4._netCDF4.Dataset a Meso-NH diachronic netCDF4 file - + Dvar_input : Dict{'var_name',('group_name','var_name')} with 'var_name' is the exact str of the netCDF4 variable name @@ -281,35 +330,40 @@ def read_BACKUPfile(theFile, Dvar_input, Dvar, get_data_only, del_empty_dim=True e.g. : {'f1':['ZS', 'WT','ni', 'level'], 'f2':[('/LES_budgets/Cartesian/Not_time_averaged/Not_normalized/cart/',MEAN_TH'),('/Budgets/RI','AVEF')] } - + get_data_only: bool, default: True if True, the function returns Dvar as masked_array (only data) if False, the function returns Dvar as netCDF4._netCDF4.Variable - + del_empty_dim: bool, default: True if get_data_only=True and del_empty_dim=True, returns Dvar as masked_array without dimensions with size 1 and 0 e.g. : an array of dimensions (time_budget, cart_level, cart_nj, cart_ni) with shape (180,1,50,1) is returned (180,50) - + Returns ------- - Dvar : Dict + Dvar : Dict Dvar['var_name'] if the group contains only one variable Dvar[('group_name','var_name')] if the group contains more than one variable - """ - # Reading date since beginning of the model run - Dvar['time'] = theFile.variables['time'][0] - Dvar['date'] = nc.num2date(Dvar['time'],units=theFile.variables['time'].units, calendar = theFile.variables['time'].calendar) + """ + # Reading date since beginning of the model run if theFile is a Meso-NH netCDF file + if 'MNH_cleanly_closed' in theFile.ncattrs(): + try: + Dvar['time'] = theFile.variables['time'][0] + Dvar['date'] = nc.num2date(Dvar['time'], units=theFile.variables['time'].units, calendar=theFile.variables['time'].calendar) + except BaseException: + pass for var in Dvar_input: - if type(var) == tuple: + if isinstance(var, tuple): Dvar = read_from_group(theFile, Dvar, var[0], var[1], get_data_only, del_empty_dim, removeHALO) else: Dvar = read_var(theFile, Dvar, var, get_data_only, del_empty_dim, removeHALO) # For all variables except scalars, change Fill_Value to NaN if get_data_only: - Dvar[var]= np.where(Dvar[var] != -99999.0, Dvar[var], np.nan) - Dvar[var]= np.where(Dvar[var] != 999.0, Dvar[var], np.nan) + Dvar[var] = np.where(Dvar[var] != -99999.0, Dvar[var], np.nan) + Dvar[var] = np.where(Dvar[var] != 999.0, Dvar[var], np.nan) return Dvar + def read_TIMESfiles_55(theFile, Dvar_input, Dvar, get_data_only, del_empty_dim=True, removeHALO=True): """Read variables from Meso-NH MASDEV >= 5.5.0 diachronic file For all variables in Dvar_input of one file, call functions to read the variable of the group+variable @@ -318,7 +372,7 @@ def read_TIMESfiles_55(theFile, Dvar_input, Dvar, get_data_only, del_empty_dim=T ---------- theFile : netCDF4._netCDF4.Dataset a Meso-NH diachronic netCDF4 file - + Dvar_input : Dict{'var_name',('group_name','var_name')} with 'var_name' is the exact str of the netCDF4 variable name @@ -326,57 +380,58 @@ def read_TIMESfiles_55(theFile, Dvar_input, Dvar, get_data_only, del_empty_dim=T e.g. : {'f1':['ZS', 'WT','ni', 'level'], 'f2':[('/LES_budgets/Cartesian/Not_time_averaged/Not_normalized/cart/',MEAN_TH'),('/Budgets/RI','AVEF')] } - + get_data_only: bool, default: True if True, the function returns Dvar as masked_array (only data) if False, the function returns Dvar as netCDF4._netCDF4.Variable - + del_empty_dim: bool, default: True if get_data_only=True and del_empty_dim=True, returns Dvar as masked_array without dimensions with size 1 and 0 e.g. : an array of dimensions (time_budget, cart_level, cart_nj, cart_ni) with shape (180,1,50,1) is returned (180,50) - + Returns ------- - Dvar : Dict + Dvar : Dict Dvar[ifile]['var_name'] if the group contains only one variable Dvar[ifile][('group_name','var_name')] if the group contains more than one variable - """ + """ for var in Dvar_input: - if type(var) == tuple: + if isinstance(var, tuple): Dvar = read_from_group(theFile, Dvar, var[0], var[1], get_data_only, del_empty_dim, removeHALO) else: Dvar = read_var(theFile, Dvar, var, get_data_only, del_empty_dim, removeHALO) return Dvar + def removetheHALO(idim, var): """Remove a NHALO=1 point [1:-1] at a given dimension idim of a variable var TODO: with NHALO=3 with WENO5 TODO: with the use of get_data_only=False (NetCDF4 output) - + Parameters ---------- idim: int the dimension over which remove the first and last point - + var: array a Meso-NH netCDF4 variable name - + Returns ------- - var : array - """ + var : array + """ if idim == 1: var = var[1:-1] elif idim == 2: - var = var[:,1:-1] + var = var[:, 1:-1] elif idim == 3: - var = var[:,:,1:-1] + var = var[:, :, 1:-1] elif idim == 4: - var = var[:,:,:,1:-1] + var = var[:, :, :, 1:-1] elif idim == 5: - var = var[:,:,:,:,1:-1] + var = var[:, :, :, :, 1:-1] elif idim == 6: - var = var[:,:,:,:,:,1:-1] + var = var[:, :, :, :, :, 1:-1] elif idim == 7: - var = var[:,:,:,:,:,:,1:-1] + var = var[:, :, :, :, :, :, 1:-1] return var