:author: 周茂
:date: 2026-01-04

.. index:: ! img2grd
.. program:: img2grd

img2grd
========

:官方文档: :doc:`gmt:supplements/img/img2grd`
:简介: 从墨卡托 img 格式文件中提取网格数据

**img2grd** 读取 img 格式的文件，提取数据并写入到网格文件中。:option:`-M` 选项用来设置
在提取的网格文件中保留球墨卡托投影，如果不使用该选项，则取消保留球墨卡托投影，结果为
地理坐标网格。


img 格式是一种二进制网格文件，数据为球墨卡托投影后的结果。球墨卡托投影原理决定其纬度间隔在转换为地理坐标
时不是等间距的，因此导致了下述两个坐标系统下的区域范围轻微的差异。球墨卡托投影原理还决定了 img 文件
在纬度方向并不是全球覆盖的，在极区有空白，即 :option:`-D` 选项设置的最大纬度 80.738 ；
在经度方向则是全覆盖的。虽然该模块可以处理任何 img
格式文件，但目前较少的产品或模型以 img 格式发布。Sandwell 和 Smith 基于卫星测高数据创建的
海洋重力场（包括重力异常、垂线偏差和垂直重力梯度）和海底地形模型采用了该格式发布。
由于测高卫星轨道倾角的限制，在不同时期，
上述模型的纬度范围不同，早期最大纬度为 72.006，当前最大纬度范围为 80.738 。

Sandwell 和 Smith 模型的大小和文件名有一定规律，用于在 GMT 中自动设置参数并处理，文件名在不同时期的
命名方式不同，但规律类似，例如：world_grav.img.7.2 和 grav_28.1.img

- grav 表示该文件为自由空气重力异常
- img 表示文件格式为 img
- 7.2 表示该模型版本号为 7，分辨率为 2 分；28.1 表示版本号为 28，分辨率为 1 分

语法
----

**gmt img2grd**
*imgfile*
:option:`-G`\ *grdfile*
:option:`-R`\ *region*
[ :option:`-D`\ [ *minlat/maxlat* ] ]
[ :option:`-E` ]
[ :option:`-F` ]
[ :option:`-I`\ *inc* ]
[ :option:`-M` ]
[ :option:`-N`\ *navg* ]
[ :option:`-S`\ [*scale*] ]
[ :option:`-T`\ *type* ]
[ :option:`-V`\ [*level*] ]
[ :option:`-W`\ *maxlon* ]
[ :option:`-n`\ *flags* ]
[ :doc:`--PAR=value </conf/overview>` ]

输入数据
----------------

*imgfile*
    输入墨卡托 img 格式文件名，例如： Sandwell 和 Smith 创建的海洋重力场
    或海底地形文件。

必须选项
------------------

.. include:: explain_grd_out.rst_

.. include:: explain_-R.rst_

可选选项
--------

.. option:: -D

**-D**\ [*minlat/maxlat*]
    设置纬度范围为 -80.738/+80.738。附加 *minlat/maxlat* 可以指定输入的 img 文件
    的纬度范围 [ 默认为 -72.006/72.006 ]。该选项通常不必使用，因为本模块可以自动确定
    输入文件的纬度范围。

.. option:: -E

**-E**
    当输出网格为地理坐标时，即未设置 :option:`-M` 选项时，使用 :option:`-E` 选项可以强制地理坐标格网符
    合纬度范围并强制经纬度间隔相等。默认情况下，本模块首先从输入文件中提取 :option:`-R` 指定范围内的
    网格，网格的坐标形式和原始文件相同，即球墨卡托投影。因此，当输出网格为地理坐标时，需将球墨
    卡托投影转换为地理坐标。在转换过程中，通常会稍微超出指定的纬度范围，且在纬度方向的间隔和经
    度方向不匹配。但是，这种强制经纬度范围等间隔的重采样会导致内插误差。因此，只有强制输出网格
    为指定范围内，且经纬度间隔相同时才建议使用。在这种情况下，:option:`-R` 设置的区域必须为经纬度间
    隔的整数倍。

.. option:: -F

**-F**
    平移 x 和 y 墨卡托坐标，使其位于左下角，即 (0,0)

.. option:: -I

**-I**\ *inc*
    指定输入的 img 文件纬度间隔 *inc* [ 默认为 2 ]。添加 **m** [ 默认 ] 或 **s**
    来指定单位。该选项通常不必使用，因为本模块可以根据输入 img 文件大小自动判断。

.. option:: -M

**-M**
    输出球墨卡托网格 [ 默认为地理 lon/lat 网格 ]。输入文件的的球墨卡托投影被保留，因此 :option:`-R` 指定的区域
    被轻微地调整，调整后的区域对应于像素的边界。设置输出网格文件的头部信息，x 和 y 轴方向的距离相对于使用 **-Jm1**
    投影下的投影中心。使用 :option:`-F` 选项可以将 x 和 y 轴的距离设置为相对于输出范围 :option:`-R` 的西南角，距离单位为用户默认
    的单位。

.. option:: -N

**-N**\ *navg*
    对输入的墨卡托像素在 *navg* x *navg* 的范围内平均，并创建平均后的像素点。
    如果和 **-T3** 同时使用，输出网格值为平均约束值，范围为 0 到 1。
    如果和 **-T2** 同时使用，输出网格值为像素平均值或 NaN，平均约束值 > 0.5 时，为像素平均值，
    平均约束值 < 0.5 时，为 NaN。
    *navg* 的单位为像素 [ 默认值为 1，表示不进行平均]。

.. option:: -S

**-S**\ [*scale*]
    对输出的网格文件进行 *scale* 比例缩放，以保证结果单位正确[ 默认是 1.0 ]。
    对于目前的 img 文件：
    海底地形，单位为米，因此设置为 :option:`-S`\ 1 得到米 ；
    自由空气重力异常，单位为 mGal\*10 ，因此设置为 :option:`-S`\ 0.1 得到 mGal ；
    垂线偏差，单位为 micro-radians\*10 ， 因此设置为 :option:`-S`\ 0.1 得到 micro-radians ；
    垂直重力梯度，单位为 Eotvos\*10 ， 因此设置为 :option:`-S`\ 0.1 得到 Eotvos 或者设置为 :option:`-S`\ 0.01 得到 mGal/km。
    如果不给定 *scale*，则将通过文件名来自动确定。

.. option:: -T

**-T**\ *type*
    *type* 处理约束信息的编码（译注：约束信息可理解为数据可靠性的一种衡量标准）。
    *type* = 0 表明没有约束信息编码到 img 文件（1995 年前的重力场模型）中，并获取所有数据。
    *type* > 0 表明有约束信息编码到 img 文件 (1995 以后以及当前版本的 img 文件)。
    **-T1** 获取所有数据；
    **-T2** 获得约束点的数据值，内插点的数值设置为 NaN ；
    **-T3** 约束点的数据值设置为 1，内插点设置为 0 [ 默认值为 1]。

.. include:: explain_-V.rst_

.. option:: -W

**-W**\ *maxlon*
    指定输入 img 文件的最大经度 *maxlon*。从 1995 年开始，Sandwell 和 Smith 的模型的
    *maxlon* = 360.0 ，但是某些早期的模型 *maxlon* = 390.0 [ 默认为 360.0 ]。

.. include:: explain_-n.rst_

.. include:: explain_help.rst_

地理坐标示例
------------

如果需要输出地理坐标网格，则不能使用 :option:`-M` 选项。从 *grav_28.1.img* 文件中提取数据并转换为地理坐标::

    gmt img2grd grav_28.1.img -Gmerc_grav.nc -R-40/40/-70/-30 -V

由于在 img 文件的墨卡托单位中，纬度间隔是相等的，转换为地理坐标时导致纬度间隔不相等，因此最终提取的网格并
不能严格匹配 :option:`-R` 选项指定的纬度范围。使用 :option:`-E` 选项可以精确的匹配 :option:`-R` 指定的范围，并且强制经度
和纬度的间隔相同::

    gmt img2grd grav_28.1.img -Gmerc_grav.nc -R-40/40/-70/-30 -E -V


墨卡托示例
----------

**数据处理示例**

由于 img 格式为墨卡托投影，因此如果最终数据处理结果为墨卡托地图，不建议使用该模块提取地理坐标，而是直接使
用墨卡托投影坐标处理。如果强行从 img 提取到地理坐标再转换为墨卡托投影会导致丢失短波信息。因此，最好使用
:option:`-M` 选项并与所需的墨卡托投影设置一定的线性比例缩放（见 GMT 例 29）。从 *world_grav.img.7.2* 中
提取 :option:`-R`\ -40/40/-70/-30 范围内的数据::

    gmt img2grd -M world_grav.img.7.2 -Gmerc_grav.nc -R-40/40/-70/-30 -V

注意到 :option:`-V` 选项告诉用户区域范围被调整为 :option:`-R`\ -40/40/-70.0004681551/-29.9945810754 。如果在脚本
编写中需要使用该范围，可以使用 **grdinfo -Ii**。此外，用户也可以从 **grdinfo** 的信息中看到，网格文件头
部信息显示区域范围调整为 :option:`-R`\ -40/40/-99.4333333333/-31.4666666667。该纬度范围为使用 **-Jm**\ 1 和
:option:`-R`\ -40/40/-70.0004681551/-29.9945810754 经过球墨卡托投影的结果。使用 ship.lonlatgrav 对 merc_grav.nc
进行采样的代码为::

    gmt set PROJ_ELLIPSOID Sphere

    gmt mapproject -R-40/40/-70.0004681551/-29.9945810754 -Jm1i -C ship.lonlatgrav | \
    gmt grdtrack -Gmerc_grav.nc | \
    gmt mapproject -R-40/40/-70.0004681551/-29.9945810754 -Jm1i -I -C > ship.lonlatgravsat

如上述代码，建议对数据投影和重投影，而不是对 img 文件进行类似的操作。因为对数据处理只进行一次插值（即 grdtrack）,
对 img 文件操作会进行两次插值（分别在转换和重采样中）。

如果想从上面的墨卡托网格中提取地理坐标网格，可以使用::

    gmt grdproject merc_grav.nc -R-40/40/-70.0004681551/-29.9945810754 -Jm1i -I -D2m -Ggrav.nc

在某些情况下，上述操作中的两个坐标的范围不能很好的对齐。这种情况下，可以使用下面的代码::

    gmt grd2xyz merc_grav.nc | \
    gmt mapproject -R-40/40/-70.0004681551/-29.994581075 -Jm1i -I | \
    gmt surface -R-40/40/-70/70 -I2m -Ggrav.nc

**绘图示例**

如果要获取上述区域的墨卡托地图，假定用户的 gmt.conf 中的 **PROJ_LENGTH_UNIT** 是 inch 。因为上面的 merc_grav.nc
是使用 **-Jm**\ 1i 投影得到的，1度对应地图 1 英寸，因此为 80 英寸宽。可以通过使用 **-Jx**\ 0.1i 投影来制作
8 英寸宽的地图以在其他程序中应用 （例如：grdcontour，grdimage，grdview），如果要叠加地理坐标图，可以使用上面
的 :option:`-R` 和 **-Jm**\ 0.1 来使两个坐标系统匹配，避免导致移位。

除了上述的方法，还有更加方便的方法。注意到输入的 img 文件的像素宽为 2 分（可以通过 ``grdinfo mercgrav.nc`` 查看）,
merc_grav.nc 有 2400 x 2039 个像素。如果图的宽度为 8 英寸，每英寸中包含 300 个像素。但实际绘图中不需要这么高的分辨率，
每英寸 100 个像素即可满足，因此可以把原始数据平均为 3 x 3 的像素。（如果想绘制等值线图的话，可能需要更多的平均保证等值线
光滑，例如：6 x 6）由于 2039 不能被 3 整除，因此实际的 :option:`-R` 会被调整::

    gmt img2grd -M world_grav.img.7.2 -Gmerc_grav_2.nc -R-40/40/-70/-30 -N3 -V

调整以后的区域为 :option:`-R`\ -40/40/-70.023256525/-29.9368261101 ，输出结果为 800 x 601 个像素。下面可以生成一个可以使用
:doc:`grdgradient` 生成一个光照文件::

    gmt grdgradient merc_grav_2.nc -Gillum.nc -A0/270 -Ne0.6

如果已经有一个 "grav.cpt" 文件，可以使用下面的命令生成阴影图::

    gmt begin
        gmt grdimage merc_grav_2.nc -Iillum.nc -Cgrav.cpt -Jx0.1i
        gmt basemap -R-40/40/-70.023256525/-29.9368261101 -Jm0.1i -Ba10
    gmt end show

如果要从 img 文件中只获取地理坐标的约束数据值，则需要同时使用 :option:`-T`\ 2 选项，使用 grd2xyz 输出数据值，通过管道和 ``grep -v NaN``
可以消除 NaN 值。


相关模块
--------

:doc:`img2google`
