I actually got your code to work by making very slight modifications:
put.var.ncdf(nc.ex, varz, arr1980[,,1])
# reading it in with
nc.open <- open.ncdf("rf.nc")
get.var.ncdf(nc.open, "Rainfall")->dummy
image(dummy) # plots a map of india
On Sat, Oct 23, 2010 at 8:41 PM, <govin...@msu.edu> wrote:
>
>
> 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
>
>
[[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
