从全球netCDF4文件提取特定国家tasmax数据及相关处理咨询

我从ESGF网站下载了1850-2014年的tasmax全球数据集（netCDF4格式），需要完成以下任务：

• 提取孟加拉国、斯里兰卡等特定国家的数据，用于Mann Kendall检验、Spearman秩相关检验、Sen's斜率等趋势分析（目前仅能提取单点数据，效率极低）
• 解决将数据导入QGIS时与世界矢量图不匹配、图层翻转的问题
• 将100km分辨率降为25km，寻求Python/QGIS/ArcGIS的解决方案

数据集详情 #

某包含1850-1899年数据的.nc文件信息：

table_id: day
table_info: Creation Date:(20 February 2019) MD5:951084b632bd52c3f6224e495b1cb65e
title: INM-CM4-8 output prepared for CMIP6
tracking_id: hdl:21.14100/819efd11-8094-46b4-86a9-1296bfe83b55
variable_id: tasmax

dict_keys(['time', 'time_bnds', 'lat', 'lat_bnds', 'lon', 'lon_bnds', 'height', 'tasmax']) 

dimensions(sizes): time(5475), lat(120), lon(180), bnds(2)

tasmax变量详情：

tasmax=ds.variables['tasmax']
print(tasmax)

输出结果：

<class 'netCDF4._netCDF4.Variable'>
float32 tasmax(time, lat, lon)
    standard_name: air_temperature
    long_name: Daily Maximum Near-Surface Air Temperature
    comment: maximum near-surface (usually, 2 meter) air 
temperature (add cell_method attribute 'time: max')
    units: K
    original_name: tasmax
    cell_methods: area: mean time: maximum (interval: 1 day)
    cell_measures: area: areacella
    history: 2019-05-28T11:29:24Z altered by CMOR: Treated 
scalar dimension: 'height'.
    coordinates: height
    missing_value: 1e+20
    _FillValue: 1e+20
unlimited dimensions: time
current shape = (5475, 120, 180)
filling on, default _FillValue of 9.969209968386869e+36 used


tasmax.dimentions
('time', 'lat', 'lon')

我曾尝试提取单点坐标数据：

print(temp[:,23.35,90])

输出：[282.4348 281.53607 283.49426 ... 280.90012 282.33377 282.48367]
但这种方式效率极低，无法覆盖整个国家范围。

一、批量提取特定国家数据 #

Python方案（推荐） #

用xarray+regionmask快速筛选国家范围数据，无需手动匹配坐标：

1. 安装依赖：

pip install xarray netCDF4 regionmask

2. 代码实现：

import xarray as xr
import regionmask

# 加载netCDF数据
ds = xr.open_dataset("your_file.nc")

# 调用内置国家边界数据集（1:110m分辨率）
countries = regionmask.defined_regions.natural_earth.countries_110

# 目标国家代码：孟加拉国(15)、斯里兰卡(142)，可通过`countries.names`查询所有国家
target_mask = countries[[15, 142]].mask(ds.lon, ds.lat)

# 筛选并保留国家范围内的网格点
ds_countries = ds.where(target_mask.notnull(), drop=True)

# 保存筛选后的数据集（可选）
ds_countries.to_netcdf("tasmax_bangladesh_srilanka.nc")

• 后续可直接基于xarray结合scipy/xarray-stats实现趋势检验，无需额外格式转换

QGIS方案 #

1. 将netCDF文件拖入QGIS，选择tasmax图层加载
2. 导入国家边界矢量数据（如自然地球1:10m国家数据集）
3. 打开处理工具箱 → 提取 → 按掩膜提取，选择tasmax栅格和国家边界作为掩膜
4. 导出提取后的栅格（可保存为GeoTIFF或netCDF格式）

ArcGIS方案 #

1. 通过多维工具 → 创建多维栅格图层加载netCDF数据
2. 加载国家边界矢量数据
3. 打开提取分析 → 按掩膜提取，选择多维栅格和国家边界完成提取
4. 可直接用多维分析工具对提取后的数据做趋势计算

二、解决QGIS图层匹配与翻转问题 #

图层翻转修复 #

• 打开tasmax图层属性 → 源，查看像素Y尺寸：若为负数，说明纬度是北到南倒序
• 在处理工具箱搜索翻转栅格，选择Y轴翻转，输出新图层即可修正

坐标系统匹配 #

• 查看世界矢量图的坐标参考系（CRS，通常为WGS84/EPSG:4326）
• 右键tasmax图层 → 设置CRS，选择与矢量图一致的CRS（若数据本身是经纬度坐标，直接设为EPSG:4326）
• 若仍不匹配，用投影变换工具将栅格转换为矢量图的CRS

三、分辨率从100km降为25km（降尺度） #

Python方案（xarray+rioxarray） #

import rioxarray as rxr
from rasterio.enums import Resampling

# 加载数据并设置CRS
ds = rxr.open_rasterio("tasmax_bangladesh_srilanka.nc").rio.set_crs("EPSG:4326")

# 缩放因子：100km转25km，分辨率提升4倍（原1°≈111km，对应0.25°≈27km，接近目标）
scale_factor = 4

# 双线性插值重采样
ds_resampled = ds.rio.reproject(
    ds.rio.crs,
    resolution=(ds.rio.resolution()[0]/scale_factor, ds.rio.resolution()[1]/scale_factor),
    resampling=Resampling.bilinear
)

# 保存结果
ds_resampled.rio.to_raster("tasmax_25km.nc")

• 重采样方法可根据需求选择：最近邻（适合分类数据）、双线性/立方卷积（适合连续气候数据）

QGIS方案 #

1. 加载需要降尺度的栅格图层
2. 打开处理工具箱 → 栅格 → 重采样
3. 设置输出分辨率为原分辨率的1/4（如原1°设为0.25°）
4. 选择重采样方法（推荐双线性插值），输出新图层

ArcGIS方案 #

1. 加载多维栅格图层
2. 打开数据管理工具 → 栅格 → 栅格处理 → 重采样
3. 设置输出分辨率为0.25°（对应约25km）
4. 选择重采样技术（如双线性插值）完成降尺度

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