era5 数据 metpy.calc.wind_direction 中的错误

Error in metpy.calc.wind_direction with era5 data

我从 era-5 下载了 Era5 U 和 V 风组件,我正在使用 xarray 读取 .nc 文件和 select 数据中的几个经纬度点。在我需要使用 metpy.calc 函数计算风速和风向之后:

ds = xr.open_mfdataset('ERA5_u10_*.nc',combine='by_coords').metpy.parse_cf()
dv = xr.open_mfdataset('ERA5_v10_*.nc',combine='by_coords').metpy.parse_cf()
#direction =mpcalc.wind_direction(ds['u10'],dv['v10'])
#dd=mpcalc.wind_direction(ds.u10,dv.v10)
locations=[]
k=0
u10 = []
v10=[]
speed=[]
direction=[]
u100 = []
v100=[]
for i,j in zip(lats,lons):
    stations_u=ds.sel(longitude=j,latitude=i,method='nearest')
    stations_v=dv.sel(longitude=j,latitude=i,method='nearest')

    u10.append(stations_u.u10)### mudar variavel
    v10.append(stations_v.v10)### mudar variavel
    speed.append(mpcalc.wind_speed(stations_u.u10,stations_v.v10))
    direction.append(mpcalc.wind_direction(stations_u.u10,stations_v.v10))

风速正常,但风向出错:

runfile('/mnt/data1/ERA5/wind/untitled0.py', wdir='/mnt/data1/ERA5/wind')
Found valid latitude/longitude coordinates, assuming latitude_longitude for projection grid_mapping variable
Found valid latitude/longitude coordinates, assuming latitude_longitude for projection grid_mapping variable
Traceback (most recent call last):

  File "/usr/local/python/anaconda3/envs/my_env/lib/python3.6/site-packages/dask/array/core.py", line 1615, in __setitem__
    y = where(key, value, self)

  File "/usr/local/python/anaconda3/envs/my_env/lib/python3.6/site-packages/dask/array/routines.py", line 1382, in where
    return elemwise(np.where, condition, x, y)

  File "/usr/local/python/anaconda3/envs/my_env/lib/python3.6/site-packages/dask/array/core.py", line 4204, in elemwise
    broadcast_shapes(*shapes)

  File "/usr/local/python/anaconda3/envs/my_env/lib/python3.6/site-packages/dask/array/core.py", line 4165, in broadcast_shapes
    "shapes {0}".format(" ".join(map(str, shapes)))

ValueError: operands could not be broadcast together with shapes (262968,) (nan,) (262968,)


The above exception was the direct cause of the following exception:

Traceback (most recent call last):

  File "/mnt/data1/ERA5/wind/untitled0.py", line 40, in <module>
    direction.append(mpcalc.wind_direction(stations_u.u10,stations_v.v10))

  File "/usr/local/python/anaconda3/envs/my_env/lib/python3.6/site-packages/metpy/xarray.py", line 1206, in wrapper
    result = func(*bound_args.args, **bound_args.kwargs)

  File "/usr/local/python/anaconda3/envs/my_env/lib/python3.6/site-packages/metpy/units.py", line 246, in wrapper
    return func(*args, **kwargs)

  File "/usr/local/python/anaconda3/envs/my_env/lib/python3.6/site-packages/metpy/calc/basic.py", line 104, in wind_direction
    wdir[mask] += 360. * units.deg

  File "/usr/local/python/anaconda3/envs/my_env/lib/python3.6/site-packages/pint/quantity.py", line 1868, in __setitem__
    self._magnitude[key] = factor.magnitude

  File "/usr/local/python/anaconda3/envs/my_env/lib/python3.6/site-packages/dask/array/core.py", line 1622, in __setitem__
    ) from e

ValueError: Boolean index assignment in Dask expects equally shaped arrays.
Example: da1[da2] = da3 where da1.shape == (4,), da2.shape == (4,) and da3.shape == (4,).

如果我尝试在 U 和 V 中仅使用 1 个元素来计算 Wdir,令人惊讶的是它有效:

mpcalc.wind_direction(stations_u.u10[0],stations_v.v10[0]).values
Out[12]: array(173.41553, dtype=float32)

metpy.calc.wind_direction 不幸的是 known not to work with Dask arrays--and in fact many places in MetPy don't work well with Dask, though we definitely want to。现在,要使用 wind_direction,您需要使用例如将 stations_u 等转换为标准的 numpy 数组.compute().

感谢您的回复!是的,这是有道理的...

我正在使用的数据数组非常大,所以我需要将它保留为数据数组,因为将它复制到普通 np 数组会使用比我拥有的更多的内存 xD 我最终尝试使用 statins_u 作为 u 和 v 进行计算,令人惊讶的是它有效:

mpcalc.wind_direction(stations_u.u10,stations_u.u10)## Heading ##
Out[7]: 
<xarray.DataArray 'sub-404afe706bcfce1f36810c70766c0647' (time: 262968)>
<Quantity(dask.array<sub, shape=(262968,), dtype=float32, chunksize=(87672,), chunktype=numpy.ndarray>, 'degree')>
Coordinates:
    longitude  float32 -7.1
    latitude   float32 43.8
  * time       (time) datetime64[ns] 1990-01-01 ... 2019-12-31T23:00:00

为什么它在这种情况下有效?在这种情况下 metpy.calc_wind_direction 不使用 dask 数组吗? 我发现属性名称和单位有问题 pint,但无法理解它是什么...