sample1d

Resample 1-D table data using splines

Synopsis

gmt sample1d [ table ] [ -A[f|p|m|r|R][+d][+l] ] [ -C[section/]master|cpt|color\(_1\),color\(_2\)[,color\(_3\),…][+h[hinge]][+idz][+u|Uunit][+sfname] ] [ -E ] [ -Fa|c|e|l|n|sp[+d1|2] ] [ -Ncol ] [ -T[min/max/]inc[+a][+i|n][+u] | [-Tfile|list]] [ -V[level] ] [ -Wcol ] [ -bbinary ] [ -dnodata[+ccol] ] [ -eregexp ] [ -fflags ] [ -ggaps ] [ -hheaders ] [ -iflags ] [ -jflags ] [ -oflags ] [ -qflags ] [ -sflags ] [ -wflags ] [ -:[i|o] ] [ --PAR=value ]

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

Description

sample1d reads a multi-column ASCII [or binary] data set from file [or standard input] and interpolates the time-series or spatial profile at locations where the user needs the values. The user must provide the column number of the independent (monotonically increasing or decreasing) variable, here called time (it may of course be any type of quantity) when that is not the first column in data set. Equidistant or arbitrary sampling can be selected. All columns are resampled based on the new sampling interval. Several interpolation schemes are available, in addition to a smoothing spline which trades off misfit for curvature. Extrapolation outside the range of the input data is not supported. By giving a CPT via -C we will add r, g, b, a columns to the output based on the last input data column.

Required Arguments

table: This is one or more ASCII [of binary, see -bi] files with one column containing the independent time variable (which must be monotonically in/de-creasing) and any number of optional columns holding other data values. If no file is provided, sample1d reads from standard input.

Optional Arguments

-A[f|p|m|r|R][+d][+l]

For track resampling (if -T…unit is set) we can select how this is to be performed via these directives:

f - Keep original points, but add intermediate points if needed; note this selection does not necessarily yield equidistant points [Default].
m - Same as +f, but first follow meridian (along y) then parallel (along x).
p - Same as +f, but first follow parallel (along x) then meridian (along y).
r - Resample at equidistant locations; input points are not necessarily included in the output.
R - Same as r, but adjust given spacing to fit the track length exactly.

A few modifiers can also be used:

+d - Delete duplicate input records (identified by having no change in the time column, and
+l - If distances should be measured along rhumb lines (loxodromes).

Note: Calculation mode for loxodromes is spherical, hence -je cannot be used in combination with +l.

-C[section/]master|cpt|color\(_1\),color\(_2\)[,color\(_3\),…][+h[hinge]][+idz][+u|Uunit][+sfname]

Determine the color components based on the z values. Append name of a master CPT, an input CPT file or a comma-separated list of colors from which to build a CPT. If no argument is given then under modern mode we select the current CPT, if it is available. Generally, the input can be many things:

A standard GMT master CPT file, e.g., earth (see Of Colors and Color Legends) and can be either addressed by master or section/master (without any .cpt extension).
File name of already custom-made cpt file (e.g, my_colors.cpt).
Build a linear continuous CPT from color\(_1\),color\(_2\)[,color\(_3\),…] automatically, where z starts at 0 and is incremented by one for each color. In this case, color\(_i\) can be a r/g/b (e.g., 255/100/75), h-s-v triplet (e.g., 180-0.8-1), a c/m/y/k quadruple (e.g., 80/50/40/30), an HTML hexadecimal color (e.g. #aabbcc), or a color name. No spaces between commas are allowed.

A few modifiers pertains to hinges and units:

+h - If given a master CPT with soft hinges then you can enable the hinge at data value hinge [0], whereas for hard-hinge CPTs you can adjust the location of the hinge [0] but not disable it.
+i - Append increment dz to quantize the grid range [Default uses the exact grid range].
+s - Append fname to save the finalized CPT in a disk file. This is useful when the CPT is generated automatically, but if used, must be at the end of the -C option.
+u - For any other master CPT, you may convert their z-values from other distance Units to meter by appending the original unit code.
+U - Likewise, you may convert their z-values from meter to other Units by appending the desired unit code.

Note: We append four extra output columns at the end of the output records. These are the red, green, and blue components color components (all in 0-255 range) and the transparency value (0-100% range, with 0 meaning opaque). If the color table is a gray table, then we output gray (0-255) instead of three identical components.

-E: If the input dataset contains records with trailing text then we will attempt to add these to output records that exactly match the input times. Output records that have no matching input record times will have no trailing text appended [Default ignores trailing text].

-Fl|a|c|e|l|n|sp[+d1|2]

Choose the desired 1-D spline interpolant from one of these directives:

a - Akima spline [Default].
e - Step curve.
c - Natural cubic spline.
l - Linear interpolant.
n - No interpolation, just select nearest point.
s - Smoothing cubic spline; append fit parameter p for an approximate fit.

You may change the GMT default interpolant; see GMT_INTERPOLANT in your gmt.conf file. You may optionally evaluate the first or second derivative of the spline by appending +d1 or +d2, respectively. Note: If you use derivatives with directives e or n then the result is always zero.

_images/GMT_splines.png — The -F option lets you choose among several interpolants, including one that is approximate (the smoothing spline). You can also specify that you actually need a derivative of the solution instead of the value.

-Ncol: Sets the column number of the independent time variable [Default is 0 (first)].

-T[min/max/]inc[+a][+i|n][+u] | [-Tfile|list]: Make evenly spaced time-steps from min to max by inc [Default uses input times]. The form -Tlist means a online list of time coordinates like for example: -T13,15,16,22.5. For details on array creation, see Generate 1-D Array. Note: For resampling of spatial (x, y) or (lon, lat) series you must give an increment with a valid distance unit; see Units for map units or use c if plain Cartesian coordinates. The first two columns must contain the spatial coordinates. From these we calculate distances in the chosen units and interpolate using this parametric series.

-V[level]: Select verbosity level [w]. (See full description) (See cookbook information).

-Wcol: Sets the column number of the weights to be used with a smoothing cubic spline. Requires -Fs.

-birecord[+b|l] (more …): Select native binary format for primary table input. [Default is 2 (or at least the number of columns implied by -T)].

-borecord[+b|l] (more …): Select native binary format for table output. [Default is same as input].

-d[i|o][+ccol]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.

-gx|y|z|d|X|Y|Dgap[u][+a][+ccol][+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][+ddivisor][+sscale|d|k][+ooffset][,…][,t[word]] (more …): Select input columns and transformations (0 is first column, t is trailing text, append word to read one word only).

-je|f|g (more …): Determine how spherical distances or coordinate transformations are calculated.

-ocols[+l][+ddivisor][+sscale|d|k][+ooffset][,…][,t[word]] (more …): Select output columns and transformations (0 is first column, t is trailing text, append word to write one word only).

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

-wy|a|w|d|h|m|s|cperiod[/phase][+ccol] (more …): Convert an input coordinate to a cyclical coordinate.

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

-^ 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.

Units

For map distance unit, append unit d for arc degree, m for arc minute, and s for arc second, or e for meter [Default unless stated otherwise], f for foot, k for km, M for statute mile, n for nautical mile, and u for US survey foot. By default we compute such distances using a spherical approximation with great circles (-jg) using the authalic radius (see PROJ_MEAN_RADIUS). You can use -jf to perform “Flat Earth” calculations (quicker but less accurate) or -je to perform exact geodesic calculations (slower but more accurate; see PROJ_GEODESIC for method used).

ASCII Format Precision

The ASCII output formats of numerical data are controlled by parameters in your gmt.conf file. Longitude and latitude are formatted according to FORMAT_GEO_OUT, absolute time is under the control of FORMAT_DATE_OUT and FORMAT_CLOCK_OUT, whereas general floating point values are formatted according to FORMAT_FLOAT_OUT. Be aware that the format in effect can lead to loss of precision in ASCII output, which can lead to various problems downstream. If you find the output is not written with enough precision, consider switching to binary output (-bo if available) or specify more decimals using the FORMAT_FLOAT_OUT setting.

Generate 1-D 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 \(\log_2\) of min and max, get their nearest integers, build an equidistant \(\log_2\)-array using inc integer increments in \(\log_2\), then undo the \(\log_2\) 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 \(\log_{10}\) 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 +i if inc is a fractional number and it is cleaner to give its reciprocal value instead. To set up times for a 24-frames per second animation lasting 1 minute, run:

gmt math -o0 -T0/60/24+i T =
0
0.0416666666667
0.0833333333333
0.125
0.166666666667
...

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, let inc be 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

Notes: (1) If you need to pass the list nodes via a dataset file yet be understood as a list (i.e., no interpolation), then you must set the file header to contain the string “LIST”. (2) Should you need to ensure that the coordinates are unique and sorted (in case the file or list are not sorted or have duplicates) then supply the +u modifier.

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

If the module allows you to set up an absolute time series, append a valid time unit from the list year, month, 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).

Notes

The smoothing spline s(t) requires a fit parameter p that allows for the trade-off between an exact interpolation (fitting the data exactly; large p) to minimizing curvature (p approaching 0). Specifically, we seek to minimize

\[F_p (s)= K (s) + p E (s), \quad p > 0,\]

where the misfit is evaluated as

\[E (s)= \sum^n_{i=1} \left [ \frac{s(t_i) - y_i}{\sigma_i} \right ]^2\]

and the curvature is given by the integral over the domain of the second derivative of the spline

\[K (s) = \int ^b _a [s''(t) ]^2 dt.\]

Trial and error may be needed to select a suitable p.

Examples

Note: Below are some examples of valid syntax for this module. The examples that use remote files (file names starting with @) can be cut and pasted into your terminal for testing. Other commands requiring input files are just dummy examples of the types of uses that are common but cannot be run verbatim as written.

To resample the file profiles.tdgmb, which contains (time,distance,gravity,magnetics,bathymetry) records, at 1 km equidistant intervals using Akima’s spline, use:

gmt sample1d profiles.tdgmb -N1 -Fa -T1 > profiles_equi_d.tdgmb

To resample the file depths.txt at positions listed in the file grav_pos.dg, using a cubic spline for the interpolation, use:

gmt sample1d depths.txt -Tgrav_pos.dg -Fc > new_depths.txt

To resample the file points.txt every 0.01 from 0-6, using a cubic spline for the interpolation, but output the first derivative instead (the slope), try:

gmt sample1d points.txt -T0/6/0.01 -Fc+d1 > slopes.txt

To resample the file track.txt which contains lon, lat, depth every 2 nautical miles, use:

gmt sample1d track.txt -T2n -AR > new_track.txt

To do approximately the same, but make sure the original points are included, use:

gmt sample1d track.txt -T2n -Af > new_track.txt

To obtain a rhumb line (loxodrome) sampled every 5 km instead, use:

gmt sample1d track.txt -T5k -AR+l > new_track.txt

To sample temperatures.txt every month from 2000 to 2018, use:

gmt sample1d temperatures.txt -T2000T/2018T/1o > monthly_temp.txt

To use a smoothing spline on a topographic profile for a given fit parameter, try:

gmt sample1d @topo_crossection.txt -T300/500/0.1 -Fs0.001 > smooth.txt