我从era-5下载了Era5 U和V风组件,我使用xarray读取.nc文件并从数据中选择几个纬度点。之后我需要用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,).
如果我尝试计算Wdir只有一个元素在U和在V中,令人惊讶的是它工作:
mpcalc.wind_direction(stations_u.u10[0],stations_v.v10[0]).values
Out[12]: array(173.41553, dtype=float32)
metpy.calc.wind_direction
很不幸地不能与Dask数组一起工作——事实上,MetPy中的许多地方都不能很好地与Dask一起工作,尽管我们绝对希望这样做。现在,要使用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数组在这种情况下?我觉得属性名称、单位和品脱有问题,但不明白是怎么回事…