grdinterpolate
- 官方文档:
- 简介:
从一个 3D 网格 (cube) 或堆叠的一系列 2D 网格 (layer) 中采样一个 3D 网格，或一系列 2D 网格，或 1D 序列。
语法
gmt grdinterpolate cube | grd1 grd2 … -Goutfile [ -D[+xxname][+yyname][+zzname][+vvname][+sscale][+ooffset][+ninvalid][+ttitle][+rremark] ] [ -Eline ] [ -Fa|c|e|l|n|sp[+d1|2] ] [ -Rregion ] [ -Sx/y|pointfile[+hheader] ] [ -T[min/max/]inc[+i|n] |-Tfile|list ] [ -V[level] ] [ -Z[levels] ] [ -bbinary ] [ -dnodata[+ccol] ] [ -eregexp ] [ -fflags ] [ -ggaps ] [ -hheaders ] [ -iflags ] [ -nflags ] [ -oflags ] [ -qflags ] [ -sflags ] [ -:[i|o] ] [ --PAR=value ]
描述
grdinterpolate 读取一个 netCDF 格式的 3D 网格（以下称之为 cube ）或一系列 2D 网格层 （以下称每个 2D 网格层为 layer，按顺序堆叠形成第三个维度）， 沿着第三个维度（以下称之为 level ）对一个或多个 level 做插值。 要求 cube 的数据维度包含两个共有的 x 和 y 维度，第三个维度通常为距离或时间。 输出的插值结果可以保存为一个 cube 或者一系列 layer。 或者，可以在一个或多个任意点位 (x/y) 沿着 level 轴插值 (-S) ， 输出的表数据包含一个或多个 level 系列；也可以沿着任意垂直切片 (-E) 做插值， 输出为 2D 网格形式。
必选选项
- cube
指定 netCDF 格式的 3D 输入网格。或者使用 -Z 时，此处可指定一系列 2D 网格层 (layer)。
可选选项
- -D[+xxname][+yyname][+zzname][+sscale][+ooffset][+ninvalid][+ttitle][+rremark]
给定网格文件头段中的基本信息：
+xxname X变量名及其单位，格式为 varname [unit]，比如 “distance [km]”
+sscale 读入网格数据后要乘以的因子，默认值为 1
+ooffset 读入数据后并乘以因子后要加的常数，默认值为 0
+ninvalid 指定值 invalid 用于表示该节点处无有效值，默认为NaN
+ttitle 网格文件的标题
+rremark 网格文件的注释信息
其它说明：
未指定的项其值保持不变
可以给一个空值以重置某一项，比如使用 +t 而不指定标题则设置标题为空
若字符串中包含空格则需要用引号括起来
若字符串中包含加号 +，需要使用 + 对其进行转义。 或者可以使用单引号套双引号的方式，例如 ‘“title with + inside”’
若字符串中使用了shell变量，且变量值中包含加号，则需要使用
${variable/+/\\+}对于地理数据（比如 -fg ），xname 和 yname 会自动设置
- -Etable|origin|line[,line,…][+aaz][+c][+d][+g][+iinc][+llength][+nnp][+oaz][+p][+rradius][+x]
创建一条或多条测线。有以下几种方式：
给定 table 文件，其中包含点位。
给定起点 origin （和方位角 azimuth ，长度 length ，和下方相关选项）。
给定一条或多条测线 line 。每条 line 的格式为 start/stop ， 其中 start 和 stop 格式为 lon/lat 或者两个字母指定地图上（或者使用 -R 指定的范围内）的位置 （例如BL,TR等，详见 锚点 ），你也可以使用 Z-,Z+ 代表网格（或多维网格的顶层，仅限 grdinterpolate ） 的全局最小值位置和最大值位置。
支持以下选项：
+a - 设置测线走向 azimuth （配合 +l 指定长度），测线始于起点 origin 。
+c - 当多条测线首尾重合时，将其连接起来合并成一条测线。
+d - 输出结果增加沿测线的距离。
+g - 当测线起始点位于同一经度/纬度时，报告沿测线的纬度/经度变化而不是大圆弧距离。
+i - 指定采样间隔 inc ，若不指定则默认为网格最小间距的一半。
+l - 指定测线长度 length ，与 +a （或 +o ）以及 origin 一起使用。
+n - 指定采样点数 np （这间接确定了采样间隔）。
+o - 与 +a 类似，但将 origin 作为测线中点，指定测线走向 az。
+p - 当测线起始点位于同一纬度时，沿着纬线采样而不是大圆弧。
+r - 创建圆形测线，将 origin 作为圆心，指定 radius 作为半径。要求使用 +n 或 +i 。
+x - 地理距离的测量沿着恒向线（rhumb lines 或 loxodromes， 即沿着该方向行进时行进方向始终与经线保持相同夹角）而不是默认的大圆弧。
注：
对于设置长度的选项，仅能设置一种单位，设置不同单位的结果错误。 地理坐标下若不指定 -E 中的长度单位默认为千米(k)。
对于地理坐标，可使用 -j 指定距离计算方式（默认计算大圆弧路径）。
- -Fa|c|e|l|n|sp[+d1|2]
选择 1-D 插值方式:
a - Akima 样条（默认）
e - 阶梯曲线
c - 自然三次样条 (Natural cubic spline)
l - 线性插值
n - 不插值，取最近的数据点作为插值后的值
s - 平滑三次样条，指定拟合因子 p 以权衡样条拟合和曲率
用户也可以更改 GMT 默认的插值方式，见 GMT_INTERPOLANT 。 可以加上 +d1 或 +d2 来计算样条的一阶导或二阶导。 注： 如果使用 e 或 n 的插值方法来计算导数，则结果恒为 0。
- -Rxmin/xmax/ymin/ymax[+r][+uunit] (more …)
指定数据范围。
- -Sx/y|pointfile[+hheader]
在给定点位 (x/y) 或 pointfile 文件中的多个点位处沿 level 轴插值，输出空间或时间序列。 如果需要由起点和终点或类似点定义的一系列点位，可以首先使用 project 制作这样的点位文件。 默认将采样 cube 的每个 level。使用 -T 可自定义 level。 在输出的插值结果中，level 值将插入点位文件的第 3 列，后跟其它输入列，最后一列为每个 level 对应的采样值。 加上 +h 可在输出的每段序列段头后增加头段注释 header 。 在输出时，输入的每个点位的尾随文本将作为每个点位输出序列的段头。 默认情况下，表数据输出到标准输出。使用 -G 可指定输出文件名。 或者如果希望每个点位的序列输出到单独文件，可在 -G 指定的文件名中包含 C 语言的整数格式 (e.g., %d)， 此时将根据点数创建单独的输出文件。
备注
当前 GMT 6.5 版本 -S 无法与 -T 同时使用（会报错），该 bug 会在下个版本修复，详见 PR #8733 。 对于 -T 选项，程序内部是使用 sample1d 模块对 level 之间进行插值， 因此可将 -S 和 -T 拆分成 grdinterpolate -S 和 sample1d -T -N2 分开执行， 以暂时回避此 bug。
- -T[min/max/]inc[+a][+i|n][+u] | [-Tfile|list]
定义要等距采样的 level 范围（从 min 到 max ）和步长 inc （默认使用输入 cube 的每个 level ）。 具体详见 生成一维数组 。
注：
- -Z[levels]
读取命令行中传入的所有 2D 网格作为一个 cube 的每个 layer（默认读取一个 cube 文件）。 可指定 levels 来控制这些 2D 网格如何组合成一个 cube，即为这些 2D 网格分配其所在level。 levels 的设置方式和 -T 一致，若不设置则默认 level 为从 0 开始的整数。
- -V[level] (more …)
设置 verbose 等级 [w]
- -:[i|o] (more …)
交换输入或输出中的第一和第二列
- -bi[ncols][type][w][+l|b] (more …)
设置二进制输入数据的格式
- -bo[ncols][type][w][+l|b] (more …)
设置二进制输出的数据格式
- -d[i|o]nodata (more …)
将输入数据中等于 nodata 的记录替换为 NaN，或将输出数据中值为 NaN 的记录替换为 nodata
- -e[~]“pattern” | -e[~]/regexp/[i] (more …)
筛选或剔除匹配指定模式的数据记录
- -f[i|o]colinfo (more …)
指定输入或输出列的数据类型
- -g[a]x|y|d|X|Y|D|[col]zgap[+n|p] (more …)
确定数据或线段的间断
- -h[i|o][n][+c][+d][+msegheader][+rremark][+ttitle] (more …)
跳过或生成指定数目的头段记录
- -icols[+l][+sscale][+ooffset][,…][,t[word]] (more …)
设置输入数据列及简单变换（0表示第一列，t 表示文本列）
- -n[b|c|l|n][+a][+bBC][+c][+tthreshold] (more …)
设置网格文件的插值方式
- -ocols[,…][,t[word]] (more …)
设置输出数据列（0表示第一列，t 表示文本列）
- -q[i|o][~]rows[+ccol][+a|f|s] (more …)
筛选输入或输出的行或数据范围
- -s[cols][+a|+r] (more …)
设置 NaN 记录的处理方式
- -^ 或 -
显示简短的帮助信息，包括模块简介和基本语法信息（Windows下只能使用 -）
- -+ 或 +
显示帮助信息，包括模块简介、基本语法以及模块特有选项的说明
- -? 或无参数
显示完整的帮助信息，包括模块简介、基本语法以及所有选项的说明
- --PAR=value
临时修改GMT参数的值，可重复多次使用。参数列表见 配置参数
生成一维数组
下面将展示如何使用 gmtmath 生成一维数组（其中大部分操作也可通过 linux 中的 seq 命令方便地实现）
不使用任何子选项
以 0.1 为步长，生成 3.1 到 4.2 等等间隔序列
gmt math -o0 -T3.1/4.2/0.1 T =
3.1
3.2
...
4.2
+a 选项
该选项不对生成的数列进行运算操作，而是将该数列以列的形式追加到输出表数据，类似 linux 中的 paste 命令。
+b 选项
以 3 和 20 分别为起点和终点，创建一个 2 的整数幂的列表，将生成 的整数幂以 1 为步长取值 -T3/20/1+b
gmt math -o0 -T3/20/1+b T =
4
8
16
+l 选项
以 7 和 135 分别为起点和终点，先创建一个 10 的整数幂的列表，当 inc 为 1 时，输出该列表；当 inc 为 2 时，同时输出该列表中每个元素的 1，2，5 倍值（不包括超出原始数据范围的部分）； 当 inc 为 3 时，同时输出 该列表中每个元素的 1，2，…，9 倍值（同样不包括超出原始数据范围的部分）。 -T7/135/2+l 将会生成如下数列
gmt math -o0 -T7/135/2+l T =
10
20
50
100
inc 为负整数时，将实现如下效果
gmt math -o0 -T1e-4/1e4/-2+l T =
0.0001
0.01
1
100
10000
+i 选项
该选项将以 1 作为默认步长，第三个数字做为每步长中的数字个数 length。假设需要在 1 分钟内每 24 秒生成一帧
gmt math -o0 -T0/60/24+i T =
0
0.0416666666667
0.0833333333333
0.125
0.166666666667
...
60
+n 选项
生成固定长度的数列。以 3.44 和 7.82 为 起点和终点，生成长度为 5 的等间隔序列
gmt math -o0 -T3.44/7.82/5+n T =
3.44
4.535
5.63
6.725
7.82
-T 后直接加文件或逗号分隔的数列
-T 后可以直接加文件，文件中即为要生成的列表；另外，可以直接使用逗号分隔， 将列表附加到 -T 选项后，如下为 Fibonacci 数列前 6 项
gmt math -o0 -T0,1,1,2,3,5 T =
0
1
1
2
3
5
注 ： 如果数列只包含一个值，必须在其后加逗号以表明仍是一个数列
+u 选项
如果文件或者逗号分隔的数列中有重复数字或未排序，可以使用 +u 选项去重并排序。
+t 选项
生成绝对时间序列。在 inc 后分别添加 y ，o ，d ，h ， m 和 s 表示时间步长的单位为年，月，日，时，分，秒。在其后附加 +t 选项， 可以进一步强调生成时间序列，也可以不加
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
生成距离序列
如果输入文件中包含两列以上的数据，可以使用前两列计算距离，并生成等距序列。 在 inc 后分别添加 d，m，s，e，f，M，n 和 u 表示距离步长的单位为度，分，秒，米，英尺，公里，英里，海里，英尺。如果 为笛卡尔坐标，使用 c 作为距离步长单位。
+e 选项
如果只给定 inc 而从数据中获取最大值和最小值，则 (max - min)/inc 可能 不是整数，GMT 讲会自动对 inc 进行一定的调整。如果不想调整 inc ，则可以 使用 +e 选项，GMT 会固定最小值，适当调整最大值。
文件顺序
如果用户在命令行提供一系列的 2D 网格，并使用 -Z 分配各自的 level， 则用户必须确保输入的 2D 网格顺序和 -Z 分配的 level 顺序相匹配。 除非文件名是按照词典顺序排列，否则应谨慎使用通配符来列举所有网格 (e.g., *.nc)。
时间坐标
GMT 可以识别 netCDF 网格文件中的时间坐标。网格文件中变量的 unit 属性会被解析以 确定网格文件中时间坐标的起算点和单位。这些时间坐标值会被进一步根据 TIME_UNIT 和 TIME_EPOCH 转换为 GMT 内部时间值。输出时，默认以相对时间的形式输出，也可以使用 -f 选项 指定以绝对时间方式输出。
创建序列
-S 选项涉及表格的读取和输出，因此 grdinterpolate 模块也可以使用和表格输入输出相关的标准选项，例如 -b 选项 ，-i 和 -o 选项 等。 由于使用 -S 提供的点位坐标不要求正好位于网格节点， 本模块内部会使用 grdtrack 模块对每个layer做采样， 故 -n 选项 也可使用。
示例
从 cube 文件 temperature.nc 中使用三次样条插值提取一个 level=3400 的新 layer:
gmt grdinterpolate temperature.nc -T3400 -Fc -Gtemp_3400.nc
从一系列的 2D 网格文件 layers_*.nc 中使用线性插值提取一个 level=3400 的新 layer， 其中每个网格文件对应的 level 值按照 z.txt 文件分配:
gmt grdinterpolate layers_*.nc -Zz.txt -T3400 -Fl -Gtemp_3400.nc
从 cube 文件 temperature.nc 中使用 Akima 样条对所有 level 重新采样，输出一个新 cube ， 其中新 level 从 1500 到 2500，间隔 50:
gmt grdinterpolate temperature.nc -T1500/2500/50 -Gtemperature_1500_2500.nc -Fa
与上面相同，但是每个 level 输出为一个单独的 2D 网格:
gmt grdinterpolate temperature.nc -T1500/2500/50 -Gtemperature_%4.0f.nc -Fa
从一系列的 2D 网格文件 deformation_*.nc 提取出一个点位的时间序列（此处 level 为时间）， 其中每个网格文件对应的时间点按照 dates.txt 文件分配，在输出序列的段头添加 “Some like it hot”:
gmt grdinterpolate deformation_*.nc -Zdates.txt -S115W/33N+h"Some like it hot" > record.txt
从 cube 文件 S362ANI_kmps.nc 中提取穿过夏威夷热点的垂直切片，采样 vs（各向同性剪切波速度）， 测线起始点位于同一纬度，测线距离以经度度数表示:
gmt grdinterpolate S362ANI_kmps.nc?vs -E180/20/220/20+i1d+g+p -T25/500/25 -Gslice.nc