Question

示例贡献者

刘珠妹（作者）、陈箫翰（修订）

视频讲解

https://www.bilibili.com/video/BV1Xh411e7fi/

---

绘制三维地形图的核心是 grdview 模块以及 -JZ 选项。使用 -JZ 选项设置 Z 方向上的缩放比例，同时也必须在 -R 选项中设置 Z 方向上的最小值和最大值。此外，-p 选项用于设置三维观察的视角。

本例将展示如何绘制三维起伏地形图，并叠加卫星图像图层与地震断层图层，用于对比研究。请点击观看视频学习代码的详细讲解。

示例中使用到的数据包括：

• @earth_relief_06m：earth_relief: 全球地形起伏数据

• @earth_day_06m：earth_day 和 earth_night: 地球昼夜卫星影像

• CN-faults.gmt：CN-faults: 中国断层数据

下面的绘图脚本将三维地形、卫星图片以及断层线绘制在不同高度：

#!/usr/bin/env bash
# 将三维地形、卫星图片以及断层线绘制在不同高度

gmt begin 3D_earth_relief_1
    gmt set MAP_FRAME_TYPE plain
    # 1. 绘制地形图形，图层默认位于 z 轴底部（此例中即 -1000 的高度）
    gmt grdview @earth_relief_06m -R101/108/35/42/-1000/20000 -JM10c -JZ10c -N-1000+ggray -Qi -I -Ba -BwsENZ -p60/25    
    gmt colorbar -C -Ba -DJTC+o0/1c -p

    # 2. 绘制断层线图层，-p 设置图层抬高到10000 
    gmt plot CN-faults.gmt -W2p,brown -p60/25/10000 -Ba -BWSen
    gmt basemap -TdjLT+w1.5c+l+o1c -p
    
    # 3. 绘制卫星图片图层，-p 设置图层抬高到20000
    gmt grdimage @earth_day_06m -p60/25/20000 -Ba -BWSen
gmt end show

下面的脚本将卫星图片和断层线直接绘制在三维地形上:

#!/usr/bin/env bash
# 将卫星图片和断层线直接绘制在三维地形上

gmt begin 3D_earth_relief_2
    gmt set MAP_FRAME_TYPE plain
    # 预处理1：生成断层的三维采样，输出数据为断层的经度、纬度以及对应的地形高程
    gmt grdtrack -R101/108/35/42 -G@earth_relief_06m CN-faults.gmt > faults.xyz

    # 预处理2：截取区域卫星图片
    # 下载 6 弧分的卫星图片数据 earth_day_06m_p 到当前目录
    gmt which -Gl @earth_day_06m_p
    # 使用 gdal_translate 截取区域的卫星图片
    # -projwin 后的四个参数指定了区域范围：左上角经度 左上角纬度 右下角经度 右下角纬度
    # gdal_translate 的详细用法见 https://gdal.org/programs/gdal_translate.html
    gdal_translate -of GTIFF -projwin 101 42 108 35 earth_day_06m_p.tif day.tif

    # 绘制卫星图片，-G 表示在 DEM 上贴卫星图片
    gmt grdview @earth_relief_06m -R101/108/35/42/-1000/5000 -JM10c -JZ4c -N-1000+ggray -Gday.tif -Qi -Ba -Bz -BwsENZ -p60/25   
    # 绘制断层线
    gmt plot3d faults.xyz -W1p,brown -p
    # 绘制色标
    gmt colorbar -C -Ba -DJTC+o0/1c -p
    # 绘制底图（含指北针）
    gmt basemap -TdjLT+w1.5c+l+o1c -p60/25/2000

    rm day.tif faults.xyz earth_day_06m_p.tif
gmt end show