Xarray：在没有维度的情况下切片坐标

Question

Xarray：在没有维度的情况下切片坐标

9

虽然这个主题看起来应该很简单，但我仍然对它感到困难。

我想使用一组纬度和经度坐标来切割xarray数据集。

这是我的数据集的样子：

In [31]: data = xr.open_mfdataset(open_file, decode_cf=True)

In [32]: data
Out[32]: 
<xarray.Dataset>
Dimensions:  (time: 108120, x: 349, y: 277)
Coordinates:
    lons     (y, x) float64 -145.5 -145.3 -145.1 -144.9 -144.8 -144.6 -144.4 ...
    lats     (y, x) float64 1.0 1.104 1.208 1.312 1.416 1.519 1.621 1.724 ...
  * time     (time) datetime64[ns] 1980-01-01 1980-01-01T03:00:00 ...
Dimensions without coordinates: x, y
Data variables:
    stp      (time, y, x) float64 0.1235 0.0867 0.07183 0.05389 0.05901 ...

以下是我进行切片的步骤：

In [48]: lat_bnd = [25,30]
    ...: lon_bnd = [-80,-75]

In [49]: r = data.sel(y=slice(*lat_bnd),x=slice(*lon_bnd))

一切看起来都很棒：

In [50]: r
Out[50]: 
    <xarray.Dataset>
    Dimensions:  (time: 108120, x: 5, y: 5)
    Coordinates:
        lons     (y, x) float64 -82.52 -82.28 -82.05 -81.81 -81.57 -82.44 -82.2 ...
        lats     (y, x) float64 13.54 13.46 13.38 13.3 13.22 13.77 13.69 13.61 ...
      * time     (time) datetime64[ns] 1980-01-01 1980-01-01T03:00:00 ...
    Dimensions without coordinates: x, y
    Data variables:
        stp      (time, y, x) float64 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ...

但是我的纬度/经度数值不正确：

In [53]: r.lats.values
Out[53]: 
array([[ 13.53542397,  13.45647916,  13.37686013,  13.296571  ,
         13.21561592],
       [ 13.76719053,  13.6878189 ,  13.60776989,  13.52704767,
         13.44565641],
       [ 13.99938176,  13.91958109,  13.83909988,  13.75794233,
         13.67611265],
       [ 14.2319952 ,  14.15176326,  14.07084762,  13.98925249,
         13.90698214],
       [ 14.46502833,  14.3843629 ,  14.30301059,  14.22097564,
         14.13826236]])

In [54]: r.lons.values
Out[54]: 
array([[-82.52229969, -82.28438922, -82.0469968 , -81.8101255 ,
        -81.57377834],
       [-82.44118948, -82.20260881, -81.96455096, -81.72701901, -81.490016  ],
       [-82.3595596 , -82.12030558, -81.8815792 , -81.64338357,
        -81.40572174],
       [-82.27740522, -82.03747469, -81.79807668, -81.55921433,
        -81.32089068],
       [-82.19472148, -81.95411126, -81.71403851, -81.47450637, -81.2355179 ]])

当然，如果我尝试使用经纬度坐标进行切片，由于维度不匹配，会出现错误。

    In [55]: r = data.sel(lats=slice(*lat_bnd),lons=slice(*lon_bnd))
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-55-7c6237be5f22> in <module>()
----> 1 r = data.sel(lats=slice(*lat_bnd),lons=slice(*lon_bnd))

/lib/anaconda2/lib/python2.7/site-packages/xarray/core/dataset.pyc in sel(self, method, tolerance, drop, **indexers)
   1204         """
   1205         pos_indexers, new_indexes = indexing.remap_label_indexers(
-> 1206             self, indexers, method=method, tolerance=tolerance
   1207         )
   1208         result = self.isel(drop=drop, **pos_indexers)

/lib/anaconda2/lib/python2.7/site-packages/xarray/core/indexing.pyc in remap_label_indexers(data_obj, indexers, method, tolerance)
    275     new_indexes = {}
    276 
--> 277     dim_indexers = get_dim_indexers(data_obj, indexers)
    278     for dim, label in iteritems(dim_indexers):
    279         try:

/lib/anaconda2/lib/python2.7/site-packages/xarray/core/indexing.pyc in get_dim_indexers(data_obj, indexers)
    243     if invalid:
    244         raise ValueError("dimensions or multi-index levels %r do not exist"
--> 245                          % invalid)
    246 
    247     level_indexers = defaultdict(dict)

ValueError: dimensions or multi-index levels ['lons', 'lats'] do not exist

我对于这个NARR数据集的理解是否有所遗漏？

- Maria Molina

3个回答

2

也许有人会感兴趣 - 在使用xarray的最近邻插值功能无法获得收益时。我遇到了类似的问题 - 一个带有2D数组作为lons/lats坐标的曲线网格。此外，我正在寻找最接近给定点的坐标。我使用以下附加函数将“任何”lon/lat对转换为（x，y）元组：

def find_nearest(lons, lats, lon0,lat0):
   idx = ((lons - lon0)**2+(lats - lat0)**2).argmin()
   value_lat =  lats.flat[idx]
   return tuple(np.squeeze(np.where(lats2d == value_lat)))

然后，将其用作：

find_nearest(tmp.longitude.values,tmp.latitude.values, -22.16,32.3)

- Marcin Kawka

0

一种解决方法可能是基于 x、y 坐标进行切片。为了检查是否符合您的领域要求，您可以查看快速绘图并调整 x、y 值以相应地进行切片。但更好的方法是将您的经纬度转换为 x、y 坐标，然后根据相应的 x、y 进行切片。

- user11683337

1

亲爱的用户1168337，您的建议非常重要。然而，如何将lat/lon坐标转换为netcdf_xarray x,y？您能否提供一个实际的例子？此致。 - Philipe Riskalla Leal

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Dan · Accepted Answer

更新于2020-04-30

如果您想根据经纬度选择数据，可以使用where()做如下操作：

data.where((data.lats > 25) & (data.lats < 30) & (data.lons > -80) & (data.lons < -75))

您可以添加drop=True以返回一个较小的数据集，而不是用NA填充非匹配值。

原始回答

在您的第一个示例中，您没有按照经度/纬度进行索引，而是通过每个x和y的数字索引进行切片。也就是说，您正在切片第25到30个y和第-80到-75个x值之间。这就解释了为什么您的输出中经纬度值没有意义。

您可以使用xr.Dataset.set_index()将您的坐标与维度相关联，例如：

data = data.set_index(y='lats')
data = data.set_index(x='lons')