Source code for epygram.fields.V1DField

#!/usr/bin/env python
# -*- coding: utf-8 -*-
# Copyright (c) Météo France (2014-)
# This software is governed by the CeCILL-C license under French law.
# http://www.cecill.info
"""
Contains the class that handle a Vertical 1D field.
"""

from __future__ import print_function, absolute_import, unicode_literals, division

import datetime
import numpy
import six

import footprints
from bronx.graphics.axes import set_figax, set_nice_time_axis

from .D3Field import D3CommonField, D3Field, D3VirtualField
from epygram import epygramError, config, util
from epygram.geometries import V1DGeometry

epylog = footprints.loggers.getLogger(__name__)


[docs]class V1DCommonField(D3CommonField): """ Vertical 1-Dimension (column) virtual or not field class. A field is defined by its identifier 'fid', its data, its geometry, and its validity. """ _collector = ('field',) _abstract = True _footprint = dict( attr=dict( structure=dict( values=set(['V1D'])), ) ) ################### # PRE-APPLICATIVE # ################### # (but useful and rather standard) ! # [so that, subject to continuation through updated versions, # including suggestions/developments by users...]
[docs] def plotfield(self, *args, **kwargs): """ Interface method to methods plotprofiles() and plotverticalhovmoller(), depending on time dimension (or not). Cf. these functions for arguments. """ if len(self.validity) == 1: return self.plotprofiles(*args, **kwargs) else: return self.plotverticalhovmoller(*args, **kwargs)
[docs] def plotprofiles(self, *args, **kwargs): """Cf. eponymous function of module for arguments.""" return plotprofiles(self, *args, **kwargs)
[docs] def plotverticalhovmoller(self, *args, **kwargs): """Cf. eponymous function of module for arguments.""" return plotverticalhovmoller(self, *args, **kwargs)
[docs] def plotanimation(self, *args, **kwargs): """Cf. eponymous function of module for arguments.""" return plotanimation(self, *args, **kwargs)
# FUNCTIONS # #############
[docs]def plotverticalhovmoller(profile, over=(None, None), fidkey=None, Ycoordinate=None, title=None, logscale=False, zoom=None, colorbar='vertical', graphicmode='colorshades', minmax=None, levelsnumber=21, center_cmap_on_0=False, colormap='jet', minmax_in_title=True, contourcolor='k', contourwidth=1, contourlabel=True, datefmt=None, showgrid=True, figsize=(6., 9.), rcparams=None): """ Makes a simple vertical Hovmöller plot of the field. :param profile: being a :class:`epygram.fields.V1DField` :param over: any existing figure and/or ax to be used for the plot, given as a tuple (fig, ax), with None for missing objects. *fig* is the frame of the matplotlib figure, containing eventually several subplots (axes); *ax* is the matplotlib axes on which the drawing is done. When given (is not None), these objects must be coherent, i.e. ax being one of the fig axes. :param fidkey: type of fid for entitling the plot with *fid[fidkey]*, if title is *None*; if *None*, labels with raw fid. :param Ycoordinate: label for the Y coordinate. :param title: title for the plot. :param logscale: to set Y logarithmic scale :param zoom: a dict containing optional limits to zoom on the plot. \n Syntax: e.g. {'ymax':500, ...}. :param colorbar: if *False*, hide colorbar the plot; else, befines the colorbar orientation, among ('horizontal', 'vertical'). Defaults to 'vertical'. :param graphicmode: among ('colorshades', 'contourlines'). :param minmax: defines the min and max values for the plot colorbar. \n Syntax: [min, max]. [0.0, max] also works. Default is min/max of the field. :param levelsnumber: number of levels for contours and colorbar. :param center_cmap_on_0: aligns the colormap center on the value 0. :param colormap: name of the **matplotlib** colormap to use. :param minmax_in_title: if True and minmax is not None, adds min and max values in title :param contourcolor: color or colormap to be used for 'contourlines' graphicmode. It can be either a legal html color name, or a colormap name. :param contourwidth: width of contours for 'contourlines' graphicmode. :param contourlabel: displays labels on contours. :param datefmt: date format to use, e.g. "%Y-%m-%d %H:%M:%S %Z" :param showgrid: True/False to show grid or not :param figsize: figure sizes in inches, e.g. (5, 8.5). If None, get the default figsize in config.plotsizes. :param rcparams: list of (*args, **kwargs) to be passed to pyplot.rc() defaults to [(('font',), dict(family='serif')),] Warning: requires **matplotlib**. """ import matplotlib import matplotlib.pyplot as plt from mpl_toolkits.axes_grid1 import make_axes_locatable if rcparams is None: rcparams = [(('font',), dict(family='serif')),] for args, kwargs in rcparams: plt.rc(*args, **kwargs) # User colormaps if colormap not in plt.colormaps(): util.load_cmap(colormap) # Figure, ax fig, ax = set_figax(*over, figsize=figsize) # coords z = numpy.zeros((len(profile.validity), len(profile.geometry.vcoordinate.levels))) for k in range(len(profile.geometry.vcoordinate.levels)): z[:, k] = profile.geometry.vcoordinate.levels[k] x = numpy.zeros((len(profile.validity), len(profile.geometry.vcoordinate.levels))) validities = {profile.validity[i].get():i for i in range(len(profile.validity))} xaxis_label = 'Validity' if len(validities) == 1 and len(profile.validity) != 1: xaxis_label = 'Basis' validities = {profile.validity[i].getbasis():i for i in len(profile.validity)} epoch = datetime.datetime(1970, 1, 1) for i in range(len(profile.validity)): d = profile.validity[i].get() if xaxis_label == 'Validity' else profile.validity[i].getbasis() timedelta = d - epoch p = (timedelta.microseconds + (timedelta.seconds + timedelta.days * 24 * 3600) * 1e6) / 1e6 x[i, :] = matplotlib.dates.epoch2num(p) data = profile.getdata() if profile.geometry.vcoordinate.typeoffirstfixedsurface in (119, 100): reverseY = True else: reverseY = False # min/max m = data.min() M = data.max() if minmax is not None: if minmax_in_title: minmax_in_title = '(min: ' + \ '{: .{precision}{type}}'.format(m, type='E', precision=3) + \ ' // max: ' + \ '{: .{precision}{type}}'.format(M, type='E', precision=3) + ')' try: m = float(minmax[0]) except Exception: m = data.min() try: M = float(minmax[1]) except Exception: M = data.max() else: minmax_in_title = '' if abs(m - M) > config.epsilon: levels = numpy.linspace(m, M, levelsnumber) else: raise epygramError("cannot plot uniform field.") if center_cmap_on_0: vmax = max(abs(m), M) vmin = -vmax else: vmin = m vmax = M L = int((levelsnumber - 1) // 15) + 1 hlevels = [levels[l] for l in range(len(levels) - L // 3) if l % L == 0] + [levels[-1]] # plot if reverseY and not ax.yaxis_inverted(): ax.invert_yaxis() if logscale: ax.set_yscale('log') ax.grid() if graphicmode == 'colorshades': pf = ax.contourf(x, z, data, levels, cmap=colormap, vmin=vmin, vmax=vmax) if colorbar: position = 'right' if colorbar == 'vertical' else 'bottom' cax = make_axes_locatable(ax).append_axes(position, size="5%", pad=0.1) cb = plt.colorbar(pf, orientation=colorbar, ticks=hlevels, cax=cax) if minmax_in_title != '': cb.set_label(minmax_in_title) elif graphicmode == 'contourlines': pf = ax.contour(x, z, data, levels=levels, colors=contourcolor, linewidths=contourwidth) if contourlabel: ax.clabel(pf, colors=contourcolor) # time set_nice_time_axis(ax, 'x', showgrid=showgrid, datefmt=datefmt) # decoration surf = z[-1, :] bottom = max(surf) if reverseY else min(surf) ax.fill_between(x[-1, :], surf, numpy.ones(len(surf)) * bottom, color='k') if Ycoordinate is None: if profile.geometry.vcoordinate.typeoffirstfixedsurface == 119: Ycoordinate = 'Level \nHybrid-Pressure \ncoordinate' elif profile.geometry.vcoordinate.typeoffirstfixedsurface == 100: Ycoordinate = 'Pressure (hPa)' elif profile.geometry.vcoordinate.typeoffirstfixedsurface == 102: Ycoordinate = 'Altitude (m)' elif profile.geometry.vcoordinate.typeoffirstfixedsurface == 103: Ycoordinate = 'Height (m)' elif profile.geometry.vcoordinate.typeoffirstfixedsurface == 118: Ycoordinate = 'Level \nHybrid-Height \ncoordinate' elif profile.geometry.vcoordinate.typeoffirstfixedsurface == 109: Ycoordinate = 'Potential \nvortex \n(PVU)' else: Ycoordinate = 'unknown \ncoordinate' ax.set_xlabel(xaxis_label) ax.set_ylabel(Ycoordinate) if zoom is not None: ykw = {} xkw = {} for pair in (('bottom', 'ymin'), ('top', 'ymax')): try: ykw[pair[0]] = zoom[pair[1]] except Exception: pass for pair in (('left', 'xmin'), ('right', 'xmax')): try: xkw[pair[0]] = zoom[pair[1]] except Exception: pass ax.set_ylim(**ykw) ax.set_xlim(**xkw) if title is None: if fidkey is None: fid = profile.fid[sorted(profile.fid.keys())[0]] else: fid = profile.fid[fidkey] title = u'Vertical Hovmöller of ' + str(fid) ax.set_title(title) return (fig, ax)
[docs]def plotprofiles(profiles, over=(None, None), labels=None, fidkey=None, Ycoordinate=None, unit='SI', title=None, logscale=False, ema=False, zoom=None, figsize=(6., 9.), rcparams=None): """ To plot a series of profiles. Returns a tuple of :mod:`matplotlib` (*Figure*, *ax*). :param profiles: a :class:`epygram.base.FieldSet` of :class:`epygram.fields.V1DField`, or a single :class:`epygram.fields.V1DField`. All profiles are supposed to have the same unit, and the same vertical coordinate. :param over: any existing figure and/or ax to be used for the plot, given as a tuple (fig, ax), with None for missing objects. *fig* is the frame of the matplotlib figure, containing eventually several subplots (axes); *ax* is the matplotlib axes on which the drawing is done. When given (is not None), these objects must be coherent, e.g. ax being one of the fig axes. :param labels: a list of labels for the profiles (same length and same order). :param fidkey: key of fid for labelling the curve with *fid[fidkey]*; if *None*, labels with raw fid. :param Ycoordinate: label for the Y coordinate. :param unit: label for X coordinate. :param title: title for the plot. :param logscale: to set Y logarithmic scale :param ema: to make emagram-like plots of Temperature :param zoom: a dict containing optional limits to zoom on the plot. \n Syntax: e.g. {'ymax':500, ...}. :param figsize: figure sizes in inches, e.g. (5, 8.5). If None, get the default figsize in config.plotsizes. :param rcparams: list of (*args, **kwargs) to be passed to pyplot.rc() defaults to [(('font',), dict(family='serif')), (('figure',), dict(autolayout=True))] """ import matplotlib.pyplot as plt if rcparams is None: rcparams = [(('font',), dict(family='serif')), (('figure',), dict(autolayout=True))] for args, kwargs in rcparams: plt.rc(*args, **kwargs) colors = ['red', 'blue', 'green', 'orange', 'magenta', 'darkolivegreen', 'yellow', 'salmon', 'black'] linestyles = ['-', '--', '-.', ':'] if isinstance(profiles, V1DField): profiles = [profiles] if isinstance(labels, six.string_types): labels = [labels] p0 = profiles[0] if p0.geometry.vcoordinate.typeoffirstfixedsurface in (119, 100): reverseY = True else: reverseY = False if p0.geometry.vcoordinate.typeoffirstfixedsurface in (118, 119): Y = p0.geometry.vcoordinate.levels if Ycoordinate is None: if p0.geometry.vcoordinate.typeoffirstfixedsurface == 119: Ycoordinate = 'Level \nHybrid-Pressure \ncoordinate' elif p0.geometry.vcoordinate.typeoffirstfixedsurface == 100: Ycoordinate = 'Pressure (hPa)' elif p0.geometry.vcoordinate.typeoffirstfixedsurface == 102: Ycoordinate = 'Altitude (m)' elif p0.geometry.vcoordinate.typeoffirstfixedsurface == 103: Ycoordinate = 'Height (m)' elif p0.geometry.vcoordinate.typeoffirstfixedsurface == 118: Ycoordinate = 'Level \nHybrid-Height \ncoordinate' elif p0.geometry.vcoordinate.typeoffirstfixedsurface == 109: Ycoordinate = 'Potential \nvortex \n(PVU)' else: Ycoordinate = 'unknown \ncoordinate' # Figure fig, ax = set_figax(*over, figsize=figsize) if logscale: ax.set_yscale('log') for i, p in enumerate(profiles): if len(p.validity) != 1: raise epygramError("plotprofiles can handle only profiles with one validity.") Y = numpy.array(p.geometry.vcoordinate.levels).flatten() if labels is not None: label = labels[i] else: if fidkey is not None: label = p.fid.get(fidkey, p.fid) else: label = str(p.fid) data = p.getdata() if ema: mindata = numpy.inf maxdata = -numpy.inf alpha = 0.75 templines = numpy.arange(round(min(data), -1) - 10, round(max(data), -1) + 10 + 50, 10) for t in templines: ax.plot([t, t + (max(data) - min(data)) * alpha], [max(Y), min(Y)], color='grey', linestyle=':') ax.set_yticks(numpy.linspace(0, 1000, 11)) data = data + abs(Y - Y[-1]) * (data.max() - data.min()) / \ (Y.max() - Y.min()) * alpha unit = 'K' mindata = min(mindata, data.min()) maxdata = max(maxdata, data.max()) plot_kwargs = {} if len(profiles) > 1: plot_kwargs['color'] = colors[i % len(colors)] plot_kwargs['linestyle'] = linestyles[i // len(colors)] ax.plot(data.flatten(), Y.flatten(), label=label, **plot_kwargs) if reverseY and not ax.yaxis_inverted(): ax.invert_yaxis() # Decoration if reverseY: ax.set_ylim(bottom=numpy.array(Y).max()) else: ax.set_ylim(bottom=numpy.array(Y).min()) if zoom is not None: if 'ymin' in zoom: ax.set_ylim(bottom=zoom['ymin']) if 'ymax' in zoom: ax.set_ylim(top=zoom['ymax']) if 'xmin' in zoom: ax.set_xlim(left=zoom['xmin']) if 'xmax' in zoom: ax.set_xlim(right=zoom['xmax']) if title is not None: ax.set_title(title) legend = ax.legend(loc='upper right', shadow=True) for label in legend.get_texts(): label.set_fontsize('medium') ax.set_xlabel(r'$' + unit + '$') ax.set_ylabel(Ycoordinate) if ema: ax.grid(axis='y') ax.set_xlim(mindata - 10, maxdata + 10) else: ax.grid() return (fig, ax)
[docs]def plotanimation(profile, title='__auto__', repeat=False, interval=1000, **kwargs): """ To plot a time-dependent profile as an animation. Returns a :class:`matplotlib.animation.FuncAnimation`. :param profile: the :class:`epygram.fields.V1DField` to plot. :param title: title for the plot. '__auto__' (default) will print the current validity of the time frame. :param repeat: to repeat animation :param interval: number of milliseconds between two validities Other kwargs passed to plotprofiles(). """ import matplotlib.animation as animation if len(profile.validity) == 1: raise epygramError("plotanimation can handle only profile with several validities.") if title is not None: if title == '__auto__': title_prefix = '' else: title_prefix = title title = title_prefix + '\n' + profile.validity[0].get().isoformat(sep=b' ') else: title_prefix = None profile0 = profile.getvalidity(0) mindata = profile.getdata().min() maxdata = profile.getdata().max() mindata -= (maxdata - mindata) / 10. maxdata += (maxdata - mindata) / 10. if kwargs.get('ema', False): epylog.warning("'ema' option not fully tested in animation: min/max may not be optimised.") zoom = kwargs.get('zoom') zoom = util.ifNone_emptydict(zoom) if 'xmax' not in zoom: zoom.update(xmax=maxdata) if 'xmin' not in zoom: zoom.update(xmin=mindata) kwargs['zoom'] = zoom fig, ax = plotprofiles(profile0, title=title, **kwargs) if kwargs.get('colorbar_over') is None: kwargs['colorbar_over'] = fig.axes[-1] # the last being created, in plotfield() kwargs['over'] = (fig, ax) def update(i, ax, myself, profilei, title_prefix, kwargs): if i < len(myself.validity): ax.clear() profilei = myself.getvalidity(i) if title_prefix is not None: title = title_prefix + '\n' + profilei.validity.get().isoformat(sep=b' ') profilei.plotfield(title=title, **kwargs) anim = animation.FuncAnimation(fig, update, fargs=[ax, profile, profile0, title_prefix, kwargs], frames=list(range(len(profile.validity) + 1)), # AM: don't really understand why but needed for the last frame to be shown interval=interval, repeat=repeat) return anim
[docs]class V1DField(V1DCommonField, D3Field): """ Vertical 1-Dimension (column) real field class. A field is defined by its identifier 'fid', its data, its geometry, and its validity. """ _collector = ('field',) _footprint = dict( attr=dict( structure=dict( values=set(['V1D'])), geometry=dict( type=V1DGeometry, access='rwx'), ) )
[docs]class V1DVirtualField(V1DCommonField, D3VirtualField): """ Vertical 1-Dimension (column) virtual field class. A field is defined by its identifier 'fid', its data, its geometry, and its validity. """ _collector = ('field',) _footprint = dict( attr=dict( structure=dict( values=set(['V1D'])), ) )