python-3.x - 重新投影 Xarray 数据集

Question

我正在尝试将 Lambert Conformal 数据集重新投影到 Plate Carree。我知道这可以很容易地使用 cartopy 在视觉上完成。但是，我正在尝试创建一个新的数据集，而不仅仅是显示重新投影的图像。以下是我制定的方法，但我无法正确地对数据集进行子集化(Python 3.5、MacOSx)。

from siphon.catalog import TDSCatalog
import xarray as xr
from xarray.backends import NetCDF4DataStore
import numpy as np
import cartopy.crs as ccrs
from scipy.interpolate import griddata
import numpy.ma as ma
from pyproj import Proj, transform
import metpy

# Declare bounding box
min_lon = -78
min_lat = 36
max_lat = 40
max_lon = -72
boundinglat = [min_lat, max_lat]
boundinglon = [min_lon, max_lon]

# Load the dataset
cat = TDSCatalog('https://thredds.ucar.edu/thredds/catalog/grib/NCEP/HRRR/CONUS_2p5km/latest.xml')
dataset_name = sorted(cat.datasets.keys())[-1]
dataset = cat.datasets[dataset_name]
ds = dataset.remote_access(service='OPENDAP')
ds = NetCDF4DataStore(ds)
ds = xr.open_dataset(ds)

# parse the temperature at 850 and @ 0z reftime
tempiso = ds.metpy.parse_cf('Temperature_isobaric')
t850 = tempiso[0][2]

# transform bounding lat/lons to src_proj
src_proj = tempiso.metpy.cartopy_crs #aka lambert conformal conical
extents = src_proj.transform_points(ccrs.PlateCarree(), np.array(boundinglon), np.array(boundinglat))

# subset the data using the indexes of the closest values to the src_proj extents
t850_subset = t850[(np.abs(tempiso.y.values - extents[1][0])).argmin():(np.abs(tempiso.y.values - extents[1][1])).argmin()][(np.abs(tempiso.x.values - extents[0][1])).argmin():(np.abs(tempiso.x.values - extents[0][0])).argmin()]

# t850_subset should be a small, reshaped dataset, but it's shape is 0x2145
# now use nplinspace, npmeshgrid & scipy interpolate to reproject

我的变换点 > 查找最近值子集不起作用。它声称最近的点在数据集的范围之外。如前所述，我计划使用 nplinspace、npmeshgrid 和 scipy interpolate 从 t850_subset 创建一个新的方形纬度/经度数据集。

是否有更简单的方法来调整和重新投影 xarray 数据集？

Answer 1

您最简单的前进道路是利用 xarray 的能力来进行类似 pandas 的数据选择；这是 IMO xarray 的最佳部分。将最后两行替换为:

# By transposing the result of transform_points, we can unpack the
# x and y coordinates into individual arrays.
x_lim, y_lim, _ = src_proj.transform_points(ccrs.PlateCarree(),
    np.array(boundinglon), np.array(boundinglat)).T
t850_subset = t850.sel(x=slice(*x_lim), y=slice(*y_lim))

您可以在 xarray 的 selection and indexing functionality 文档中找到更多信息。 .您可能还会对 xarray 的 built-in support for interpolation 感兴趣.如果对 SciPy 之外的插值方法感兴趣，MetPy 也有一套 other interpolation methods .