python netCDF4不同维度的变量乘法

Question

我有第一个变量 h 具有维度 (111,141) 另一个变量 cs_w 具有维度 (51,) .基本上我的数据是 ROMS 历史输出数据。现在我想将 h 与 cs_w 相乘，最终结果的维度应该是 (51,111,141)。在这里，我被卡住了，无法继续。下面是我的代码

import numpy as np

import matplotlib.pyplot as plt

import netCDF4 as nc

f_in = nc.Dataset('ocean_his_0010.nc', "r")



h = f_in.variables['h']

cs_w = f_in.variables['Cs_w']


z=[[],[],[]]
for i in range(len(h[0])):
    for j in range(len(h[1])):
        for k in range(len(cs_w)):
            z[i][j][k] = h[i][j]*cs_w[k]

这是我要使用的两个变量的描述。

输出[88]： float64 Cs_w(s_w) long_name：在 W 点的 S 坐标拉伸曲线 有效最小值：-1.0 有效最大值：0.0 字段：Cs_w，标量 无限维度： 当前形状 = (51,) 填充，默认 _FillValue 使用 9.969209968386869e+36

h 出[89]： float64 h(eta_rho, xi_rho) long_name：RHO 点的水深测量 单位：米 网格：网格 地点：脸 坐标：lon_rho lat_rho 场：浴，标量 无限维度： 当前形状 = (111, 141) 填充，默认 _FillValue 使用 9.969209968386869e+36

Below is the ncdump ocean_his_0010.nc

netcdf ocean_his_0010 {
dimensions:
        xi_rho = 141 ;
        xi_u = 140 ;
        xi_v = 141 ;
        xi_psi = 140 ;
        eta_rho = 111 ;
        eta_u = 111 ;
        eta_v = 110 ;
        eta_psi = 110 ;
        N = 50 ;
        s_rho = 50 ;
        s_w = 51 ;
        tracer = 2 ;
        boundary = 4 ;
        ocean_time = UNLIMITED ; // (360 currently)

        double Cs_w(s_w) ;
                Cs_w:long_name = "S-coordinate stretching curves at W-points" ;
                Cs_w:valid_min = -1. ;
                Cs_w:valid_max = 0. ;
                Cs_w:field = "Cs_w, scalar" ;
        double h(eta_rho, xi_rho) ;
                h:long_name = "bathymetry at RHO-points" ;
                h:units = "meter" ;
                h:grid = "grid" ;
                h:location = "face" ;
                h:coordinates = "lon_rho lat_rho" ;
                h:field = "bath, scalar" ;

Answer 1

您不能只将（多维）列表定义为z=[[],[],[]] 并按照您尝试的方式开始填充它，首先需要适当地调整其大小。参见例如this question/answer 处理同样的问题。

Numpy 处理 nD 数组通常更方便，你的 z 数组可以简单地定义为：

z = np.zeros((51,111,141))

并使用诸如嵌套循环之类的内容或使用vectorized instructions 填充，例如：

for k in range(51):
    z[k,:,:] = cs_w[k] * h[:,:]

甚至全自动（甚至无需事先定义z）：

import numpy as np

h    = np.zeros((111,141))
cs_w = np.zeros(51)

z = cs_w[:,np.newaxis,np.newaxis] * h

使用这些矢量化操作通常比手动写出循环快很多。