grdcut・grdgradient・grdimageの使い方

grdcut — GMT 6.2.0 documentation
grdgradient — GMT 6.2.0 documentation
grdimage — GMT 6.2.0 documentation

grdcut・grdgradient・grdimage はDEMを使用する時によく使うモジュールです。
DEM3兄弟 とか呼んであげてください。

DEMとは

DEM(数値標高モデル) | 国土地理院

GMTではオンラインでDEMを取得することができます。
取得できるものは以下の通り。

Remote Data Sets — GMT 6.2.0 documentation
CodeDimensionsRegSizeDescription
01d360 x 180g,p128 KB1 arc degree global relief (SRTM15+V2.1 @ 111 km)
30m720 x 360g,p435 KB30 arc minute global relief (SRTM15+V2.1 @ 55 km)
20m1080 x 540g,p918 KB20 arc minute global relief (SRTM15+V2.1 @ 37 km)
15m1440 x 720g,p1.6 MB15 arc minute global relief (SRTM15+V2.1 @ 28 km)
10m2160 x 1080g,p3.4 MB10 arc minute global relief (SRTM15+V2.1 @ 18 km)
06m3600 x 1800g,p8.8 MB6 arc minute global relief (SRTM15+V2.1 @ 10 km)
05m4320 x 2160g,p13 MB5 arc minute global relief (SRTM15+V2.1 @ 9 km)
04m5400 x 2700g,p19 MB4 arc minute global relief (SRTM15+V2.1 @ 7.5 km)
03m7200 x 3600g,p33 MB3 arc minute global relief (SRTM15+V2.1 @ 5.6 km)
02m10800 x 5400g,p71 MB2 arc minute global relief (SRTM15+V2.1 @ 3.7 km)
01m21600 x 10800g,p258 MB1 arc minute global relief (SRTM15+V2.1 @ 1.9 km)
30s43200 x 21600g,p935 MB30 arc second global relief (SRTM15+V2.1 @ 1.0 km)
15s86400 x 43200p3.2 GB15 arc second global relief (SRTM15+V2.1)
03s432000 x 216000g6.8 GB3 arc second global relief (SRTM3S)
01s1296000 x 432000g41 GB1 arc second global relief (SRTM1S)

おすすめは、15s 03s 01s
かなり広い範囲を描画する場合は30s 01mも検討してみてください。

GMT上では @earthrelief【code】と指定してください。今回は15sを選択します。

DEMを切り出す (grdcut)

gmt grdcut [使用するDEM] -R[切り出し範囲] -G[出力ファイル名]

DEMを使うためには、まず自分が描きたい範囲のDEMを切り出す必要があります。

たとえば、

#!/bin/bash
gmt grdcut @earth_relief_15s -R120/150/20/50 -Gdem.nc

これを行うと、dem.ncというDEMが出力されます。

DEMの勾配計算 (grdgradient)

gmt grdgradient [使用するDEM] -G[出力ファイル名] -A[光を当てる角度] -N[正規化値] -G[出力ファイル名]

かっこいい地図を描く場合は勾配計算を行いましょう。
光を当てる角度は45度。正規化値はe0.8がおすすめです。

#!/bin/bash
gmt grdcut @earth_relief_15s -R120/150/20/50 -Gdem.nc
gmt grdgradient dem.nc -Ggrad.grd -A45 -Ne0.8

これを行うと、grad.grdという勾配ファイルが出力されます。

DEMを描画 (grdimage)

gmt grdimage [DEMファイル] -I[勾配ファイル] -C

最後にDEMを地図に描画します。

makecptでDEM用のカラーパレットを作成することを忘れずに。
geo とか globe がおすすめ。

#!/bin/bash
gmt begin hogehoge jpg
    gmt basemap -JM12 -R120/150/20/50 -Bafg -BWSNE
    gmt coast -Df -W0.25 -Lx11/-1.5+c20+w1000+f
    gmt makecpt -Cgeo -T-8000/8000/200 -Z
    gmt grdcut @earth_relief_15s -R120/150/20/50 -Gdem.nc
    gmt grdgradient dem.nc -Ggrad.grd -A45 -Ne0.8
    gmt grdimage dem.nc -Igrad.grd -C
gmt end
rm *.nc *.grd

この結果が、

となる。

コメント

タイトルとURLをコピーしました