Skew-T

pyMeteo documentation

Contents:

Skew-T/Log-P plotting, meteorology and interfacing with CM1

This package provides routines for plotting Skew-T/Log-P diagrams and working with CM1 output. For contributions, please submit pull requests to https://github.com/cwebster2/pymeteo. Please report bugs using the github issue tracker at https://github.com/cwebster2/pymeteo/issues.

image.png

EXAMPLE:Plotting Skew-T diagrams in Python

---

"""
===========================================================
SkewT-logP diagram: using transforms and custom projections
===========================================================

This serves as an intensive exercise of matplotlib's transforms and custom
projection API. This example produces a so-called SkewT-logP diagram, which is
a common plot in meteorology for displaying vertical profiles of temperature.
As far as matplotlib is concerned, the complexity comes from having X and Y
axes that are not orthogonal. This is handled by including a skew component to
the basic Axes transforms. Additional complexity comes in handling the fact
that the upper and lower X-axes have different data ranges, which necessitates
a bunch of custom classes for ticks,spines, and the axis to handle this.

"""

from matplotlib.axes import Axes
import matplotlib.transforms as transforms
import matplotlib.axis as maxis
import matplotlib.spines as mspines
from matplotlib.projections import register_projection


# The sole purpose of this class is to look at the upper, lower, or total
# interval as appropriate and see what parts of the tick to draw, if any.
class SkewXTick(maxis.XTick):
    def update_position(self, loc):
        # This ensures that the new value of the location is set before
        # any other updates take place
        self._loc = loc
        super(SkewXTick, self).update_position(loc)

    def _has_default_loc(self):
        return self.get_loc() is None

    def _need_lower(self):
        return (self._has_default_loc() or
                transforms.interval_contains(self.axes.lower_xlim,
                                             self.get_loc()))

    def _need_upper(self):
        return (self._has_default_loc() or
                transforms.interval_contains(self.axes.upper_xlim,
                                             self.get_loc()))

    @property
    def gridOn(self):
        return (self._gridOn and (self._has_default_loc() or
                transforms.interval_contains(self.get_view_interval(),
                                             self.get_loc())))

    @gridOn.setter
    def gridOn(self, value):
        self._gridOn = value

    @property
    def tick1On(self):
        return self._tick1On and self._need_lower()

    @tick1On.setter
    def tick1On(self, value):
        self._tick1On = value

    @property
    def label1On(self):
        return self._label1On and self._need_lower()

    @label1On.setter
    def label1On(self, value):
        self._label1On = value

    @property
    def tick2On(self):
        return self._tick2On and self._need_upper()

    @tick2On.setter
    def tick2On(self, value):
        self._tick2On = value

    @property
    def label2On(self):
        return self._label2On and self._need_upper()

    @label2On.setter
    def label2On(self, value):
        self._label2On = value

    def get_view_interval(self):
        return self.axes.xaxis.get_view_interval()


# This class exists to provide two separate sets of intervals to the tick,
# as well as create instances of the custom tick
class SkewXAxis(maxis.XAxis):
    def _get_tick(self, major):
        return SkewXTick(self.axes, None, '', major=major)

    def get_view_interval(self):
        return self.axes.upper_xlim[0], self.axes.lower_xlim[1]


# This class exists to calculate the separate data range of the
# upper X-axis and draw the spine there. It also provides this range
# to the X-axis artist for ticking and gridlines
class SkewSpine(mspines.Spine):
    def _adjust_location(self):
        pts = self._path.vertices
        if self.spine_type == 'top':
            pts[:, 0] = self.axes.upper_xlim
        else:
            pts[:, 0] = self.axes.lower_xlim


# This class handles registration of the skew-xaxes as a projection as well
# as setting up the appropriate transformations. It also overrides standard
# spines and axes instances as appropriate.
class SkewXAxes(Axes):
    # The projection must specify a name.  This will be used be the
    # user to select the projection, i.e. ``subplot(111,
    # projection='skewx')``.
    name = 'skewx'

    def _init_axis(self):
        # Taken from Axes and modified to use our modified X-axis
        self.xaxis = SkewXAxis(self)
        self.spines['top'].register_axis(self.xaxis)
        self.spines['bottom'].register_axis(self.xaxis)
        self.yaxis = maxis.YAxis(self)
        self.spines['left'].register_axis(self.yaxis)
        self.spines['right'].register_axis(self.yaxis)

    def _gen_axes_spines(self):
        spines = {'top': SkewSpine.linear_spine(self, 'top'),
                  'bottom': mspines.Spine.linear_spine(self, 'bottom'),
                  'left': mspines.Spine.linear_spine(self, 'left'),
                  'right': mspines.Spine.linear_spine(self, 'right')}
        return spines

    def _set_lim_and_transforms(self):
        """
        This is called once when the plot is created to set up all the
        transforms for the data, text and grids.
        """
        rot = 30

        # Get the standard transform setup from the Axes base class
        Axes._set_lim_and_transforms(self)

        # Need to put the skew in the middle, after the scale and limits,
        # but before the transAxes. This way, the skew is done in Axes
        # coordinates thus performing the transform around the proper origin
        # We keep the pre-transAxes transform around for other users, like the
        # spines for finding bounds
        self.transDataToAxes = self.transScale + \
            self.transLimits + transforms.Affine2D().skew_deg(rot, 0)

        # Create the full transform from Data to Pixels
        self.transData = self.transDataToAxes + self.transAxes

        # Blended transforms like this need to have the skewing applied using
        # both axes, in axes coords like before.
        self._xaxis_transform = (transforms.blended_transform_factory(
            self.transScale + self.transLimits,
            transforms.IdentityTransform()) +
            transforms.Affine2D().skew_deg(rot, 0)) + self.transAxes

    @property
    def lower_xlim(self):
        return self.axes.viewLim.intervalx

    @property
    def upper_xlim(self):
        pts = [[0., 1.], [1., 1.]]
        return self.transDataToAxes.inverted().transform(pts)[:, 0]


# Now register the projection with matplotlib so the user can select
# it.
register_projection(SkewXAxes)

if __name__ == '__main__':
    # Now make a simple example using the custom projection.
    from matplotlib.ticker import (MultipleLocator, NullFormatter,
                                   ScalarFormatter)
    import matplotlib.pyplot as plt
    from six import StringIO
    import numpy as np

    # Some examples data
    data_txt = '''
    841.08 1514.76 20.36 -9.33 12.57 NaN 115.43 4.38 NaN NaN NaN
    837.31 1553.84 19.82 -9.70 12.62 NaN 119.81 5.82 NaN NaN NaN
    832.20 1606.34 19.30 -9.92 12.81 NaN 122.46 6.20 NaN NaN NaN
    826.47 1665.29 18.76 -10.13 13.03 NaN 125.07 6.34 NaN NaN NaN
    820.16 1731.04 18.18 -10.36 13.28 NaN 127.81 6.41 NaN NaN NaN
    813.13 1804.35 17.53 -10.66 13.51 NaN 131.46 6.46 NaN NaN NaN
    805.34 1886.10 16.83 -11.08 13.66 NaN 136.52 6.78 NaN NaN NaN
    796.74 1977.28 16.12 -11.79 13.51 NaN 140.35 7.39 NaN NaN NaN
    787.22 2079.08 15.80 -12.47 13.05 NaN 137.56 6.74 NaN NaN NaN
    776.71 2192.77 15.28 -12.90 13.03 NaN 129.62 5.61 NaN NaN NaN
    765.12 2319.63 14.31 -13.95 12.74 NaN 119.53 4.54 NaN NaN NaN
    752.34 2461.01 12.79 -14.52 13.44 NaN 126.77 4.66 NaN NaN NaN
    738.30 2618.42 11.08 -15.43 13.96 NaN 140.59 4.80 NaN NaN NaN
    722.87 2793.92 10.10 -17.27 12.78 NaN 161.82 4.53 NaN NaN NaN
    705.98 2990.00 9.32 -18.02 12.64 NaN 172.01 3.65 NaN NaN NaN
    687.49 3209.10 8.22 -19.32 12.20 NaN 189.58 2.72 NaN NaN NaN
    667.34 3453.58 6.52 -21.23 11.62 NaN 212.09 2.50 NaN NaN NaN
    645.43 3726.18 4.55 -23.78 10.64 NaN 247.18 2.35 NaN NaN NaN
    621.67 4030.01 2.23 -27.22 9.19 NaN 284.45 2.96 NaN NaN NaN
    596.00 4368.56 -0.37 -31.88 7.15 NaN 311.85 4.51 NaN NaN NaN
    568.38 4745.81 -3.17 -37.99 4.80 NaN 329.14 7.25 NaN NaN NaN
    538.80 5166.03 -6.28 -42.48 3.80 NaN 338.22 11.46 NaN NaN NaN
    507.31 5633.85 -9.77 -42.60 4.93 NaN 339.77 17.20 NaN NaN NaN
    473.98 6154.49 -13.71 -39.63 9.23 NaN 335.29 23.77 NaN NaN NaN
    438.95 6733.42 -18.19 -37.18 17.24 NaN 327.79 29.99 NaN NaN NaN
    402.44 7376.36 -23.27 -37.34 26.36 NaN 319.57 34.03 NaN NaN NaN
    364.73 8089.49 -28.91 -41.51 28.65 NaN 313.96 34.62 NaN NaN NaN
    326.20 8879.15 -35.21 -47.90 26.15 NaN 316.09 32.98 NaN NaN NaN
    287.40 9750.00 -42.60 -53.14 30.56 NaN 323.32 32.35 NaN NaN NaN
    248.96 10704.66 -50.76 -57.78 42.99 NaN 322.31 31.86 NaN NaN NaN
    214.00 11676.93 -55.91 -64.39 33.73 NaN 310.77 28.91 NaN NaN NaN
    184.74 12608.61 -56.90 -69.48 18.87 NaN 291.23 29.02 NaN NaN NaN
    160.04 13517.79 -56.10 -75.10 7.50 NaN 278.57 33.69 NaN NaN NaN
    139.10 14406.03 -56.83 -80.45 1.57 NaN 267.59 35.76 NaN NaN NaN
    121.48 15255.93 -60.15 -80.45 1.69 NaN 258.59 33.93 NaN NaN NaN
    106.72 16055.47 -63.59 -80.45 3.69 NaN 250.17 30.95 NaN NaN NaN
    '''

    # Parse the data
    sound_data = StringIO(data_txt)
    p, h, T, Td = np.loadtxt(sound_data, usecols=range(0, 4), unpack=True)

    # Create a new figure. The dimensions here give a good aspect ratio
    fig = plt.figure(figsize=(6.5875, 6.2125))
    ax = fig.add_subplot(111, projection='skewx')

    plt.grid(True)

    # Plot the data using normal plotting functions, in this case using
    # log scaling in Y, as dictated by the typical meteorological plot
    ax.semilogy(T, p, color='C3')
    ax.semilogy(Td, p, color='C2')

    # An example of a slanted line at constant X
    l = ax.axvline(0, color='C0')

    # Disables the log-formatting that comes with semilogy
    ax.yaxis.set_major_formatter(ScalarFormatter())
    ax.yaxis.set_minor_formatter(NullFormatter())
    ax.set_yticks(np.linspace(100, 1000, 10))
    ax.set_ylim(1050, 100)
    ax.xaxis.set_major_locator(MultipleLocator(10))
    ax.set_xlim(-50, 50)
    plt.show()

image.png