proto for fasterize cell abstraction

[mdsumner]

## classic fasterize

library(fasterize)
library(raster)
library(stars)
r0 <- raster(extent(-180, 180, -90, 90), nrows = 360, ncols = 720, crs = "+init=epsg:4326")
s0 <- st_as_stars(st_bbox(c(xmin = -180, ymin = -90, xmax = 180, ymax = 90), crs = 4326), 
                  nx = 720, ny = 360) 

pols <- rnaturalearth::ne_countries(returnclass = "sf")
pols$row_index <- seq_len(nrow(pols))
## fudge this happy future day of cell abstraction
fudge <- fasterize(pols, r0, field = "row_index")

nc <- ncol(fudge)
l <- vector("list", nrow(fudge))
set_na_Inf <- function(x) {x[is.na(x)] <- Inf; x}
for (rownum in seq_along(l)) {
  rl <- rle(set_na_Inf(raster::extract(fudge, cellFromRow(fudge, rownum))))
  ## clean up NA
 # bad <- is.na(rl$values)
  ## catch start index of non-NA values
#  v <- rl$values
#  v[is.na(v)] <- 1
#  vstart <- cumsum(v)[!is.na(rl$values)]

    l[[rownum]] <- tibble::tibble(row = rownum, value = rl$values[!bad], length = rl$lengths[!bad])
}

library(dplyr)
index <- bind_rows(l)
pryr::object_size(index)
pryr::object_size(fudge)

set_inf_NA <- function(x) {x[!is.finite(x)] <- NA; x}
unpack_rle_stars <- function(rl, template) {
  template[[1]][] <- set_inf_NA(rep(rl$value, rl$length))
  template
}

unpack_rle_raster <- function(rl, template) {
  setValues(template, rep(rl$value, rl$length))
}

plot(unpack_rle_stars(index, s0))
plot(unpack_rle_raster(index, r0))

[mdsumner]

This way, not the above

Next to write a pack_raster() function (for this kind of demo, and for quantifying this encoding), and an unpack_raster() or something ...

library(fasterize)
#> 
#> Attaching package: 'fasterize'
#> The following object is masked from 'package:graphics':
#> 
#>     plot
library(raster)
#> Loading required package: sp
library(stars)
#> Loading required package: abind
#> Loading required package: sf
#> Linking to GEOS 3.7.0, GDAL 2.4.0, PROJ 5.2.0
r0 <- raster(extent(-180, 180, -90, 90), nrows = 360, ncols = 720, crs = "+init=epsg:4326")
s0 <- st_as_stars(st_bbox(c(xmin = -180, ymin = -90, xmax = 180, ymax = 90), crs = 4326), 
                  nx = 720, ny = 360) 

pols <- rnaturalearth::ne_countries(returnclass = "sf")
pols$row_index <- seq_len(nrow(pols))
## fudge this happy future day of cell abstraction
dummy <- fasterize(pols, r0, field = "row_index")

nc <- ncol(r0)
l <- vector("list", nrow(r0))
set_na_Inf <- function(x) {x[is.na(x)] <- Inf; x}

rownum <- 100

for (rownum in seq_along(l)) {
  rl <- rle(raster::extract(dummy, cellFromRow(dummy, rownum)))
  runs <- tibble::tibble(row = rownum, 
                         value = rl$values,
                         length = rl$lengths,
                         start = cumsum(length) - length + 1) %>% 
    dplyr::filter(!is.na(value))

    l[[rownum]] <- runs
}

library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:raster':
#> 
#>     intersect, select, union
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union
index <- bind_rows(l)
index$cell <- cellFromRowCol(dummy, index$row, index$start)
index$row <- index$start <- NULL
cells <- tibble::tibble(cell = 
                          unlist(purrr::map(purrr::transpose(index), ~seq(.x$cell, length.out = .x$length))), 
                        value = 
                          unlist(purrr::map(purrr::transpose(index), ~rep_len(.x$value, length.out = .x$length))))
rout <- dummy * NA_real_
rout[cells$cell] <- cells$value
plot(rout)

  (index_size <- pryr::object_size(index))
#> Registered S3 method overwritten by 'pryr':
#>   method      from
#>   print.bytes Rcpp
#> 81.3 kB
  (raster_size <- pryr::object_size(dummy))
#> 2.08 MB
  raster_size/index_size
#> 25.6 B

^{Created on 2019-06-05 by the reprex package (v0.3.0)}