続・Python による気象解析入門

地図をかく (2)

[その (1)] につづいて、cartopy によって地図をかく。 ここでは、正距円筒図法 (plate carree) の世界地図をつくるばあいにしぼって、 プログラムの組み立てを見ていく。

このページの例文は、つぎのウェブサイトで提供された例文をまねしたところがある。

• 山下 陽介, 気象データ解析のためのmatplotlibの使い方 (2022年3月13日版) https://yyousuke.github.io/matplotlib/index.html
  • 7章 cartopy の利用 (PDFファイル) https://yyousuke.github.io/matplotlib/matplotlib-7.pdf (2022-06-24 閲覧)

投影の指定

matplotlib と cartopy を利用する。

matplotlib で図をかくための手順にはいくつかの方法があるが、 山下さんの例文にあわせて、ひとつの figure の中に複数の subplot をもつ形を標準とし、 ひとつだけの図をかくときも add_subplot を 「fig を 1×1に分割してそのうち1番めの場所に作図する」ような形で 呼ぶ手順をとることにする。

例文では、プロット領域 (figure) として「fig」を確保し、 その第1の (このばあいには唯一の) subplot で、PlateCarree という地図投影をつかうこと、 その中心経線は 東経150度 (この数値は増田が変更した) とすることを指定している。

最後に plt.show() で画面に図を表示するが、 この例文のままでは、表示されるのは長方形の枠だけである。

[22a_mapworld0.py]

import matplotlib.pyplot as plt
import cartopy.crs as ccrs

# プロット領域の作成
fig = plt.figure()

# cartopy 呼び出し
ax = fig.add_subplot(1,1,1,projection=ccrs.PlateCarree(central_longitude=150))

# 画面に表示
plt.show()

海岸線をひく

「cartopy 呼び出し」の部分と「画面に表示」の部分のあいだに、つぎのコードを入れる。

# 海岸線を描く
ax.coastlines()

経線・緯線をひく

つぎのように ax.gridlines() をつかって、経線・緯線をひくことができる。 線をひく経度・緯度の値をきめるために、matplotlib.ticker の FixedLocator という関数をつかう。

# 経度・緯度線をひく
gl = ax.gridlines(crs=ccrs.PlateCarree())
gl.xlocator = mticker.FixedLocator(np.arange(-180, 180, 30)) # 経度線
gl.ylocator = mticker.FixedLocator(np.arange( -90,  90, 30)) # 緯度線

その準備として、つぎの手続きが必要になる。

import matplotlib.ticker as mticker

これを組みこんだプログラムは、つぎのようになる。

[22c_mapworld2.py]

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.ticker as mticker
import cartopy.crs as ccrs

# プロット領域の作成
fig = plt.figure()

# cartopy 呼び出し
ax = fig.add_subplot(1,1,1,projection=ccrs.PlateCarree(central_longitude=150))

# 海岸線を描く
ax.coastlines()

# 経度・緯度線をひく
gl = ax.gridlines(crs=ccrs.PlateCarree())
gl.xlocator = mticker.FixedLocator(np.arange(-180, 180, 30)) # 経度線
gl.ylocator = mticker.FixedLocator(np.arange( -90,  90, 30)) # 緯度線

# 画面に表示
plt.show()

なお、ax.gridlines で、線の種類を指定することができる。 たとえば、つぎのように変えると、赤の細めの破線になる。 (「\」 (バックスラッシュ、端末のフォントによっては円マーク) は、Pythonの文法上の約束で、 プログラムの行の継続をしめす。)

gl = ax.gridlines(crs=ccrs.PlateCarree(),\
    linewidth=0.5, linestyle='--', color='r', alpha=0.8)

経線・緯線の経度・緯度のラベルをつける

ax.gridlines で「draw_labels=True」というオプションをつかうと、 緯線・経線の緯度・経度の数値をしめすラベルが追加される。

gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=True)

ただし、この例の結果はあまりうまくない。 ラベルの文字が省略されてしまっていることがある。文字がかさなるのをさけるためらしい。 ax.gridlines でかかれるラベルの文字の大きさを制御する方法がよくわからないが、 plt.figure の呼び出しをつぎのように変更し、図をかく紙にあたるものを大きくすれば、 (かさなる心配がなくなるので) 省略されなくなる。

fig = plt.figure(figsize=(8, 6))

(図の右端の経度のラベルが省略されることはあいかわらずである。 右端が地球上で左端と同じ場所だからだと思うが、対策はわからない。)

また、図の上側と下側のラベルは一方でよいことが多い。 つぎのようにすれば、上側のラベルを省略できる。

gl.top_labels = False

右側のレベルは表示されていないが、 本来は表示されるはずのものである (中心経線を指定しないときは表示される)。 つぎのようにすれば、右側のラベルを省略することを指示できる。

gl.right_labels = False

以上の変更をしたプログラムが、つぎのものである。

[22d_mapworld3.py]

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.ticker as mticker
import cartopy.crs as ccrs

# プロット領域の作成
#fig = plt.figure()
fig = plt.figure(figsize=(8, 6)) # 図のサイズを大きくする

# cartopy 呼び出し
ax = fig.add_subplot(1,1,1,projection=ccrs.PlateCarree(central_longitude=150))

# 海岸線を描く
ax.coastlines()

# 経度・緯度線をひく
gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=True)
gl.xlocator = mticker.FixedLocator(np.arange(-180, 180, 30)) # 経度線
gl.ylocator = mticker.FixedLocator(np.arange( -90,  90, 30)) # 緯度線
gl.top_labels   = False # 上のラベルを消す
gl.right_labels = False # 右のラベルを消す

# 画面に表示
plt.show()

山下さんは、さらに、経度・緯度線とラベルとを別々にかくことによって、 ラベルの文字の大きさや、書きかた (たとえば経度に「E」をつけるかどうか) を 制御することをこころみているが、ここではそれを追いかけないことにする。

領域を限定した地図

[その (1)] でかいた「東アジア」の地図のように 経度・緯度で領域を限定した地図をかくことができる。

「cartopy 呼び出し」のすぐあとで、つぎのように ax.set_extent で 経度の下限・上限、緯度の下限・上限をまとめたリストで指定する。 符号は、東経・北緯が正、西経・南緯が負である。 (set_extent での 投影方法の指定は、位置座標の解釈のためのものなので、 作図する地図の投影図法が正距円筒図法でなくても PlateCarree を指定する。)

# 経度・緯度範囲を限定する
ax.set_extent([110, 170, 10, 60], ccrs.PlateCarree())

経度線・緯度線の間隔を図の範囲にあわせて設定しなおす。ここでは 10度にしてみた。

gl.xlocator = mticker.FixedLocator(np.arange(-180, 180, 10)) # 経度線
gl.ylocator = mticker.FixedLocator(np.arange( -90,  90, 10)) # 緯度線

海岸線のデータが世界地図用のものでは あらいので、もうすこしこまかいものにしてみる。 「50m」の m は (メートルではなく) million で、縮尺が 5千万分の1の地図に適切なくわしさ、という意味だそうだ。

# 海岸線を描く
ax.coastlines(resolution='50m')

以上の変更をしたプログラムはつぎのようになる。

[22e_map_eastasia.py]

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.ticker as mticker
import cartopy.crs as ccrs

# プロット領域の作成
fig = plt.figure()

# cartopy 呼び出し
ax = fig.add_subplot(1,1,1,projection=ccrs.PlateCarree(central_longitude=150))

# 経度・緯度範囲を限定する
ax.set_extent([110, 170, 10, 60], ccrs.PlateCarree())

# 海岸線を描く
ax.coastlines(resolution='50m')

# 経度・緯度線をひく
gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=True)
gl.xlocator = mticker.FixedLocator(np.arange(-180, 180, 10)) # 経度線
gl.ylocator = mticker.FixedLocator(np.arange( -90,  90, 10)) # 緯度線
gl.top_labels   = False # 上のラベルを消す
gl.right_labels = False # 右のラベルを消す

# 画面に表示
plt.show()