 登録日: 2004-7-29 居住地: 地球 投稿: 303 |
COWL Index の計算 (Gaussian Grid の場合のデータの取り扱い) NCEP/NCAR の再解析データを扱う上で、時々出てくるGaussian Grid 。面倒だな、と、思っていたら、何だ、自分で既に準備していたんじゃないか。 Gaussian Grid の場合のデータの取り扱いこれを元に、COWL Index を計算する。ただし、データは、/DATA/DL/NCEP_NCAR 以下に年毎にわかれて置かれているものとします。 1. 陸と海を別けるためのデータを準備する wgrib -d all -nh -o ~/landmask.bin 02/data/fixed/land.sfc
2. 各年のデータを準備する。 cd /DATA/DL/NCEP_NCAR/
for dir in ??
do
wgrib -nh -d all -o ~/${dir}tmpsfc.bin ${dir}/data/monthly/sfc/tmp.sfc
done
3. R で読み込む 3.1 ランドマスクを読み込み、cos(緯度)の重みをかける。 zz <- file("landmask.bin","rb")
landmask <- readBin(zz, "double", 192*94, 4 )
close(zz)
dim(landmask) <- c(192, 94)
landmask <- landmask[,seq(94,1,-1)]
lon <- seq(0, 360-360/192, 360/192)
lat <- c( -88.542,-86.653,-84.753,-82.851,-80.947,-79.043,-77.139,-75.235,-73.331,-71.426,
-69.522,-67.617,-65.713,-63.808,-61.903,-59.999,-58.094,-56.189,-54.285,-52.380,
-50.475,-48.571,-46.666,-44.761,-42.856,-40.952,-39.047,-37.142,-35.238,-33.333,
-31.428,-29.523,-27.619,-25.714,-23.809,-21.904,-20.000,-18.095,-16.190,-14.286,
-12.381,-10.476,-8.571,-6.667,-4.762,-2.857,-0.952,0.952,2.857,4.762,6.667,8.571,
10.476,12.381,14.286,16.190,18.095,20.000,21.904,23.809,25.714,27.619,29.523,
31.428,33.333,35.238,37.142,39.047,40.952,42.856,44.761,46.666,48.571,50.475,
52.380,54.285,56.189,58.094,59.999,61.903,63.808,65.713,67.617,69.522,71.426,
73.331,75.235,77.139,79.043,80.947,82.851,84.753,86.653,88.542)
landmaskweight <- t( matrix(cos(lat/180*pi), 94,1) %*% matrix( rep(1,192), 1, 192 ) ) * landmask
filled.contour(landmaskweight)
rm(landmask)
3.2 20 度以北で平均気温を計算する。 for( f in list.files(pattern="tmpsfc.bin") ){
zz <- file(f, "rb")
data <- readBin(zz,"double", 192*94*12, 4)
dim(data) <- c(192, 94, 12 )
data <- data[,seq(94,1,-1),]
close(zz)
# filled.contour(lon, lat, data[,,7] * landmaskweight )
for( month in 1:12 ){
write(
sum( data[,lat>20,month] * landmaskweight[,lat>20 ] ) /
sum( landmaskweight[, lat> 20] ) - 273.15,
file="COWL.txt", append=T)
}
}
※ これでは、00, 01, 02 の三年間がはじめに来てしまうことに注意。 3.3 北半球全部で計算する。 for( f in list.files(pattern="tmpsfc.bin") ){
zz <- file(f, "rb")
data <- readBin(zz,"double", 192*94*12, 4)
dim(data) <- c(192, 94, 12 )
data <- data[,seq(94,1,-1),]
close(zz)
# filled.contour(lon, lat, data[,,7] * landmaskweight )
for( month in 1:12 ){
write(
sum( data[,lat>=0,month] * landmaskweight[,lat>=0 ] ) /
sum( landmaskweight[, lat>=0] ) - 273.15,
file="COWLNH.txt", append=T)
}
}
※ これも、00, 01, 02 の三年間がはじめに来てしまうことに注意。 以上でおしまい。 |
メイン