講義のページ/Seminar のバックアップ(No.11) - 講義のページ

[ トップ ] [ 新規 | | 検索 | 最終更新 | ヘルプ | ログイン ]

Top/講義のページ/Seminar

専攻演習(ゼミ) †

気象庁予報データ †

データの取得 †

京都大学生存圏研究所のデータアーカイブ

気象庁データの数値予報GPV
※ 予報モデルGPVのファイル名規則：気象業務支援センター/オンライン気象情報/ファイル形式データのファイル形式データ一覧表
※ ここの気象庁予報データは、研究教育の目的のために利用できます。

ファイル名の例

ファイル名の末尾	意味
GSM_GPV_Rjp_L-pall_FD0000-0312_grib2.bin	全球数値予報モデルGPV ／日本領域／気圧面についてのデータ／予報時間
MSM_GPV_Rjp_L-pall_FH00-15_grib2.bin	メソ数値予報モデルのGPV ／日本領域／気圧面についてのデータ／予報時間 00～15時間

扱い方 †

準備
```
※ 個人的にはLinuxを使うが、授業なので、Windows を使うことを前提にして書く。
```
- wgrib2
  気象庁予報データは、grib2 という形式のデータで提供されている。これを netCDF 形式に変換するために、wgrib2 コマンドが必要になる。
  - ダウンロード
    https://www.ftp.cpc.ncep.noaa.gov/wd51we/wgrib2/Windows10/v3.1.1/からwgirb2.exe と、dll で終わるファイルをすべてダウンロードする。
    このとき、セキュリティーの警告が出る。すべて、「保存」を選択して保存する。
- Rstudio
  - ダウンロードとインストール
    【2021年11月版】WindowsにR・RStudioをインストールするなどを参照しながらインストールする。
    RStudio だけでなく、R 自身もインストールすること。
  - ncdf4 ライブラリのインストール
    Rstudio を起動し、左下のコマンドプロンプトで
```
install.packages("ncdf4")
```
    と入力する。
  - 作業ディレクトリの特定
    Rstudio を起動すると、初期設定として、「ドキュメント」フォルダが利用される。他のファイルと混ざらないように、work ディレクトリを作成して、そこで作業する。
    Rstudio の右下の画面で、「New Folder」をクリックして、「work」と入力してOKボタンを押す。「work」フォルダが現れるので、そこをダブルクリックする。
    「New Folder」の右側に「More▼」があるので、これをクリックして、「Set As Working Directory」をクリックする。
  - wgrib2 の移動
    ダウンロードした wgrib2.exe と dll ファイルを、上述のドキュメント以下の work ディレクトリに移動する。
- データの取得と変換
  京都大学生存圏研究所では、気象庁の予報データなどのファイルを保存し、公開してくれている。教育・研究目的に利用できる。これをありがたく使わせていただく。
  - 気象庁の予報データ
    http://database.rish.kyoto-u.ac.jp/arch/jmadata/data/gpv/original/から特定の日時を選択してダウンロードする。ファイルは、先程の work ディレクトリに移動しておく。ファイル名の意味については、…後日
  - データ形式の変換
    いろいろと面倒なので、Rstudio を立ち上げた状態で実行する。
```
system("wgrib2 ダウンロードしたファイル.bin -netcdf 適当なファイル名.ncdf")
```
    いかにもそれっぽい表示（日付とか、データの変数名とかの表示）がされたら、変換が成功した証である。

データの扱い
R で具体的に処理を始める。

library(ncdf4)
nc <- nc_open("適当なファイル名.ncdf")
print(nc)

変数名(var) の一覧が表示されるので、その中から利用する変数を選ぶ。そのデータを読み込む。

gph1000 <- ncvar_get(nc, "HGT_1000mb")
lon <- ncvar_get(nc, "longitude")
lat <- ncvar_get(nc, "latitude")
tim <- as.POSIXct(ncvar_get(nc, "time"),origin="1970-01-01",tz="Asia/Tokyo")

入っている変数

MSM_GPV_Rjp_Lsurf

変数名	意味	単位
PRMSL_meansealevel	海面気圧	Pa
PRES_surface	現地気圧	Pa
UGRD_10maboveground	10m での東西風	m/s
VGRD_10maboveground	10m での南北風	m/s
TMP_1D5maboveground	1.5m での気温	K
RH_1D5maboveground	1.5m での相対湿度	%
LCDC_surface	下層雲の被覆率	%
MCDC_surface	中層雲の被覆率	%
HCDC_surface	上層雲の被覆率	%
TCDC_surface	全雲量	%
APCP_surface	降水量	kg/m^2
DSWRF_surface	下向き短波放射フラックス	W/m^2

MSM_GPV_Rjp_L-pall

変数名	意味	単位
HGT_500mb	500hPa ジオポテンシャル高度	m
UGRD_500mb	500hPa 東向き風速	m/s
VGRD_500mb	500hPa 北向き風速	m/s
TMP_500mb	500hPa 気温	K
VVEL_500mb	500hPa 鉛直速度	Pa/s
RH_500mb	500hPa 相対湿度	%

気圧面
1000, 975, 950, 925, 900, 850, 800, 700, 600, 500, 400, 300, 250, 200, 150, 100

作図例 †

水平面内の等値線図
```
filled.contour( lon, lat, gph1000[,,1] )
```
- 海岸線の追加
  水平面内で図を描く場合には、地図を入れると見通しが良くなる。
  - 準備
    次のコマンドを実行する。RStudio にパッケージがインストールされる。1回だけ行えば良い。
```
install.packages("devtools")
```
    次のコマンドを実行する。RStudio を起動するたびに1回行えば良い。
```
library(devtools)
source_url("https://www2.obirin.ac.jp/moriat/OpenData/addmap.R")
```
  - 水平面内の図に地図を入れる例
```
filled.contour( lon, lat, gph1000[,,1], plot.axes={axis(1); axis(2); addmap()} )
```
特定の地点の値の時間変化
経度lon は、0.125°ごと、緯度lat は 0.1° ごとに入っている。緯度と経度を指定して、値の時間変化を図にする。
```
plot( tim, gph1000[lon==140.125,lat==35.7,], type='b')
```
特定の経度の、緯度と時刻による値の変化
```
filled.contour( tim, lat, t(gph1000[lon==140.125,,]) )
```

データの加工 †

多次元配列の扱い

連続したデータを配列にする

sampledata <- 1:30
dim(sampledata) <- c(5,3,2)
sampledata[,,1]
sampledata[,,2]

多次元配列の組み換え（引数の場所の変更）

str(sampledata)
sampledata <- aperm(sampledata, c(1,3,2))
str(sampledata)
sampledata[,1,]
sampledata[,2,]

多次元配列を単なる1次元の配列に戻す
```
c(sampledata)
```

温位のデータを作る

ncvar_gettheta <- function(nc){
  levs <- c(1000, 975, 950, 925, 900, 850, 800, 700, 600, 500, 400, 300, 250, 200, 150, 100)

  lon <- ncvar_get(nc, "longitude")
  lat <- ncvar_get(nc, "latitude")
  tim <- as.POSIXct(ncvar_get(nc, "time"),origin="1970-01-01",tz="Asia/Tokyo")

  theta     <- c()
  validlevs <- c()
  for( lev in levs ){
    tmpdata <- c()
    try( tmpdata <- c(tmpdata, ncvar_get(nc,paste("TMP_",lev,"mb",sep=""))) )
    if( length(tmpdata) != 0 ){
      theta <- c( theta, c(tmpdata * (1000/lev)^(2/7)))
      validlevs <- c( validlevs, lev)
    }
  }
  dim(theta) <- c(length(lon), length(lat), length(tim), length(theta)/length(lon)/length(lat)/length(tim) )
  theta <- aperm(theta, c(1,2,4,3))
 
  return( list(lon=lon, lat=lat, levs=validlevs, tim=tim, theta=theta))
}

次のようにしてこの関数を取り込める。

library(devtools)
source_url("https://www2.obirin.ac.jp/moriat/OpenData/ncvar_gettheta.R")

緯度高度断面の等値線図を描く。返り値がリストなので、$ をつけて各変数を参照する。

result <- ncvar_gettheta(nc)
filled.contour( result$lat, -result$lev, result$theta[1,,,1]

※ filled.contour では、x, y 座標の値は増えていく順（昇順）でなければならない。気圧レベルは下ほど値が大きいので、負号をつけて上下を逆転させている。

関連文献 †

連成振動 †

　　　　　▼(ひらく)...

群速度 †

　　　　　▼(ひらく)...

2021年度受講生用ページ †

力学 †

トピックス †

円運動と曲率半径
曲率円に関するページ(外部へのリンク)

国家公務員試験解答例 †

過去の専攻演習(ゼミ)のページ †