compare tiff_predictor

[mdsumner]

compare with this, a rich example for comparing implemntation (and a pending potential GDAL PR to add tiff predictor to ZARR driver).

virtual-zarr/virtual-tiff#80

I'm caught between zaro and a 'wholecat' experiment, but this was yesterday's workflow, ghrsstt-refs.parquet is generated by rustycogs::tiff_refs()` from the ausantarctic source.coop COGs.

source("zaro-enriched-store.R")
source("enrich-rustycogs-parquet.R")
zstd_decompress <- function(raw, n_bytes) {
  buf <- arrow::buffer(raw)
  reader <- arrow::BufferReader$create(buf)
  stream <- arrow::CompressedInputStream$create(reader, "zstd")
  stream$Read(n_bytes)$data()
}

decode_tile <- function(path, offset, length, tile_w = 512L, tile_h = 512L, size = 2L) {
  raw <- fetch_bytes(path, offset, length)
  decoded <- zstd_decompress(raw, tile_w * tile_h * size)
  vals <- readBin(decoded, integer(), n = tile_w * tile_h, size = size,
                  signed = TRUE, endian = "little")
  mat <- matrix(vals, nrow = tile_h, ncol = tile_w, byrow = TRUE)
  ## tiff_predictor == 2
  for (i in seq_len(tile_h)) mat[i, ] <- cumsum(mat[i, ])
  mat
}
refs <- arrow::read_parquet("ghrsst-refs.parquet")
a <- enrich_ghrsst_mur(refs)
arrow::write_parquet(a, "a.parquet")

store <- EnrichedParquetStore("a.parquet")
path1 <- refs$path[1]
store@refs |>  dplyr::collect() |> dplyr::mutate(id = dplyr::row_number()) |>  dplyr::filter(tile_col == 4, tile_row == 13) |>  dplyr::filter(path == path1, ifd == 0)
row <- store@refs[928, ] |> dplyr::collect()
tt <- decode_tile(row$path, row$offset, row$length)
tt[tt > 20000] <- NA
#tt <- tt[nrow(tt):1, ]
library(terra)
info <- vapour::vapour_raster_info(sprintf("/vsicurl/%s", path))
g <- grout::tile_index(grout::grout(info$dimension, info$extent, c(512, 512)))
g |> dplyr::filter(tile_col == row$tile_col + 1, tile_row == row$tile_row + 1)
r <- rast(sprintf("vrt:///vsicurl/%s?ovr=2", path))
par(mfrow = c(2, 2))
ex <- c(-158, -153, 18, 24)
plot(crop(r, ext(ex)))
ximage::xcontour(tt, unlist(row |> dplyr::select(bbox_xmin, bbox_xmax, bbox_ymin, bbox_ymax)), add = TRUE)
plot(crop(r, ext(ex)))
plot(crop(r, ext(ex)))
ximage::ximage(tt, unlist(row |> dplyr::select(bbox_xmin, bbox_xmax, bbox_ymin, bbox_ymax)), add = TRUE, col = hcl.colors(24))

[mdsumner]

draft PR for GDAL ZARR driver (it works), just needs testing and some resolution for the fact that kerchunk/parquet zarr and v2 zarr seemingly has no support, so while GDAL could add it it's a weird cowboy format at this stage OSGeo/gdal@master...mdsumner:gdal:zarr-tiff-predictor-DRAFT