Question

官方文档

mapproject

简介

地图投影正变换及逆变换

译文作者

@cugliming

mapproject 从表数据中读取（经度，纬度）位置，并使用指定的地图投影和比例尺计算（x，y）坐标。

也可以读取（x，y）位置并通过逆变换计算（经度，纬度）值。通过这种方式，我们可以将已知投影类型的地图数字化为线性（x，y）点，再将其变换为地理坐标。

也可以计算沿轨道到固定点的距离，或者最接近直线的距离。或用于执行各种基准转换。

读入的数据字段必须具有（经度，纬度）或（x，y）；关于如何读取（纬度，经度）文件，请参阅选项 -:。

mapproject还可以根据构成轨迹的输入坐标来计算各种辅助输出数据。方位角、距离、到其它直线的距离以及沿着直线的行进时间等，都可以通过使用 -A，-G，-L，-Z选项中的一个或多个来计算。

语法

gmt mapproject [ tables ]

-Jprojection (more …)

设置地图投影方式

-Rxmin/xmax/ymin/ymax[+r][+uunit] (more …)

指定数据范围

UTM投影比较特殊：如果使用

table

一个或多个ASCII或二进制表数据。若不提供表数据，则会从标准输入中读取。

-Ab|B|f|F|o|O[lon0/lat0][+v]

计算沿着某轨迹，或到某固定点（lon/lat）的方位角。-Af计算每个数据点的（前向）方位角，使用-Ab计算数据点到固定点的反方位角。使用-Ao可获得（-90/90）方向，而非方位角（0/360）。大写的F，B，O，表示从大地纬度转换为地心纬度，并估计测地线的方位角（假设当前的椭球体不是球体）。如果没有给出固定点，则从上一个点开始计算方位角（或反方位角）。或者，通过输入文件中的3-4列，利用+v选项以获得可变的第二个点（lon0 / lat0）。有关-A如何影响输出记录的信息，请参见输出顺序。

-C[dx/dy]

将投影坐标的中心设置为地图投影中心(默认为左下角)。也可以通过指定dx/dy，在投影结果中加上一定的偏移量(指定-I后，则减去该偏移量)，如特定投影类型中的东、北伪偏移量。偏移量的单位是有效的绘图距离单位（请参阅 PROJ_LENGTH_UNIT），但如果指定了-F，偏移量单位为米。

-Dc|i|p

临时覆盖 PROJ_LENGTH_UNIT，并改用c(厘米)，i(英寸)或p(像素)为单位。该选项不能与-F一起使用。

-E[datum]

将大地坐标（lon，lat，height）转换为以地球为中心的地球固定坐标（ECEF）（x，y，z）（逆向转换，添加-I）。

-F[unit]

强制使用1:1比例尺，即输出数据（或输入，参见-I）为真实投影下的米数。如果要指定其他单位，请附加所需的单位（请参阅units）。如果没有-F，则输出（或输入，请参见-I）以 PROJ_LENGTH_UNIT 指定的单位（或者采用-D）。

-G[lon0/lat0][+a][+i][+uunit][+v]

计算沿某轨迹的距离，或到某固定点（用-Glon0 / lat0指定）的距离。

在距离后加上+u以指定单位（参阅Units，以了解可用单位及如何计算距离[真实半径的大圆]）；例如，c表示使用输入坐标的笛卡尔距离，C使用投影坐标的笛卡尔距离。C单位需要设置-R和-J。如果没有给出固定点，将计算沿输入点定义的轨迹的累积距离[指定+a，也表示计算累积距离]; 附加+i可获取连续点之间的增量距离，或者同时附加+a +i，以同时获得累积距离及增量距离。或者，通过输入文件中的3-4列，附加+v以获得可变的第二个点（lon0 / lat0）。有关-G如何影响输出记录的信息，参见输出顺序。

-I

进行逆变换，即从（x，y）数据获取（经度，纬度）。

-Lline.xy[+uunit][+p]

判断从输入数据点到ASCII多段文件line.xy中指定的行的最短距离。

距离和最近点的坐标将作为三个新列添加到输出中。追加距离单位（参阅Unites以获取可用单位以及如何计算距离[使用真实半径的大圆]），例如c(使用输入坐标的笛卡尔距离)，C(使用投影坐标的笛卡尔距离)。C单位需要设置-R和-J。附加+p以报告线段ID和小数点编号，而不是最近点的lon/lat。 有关-L如何影响输出记录的信息，请参见输出顺序。

-N[a|c|g|m]

从大地纬度（使用当前的椭球，参阅 PROJ_ELLIPSOID）转换为四个不同辅助纬度中的一个，经度不受影响。从authalic（等积？），conformal（等角？），geocentric（地心坐标系？）和meridional latitudes（子午线纬度？）[地心坐标系]中选择。使用-I将辅助纬度转换为大地纬度。

-Q[d|e]

列出所有投影参数。

-Qd，仅列出基准；-Qe仅列出椭圆体。

-S

不转换区域之外的点。

-T[h]from[/to]

使用标准的Molodensky转换对datums(基准面？)之间进行坐标转换。如果第3个输入列的高度高于椭圆形，则使用-Th [默认假定height = 0，即在椭圆形上]。使用基准面的ID（参见-Qd）指定基准或给出椭球：dx，dy，dz，其中椭球可以是椭球ID（参见-Qe），也可以指定为a[，inv_f]，其中a是半长轴，inv_f是反展平（如果省略，则为0）。

-V[level] (more …)

设置 verbose 等级 [w]

-W[g|h|j|n|w|x]

在标准输出上打印地图的宽度和高度。 没有读取输入文件。要仅输出宽度或高度，请分别附加w或h。要输出地图点的绘图坐标，请输入glon /lat。报告的绘图尺寸的单位可以通过-D进行更改。要输出参考点的地图坐标，请选择jcode（使用标准的两个字符对齐代码），nrx / ry（其中参考点以0-1范围内的标准化位置给出）或xpx / py，直接给出点。[默认值返回地图的宽度和高度]。

-Z[speed][+a][+i][+f][+tepoch]

计算沿-G指定轨迹的行进时间；附加一个恒定的单位，如果缺失，则会从第3列读取。

速度为单位时间(单位由 TIME_UNIT [m/s]指定)内通过的距离（单位由-G指定）。附加+i可以输出连续点之间的增量行进时间，+a可以获取累积的行进时间，或者同时加上以获取两种时间信息。使用+f将累积的（经过的）旅行时间格式化为ISO 8601规定格式。至于用于表示秒的小数位数，请查阅 FORMAT_CLOCK_OUT设置。附加+tepoch以报告连续点的绝对时间（ETA）。最后，由于需要增加距离，因此需要带有+i修饰符的-G选项。有关-Z如何影响输出记录的信息，请参见输出顺序。

-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 表示文本列）

-je|f|g (more …)

设置球面距离的计算方式

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

设置输出数据列（0表示第一列，t 表示文本列）

-p[x|y|z]azim[/elev[/zlevel]][+wlon0/lat0[/z0]][+vx0/y0] (more …)

设置3D透视视角

-s[cols][+a|+r] (more …)

设置 NaN 记录的处理方式

-:[i|o] (more …)

交换输入或输出中的第一和第二列

-^ 或 -

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

-+ 或 +

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

-? 或无参数

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

--PAR=value

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

输出顺序

选项 -A、-G、-L 和 -Z 生成的数据一律按这几个选项的字母顺序排序。因此，这些选项在命令中出现的顺序无关紧要。当然，输出顺序也可以用 -o 选项调整。

示例

将经纬度坐标数据转换为Mercator网格上以cm为单位的（x,y）

本示例中，Mercator每度表示为2cm:

$ gmt mapproject  -R105/108/33/36 -Jm2c << EOF
105 33
106 34
107 35
108 36
EOF

输出结果:

0       0
2       2.38727963335
4       4.80309472214
6       7.24888693724

UTM转化为地理坐标

• UTM坐标单位为m

• 需要知道该坐标的UTM分带号，本示例坐标所在带号为51

# 将UTM坐标转换为经纬度
$ echo 311072.4 3542183.5 | gmt mapproject -Ju51/1:1 -C -F -I
# 输出结果：
121.000000401   32.0000000858

# 将Gauss-Kruger坐标转换为经纬度
$ echo 310996.760635 3543600.93152 | gmt mapproject -Ju51/1:1 -C -F -I --PROJ_SCALE_FACTOR=1
# 输出结果：
121 31.9999999997

注意： Gauss-Kruger坐标可直接使用UTM的-Ju投影，曲别在于Gauss-Kruger的比例因子为1，而UTM比例因子为0.9996.

经纬度转化为UTM

# 经纬度转换为Gauss-Kruger
$ echo 121 32 | gmt mapproject -Ju51/1:1 -C -F --PROJ_SCALE_FACTOR=1
# 输出结果
310996.760635 3543600.93152

# 经纬度转换为UTM
$ echo 121 32 | gmt mapproject -Ju51/1:1 -C -F
# 输出结果
311072.361931 3542183.49115

在经纬度坐标下，计算点 (104.92, 29) 相对于点 (104.92, 31) 的方位角和大圆距离（单位为 km）。下面的例子使用了 -je 选项以获得高精度的测地线距离:

$ echo 104.92 29 | gmt mapproject -Af104.92/31 -G104.92/31+uk -je
104.92  29      180     221.704999623