Pythonモジュール『wxparams』の使い方

はじめに

2020年春、私たちの生活を一変させた某世界的イベントの影響を受けて暇になったので(汗)、気象パラメータを計算できるPythonモジュール『wxparams』を作ってGitHubで公開しました。

Yoshiki443/weather_parameters

このnoteは『wxparams』のチュートリアル的な位置付けとして書いています。チュートリアルだなんて大袈裟な…と我ながら思うのですが、幸運にも人に紹介する機会があったりするので、まとめておくと役に立つかなと思って書きました。

『wxparams』のマニュアルとしてはQiitaにまとめていますので、noteでは使い方の例をご紹介いたします。

気象要素を計算するPythonモジュールを公開します - Qiita

基本的な使い方

wxparamsのインストールについては、上のQiita記事をご参照ください。ここではインストールは済んでいるものとして進めます。

まずはwxparamsのインポートをします。ここではwxという省略名を付けてインポートしています。また一緒に使うnumpy・pandasもインポートします。

import wxparams as wx
import numpy as np
import pandas as pd

ではさっそく使い方です。

wxparamsでは、ある気象要素から別の気象要素を計算する関数を提供しています。まず代表的な関数として、気温・相対湿度から露点温度を計算するRH_to_Td関数を使ってみます。

以下のコードを実行してみて下さい。気温20度・相対湿度75%の時の露点温度を計算しようとしていますが、このコードではエラーになってしまいます。

'''
エラーになります
'''
TEMP = 20.0
RH = 75.0
Td = wx.RH_to_Td(TEMP, RH)
print(Td)

wxparamsはインプットとしてnumpy.ndarray型を想定していますので、例えば以下のようなコードにして実行します。

'''
numpy.ndarray型でインプット
'''
TEMP = np.array([20.0])
RH = np.array([75.0])
Td = wx.RH_to_Td(TEMP, RH)
print(Td)

'''
(出力)
[15.4380164] 
'''​

戻り値もnumpy.ndarray型になっています。

またpandas.Series型でも動作します。上で使ったデータをデータフレームに格納して、以下のように実行します。

'''
pandas.DataFrameを作成
'''
df = pd.DataFrame({"TEMP":TEMP, "RH":RH})
df["Td"] = wx.RH_to_Td(df["TEMP"], df["RH"])
print(df)

'''
(出力)
   TEMP    RH         Td
0  15.8  82.3  15.438016
'''​

データフレームで気象データを扱う

上述の通り、pandas.Series型で引数として渡すことができるので、いわゆるテーブルデータの形で気象データを扱う時に、wxparamsを利用できます。

ここでは練習用データフレームを以下のように作成します。データの内容は、気温[K]・相対湿度[%]・気圧[Pa]・東西風[m/s]・南北風[m/s]です。

'''
気温[K], 相対湿度[%], 気圧[Pa], 風[m/s]
'''
df = pd.DataFrame({"気温":[270.65, 292.45, 303.95],
                   "相対湿度":[85.2, 98.3, 48.1],
                   "気圧":[50000, 92500, 101325],
                   "東西風":[2.9, -5.2, 13.5],
                   "南北風":[-8.1, -2.9, 5.8]})

ここで単位に注意して下さい。あえて単位を[K]や[Pa]にしています。
数値予報データは、通常上記のような単位でデータを持っていますが、wxparamsでインプットとして想定している単位とは異なる場合があります。どのような単位を想定しているかはQiita記事にまとめてありますので、ご参照ください。

それでは上で定義したデータフレームを使い、露点温度・相当温位・風向風速の計算をしてみます。

'''
露点温度[C]を計算
気温[K]は[C]に変換する
'''
df["気温"] = df["気温"] - 273.15
df["露点温度"] = wx.RH_to_Td(df["気温"], df["相対湿度"])

'''
相当温位[K]を計算する
気圧は100で割って[hPa]に変換する
'''
df["気圧"] = df["気圧"] / 100
df["相当温位"] = wx.Theta_e(df["気温"], df["露点温度"], df["気圧"])

'''
風向風速
'''
df["風速"], df["風向"] = wx.UV_to_SpdDir(df["東西風"], df["南北風"])
print(df)

'''
(出力)
    気温   相対湿度      気圧    東西風    南北風
0  -2.5     85.2    500.00     2.9     -8.1
1  19.3     98.3    925.00    -5.2     -2.9
2  30.8     48.1   1013.25     13.5     5.8

      露点温度      相当温位         風速         風向
0   -4.642878   348.350507    8.603488   340.301379
1   19.025067   343.365456    5.953990    60.851928
2   18.571997   342.789011   14.693196   246.750207
'''​

露点温度を計算するRH_to_Td関数は、引数として気温[C]・相対湿度[%]を取り、露点温度[C]を返します。このため、気温が[K]単位であれば[C]へ変換します。

相当温位を計算するTheta_e関数は、引数として気温[C]・露点温度[C]・気圧[hPa]を取り、相当温位[K]を返します。このため、気圧は[Pa]を[hPa]へ変換します。

風向風速を計算するUV_to_SpdDir関数は、風のUV成分から風向・風速を計算します。風速の単位は[m/s]でも[kt]でも、なんでもOKです。なおこの関数は戻り値が風速・風向の順番で2種類ある点に注意が必要です。

多次元配列で気象データを扱う

数値予報などのGPVデータをPythonで読み込んだ場合、気象要素ごとに2次元配列になっているか、時間軸も併せて3次元配列になっていると思います。wxparamsでは、Numpy配列データをそのまま引数としてインプットできます。

では練習用として、乱数を発生させてデータを作成します。ここでは(3, 3)のNumpy配列を作成します。

np.random.seed(0)
T850 = np.random.normal(15, 2, (3,3))
RH850 = np.clip(0, 100, np.random.normal(80, 10, (3,3)))
T500 = np.random.normal(-5, 1, (3,3))
Usfc = np.random.normal(0, 5, (3,3))
Vsfc = np.random.normal(0, 5, (3,3))

print("850hPa気温")
print(T850)
print("850hPa相対湿度")
print(RH850)
print("500hPa気温")
print(T500)
print("地上風U成分")
print(Usfc)
print("地上風V成分")
print(Vsfc)

'''
(出力)
850hPa気温
[[18.52810469 15.80031442 16.95747597]
[19.4817864  18.73511598 13.04544424]
[16.90017684 14.69728558 14.7935623 ]]
850hPa相対湿度
[[84.10598502 81.44043571 94.54273507]
[87.61037725 81.21675016 84.43863233]
[83.33674327 94.94079073 77.94841736]]
500hPa気温
[[-4.6869323  -5.85409574 -7.55298982]
[-4.3463814  -4.1355638  -5.74216502]
[-2.73024538 -6.45436567 -4.95424148]]
地上風U成分
[[-0.93591925  7.66389607  7.34679385]
[ 0.77473713  1.8908126  -4.43892874]
[-9.90398234 -1.73956075  0.78174485]]
地上風V成分
[[ 6.1514534   6.01189924 -1.93663409]
[-1.51151375 -5.24276483 -7.10008969]
[-8.53135095  9.75387698 -2.54826091]]
'''

乱数シードを固定することで、上記出力と同じ値が得られると思います。また相対湿度は0〜100に値が収まるよう、念のためクリッピング処理を入れています。

それでは作成した850hPa気温と、850hPa相対湿度を使って、850hPa露点温度を計算してみます。

Td850 = wx.RH_to_Td(T850, RH850)

print("850hPa露点温度")
print(Td850)

'''
(出力)
850hPa露点温度
[[15.79482329 12.63136282 16.07567762]
[17.37275694 15.45160182 10.48394266]
[14.05883974 13.89537795 10.98765309]]
'''​

次に地上風のUV成分から風向・風速を計算してみます。

Wspd, Wdir = wx.UV_to_SpdDir(Usfc, Vsfc)

print("地上風速")
print(Wspd)
print("地上風向")
print(Wdir)

'''
地上風速
[[ 6.22224428  9.74054596  7.59775832]
[ 1.69849682  5.57330739  8.37349162]
[13.0718329   9.9077842   2.66547532]]
地上風向
[[171.34901931 231.88779301 284.76739064]
[332.86226872 340.16806954  32.01334813]
[ 49.25818908 169.88786969 342.94525245]]
'''

最後にSSI関数を使ってSSIを計算してみます。

SSI関数では引数として気圧・気温・露点温度を取ります。
850hPa気温・850hPa露点温度・500hPa気温のサンプルデータは作ってありますので、あとは気圧の値が必要です。

気圧は850hPaと500hPaそのものなのですが、引数間でNumpy配列のshapeを合わせておく必要があります。そこで次のようにコードを記述します。

P850 = np.full_like(T850, 850)
P500 = np.full_like(T850, 500)

print("850hPa気圧")
print(P850)
print("500hPa気圧")
print(P500)

'''
(出力)
850hPa気圧
[[850. 850. 850.]
[850. 850. 850.]
[850. 850. 850.]]
500hPa気圧
[[500. 500. 500.]
[500. 500. 500.]
[500. 500. 500.]]
'''​

これで準備は整いましたので、SSIを計算してみましょう。

SSI = wx.SSI(P850, P500, T850, T500, Td850)

print("SSI")
print(SSI)

'''
(出力)
SSI
[[-0.17765495  3.04421969 -2.58777985]
[-1.83295367  0.61266618  6.66078908]
[ 4.24912474  1.80994707  5.99943283]]
'''

最後に

今回は露点温度・相当温位・風向風速・SSIを計算する関数を紹介しました。これ以外にも水蒸気量や絶対湿度などを計算する関数を用意してますので、ぜひ様々な用途でご活用いただければ幸いです。

最後までお読みいただき、ありがとうございました。

仕事のご依頼･お問合せ | Weather Data Science

(※)
このnoteは、2020年に私自身が作成に参加し、現在気象庁ウェブサイトで公開されている『気象データ利活用素材集』の一部を利用し、コンテンツを追加して執筆しました。

気象庁 | 「気象データアナリスト育成講座」の認定制度について