Python修改NetCDF文件维度与数值并保存的实现方法

注意：你原有代码中new_values的初始化写法有误，np.full的第一个参数需要传入形状元组，第二个参数为填充值，正确初始化语句应为new_values = np.full((len(new_lev), len(lat), len(lon)), np.nan)，否则运行时会触发参数数量报错。

保存时建议输出到新文件，不要直接修改原始文件，避免原始数据损坏。完整代码如下：

import numpy as np
import netCDF4
from scipy.interpolate import interp1d

def interpolation(a,b,c):
    f = interp1d(a,b,kind='linear', bounds_error=False)
    return f(c)

# 第一步：读取原始文件的所有需要保留的数据，读完立刻关闭原文件
file = 'directory/test.nc'
data = netCDF4.Dataset(file)
# 读取经纬度数值和属性
lon = data.variables['lon'][:]
lat = data.variables['lat'][:]
lon_attrs = {k: data.variables['lon'].getncattr(k) for k in data.variables['lon'].ncattrs()}
lat_attrs = {k: data.variables['lat'].getncattr(k) for k in data.variables['lat'].ncattrs()}
# 读取原始垂直层和待插值变量
lev = data.variables['lev'][:]
values_attrs = {k: data.variables['values'].getncattr(k) for k in data.variables['values'].ncattrs()}
values = data.variables['values'][:]
# 关闭原文件，释放文件锁
data.close()

# 第二步：做插值计算
new_lev = np.linspace(0,1,200)
# 修正初始化错误
new_values = np.full((len(new_lev), len(lat), len(lon)), np.nan)
for loop_lat in range(len(lat)):
    for loop_lon in range(len(lon)):
        new_values[:, loop_lat, loop_lon] = interpolation(lev, values[:,loop_lat,loop_lon], new_lev)

# 第三步：写入新的NetCDF文件
out_file = 'directory/test_interp.nc'
with netCDF4.Dataset(out_file, 'w', format='NETCDF4') as nc_out:
    # 定义维度
    nc_out.createDimension('lon', len(lon))
    nc_out.createDimension('lat', len(lat))
    nc_out.createDimension('lev', len(new_lev))
    
    # 创建维度变量并写入值、保留原属性
    lon_var = nc_out.createVariable('lon', lon.dtype, ('lon',))
    lon_var[:] = lon
    for attr_name, attr_val in lon_attrs.items():
        lon_var.setncattr(attr_name, attr_val)
    
    lat_var = nc_out.createVariable('lat', lat.dtype, ('lat',))
    lat_var[:] = lat
    for attr_name, attr_val in lat_attrs.items():
        lat_var.setncattr(attr_name, attr_val)
    
    lev_var = nc_out.createVariable('lev', new_lev.dtype, ('lev',))
    lev_var[:] = new_lev
    # 可根据需要给新lev设置属性，比如单位
    lev_var.setncattr('units', 'level')
    lev_var.setncattr('long_name', 'vertical level')
    
    # 创建插值后的变量并写入
    values_var = nc_out.createVariable('values', new_values.dtype, ('lev', 'lat', 'lon'))
    values_var[:] = new_values
    # 保留原变量的属性
    for attr_name, attr_val in values_attrs.items():
        values_var.setncattr(attr_name, attr_val)

• 插值函数里加了bounds_error=False参数，避免新垂直层范围超出原始层范围时直接报错，超出范围的位置会返回填充值NaN，你也可以根据需求加fill_value参数指定外插填充值。
• 代码中保留了原始经纬度和目标变量的所有属性，不会丢失单位、坐标系统等元信息。
• 如果原文件还有其他依赖垂直维度的变量，需要全部完成插值后再写入，否则会出现维度不匹配的问题。
• 用with上下文管理器写文件可以自动在操作完成后关闭文件，避免文件损坏。

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