Python wind speed and direction plot#
Software requirements:
Python 3
xarray
numpy
matplotlib
Example script#
wind_dir_speed_circular_plot.py
#!/usr/bin/env python
# coding: utf-8
'''
DKRZ example
Wind rose
This example demonstrate how to plot the 10-days wind direction and speed
average for Germany.
Bar length: wind speed
-------------------------------------------------------------------------------
2022 copyright DKRZ licensed under CC BY-NC-SA 4.0
(https://creativecommons.org/licenses/by-nc-sa/4.0/deed.en)
-------------------------------------------------------------------------------
'''
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 = '../../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.switch_backend('agg')
#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()