0

I'm trying to combine a large number of raster tiles to a single mosaic using R codes as follows. The error that appears is:

Error in if (xn == xx) { : missing value where TRUE/FALSE needed

The error appears after the for loop.

I will highly appreciate your suggestion.

require(raster)
rasters1 <- list.files("D:/lidar_grid_metrics/ElevMax", 
                   pattern="*.asc$", full.names=TRUE, recursive=TRUE)
rast.list <- list()
   for(i in 1:length(rasters1)) { rast.list[i] <- raster(rasters1[i]) }

rast.list$fun <- mean
    rast.mosaic <- do.call(mosaic,rast.list)
plot(rast.mosaic)
Ferdinando
  • 964
  • 1
  • 12
  • 23
Ram
  • 19
  • 4

1 Answers1

1

First a better way to write what you do (use lapply)

library(raster)
ff <- list.files("D:/lidar_grid_metrics/ElevMax", 
                   pattern="\\.asc$", full.names=TRUE, recursive=TRUE)
rast.list <- lapply(ff, raster)

rast.list$fun <- mean
rast.mosaic <- do.call(mosaic,rast.list)

Now, to the error your get. It is useful to show the results of traceback() after the error occurs. But from the error message you get, I infer that one of the RasterLayers has an extent with an NA value. That makes it invalid. You can check if that is true (and if so figure out what is going on) by doing 

t(sapply(rast.list, function(i) as.vector(extent(i))))

EDIT

With the files Ram send me I figured out what was going on. There was a bug when creating a RasterLayer from an ascii file with the native driver if the file specifies "xllcenter" rather than "xllcorner".

This is now fixed on the development version (2.9-1) available on github.

The problem can also be avoided by installing rgdal because if rgdal is available, the native driver won't be used.

Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63
  • Hi Rob, Thank you for your inputs. The same code works fine on another machine but not on the machine on which I work. The results of traceback() after the error are 19: .local(x, ...) 18: extent(object) 17: extent(object) 16: methods::validObject(extent(object)) 15: validityMethod(as(object, superClass)) 14: anyStrings(validityMethod(as(object, superClass))) 13: validObject(.Object) 12: initialize(value, ...) 10: methods::new("RasterLayer", extent = ext, nrows = nrows, ncols = ncols, crs = crs) 9: .local(...) – Ram Jan 28 '19 at 21:01
  • t(sapply(rast.list, function(i) as.vector(extent(i)))) returns [1, ] – Ram Jan 28 '19 at 21:03
  • Adding "rgdal" package to solved the problem. Thanks, – Ram Jan 28 '19 at 22:49
  • Ah, so these are ascii files that gdal can read, but the native raster code fails. --- it would be helpful if you could point me to (or send m)e one of these files – Robert Hijmans Jan 29 '19 at 04:16
  • Hi Rob, Sure I will send you one of the ascii files. How do you want me to send the file (via email)? – Ram Jan 30 '19 at 16:02
  • Yes, if the file is small enough. Else perhaps add a link to your post – Robert Hijmans Jan 31 '19 at 05:41