気象庁二重偏波気象レーダー極座標データのGRIBをnetCDF(CF/Radial 規約)に変換する

# TL;DR

変換プログラム（作りかけ）を[GitHub](https://github.com/wm-ytakano/jmardr_cfradial)に置いておくので自己責任でどうぞ。

# はじめに

[気象庁レーダー毎極座標GPVのGRIBをnetCDF(CF/Radial 規約)に変換する](https://qiita.com/wm-ytakano/items/c1f460ec257d4ad2a325)の続きです。CF/Radial規約に関する説明はこの記事を参照ください。

前回扱った気象業務支援センター経由で配信されている[気象庁レーダー毎極座標データ](https://www.data.jma.go.jp/add/suishin/shiyou/pdf/no13702)は、リアルタイムで提供されているものの、レーダーが二重偏波化したことにより得られる情報（偏波パラメータ）は含まれていません。この偏波パラメータを含んだ東京の二重偏波レーダーの過去データは、「利用ニーズ把握のために提供する気象データ」という形で、[気象過去データの利用環境（令和2年度）](https://www.data.jma.go.jp/developer/past_data/index.html)という枠組みでAWS上で公開されています。今回はこのデータをNetCDF([CF/Radial規約](https://github.com/NCAR/CfRadial))に変換することにします。

# データについて

## パラメータ(変数)

[二重偏波気象レーダー極座標データ](https://www.data.jma.go.jp/developer/past_data/app/2hempa_format_20200914.pdf)は、次の8変数が提供されています：

| パラメータ名                     | 略称 | パラメータ番号(実際) | パラメータ番号(仕様) |
| -------------------------- | ---- | -------------------: | -------------------: |
| 水平偏波反射強度           | ref  |                  195 |                  195 |
| ドップラー速度             | vel  |                  228 |                    2 |
| ドップラー速度幅           | vsw  |                  230 |                    0 |
| 反射因子差 (Zdr)           | zdr  |                  197 |                  197 |
| 受信信号偏波間位相差 (ψdp) | psd  |                  198 |                  198 |
| 偏波間位相差変化率 (Kdp)   | kdp  |                  202 |                  202 |
| 偏波間相関係数 (ρhv)       | rhv  |                  199 |                  199 |
| 品質管理情報               | qci  |                  215 |                  206 |

仕様書にあるパラメータ番号と、実際にGRIBに格納されているパラメータ番号が異なるものがあるので、併記しています。仕様書にはパラメータ番号207-254は保留と書かれているので、単に記載されていないだけだと思われます。

反射強度について、GRIBファイル名では ref (基底反射率) が使われていますが、実際にはパラメータ番号195の zhh (水平偏波反射強度) が格納されている点に少し注意が必要です。実用上は両者は同じ変数だとみなせるので、こちらは大きな問題ではありません。

変数の意味についてはこれは技術ノートなので解説しませんが、日本語であれば 深尾・浜津(2005)『気象と大気のレーダーリモートセンシング』京大出版会 を読むとよいでしょう。京都大学学術情報リポジトリから[全文PDF](https://repository.kulib.kyoto-u.ac.jp/dspace/bitstream/2433/49766/1/fukao2005a.pdf)をダウンロード出来ます。

## スキャンシーケンス

東京レーダーのスキャンシーケンスは次の表の通りです。

ただし、走査番号20の仰角90°(vertical pointing)については、偏波レーダー用の更正用の走査ですが、時刻(reference time)の分の下一桁が0のものについてのみ行われ、下一桁が5の場合には行われません。つまりボリュームスキャン2回に1回、10分おきということになります。

したがって、走査番号20を除けば、レンジ間隔は250mで一定であるということになります。

| 走査番号 | 仰角(°) | 最大レンジ(km) | レンジ間隔(m) | 方位角間隔(°) |
| -------: | ------: | -------------: | ------------: | ------------: |
| 1        | 25.0    | 64             | 250           | 0.70          |
| 2        | 17.8    | 64             | 250           | 0.70          |
| 3        | 12.8    | 64             | 250           | 0.70          |
| 4        | 9.3     | 120            | 250           | 0.70          |
| 5        | 6.7     | 120            | 250           | 0.70          |
| 6        | 4.8     | 180            | 250           | 0.70          |
| 7        | 3.5     | 180            | 250           | 0.70          |
| 8        | 2.5     | 400            | 250           | 0.35          |
| 9        | 1.7     | 400            | 250           | 0.35          |
| 10       | 1.1     | 400            | 250           | 0.35          |
| 11       | 0.7     | 400            | 250           | 0.35          |
| 12       | 0.3     | 400            | 250           | 0.35          |
| 13       | 0.0     | 400            | 250           | 0.35          |
| 14       | 0.3     | 250            | 250           | 0.70          |
| 15       | 0.3     | 250            | 250           | 0.35          |
| 16       | 0.7     | 250            | 250           | 0.35          |
| 17       | 0.7     | 250            | 250           | 0.70          |
| 18       | 1.3     | 150            | 250           | 0.70          |
| 19       | 2.2     | 150            | 250           | 0.70          |
| 20       | 90.0    | 64             | 125           | 0.70          |

## 従来のレーダー極座標データからの変更点

[二重偏波気象レーダー極座標データ](https://www.data.jma.go.jp/developer/past_data/app/2hempa_format_20200914.pdf)は、前回扱った[気象庁レーダー毎極座標データ](https://www.data.jma.go.jp/add/suishin/shiyou/pdf/no13702)とGRIB形式である点は同じですが、いくつか異なる点があります：

- 使われているテンプレートが3.50120から3.50121、4.51022から4.51123に変更になった。
    - これにより例えば、仰角の定義が第4節（プロダクト定義節）から第3節（格子系定義節）へと記述される場所が変更になった

- データの圧縮がランレングス圧縮（テンプレート5.200および7.200）から、GRIB全体をgzip圧縮することを前提として、（圧縮率の低い）単純圧縮（テンプレート5.0および7.0）に変更された

- ボリュームスキャンにかかる時間が10分から5分になったことに対応して、データも5分ごとに別ファイルに分かれるようになった

- 距離(range)と方位角(azimuth)の間隔は一定だったのが、(250m, 0.7°), (250m, 0.35°) および (125m, 0.35°)の3種類になった（詳細はスキャンシーケンスの表を参照）

- これまでは変数別で1つのファイルに全ての走査(仰角)が格納されていたが、走査(仰角)ごと別ファイルに格納されることになった

- ドップラー速度は以前は観測される仰角のうち一部の仰角でしか存在しなかったが、反射強度と同じく全ての仰角で格納されるようになった



# 変換プログラムの作成
前回とやることは同じで、[二重偏波気象レーダー極座標データフォーマット (PPI・RHI観測共通)](https://www.data.jma.go.jp/developer/past_data/app/2hempa_format_20200914.pdf)と[CfRadial Data File Format Version 1.4](https://github.com/NCAR/CfRadial/blob/master/docs/CfRadialDoc.v1.4.20160801.pdf)を見比べて変換するコードを書くだけです。前回と変わった点は次の通りです：

- 単純圧縮になったため、wgrib2を使わずPythonで資料節を直接デコードするようにした

現在、いくつかの課題が残っています：

- ボリュームスキャンを1つのNetCDFにまとめて格納するように変更する

- 8bitで表現される品質管理情報をそのまま格納するか、単偏波系データ品質・二重偏波系データ品質・選択的MTIによる結果に変数を分けて格納するかを決める

# 謝辞

- 二重偏波気象レーダー極座標データは[気象過去データの利用環境（令和2年度）](https://www.data.jma.go.jp/developer/past_data/index.html)の枠組みで取得しました。データ公開の要望を実現して頂いた気象庁の方々に御礼申し上げます。