根据特定索引屏蔽 xarray 或数据集

mask an xarray or dataset based on specific indices

我有一个网格化数据集(例如气温数据),我想在除海岸线以外的所有地方进行屏蔽(例如,海上 3 个网格单元,陆上 3 个网格单元,包含海岸线的 6 个网格单元缓冲区) .解决此问题的最佳方法(我认为)是使用包含 land/ocean 信息的掩码数组。

我的land/ocean面具在这里。不幸的是,我不知道如何通过代码重新创建这个文件(我已经下载了文件并在此处提供了文件详细信息)。

>>> mask.shape
(1, 81, 1440)
>>> mask
masked_array(
  data=[[[0.0000000e+00, 0.0000000e+00, 0.0000000e+00, ...,
          0.0000000e+00, 0.0000000e+00, 0.0000000e+00],
         [0.0000000e+00, 0.0000000e+00, 0.0000000e+00, ...,
          0.0000000e+00, 0.0000000e+00, 0.0000000e+00],
         [0.0000000e+00, 0.0000000e+00, 0.0000000e+00, ...,
          0.0000000e+00, 0.0000000e+00, 0.0000000e+00],
         ...,
         [9.4939685e-01, 7.2325826e-01, 4.7175312e-01, ...,
          9.9737060e-01, 9.8022592e-01, 9.5705426e-01],
         [7.6317883e-01, 5.6850266e-01, 2.9611492e-01, ...,
          9.9931705e-01, 9.8200846e-01, 9.4740820e-01],
         [3.4933090e-01, 1.7199135e-01, 3.2305717e-05, ...,
          9.5507932e-01, 7.9238594e-01, 5.7203436e-01]]],
  mask=False,
  fill_value=1e+20,
  dtype=float32)

我找到了如何获取掩码介于 0(海洋)和 1(陆地)之间的索引:

indices_btwn_0_1 = np.where(np.logical_and(mask>0,mask<1))
lat_indices_btwn_0_1 = indices_btwn_0_1[1] #
lon_indices_btwn_0_1 = indices_btwn_0_1[2]

现在我想我必须以某种方式使用这些指数来创建下面的数据集的子集,这样我就有了海岸线周围温度的子集。

ds = xr.open_dataset(ifile_climate_var_data)
>>> ds
<xarray.Dataset>
Dimensions:    (latitude: 81, longitude: 1440, time: 8760)
Coordinates:
  * longitude  (longitude) float32 -180.0 -179.75 -179.5 -179.25 -179.0 ...
  * latitude   (latitude) float32 85.0 84.75 84.5 84.25 84.0 83.75 83.5 ...
  * time       (time) datetime64[ns] 2019-01-01 2019-01-01T01:00:00 ...
Data variables:
    t2m        (time, latitude, longitude) float32 253.0178 252.98686 ...
Attributes:
    Conventions:  CF-1.6
    history:      2021-02-15 19:05:54 GMT by grib_to_netcdf-2.16.0: /opt/ecmw...

temp_2m = ds.t2m

我不明白的是如何从这里开始。我已经查看了高级索引的 xarray 文档,但我仍然不确定如何使用它来解决我的问题,我们将不胜感激。

我认为您在使用掩码方面走在正确的轨道上,但我认为 xarray 的高级索引不是您正在寻找的答案。

我建议看看 scipy.ndimage 中的 binary_dilationbinary_erosionhttps://docs.scipy.org/doc/scipy/reference/generated/scipy.ndimage.binary_dilation.html https://docs.scipy.org/doc/scipy/reference/generated/scipy.ndimage.binary_erosion.html

如果你用1s标记陆地,用0s标记海洋(或布尔值),你可以看到与scipy示例的对应关系:

from scipy import ndimage
a = np.zeros((7,7), dtype=int)
a[1:6, 2:5] = 1
a

所以原始数组看起来像:

array([[0, 0, 0, 0, 0, 0, 0],
       [0, 0, 1, 1, 1, 0, 0],
       [0, 0, 1, 1, 1, 0, 0],
       [0, 0, 1, 1, 1, 0, 0],
       [0, 0, 1, 1, 1, 0, 0],
       [0, 0, 1, 1, 1, 0, 0],
       [0, 0, 0, 0, 0, 0, 0]])

侵蚀的单次迭代:

ndimage.binary_erosion(a).astype(a.dtype)
array([[0, 0, 0, 0, 0, 0, 0],
       [0, 0, 0, 0, 0, 0, 0],
       [0, 0, 0, 1, 0, 0, 0],
       [0, 0, 0, 1, 0, 0, 0],
       [0, 0, 0, 1, 0, 0, 0],
       [0, 0, 0, 0, 0, 0, 0],
       [0, 0, 0, 0, 0, 0, 0]])

所以在你的例子中,你可以这样做:

is_land = (mask > 0)  # this might not be the right condition
eroded = is_land.copy(data=ndimage.binary_erosion(is_land.values, iterations=3)
dilated = is_land.copy(data=ndimage.binary_dilation(is_land.values, iterations=3)
coastal_zone = eroded | dilated