python netCDF4 variable multiplication of different dimension

Question

I have first variable h having dimension (111,141) another variable cs_w having dimension (51,). Basically my data is ROMS history output data. Now I want to multiply h with cs_w and the final result should be of dimension (51,111,141). Here I got stuck and unabble to proceed. below is my code

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]

Here is the description of the two variable which I want to use.

Out[88]: float64 Cs_w(s_w) long_name: S-coordinate stretching curves at W-points valid_min: -1.0 valid_max: 0.0 field: Cs_w, scalar unlimited dimensions: current shape = (51,) filling on, default _FillValue of 9.969209968386869e+36 used

h Out[89]: float64 h(eta_rho, xi_rho) long_name: bathymetry at RHO-points units: meter grid: grid location: face coordinates: lon_rho lat_rho field: bath, scalar unlimited dimensions: current shape = (111, 141) filling on, default _FillValue of 9.969209968386869e+36 used

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" ;

Bart · Accepted Answer · 2019-01-12T07:52:43.093

1

You can't just define a (multi-dimensional) list as z=[[],[],[]] and start filling it the way you are attempting, it needs to be properly sized first. See e.g. this question/answer which deals with the same problem.

Numpy is usually more convenient for dealing with nD-arrays, your z array can simply be defined as:

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

And be filled using either something like your nested loop, or using vectorized instructions, for example:

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

Or even fully automatic (without even having to define z beforehand):

import numpy as np

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

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

Using these vectorized operations is usually a lot faster than writing out the loop manually.

edited Jan 12 '19 at 07:52

answered Jan 11 '19 at 10:23

Bart

9,825
5
47
73

Thanks for your suggestion. The first comment is working properly, but the last I tried whic is giving below error z = cs_w[:,np.newaxis, np.newaxis]*h IndexError: only integers, slices (`:`), ellipsis (`...`), and 1-d integer or boolean arrays are valid indices – Navin Chandra Jan 12 '19 at 05:38
That's strange, just to be sure I added a full minimal example, can you try that? It runs fine on my system. – Bart Jan 12 '19 at 07:53
p.s.: I still use the `np.newaxis` approach, but I believe that [`np.broadcast_to`](https://docs.scipy.org/doc/numpy-1.10.0/reference/generated/numpy.broadcast_to.html) is now the preferred method. – Bart Jan 12 '19 at 08:00

python netCDF4 variable multiplication of different dimension

1 Answers1