はじめに
juyter-notebookでpygribを使用し，気象庁のＬＦＭ格子点データを読んだ．
ファイルをダウンロード
京都大学生存圏研究所のデータベースからデータをダウンロードする．
ダウンロードが成功すれば下記の「Out」のように出力される．
In
import urllib.request
from bs4 import BeautifulSoup
url = 'http://database.rish.kyoto-u.ac.jp/arch/jmadata/data/gpv/original/2020/08/01/'
# urlの情報を取得
response = urllib.request.urlopen(url)
soup = BeautifulSoup(response)
response.close()
## a href にファイル名が入っていそうなので，１行ずつリストに入れる
file_list = []
for row in soup.findAll('a'):
file_list.append(row.get('href'))
## binファイルのみ抽出
file_list = [file for file in file_list if file.endswith('.bin')==True]
# ダウンロードしたファイルを置く場所のPATH
DL_PATH = 'F:/'
# file_list中の始めのファイルのみダウンロード
urllib.request.urlretrieve(url+fname, DL_PATH+file_list[0])
Out
('F:/Z__C_RJTD_20200801000000_CWM_GPV_Rjp_Gll0p05deg_FD0000-0300_grib2.bin',
<http.client.HTTPMessage at 0x1d09ef42d30>)
ファイルを開く
ファイルを開き，read()することでgribmessageのリストとして中身を確認できる．
In
import pygrib
grbs = pygrib.open('F:/Z__C_RJTD_20171205010000_LFM_GPV_Rjp_L-pall_FH0100_grib2.bin')
grbs.read()
Out
[1:Geopotential height:gpm (instant):regular_ll:isobaricInhPa:level 100000.0 Pa:fcst time 1 mins:from 201712050100,
2:U component of wind:m s**-1 (instant):regular_ll:isobaricInhPa:level 100000.0 Pa:fcst time 1 mins:from 201712050100,
3:V component of wind:m s**-1 (instant):regular_ll:isobaricInhPa:level 100000.0 Pa:fcst time 1 mins:from 201712050100,
...
90:V component of wind:m s**-1 (instant):regular_ll:isobaricInhPa:level 10000.0 Pa:fcst time 1 mins:from 201712050100,
91:Temperature:K (instant):regular_ll:isobaricInhPa:level 10000.0 Pa:fcst time 1 mins:from 201712050100,
92:Vertical velocity:Pa s**-1 (instant):regular_ll:isobaricInhPa:level 10000.0 Pa:fcst time 1 mins:from 201712050100]
selectメソッド
引数の条件でgribmessageを検索できる．
| 引数 | 入力型 |
|---|---|
| name | str |
| shortName | str |
| units | str |
| forecastTime | int |
In
grbs.select(name='Geopotential height')
grbs(name='Geopotential height') # selectメソッドと同じことができる
Out
[1:Geopotential height:gpm (instant):regular_ll:isobaricInhPa:level 100000.0 Pa:fcst time 1 mins:from 201712050100,
7:Geopotential height:gpm (instant):regular_ll:isobaricInhPa:level 97500.0 Pa:fcst time 1 mins:from 201712050100,
13:Geopotential height:gpm (instant):regular_ll:isobaricInhPa:level 95000.0 Pa:fcst time 1 mins:from 201712050100,
...
78:Geopotential height:gpm (instant):regular_ll:isobaricInhPa:level 20000.0 Pa:fcst time 1 mins:from 201712050100,
83:Geopotential height:gpm (instant):regular_ll:isobaricInhPa:level 15000.0 Pa:fcst time 1 mins:from 201712050100,
88:Geopotential height:gpm (instant):regular_ll:isobaricInhPa:level 10000.0 Pa:fcst time 1 mins:from 201712050100]
変数の読み込み
selectメソッドで選んだ変数には様々な情報が格納されている．
In
grb = grbs.select(name='Geopotential height')[0]
# keyを取得し，その中に格納されている値を読む
for grb_key in grb.keys():
if grb.valid_key(grb_key)==True: # keyの欠落がなく，値が読み取れるものを出力
print(grb_key, ' : ',grb[grb_key])
else: # 読めなかったkeyは赤字で出力
print('\033[31m'+grb_key+' : Can not read'+'\033[0m')
Out
globalDomain : g
GRIBEditionNumber : 2
tablesVersionLatestOfficial : 29
tablesVersionLatest : 29
...
productType : unknown
md5Section7 : 4caa43e9fc2f9bd2ba637aed3af7dc72
section8Length : 4
valuesメソッド
数値データを含むmasked_arrayを返す．
In
grb = grbs(name='Temperature')[0]
grb.values
Out
masked_array(
data=[[--, --, --, ..., --, --, --],
[--, --, --, ..., --, --, --],
[--, --, --, ..., --, --, --],
...,
[291.39064025878906, 291.29689025878906, 291.23439025878906, ...,
--, --, --],
[291.48439025878906, 291.40626525878906, 291.34376525878906, ...,
--, --, --],
[291.59376525878906, 291.51564025878906, 291.48439025878906, ...,
--, --, --]],
mask=[[ True, True, True, ..., True, True, True],
[ True, True, True, ..., True, True, True],
[ True, True, True, ..., True, True, True],
...,
[False, False, False, ..., True, True, True],
[False, False, False, ..., True, True, True],
[False, False, False, ..., True, True, True]],
fill_value=9999.0)
lationsメソッド
緯度，経度データを返す．
In
grb = grbs(name='Temperature')[0]
grb.latlons()
Out
(array([[47.6 , 47.6 , 47.6 , ..., 47.6 , 47.6 , 47.6 ],
[47.56, 47.56, 47.56, ..., 47.56, 47.56, 47.56],
[47.52, 47.52, 47.52, ..., 47.52, 47.52, 47.52],
...,
[22.48, 22.48, 22.48, ..., 22.48, 22.48, 22.48],
[22.44, 22.44, 22.44, ..., 22.44, 22.44, 22.44],
[22.4 , 22.4 , 22.4 , ..., 22.4 , 22.4 , 22.4 ]]),
array([[120. , 120.05, 120.1 , ..., 149.9 , 149.95, 150. ],
[120. , 120.05, 120.1 , ..., 149.9 , 149.95, 150. ],
[120. , 120.05, 120.1 , ..., 149.9 , 149.95, 150. ],
...,
[120. , 120.05, 120.1 , ..., 149.9 , 149.95, 150. ],
[120. , 120.05, 120.1 , ..., 149.9 , 149.95, 150. ],
[120. , 120.05, 120.1 , ..., 149.9 , 149.95, 150. ]]))
dataメソッド
数値データ，緯度・経度データを返す．
出力結果はvaluesメソッドとlationsメソッドの結果を合わせたような感じ．
In
grb = grbs(name='Temperature')[0]
grb.data()
その他のメソッド
keyに入っている情報をメソッドとして取り出せるのでよく使いそうなものをまとめた．
| メソッド | 戻り値の内容 |
|---|---|
| .name | 変数名 |
| .shortName | 変数名の省略形 |
| .units | 単位 |
| .maximum | 最大値 |
| .minimum | 最小値 |
| .average | 平均値 |
| .validDate | 予報時刻(UTC) |
| .analDate | 解析基準時刻(UTC) |
| .forecastTime | 予報時間 |
| .missingValue | 欠損値 |
おまけ
ファイルの中にどのような変数があるのかを簡単に見るために，気になる要素を抽出して出力するようにした．
In
import pygrib
import pandas as pd
grbs = pygrib.open('F:/Z__C_RJTD_20171205010000_LFM_GPV_Rjp_L-pall_FH0100_grib2.bin')
name = []
s_name = []
unit = []
level = []
level_type = []
# 各リストに要素を追加
for grb in grbs:
name.append(grb.name)
s_name.append(grb.shortName)
unit.append(grb.units)
level.append(grb.level)
level_type.append(grb.typeOfLevel)
# listをDataFrameに変換
df = pd.DataFrame([s_name, name, unit, level, level_type])
df = df.T
df.columns=['shortName', 'name', 'units', 'level', 'typeOfLevel']
# 行と列をすべて表示させる
pd.set_option('display.max_rows', None)
pd.set_option('display.max_columns', None)
print(df)
Out
shortName name units level typeOfLevel
0 gh Geopotential height gpm 1000 isobaricInhPa
1 u U component of wind m s**-1 1000 isobaricInhPa
2 v V component of wind m s**-1 1000 isobaricInhPa
3 t Temperature K 1000 isobaricInhPa
4 w Vertical velocity Pa s**-1 1000 isobaricInhPa
...
87 gh Geopotential height gpm 100 isobaricInhPa
88 u U component of wind m s**-1 100 isobaricInhPa
89 v V component of wind m s**-1 100 isobaricInhPa
90 t Temperature K 100 isobaricInhPa
91 w Vertical velocity Pa s**-1 100 isobaricInhPa