grdinterpolate¶

Interpolate 2-D grids or 1-D series from a 3-D data cube

Synopsis¶

Note: No space is allowed between the option flag and the associated arguments.

Description¶

grdinterpolate reads a single 3-D netCDF data cube (or a set of 2-D layers) and interpolates along the 3rd dimension for one or more output levels. The data cube must be organized with one or more layers representing the x and y dimensions while the 3rd dimension may represent distance or time; we refer to this dimension as the level. The output layers may be written as a single 3-D cube or as a set of 2-D layers. Alternatively, we interpolate the cube along the level-axis at one or more arbitrary (x/y) coordinates (-S), resulting in a data table with one or more level-series, or we slice the 3-D cube along an arbitrary vertical slice and write that 2-D slice to a grid file (-E).

Required Arguments¶

3Dgrid: Name of a 3-D netCDF data cube to be interpolated. Alternatively, with -Zi, you can specify a set of 2-D grid layers instead.

-Goutfile: This is the output 3D data cube file. If -T only selects a single layer then the data cube collapses to a regular 2-D grid file. If -Zo is used then outfile must contain a C-format statement for a floating point number. Also see -S for a similar use to write individual level-series tables.

-T[min/max/]inc[+n] | -Tfile|list: Make evenly spaced time-steps from min to max by inc [Default uses input times]. For details on array creation, see Generate 1D Array.

Optional Arguments¶

-Etable|line[,line,…][+aaz][+g][+iinc][+llength][+nnp][+oaz][+p][+rradius][+x]: Specify a crossectinonal profile via a file or from specified line coordinates and modifiers. If a file, it must be contain a single segment with either lon lat or lon lat dist records. These must be equidistant. Alternatively, the format of each line is start/stop, where start or stop are lon/lat (x/y for Cartesian data). You may append +iinc to set the sampling interval; if not given then we default to half the minimum grid interval. If your line starts and ends at the same latitude you can force sampling along the parallel with +p [great circle]. For a line along parallels or meridians you can add +g to report degrees of longitude or latitude instead of great circle distances starting at zero. Append +x to compute distances along a loxodrome (rhumbline) instead of great circle. Instead of two coordinates you can specify an origin and one of +a, +o, or +r. The +a sets the azimuth of a profile of given length starting at the given origin, while +o centers the profile on the origin; both require +l. For circular sampling specify +r to define a circle of given radius centered on the origin; this option requires either +n or +i. The +nnp modifier sets the desired number of points, while +llength gives the total length of the profile. Also note that only one distance unit can be chosen. Giving different units will result in an error. If no units are specified we default to great circle distances in km (if geographic). If working with geographic data you can use -j to control distance calculation mode [Great Circle]. Use -G to set the output grid file name.

-Fl|a|c|n[+1|2]: Choose from l (Linear), a (Akima spline), c (natural cubic spline), and n (no interpolation: nearest point) [Default is -Fa]. You may change the default interpolant; see GMT_INTERPOLANT in your gmt.conf file. You may optionally evaluate the first or second derivative of the spline by appending +1 or +2, respectively.

-Rxmin/xmax/ymin/ymax[+r][+uunit] (more …): Specify the region of interest. Using the -R option will select a subsection of the grid. If this subsection exceeds the boundaries of the grid, only the common region will be output.

-Sx/y|pointfile[+hheader]: Rather that compute gridded output, create tile/spatial series through the stacked grids at the given point (x/y) or the list of points in pointfile. If you need a series of points defined by an origin and an end point or similar, you can make such a file first with project. By default we simply sample the cube at each level. Use -T to interpolate the series. The grid level (e.g., depth or time) will be appended as the last numerical value in the series records. Use the optional +h modifier to append header to the trailing text of these input points. On output the trailing text will become the segment header for the series that originate from each point. By default, the table output is written to standard output. Use -G to specify a file name. Alternatively, if you wish each series to be written to its own data file, let the filename in -G have a C-format integer specifier (e.g., %d) and we will use the running point number to create unique file names.

-V[level] (more …): Select verbosity level [w].

-Zilevels|o: Use -Zi to obtain the levels and then we read the corresponding number of 2-D input grids given on the command line [Default is a single 3-D data cube]. The levels are specified the same way as in -T. Use -Zo to write the 3-D data cube as a series of 2-D grids instead. If used, then the outgrid name given by -G must contain a C-language format statement for a floating point number (for instance, one can try layer_%6.6f.grd) which will contain the level for each grid [Default is a 3-D data cube, unless only one layer is implied by -T].
-:: Toggles between (longitude,latitude) and (latitude,longitude) input/output in -Stable. [Default is (longitude,latitude)].

-bi[ncols][t] (more …): Select native binary format for primary input. [Default is 2 input columns].

-bo[ncols][type] (more …): Select native binary output. [Default is one more than input].

-d[i|o]nodata (more …): Replace input columns that equal nodata with NaN and do the reverse on output.

-e[~]“pattern” | -e[~]/regexp/[i] (more …): Only accept data records that match the given pattern.

-f[i|o]colinfo (more …): Specify data types of input and/or output columns.

-g[a]x|y|d|X|Y|D|[col]zgap[+n|p] (more …): Determine data gaps and line breaks.

-h[i|o][n][+c][+d][+msegheader][+rremark][+ttitle] (more …): Skip or produce header record(s).

-icols[+l][+sscale][+ooffset][,…][,t[word]] (more …): Select input columns and transformations (0 is first column, t is trailing text, append word to read one word only).

-n[b|c|l|n][+a][+bBC][+c][+tthreshold] (more …): Select interpolation mode for grids.

-ocols[,…][,t[word]] (more …): Select output columns (0 is first column; t is trailing text, append word to write one word only).

-q[i|o][~]rows[+ccol][+a|f|s] (more …): Select input or output rows or data range(s) [all].

-s[cols][+a][+r] (more …): Set handling of NaN records.

-^ or just -: Print a short message about the syntax of the command, then exit (NOTE: on Windows just use -).
-+ or just +: Print an extensive usage (help) message, including the explanation of any module-specific option (but not the GMT common options), then exit.
-? or no arguments: Print a complete usage (help) message, including the explanation of all options, then exit.
--PAR=value: Temporarily override a GMT default setting; repeatable. See gmt.conf for parameters.

Generate 1D Array¶

We will demonstrate the use of options for creating 1-D arrays via gmtmath. Make an evenly spaced coordinate array from min to max in steps of inc, e.g.,:

gmt math -o0 -T3.1/4.2/0.1 T =
1
2
3
4
5
6
7

Append +b if we should take log2 of min and max, get their nearest integers, build an equidistant log2-array using inc integer increments in log2, then undo the log2 conversion. E.g., -T3/20/1+b will produce this sequence:

gmt math -o0 -T3/20/1+b T =
4
8
16

Append +l if we should take log10 of min and max and build an array where inc can be 1 (every magnitude), 2, (1, 2, 5 times magnitude) or 3 (1-9 times magnitude). E.g., -T7/135/2+l will produce this sequence:

gmt math -o0 -T7/135/2+l T =
10
20
50
100

For output values less frequently than every magnitude, use a negative integer inc:

gmt math -o0 -T1e-4/1e4/-2+l T =
0.0001
0.01
1
100
10000

Append +n if inc is meant to indicate the number of equidistant coordinates instead. To have exactly 5 equidistant values from 3.44 and 7.82, run:

gmt math -o0 -T3.44/7.82/5+n T =
44
535
63
725
82

Alternatively, give a file with output coordinates in the first column, or provide a comma-separated list of specific coordinates, such as the first 6 Fibonacci numbers:

gmt math -o0 -T0,1,1,2,3,5 T =
0
1
1
2
3
5

If you only want a single value then you must append a comma to distinguish the list from the setting of inc.

If the module allows you to set up an absolute time series, append a valid time unit from the list year, month, week, day, hour, minute, and second to the given increment; add +t to ensure time column (or use -f). Note: The internal time unit is still controlled independently by TIME_UNIT. The first 7 days of March 2020:

gmt math -o0 -T2020-03-01T/2020-03-07T/1d T =
2020-03-01T00:00:00
2020-03-02T00:00:00
2020-03-03T00:00:00
2020-03-04T00:00:00
2020-03-05T00:00:00
2020-03-06T00:00:00
2020-03-07T00:00:00

A few modules allow for +a which will paste the coordinate array to the output table.

Likewise, if the module allows you to set up a spatial distance series (with distances computed from the first two data columns), specify a new increment as inc with a geospatial distance unit from the list degree (arc), minute (arc), second (arc), meter, foot, kilometer, Miles (statute), nautical miles, or survey foot; see -j for calculation mode. To interpolate Cartesian distances instead, you must use the special unit c.

Finally, if you are only providing an increment and will obtain min and max from the data, then it is possible (max - min)/inc is not an integer, as required. If so, then inc will be adjusted to fit the range. Alternatively, append +e to keep inc exact and adjust max instead (keeping min fixed).

File Order¶

If you provide a series of 2-D files and thus separately assigning the level via -Zi, then you must make sure that the order the grids are given on the command line matches the levels you provide via -Zi. Unless your files are named in lexical order you must be careful with using wildcards to list all the grids (e.g., *.nc).

Time Coordinates¶

Time coordinates in netCDF grids, be it the x, y, or z coordinate, will be recognized as such. The variable’s unit attribute is parsed to determine the unit and epoch of the time coordinate in the grid. Values are then converted to the internal time system specified by TIME_UNIT and TIME_EPOCH in the gmt.conf file or on the command line. The default output is relative time in that time system, or absolute time when using the option -f0T, -f1T, or -f2T for x, y, or z coordinate, respectively.

Series creation¶

The (optional) table-reading and table-producing -S option may require some of the standard common options associated with table i/o, such as -b, -i, o, etc., thus these options are available to grdinterpolate as well. Because the coordinates given via -S are not required to equal the coordinates of the grid nodes, we are resampling each 2-D layer at the given points via grdtrack, hence the availability of the -n option.

Examples¶

To extract a single, new 2-D layer from the temperature.nc 3-D cube for level 3400 using a cubic spline, try:

gmt grdinterpolate temperature.nc -T3400 -Fc -Gtemp_3400.nc

To extract a single, new 2-D layer from the 3-D cube implied by the individual grids layers_*.nc, with individual layer values given via z.txt, for level 3400 using a linear spline, try:

gmt grdinterpolate layers_*.nc -Ziz.txt -T3400 -Fl -Gtemp_3400.nc

To resample the the temperature.nc 3-D cube for all levels from 1500 to 2500 in steps of 50, using an Akima spline, try:

gmt grdinterpolate temperature.nc -T1500/2500/50 -Gtemperature_1500_2500.nc -Fa

The same, but this time write individual 2-D grids per layer, try:

gmt grdinterpolate temperature.nc -T1500/2500/50 -Gtemperature_%4.0f.nc -Fa -Zo

To extract a time-series through the grids deformation_*.nc at the location (115W, 33N), with the times of each grid provided by the file dates.txt, and append the string “Some like it hot” to the segment header for the series, try:

gmt grdinterpolate deformation_*.nc -Zidates.txt -S115W/33N+h"Some like it hot" > record.txt

To extract a vertical slice of the 3-D grid S362ANI_kmps.nc with seismic velocities that goes through the Hawaii hotspot, selecting cube vs (Isotropic Shear Velocity) and letting the distances be longitude degrees along the parallel, try:

gmt grdinterpolate S362ANI_kmps.nc?vs -E180/20/220/20+i1d+g+p -T25/500/25 -Gslice.nc