grdpmodeler

贡献者:

何星辰

最近更新日期:

2025-10-18


官方文档:

grdpmodeler

简介:

评估地理网格上的板块运动模型

grdpmodeler 读取一个地质年龄网格和一个板块运动模型,并计算若干种模型预测结果。 用户还可以选择提供一个多段格式的裁剪多边形文件;在这种情况下, 仅使用多边形内部的网格部分来计算模型预测,外部区域的网格值将被设置为 NaN。

语法

gmt grdpmodeler [ agegrdfile ] -Erot_file|ID1-ID2|lon/lat/angle[+i] -Sflags [ -Fpolygonfile ] [ -Goutgrid ] [ -Iincrement ] [ -Rregion ] [ -Tage ] [ -V[level] ] [ -bbinary ] [ -dnodata[+ccol] ] [ -hheaders ] [ -oflags ] [ -rreg ] [ -:[i|o] ] [ --PAR=value ]

注意:选项标志和相关参数之间不允许有任何空格。

必选选项

-Erot_file|ID1-ID2|lon/lat/angle[+i]

可通过以下三种方式之一指定旋转参数:

  • 提供一个包含旋转参数的文件。该文件必须为每个旋转包含一条记录, 每条记录的格式如下:

    lon lat tstart [tstop] angle [ khat a b c d e f g df ]

    其中 tstarttstop 以百万年(Myr)为单位,lon lat angle 均以度为单位。 tstarttstop 分别为阶段的起始与结束年龄。 若记录中未提供 tstop,则表示为一次总重建旋转,此时 tstop 默认设为 0, 且文件中其他记录也不应指定 tstop。 若旋转的协方差矩阵 C 可用,则必须按照括号中列出的九个可选项提供: 此处 C = (g/khat)*[ a b d; b c e; d e f ], 即 C 由三个行向量组成。 若拟合旋转的自由度 (df) 为 0 或未给出,则默认设置为 10000。 文件中空行以及以 # 开头的记录将被忽略。

  • 提供由两个板块 ID 组成的文件名(以连字符连接,例如 PAC-MBL), 程序将从 GPlates 旋转数据库中提取相应的旋转参数。 若未找到对应旋转,将返回错误。

  • 直接指定 lon/lat/angle,即旋转极的经度、纬度与旋转角度(均以度为单位), 用于单一的总重建旋转。

    无论采用哪种方式,都可以在参数后附加 +i, 表示希望对旋转进行反转(invert rotation)。

-Sflags

指定要计算的模型预测类型,可附加一个或多个标识符: - a:板块运动方位角(azimuth) - d:当前位置与洋中脊生成点的球面距离(单位:km) - s:板块运动阶段 ID(1 表示最年轻阶段) - v:板块运动速度(mm/yr) - w:板块旋转速率(度/Myr) - x:相对于生成点的经度变化量 - y:相对于生成点的纬度变化量 - X:地壳生成时的经度 - Y:地壳生成时的纬度

如果未指定任何参数,则默认计算全部 [adsvwxyXY]。

可选选项

ingrid

输入的地理网格文件名(经度、纬度坐标),其数值为地壳年龄(单位:Myr)。 如果未提供网格文件,可使用 -R-I 和可选的 -r 定义计算域。

-Fpolygonfile

指定一个多段闭合多边形文件,描述模型计算的区域内部。 多边形外部的网格将被设为 NaN。 [默认对整个网格进行模型计算。]

-Goutgrid[=ID][+ddivisor][+ninvalid][+ooffset|a][+sscale|a][:driver[dataType][+coptions]]

输出网格文件名。通过追加 =ID 可指定 网格格式。 对网格追加子选项可进行额外设置:

  • +d 将网格除以一个数 divisor [默认为 1]

  • +n 将网格中的 invalid 替换为 NaN

  • +o 将网格中的值进行一定的偏移,即加一个数 offset,或使用 a 自动对值 进行调整以保证整数网格的精度 [默认为 0]

  • +s 将网格中的值缩放,即乘以 scale,或使用 a 自动对网格缩放以保证整数 网格的精度 [默认为 1]

offset 操作位于 scale 之前;+sa 将会同时设置 +oa。如果需 使用 GDAL 指定网格格式,则 ID 应设置为 gd,并指定 :driver 和可选的数据 类型 dataType,以及 +coptions 选项传递给 GDAL 的 -co 选项。

-Ixinc[+e|n][/yinc[+e|n]]

指定X和Y方向的网格间隔

  • xincyinc 为 X 和 Y 方向的网格间隔。对于地理坐标,可以指定网格间隔单位 [默认单位为度]

  • +e 微调X和Y方向范围的最大值,使得其是网格间隔的整数倍(默认会微调 网格间隔以适应给定的数据范围)

  • +n 表明 xincyinc 不是网格间隔,而是X和Y方向的节点数。 此时会根据节点数、网格区域范围以及网格配准方式重新计算网格间隔。

注意:

  • yinc 设置为0,则表示其与 xinc 相同

  • 若使用 -Rgrdfile 选项,则网格间隔和配准方式已经根据网格文件自动初始化, 此时依然可以使用 -I-r 覆盖相应的值

-Rwest/east/south/north[/zmin/zmax][+r][+uunit]

指定感兴趣的区域。 |Add_-Rgeo|

该选项可以通过多种方式来定义区域:

  1. -Rwest/east/south/north 这是在使用经纬线呈直线的地图投影时指定地理区域的标准方式。 坐标可用十进制度表示,也可以使用 [±]dd:mm[:ss.xxx][W|E|S|N] 格式表示。

  2. -Rwest/south/east/north+r 当使用倾斜(斜轴)投影时,经线和纬线不再适合作为地图边界。 在这种情况下,可以通过给出矩形的左下角与右上角的地理坐标,并添加 +r 修饰符, 来保证地图输出区域为矩形,即使经纬线不是直线。

  3. -Rg-Rd 用于快速指定全球范围。-Rg 表示经度 0 到 360,纬度 -90 到 +90; -Rd 表示经度 -180 到 +180,纬度 -90 到 +90。

  4. -Rcode1,code2,…[+e|r|Rincs] 通过查阅 DCW(数字世界地图)数据库,间接地根据一个或多个国家代码来确定边界区域。 可以使用两个字符的 ISO 3166-1 alpha-2 国家代码(例如 NO)或完整国名(例如 Norway)。 如果要选择国家内的某个州(若可用),请追加 .state(如 US.TX),或使用完整州名(如 Texas)。 若要选择整个大陆,请直接拼写完整名称(如 -RAfrica )。 还可指定 DCW数据 的缩写或完整名称。所有名称大小写不敏感。 可附加以下修饰符:

    • +r:将区域边界调整为 incxinc/yincwinc/einc/sinc/ninc 的整数倍(默认不调整)。 例如,-RFR+r1 将法国的边界取整到最接近的整数度。inc 为正数表示扩大范围,负数表示缩小。

    • +R:在区域边界的基础上增加或减少 incxinc/yincwinc/einc/sinc/ninc 的数值(默认不扩展)。 inc 为正数表示扩大区域,负数表示缩小。

    • +e:将区域边界调整为 incxinc/yincwinc/einc/sinc/ninc 的整数倍, 并确保边界至少调整 inc 的 0.25 倍(默认不调整)。

  5. -Rxmin/xmax/ymin/ymax[+uunit] 在投影单位(如 UTM 米制)下指定区域。 其中 xmin/xmax/ymin/ymax 为与所选投影 (-J) 兼容的笛卡尔坐标。 unit 为允许的 单位-j 选项 (默认是 e)。 该选项会反算出实际的矩形地理区域。

    对于以 (0,0) 为中心的投影区域,可使用简写形式 -Rhalfwidth[/halfheight]+uunit, 其中 halfheight 默认为 halfwidth。此简写形式必须带 +u 修饰符。

  6. -Rjustifylon0/lat0/nx/ny 其中 justify 为两字符组合: L|C|R (左,中,右) and T|M|B (上,中,下),(例如, BL 为左下角) justify 指明 lon0/lat0 是矩形区域的哪个点, 而 nxny 与网格间距(通过 -I 设置)共同决定区域范围。 该形式常用于创建网格。 例如:-RCM25/25/50/50 表示一个以 (25,25) 为中心、尺寸为 50×50 的网格区域。

  7. -Rgridfile 从指定网格文件中复制区域范围设置。 根据调用模块的不同,此方式可能同时设置网格间距与网格配准方式 (参见 网格配准 )。

  8. -Ra[uto] 或 -Re[xact] 仅在现代模式下的绘图模块可用。自动从输入数据中确定区域:

    • -Re:精确匹配数据的范围(默认若未指定 -R)。

    • -Ra:在数据范围基础上略微扩大,使区域边界为数据范围的合理倍数。

-Tage

使用固定年龄进行模型计算(即忽略输入年龄网格中的年龄值)。 此选项可用于生成某一时刻的模型快照。 若未提供年龄网格,则必须指定该选项。

-V[level] (more …)

设置 verbose 等级 [w]

-bi[ncols][type][w][+l|b] (more …)

设置二进制输入数据的格式

-d[i|o]nodata (more …)

将输入数据中等于 nodata 的记录替换为 NaN,或将输出数据中值为 NaN 的记录替换为 nodata

-h[i|o][n][+c][+d][+msegheader][+rremark][+ttitle] (more …)

跳过或生成指定数目的头段记录

-ocols[,…][,t[word]] (more …)

设置输出数据列(0表示第一列,t 表示文本列)

-r[g|p] (more …)

设置网格配置方式 [默认为网格线配准]

-^-

显示简短的帮助信息,包括模块简介和基本语法信息(Windows下只能使用 -

-++

显示帮助信息,包括模块简介、基本语法以及模块特有选项的说明

-? 或无参数

显示完整的帮助信息,包括模块简介、基本语法以及所有选项的说明

--PAR=value

临时修改GMT参数的值,可重复多次使用。参数列表见 配置参数

大地坐标与地心坐标

所有球面旋转均作用于地心坐标。 这意味着输入的数据点和网格被视为采用大地坐标表示, 在进行旋转之前,必须先将其转换为地心坐标。 旋转操作完成后,最终重建得到的点将再转换回大地坐标。 如果将椭球设置参数 PROJ_ELLIPSOID 更改为 Sphere, 则可跳过此默认的转换过程。

点位于多边形内/外

为了确定点在多边形内,外或在边界上,GMT 会平衡数据类型以及多边形形状等因素来 确定算法。对于笛卡尔坐标,GMT 使用 non-zero winding 算法,该算法非常快。对于 地理坐标,如果多边形不包括两极点且多边形的经度范围不超过 360 度,同样使用该 算法。否则,GMT 会采用 full spherical ray-shooting 方法。

示例

我们将使用一个太平洋地壳年龄的网格 (pac_age.nc)、一个板块运动 模型 (Pac_APM.txt) 以及一个包含当前太平洋板块轮廓的多边形 (pac_clip_path.txt)。 要评估太平洋板块当前的板块运动方位角,请尝试:

gmt grdpmodeler pac_age.nc -EPac_APM.txt -V -Fpac_clip_path.txt \
            -Gpac_dir_0.nc -Sa -T0

要确定自地壳形成以来整个太平洋板块的纬度变化,请尝试:

gmt grdpmodeler pac_age.nc -EPac_APM.txt -V -Fpac_clip_path.txt \
     -Gpac_dlat.nc -Sy

要确定太平洋地壳形成时的有效板块运动速度,请尝试:

gmt grdpmodeler pac_age.nc -EPac_APM.txt -V -Fpac_clip_path.txt \
     -Gpac_vel.nc -Sv

要确定地壳自形成以来已经移动了多远,请尝试:

gmt grdpmodeler pac_age.nc -EPac_APM.txt -V -Fpac_clip_path.txt \
     -Gpac_dist.nc -Sd

要将地壳形成时的坐标保存到单独的网格中,请尝试:

gmt grdpmodeler pac_age.nc -EPac_APM.txt -V -Fpac_clip_path.txt \
     -Gpac_origin_%s.nc -SXY

要重复相同的练习,但将输出的 lon, lat, age, xorigin, yorigin 保存到一个表格中,请使用:

gmt grdpmodeler pac_age.nc -EPac_APM.txt -V -Fpac_clip_path.txt -SXY > origin.txt

备注

GMT 发布了 EarthByte 旋转模型 Global_EarthByte_230-0Ma_GK07_AREPS.rot。 如需使用替代的旋转文件,请创建一个名为 GPLATES_ROTATIONS 的环境变量,使其指向该替代的C旋转文件。

相关模块

backtracker, pmodeler, grdrotater, grdspotter, hotspotter, originater, rotconverter