メソ数値予報モデル(MSM)のデータ領域可視化(その1)

メソ数値予報モデル(MSM)のデータ領域を可視化してみました。まずは結果から貼り付けておきます。

MSMのデータ領域(ランベルト座標、等緯度経度座標)
必要なパラメータさえ分かればMSMのデータ領域は自分で計算することはできますが、ここでは、公開されているMSMサンプルデータ(GRIB2)をpygribを用いて読み込み、緯度経度の情報を取り出して図示しています。必要なものは以下のとおりです。

Pythonライブラリについては、以下のimportが成功すれば準備OKです。

import pygrib
import numpy as np
import matplotlib.pyplot as plt

Debian/Ubuntu系であれば以下のようにaptでインストールできるでしょう。

sudo apt install -y python3-grib python3-numpy python3-matplotlib

コードは以下のとおりです。極めて単純です。気をつけないといけないのは、経度方向にnx, 緯度方向にny格子点があるとき、latlons() の返り値 lats, lons は (ny, nx) の2次元配列になっていることですかね(慣れない)。

# -*- coding: utf-8 -*-
# MSMサンプルデータをpygribで読み込み、データ領域を描画する
import pygrib
import numpy as np
import matplotlib.pyplot as plt

def extract_latlons(filename):
    data = pygrib.open(filename)
    lats,lons = data.message(1).latlons() # 読み込む物理量はなんでもよい
    nlats = np.concatenate([lats[0,:],lats[:,-1],lats[-1,::-1],lats[::-1,0]])
    nlons = np.concatenate([lons[0,:],lons[:,-1],lons[-1,::-1],lons[::-1,0]])
    return nlats,nlons

if __name__ == '__main__':
    file1 = "Z__C_RJTD_20240228000000_MSM_GPV_Rjp_Glm5km_Lm1-39_Pqc_FH00-39_grib2.bin"
    file2 = "Z__C_RJTD_20171205000000_MSM_GPV_Rjp_Lsurf_FH00-15_grib2.bin"
    lats1,lons1 = extract_latlons(file1)
    lats2,lons2 = extract_latlons(file2)
    
    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.set_aspect("equal")
    ax.set_xlim(100,160)
    ax.set_ylim(15,60)

    ax.plot(lons1,lats1,color="red",label="Lambert (817x661)")
    ax.plot(lons2,lats2,color="blue",linestyle=":",label="LatLon (481x505)")

    ax.text(0.01,-0.1,__file__,color="blue",fontsize="smaller",
            transform=ax.transAxes)
    
    plt.legend()
    plt.title("Meso-Scale Model's domain")
    plt.savefig("figure.png",dpi=150)
    plt.close

なお、今回の図示は簡素化したものであって、海岸線とかを一切描画していません。続編として今度はCartopyで凝った図を描いてみたいと思いますが今回はこのへんで。