【问题标题】:Calculating vorticity for multiple vertical levels in MetPy在 MetPy 中计算多个垂直水平的涡度
【发布时间】:2019-03-05 16:45:42
【问题描述】:

我正在尝试在 MetPy 中计算多个(连续)垂直级别的涡度。当我尝试计算单个级别时,一切正常。

这是代码;我使用了来自https://unidata.github.io/MetPy/latest/examples/cross_section.html#sphx-glr-examples-cross-section-py 的横截面示例。

import cartopy.crs as ccrs
import cartopy.feature as cfeature
import matplotlib.pyplot as plt
import numpy as np
import xarray as xr

import metpy.calc as mpcalc
from metpy.cbook import get_test_data
from metpy.interpolate import cross_section
from metpy.units import units

data = xr.open_dataset(get_test_data('narr_example.nc', False))
data = data.metpy.parse_cf().squeeze()

data_crs = data['Temperature'].metpy.cartopy_crs
lat = data['lat']
lon = data['lon']
f = mpcalc.coriolis_parameter(lat)
dx, dy = mpcalc.lat_lon_grid_deltas(lon, lat, initstring=data_crs.proj4_init)

然后进行涡度的计算。

vort = mpcalc.vorticity(data['u_wind'], data['v_wind'], dx, dy)

回溯:

Traceback (most recent call last):
  File "E:\Временные файлы\cross_section (1).py", line 63, in <module>
    vort = mpcalc.vorticity(data['u_wind'], data['v_wind'], dx, dy)
  File "C:\ProgramData\Miniconda3\lib\site-packages\metpy\xarray.py", line 436, in wrapper
    return func(*args, **kwargs)
  File "C:\ProgramData\Miniconda3\lib\site-packages\metpy\calc\kinematics.py", line 60, in wrapper
    ret = func(*args, **kwargs)
  File "C:\ProgramData\Miniconda3\lib\site-packages\metpy\calc\kinematics.py", line 121, in vorticity
    dudy = first_derivative(u, delta=dy, axis=-2)
  File "C:\ProgramData\Miniconda3\lib\site-packages\metpy\calc\tools.py", line 920, in wrapper
    return preprocess_xarray(func)(f, **kwargs)
  File "C:\ProgramData\Miniconda3\lib\site-packages\metpy\xarray.py", line 436, in wrapper
    return func(*args, **kwargs)
  File "C:\ProgramData\Miniconda3\lib\site-packages\metpy\calc\tools.py", line 1014, in first_derivative
    combined_delta = delta[tuple(delta_slice0)] + delta[tuple(delta_slice1)]
  File "C:\ProgramData\Miniconda3\lib\site-packages\pint\quantity.py", line 1400, in __getitem__
    value = self._magnitude[key]
IndexError: too many indices for array

我完全被困住了。搜索“metpy 多级计算”(没有实际报价)没有相关结果。 医生说:

metpy.calc.vorticity(u, v, dx, dy)[source]
Calculate the vertical vorticity of the horizontal wind.

Parameters: 
u ((M, N) ndarray) – x component of the wind
v ((M, N) ndarray) – y component of the wind
dx (float or ndarray) – The grid spacing(s) in the x-direction. If an array, there should be one item less than the size of u along the applicable axis.
dy (float or ndarray) – The grid spacing(s) in the y-direction. If an array, there should be one item less than the size of u along the applicable axis.
dim_order (str or None, optional) – The ordering of dimensions in passed in arrays. Can be one of None, 'xy', or 'yx'. 'xy' indicates that the dimension corresponding to x is the leading dimension, followed by y. 'yx' indicates that x is the last dimension, preceded by y. None indicates that the default ordering should be assumed, which is ‘yx’. Can only be passed as a keyword argument, i.e. func(…, dim_order=’xy’).
Returns:    
(M, N) ndarray – vertical vorticity

我的结论是输入可以有超过 2 维,但 3 维输入(在我的情况下)会出错。可以做些什么来修复它们?

我对 Python 完全陌生,所以我可能会犯一个愚蠢的错误。

【问题讨论】:

  • 你尝试过简单的二维切片数组吗?
  • @Alessandro,我使用了 MetPy 的 loc 命令。无论如何,@JonThielen 的回答可以完成这项工作。

标签: python metpy


【解决方案1】:

不幸的是,如果您不知道要查找什么,出现的错误消息在这种情况下并没有太大帮助!

在您的示例中,vorticity 函数调用的问题是您的输入变量的维度不匹配。 data['u_wind']data['v_wind'] 是形状为 (29, 118, 292) 的 3D 数组,但 dxdy,因为它们是从 lat_lon_grid_deltas 计算出来的,所以分别是形状为 (118, 291)(117, 292) 的二维数组。因此,我们需要获取适当广播的数组...您可以通过多种不同的方式来做到这一点,但这里有两个我推荐的选项:

选项 1:手动广播

由于dxdy 缺少的“额外”维度是第一个维度(垂直方向),我们可以通过插入大小将dxdy 制作成正确对齐的3D 数组-一个主要维度:

dx, dy = mpcalc.lat_lon_grid_deltas(lon, lat, initstring=data_crs.proj4_init)
dx = dx[None, :]
dy = dy[None, :]

vort = mpcalc.vorticity(data['u_wind'], data['v_wind'], dx, dy)

选项 2:使用 grid_deltas_from_dataarray() 辅助函数

MetPy 还有一个辅助函数,可以轻松地从 xarray DataArray 中提取网格增量。它还确保广播正确发生,因此您不必自己做。在您的示例中使用它,它将是:

dx, dy = mpcalc.grid_deltas_from_dataarray(data['u_wind'])

vort = mpcalc.vorticity(data['u_wind'], data['v_wind'], dx, dy)

【讨论】:

    猜你喜欢
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2021-02-08
    • 2012-03-04
    • 2021-08-16
    • 2017-03-14
    • 2015-12-01
    相关资源
    最近更新 更多