MetPyでSSIを計算する

以前のnoteで紹介したが、MSMを使ってエマグラムを予想図として描いていた。

その後さらにMetPyが提供しているいくつかのAPIを使ってSSIを計算し、エマグラム上に表示させていた。

MetPyのReference GuideのcalcのAPIを見ていたら、直接SSIの計算ができるAPI（showalter_index）を発見。リリースノートを確認すると、Ver.1.1でリリースされたようだ。

Anacondaの環境にインストールしているMetPyのバージョンが1.0.1だったので、1.1へバージョンアップして、このAPIを使ってSSIを計算するコードに変更してみた。

Reference Guideを見ただけでは使い方がよく分からなかったのだが、GitHubにあるMetPyのソースコードを見ていたら、テストコードがあったのでこれを参考にした。

（GitHubのMetPyのテストコードより抜粋）

def test_showalter_index():
   """Test the Showalter index calculation."""
   pressure = units.Quantity(np.array([931.0, 925.0, 911.0, 891.0, 886.9, 855.0, 850.0, 825.6,
                                       796.3, 783.0, 768.0, 759.0, 745.0, 740.4, 733.0, 715.0,
                                       700.0, 695.0, 687.2, 684.0, 681.0, 677.0, 674.0, 661.9,
                                       657.0, 639.0, 637.6, 614.0, 592.0, 568.9, 547.4, 526.8,
                                       500.0, 487.5, 485.0]), 'hPa')
   temps = units.Quantity(np.array([18.4, 19.8, 20.0, 19.6, 19.3, 16.8, 16.4, 15.1, 13.4,
                                    12.6, 11.2, 10.4, 8.6, 8.3, 7.8, 5.8, 4.6, 4.2, 3.4, 3.0,
                                    3.0, 4.4, 5.0, 5.1, 5.2, 3.4, 3.3, 2.4, 1.4, -0.4, -2.2,
                                    -3.9, -6.3, -7.6, -7.9]), 'degC')
   dewp = units.Quantity(np.array([9.4, 8.8, 6.0, 8.6, 8.4, 6.8, 6.4, 4.0, 1.0, -0.4, -1.1,
                                   -1.6, 1.6, -0.2, -3.2, -3.2, -4.4, -2.8, -3.6, -4.0, -6.0,
                                   -17.6, -25.0, -31.2, -33.8, -29.6, -30.1, -39.0, -47.6,
                                   -48.9, -50.2, -51.5, -53.3, -55.5, -55.9]), 'degC')
   result = showalter_index(pressure, temps, dewp)

気圧、温度、露点温度（温度と湿度から計算）の各データを下層から順にndarray配列して引数として渡せば計算できる。実際には850hPaと500hPaの温度、露点温度があれば計算できるのではないだろうか。

SSIのAPIを使わずに計算していた以前の値と比較しても齟齬はなかった。

SSIの天気図が描ければ興味深いのだが、天気図にすると実行時間が結構かかりそうな気がする。

MetPyには多くのAPIが提供されているので今後もいろいろと試してみたい。

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