MSCPROGS post-processing tools

MSCPROGS is a collection of Fortran (and C) programs for post-processing and analysing HYCOM model output. It is a general-purpose toolkit covering grid projection, field extraction, temporal averaging, section transports, Lagrangian ice drift, grid utilities, and more. The tools live in hycom/MSCPROGS/.

All of these post-processing and analysis steps can alternatively be done with xhycom, a Python-based toolkit for working with HYCOM output.

For compilation instructions, see Compile MSCPROGS.

Important

There are two ways to post-process and analyse HYCOM output — pick one:

  • MSCPROGS (this page) — the original Fortran/C toolkit, described below.

  • xhycom — a Python-based toolkit covering the same post-processing and analysis tasks.

You don’t need both. See xhycom if you’d rather work in Python.

Running the tools

Environment setup

Before running any MSCPROGS tool, load the HPC modules and set $MSCPROGS. See the HPC environment section for details.

Source HPC environment — Betzy (NRIS/Sigma2)
source ${HOME}/NERSC-HYCOM-CICE/environment/betzy_env.sh

source $WORK/<CONFIGNAME>/REGION.src    # sets $MSCPROGS and other paths

REGION.src also adds $MSCPROGS/bin and $MSCPROGS/bin_setup to your PATH, so all MSCPROGS executables are immediately available.

Working directory

The natural place to run MSCPROGS tools is the experiment data directory:

cd $WORK/<CONFIGNAME>/expt_<EXPT_ID>/data/

expt_postprocess.sh copies regional.grid.* and regional.depth.* here automatically after each model run, so the grid files are already in place.

You still need to copy the relevant input files into the working directory. Example input files are provided in $MSCPROGS/Input/.

Tip

You can also run MSCPROGS tools directly from the scratch directory (expt_<EXPT_ID>/SCRATCH/) to analyse output before the run is finalised. expt_preprocess.sh already copies regional.grid.* and regional.depth.* there, so only the input files need to be copied in — the same step as above.

Typical workflow

The general pattern for any MSCPROGS tool is:

  1. Navigate to the data directory and copy the needed input files (above).

  2. Run the tool — most print usage information when called with no arguments:

    hyc2proj           # prints usage
m2nc               # prints usage
m2section          # prints usage

  3. Process the output NetCDF files with your analysis scripts.

hyc2proj and hyc2stations — projection and interpolation

hyc2proj projects HYCOM isopycnal-layer output onto a regular grid and interpolates vertically to fixed depth levels, writing CF-compliant NetCDF. Four output projections are supported:

Projection

proj.in keyword

Notes

Regular lon/lat

regular

Uniform degree spacing; grid defined by start/end/spacing in degrees

Mercator

mercator

Uniform spacing in projected coordinates (non-uniform in degrees); defined by start/end bounds + number of points; conformal

Polar stereographic

polar_stereographic

Defined by projection coordinates + number of points; optionally rotate vectors onto the PS grid

Model native

native

No horizontal interpolation; extract sub-domain by index range

hyc2proj requires the following files in the working directory:

  • regional.grid.a/.b and regional.depth.a/.b — already present in data/ and SCRATCH/

  • grid.info — conformal mapping parameters for the model grid; copy from topo/:

    cp $WORK/<CONFIGNAME>/topo/grid.info .

  • proj.in — target projection and grid (see Input files)

  • depthlevels.in — vertical depth levels to interpolate to

  • an extract.* file — fields to extract; the tool auto-selects the file by name based on the input file type (e.g. extract.archm for archm files), so the name must not be changed

Horizontal interpolation is bilinear; vertical interpolation defaults to cubic spline. staircase and linear options are also available (faster but lower quality).

cp $WORK/<CONFIGNAME>/topo/grid.info .
cp $MSCPROGS/Input/proj.in.regular_grid proj.in     # choose and edit a sample proj.in
cp $MSCPROGS/Input/depthlevels.in .
cp $MSCPROGS/Input/extract.archm .
hyc2proj archm.*.a                                   # process one or more .a files

One NetCDF output file is created per input file. The name is derived from the input filename with the extension changed to .nc and the date embedded — for example archm_20100101_00.nc for an archm file. The original .a/.b files are not modified.

Note

hyc2proj can be slow on a login node for non-native projections. As a reference, a single TP2 archm file with 4 three-dimensional fields and 40 target depth levels takes roughly 2 minutes on a regular 0.125° grid covering the Arctic (−38°–20°E, 65°–88°N). For many files, consider running in an interactive job or a short batch script.

hyc2stations is a companion routine that interpolates HYCOM fields onto specified spatial positions (station locations or cross-section points) rather than a grid. It reads stations.in instead of proj.in, together with depthlevels.in and an extract.* file. Output is one NetCDF file per station group, named <inputbase>_group<groupname>.nc. A list of all generated files is written to hyc2stations.filelist. A helper script setupstations.py generates stations.in entries from start/end coordinates along a rhumb line.

Batch post-processing

For a single run or a few files, running hyc2proj interactively (above) is fine. For many files — e.g. a multi-year reanalysis — submit it as a SLURM job. A convenient pattern is to parallelise by year: stage each year’s archives into its own Y<year>/ directory together with a copy of the required input files, then run one hyc2proj task per year.

The required files (see hyc2proj) must be present in each per-year directory: regional.grid.*, regional.depth.*, grid.info, proj.in, depthlevels.in, the matching extract.* file, and the compiled hyc2proj executable.

Note

The parallel-by-year wrapper does not work with detvflux — leave it off in extract.archm/extract.archv.

Example: parallel-by-year SLURM script (Olivia)

The header (#SBATCH options) and module block are machine-specific. On Olivia, use a CPU partition and load the modules established when compiling MSCPROGS.

#!/bin/bash
#SBATCH --account=nn2993k
#SBATCH --job-name=hyc2proj
#SBATCH --nodes=1
#SBATCH --ntasks=32            # one task per year
#SBATCH --mem-per-cpu=10G
#SBATCH --time=05:00:00        # ~3 h per model year, as a guide

# --- Olivia HPC modules ---
source $HOME/NERSC-HYCOM-CICE/environment/olivia_env.sh
ulimit -s 2000000

yr_s=1993
yr_e=2024
rid=archm.

# Stage each year's archives + input files into Y<year>/
for year in $(seq $yr_s $yr_e); do
    dir="Y$year"; mkdir -p "$dir"
    cp ./regional.* ./grid.info ./extract.archm ./depthlevels.in ./proj.in ./hyc2proj "$dir"/
    mv ${rid}${year}_*.a ${rid}${year}_*.b "$dir"/ 2>/dev/null
done

# Run one hyc2proj per year in parallel
for year in $(seq $yr_s $yr_e); do
    ( cd "Y$year" && srun --exclusive -n1 -c1 ./hyc2proj --vertint linear ${rid}${year}_*.a ) &
done
wait

m2nc / m2t — 2D field extraction to NetCDF

Reads HYCOM .ab files and writes selected 2D fields to tmp1.nc. Fields to extract are specified in the relevant extract.* file. Multiple input files produce multiple time records in the output. m2t additionally writes a tmp1.dat file in Tecplot ASCII format (a commercial CFD visualization tool; the format is also readable by ParaView and similar tools).

cp $MSCPROGS/Input/extract.archm .
m2nc archm.*.a

Note

No vertical interpolation — output stays on isopycnal layers.

h2nc (ExtractNC3D) — 3D field extraction to NetCDF

h2nc <file(s)>

Converts HYCOM .ab files to NetCDF with 3D fields. Like m2nc it uses the extract.* files to select fields, but retains the full vertical structure. Output is written to tmp1.nc.

hycave — temporal averaging

Computes a thickness-weighted time mean across the supplied files. filetype is one of restart, nersc_daily, nersc_weekly, archv, or archm. Requires regional.grid and regional.depth in the working directory. Output is a new .ab file pair with a fixed placeholder date — for example AVEDAILY_9999_999_9999_999.ab for nersc_daily and AVE.archv.9999_999_99.ab for archv. Input files are not modified.

hycave nersc_daily TP2daily*.a

m2section, m2transport, m2transport2 — section transports

Extracts data along user-defined sections and computes transports. Sections are defined in sections.in (see Input files). All three scripts call section_intersect first to identify model grid points along each section using a great-circle algorithm.

cp $MSCPROGS/Input/sections.in .
cp $MSCPROGS/Input/extract.archm .              # auto-selected by file type
m2section archm.*.a                            # output: section001.nc, section002.nc, ...
m2transport archm.*.a                          # output: transport_net001.dat, ...

  • m2section — extracts fields along sections to section001.nc, section002.nc, … Uses extract.* files to select fields.

  • m2transport — fast barotropic volume and ice transport across sections. Output to ASCII files transport_net001.dat, transport_net002.dat, …

  • m2transport2 — full-depth transport with finer control. Transport criteria (depth range, water mass temperature range, etc.) are set in transport.in. Also writes transports2.nc.

m2transport2 additionally supports scalar (e.g. heat) transports when scalartransport.in is present. Each line specifies:

'<variable>'  <offset>  <scale_factor>

For example, heat transport relative to 6 °C:

'temp'  6.  3987000

where 3987000 ≈ cₚ × ρ.

Note

Positive transport direction: stand at the first section point and look toward the second — positive is to the right. Sections follow great circles, not rhumb lines.

icedrift2 — Lagrangian ice drift

icedrift2            # reads Input/icedrift.in
icedrift2 <args>     # CERSAT mode: takes CERSAT ice drift positions as input

Integrates ice particle trajectories forward in time using Runge–Kutta (2nd order), reading drift velocities from HYCOM daily output. Useful for validation against IABP buoy data or as input for data assimilation. Output is written to drift.uf and drift_diag.uf (in CERSAT mode the names follow the input drift file, e.g. <driftfile>.uf and <driftfile>_diag.uf). A text summary is written to tsdrift.asc.

nestbat — nesting bathymetry

Adjusts a local-domain depth file so bathymetry at the nesting boundary is consistent with the parent (global) grid, smoothly transitioning across a user-specified zone. The nestbat-2.2 version is used by the nestbat.sh setup script. Run from an empty working directory containing the required global and local grid files; the program is interactive.

Tides_FES2014 — tidal forcing from FES atlas

Generates tidal boundary conditions for HYCOM from the FES2014 tidal atlas. Requires the external FES2014 C library and the GNU C compiler; see src/Tides_FES2014/README.md for build instructions.

Other tools

Tool

Executable

Purpose

Barstrf

hyc2bastrf

Barotropic streamfunction

Conf_grid

—

Generate conformal-mapping model grids

ConfmapRoutines

—

Convert between grid indices and lon/lat

Curviint

curviint

Interpolate restart files from one grid to another

DateTools

—

Convert between ordinal days and calendar dates

DProfile

dprofile

Extract vertical profiles from model files

Ensstat

—

Ensemble and time-domain statistics

FindLayer

findlayer

Layer thickness between variable thresholds (e.g. salinity 35–36 psu)

GP / GPdens

—

Extract and analyse grid-point time series (activated via infile.in)

GridToLL

gridtoll

Interpolate fields from tmp1.nc to a regular lon/lat grid

MkEnsemble

mkensemble

Create ensemble by perturbing restart files

Relax

rmunew

Generate relaxation masks

RelaxToNetCDF

—

Interpolate climatology files to projection/depth levels

River_Forcing

rivers

Spread point-source river discharge along the model coast

SSHFromState

ssh_from_restart

Compute SSH from restart files

Tides_CSR

csr2mod_GE

Tidal boundary forcing from CSR tidal atlas

TRIP

trip_*

River forcing from the TRIP database + ERA40/ERA-i runoff

ZONAL

zonal, mosf

Zonal averages and meridional overturning streamfunction

Post-processing runs from NIRD

NIRD access differs between machines:

Login node

SVC node

Compute node

Betzy

read-write

—

not mounted

Olivia

not mounted

read-write

read-only

The recommended workflow is to copy (not symlink) archives onto the scratch filesystem before processing — the parallel-by-year hyc2proj script moves files into per-year directories, which fails on a read-only source, and streaming large 3-D archives live over NIRD is slow.

Stage the archive files (archm.* / archv.*, both .a and .b) together with the run’s regional.* and grid.info into a scratch directory:

  • Betzy: copy from the login node to /cluster/work/users/$USER/postproc/<run>/

  • Olivia: copy from a service (SVC) node to /cluster/work/projects/nn2993k/$USER/postproc/<run>/

Then run the post-processing tool from the staging directory. On Olivia, tools that only read archives — m2nc, hycave, m2section — can also read directly from NIRD inside a batch job, skipping the staging step.

Move the resulting .nc files to /cluster/projects/nn2993k/$USER/... (or back to NIRD) before the 21-day scratch purge, then delete the staged archives.

Tip

For a large multi-year run, stage and process one year at a time to keep the scratch footprint bounded.

Input files

Most tools share a common set of input files, stored as examples in $MSCPROGS/Input/. Copy the relevant file to the working directory and edit as needed.

extract.* files

Controls which fields are extracted. The tool auto-detects the file type from the input filename and opens the corresponding extract file by name — so the file must be present in the working directory under exactly that name (do not rename it). There is a separate template for each file type:

Template

File type

Notes

extract.archm

archm — time-averaged archive

Most commonly used; averaging period set by meanfq in blkdat.input

extract.archv

archv — instantaneous snapshot archive

Written at diagfq frequency

extract.archv_wav

archv_wav — wave archive

Input/ also contains extract.daily, extract.weekly, extract.restart, and extract.sec — these are for older TOPAZ-format output files and are not used in standard NERSC HYCOM-CICE runs.

The extract.sec file in Input/ is an additional example for section tools; like all other tools, m2section auto-selects extract.archm or extract.archv based on file type — there is no special extract.sec lookup.

4                      # version number — do not change
pres     intf          # field used to locate layer interfaces — do not change
2 F F T                # flags: (unused) sphere rotate index-velocities — do not change
1 270 1 210            # section-plot index range (ia ib ja jb) — do not change
1 1   XX XX   X        # 3-D layer range header — do not change
saln    01 01   T      # 3-D field: extract all layers
temp    01 01   T
pres    01 01   T
utot    01 01   T
vtot    01 01   T
ssh     00 00   T      # 2-D field (layer indices are 00 00 for 2-D fields)
fice    00 00   T
hice    00 00   F      # F = skip this field

Each field line has four columns:

Column

Description

fieldname

Must match the field name in the .b file exactly. Run cat *.b to see what is available in your output.

layer_start

First isopycnal layer to extract. Use 00 for 2-D fields.

layer_end

Last isopycnal layer to extract. Use 00 for 2-D fields. layer_start layer_end is a literal range: 01 01 extracts only layer 1, 01 26 extracts layers 1–26.

T / F

T to extract this field, F to skip it.

The layer range is respected by m2nc and h2nc. For hyc2proj the range is ignored — it always reads all isopycnal layers and then interpolates them onto the z-levels in depthlevels.in, so 01 01 is conventional placeholder.

Note

Field names differ between file types — saln / temp in daily output, salin / thknss in archv. Always copy the extract file that matches your input file type.

proj.in

Defines the target projection for hyc2proj. Example files for all four projections are in Input/.

Regular lon/lat (proj.in.regular_grid):

regular    # projection type
-38.0      # longitude of SW corner (degrees E)
 20.0      # longitude of NE corner (degrees E)
0.125      # longitude grid spacing (degrees)
 65.0      # latitude of SW corner (degrees N)
 88.0      # latitude of NE corner (degrees N)
0.125      # latitude grid spacing (degrees)

Mercator (proj.in.mercator):

mercator
-100.0     # longitude W edge (degrees E)
  20.0     # longitude E edge (degrees E)
   800     # number of grid points in E-W direction
 -20.0     # latitude S edge (degrees N)
  70.0     # latitude N edge (degrees N)
   400     # number of grid points in S-N direction
    0.     # central longitude

Polar stereographic (proj.in.stereographic):

polar_stereographic
-3800.     # left  edge of projection area (km from pole)
 3800.     # right edge of projection area (km)
 12.5      # grid spacing in x-direction (km)
-5500.     # bottom edge of projection area (km)
 5500.     # top    edge of projection area (km)
 12.5      # grid spacing in y-direction (km)
-45.0      # central longitude of projection
 90.0      # central latitude  of projection
 0.01      # conversion factor (km → m if 0.001, km → km if 1.0)
 F         # rotate vectors to output grid (F = keep East/North components)

Model native (proj.in.native):

native
    1      # first i-index
  100      # last  i-index
    1      # i-index stride
    1      # first j-index
  100      # last  j-index
    1      # j-index stride

Extracts a sub-domain without horizontal interpolation.

depthlevels.in

Specifies the target vertical depth levels (in metres) for hyc2proj and hyc2stations. The first line is the number of levels; one depth per line follows.

40        # number of depth levels
0.0
2.0
4.0
...
4000.0

sections.in

Defines named sections as lon/lat endpoint pairs for m2section, m2transport, and m2transport2.

1.2 F          # file version; F = quiet output (T = verbose)
#####
FRAM_STRAIT    # section name (used in output filenames)
-20. 79.       # start point: lon, lat (degrees)
15.  79.       # end  point:  lon, lat (degrees)
#
DENMARK_STRAIT
-37.  66.1
-22.68 66.

Sections follow great circles between the two endpoints. The transport sign convention is: positive to the right when standing at the start point and facing the end point.

transport.in

Used by m2transport2 to filter transport by depth range, salinity range, or temperature range. Each line defines one transport component:

SVINOY-ATL-N   SVINOY   SALINITY   n   35.00   39.00
SVINOY         SVINOY   DEPTH      n    0.00 1000.00

Column

Description

Label

Output name used in result files

Section

Section name — must match an entry in sections.in

Criterion

DEPTH, SALINITY, or TEMPERATURE

Direction

n (no direction filter)

Min

Lower bound of the criterion range

Max

Upper bound of the criterion range

stations.in

Defines groups of station positions for hyc2stations. Each group starts with #GROUPNAME and is followed by lon/lat pairs, one per line.

#FRAM_STRAIT
-20.00  79.00
-10.00  79.25
  0.00  79.50
 15.00  79.00
#DENMARK_STRAIT
-37.00  66.10
...

The helper script setupstations.py generates these point lists by interpolating at regular spacing along a rhumb line between two endpoints.