如何使用每日 NetCDF 数据识别或确定极端事件的连续网格单元

Question

我一直在处理我域中维度为 time:1322、lat:68、lon:136 的每日降水 netcdf 数据。对于每个时间步长（日期），我想识别/确定连续的网格单元，即日降水量超过我的阈值（第 95 个百分位数值）的日期中的网格单元，该阈值进一步连接到另一个更大的网格单元比阈值使用洪水填充算法功能

我已经编写并部署了 Flood Fill Algorithm 函数来访问每个单元格，如果单元格值高于阈值，但我不知道如何确定单元格是否 connected 到上，下，左或右的另一个单元格大于 threshold 的方向并指定值 1。同样，如果它小于 threshold，则用 NAN 或 0

替换单元格的值

我需要帮助/建议来完成此代码以执行上述任务并将结果作为新变量附加到加载的数据以进行进一步分析或其他合适的方法来执行此任务

这是我的代码：

import numpy as np
import xarray as xr

# Load the daily rainfall data in netcdf format
path_gpcc = '/path_to_my_file'
ds = xr.open_dataset(path_gpcc)

# Extract the rainfall variable
rainfall = ds.precip

# Define the Flood Fill Algorithm function
def flood_fill(t, x, y):
    """
    Recursive function to perform Flood Fill Algorithm
    """
    # Set current cell to visited
    visited[t, x, y] = True
    
    # Check if current cell exceeds the threshold and has not been visited
    if rainfall[t, x, y] >= threshold[x, y] and not visited[t, x, y]:
        # Check neighboring cells
        if x > 0:
            flood_fill(t, x-1, y)
        if x < rainfall.shape[1]-1:
            flood_fill(t, x+1, y)
        if y > 0:
            flood_fill(t, x, y-1)
        if y < rainfall.shape[2]-1:
            flood_fill(t, x, y+1)
            
# Calculate the 90th percentile value of the rainfall for each grid cell at all time steps
threshold = np.percentile(rainfall, 90, axis=0)

# Create an array to keep track of visited cells for all time steps
visited = np.zeros_like(rainfall, dtype=bool)

# List to store the contiguous grid cells
contiguous_cells = []

# Loop through each time step in the data period
for t in range(rainfall.shape[0]):
    # Loop through each cell in the grid for the current time step
    for x in range(rainfall.shape[1]):
        for y in range(rainfall.shape[2]):
            # Check if cell has been visited
            if not visited[t, x, y]:
                # Start Flood Fill Algorithm
                flood_fill(t, x, y)
                # Add contiguous grid cells to list
                contiguous_cells.append(np.where(visited))
                # Reset visited array
                visited = np.zeros_like(rainfall, dtype=bool)

# Create a new variable with the contiguous grid cells for all time steps
contiguous_cells_var = xr.Variable(dims=('time', 'lat', 'lon'), data=np.zeros_like(rainfall, dtype=int))
for cells in contiguous_cells:
    contiguous_cells_var[cells] = 1

PS：hilorywilsmart的方法和@Michael Delgado提供的解决方案在极端天气事件的连续网格单元分析中没有完全满足条件邻域