見出し画像

PythonでGMTを使わずにGSMのGPVデータから天気図を描く

はじめに

ここ一ヶ月余り、週末にはPythonで天気図を描くことをテーマに、試行錯誤してきた。

とりあえず、先週はエマグラムを書いてみたりした。

MetPyExample Galleryのコードを見ていたら、GSMやMSMのGPVデータから、PyGMTを使わずとも簡単に天気図を描けそうな気がしてきた。

500hPa高度場、温度場を描く

手始めに、AUPQ35の下図の500hPa高層天気図相当の天気図を、GSMのGPVデータを使って描いてみる。

基本的にnoteにはソースコードをそのまま載せないことにしているのだが、今回はそのまま掲載する。

#!/usr/bin/env python
# coding: utf-8

import cartopy.crs as ccrs
import cartopy.feature as cfeature
import matplotlib.pyplot as plt
import matplotlib.ticker as mticker
from metpy.units import units
import numpy as np
from scipy.ndimage import gaussian_filter
import pygrib as grib

gpv = grib.open('Z__C_RJTD_20200821120000_GSM_GPV_Rgl_FD0000_grib2.bin')


height, lat, lon = gpv.select(shortName = 'gh', level = 500)[0].data()
temp, _, _ = gpv.select(shortName = 't', level = 500)[0].data()


height_500 = gaussian_filter(height, sigma = 3.0)
temp_500 = (gaussian_filter(temp, sigma = 3.0) * units.kelvin).to(units.celsius)

mapcrs = ccrs.LambertConformal(central_longitude=140,
                              central_latitude=35,
                              standard_parallels=(30, 60))

datacrs = ccrs.PlateCarree()

fig = plt.figure(1, figsize = (294 / 25.4, 210 / 25.4))
ax = plt.subplot(111, projection = mapcrs)
ax.set_extent([110,170,10,60], ccrs.PlateCarree())

ax.add_feature(cfeature.COASTLINE.with_scale('50m'))
ax.add_feature(cfeature.BORDERS.with_scale('50m'))

ax.gridlines(xlocs=mticker.MultipleLocator(10), 
            ylocs=mticker.MultipleLocator(10),
            linestyle='-', color='gray')

clevs_500_hght = np.arange(0, 8000, 60)
cs = ax.contour(lon, lat, height_500, clevs_500_hght, colors='black', transform = datacrs)
plt.clabel(cs, fmt='%d')

clevs_500_temp= np.arange(-60, 15, 3)
cs = ax.contour(lon, lat, temp_500, clevs_500_temp, colors = 'black', transform = datacrs)
plt.clabel(cs, fmt='%d')

plt.title('500-hPa Geopotential Heights (m) and Temperature(C)', loc = 'left')

plt.subplots_adjust(bottom=0, top=1)​

下記ができあがった500hPaの高層天気図。

画像2

ついでに、300hPaの高度と風速を描いてみると下記のようになる。

画像2

さいごに

現在、BashプログラムでGMT6.xを使って自動作成している色塗り天気図システム。

全ての天気図をこの方法で置き換えられるか検証が必要だが、ここまでmatplotlib、cartopy、scipy、pygrib、そしてMetPyを使ってPythonで天気図が書けてしまうとGMTやPyGMTに固執する意味はないだろう。

今後、他の天気図を含めてPythonへの置き換えができるのか、確認してみることにする。

天気図の作図にあたっては、京都大学生存圏研究所の生存圏データベースのGPVデータを使用させて頂いています。

よろしければ、サポートをお願いします。 より多くの方に役立つnoteを書けるよう頑張ります!!