概要
マップ上に図形等を描画するのは GIS が得意とするところですが、大量のファイルの入力が必要であったり、描画した画像を大量に出力したりする必要があるような時には Python のコードで描画の処理を行いたいところです。
ここでは、GIS で行うような基本的なマップや図形の描画を Python で行う方法を紹介します。
環境
主に以下のライブラリを利用します。
- Python 3.11.6
- cartopy 0.22.0
- contextily 1.4.0
cartopy
は matplotlib
のグラフで地理空間を扱うためのライブラリです。
contextily
は地図タイルを取得するライブラリで、取得されたタイルを cartopy
で地理空間に拡張されたグラフに投影します。
ここでは、conda
を用いて以下のコマンドで環境を作成しています。
conda create -n geo python=3.11
conda activate geo
conda install contextily cartopy
実装
ここでは、東京駅周辺のマップと合わせて東京駅の位置を示すポイントを描画することとして、以下のように実装しました。
import cartopy.crs as ccrs
import contextily as cx
import matplotlib.pyplot as plt
import pyproj
from matplotlib.offsetbox import AnchoredText
def main():
# 日本語フォントの設定
plt.rcParams["font.family"] = "Meiryo"
# 対象地域の適当な座標系を設定
crs = ccrs.epsg(6677) # 平面直角座標系 9系(JGD2011)
# 座標系に対応したオブジェクト
fig = plt.figure(figsize=(8, 8)) # 実際の出力サイズは描画領域で決まるため、ここでは最大サイズとして指定。
ax = fig.add_subplot(1, 1, 1, projection=crs) # projectionへのcrsの指定により地理空間に投影
# 描画領域の設定(東京駅周辺の座標を確認して入力)
ax.set_xlim(-8000.0, -4000.0)
ax.set_ylim(-37000.0, -34000.0)
# ベースマップ(OpenStreetMap)の追加
cx.add_basemap(ax, crs=crs, zoom=16, source="http://{s}.tile.openstreetmap.org/{z}/{x}/{y}.png")
# 図形(ポイント)の描画
lat, lon = 35.6810912, 139.7671861 # 東京駅の緯度経度
transformer = pyproj.Transformer.from_crs("EPSG:4326", crs)
y, x = transformer.transform(lat, lon)
plt.scatter(x, y, s=200, c="red", alpha=0.5, label="東京駅")
# 凡例の描画
ax.legend(loc="upper right", prop={"size": 20})
# OpenStreetMapのクレジットを右下隅に表示
text_box = AnchoredText("© OpenStreetMap contributors", loc="lower right", borderpad=0, frameon=True)
plt.setp(text_box.patch, facecolor="white", alpha=0.6, linewidth=0)
ax.add_artist(text_box)
# プロットエリアの枠を非表示
ax.set_axis_off()
# プロットエリア以外の部分をなくし、地図部分のみをファイルに出力
fig.tight_layout(pad=0)
fig.savefig("output.png", bbox_inches="tight", pad_inches=0)
if __name__ == "__main__":
main()
以下のようにベースマップやポイント、凡例やOpenStreetMapのクレジットの表示ができました。
参考文献