一、数据下载
#中国1km逐月平均气温数据集（1901-2025）https://www.tpdc.ac.cn/zh-hans/data/71ab4677-b66c-4fd1-a004-b2a541c4d5bf
#中国1km逐月降水量数据集（1901-2025）https://www.tpdc.ac.cn/zh-hans/data/faae7605-a0f2-4d18-b28f-5cee413766a2
二、Panoply查看NC信息
三、Python提取逐月逐年栅格
基于python提取逐月栅格，并合成年数据，以下是代码和空间分布草图。
import os, globimport numpy as npimport xarray as xrimport geopandas as gpdimport rioxarrayout_dir = r"E:\tmppeng\TIF"os.makedirs(out_dir, exist_ok=True)china = gpd.read_file(r"E:\研究区\中国标准行政区划数据GS（2024）0650号\shp格式\中国边界.shp")if china.crs.to_string() != "EPSG:4326": china = china.to_crs("EPSG:4326")mask = china.geometry.values.tolist()for f in sorted(glob.glob(r"E:\tmppeng\tmp_*.nc")): year = os.path.basename(f).split("_")[1].split(".")[0] print(f"处理 {year} ...") ds = xr.open_dataset(f) data = ds["tmp"].values * 0.1 lon, lat = ds["lon"].values, ds["lat"].values ds.close() for m in range(12): da = xr.DataArray(data[m], dims=["y", "x"], coords={"y": lat, "x": lon}) da.rio.write_crs("EPSG:4326").rio.clip(mask, crs="EPSG:4326", drop=True).rio.to_raster( os.path.join(out_dir, f"tmp_{year}_{m+1:02d}.tif"), nodata=-32768) annual = np.nanmean(data, axis=0) da = xr.DataArray(annual, dims=["y", "x"], coords={"y": lat, "x": lon}) da.rio.write_crs("EPSG:4326").rio.clip(mask, crs="EPSG:4326", drop=True).rio.to_raster( os.path.join(out_dir, f"tmp_{year}_annual.tif"), nodata=-32768) print(f" {year} 完成")print("全部完成。")
import os, globimport numpy as npimport xarray as xrimport geopandas as gpdimport rioxarrayout_dir = r"E:\Pre_peng\TIF"os.makedirs(out_dir, exist_ok=True)china = gpd.read_file(r"E:\研究区\中国标准行政区划数据GS（2024）0650号\shp格式\中国边界.shp")if china.crs.to_string() != "EPSG:4326": china = china.to_crs("EPSG:4326")mask = china.geometry.values.tolist()for f in sorted(glob.glob(r"E:\Pre_peng\pre_*.nc")): year = os.path.basename(f).split("_")[1].split(".")[0] print(f"处理 {year} ...") ds = xr.open_dataset(f) data = ds["pre"].values * 0.1 lon, lat = ds["lon"].values, ds["lat"].values ds.close() for m in range(12): da = xr.DataArray(data[m], dims=["y", "x"], coords={"y": lat, "x": lon}) da.rio.write_crs("EPSG:4326").rio.clip(mask, crs="EPSG:4326", drop=True).rio.to_raster( os.path.join(out_dir, f"pre_{year}_{m+1:02d}.tif"), nodata=-32768) annual = np.nansum(data, axis=0) annual[np.all(np.isnan(data), axis=0)] = np.nan da = xr.DataArray(annual, dims=["y", "x"], coords={"y": lat, "x": lon}) da.rio.write_crs("EPSG:4326").rio.clip(mask, crs="EPSG:4326", drop=True).rio.to_raster( os.path.join(out_dir, f"pre_{year}_annual.tif"), nodata=-32768) print(f" {year} 完成")print("全部完成。")
