個人的ノート集トップ
PyGMT
第12章　グラフ紙面上の処理

第12章　グラフ紙面上の処理

ここではグラフ内の座標系ではなく，紙面上での図の処理について扱う．

描画位置の原点の移動

プロット時の原点の位置を移動する． PyGMTの図は左下の角が原点となるようにされる．

fig.shift_origin(
    xshift = <xshift>,
    yshift = <yshift>,
)

それぞれ右，上が正． 例えば"3c"とすれば3cmプロット時の原点を移動する．

地図の紙面上サイズの取得

例えばprojection=“M6c"のようにしたとき， 地図のサイズは横方向で指定され，縦方向は自動的に決定される．

このようなときに地図のサイズが取得できると， メインの地図の横や上下に断面図などを追加する際に， 位置やサイズの調整のために便利である．

GMTのmapprojectメソッドで地図のサイズを取得できるが， PyGMTではmapprojectがまだ未実装なようである．

したがって手順は増えるが下記のようにして地図サイズを取得する．

参考: https://docs.generic-mapping-tools.org/dev/mapproject.html 
(figを先に作成しておく)
with pygmt.clib.Session() as session:
    with pygmt.helpers.GMTTempFile() as tmp:
        session.call_module("mapproject",f"-W<オプション> ->{tmp.name}")
        size = tmp.read()

このコードで得られる「size」はcm単位の文字列． -Wに続くオプションによって，縦横いずれのサイズを出力するかを指定する：

  • h: 高さのみ出力
  • w: 横幅のみ出力
  • なし: 高さ，横幅の両方を出力

これで得られた画像サイズを shift_origin(->第*※相互参照無効*節)と組み合わせれば， 適切な位置とサイズで図を追加できる．

例）深さ断面図の追加

import pygmt

rgn = [129,135,30,35]
prj = 'M6c'
fig = pygmt.Figure()

# =======================
# MAIN
# =======================
fig.coast(
    projection = prj,
    region = rgn,
    frame = ['agf','NWes'],
    shorelines = '0.1p',
    resolution = 'h',
    land = 'grey60',
)

# =======================
# getting size
# =======================
with pygmt.clib.Session() as session:
    with pygmt.helpers.GMTTempFile() as tmpfile:
        session.call_module("mapproject",f"-W ->{tmpfile.name}")
        size = tmpfile.read().split()
        size = [float(s) for s in size]
        w = size[0]
        h = size[1]

# =======================
# cross section
# =======================
with pygmt.config(
    FONT_LABEL = '8p'
):
    # lat vs z -----------
    fig.shift_origin(
        xshift = w*1.1,
    )
    fig.basemap(
        frame = ['xagf+lDepth[km]','yagf+u@.','NEws'],
        projection = f'X{w*0.5}c/{h}c',
        region = [0,60,fig.region[2],fig.region[3]],
    )

    # lon vs z -----------
    fig.shift_origin(
        xshift = -1*w*1.1,
        yshift = -w*(0.5+0.1),
    )
    fig.basemap(
        frame = ['xagf+u@.','yagf+lDepth[km]','SWne'],
        projection = f'X{w}c/{-w*0.5}c',
        region = [fig.region[0],fig.region[1],0,60],
    )

fig.show()
メイン地図のサイズを取得することでプロットした二方向の深さ方向断面図．ただし図の大きさや位置は揃うものの，グリッドがCartesian座標であるため，緯度はメインの地図からわずかにずれが生じてしまう．

メイン地図のサイズを取得することでプロットした二方向の深さ方向断面図．ただし図の大きさや位置は揃うものの，グリッドがCartesian座標であるため，緯度はメインの地図からわずかにずれが生じてしまう．