Hfrhyu

I am using mosaic_rasters from gdalUtils to combine different raster files. More specifically, I am using the 12 tiles that cover Southern Africa from the well-known Hansen et al. (2013) forest map (1.1) that can be downloaded here: http://earthenginepartners.appspot.com/science-2013-global-forest/download_v1.1.html. I am using the tree cover, forest gain and forest loss data. All is working fine for the forest gain and forest loss files but I receive an error when combing the tree cover files (which are the largest).

Below you will find my code. Note that "tree_cover_files" is a vector that points towards the 12 tiles stored in a local folder. I am receiving the following error after some time:

Checking gdal_installation...
Scanning for GDAL installations...
Checking the gdalUtils_gdalPath option...
GDAL version 2.2.4
GDAL command being used: "C:OSGeo4W64bingdalbuildvrt.exe" -input_file_list "c:TempRtmpkRd30lfile4687024480a.txt" "c:TempRtmpkRd30lfile4684596162f.vrt" 
Checking gdal_installation...
Scanning for GDAL installations...
Checking the gdalUtils_gdalPath option...
GDAL version 2.2.4
GDAL command being used: "C:OSGeo4W64bingdal_translate.exe" -of "GTiff" "c:TempRtmpkRd30lfile4684596162f.vrt" "P:/globiom/Projects/ISWEL/data/forest/combined_tiles/tree_cover.tif"
Input file size is 120000, 1600000...10...20...30...40...
ERROR 1: TIFFFillStrip:Read error at scanline 39921; got 3204 bytes, expected 10047
ERROR 1: TIFFReadEncodedStrip() failed.
ERROR 1: P:/globiom/Projects/ISWEL/data/forest/tree_cover/Hansen_GFC2014_treecover2000_10S_030E.tif, band 1: IReadBlock failed at X offset 0, Y offset 39922
ERROR 1: GetBlockRef failed at X block offset 0, Y block offset 39922

I understand that this might be thread error, related to how the files are read and intermediate vrt file is constructed (https://github.com/mapnik/node-mapnik/issues/437 ) and could be solved by setting "VRT_SHARED_SOURCE" to 0. I tried to do this using setCPLConfigOption("VRT_SHARED_SOURCE", "0") in R but still receive the same error. Is this really the way to solve this or is the file perhaps corrupt (unlikely as this dataset is used by many people - I also downloaded them twice)? I hope somebody can give me advice on how to make this work.

library(purrr)
library(curl)

### FUNCTIONS
# Function to download data
# NB: need to add mode = "wb" for binary files otherwise the tif will be corrupted!!!!
download_f <- function(url, folder)
 print(url)
 dstfile <- file.path(proj_path, paste0(folder, basename(url)))
 if(!file.exists(dstfile)) download(url, mode = "wb", destfile = dstfile)


### DOWNLOAD DATA
#source: Hansen et al. (2013), http://science.sciencemag.org/content/342/6160/850
#Website: http://earthenginepartners.appspot.com/science-2013-global-forest/download_v1.1.html

# The files are at a 30m resolution and therefore very large. We only download the tiles for Southern Africa using the following
# NS and EWE combinations
NS <- c("00N", "10S", "20S", "30S")
EW <- c("010E", "020E", "030E")

NS_EW <- expand.grid(NS = NS, EW = EW, stringsAsFactors = F)
NS_EW <- paste(NS_EW$NS, NS_EW$EW, sep = "_")
NS_EW <- paste(NS_EW, collapse = "|")

## Treecover
# list of urls for treecover tiles
url <- "http://commondatastorage.googleapis.com/earthenginepartners-hansen/GFC2014/treecover2000.txt"

# Select relevant tiles
filenames = getURL(url)
filenames = strsplit(filenames, "n")[[1]]
tiles = grep(NS_EW, filenames, value = T)

# Set location
# Create "data" subfolder manually 
proj_path <- "C:/temp"

# download
walk(tiles, download_f, "data/")


### MOSAIC 
e <- extent(20, 40, -20, -10)
template <- raster(e)
proj4string(template) <- CRS("+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0") 
writeRaster(template, file = file.path(proj_path, "data/tree_cover.tif"), format="GTiff", overwrite = T)

# Create list of tif files
tree_cover_files <- list.files(file.path(proj_path, "data"), pattern = glob2rx("*.tif"), full.names = T)

setCPLConfigOption("VRT_SHARED_SOURCE", "0")
mosaic_rasters(gdalfile = tree_cover_files,
 dst_dataset = file.path(proj_path, "data/tree_cover.tif"), of="GTiff",
 verbose = T)