Python matplotlib example wind speed and direction plot#

Software requirements:

  • Python 3

  • xarray

  • numpy

  • matplotlib

Script wind_dir_speed_circular_plot.py:

#!/usr/bin/env python
# coding: utf-8
#
# DKRZ example scripr:      wind_dir_circular_plot.py
#
# Plot the 10-days wind direction and speed average for Germany.
#           bar length:  wind speed
#
# 29.04.22  kmf

import xarray as xr
import numpy as np

import matplotlib.pyplot as plt
from matplotlib import cm
from matplotlib import colors
import matplotlib as mpl

def main():
    #-- Read the data set from file
    fname = '/Users/k204045/NCL/NCL_User_Guide/Version_1.0/data/rectilinear_grid_2D.nc'
    ds = xr.open_dataset(fname)
    #print(ds.info())

    #-- Next, we want to extract a smaller region (northern Germany) to take a closer
    #-- look at the data there.
    #-- First, we create the mask using the dataset coordinates.
    mask = ((ds.coords['lat'] > 47.25)
          & (ds.coords['lat'] < 55.5)
          & (ds.coords['lon'] > 5.8)
          & (ds.coords['lon'] < 15))

    #-- Now, use the mask to set all values outside the region to missing value.
    u = xr.where(mask, ds.u10, np.nan)
    v = xr.where(mask, ds.v10, np.nan)

    #-- Compute the wind direction
    wind_dir = np.arctan2(v,u) * (180/np.pi) + 180.

    #-- Compute the wind speed (magnitude)
    wind_speed = np.sqrt(u**2 + v**2)

    #-- Compute the averages
    wind_diravg   = np.mean(wind_dir,   axis=(0,1))
    wind_speedavg = np.mean(wind_speed, axis=(0,1))

    print('wind_diravg:         ', wind_diravg)
    print('')
    print('wind_speedavg:       ', wind_speedavg)
    print('wind_diravg.min():   ', wind_diravg.min().values)
    print('wind_diravg.max():   ', wind_diravg.max().values)
    print('')
    print('wind_speedavg.min(): ', wind_speedavg.min().values)
    print('wind_speedavg.max(): ', wind_speedavg.max().values)

    #-- Set levels
    levels_min  =  0.
    levels_max  = 12.
    levels_step =  0.25
    levels      = np.arange(levels_min, levels_max, levels_step)
    nlevels  = levels.size

    print('levels_min, levels_max, levels_step, nlevels: ',
           levels_min, levels_max, levels_step, nlevels)

    #-- Count values by wind_speed
    counts, bin_edges = np.histogram(wind_speedavg, bins=levels)

    max_counts = counts.max()

    print('counts       = ', counts)
    print('counts.max() = ', max_counts)
    print('counts.size  = ', counts.size)
    print('bin_edges    = ', bin_edges)

    #-- Decrease the number of colors to be used.
    colors = plt.cm.RdBu(wind_speedavg / levels_max)
    colors = plt.cm.gnuplot(wind_speedavg / levels_max)
    colors = plt.cm.Greys(wind_speedavg / levels_max)

    #-- Define angles array
    theta_step = 5
    theta_num  = int(360/theta_step)
    theta_arr  = np.radians(wind_diravg)

    #-- bar width
    bwidth = (2*np.pi) / theta_num

    # Create the circular bar plot
    #
    # The plot function `bar()` has to be used with the **polar projection** which
    # is set in the `subplots()` call. After drawing the bars, the color of each
    # bar is changed according to the wind speed.
    # Per default the angles are drawn counterclockwise and we have to set the 0
    # degree to top and the theta direction to clockwise.

    #plt.style.use('dark_background')

    fig, ax = plt.subplots(subplot_kw={'projection': 'polar'})

    bars = ax.bar(theta_arr, wind_speedavg,
                  color=colors,
                  width=bwidth,
                  bottom=0.,
                  alpha=0.9)

    ax.set_theta_zero_location('N')    # theta=0 at the top
    ax.set_theta_direction(-1)         # theta increasing clockwise

    #-- add geographic directions
    plt.gcf().text(0.5,   0.956,  'N', fontsize=14)
    plt.gcf().text(0.5,   0.016,  'S', fontsize=14)
    plt.gcf().text(0.853, 0.481,  'O', fontsize=14)
    plt.gcf().text(0.131, 0.481,  'W', fontsize=14)

    #-- add some text
    plt.gcf().text(0.02, 0.95,'Wind direction and speed', fontsize=12)
    plt.gcf().text(0.02, 0.9, 'Germany, 10-days average', fontsize=9)
    plt.gcf().text(0.02, 0.87,'2001-01-1 to 2001-01-10 (6 hourly data)', fontsize=6)

    #-- add units to plot
    plt.figtext(0.65, 0.8, ds.u10.attrs['units'], ha="right", fontsize=6)

    #-- add copyright
    plt.figtext(0.99, 0.01, "© 2021 DKRZ", ha="right", fontsize=5)

    #-- save plot to PNG file
    plt.savefig('plot_wind_dir_speed_circular.png', bbox_inches='tight', dpi=200)


if __name__ == '__main__':
    main()

Plot result:

image0