MSM-S 解析値を、GrADS で単純バイナリに変換し、Python で作図する (2022年8月の事例)

気象庁 MSM 解析値の概略は、、 [気象庁のMSM解析値 (京大生存研のNetCDF版) を GrADS で 読む] のページに書いた。

気象庁から発信されている形式は GRIB だが、 ここでは、京大生存研で netCDF形式に変換されたものをつかう。 Python には netCDF を読み書きするパッケージがあるので、 それをつかうのがよさそうだ。 しかし、わたしは、まだそのつかいかたを修得していない。

ここでは、 GrADS で netCDF を読んで単純バイナリ形式で書きだし、 Python で単純バイナリ形式ののデータを読んで作図する手順をしめす。

この手順には、つぎのような欠点がある。

• 単純バイナリ形式のファイルの内容には、 気象要素の名まえ、物理量の単位、日時、格子点の位置や配置順序などがふくまれない。 ファイル名でいくらかの情報をつたえることができるが、そのほかは人が別途管理するしかない。 netCDF ならば、そのような情報をファイル内に「メタデータ」あるいは追加の変数値として ふくめておくことができる。
• GrADS の fwrite で書きだせる単純バイナリデータは、2次元の格子点値である。 3次元の格子点値を書きだすには、層ごとにわける必要があるらしい。 層ごとの水平分布を見るにはそれでもよいが、断面を見たいときなどには不便だ。

したがって、この方法は MSM解析値のうち水平2次元の量をふくむサブセット「MSM-S」にかぎってつかうことにする。 気圧座標の3次元の量をふくむ「MSM-P」については、 Python で netCDF を読む方法を修得してからあつかうことにしたい。

事例

2022年 8月21日 から 9月17日までを対象とする。 (8月24-25日の観測実習の期間をふくむように設定し、その後に延長した。)

京大生存研の http://database.rish.kyoto-u.ac.jp/arch/jmadata/data/gpv/netcdf/MSM-S/ の下の 2022 のディレクトリにある各日のファイルをダウンロードした。 ファイルは 1日ごと (UTCの0時から23時まで) になっており、 8月21日ならば「0821.nc」のような名まえになっている。

GrADSで読む

この netCDF ファイルは、GrADS で「sdfopen 0821.nc」のようにすれば 開くことができる。 しかし、時系列としてつながるファイルをつづけて読むのにはこの方法では 不便なので、つぎのような ctl ファイルをつくり、 GrADS で「xdfopen msm_s_hourly202208.ctl」の形で開くことにする。 (この ctl ファイルは xdfopen 用であり、open 用でないことに注意。)

• msm_s_hourly202208.ctl

つぎの変数をあつかうことができる。

• psea: 海面気圧 [Pa]
• sp: 地表面気圧 [Pa] [注]
• u: 風の東西成分 (西風が正) [m/s]
• v: 風の南北成分 (南風が正) [m/s]
• temp: 気温 [K]
• rh: 相対湿度 [無次元]
• このほかに、降水量 (予報値) と 雲量 (全雲量、上層、中層、下層) があるが、ここでは省略する。

[注]

• これは海面気圧ではなく、 MSMモデルの地表面 (海では海面、陸では地面) の高さでの気圧である。 観測点での現地気圧と比較する際には、モデルの地表面の標高と観測点の標高の差を考慮する必要がある。
• 地表面気圧の変数名が、この教材ページ内で、不統一になってしまった。 増田は、つぎの段階で書きだす単純バイナリファイル名では「psfc」とした。 気圧をあらわす変数名を p ではじまるものにそろえ、海面気圧との区別を明確にしたかったからである。 そのつぎの段階で書きだす画像ファイル名では「psl」とした。 これは、まえにつくった教材 (JRA-55データにもとづく作図) にあわせたものである。

GrADSで単純バイナリ形式で書きだす

上にあげた 6つの変数の値を、それぞれ単純バイナリ形式で書きだす。 ただし、毎時ではなく 3時間おき (UTCの 0, 3, 6, 9, 12, 15, 18, 21時) の値を 書きだすことにした。 これが解析値 (観測値をとりこんだ値) で、他の時刻の値は予報値だからである。 上の ctl ファイルがあることを前提として、 つぎのスクリプトを GrADS の中で実行すればよい。 (結果のファイルが置かれるディレクトリなどの準備が必要だが、質問があれば説明する。)

• fwrite202208.gs

Python による作図: 等圧線

東アジア (日本付近) の 北緯 20～50度、東経 120～150度 の範囲の地図を Cartopy をつかって描き、 対象期間の3時間きざみの各時刻 の 海面気圧の等値線を contour でひいてみた。 プログラムはつぎのものである。

• contour_msm_psea202208.py (文字コード UTF-8)

[格子点型の気象データ (JRA-55 海面気圧) の地図上表示 (セミナー教材 2022-07-05)]とは、 材料としたデータがちがうだけである。 そこでつかった JRA-55 は緯度経度 1.25度格子、 MSM-S は 緯度 0.05度、経度 0.0625度 (約 5 km) の格子なので、 等圧線の空間的くわしさがかなりちがう。

例として、2022年8月24日 06 UTC (15 JST) の MSM-S による等圧線図をしめす。

Python による作図: 地上風

関東地方とその周辺 (仮称 the Greater Kanto なので gkanto) の 北緯 35～38度、東経 137.5～141.5度 の範囲の地図を Cartopy をつかって描き、 対象期間の3時間きざみの各時刻 の 地上風 の風速ベクトルの大きさによって contourf でぬりわけ、 風速ベクトルを quiver で矢印表示してみた。 矢印の寸法や形、ぬりわけの色の選択などが、うまくないが、ひとまずできたものをしめしておく。 プログラムはつぎのものである。

• map_msm_wind_gkanto202208.py (文字コード UTF-8)

例として、2022年8月24日 06 UTC (15 JST) の MSM-S による地上風の図をしめす。

海上のほうが陸上よりも風が強い。実際にそうなのだろうが、モデルでは誇張されているかもしれない。 海上で矢印が明確に見えるようにすると陸上では風向がよくわからず、 陸上の風向がわかるようにすると海上で矢印がぶつかりあってしまう。 どちらを見たいかによって、矢印の寸法のきめかたをつかいわけるしかなさそうだ。