見出し画像

JRA-55を使った力学的圏界面の天気図

黒良 龍太

今回は、圏界面付近の気象状況を把握できる等渦位面天気図を作図します。

力学的圏界面について

力学的圏界面(渦位1.5PVUや2PVU面)における温位や風、気圧を描画した天気図を作成します。この天気図からJETや力学的圏界面の高度や、圏界面の折り込みなどを確認することで、対流圏における中緯度の上層擾乱などの構造把握の参考となります。下層擾乱との相互作用としての温帯低気圧の発達においては、圏界面下降などの特徴をこの天気図からみることができます。

等渦位面天気図に関する見方などについては、気象庁の総観気象学 応用編 第3章にある「力学的圏界面と等温位面解析」、や気象庁の量的予報研修テキストにある「予報作業における渦位の利用について」が参考になるでしょう。

MetPyの metpy.interpolate では、等圧面のデータから、等渦位面や等温位面のデータを簡単に計算してくれる関数が用意されています。JRA-55のデータを活用して、これらについて紹介します。なお、JRA-55のデータ読み込みなどの部分については、以前の記事で説明していますので、省いています。

画像1
今回作成する力学的圏界面(1.5PVU)の天気図 (2010年5月7日9時)
温位分布(シェード)、等風速線(白実線)、風(矢羽)、等圧線(赤実線)

等渦位面データ作成の準備

JRA-55の等圧面データを読み込んだ後、等渦位面データセットを作成するために等圧面のデータセットを準備します。このデータセットの中に、等圧面ではありますが、等渦位面の気圧を求めるために気圧の３次元配列も便宜上準備します。もちろんこの値は、等圧面と同じ気圧の一定値を代入することになります。

下のコードでは、最初の部分にその配列を作成し、等圧面データのデータセットを作成しています。

# 等渦位面の気圧を求めるために、気圧の3次元データ作成
cube_pre = np.zeros(dataHgt.values.shape)
for i in range(len(dataHgt['level'].values)):
   cube_pre[i,:,:] = dataHgt['level'].values[i] 
#
### Dataset作成                                                                                                                          
ds = xr.Dataset(
 {
   "Geopotential_height_isobaric": (["level", "lat", "lon"], dataHgt.values),
   "u-component_of_wind_isobaric": (["level", "lat", "lon"], dataUgrd.values),
   "v-component_of_wind_isobaric": (["level", "lat", "lon"], dataVgrd.values),
   "Temperature_isobaric": (["level", "lat", "lon"], dataTmp.values),
   "Pressure": (["level", "lat", "lon"], cube_pre),
 },
 coords={
   "time": dataHgt["time"].values,
   "level": dataHgt["level"].values,
   "lat": dataHgt["lat"].values,
   "lon": dataHgt["lon"].values,
 },
)
## 単位設定                                                                      
ds['Geopotential_height_isobaric'].attrs['units'] = 'meter'
ds['u-component_of_wind_isobaric'].attrs['units'] = 'm/s'
ds['v-component_of_wind_isobaric'].attrs['units'] = 'm/s'
ds['Temperature_isobaric'].attrs['units'] = 'K'
ds['Pressure'].attrs['units'] = 'hPa'
#ds['time'].attrs['units'] = ''                                                  
ds['level'].attrs['units'] = 'hPa'
ds['lat'].attrs['units'] = 'degrees_north'
ds['lon'].attrs['units'] = 'degrees_east'

次に、温位と渦位を計算し、この等圧面データセットに追加ます。計算にはmetpy.calc にある関数を利用しています。

等渦位面の気圧、温位、風を得るために、metpy.interpolate　にある、interpolate_to_isosurface関数を利用します。この関数は、指定した物理量が一定となる面における２次元データを、与えた任意の３次元データを内挿して計算してくれます。この関数には、２つの物理量(ここでは渦位と、等渦位面上で求める物理量)の３次元データを与え、作図したい面の物理量の定数（ここでは渦位1.5や2.0PVUなど)を指定します。面の高度を求めるときに地上から上にサーチするか、その逆かも指定できます。この関数の引数は具体的には次の通り指定します。
　①　作図する等値面の物理用3次元の値
　　　　等渦位面を扱う場合は渦位の3次元データを指定
　②　①の値が一定となる面上で求める物理量（３次元のデータセット）
　　　　等渦位面で求める物理量を指定（気圧や温位など）
　③   作図したい①の値が一定となる面のその値（1次元の値）
              力学的圏界面であれば、①は渦位で、③は2PVUなど。
　④   サーチの方向(省略可能)
　　　　Defaultは地上からサーチ、Falseとすると、Topからサーチ。

そして、得られた等渦位面の配列データをデータセット(dsDT)を作成しています。このデータセットを作ることができれば、これまでの記事の通り、天気図を比較的簡単に作成できます。描画のコードは省略します。

### 計算
# 温位
ds['PT'] = mpcalc.potential_temperature(ds['level'],ds['Temperature_isobaric'])
# 渦位
ds['PV'] = mpcalc.potential_vorticity_baroclinic(ds['PT'],ds['level'],
                                                ds['u-component_of_wind_isobaric'],
                                                ds['v-component_of_wind_isobaric'])
### 等渦位面データの作成
# 気圧
Pr_dt = minterpolate.interpolate_to_isosurface(ds['PV'].values,ds['Pressure'].values,const_pv)
# 温位
PT_dt = minterpolate.interpolate_to_isosurface(ds['PV'].values,ds['PT'].values,const_pv)
# 風 東西成分
U_dt = minterpolate.interpolate_to_isosurface(ds['PV'].values,
                                             ds['u-component_of_wind_isobaric'].values,const_pv)
# 風 南北成分
V_dt = minterpolate.interpolate_to_isosurface(ds['PV'].values,
                                             ds['v-component_of_wind_isobaric'].values,const_pv)
#
## Data Set of Dynamic Tropapause
dsDT = xr.Dataset(
 {
   "u-component_of_wind": (["lat", "lon"], U_dt),
   "v-component_of_wind": (["lat", "lon"], V_dt),
   "PotentialTemperature": (["lat", "lon"], PT_dt),
   "Pressure": (["lat", "lon"], Pr_dt),
 },
 coords={
   "time": dataHgt["time"].values,
   "level": const_pv,
   "lat": dataHgt["lat"].values,
   "lon": dataHgt["lon"].values,
 },
)
## 単位設定                                                                      
dsDT['u-component_of_wind'].attrs['units'] = 'm/s'
dsDT['v-component_of_wind'].attrs['units'] = 'm/s'
dsDT['PotentialTemperature'].attrs['units'] = 'K'
dsDT['Pressure'].attrs['units'] = 'hPa'
dsDT['level'].attrs['units'] = 'kelvin * meter ** 2 / kilogram / second'
dsDT['lat'].attrs['units'] = 'degrees_north'
dsDT['lon'].attrs['units'] = 'degrees_east'
#
## knotsへ変換
dsDT['u-component_of_wind'] = (dsDT['u-component_of_wind']).metpy.convert_units('knots')
dsDT['v-component_of_wind'] = (dsDT['v-component_of_wind']).metpy.convert_units('knots')
dsDT['wind_speed'] = mpcalc.wind_speed(dsDT['u-component_of_wind'], dsDT['v-component_of_wind'])

サンプルコード

下に、サンプルコード(Jupyter Notebook用)を掲載します。JRA-55のファイルを準備し、サンプルコードの、netCDF形式の保存先フォルダーや、年月日などの部分を書き換えることで、紹介した力学的圏界面天気図を作成することができます。是非、試してみてください。

次回からは、実況を把握するための天気図などを紹介する予定です。

この記事が気に入ったらサポートをしてみませんか？