11. Utilities¶
11.1. MITgcmutils¶
This Python package includes a number of helpful functions and scripts for dealing with MITgcm output. You can install it from the model repository (in directory utils/python/MITgcmutils) or from the Python Package Index:
pip install --user MITgcmutils
The following functions are exposed at the package level:
from module mnc:
rdmnc()
andmnc_files()
from module ptracers:
iolabel()
and:iolabel2num()
from module diagnostics:
readstats()
The package also includes a standalone script for joining tiled mnc files: gluemncbig.
For more functions, see the individual modules:
11.1.1. mds¶
- MITgcmutils.mds.parsemeta(metafile)[source]¶
parses metafile (file object or filename) into a dictionary of lists of floats, ints or strings
- MITgcmutils.mds.rdmds(fnamearg, itrs=-1, machineformat='b', rec=None, fill_value=0, returnmeta=False, astype=<class 'float'>, region=None, lev=(), usememmap=False, mm=False, squeeze=True, verbose=False)[source]¶
Read meta-data files as written by MITgcm.
Call signatures:
a = rdmds(fname,…)
a,its,meta = rdmds(fname,…,returnmeta=True)
- Parameters
fname (string) –
name of file to read, without the ‘.data’ or ‘.meta’ suffix. If itrs is given, the iteration number is added to fname as well. fname may contain shell wildcards, which is useful for tile files organized into directories, e.g.,
T = rdmds(‘prefix*/T’, 2880)
will read prefix0000/T.0000002880.*, prefix0001/T.0000002880.*, … (and any others that match the wildcard, so be careful how you name things!)
itrs (int or list of ints or np.NaN or np.Inf) –
Iteration number(s). With itrs=-1, will try to read
fname.meta or fname.001.001.meta, …
If itrs is a list of integers of an integer, it will read the corresponding
fname.000000iter.meta, …
If itrs is np.NaN, it will read all iterations for which files are found. If itrs is np.Inf, it will read the highest iteration found.
machineformat (int) – endianness (‘b’ or ‘l’, default ‘b’)
rec (list of int or None) – list of records to read (default all) useful for pickups and multi-field diagnostics files
fill_value (float) – fill value for missing (blank) tiles (default 0)
astype (data type) – data type to return (default: double precision) None: keep data type/precision of file
region (tuple of int) – (x0,x1,y0,y1) read only this region (default (0,nx,0,ny))
lev (list of int or tuple of lists of int) – list of levels to read, or, for multiple dimensions (excluding x,y), tuple(!) of lists (see examples below)
usememmap (bool) – if True, use a memory map for reading data (default False) recommended when using lev, or region with global files to save memory and, possibly, time
- Returns
a (array_like) – numpy array of the data read
its (list of int) – list of iteration numbers read (only if returnmeta=True)
meta (dict) – dictionary of metadata (only if returnmeta=True)
Examples
>>> XC = rdmds('XC') >>> XC = rdmds('res_*/XC') >>> T = rdmds('T.0000002880') >>> T = rdmds('T',2880) >>> T2 = rdmds('T',[2880,5760]) >>> T,its = rdmds('T',numpy.Inf) >>> VVEL = rdmds('pickup',2880,rec=range(50,100)) >>> a5 = rdmds('diags',2880,rec=0,lev=[5]) >>> a = rdmds('diags',2880,rec=0,lev=([0],[0,1,5,6,7])) >>> from numpy import r_ >>> a = rdmds('diags',2880,rec=0,lev=([0],r_[:2,5:8])) # same as previous >>> a = rdmds('diags',2880,rec=0)[0, [0,1,5,6,7], ...] # same, but less efficient >>> a = rdmds('diags',2880)[0, 0, [0,1,5,6,7], ...] # even less efficient
- MITgcmutils.mds.scanforfiles(fname)[source]¶
return list of iteration numbers for which metafiles with base fname exist
- MITgcmutils.mds.wrmds(fbase, arr, itr=None, dataprec='float32', ndims=None, nrecords=None, times=None, fields=None, simulation=None, machineformat='b', deltat=None, dimlist=None)[source]¶
Write an array to an mds meta/data file set.
If itr is given, the files will be named fbase.0000000itr.data and fbase.0000000itr.meta, otherwise just fbase.data and fbase.meta.
- Parameters
fbase (string) – Name of file to write, without the ‘.data’ or ‘.meta’ suffixes, and without the iteration number if itr is give
arr (array_like) – Numpy array to write
itr (int or None) – If given, this iteration number will be appended to the file name
dataprec (string) – precision of resulting file (‘float32’ or ‘float64’)
ndims (int) – number of non-record dimensions; extra (leading) dimensions will be folded into 1 record dimension
nrecords (int) – number of records; will fold as many leading dimensions as necessary (has to match shape!)
times (float or list of floats) – times to write into meta file. Either a single float or a list of two for a time interval
fields (list of strings) – list of fields
simulation (string) – string describing the simulation
machineformat (string) – ‘b’ or ‘l’ for big or little endian
deltat (float) – time step; provide in place of either times or itr to have one computed from the other
dimlist (tuple) – dimensions as will be stored in file (only useful when passing meta data from an existing file to wrmds as keyword args)
11.1.2. mnc¶
- class MITgcmutils.mnc.MNC(fpatt, layout=None, multitime=False)[source]¶
A file object for MNC (tiled NetCDF) data.
Should behave mostly like scipy.io.netcdf.netcdf_file in ‘r’ mode.
- Parameters
fpatt (string) – glob pattern for tile files
layout (string) –
which global layout to use:
- ’model’
use layout implied by Nx, Ny
- ’exch2’
use exch2 global layout
- ’faces’
variables are lists of exch2 faces
default is to use exch2 layout if present, model otherwise
Example
>>> nc = mnc_files('mnc_*/state.0000000000.t*.nc') >>> temp = nc.variables['Temp'][:] >>> salt = nv.variables['S'][:] >>> nc.close() temp and salt are now assembled (global) arrays of shape (Nt, Nr, Ny, Nx) where Nt is the number iterations found in the file (in this case probably 1).
Notes
The multitime option is not implemented, i.e., MNC cannot read files split in time.
- MITgcmutils.mnc.mnc_files(fpatt, layout=None)[source]¶
A file object for MNC (tiled NetCDF) data.
Should behave mostly like scipy.io.netcdf.netcdf_file in ‘r’ mode.
- Parameters
fpatt (string) – glob pattern for tile files
layout (string) –
which global layout to use:
- ’model’
use layout implied by Nx, Ny
- ’exch2’
use exch2 global layout
- ’faces’
variables are lists of exch2 faces
default is to use exch2 layout if present, model otherwise
Example
>>> nc = mnc_files('mnc_*/state.0000000000.t*.nc') >>> temp = nc.variables['Temp'][:] >>> salt = nv.variables['S'][:] >>> nc.close() temp and salt are now assembled (global) arrays of shape (Nt, Nr, Ny, Nx) where Nt is the number iterations found in the file (in this case probably 1).
Notes
The multitime option is not implemented, i.e., MNC cannot read files split in time.
- MITgcmutils.mnc.rdmnc(fpatt, varnames=None, iters=None, slices=Ellipsis, layout=None)[source]¶
Read one or more variables from an mnc file set.
- Parameters
fpatt (string) – glob pattern for netcdf files comprising the set
varnames (list of strings, optional) – list of variables to read (default all)
iters (list of int, optional) – list of iterations (not time) to read
slices (tuple of slice objects) – tuple of slices to read from each variable (typically given as numpy.s_[…])
- Returns
dictionary of variable arrays
- Return type
dict of numpy arrays
Example
>>> S = rdmnc("mnc_*/state.0000000000.*', ['U', 'V'], slices=numpy.s_[..., 10:-10, 10:-10]) >>> u = S['U'] >>> v = S['V']
Notes
Can currently read only one file set (i.e., 1 file per tile), not several files split in time.
Consider using mnc_files for more control (and similar convenience). The same restriction about multiple files applies, however.
11.1.3. diagnostics¶
- MITgcmutils.diagnostics.readstats(fname)[source]¶
locals,totals,itrs = readstats(fname)
Read a diagstats text file into record arrays (or dictionaries).
- Parameters
fname (string) – name of diagstats file to read
- Returns
locals (record array or dict of arrays) – local statistics, shape (len(itrs), Nr, 5)
totals (record array or dict of arrays) – column integrals, shape (len(itrs), 5)
itrs (list of int) – iteration numbers found in the file
Notes
The 5 columns of the resulting arrays are average, std.dev, min, max and total volume.
There is a record (or dictionary key) for each field found in the file.
11.1.4. ptracers¶
11.1.5. jmd95¶
Density of Sea Water using the Jackett and McDougall 1995 (JAOT 12) polynomial
- MITgcmutils.jmd95.dens(s, theta, p)¶
Computes in-situ density of sea water
Density of Sea Water using Jackett and McDougall 1995 (JAOT 12) polynomial (modified UNESCO polynomial).
- Parameters
s (array_like) – salinity [psu (PSS-78)]
theta (array_like) – potential temperature [degree C (IPTS-68)]; same shape as s
p (array_like) – pressure [dbar]; broadcastable to shape of s
- Returns
dens – density [kg/m^3]
- Return type
array
Example
>>> densjmd95(35.5, 3., 3000.) 1041.83267
Notes
AUTHOR: Martin Losch 2002-08-09 (mlosch@mit.edu)
Jackett and McDougall, 1995, JAOT 12(4), pp. 381-388
- MITgcmutils.jmd95.densjmd95(s, theta, p)[source]¶
Computes in-situ density of sea water
Density of Sea Water using Jackett and McDougall 1995 (JAOT 12) polynomial (modified UNESCO polynomial).
- Parameters
s (array_like) – salinity [psu (PSS-78)]
theta (array_like) – potential temperature [degree C (IPTS-68)]; same shape as s
p (array_like) – pressure [dbar]; broadcastable to shape of s
- Returns
dens – density [kg/m^3]
- Return type
array
Example
>>> densjmd95(35.5, 3., 3000.) 1041.83267
Notes
AUTHOR: Martin Losch 2002-08-09 (mlosch@mit.edu)
Jackett and McDougall, 1995, JAOT 12(4), pp. 381-388
11.1.6. mdjwf¶
Density of Sea Water using McDougall et al. 2003 (JAOT 20) polynomial
- MITgcmutils.mdjwf.dens(s, theta, p)¶
Computes in-situ density of sea water
Density of Sea Water using McDougall et al. 2003 (JAOT 20) polynomial (Gibbs Potential).
- Parameters
s (array_like) – salinity [psu (PSS-78)]
theta (array_like) – potential temperature [degree C (IPTS-68)]; same shape as s
p (array_like) – pressure [dbar]; broadcastable to shape of s
- Returns
dens – density [kg/m^3]
- Return type
array
Example
>>> densmdjwf(35., 25., 2000.) 1031.654229
Notes
AUTHOR: Martin Losch 2002-08-09 (Martin.Losch@awi.de)
McDougall et al., 2003, JAOT 20(5), pp. 730-741
- MITgcmutils.mdjwf.densmdjwf(s, theta, p)[source]¶
Computes in-situ density of sea water
Density of Sea Water using McDougall et al. 2003 (JAOT 20) polynomial (Gibbs Potential).
- Parameters
s (array_like) – salinity [psu (PSS-78)]
theta (array_like) – potential temperature [degree C (IPTS-68)]; same shape as s
p (array_like) – pressure [dbar]; broadcastable to shape of s
- Returns
dens – density [kg/m^3]
- Return type
array
Example
>>> densmdjwf(35., 25., 2000.) 1031.654229
Notes
AUTHOR: Martin Losch 2002-08-09 (Martin.Losch@awi.de)
McDougall et al., 2003, JAOT 20(5), pp. 730-741
11.1.7. cs¶
- MITgcmutils.cs.pcol(x, y, data, projection=None, vmin=None, vmax=None, **kwargs)[source]¶
Plots 2D scalar fields on the MITgcm cubed sphere grid with pcolormesh.
- Parameters
x (array_like) – ‘xg’, that is, x coordinate of the points one half grid cell to the left and bottom, that is vorticity points for tracers, etc.
y (array_like) – ‘yg’, that is, y coordinate of same points
data (array_like) – scalar field at tracer points
projection (Basemap instance, optional) – used to transform if present. Unfortunatly, cylindrical and conic maps are limited to the [-180 180] range. projection = ‘sphere’ results in a 3D visualization on the sphere without any specific projection. Good for debugging.
Example
>>> from mpl_toolkits.basemap import Basemap >>> import MITgcmutils as mit >>> import matplotlib.pyplot as plt >>> from sq import sq >>> >>> x=mit.rdmds('XG'); y=mit.rdmds('YG'); e=mit.rdmds('Eta',np.Inf) >>> fig = plt.figure(); >>> mp = Basemap(projection='moll',lon_0 = 0., >>> resolution = 'l', area_thresh = 1000.) >>> plt.clf() >>> h = mit.cs.pcol(x,y,sq(e), projection = mp) >>> mp.fillcontinents(color = 'grey') >>> mp.drawmapboundary() >>> mp.drawmeridians(np.arange(0, 360, 30)) >>> mp.drawparallels(np.arange(-90, 90, 30)) >>> plt.show()
11.1.8. llc¶
- MITgcmutils.llc.contour(*arguments, **kwargs)[source]¶
Create a contour plot of a 2-D llc array (with tricontour).
Call signatures:
contour(X, Y, C, N, **kwargs) contour(X, Y, C, V, **kwargs)
- Parameters
X (array-like) – x coordinates of the grid points
Y (array-like) – y coordinates of the grid points
C (array-like) – array of color values.
N (int) – number of levels
V (list of float) – list of levels
kwargs – passed to tricontour.
- MITgcmutils.llc.contourf(*arguments, **kwargs)[source]¶
Create a contourf plot of a 2-D llc array (with tricontour).
Call signatures:
contourf(X, Y, C, N, **kwargs) contourf(X, Y, C, V, **kwargs)
- Parameters
X (array-like) – x coordinates of the grid points
Y (array-like) – y coordinates of the grid points
C (array-like) – array of color values.
N (int) – number of levels
V (list of float) – list of levels
kwargs – passed to tricontour.
- MITgcmutils.llc.div(u, v, dxg=None, dyg=None, rac=None, hfw=None, hfs=None)[source]¶
Compute divergence of vector field (U,V) on llc grid
Call signatures:
divergence = div(U, V, DXG, DYG, RAC, HFW, HFS) divergence = div(U, V) divergence = div(U, V, DXG, DYG) divergence = div(U, V, DXG, DYG, RAC) divergence = div(U, V, DXG, DYG, hfw=HFW, hfs=HFS)
- Parameters
u (array-like (timelevel,depthlevel,jpoint,ipoint)) – x-component of vector field at u-point
v (array-like (timelevel,depthlevel,jpoint,ipoint)) – y-component of vector field at v-point
dxg (array-like (jpoint,ipoint), optional) – grid spacing in x across v-point, defaults to one
dyg (array-like (jpoint,ipoint), optional) – grid spacing in y across u-point, defaults to one
rac (array-like (jpoint,ipoint), optional) – grid cell area, defaults to dxg*dyg
hfw (array-like (depthlevel,jpoint,ipoint), optional) – hFac at u-point, defaults to one
hfs (array-like (depthlevel,jpoint,ipoint), optional) – hFac at v-point, defaults to one
- MITgcmutils.llc.faces2mds(ff)[source]¶
convert 6 faces to mds 2D data, inverse opertation of llc.faces
- MITgcmutils.llc.flat(fld, **kwargs)[source]¶
convert mds data into global 2D field only fields with 2 to 5 dimensions are allowed
- MITgcmutils.llc.grad(X, dxc=None, dyc=None, hfw=None, hfs=None)[source]¶
Compute horizontal gradient of scalar field X on llc grid
Call signatures:
dXdx, dXdy = div(X, DXC, DYC, HFW, HFS) dXdx, dXdy = div(X) dXdx, dXdy = div(X, DXC, DYC) dXdx, dXdy = div(X, hfw=HFW, hfs=HFS)
- Parameters
X (array-like (timelevel,depthlevel,jpoint,ipoint)) – scalar field at c-point
dxc (array-like (jpoint,ipoint), optional) – grid spacing in x across u-point, defaults to one
dyc (array-like (jpoint,ipoint), optional) – grid spacing in y across v-point, defaults to one
hfw (array-like (depthlevel,jpoint,ipoint), optional) – hFac at u-point, defaults to one
hfs (array-like (depthlevel,jpoint,ipoint), optional) – hFac at v-point, defaults to one
- MITgcmutils.llc.mds(fld, center='Atlantic')[source]¶
convert global ‘flat’ field into mds data; only fields with 2 to 5 dimensions are allowed
- MITgcmutils.llc.pcol(*arguments, **kwargs)[source]¶
Create a pseudo-color plot of a 2-D llc array (with plt.pcolormesh).
Call signatures:
pcol(X, Y, C, **kwargs) pcol(X, Y, C, m, **kwargs)
- Parameters
X (array-like) – x coordinates of the grid point corners (G-points)
Y (array-like) – y coordinates of the grid point corners (G-points)
C (array-like) – array of color values.
m (Basemap instance, optional) – map projection to use. NOTE: currently not all projections work
kwargs – passed to plt.pcolormesh.
- MITgcmutils.llc.uv2c(u, v)[source]¶
Average vector component (u,v) to center points on llc grid
Call signatures:
uc,vc = uv2c(U,V)
- Parameters
U (array-like (timelevel,depthlevel,jpoint,ipoint)) – x-component of vector field at u-point
V (array-like (timelevel,depthlevel,jpoint,ipoint)) – y-component of vector field at v-point
11.1.9. gluemncbig¶
This command line script is part of MITgcmutils and provides a convenient method for stitching together NetCDF files into a single file covering the model domain. Be careful though - the resulting files can get very large.
Usage: gluemncbig [-2] [-q] [--verbose] [--help] [--many] [-v <vars>] -o <outfile> <files>
-v <vars> comma-separated list of variable names or glob patterns
-2 write a NetCDF version 2 (64-Bit Offset) file allowing for large records
--many many tiles: assemble only along x in memory; less efficient
on some filesystems, but opens fewer files simultaneously and
uses less memory
-q suppress progress messages
--verbose report variables
--help show this help text
All files must have the same variables.
Each variable (or 1 record of it) must fit in memory.
With --many, only a row of tiles along x must fit in memory.
Examples:
gluemncbig -o ptr.nc mnc_*/ptr_tave.*.nc
gluemncbig -o BIO.nc -v 'BIO_*' mnc_*/ptr_tave.*.nc