[R-sig-Geo] reproject a multi-layered netcdf file using projectRaster with "nearest neighbor method"

[Yerguner]

I'm getting a curious error when trying to reproject a multi-layered netcdf
file using projectRaster() with the "nearest neighbor method" of the raster
package. But, interestingly, projectRasters() works fine when I set the
default "bilinear" method.
I need to use 'ngb' method since I want to keep the values same. And really
hope it can also reduce the process time since I've hundreds of files to
repeat this procedure and planning to use Rmpi.

Actually, another thing that I wonder is the performance of projectRaster
compared to gdalwarp.

Any help would be greatly appreciated.
Thanks,

Yasemin

library(raster)
library(rgdal)
library(ncdf)
rasterOptions(tmpdir=“/some_tmp_path/“)
input_file <- “cold2002_s1.nc"
cold02s1 <- brick(input_file)
newproj <- "+proj=laea +lat_0=45 +lon_0=-100 +x_0=0 +y_0=0 +a=6370997
+b=6370997 +units=m +no_defs"
prj1 <- projectExtent(cold02s1, newproj)
res(prj1) <- c(231.6564, 231.6563)
prj1 <- projectRaster(cold02s1, prj1, method='ngb')
Error in getfun(nc, varid = zvar, start = start, count = count) :
NAs in foreign function call (arg 3)
> traceback()
6: .C("R_nc_get_vara_double", as.integer(nc$id), as.integer(varid -
1), as.integer(c.start), as.integer(c.count), data = as.double(rv$data),
error = as.integer(rv$error), PACKAGE = "ncdf", DUP = FALSE)
5: getfun(nc, varid = zvar, start = start, count = count)
4: .readBrickCellsNetCDF(x, cells, layer, nl)
3: .cellValues(object, cells, layer = layer, nl = nl)
2: .xyValues(from, xy, method = method)
1: projectRaster(cold02s1, prj1, method = "ngb")

## I also tried: cold02s1 <- brick(input_file, varname="tmin") and got the
same error

##### here some info about the object and the brick
> prj1
class       : RasterLayer
dimensions  : 2044, 4457, 9110108  (nrow, ncol, ncell)
resolution  : 231.6564, 231.6563  (x, y)
extent      : -52180.77, 980311.8, 242137.5, 715643  (xmin, xmax, ymin,
ymax)
coord. ref. : +proj=laea +lat_0=45 +lon_0=-100 +x_0=0 +y_0=0 +a=6370997
+b=6370997 +units=m +no_defs

> cold02s1
class       : RasterBrick
dimensions  : 450, 990, 445500, 365  (nrow, ncol, ncell, nlayers)
resolution  : 1000, 1000  (x, y)
extent      : -50000, 940000, 496000, 946000  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=lcc +x_0=0 +y_0=0 +lon_0=-100 +lat_0=42.5 +lat_1=25
+lat_2=60 +ellps=WGS84
data source : /mypath/cold2002_s1.nc
names       : X2002.01.01, X2002.01.02, X2002.01.03, X2002.01.04,
X2002.01.05, X2002.01.06, X2002.01.07, X2002.01.08, X2002.01.09,
X2002.01.10, X2002.01.11, X2002.01.12, X2002.01.13, X2002.01.14,
X2002.01.15, ...
Date        : 2002-01-01, 2002-12-31 (min, max)
varname     : tmin

### I tried these and got the same error
prj1n <- projectRaster(cold02s1, prj1, method = "ngb")
Error in getfun(nc, varid = zvar, start = start, count = count) :
NAs in foreign function call (arg 3)
> prj1m <- projectRaster(cold02s1, crs = newproj, res= c(231.6564,
> 231.6563), method = "ngb")
Error in getfun(nc, varid = zvar, start = start, count = count) :
NAs in foreign function call (arg 3)

#### Interestingly when I try with the “bilinear” method which is default,
it works fine (it takes ~2 hrs).
prj1b <- projectRaster(cold02s1, prj1)
prj1b
class       : RasterBrick
dimensions  : 2044, 4457, 9110108, 365  (nrow, ncol, ncell, nlayers)
resolution  : 231.6564, 231.6563  (x, y)
extent      : -52180.77, 980311.8, 242137.5, 715643  (xmin, xmax, ymin,
ymax)
coord. ref. : +proj=laea +lat_0=45 +lon_0=-100 +x_0=0 +y_0=0 +a=6370997
+b=6370997 +units=m +no_defs
data source :
/pegasus/yb3/modvar2002/raster_tmp_2014-05-22_121316_13547_40735.grd
names       : X2002.01.01, X2002.01.02, X2002.01.03, X2002.01.04,
X2002.01.05, X2002.01.06, X2002.01.07, X2002.01.08, X2002.01.09,
X2002.01.10, X2002.01.11, X2002.01.12, X2002.01.13, X2002.01.14,
X2002.01.15, ...
min values  :    -29.5000,    -59.0000,    -86.0000,   -106.3880,  
-126.8480,   -150.4834,   -177.1940,   -198.6912,   -210.1884,   -218.1140,  
-231.6112,   -246.9922,   -278.9922,   -313.4882,   -344.3330, ...
max values  :    -13.5000,    -28.5000,    -41.5000,    -52.0000,   
-62.0000,    -76.0000,    -92.5000,   -104.0000,   -109.5000,   -114.5000,  
-121.5000,   -130.5000,   -143.5000,   -157.0000,   -171.0000, ...

#### my environment
> sessionInfo()
R version 3.0.2 (2013-09-25)
Platform: x86_64-redhat-linux-gnu (64-bit)
locale:
[1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C
[3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8
[5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8
[7] LC_PAPER=en_US.UTF-8       LC_NAME=C
[9] LC_ADDRESS=C               LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base
other attached packages:
[1] ncdf_1.6.6    rgdal_0.8-16  raster_2.2-31 sp_1.0-15
loaded via a namespace (and not attached):
[1] grid_3.0.2      lattice_0.20-23



--
View this message in context: http://r-sig-geo.2731867.n2.nabble.com/reproject-a-multi-layered-netcdf-file-using-projectRaster-with-nearest-neighbor-method-tp7586498.html
Sent from the R-sig-geo mailing list archive at Nabble.com.