Global cropland extent (fractions) annual 2000-2022 at 250 m and 1 km

Global cropland extent annual for 2000-2022 based on the Potapov et al. (2021). Cropland defined as: land used for annual and perennial herbaceous crops for human consumption, forage (including hay), and biofuel. Perennial woody crops, permanent pastures, and shifting cultivation are excluded from the definition. The original 30-m resolution data (0/1 values) was interpolated from time-series 2003, 2007, 2011, 2015, 2019 to annual values 2000 to 2022 using linear interpolation. All values shown are in principle fractions 0-100%. The 30-m and 100-m resoluton images (COGs) are too large for Zenodo but you can access them from URLs in the filenames_openlandmap_cropland.txt file. See for example (drop the URL in QGIS):

• https://s3.eu-central-1.wasabisys.com/openlandmap/layers30m/cropland_glad.potapov.et.al_p_30m_s_20030101_20031231_go_epsg.4326_v20240624.tif (3.2GB)
• https://s3.eu-central-1.wasabisys.com/openlandmap/layers100m/cropland_glad.potapov.et.al_p_100m_s_20030101_20031231_go_epsg.4326_v20240624.tif (2.3GB)

Disclaimer: linear interpolation has limited accuracy and is basically only used to gap-fill the missing years. The remaining missing values in the maps can be ALL consider to be 0 value for cropland. A more detailed up-to-date cropland map of the world is provided by van Tricht et al., (2023), however only single year (2021) has been mapped at 10-m resolution within the WorldCereal project.

The temporal interpolation was implemented using terra package ii.e. using the following fuction:

library(terra)
y.l = c(2003, 2007, 2011, 2015, 2019)
out.years = 2000:2022
i = parallel::mclapply(y.l, function(x){system(paste0('gdal_translate Global_cropland_', x, '.vrt Global_cropland_', x, '.tif -co TILED=YES -co BIGTIFF=YES -co COMPRESS=DEFLATE -co ZLEVEL=9 -co BLOCKXSIZE=1024 -co BLOCKYSIZE=1024 -co NUM_THREADS=8 -co SPARSE_OK=TRUE -a_nodata 255 -scale 0 1 0 100 -ot Byte'))}, mc.cores = length(y.l))
## land mask at 1 deg (100x100km) ----
x = parallel::mclapply(y.l, function(x){system(paste0("gdal_translate Global_cropland_", x, ".tif Global_cropland_", x, "_1d.tif -tr 1 1 -r average -co BIGTIFF=YES -ot Byte -co NUM_THREADS=10"))}, mc.cores = length(y.l))
## 2 hrs
g1 = terra::rast(paste0("Global_cropland_", y.l, "_1d.tif"))
gs = sum(g1, na.rm=TRUE)
plot(gs)
gs.p <- as.polygons(gs, values = TRUE, extent=FALSE, dissolve=FALSE, na.rm=TRUE)
## Input layers:
r = terra::rast(paste0("Global_cropland_", y.l, ".tif"))
int.mc = function(r, tile, y.l, out.years=2000:2022){
  bb = paste(as.vector(ext(tile)), collapse = ".")
  if(any(!file.exists(paste0("./tmp/", out.years, "/Global_cropland_", out.years, "_", bb, ".tif")))){
    r.t = terra::crop(r, ext(tile))
    ## each tile is 16M pixels
    r.x = as.data.frame(r.t, xy=TRUE, na.rm=FALSE)
    rs = rowSums(r.x[,-c(1:2)], na.rm=TRUE)
    ## if sum is == 0 means no cropland throughout the time-series
    sel = which(rs>0)
    ## extract complete values:
    r.x0 = r.x[sel,-c(1:2)]
    r.x0[is.na(r.x0)] = 0
    ## interpolate between values:
    t1s = as.data.frame(t(apply(r.x0, 1, function(y){ try( approx(y.l, as.vector(y), xout=out.years, rule=2)$y ) })))
    ## write to GeoTIFFs
    t1s$x <- r.x$x[sel];   t1s$y <- r.x$y[sel]
    ## convert to RasterLayer:
    r.x = rast(t1s[,c("x","y",paste0("V", 1:length(out.years)))], type="xyz", crs="+proj=longlat +datum=WGS84 +no_defs")
    for(j in 1:length(out.years)){
      writeRaster(r.x[[j]], filename=paste0("./tmp/", out.years[j], "/Global_cropland_", out.years[j], "_", bb, ".tif"), gdal=c("COMPRESS=DEFLATE"), datatype='INT1U', NAflag=0, overwrite=FALSE)
    }
  }
}
## test it:
#int.mc(r, tile=gs.p[1000], y.l)
## run in parallel ----
## takes 12 hrs... 1TB RAM
i = parallel::mclapply(sample(1:length(gs.p)), function(x){try( int.mc(r, tile=gs.p[x], y.l) )}, mc.cores = 70)