見出し画像

JRA-55を使った地上天気図

JRA-55(NetCDF形式)を利用して、地上気圧や風、温位、降水量を天気図として作図するコードについて紹介します。このコードで作図した天気図を下に示します。1991年台風第19号(通称リンゴ台風)が長崎県に上陸間近の9月27日15時の天気図です。この時の中心気圧は935hPaでした!

この記事の最後に、サンプルコードを掲載しています。ご利用ください。

画像1
1991年9月27日15時のJRA-55を使った地上天気図

地上天気図

ご存知の通り、地上天気図は最も一般的な天気図です。高気圧や低気圧の位置や、等圧線の間隔や曲率などから天気や風をある程度推測できます。

天気図には、観測地点の観測値の表記や等圧線、高気圧、低気圧、台風、低圧部、前線が示されています。これに加えて、等温度線やドライライン、シアーライン、海上警報に関する情報などが解析されている天気図もあります。

気象庁の現業で行われている地上天気図については、この気象庁HPにある「気象庁予報現業で作成する天気図について」気象庁が作成する実況天気図の解説が参考になります。また、アメリカの天気図解析マニュアルの最後のページにドライラインなどの表記が示されています。

数値予報の予想地上天気図になると、風の予測が格子点上に示され、予想降水量分布など追加されます。アメリカの予想地上天気図は、降水種別がわかるように、雨、雪、みぞれを区別して分布が示されています。

JRA-55のデータについて

JRA-55には地上付近のデータとして、海面校正気圧や気温(地上2m高)、風(地上10m高)などの解析値があります。降水量は降水量を同化しておらず予報値です。実際に観測された値と大きく異なることがあることに留意しなければなりません。もちろん解析値にも誤差が含まれていますが、予報値よりは観測値との整合性がよくなっています。

ここで作成する天気図は、JRA-55のデータから、等圧線、温位(地上2m高)、高気圧と低気圧の位置、予想降水量分布、風を表示します。温位は前線解析する際の参考とするために表示させます。

JRA-55データの読み込み

最初にJRA-55のデータ読み込みを紹介したコンテンツ「JRA-55から天気図を作成する」では、xarrayのopen_dataset()を利用しましたが、今回は netCDF4クラスを利用して読み込みます。どちらを利用するかは、好みかと思います。私は多くの要素を読み込む時は、今回紹介するようなコードを書くと思います。

京都大学のサイトからダウンロードしているJRA-55(NetCDF形式)は、高層データでは月別ですが、地上付近のデータは年別になっています。目的の日時のデータの位置を調べるために、1月1日から何日目か計算し、6時間毎にデータがあることを考慮して算出しています。

import datetime
#
## 読み込みデータの指定
# JRA55の読み込む年月日時をUTCで与えます。
i_year =1964
i_month = 7
i_day = 16
i_hourZ = 18
#
## データ読み込み、データセット作成
# time番号  地上は年単位のNetCDF形式
d_dt =  datetime.datetime(i_year,i_month,i_day) - datetime.datetime(i_year,1,1)
time_no= d_dt.days * 4 + i_hourZ // 6

準備としてNetCDF形式のファイル名や要素名を指定して、利用する温位などの物理量と、座標である緯度・経度・時間を配列として読み込んでいきます。

import netCDF4

# データ格納先フォルダー名
#  PRMSL_msl_yyyy.nc  などの保存先
DataFd="./data/Jra55/"

# JRA55 faile name
## 平均海面校正気圧
nc_pre_fname= DataFd + "PRMSL_msl_{0:4d}.nc".format(i_year)
elem_pre="PRMSL_msl"
## 地上(2m)  温位
nc_pt_fname= DataFd + "POT_sfc_{0:4d}.nc".format(i_year)
elem_pt="POT_sfc"
## 地上風(10m) UGRD,VGRD
ncWu_fname= DataFd + "UGRD_fhg_{0:4d}.nc".format(i_year)
ncWv_fname= DataFd + "VGRD_fhg_{0:4d}.nc".format(i_year)
elem_wu="UGRD_fhg"
elem_wv="VGRD_fhg"
## 降水量
nc_prc_fname= DataFd + "APCP_sfc_{0:4d}.nc".format(i_year)
elem_prc="APCP_sfc"

## netcdfを読む
# 気圧
nc_pre  = netCDF4.Dataset(nc_pre_fname,"r")
dat_pre = nc_pre.variables[elem_pre][time_no,:,:] / 100.0 # Pa => hPa             
lat_pre = nc_pre.variables["lat"][:]
lon_pre = nc_pre.variables["lon"][:]
tim_pre = nc_pre.variables["time"][:]
nc_pre.close()
# 温位
nc_pt  = netCDF4.Dataset(nc_pt_fname,"r")
dat_pt = nc_pt.variables[elem_pt][time_no,:,:] # Kelvin                         
lat_pt = nc_pt.variables["lat"][:]
lon_pt = nc_pt.variables["lon"][:]
tim_pt = nc_pt.variables["time"][:]
nc_pt.close()
# 風
ncWu   = netCDF4.Dataset(ncWu_fname,"r")
ncWv   = netCDF4.Dataset(ncWv_fname,"r")
dat_u  = ncWu.variables[elem_wu][time_no,:,:] / 0.514444 # m/s => kt              
dat_v  = ncWv.variables[elem_wv][time_no,:,:] / 0.514444 # m/s => kt              
lat_w  = ncWu.variables["lat"][:]
lon_w  = ncWu.variables["lon"][:]
tim_w  = ncWu.variables["time"][:]
ncWu.close()
ncWv.close()

予報値の降水量は3時間毎にデータがあります。このため、データの位置は前述の温位などと異なるし、表示する降水量の期間が6時間以上の場合は積算することになります。ここでは6時間雨量を表示することとしています。

# 降水量の期間
prc_hours_target = 6    # R6を出力  3の倍数とすること!   

# 降水量
nc_prc  = netCDF4.Dataset(nc_prc_fname,"r")
# mm / day => mm / 3hour                                                          
dat_prc = nc_prc.variables[elem_prc][time_no * 2,:,:] / 8
# ループして3時間降水量から足す
prc_hours = 3
for i in np.arange(1, prc_hours_target / 3, 1):
 ii = time_no * 2 - i
 if ii < 0:
   break
 dat_prc = dat_prc + nc_prc.variables[elem_prc][ii,:,:] / 8
 prc_hours = prc_hours + 3
lat_prc = ncWu.variables["lat"][:]
lon_prc = ncWu.variables["lon"][:]
tim_prc = ncWu.variables["time"][:]
nc_prc.close()

降水量分布のハッチについて

これまでのコード説明と少し異なるところは、降水量のハッチの色指定でしょう。その部分を抜き出したコードを以下に抜き出します。

import cartopy.crs as ccrs
import matplotlib.pyplot as plt

## 図のSIZE指定inch 
fig_size = (10,8)
# 地図の描画範囲指定
i_area = [115, 151, 20, 50]  # 日本付近

## 図のSIZE指定
fig = plt.figure(figsize= fig_size)
## 余白設定                                                                                
plt.subplots_adjust(left=0, right=1, bottom=0.06, top=0.98)
## 図法指定
proj = ccrs.Stereographic(central_latitude=60, central_longitude=140)
latlon_proj = ccrs.PlateCarree()
## 図に関する設定
plt.rcParams["contour.negative_linestyle"] = 'solid'
## 作図開始                                                                                    
ax = fig.add_subplot(1, 1, 1, projection=proj)
ax.set_extent(i_area, latlon_proj)

## 降水量
# 降水量のハッチ
#   0-1mm:white, -10:cyan, -20:palegreen -30:lime -40:yellow
#   -60:orange -80:darkorange -100:red  100-:red
prc_hatchf = ax.contourf(lon_prc, lat_prc,
                        dat_prc, [0,1,10,20,30,40,60,80,100],
                        colors=['white','cyan','palegreen','lime','yellow',
                                'orange','darkorange','red'],
                        extend='both',
                        alpha=0.5, transform=latlon_proj)
# 10mm毎に等値線をかく
prc_hatch = ax.contour(lon_prc, lat_prc, dat_prc, np.arange(10,200,10),
                      colors='red', linewidths=1, linestyle='solid',
                      transform=latlon_proj)
ax.clabel(prc_hatch, fontsize=8, inline=True, colors='red',
         inline_spacing=5, fmt='%i', rightside_up=True)
# colorbarの位置と大きさ指定                                                     
#  add_axes([左端の距離, 下端からの距離, 横幅, 縦幅])                            
#ax_reld = fig.add_axes([0.1, 0.0, 0.8, 0.02])  # 図の下
ax_prc = fig.add_axes([0.1, 0.1, 0.8, 0.02])  # 図の中
cb_prc = fig.colorbar(prc_hatchf, orientation='horizontal', shrink=0.74,
                      aspect=40, pad=0.01, cax=ax_prc)

サンプルコード

地上天気図を作図するサンプルコード(Jupyter Notebook用)を次からダウンロードできます。当然のことながら、JRA-55のNetCDFデータが必要です。

次回は、事例解析を取り上げる予定です。






いいなと思ったら応援しよう!