CODEWITHSUNDEEP
rrasternetcdfnetcdf4
89_Simple
89_Simple

Reputation: 3805

Processing netcdf file in R with multiple variables

EDIT

I have edited my question since I made some progress.

I am trying to process the berkely earth data which is present in this website. http://berkeleyearth.org/data/. To get a sample data:

Scroll down to:

Daily Land (Experimental; 1880 - Recent) > Average High Temperature (TMAX) > 1º x 1º Latitude-Longitude Grid (each decade 200-450 MB) > 1990-1999

This netCDF file describes contains mean tmax and anomaly for each day rom 1990-1999 as a function of longitude and latitude. I also have a different set of lat lon for which I want to extract the daily mean tmax and anomaly data. This is how I am going about it. 

dat <- structure(list(locatioID = paste0('ID', 1:16), lon = c(73.73, 86, 73.45, 86.41, 85.36, 81.95, 82.57, 75.66, 82.03, 81.73, 85.66, 85.31, 81.03, 81.70, 87.03, 73.38), lat = c(24.59, 20.08, 22.61, 23.33, 23.99, 19.09, 18.85, 15.25, 26.78, 16.63, 25.98, 23.28, 24.5, 21.23, 25.08, 21.11)), 
             row.names = c(1L, 3L, 5L, 8L, 11L, 14L, 17L, 18L, 19L, 21L, 
                           23L, 26L, 29L, 32L, 33L, 35L), class = "data.frame")

xy <- dat[,2:3]
spts <- SpatialPoints(xy, proj4string=CRS("+proj=longlat +datum=WGS84"))

# read the netcdf file 
 ncname <- file.path(dirList$inputDir, 'climate','tmax','Complete_TMAX_Daily_LatLong1_1990.nc')
 ncin <- nc_open(ncname)
 print(ncin)

 # first I find the length of the datapoints in the netcdf 
 lengthvar <- ncin$var[[1]][['varsize']]
 # there are 3652 days in the netcdf

 # define the variables I need 
 dname1 <- 'temperature' # anomaly 
 dname2 <- 'climatology' # climatology 

 # load as a stack separately for climatology and anomaly 
  climtest <- brick(ncname, varname = 'climatology')
  anotest <- brick(ncname, varname = 'temperature')

# Now iterate through each day and extract the climatology and anomaly   

 for(i in 1:lengthvar){

    tempClim <- climtest[[i]]
    tempAno <- anotest[[i]]

    # Extract data by spatial point
     temp2 <- extract(tempClim, spts)
     temp2 <- extract(tempAno, spts)

      # Error in if (x@file@nodatavalue < 0) { : missing value where TRUE/FALSE needed
    }

I am not sure what does this error mean and what I am doing wrong here leading to this error.

Upvotes: 1

Views: 1170

Answers (1)

Spacedman
Spacedman

Reputation: 94202

If a raster layer is being kept on disk because of size and not stored in memory:

> inMemory(tempClim)
[1] FALSE

then extract fails:

> extract(tempClim, spts)
Error in if (x@file@nodatavalue < 0) { : 
  missing value where TRUE/FALSE needed

Ask R to read it in and it work:

> extract(readAll(tempClim), spts)
 [1] 23.83663 27.29292 28.61940 24.58966 23.86629 27.96306 27.10444 30.01295
 [9] 22.47180 29.46369 22.79476 23.86629 23.55748 27.75946 22.86436 30.12937

Although tempClim isn't large (180x360) its come from a large structure:

> dim(climtest)
[1] 180 360 365

and R is putting off reading from disk until you tell it to. By wrapping in readAll it does that. I don't know why extract doesn't print a better error message but I'm also getting deja vu here...

Upvotes: 1

Related Questions