且构网

分享程序员开发的那些事...
且构网 - 分享程序员编程开发的那些事

如何摆脱Windrose pcolormap图python上的不连续性

更新时间:2023-02-26 22:31:58

您还可以为 wd_rad 提供小于 0 且大于 2* 的值pi 使用 np.where ,为小值添加 2 * pi ,为大值减去 2 * pi . np.tile(ws,2) np.tile(conc,2)然后将 wd_rad 的扩展版本与相同的浓度值相关联.在 interpolate.griddata 中也使用这些扩展值可确保浓度值环绕在 02*pi 处.

顺便说一句,请注意jet"是一个看起来不错的颜色图,但

如果你不想插值,但想绘制 10 度宽的线段来表示每个区域,

I am attempting to plot a windrose with binned concentration values. Following the advice from this post, and some modifications, I have created a plot. However, there is a discontinuity at around 0 deg. Any help would be greatly appreciated!

Here is my code:

wd = list(merge_all_apr['Wind Dir (10s deg)'])
conc = list(merge_all_apr['Mean_CO2'])
ws = list(merge_all_apr['Wind Spd (km/h)'])

wd_rad = np.radians(np.array(wd))
conc = np.array(conc, dtype=np.float)

wind_speed = np.linspace(min(ws), max(ws), 16)

WD, WS = np.meshgrid(np.linspace(0, 2*np.pi, 36), wind_speed)

print (WS)

Z = interpolate.griddata((wd_rad, ws), oz, (WD, WS), method='linear')

fig, ax = plt.subplots(subplot_kw={"projection": "polar"})
cmap = plt.get_cmap('jet')
#cmap.set_under('none')
img = ax.pcolormesh(WD, WS, Z, cmap=cmap, vmin=0, vmax=40, alpha = 0.70)
ax.set_theta_zero_location('N')
ax.set_theta_direction(-1) 
plt.colorbar(img)
plt.show()

And the result is:

As a scatter, it works fine and looks like

I'm not sure how to provide the data in a concise way, but any help would be greatly appreciated!!

You could also provide values for wd_rad just lower than 0 and larger than 2*pi using np.where, adding 2*pi for small values and subtracting 2*pi for large values. np.tile(ws, 2) and np.tile(conc, 2) then associate the extended version of wd_rad with the same concentration values. Using also these extended values in interpolate.griddata makes sure the concentration values wrap around at 0 and 2*pi.

As an aside, note that 'jet' is a colormap that looks nice, but is very misleading as it creates yellow highlights at the wrong spots. (Also, converting pandas columns to lists is quite slow and memory consuming, better leave them in their numpy array format.)

The code below supposes that oz in the question is the same array as conc.

import matplotlib.pyplot as plt
import numpy as np
from scipy import interpolate

# wd = merge_all_apr['Wind Dir (10s deg)']
# conc = merge_all_apr['Mean_CO2']
# ws = merge_all_apr['Wind Spd (km/h)']

N = 100
wd = np.random.uniform(0, 360, N)
conc = np.random.uniform(0, 40, N)
ws = np.random.uniform(0, 45, N)

wd_rad = np.radians(np.array(wd))
conc = np.array(conc, dtype=np.float)

wd_rad_ext = np.where(wd_rad < np.pi, wd_rad + 2 * np.pi, wd_rad - 2 * np.pi)

wind_speed = np.linspace(min(ws), max(ws), 16)

WD, WS = np.meshgrid(np.linspace(0, 2 * np.pi, 37), wind_speed)

Z = interpolate.griddata((np.hstack([wd_rad, wd_rad_ext]), np.tile(ws, 2)), np.tile(conc, 2),
                         (WD, WS), method='linear', rescale=True)

fig, axes = plt.subplots(ncols=2, figsize=(10, 4), subplot_kw={"projection": "polar"})
cmap = plt.get_cmap('magma')
for ax in axes:
    if ax == axes[0]:
        img = ax.pcolormesh(WD, WS, Z, cmap=cmap, vmin=0, vmax=40)
    else:
        img = ax.scatter(wd_rad, ws, c=conc, cmap=cmap, vmin=0, vmax=40)
    ax.set_theta_zero_location('N')
    ax.set_theta_direction(-1)
    plt.colorbar(img, ax=ax, pad=0.12)
plt.show()

If you don't want an interpolation, but want to draw segments of 10 degrees wide to represent each area, plt.hist2d can be employed.

Special parameters to hist2d:

  • bins=(np.linspace(0, 2 * np.pi, 37), np.linspace(min(ws), max(ws), 17)): the wind directions will be divided into 36 regions (37 boundaries); the speeds will be divided into 16 regions
  • weights=conc: instead of a usual histogram that counts the number of values into each little region, use the concentrations; when multiple concentrations are measured in the same little region, they are averaged
  • cmin=0.001: leave regions blank when their concentration value is less than 0.001
  • cmap = 'magma_r': use the reversed 'magma' colormap, so the high values get dark and the low get light colors (see the docs about other colormaps that might be more suitable to illustrate the data, but try not to use 'jet')

The return values of hist2d are a matrix of histogram values, the bin boundaries (x and y) and a collection of colored patches (which can be used as input for the colormap).

import matplotlib.pyplot as plt
import numpy as np
from scipy import interpolate

N = 100
wd = np.random.uniform(0, 360, N)
conc = np.random.uniform(0, 40, N)
ws = np.random.uniform(0, 45, N)

wd_rad = np.radians(np.array(wd))
conc = np.array(conc, dtype=np.float)

fig, axes = plt.subplots(ncols=2, figsize=(10, 4), subplot_kw={"projection": "polar"})
cmap = 'magma_r'
for ax in axes:
    if ax == axes[0]:
        _, _, _, img = ax.hist2d(wd_rad, ws, bins=(np.linspace(0, 2 * np.pi, 37), np.linspace(min(ws), max(ws), 17)),
                                 weights=conc, cmin=0.001, cmap=cmap, vmin=0, vmax=40)
    else:
        img = ax.scatter(wd_rad, ws, c=conc, cmap=cmap, vmin=0, vmax=40)
    ax.set_theta_zero_location('N')
    ax.set_theta_direction(-1)
    plt.colorbar(img, ax=ax, pad=0.12)
plt.show()