如何用icclim/xclim计算不足365天时段的tn10p气候指标

我的目标是基于icclim包计算气候指标tn10p（即日最低气温Tmin低于第10百分位的天数占比），同时也尝试调用xclim包提供的同名指标进行计算。我需要计算跨自然年、时长不超过12个月的特定时段（例如'1960-12-01'至'1961-01-31'）的该指标值。

• 步骤1：读取时间分辨率为3小时的xarray数据集

t2m = xa.open_dataset('filepath.nc', decode_cf = True, decode_coords = "all").sel(time=slice('1960-12-01', '1961-01-31'))

• 步骤2：计算日最低气温值

t2m_min = t2m.t2m.resample(time='1D').min(keep_attrs = True)

• 步骤3.1：icclim实现代码

icclim_tn10p = icclim._generated_api.tn10p(in_files=t2m_min, slice_mode=['season',([12,1])])

• 步骤3.2：xClim实现代码

t2m_min_q10 = percentile_doy(arr = t2m_min, window=5, per=10).sel(percentiles=10)
xclim_tn10p = xclim.indicators.atmos.tn10p(tasmin = t2m_min, t10 = t2m_min_q10)

两种实现方式均抛出如下ValueError：

ValueError: conflicting sizes for dimension 'dayofyear': length 61 on <this-array> and length 365 on {'longitude': 'longitude', 'latitude': 'latitude', 'dayofyear': 'dayofyear', 'percentiles': 'percentiles'}

报错根源是percentile_doy函数需要输入包含完整365或366个日历日日序的数据计算百分位阈值，直接传入仅包含61天的短时段数据时，短序列的日序长度和函数要求的全年日序长度不匹配，该问题也对应xclim官方仓库的已知相关issue。

核心逻辑是不要直接截取目标短时段数据做百分位计算：百分位基准期计算必须使用覆盖完整日序（365/366天）的长时序数据，再将结果匹配到目标短时段上，具体操作如下：

• 第一步：读取完整基准期的长时序气温数据（至少覆盖要计算的百分位窗口对应的所有日序，通常建议用30年标准气候基准期数据，比如1991-2020年），不要一开始就截取1960-12到1961-01的短片段

# 读取完整长时序数据，不提前做时间切片
t2m_full = xa.open_dataset('filepath_full.nc', decode_cf=True, decode_coords="all")
# 基于完整长时序计算日最低温
t2m_min_full = t2m_full.t2m.resample(time='1D').min(keep_attrs=True)

• 第二步：在完整日最低温数据集上计算第10百分位阈值，此时percentile_doy拿到的是包含完整365/366天日序的数据，不会出现维度长度不匹配的问题

# xclim方案的百分位计算
t2m_min_q10_full = percentile_doy(arr=t2m_min_full, window=5, per=10).sel(percentiles=10)

# 如果用icclim，直接传入完整的日最低温数据，指定目标时段即可，不需要提前切片
icclim_tn10p = icclim.tn10p(
    in_files=t2m_min_full,
    time_range=['1960-12-01', '1961-01-31'],
    slice_mode=['season',([12,1])]
)

• 第三步：如果用xclim计算，从完整的日最低温和百分位阈值数据中，截取目标时段的数据再传入指标计算函数即可

# 截取目标时段的日最低温
t2m_min_target = t2m_min_full.sel(time=slice('1960-12-01', '1961-01-31'))
# 提取目标时段对应的dayofyear，匹配百分位数据的dayofyear维度
target_doys = t2m_min_target.time.dt.dayofyear
t2m_min_q10_target = t2m_min_q10_full.sel(dayofyear=target_doys)
# 计算tn10p
xclim_tn10p = xclim.indicators.atmos.tn10p(tasmin=t2m_min_target, t10=t2m_min_q10_target)

注意：如果原始数据本身只有1960-12到1961-01这两个月的片段，没有更长时间的基准期数据，无法正确计算tn10p。这类气候百分位指标的定义本身就要求用基准期全序列统计百分位阈值，单靠两个月短片段算出的阈值不符合指标定义，也无法绕过percentile_doy的日序维度校验。

内容的提问来源于stack exchange，提问作者ProgMau