GPM逐日降水数据文件合并

import xarray as xr
import glob
import os
from datetime import datetime

# 配置路径参数
input_dir = r"F:\GPM-processed\merged"
output_path = r"F:\GPM-processed\GPM_2001-2020_summer_daily.nc"
file_pattern = "*.nc"

# 获取并排序文件列表
def extract_date(filename):
    """修正后的日期提取函数"""
    basename = os.path.basename(filename)
    #parts = basename.split(".")
    #date_segment = parts[4]  # 正确索引位置
    #return datetime.strptime(date_segment.split("-")[0], "%Y%m%d")
    year_str = basename.split("_")[1]  # 根据实际文件名结构调整索引
    return int(year_str)

files = sorted(glob.glob(os.path.join(input_dir, file_pattern)),
               key=lambda x: extract_date(x))

# 验证文件数量
#print(f"发现 {len(files)} 个文件，预期92个")
#assert len(files) == 92, "文件数量不符合夏季天数"

print(f"发现 {len(files)} 个文件，预期20个")
assert len(files) == 20, "文件数量不符合夏季天数"

# 自定义预处理函数
def preprocess(ds):
    #filename = ds.encoding["source"]
    #date_str = os.path.basename(filename).split(".")[4].split("-")[0]
    #new_time = datetime.strptime(date_str, "%Y%m%d")
    #return ds.assign_coords(time=[new_time]).isel(time=0)


    year = ds.time.dt.year.item()  # 假设每个文件包含单一年度数据
    start_date = f"{year}-06-01"
    end_date = f"{year}-08-31"
    new_time = xr.date_range(start_date, end_date, freq="D")
    return ds.assign_coords(time=new_time)

# 合并数据集
combined = xr.open_mfdataset(
    files,
    combine="nested",
    concat_dim="time",
    preprocess=preprocess,
    parallel=True,
    chunks={"time": 10}
)

# 保存合并结果
combined.to_netcdf(output_path)
print(f"合并完成！文件已保存至：{output_path}")

# 获取并排序文件列表
def extract_year(filename):
    """从文件名中提取年份"""
    basename = os.path.basename(filename)
    # 假设文件名类似 GPM_2001_summer_daily.nc
    year_str = basename.split("_")[1]  # 根据实际文件名结构调整索引
    return int(year_str)


files = sorted(glob.glob(os.path.join(input_dir, file_pattern)),
               key=lambda x: extract_year(x))

# 验证文件数量
print(f"发现 {len(files)} 个文件，预期20个")
assert len(files) == 20, "文件数量不符合年度数量"


# 自定义预处理函数
def preprocess(ds):
    """统一时间坐标处理"""
    # 从文件名提取年份
    filename = ds.encoding["source"]
    year = extract_year(filename)

    # 重建时间坐标
    start_date = f"{year}-06-01"
    end_date = f"{year}-08-31"
    new_time = xr.date_range(start_date, end_date, freq="D")

    # 确保时间维度与数据对齐
    if len(new_time) != len(ds.time):
        raise ValueError(f"时间维度不匹配：文件 {filename} 的时间长度为 {len(ds.time)}，预期为 {len(new_time)}")

    return ds.assign_coords(time=new_time)


# 合并数据集
combined = xr.open_mfdataset(
    files,
    combine="nested",
    concat_dim="time",
    preprocess=preprocess,
    parallel=True,
    chunks={"time": 366}  # 按年分块
)

# 优化元数据
combined.time.attrs.update({
    "long_name": "Time",
    "description": "Daily data from June to August"
})

# 配置输出编码
encoding = {
    "precipitation": {
        "zlib": True,
        "complevel": 5,
        "dtype": "float32",
        "_FillValue": -9999.0
    },
    "time": {"dtype": "int32", "units": "days since 2001-01-01"}
}

# 保存合并结果
combined.to_netcdf(
    output_path,
    encoding=encoding,
    unlimited_dims=["time"]
)

print(f"合并完成！文件已保存至：{output_path}")

# 验证结果
print("\n合并后数据集摘要：")
print(xr.open_dataset(output_path))

import xarray as xr

gpm = xr.open_dataset("C:/Users/lenovo/Desktop/GPM_2001-2020_summer_daily.nc")
print(gpm)