topobath
topobath adds, removes, or replaces topography and/or
bathymetry in a three-dimensional multiparameter NetCDF cube. It can shift an
entire model column, construct a one-dimensional or linear land extension,
remove a water column, or rebuild water columns above a newly created seafloor.
Synopsis
gmt topobath model.nc [new_surface.nc] -Goutput.nc -Ooperation [-Mmethod] [options]
Coordinates
Source scaling is applied before classification, sampling, or any surface
operation. The transformed model z axis must increase strictly from the top
of the model to its bottom. Thus, elevations above sea level are negative and
depths below sea level are positive. Old and new surface grids must use the
same units and signs. Use the model +z modifier and the surface +v or
+z modifier to establish this convention when necessary. Scaling changes
coordinate values without reordering model layers.
-T<zmin>/<zmax>/<dz> defines the working output axis, where zmin is
the top, zmax is the bottom, and zmin < zmax. Without -T, the
module derives increasing bounds at the smallest transformed input spacing.
-Z is applied after every topography operation. Its +x, +y, and
+z modifiers scale final coordinates, while +X, +Y, and +Z
set their unit metadata. The +v and +V lists independently scale and
name the units of selected output fields in model-selector order. Scaling is
performed in place without reordering coordinates, model layers, fields, or
ancillary variables. A negative axis scale therefore writes a decreasing
coordinate. For example,
-Z+x0.001+Xkm+y0.001+Ykm+z-1000+Zm+v0.001,0.001+Vkm/s,km/s
can convert horizontal coordinates from metres to kilometres, restore a
positive-up vertical coordinate in metres, and convert two velocity fields
from metres per second to kilometres per second.
Operations and scopes
Every run explicitly selects an operation with -O:
-Oaadds a requested surface to a model that normally begins at sea level. Nonzero old surfaces are accepted with a warning.-Orremoves the selected existing surface. Land is flattened to zero. Targeted seafloors move to zero and their water columns are removed.-Oxreplaces the selected old surface with the supplied new surface.
Append +t to operate only on dry-land topography or +b to operate
only on wet-region bathymetry. With no modifier, both classes are selected. A
column outside the selected scope is unchanged. The wet/dry classification,
rather than the sign of a surface value, controls this selection, so dry land
below sea level is valid (as is the case in Death Valley or the Dead Sea).
If interpolated relief lies above sea level in a cell classified as wet, the
requested seafloor is clamped to sea level with a warning. This accommodates
small coastline differences between a categorical mask and a continuous relief
grid without changing the authoritative class.
Adding and replacing require one construction method:
-Mpshifts the complete column to the new surface.-Memaps the old dry-land surface to zero and extends its value vertically to the new elevation.-Mlmaps the old dry-land surface to zero and linearly grades its value to the field-specific surface minimum supplied with-L.
Wet columns always shift the solid model to the requested seafloor and fill
the newly created water column with values set with -W. Dry columns whose
requested surface is below sea level also shift as complete columns. -M is
therefore required by add(a) and replace(x) operations, but is not used by
remove(r).
Wet/dry classification and inference
-Cg gives GMT shoreline classification priority and is the default.
-Cm lets model evidence override that prior where the model resolves a
class. -Cl and -Cw classify the complete domain as land or wet,
respectively. For geographic models, -D selects the GSHHG resolution and
-A controls which hierarchy levels are wet or dry. Selected oceans, lakes,
and ponds are wet; land and islands in lakes are dry. A user -K grid with
wet = 0 and land = 1 is authoritative and replaces both GMT and model
classification.
Supplying the old surface with -E is strongly recommended. Without it,
topobath finds the shallowest valid boundary among all selected fields.
Different field boundaries produce a warning, and the shallowest boundary is
used. In wet columns, field-specific -W values and tolerances identify the
water column and seafloor. Model classification can only recognize wet
regions when those signatures are available.
A parameter missing exactly at a surface identified by another field remains
missing. Values are not pulled upward merely to hide that disagreement.
Consequently, -E is particularly important when a legitimate surface
parameter may be NaN, because a surface inferred from that parameter alone
could mistake the missing value for air. -Q writes the inferredold surface
and, with +c, classification codes 0 (unresolved), 1 (land), and 2 (wet).
Air and water
-F<field>/<air> sets a field-specific air value. Unspecified air values
are NaN. The air value applies strictly above the identified free surface; a
NaN exactly on that surface remains NaN.
-W<field>/<water>[+t<tolerance>] serves two roles: it recognizes old
water during surface inference and supplies values for newly created water.
It is required for every selected field when a wet surface must be inferred,
or when add/replace creates or rebuilds a water column. It is not required
merely because wet cells exist when -E supplies the old surface and the
operation either removes their water or leaves the wet columns unchanged.
Sampling and missing values
-H optionally fills strictly internal holes in each native x-y model
layer. Nearest-neighbor, linear Delaunay, local weighted-average, spline, and
minimum-curvature methods are available. Boundary-connected missing regions
and original valid nodes remain unchanged, and +m<maxgap> can restrict
filling by the x and y spans of each hole. Without -H, native horizontal
holes remain missing.
-R must lie within the transformed model domain. -I and common
-n resample the selected fields, inferred surface, and relief onto the
horizontal output lattice. This resampling is distinct from -H and does
not deliberately fill data holes. -T defines the output vertical
lattice. -S selects the GMT vertical interpolant: Akima (a), cubic
(c), step-up (e), linear (l), nearest (n), or smoothing spline
(s<p>) with non-negative fit parameter p. Linear is the default. Its
+g modifier bridges internal missing layers, with an optional maximum gap.
Exterior missing regions and air are not extrapolated. Packed values and
declared missing sentinels are decoded on input, and output missing values are
set to NaN.
The coordinate modifiers of -Z also transform the optional -Q grids:
+x and +y scale their coordinates, and +z scales the inferred old
surface. Classification codes cannot be scaled. Field +v and +V
modifiers apply only to the selected three-dimensional model fields.
The complete order is input scaling, optional native-layer -H filling,
-n horizontal resampling, wet/dry classification, old-surface selection
or inference, the requested operation, -S vertical interpolation, and
final -Z output scaling. Note that a sufficiently large pull-up operation
(i.e., upward column shift) can move the model base outside the source cube.
Such newly exposed base values are set to NaN, and the module reports a warning.
Field selectors limit the processed 3-D variables. Without one, every numeric
3-D variable sharing the first compatible x/y/z dimensions is selected.
Compatible global, coordinate, scalar, and horizontal ancillary metadata are
unchanged. Coordinate-dependent ancillary variables are omitted when -R or
-I changes the horizontal lattice.
Usage
gmt topobath [gq] 6.5.0 [64-bit] - Add, remove, or replace topography and/or bathymetry in three-
dimensional multiparameter NetCDF cubes
usage: gmt topobath <model.nc>[?fields][+n<missing>][+x<sx>][+X<unit>][+y<sy>][+Y<unit>][+z<sz>]⏎
…[+Z<unit>][+v<scales>][+V<units>] [<new_surface>[?field][+x<sx>][+X<unit>][+y<sy>][+Y<unit>][+v|⏎
…z<scale>][+V|Z<unit>]] -G<output.nc> -Oa|r|x[+t|b] [-A<min_area>[/<min_level>/<max_level>]]
[-Cg|m|l|w] [-D<a|f|h|i|l|c|n>] [-E<old_surface>] [-F<field>/<air>] [-H[n|l|a|s|m[<arg>]]⏎
…[+m<maxgap>]] [-K<landmask>[?field][+x<sx>][+X<unit>][+y<sy>][+Y<unit>]] [-L<field>/<minimum>]
[-Mp|e|l] [-Q<surface.nc>[+c<classification.nc>]] [-Sa|c|e|l|n|s<p>[+g[<maxgap>]]] [-T<zmin>/⏎
…<zmax>/<dz>] [-W<field>/<value>[+t<tolerance>]] [-Z[+x<sx>][+X<xunit>][+y<sy>][+Y<yunit>]⏎
…[+z<sz>][+Z<zunit>][+v<scales>][+V<units>]] [-R<west>/<east>/<south>/<north>[+r]] [-I<dx>[/<dy>]]
[-V[q|e|w|t|i|c|d]] [-di<nodata>[+c<col>]] [-n[b|c|l|n][+a][+b<BC>][+c][+t<threshold>]]
REQUIRED ARGUMENTS:
<model.nc>[?field1,field2,...][+n<missing>][+x<sx>][+X<unit>][+y<sy>][+Y<unit>][+z<sz>]⏎
…[+Z<unit>][+v<scales>][+V<units>]
Read selected numeric 3-D variables sharing x, y, and z. Without a selector, process every
compatible 3-D variable. Lowercase coordinate and value modifiers scale the unpacked input.
Uppercase modifiers set output unit metadata. One +v/+V entry is broadcast. Otherwise entries
follow the selector field order. Use +z to standardize the model vertical coordinate. After
scaling, z must increase from negative elevation above sea level to positive depth below sea
level. Scaling changes coordinate values but does not reorder model layers. For example, +z-
0.001+Zkm converts a stored z axis of 4000 ... -16000 m to -4 ... 16 km.
<new_surface>[?field][+x<sx>][+X<unit>][+y<sy>][+Y<unit>][+v|z<scale>][+V|Z<unit>]
Supply the requested topographic and/or bathymetric surface for -Oa and -Ox. Its transformed
values must use the model's working units and convention: negative topographic elevation and
positive bathymetric depth. It is sampled after its coordinate and value scales are applied.
For a surface grid, +z/+Z are aliases for +v/+V. +s<scale> is also accepted for compatibility
with older commands. New commands should use +v<scale> or +z<scale>. For example, +z-
0.001+Zkm converts positive-up relief in metres to the required signed values in kilometres.
Omit this grid with -Or.
-G<output.nc>
Write selected variables to a NetCDF model ordered (z,y,x).
-Oa|r|x[+t|b]
Select add (a), remove (r), or replace (x). Add maps a nominally sea-level surface to the
requested surface and warns when selected columns already have nonzero relief. Remove shifts
the selected free surface or seafloor to zero. Removing bathymetry removes the water column.
Replace removes the old selected surface and applies the requested one. Append +t to operate
on dry-land topography only or +b to operate on wet-region bathymetry only. Otherwise operate
on both. The Wet/dry classification rather than the sign of a surface value, determines which
columns each scope selects.
OPTIONAL ARGUMENTS:
-R<west>/<east>/<south>/<north>[+r]
Specify the min/max coordinates of your data region in user units. Use dd:mm[:ss] for regions
given in arc degrees, minutes [and seconds]. Use -R<xmin>/<xmax>/<ymin>/<ymax>[+u<unit>] for
regions given in projected coordinates, with <unit> selected from e|f|k|M|n|u [Default: e]. If
+u is set, projected regions centered on (0,0) may be set via -R<halfwidth>[/⏎
…<halfheight>]+u<unit>, where <halfheight> defaults to <halfwidth> if not given. Use
[yyyy[-mm[-dd]]]T[hh[:mm[:ss[.xxx]]]] format for time axes. Append +r if -R specifies the
coordinates of the lower left and upper right corners of a rectangular area.
Alternatively, use -R<code><x0>/<y0>/<n_columns>/<n_rows> for origin and grid dimensions,
where <code> is a 2-char combo from [T|M|B][L|C|R] (top/middle/bottom/left/center/right) and
grid spacing must be specified via -I<dx>[/<dy>] (also see -r).
Set an output horizontal region contained within the transformed model domain. The default
uses the complete model region.
-I<dx>[/<dy>]
Set positive output horizontal increments. By default, retain the model increments. Selected
model fields, the inferred surface, and relief are sampled onto this lattice. Additional
variables that depend on x or y are omitted when -R or -I changes the horizontal grid because
topobath does not resample them. Single-value metadata, such as map-projection information, is
retained.
-A<min_area>[/<min_level>/<max_level>][+a<antarctica>][+l|r][+p<percent>]
Select GSHHG features used by automatic shoreline classification. Features below <min_area>
km^2 and hierarchy levels outside <min_level>/<max_level> are skipped. The default is 0/0/1,
which uses ocean and land only. Levels are 0 ocean, 1 land, 2 lake, 3 island in lake, and 4
pond. GMT's +a, +l, +r, and +p modifiers follow grdlandmask.
-Cg|m|l|w
Choose classification precedence. Use g for GMT shoreline priority (the default), m to let
model evidence override GMT where model evidence resolves a class, l to classify the entire
domain as land, or w to classify it as wet. Model wet evidence requires matching -W
signatures. A user -K mask is authoritative and cannot be combined with -Cm, -Cl, or -Cw.
-D<a|f|h|i|l|c|n>
Set the GMT shoreline resolution used by -Cg or as the prior for -Cm. The default is low
resolution (l). Use n to disable shoreline classification. This option has no effect for
Cartesian models, -K, -Cl, or -Cw.
-E<old_surface>[?field][+x<sx>][+X<unit>][+y<sy>][+Y<unit>][+v|z<scale>][+V|Z<unit>]
Supply the old free surface or seafloor in the transformed model's signed vertical convention.
Use +v or +z to scale its values. Without -E, infer the shallowest boundary identified across
the selected fields. -W signatures locate wet-region seafloors. If fields disagree, use the
shallowest boundary and warn. Supplying -E is highly recommended when a valid parameter can be
missing at the free surface, because an unrecognized NaN may otherwise be treated as air.
+s<scale> is also accepted for compatibility with older commands. New commands should use
+v<scale> or +z<scale>.
-F<field>/<air>
Set the value written above the free surface for one selected field. Repeat as needed. The
default air value is NaN for every field. A field that is NaN exactly at a known surface
remains NaN. the air value applies strictly above that surface.
-H[n|l|a|s|m[<arg>]][+m<maxgap>]
Fill strictly internal horizontal missing-data holes in every native x-y model layer before
horizontal resampling, surface inference, and the requested topography operation. Original
non-missing nodes and boundary-connected missing regions are preserved. Without -H, native
horizontal holes are not filled. Use linear Delaunay interpolation when -H is given without a
method. Available methods are:
Nearest neighbor (n). Optionally append a search radius in grid nodes.
Linear Delaunay interpolation (l). This is the default.
Local weighted average (a). Optionally append radius[/sectors] in grid nodes. The default is
3/4.
Spline interpolation (s). Optionally append tension in the range 0-1. The default is 0.
Minimum-curvature interpolation (m). Optionally append tension in the range 0-1. The default
is 0.
+m Only fill holes whose x and y spans are both no larger than <maxgap> grid nodes. The
default is to fill all internal holes.
-K<landmask>[?field][+x<sx>][+X<unit>][+y<sy>][+Y<unit>]
Supply an authoritative classification grid with wet=0 and land=1. It replaces GMT and model
classification rather than serving only as a prior. Use ?field to select the mask variable in
a multi-variable file. The +x and +y modifiers scale its coordinates, while +X and +Y set
their units. Mask values cannot be scaled.
-Mp|e|l
Choose how -Oa or -Ox constructs selected dry-land topography. Pull-up/push-down (p) shifts
the whole column. Constant 1-D extension (e) first maps the old surface to zero, then extends
its surface value to the requested elevation. Linear extension (l) instead grades from that
surface value to -L at the requested elevation. Wet columns and dry columns whose requested
surface is below sea level always use pull-up/push-down. Newly created water uses -W. Option
-M is required for add/replace and is not used for remove.
-L<field>/<minimum>
Set the value reached at the requested dry-land surface by -Ml. Repeat for every selected
field when linear extension is used.
-Q<surface.nc>[+c<classification.nc>]
Write the old surface coordinate used by the operation and, optionally, classification codes:
0 unresolved, 1 land, and 2 wet. Unresolved columns outside the model coverage remain NaN.
Output -Z coordinate scales and units apply to these grids. +z scales the old surface, but
classification codes cannot be scaled.
-Sa|c|e|l|n|s<p>[+g[<maxgap>]]
Choose GMT vertical interpolation: Akima (a), cubic (c), step-up (e), linear (l), nearest
(n), or smoothing spline (s<p>) with non-negative fit parameter p. Linear (l) is the default.
Append +g to bridge internal missing layers, optionally only when the bracketing z-
coordinate distance does not exceed maxgap. This is vertical interpolation. Option -H fills
enclosed holes in native x-y layers, whereas common -n controls horizontal resampling onto the
-R/-I lattice.
-T<zmin>/<zmax>/<dz>
Set the working output vertical axis, where zmin is the top, zmax is the bottom, zmin < zmax,
and dz > 0. Values above sea level are negative and values below it are positive. Without -T,
retain the smallest transformed input spacing and derive increasing bounds. For example, -T-
5/20/0.1 spans 5 km elevation to 20 km depth.
-W<field>/<value>[+t<tolerance>]
Set one field's water value and the tolerance used to recognize old water during inference.
Repeat for every selected field when an old wet surface must be inferred, or when add/replace
creates or rebuilds water. It is not required merely because wet cells exist when -E supplies
the old surface and the operation removes or retains their water unchanged. NaN is accepted as
a water value.
Example: for model.nc?vp,vs,rho in km/s, km/s, and g/cm^3, -Wvp/1.5 -Wvs/0 -Wrho/1.03
assigns those values to newly created water and uses them to recognize an existing water
column.
Example: -Wvp/1.5+t0.05 treats existing vp values from 1.45 through 1.55 as water during
surface inference. Newly created water is written as exactly 1.5.
-Z[+x<sx>][+X<xunit>][+y<sy>][+Y<yunit>][+z<sz>][+Z<zunit>][+v<scales>][+V<units>]
Transform output coordinates and selected fields after all topography operations. +x, +y, and
+z scale output coordinates. +X, +Y, and +Z set coordinate units. +v supplies one broadcast
scale or one scale per selected field. +V sets field units. Field lists follow the model
selector order. Scaling occurs in place and does not reorder coordinates, model layers,
fields, or ancillary variables. A negative axis scale therefore produces a decreasing output
axis.
Example: for model.nc?vp,vs, -Z+x0.001+Xkm+y0.001+Ykm+z-1000+Zm+v0.001,0.001+Vkm/s,km/s
converts x and y from m to km, restores a positive-up z axis in m, and converts vp and vs
from m/s to km/s without reordering the cube.
-V[q|e|w|t|i|c|d]
Change the verbosity level (currently w). Choose among 7 levels; each level adds more detailed
messages:
q: Quiet, not even fatal error messages.
e: Error messages only.
w: Warnings [Default when no -V is given].
t: Timings (time-intensive operations only).
i: Informational messages (or just -V with no level appended).
c: Compatibility warnings.
d: Debugging messages.
-di<nodata>[+c<col>] Replace any <nodata> values in input data with NaN.
+c Append first column to be affected [2].
-n[b|c|l|n][+a][+b<BC>][+c][+t<threshold>]
Specify a grid interpolation directive:
b: B-spline.
c: Bicubic spline [Default].
l: Bilinear spline.
n: Nearest-neighbor value, i.e., no spline at all.
Optional modifiers:
+a Switch off anti-aliasing (except for linear) [Default is ON].
+b Change boundary conditions, where <BC> can be either g for geographic, p for periodic, or
n for natural boundary conditions. For p and n you may optionally append x or y [Default
is both]: x for periodic boundary conditions on x, y for periodic boundary conditions on
y [Default is Natural conditions, unless grid is known to be geographic].
+c Clip interpolated grid to input z-min/max [Default may exceed limits].
+t Change the minimum weight in vicinity of NaNs. A <threshold> of 1.0 requires all nodes
involved in interpolation to be non-NaN; 0.5 will interpolate about half way from a non-
NaN to a NaN node [Default is 0.5].
-^ (or -) Print short synopsis message.
-+ (or +) Print longer synopsis message.
-? (or no arguments) Print this usage message.
--PAR=<value> Temporarily override GMT default setting(s) (repeatable).
(See gmt.conf documentation for GMT default parameters).
OPERATION ORDER AND CAVEATS:
Input coordinate and field scaling occurs first. The transformed z axis must increase from top
to bottom. Scale a contrary axis with the model's +z modifier. Next -H optionally fills
strictly internal holes in each native x-y model layer, and -n resamples the requested -R/-I
model lattice. Topobath then classifies wet and dry columns, obtains the old surface from -E
or model inference, and applies -O/-M. Vertical -S interpolation and optional gap bridging are
used while sampling the moved columns. Output -Z coordinate and field scaling occurs last and
does not reorder coordinates, layers, fields, or ancillary variables.
Classification is independent of surface sign, so dry land below sea level is valid (e.g.,
Death Valley or Dead Sea). A requested surface above sea level in a classified wet cell is
physically inconsistent and is clamped to sea level with a warning. This accommodates small
coastline differences between masks and interpolated relief. Below each identified surface,
the existing column moves with the adjustment. Missing values in one field are preserved even
when another field identifies the surface. A large pull-up can move a column far enough that
its output base lies outside the source model. Those exposed values are written as NaN and a
warning is reported.
Examples
See topobath examples.