✍️ 周茂 • 📅 2022-06-21

# gravprisms

官方文档:

  gravprisms

简介:

  计算三维垂直棱柱产生的位异常

gravprisms 计算矩形棱柱产生的位场。输入文件/标准输入最多包涵 7 列： 前 4 列为 x y z_low z_high ，即棱柱中心的 x 和 y 坐标以及垂直方向的 高度；接下来两列为每个棱柱的增量 dx 和 dy ，见 `-E` ；最后一列为 棱柱的密度，如果 `-D` 指定了密度，则被固定密度覆盖。或者，也可以使用 `-C` 选项来创建棱柱。如果单个棱柱密度可变（ `-H` ），则每个棱柱被分解为 密度恒定，堆叠的子棱柱，垂直增量为 dz 。该模块可以计算等距网格上的异常 （通过 `-R` 和 `-I` 选项指定）或通过 `-N` 指定任意的计算点。可计算的异常 包括自由空气重力异常、垂直重力梯度异常或大地水准面异常。

`Source Code`

上图展示了三种不同密度的海山，海山都由棱柱构成。左图为常密度（`-D`）； 中图为垂向等密度（`-D`）；右图为随 r 和 z 变化的密度（`-H`）。

## 语法

gmt gravprisms [ table ] [ `-A` ] [ `-C`[+q][+w file][+z dz] ] [ `-D` density ] [ `-E` dx[/dy] ] [ `-F` f|n[lat]|v ] [ `-G` outfile ] [ `-H` H/rho_l/rho_h[+d densify][+p power] ] [ `-I` increment ] [ `-L` base ] [ `-M`[h][v] ] [ `-N` trackfile ] [ `-R` region ] [ `-S` shapegrid ] [ `-T` top ] [ `-V`[level] ] [ `-W` avedens ] [ `-Z` level|obsgrid ] [ `-bo` binary ] [ `-d` nodata[+c col] ] [ `-e` regexp ] [ `-f` flags ] [ `-h` headers ] [ `-i` flags ] [ `-o` flags ] [ `-r` reg ] [ `-x`[[-]n] ] [ --PAR=value ]

## 输入数据

table

  棱柱文件，格式为 x y z z_low z_high [dx dy] [rho] 。其中的可选项可以 通过 `-E` 和 `-D` 控制。密度的单位可以是 kg/m^3 或者 g/cm^3 。如果使用 `-C` 选项，则无需输入文件。

## 可选选项

-A

-A

  设置 z 轴向上为正 [默认向下为正]

-C

-C[+q][+w file][+z dz] ]

  与其他选项组合创建棱柱，1）与 `-S` height 组合创建棱柱；2）与 `-L` base 和 `-T` top 指定的两个面创建棱柱。对于第二种情况， base 和 top 可以为常数或者网格。对于第一种情况，则创建从 0 到 height 的棱柱。如果使用 `-H` 设置可变密度，则将每个棱柱分为多个堆叠 的具有单独密度的子棱柱，因此需要为这些子棱柱设置高度增量 dz ，第一个 和最后一个棱柱可能会对高度进行适当的调节。如果不使用 `-H` 选项，则只 创建统一密度的棱柱，通过 `-D` 可以使每个棱柱具有不同的密度。 +w file 将会将棱柱保存到文件中，若追加 +q 则保存结果后直接 退出，不进行位计算。

-D

-D density

  设置一个固定的密度差，可以覆盖输入文件中的密度设置，单位可以为 kg/m^3 或者 g/cm^3。或者可以指定一个网格，其中网格的每个格点的值为棱柱的密度， 此时，需要 `-C` 选项，并必须和 `-S` （或 `-L` 和 `-T` ）提供的网格具有 相同的大小和配准方式。

-E

-E dx[/dy]

  如果输入文件中所有棱柱都具有相同的 x/y 方向的尺度，则可以在这里设置其 增量。在这种情况下，输入文件值包含棱柱的中心坐标和高度以及可选的密度。 如果 `-E` 选项只给出 dx ，则 dy = dx 。对于地理坐标，dx 是经度 增量，在数值上保持不变，但实际宽度会随着纬度的增大而减小。

-F

-F f|n[lat]|v

  指定计算的位异常。

  - f 自由空气重力异常 [默认]

  - n 大地水准面，可指定一个平均纬度，默认使用网格或者 `-N` 指定的 纬度的中值

  - v 垂直重力梯度

-G

-G outgrid[=ID][+d divisor][+n invalid][+o offset|a][+s scale|a][:driver[dataType][+c options]]

> 输出网格文件名。通过追加 =ID 可指定 网格格式。 参数详细介绍请参考 写 netCDF 文件。

-H

-H H/rho_l/rho_h[+d densify][+p power]

  为海山设置随深度变化的密度。给定低密度和高密度分别为 rho_l 和 rho_h ，以及固定的参考高度 H 。使用 +d 和 +p 分别设置水压驱动的侧翼 密度的增量 [0] 以及幂指数变化 [1，线性变化]。该选项需要 `-S` 选项。细节 请参考 grdseamount 。

-I

-I xinc[+e|n][/yinc[+e|n]]

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

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

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

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

  注意：

  - 若 yinc 设置为0，则表示其与 xinc 相同

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

-L

-L base

  设置棱柱基底/下界面，可以为一个网格，也可以为一个常数 [0]

-M

-M[h][v]

  设置距离单位。h 表示水平方向单位为 km [m]，z 表示垂直方向距离 单位为 km [m]。

-N

-N trackfile

  指定计算异常所处的测线，结果输出到标准输出。如果没设置 `-Z` 选项，则 trackfile 必须有 3 列，最后一列为 z 值，否则使用 `-Z` 设置的恒定的 z 值。

-R

-R xmin/xmax/ymin/ymax[+r][+u unit]

  指定数据范围。 (参数详细介绍)

-S

-S height

  指定完整的海山网格名称，用于创建棱柱或者 `-H` 所需

-T

-T top

  设置棱柱的上界面，可以为一个网格，也可以为常数

-W

-W avedens

  指定由 `-C` 和 `-H` 创建的随空间变化的密度网格名，作为输出

-Z

-Z level|obsgrid

  指定观测高度，可以为常数，也可以为可变的值，若为可变的值，则通过网格给定， 该网格即确定了输出的网格区域

-V

-V[level]

  设置 verbose 等级 [w]。 (参数详细介绍)

-bo

-bo[ncols][type][w][+l|b]

  控制二进制文件的输出格式。 (参数详细介绍)

-d

-d[i|o]nodata

  将某些特定值当作 NaN。 (参数详细介绍)

-e

-e[~]"pattern" | -e[~]/regexp/[i]

  筛选或剔除匹配指定模式的数据记录。 (参数详细介绍)

-f

-f[i|o]colinfo

  显式指定当前输入或输出数据中每一列的数据类型。 (参数详细介绍)

-h

-h[i|o][n][+c][+d][+m segheader][+r remark][+t title]

  在读/写数据时跳过文件开头的若干个记录。 (参数详细介绍)

-i

-i cols[+l][+s scale][+o offset][,...][,t[word]]

  对输入的数据进行列选择以及简单的代数运算。 (参数详细介绍)

-o

-o cols[,...][,t[word]]

  对输出的数据进行列选择以及简单的代数运算。 (参数详细介绍)

-r

-r[g|p]

  设置网格配置方式 [默认为网格线配准]。 (参数详细介绍)

-x

-x[[-]n]

  限制多核算法中能使用的核数（需要GMT开启OpenMP支持）。 (参数详细介绍)

-:

-:[i|o]

  交换输入或输出数据的前两列。 (参数详细介绍)

-^ 或 -

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

-+ 或 +

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

-? 或无参数

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

--PAR=value

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

## 距离单位

GMT支持多种不同的距离单位，以及三种不同的球面距离计算方式。 详情见 单位 和 -j 选项。

## 示例

假定某棱柱构成的高斯海山高度为 6000 m，底部半径为 30 km，总共包括 2828 个棱柱， 通过 plot3d 可绘制该海山

    gmt plot3d -R-30/30/-30/30/0/7000 -JX12c -JZ3c -Ggray -So1q+b @prisms.txt -B \
        -Wfaint -p200/20 -pdf smt

计算该海山在给定网格上产生的重力异常，密度差异设置为 1700 kg/m^3，水平距离以 km 为单位，垂直距离以 m 为单位，观测高度为 7000 m

    gmt gravprisms -R-40/40/-40/40 -I1 -Mh -G3dgrav.nc @prisms.txt -D1700 -Ff -Z7000
    gmt grdimage 3dgrav.nc -B -pdf 3dgrav

计算该海山在 `crossing.txt` 中位置产生的垂直重力梯度异常

    gmt math -T-30/30/0.1 T 0 MUL = crossing.txt
    gmt gravprisms -Ncrossing.txt -Mh @prisms.txt -D1700 -Fv -Z7000 > vgg_crossing.txt
    gmt plot vgg_crossing.txt -R-30/30/-50/400 -i0,3 -W1p -B -pdf vgg_crossing

重新计算重力，但使用棱柱文件中的密度

    gmt gravprisms -Ncrossing.txt -Mh @prisms.txt -Ff -Z7000 > faa_crossing.txt
    gmt plot faa_crossing.txt -R-30/30/0/350 -i0,3 -W1p -B -pdf faa_crossing

## 注意事项

垂直棱柱造成的大地水准面的解析表达式（Nagy et al., 2000）比较复杂，包含 48 项。 由于各种抵消，最终结果比重力异常和重力异常垂直梯度更不稳定，因此，最终结果的 有效数字可能比较少。

## 参考文献

Grant, F. S. and West, G. F., 1965, Interpretation Theory in Applied Geophysics, 583 pp., McGraw-Hill.

Kim, S.-S., and P. Wessel, 2016, New analytic solutions for modeling vertical gravity gradient anomalies, Geochem. Geophys. Geosyst., 17, https://doi.org/10.1002/2016GC006263.

Nagy D., Papp G., Benedek J., 2000, The gravitational potential and its derivatives for the prism, J. Geod., 74, 552–560, https://doi.org/10.1007/s001900000116.

## 相关模块

grdmath, gravfft, grdseamount, gravmag3d, grdgravmag3d, talwani3d, talwani2d