查找 NetCDF 中多个经纬度中心的半径内的值

Question

我有一个 netCDF 文件，其中包含特定时间南半球的多个气旋位置（纬度、经度）和气温。

\n

我想要的是提取距离每个气旋位置中心 10 测地度（约 1110 公里）半径内的温度值 \xe2\x80\x8b\xe2\x80\x8b 。这个想法是确定与每个气旋相关的温度值 \xe2\x80\x8b\xe2\x80\x8b（假设距气旋中心的最大径向距离为 10\xc2\xba），并仅用这些温度绘制一个全局轮廓图价值观。

\n

我在这里搜索了很多，但我只找到了适用于距一个特定纬度中心的距离的代码（例如：如何查找从纬度和经度值的中心位置开始的半径内的值）。

\n

我陷入了如何同时将半正矢公式应用于多个中心的困境。

\n

import xarray as xr\nimport numpy as np\nimport matplotlib.pyplot as plt\n\nd = xr.open_dataset(\'cyc_temp.nc\')\nlat = d[\'lat\']\nlon = d[\'lon\']\ncyc_pos = d[\'id\'][:,:]\ntemp = d[\'temp\'][:,:]\n\n# Haversine formula\n\ndef haversine(lon1, lat1, lon2, lat2):\n    # convert decimal degrees to radians\n    lon1 = np.deg2rad(lon1)\n    lon2 = np.deg2rad(lon2)\n    lat1 = np.deg2rad(lat1)\n    lat2 = np.deg2rad(lat2)\n\n    # haversine formula\n    dlon = lon2 - lon1\n    dlat = lat2 - lat1\n    a = np.sin(dlat/2)**2 + np.cos(lat1) * np.cos(lat2) * np.sin(dlon/2)**2\n    c = 2 * np.arcsin(np.sqrt(a))\n    r = 6371\n    return c * r\n

\n

如果有人可以帮助我，我将不胜感激。

\n

Answer 1

这是一个有趣的问题；xarray 的自动广播使这变得非常干净。

我不确定旋风分离器位置数组的结构如何，但我假设它的结构如下（或者至少可以操纵为这种形式）：

centers = np.array([[12.0, -62.0], [40.0, -80.0], [67.0, -55.0]])
cyc_pos = xr.DataArray(centers, coords={"coord": ["lon", "lat"]}, dims=["cyclone", "coord"])

换句话说，每一行代表每个气旋的经度和纬度值。

通过cyc_pos这种方式定义，使用该函数获取经纬度网格中每个点到每个气旋中心的距离haversine相当简单，从那里获取所需的掩模只需多一行。

distances = haversine(cyc_pos.sel(coord="lon"), cyc_pos.sel(coord="lat"), lon, lat)

如果您想要针对特定​​风暴的面罩，您可以使用：

storm_id = 0
mask = (distances <= 1110.0).isel(cyclone=storm_id)

或者，如果您想要一个可以应对所有风暴的面罩：

mask = (distances <= 1110.0).any("cyclone")