Hi all,
I have the following issue with creating a ncdf file. Can someone suggest me
how to solve?
I have *.grd files for each year starting from 1980 to 2009. Each
file has daily rainfall values for lon: 66.5 N to 100.5 E, lat: 6.5
N to 38.5 N. Leap years have 366 days in respective files. In order
to proceed with my analysis, I am preparing a netcdf file which would
consist of daily dataset for the above said years (if you think there
is another simple way to compile every file into single file, let me
know of it too!). I have written the following R code (though not
complete yet), but the issue is that when I print variable written in
the netcdf file, I get "9.96921e+36" values in place of both "NA" and
some valid values (which I knew by comparing with the original file).
Here is the code...!
#####################################################
# define the netcdf coordinate variables -- note these have values! #
#####################################################
library(ncdf)
x <- dim.def.ncdf( "Lon", "degreesE",
seq(66.5,len=69, by=0.5), create_dimvar=TRUE)
y <- dim.def.ncdf( "Lat", "degreesN",
seq(6.5,len=65, by=0.5), create_dimvar=TRUE)
t <- dim.def.ncdf( "Time", "days since 1980-01-01", 1:366,
unlim=TRUE) #although I have not used "t" dimension now, as it
was kinda complex to start with
# define the EMPTY (RF) netcdf variable
origMissVal <- -1
varz = var.def.ncdf("Rainfall","millimeters", list(x,y), origMissVal,
longname="Daily Rainfall Data, 0.5 x 0.5 reso. (India)",
prec="double")
nc.ex = create.ncdf("rf.nc", varz, verbose=FALSE)
######################
# ReadTime function #
######################
readTime <- function(f,tm,x=69,y=65,size=4,type="double",NAvalue=-999){
c = file(f,open="rb")
seek(c,x*y*(tm-1)*size)
m=matrix(readBin(c,type,x*y,size),x,y)
#print(m)
close(c)
m[abs(m-NAvalue)<.Machine$double.eps] = NA
m
#print(m)
}
#############
# 1980.grd # (reading just the 1st day)
#############
tm <- 1
file <- "rf0.5_1980.grd"
arr1980 <- array(data=0, c(69, 65, tm))
for(i in 1:tm){
arr1980[,,i] <- matrix(readTime(file,i), 69, 65)
}
rm(tm)
#######################
# writing data
########################
nc.ex = create.ncdf("rf.nc", varz, verbose=FALSE)
put.var.ncdf(nc.ex, varz, arr1980)
close.ncdf(nc.ex)
########################################################
# checking the values of the nc file and "arr1980" from the original
file #
########################################################
nc.open <- open.ncdf("rf.nc")
v2 <- nc.open$var[[1]]
d.arr <- array(data=0, c(69, 65))
d.arr[,] <- get.var.ncdf(nc.open, v2)
close.ncdf(nc.open)
###################
# End
###################
Thank you for your time and consideration.
data file is here ...
http://www.4shared.com/file/5ftWLouJ/rf05_1980.html
Thanks again!
--
Regards,
Mahalakshmi
Graduate Student
[[alternative HTML version deleted]]
_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
