r raster cover function - error with landsat stacks

Question

I'm running into an unusual issue with the cover function in R. We are trying to fill cloudy pixels with values from another layer. I can make it work just fine with stacks as follows -

library(raster)
r1 <- raster(ncols=36, nrows=18)
r1[] <- 1:ncell(r1)
r1b <- r1a <- r1
r1_stack <- stack(r1, r1a, r1b)

r2 <- setValues(r1, runif(ncell(r1)))
r2b <- r2a <- r2
r_stack <- stack(r2, r2a, r2b)

r_stack[r_stack < 0.5] <- NA

r3 <- cover(r_stack, r1_stack)

But then i try to do the same thing with a raster stack and i get the error:

Error in as.character(x) : 
 cannot coerce type 'closure' to vector of type 'character'

The code:

# get all tifs
LS5_032_032_2008_09_21 <- list.files("LT050340302008090301T1-SC20170526100900/", 
                                     pattern = glob2rx("*band*.tif$"), full.names = T) 


# stack bands
cloudy_scene <- stack(LS5_032_032_2008_09_21)
# import cloud mask
cloud_mask <- raster('LT050340302008090301T1-SC20170526100900/LT05_L1TP_034030_20080903_20160905_01_T1_sr_cloud_qa.tif')

# mask data 
masked_data <- mask(cloudy_scene, mask = cloud_mask, maskvalue=0, inverse=TRUE)

####### get cloud free data
# get files
LS5_2008_09_19 <- list.files("LT050340302008091901T1-SC20170526101124/", 
                             pattern = glob2rx("*band*.tif$"), full.names = T) 

# subset and stack cloud free bands
cloud_free_data <- stack(LS5_2008_09_19)

# use cover function to assign NA pixels to corresponding pixels in other scene
cover <- cover(masked_data, cloud_free_data)

TRACEBACK() output:

9: toupper(format)
8: .defaultExtension(format)
7: .getExtension(filename, filetype)
6: .local(x, filename, ...)
5: writeStart(outRaster, filename = filename, format = format, datatype = datatype, 
       overwrite = overwrite)
4: writeStart(outRaster, filename = filename, format = format, datatype = datatype, 
       overwrite = overwrite)
3: .local(x, y, ...)
2: cover(masked_data, cloud_free_data)
1: cover(masked_data, cloud_free_data)

UPDATE: I tried to resample the data - still doesn't work

cloud_free_resam <- resample(cloud_free_data, masked_data)
cover <- cover(masked_data, cloud_free_resam)

ERROR:

Error in as.character(x) : 
  cannot coerce type 'closure' to vector of type 'character'

I also tried to crop both layers - same error

# find intersection boundary
crop_extent <- intersect(extent(cloud_free_data), extent(masked_data))
cloud_free_data <- crop(cloud_free_data, crop_extent)
masked_data <- crop(masked_data, crop_extent)

# use cover function to assign NA pixels to corresponding pixels in other scene
cover <- cover(masked_data, cloud_free_data)

GET THE DATA: (WARNING: 317mb download - unpacks to ~1gb) https://ndownloader.figshare.com/files/8561230

Any ideas what might be causing this error with this particular dataset? I'm sure we're missing something quite basic but...what? Thank you in advance.

Leah

Answer 1

This is a bug. Cover with two multi-layered Raster* objects cannot write to disk. This can be seen in the simple example by setting

rasterOptions(todisk=TRUE)

I have fixed this in version 2-6.1 (forthcoming)

Answer 2

I think it is to do with when the object is converted to a rasterBrick and the raster written to a temporary file. i.e., masked_data <- mask(cloudy_scene, mask = cloud_mask)

Using Spacedman's crop creates a 'in memory' rasterBrick objects that do not rely on accessing the files; in that case the example works with no error.

But using the full extent (or cropping at end of process), raster writes and accesses (temporary) files and the error occurs.

Maybe as a temporary fix is to explicitly split up the images into memory sized chunks, mask/cover and then restitch with mosaic.

The extents of the two objects are also slightly different, so that should be fixed. Also generally a good idea to avoid issues by explicitly setting band names, min-max values and values < 0 to NA.