diff --git a/.Rbuildignore b/.Rbuildignore index 7aac718..8d4141d 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -27,3 +27,4 @@ tests/testthat/Rplots.pdf ^CLAUDE\.md$ ^AUDIT\.md$ ^VISIBILITY\.md$ +^src/Makevars$ diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml index b923dd4..279453c 100644 --- a/.github/workflows/R-CMD-check.yaml +++ b/.github/workflows/R-CMD-check.yaml @@ -31,6 +31,7 @@ jobs: R_REMOTES_NO_ERRORS_FROM_WARNINGS: true GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }} KMP_DUPLICATE_LIB_OK: TRUE + _R_CHECK_FORCE_SUGGESTS_: false steps: - uses: actions/checkout@v4 @@ -42,39 +43,25 @@ jobs: - uses: r-lib/actions/setup-pandoc@v2 - - name: Query dependencies - run: | - install.packages("pak", repos = "https://r-lib.github.io/p/pak/dev/") - saveRDS(pak::pkg_deps("local::.", dependencies = TRUE), ".github/r-depends.rds") - shell: Rscript {0} - - - name: Cache R packages - if: runner.os != 'Windows' - uses: actions/cache@v4 - with: - path: ${{ env.R_LIBS_USER }} - key: ${{ matrix.config.os }}-${{ steps.install-r.outputs.installed-r-version }}-1-${{ hashFiles('.github/r-depends.rds') }} - restore-keys: ${{ matrix.config.os }}-${{ steps.install-r.outputs.installed-r-version }}-1- - - - name: Install curl - if: runner.os == 'Linux' - run: sudo apt-get install -y curl - - name: Install system dependencies + # libjpeg, libpng and libtiff are what the compiled image decoder links against; + # without them configure degrades that format to the R reader and the decoder is untested + - name: Install system dependencies (Linux) if: runner.os == 'Linux' run: | - pak::local_system_requirements(execute = TRUE) - pak::pkg_system_requirements("rcmdcheck", execute = TRUE) - shell: Rscript {0} + sudo apt-get update -y + sudo apt-get install -y libcurl4-openssl-dev libssl-dev libxml2-dev \ + libjpeg-dev libpng-dev libtiff5-dev \ + libfontconfig1-dev libharfbuzz-dev libfribidi-dev \ + libfreetype6-dev libgit2-dev - - name: Install dependencies - run: | - pak::local_install_dev_deps(upgrade = TRUE) - pak::pkg_install("rcmdcheck") - shell: Rscript {0} + - name: Install system dependencies (macOS) + if: runner.os == 'macOS' + run: brew install jpeg libpng libtiff pkg-config - - name: Install Pkg - run: | - Rscript -e "pak::local_install()" + - uses: r-lib/actions/setup-r-dependencies@v2 + with: + extra-packages: any::rcmdcheck + needs: check - name: Check I env: diff --git a/.gitignore b/.gitignore index ef31322..345c4e6 100644 --- a/.gitignore +++ b/.gitignore @@ -22,3 +22,6 @@ HISTORY.md INSIGHTS.md AUDIT.md VISIBILITY.md +src/Makevars +src/*.o +src/*.so diff --git a/DESCRIPTION b/DESCRIPTION index a643f6a..1806c2a 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -26,11 +26,13 @@ Description: The 'cito' package provides a user-friendly interface for training Encoding: UTF-8 LazyData: true Roxygen: list(markdown = TRUE) -RoxygenNote: 7.3.3 Depends: R (>= 4.1.0) Imports: + Rcpp, coro, + jpeg, + png, checkmate, torch, gridExtra, @@ -55,10 +57,12 @@ Suggests: ggplot2, knitr, tiff, - jpeg, - png, raster +LinkingTo: + Rcpp +SystemRequirements: libjpeg, libpng, libtiff (all optional; without them cito falls back to the R image readers) VignetteBuilder: knitr BugReports: https://github.com/citoverse/cito/issues URL: https://citoverse.github.io/cito/ Language: en-US +Config/roxygen2/version: 8.1.0 diff --git a/NAMESPACE b/NAMESPACE index 059f643..87d4848 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -62,3 +62,7 @@ export(transfer) export(tune) import(checkmate) import(tibble) +importFrom(Rcpp,sourceCpp) +importFrom(jpeg,readJPEG) +importFrom(png,readPNG) +useDynLib(cito, .registration = TRUE) diff --git a/NEWS.md b/NEWS.md index a8c5414..0644948 100644 --- a/NEWS.md +++ b/NEWS.md @@ -65,10 +65,16 @@ * conditionalEffects(data = ) works when the new data contains NAs. Those rows are dropped from the model frame, but the row indices were still drawn from the full data, which aborted the analytical method inside torch ("index out of bounds") and made the finite method fail on a length mismatch * the mvp likelihood returns one value per observation like every other likelihood, instead of a single reduced value * the development notes at the repository root (and E-CNN_and_MMN.md) are excluded from the build, they were shipped in the tarball before -* multi-band TIFFs that store their bands as separate planes (GDAL's INTERLEAVE=BAND) are no longer read incorrectly. tiff::readTIFF() returns wrong pixel values for that layout without raising an error, so models were trained on corrupted inputs. Such files are now refused with an explanatory message instead of being read silently * images read from a folder are no longer scaled down by a factor of 255. base_loader() returns values in [0, 1] and cito divided by 255 on top of that, so JPEG and PNG inputs reached the network in [0, 1/255] while TIFF inputs, which were not divided, reached it in [0, 1] -* folder inputs with .JPG, .PNG, .JPEG or .tif extensions are read instead of being rejected. The file type was matched case-sensitively, .tif was missing entirely, and the pattern was unanchored with an unescaped dot, so a file called "x.png.bak" was passed to the image reader while "photo.JPG" was not +* multi-band TIFFs that store their bands as separate planes (GDAL's INTERLEAVE=BAND) are no longer read incorrectly. tiff::readTIFF() returns wrong pixel values for that layout without raising an error, so models were trained on corrupted inputs. cito now reads such files through its compiled decoder, and the fallback reader refuses them with an explanatory message instead of continuing silently +* cloud-optimised GeoTIFFs can be read at all. A COG stores its overviews, and usually an internal mask, as further TIFF directories; cito passed every directory to readTIFF() and tried to stack images of different sizes and bit depths, which cannot succeed. Only full-resolution directories are used now +* palette and white-is-zero TIFFs are read by the compiled decoder instead of falling back; a palette image is expanded through its colormap, and white-is-zero returns the stored samples, which is what GDAL and tiff::readTIFF both do +* 32-bit float and signed 16-bit TIFFs are supported, which is what per-pixel embedding products and elevation models are distributed as. Their values pass through unscaled, since only unsigned integer samples carry an implied [0, 1] range +* folder images are decoded by a compiled reader when libjpeg, libpng and libtiff are available at install time, which is 9 to 14 times faster than the previous path and never materialises an R double array. Without those libraries cito falls back to the R readers, one format at a time * GeoTIFF tags no longer raise four warnings per file +* folder inputs with .JPG, .PNG, .JPEG or .tif extensions are read instead of being rejected. The file type was matched case-sensitively, .tif was missing entirely, and the pattern was unanchored with an unescaped dot, so a file called "x.png.bak" was passed to the image reader while "photo.JPG" was not +* RGBA images keep their alpha channel on every decode path +* cnn() gains num_workers for folder inputs. The dataset no longer stores tensors, which previously made torch's parallel dataloader fail with "external pointer is not valid" * weighted cross-validation in the hyperparameter tuning subsets matrix weights by row * print() of a citoarchitecture without an input_shape no longer errors * models fitted with cito <= 1.1 are rejected with an informative message instead of an internal error diff --git a/R/RcppExports.R b/R/RcppExports.R new file mode 100644 index 0000000..483c9a4 --- /dev/null +++ b/R/RcppExports.R @@ -0,0 +1,11 @@ +# Generated by using Rcpp::compileAttributes() -> do not edit by hand +# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393 + +cito_decoder_available <- function() { + .Call(`_cito_cito_decoder_available`) +} + +cito_decode_batch <- function(paths, nthreads = 1L) { + .Call(`_cito_cito_decode_batch`, paths, nthreads) +} + diff --git a/R/cito.R b/R/cito.R index cbff10a..a6758db 100644 --- a/R/cito.R +++ b/R/cito.R @@ -29,6 +29,10 @@ #' #' @example /inst/examples/dnn-example.R #' +#' @importFrom jpeg readJPEG +#' @importFrom png readPNG +#' @useDynLib cito, .registration = TRUE +#' @importFrom Rcpp sourceCpp #' @aliases cito-package #' @name cito "_PACKAGE" diff --git a/R/cnn.R b/R/cnn.R index fda254f..c32ebb9 100644 --- a/R/cnn.R +++ b/R/cnn.R @@ -16,6 +16,8 @@ #' @param validation Proportion of the data to be used for validation. Alternatively, a vector containing the indices of the validation samples can be provided. Default is 0.0. #' @param batchsize Batch size for training. If NULL, batchsize is 10% of the training data. Default is NULL. #' @param shuffle Whether to shuffle the data before each epoch. Default is TRUE. +#' @param num_workers (integer) Number of worker processes for the dataloader. Only effective for folder inputs. 0 decodes in the main process. +#' @param num_threads (integer) Number of threads the compiled image decoder uses per batch. Only effective for folder inputs. #' @param data_augmentation A list of functions used for data augmentation. Elements must be either functions or strings corresponding to inbuilt data augmentation functions. See details for more information. #' @param epochs Number of epochs to train the model. Default is 100. #' @param weights weights or other values (can be also a matrix) that the likelihood has access to @@ -108,6 +110,8 @@ cnn <- function(X, batchsize = NULL, shuffle = TRUE, data_augmentation = NULL, + num_workers = 0L, + num_threads = 4L, epochs = 100, weights = NULL, early_stopping = Inf, @@ -139,6 +143,8 @@ cnn <- function(X, checkmate::qassert(baseloss, c("0", "N1", "B1")) checkmate::qassert(device, "S+[3,)") checkmate::qassert(plot, "B1") + checkmate::qassert(num_workers, "X1[0,)") + checkmate::qassert(num_threads, "X1[1,)") checkmate::qassert(verbose, "B1") # Only return the model properties if no Y specified (Used in mmn()) @@ -223,7 +229,8 @@ cnn <- function(X, weights_torch <- if(is.null(weights)) NULL else torch::torch_tensor(weights) dls <- build_loaders(list(X), Y, weights = weights_torch, valid = valid, batch_size = batchsize, shuffle = shuffle, - from_folder = from_folder, data_augmentation = data_augmentation) + from_folder = from_folder, data_augmentation = data_augmentation, + num_workers = num_workers, num_threads = num_threads) train_dl <- dls$train_dl valid_dl <- dls$valid_dl @@ -249,6 +256,8 @@ cnn <- function(X, batchsize = batchsize, shuffle = shuffle, data_augmentation = data_augmentation, + num_workers = num_workers, + num_threads = num_threads, epochs = epochs, #redundant? early_stopping = early_stopping, burnin = burnin, @@ -345,7 +354,9 @@ predict.citocnn <- function(object, stop(paste0("Wrong dimension of newdata: [", paste(dim(newdata), collapse = ", "), "] Correct input dimension: [", paste(c("N", dim(object$data$X)[-1]), collapse = ", "), "]")) } - dl <- get_data_loader(newdata, batch_size = batchsize, shuffle = FALSE, from_folder = from_folder) + nt <- if(is.null(object$training_properties$num_threads)) 4L else object$training_properties$num_threads + dl <- get_data_loader(newdata, batch_size = batchsize, shuffle = FALSE, from_folder = from_folder, + num_threads = nt) pred <- NULL coro::loop(for(b in dl) { diff --git a/R/continue_training.R b/R/continue_training.R index ed46b19..df61a6f 100644 --- a/R/continue_training.R +++ b/R/continue_training.R @@ -269,8 +269,10 @@ continue_training.citocnn <- function(model, } } + nw <- if(is.null(model$training_properties$num_workers)) 0L else model$training_properties$num_workers + nt <- if(is.null(model$training_properties$num_threads)) 4L else model$training_properties$num_threads dls <- build_loaders(list(X), Y, weights = weights, valid = valid, batch_size = bs, shuffle = sh, - from_folder = from_folder) + from_folder = from_folder, num_workers = nw, num_threads = nt) train_dl <- dls$train_dl valid_dl <- dls$valid_dl @@ -348,8 +350,11 @@ continue_training.citommn <- function(model, } } + nw <- if(is.null(model$training_properties$num_workers)) 0L else model$training_properties$num_workers + nt <- if(is.null(model$training_properties$num_threads)) 4L else model$training_properties$num_threads dls <- build_loaders(X, Y, valid = valid, batch_size = model$training_properties$batchsize, - shuffle = model$training_properties$shuffle, from_folder = from_folder) + shuffle = model$training_properties$shuffle, from_folder = from_folder, + num_workers = nw, num_threads = nt) train_dl <- dls$train_dl valid_dl <- dls$valid_dl diff --git a/R/mmn.R b/R/mmn.R index d1424a5..93ded74 100644 --- a/R/mmn.R +++ b/R/mmn.R @@ -302,7 +302,9 @@ predict.citommn <- function(object, return(x) }) - dl <- do.call(get_data_loader, append(newdata, list(batch_size = batchsize, shuffle = FALSE, from_folder=from_folder))) + nt <- if(is.null(object$training_properties$num_threads)) 4L else object$training_properties$num_threads + dl <- do.call(get_data_loader, append(newdata, list(batch_size = batchsize, shuffle = FALSE, + from_folder = from_folder, num_threads = nt))) pred <- NULL coro::loop(for(b in dl) { diff --git a/R/utils.R b/R/utils.R index 822e9be..2fcc7a8 100644 --- a/R/utils.R +++ b/R/utils.R @@ -3,12 +3,21 @@ get_data_loader = function(..., shuffle=TRUE, from_folder = FALSE, data_augmentation = NULL, - transformations = NULL) { + transformations = NULL, + num_workers = 0L, + num_threads = 4L) { - if(from_folder | !is.null(data_augmentation)) ds = cito_dataset(..., data_augmentation = data_augmentation) - else ds <- torch::tensor_dataset(...) + if(from_folder | !is.null(data_augmentation)) { + ds = cito_dataset(..., data_augmentation = data_augmentation, worker_safe = num_workers > 0L, + num_threads = num_threads) + } else { + # a tensor_dataset holds tensors, which cannot be shipped to a worker process + num_workers = 0L + ds <- torch::tensor_dataset(...) + } - dl <- torch::dataloader(ds, batch_size = batch_size, shuffle = shuffle, pin_memory = TRUE) + dl <- torch::dataloader(ds, batch_size = batch_size, shuffle = shuffle, pin_memory = TRUE, + num_workers = num_workers) return(dl) } @@ -23,7 +32,7 @@ get_data_loader = function(..., # only to the training loader, never to the validation loader. build_loaders <- function(inputs, target, weights = NULL, valid = NULL, batch_size, shuffle, from_folder = FALSE, - data_augmentation = NULL) { + data_augmentation = NULL, num_workers = 0L, num_threads = 4L) { assemble <- function(idx, augment) { pick <- function(x) if(is.null(idx)) x else x[idx, drop = FALSE] @@ -32,7 +41,8 @@ build_loaders <- function(inputs, target, weights = NULL, valid = NULL, parts <- c(parts, list(pick(target))) do.call(get_data_loader, c(parts, list(batch_size = batch_size, shuffle = shuffle, - from_folder = from_folder, data_augmentation = augment))) + from_folder = from_folder, data_augmentation = augment, + num_workers = num_workers, num_threads = num_threads))) } if(is.null(valid)) { @@ -45,6 +55,21 @@ build_loaders <- function(inputs, target, weights = NULL, valid = NULL, valid_dl = assemble(valid, NULL)) } +# The compiled decoder dispatches on magic bytes, so the R readers must too: deciding by +# extension made an extensionless PNG readable on one path and not the other, which put the +# readability of a file at the mercy of which system library happened to be installed. +sniff_format = function(path) { + b = readBin(path, "raw", 4L) + if(length(b) < 4L) return("unknown") + if(identical(b[1:3], as.raw(c(0xFF, 0xD8, 0xFF)))) return("jpeg") + if(identical(b[1:4], as.raw(c(0x89, 0x50, 0x4E, 0x47)))) return("png") + if(identical(b[1:2], as.raw(c(0x49, 0x49))) && b[4] == as.raw(0x00) && + (b[3] == as.raw(42L) || b[3] == as.raw(43L))) return("tiff") + if(identical(b[1:2], as.raw(c(0x4D, 0x4D))) && b[3] == as.raw(0x00) && + (b[4] == as.raw(42L) || b[4] == as.raw(43L))) return("tiff") + "unknown" +} + # GeoTIFF tags make libtiff warn four times per file; those are not defects and would otherwise # drown a run over a few thousand tiles. Every other warning is left alone. quiet_tiff = function(expr) { @@ -54,8 +79,14 @@ quiet_tiff = function(expr) { } cito_dataset = torch::dataset( - initialize = function(..., data_augmentation) { - self$inputs = list(...) + initialize = function(..., data_augmentation, worker_safe = FALSE, num_threads = 4L) { + inputs = list(...) + # a torch tensor is an external pointer and cannot be serialized into a callr worker, + # so for num_workers > 0 the tensors are kept as R arrays and rebuilt per batch + if(worker_safe) inputs = lapply(inputs, function(x) if(inherits(x, "torch_tensor")) as.array(x$cpu()) else x) + self$inputs = inputs + self$decoder = cito_decoder_available() + self$num_threads = num_threads self$data_augmentation = data_augmentation }, .getbatch = function(index) { @@ -63,13 +94,16 @@ cito_dataset = torch::dataset( batch = lapply(self$inputs, function(x) { if(inherits(x, "torch_tensor")) { return(x[index, drop=FALSE]) + } else if(!is.character(x)) { + return(torch::torch_tensor(abind::asub(x, index, dims = 1, drop = FALSE))) } else { + batch <- self$decode_batch(x[index]) + if(!is.null(batch)) return(batch) X <- lapply(x[index], function(path) { - # anchored and case-insensitive: the old pattern missed .JPG/.PNG/.tif, which most - # cameras produce, and matched ".png" anywhere in the path, so "x.png.bak" got through - if(grepl("\\.(png|jpe?g)$", path, ignore.case = TRUE)) return(self$load_image(path)) - else if(grepl("\\.tiff?$", path, ignore.case = TRUE)) return(self$load_tiff(path)) - else stop(paste0("File format not supported: ", path)) + fmt <- sniff_format(path) + if(fmt == "tiff") self$load_tiff(path) + else if(fmt == "unknown") stop(paste0("File format not supported: ", path)) + else self$load_image(path, fmt) }) return(torch::torch_stack(X)) } @@ -92,27 +126,78 @@ cito_dataset = torch::dataset( return(batch) }, .length = function() { - if(inherits(self$inputs[[1]], "torch_tensor")) return(dim(self$inputs[[1]])[1]) - else return(length(self$inputs[[1]])) + x = self$inputs[[1]] + if(is.character(x)) return(length(x)) + return(dim(x)[1]) }, - load_image = function(path) { - img <- torch::torch_tensor(torchvision::base_loader(path), torch::torch_float32()) - if(img$ndim == 2) img <- img$unsqueeze(3) - img <- img$permute(c(3, 1, 2)) - # jpeg/png/tiff readers all return [0, 1] already, so dividing by 255 again scaled every - # folder input down by that factor and left TIFF inconsistent with JPEG/PNG - return(img) + # whole batch through the compiled decoder; NULL means "not applicable, use the R readers" + decode_batch = function(paths) { + if(!any(self$decoder)) return(NULL) + res <- cito_decode_batch(paths, self$num_threads) + # NULL means "the compiled path cannot read this batch", and the R readers take over. + # A genuine read failure is reported here rather than being rediscovered by the R reader. + if(!isTRUE(res$ok)) { + if(isTRUE(res$fallback)) return(NULL) + bad <- which(nzchar(res$errors)) + if(length(bad)) stop(paste0("cannot read ", paths[bad[1]], ": ", res$errors[bad[1]])) + return(NULL) + } + bad <- which(nzchar(res$errors)) + if(length(bad)) stop(paste0("cannot read ", paths[bad[1]], ": ", res$errors[bad[1]])) + dtype <- if(res$sformat == 2L) "float" else if(res$bits == 8L) "uint8" else "int16" + # torch_tensor_from_buffer() borrows the raw vector. $to(dtype) copies only when the dtype + # actually changes, so for a float32 file it is a no-op and the tensor would keep pointing + # at memory R is free to collect. clone() is the only guarantee. + t <- torch::torch_tensor_from_buffer(res$data, shape = res$dim, dtype = dtype)$clone() + t <- t$to(dtype = torch::torch_float32()) + # torch has no uint16, so unsigned 16-bit arrives signed and the top half wraps negative + if(res$sformat == 0L && res$bits == 16L) t <- torch::torch_where(t < 0, t + 65536, t) + t <- t$permute(c(1, 4, 2, 3)) + # only unsigned integer samples carry an implied [0, 1] scale; float embeddings and signed + # rasters such as elevation are physical quantities and must pass through untouched + if(res$sformat == 0L) t <- t/(2^res$bits - 1) + return(t) + }, + # the reference decode. A third path through magick was tried and removed: it quantised + # 16-bit PNG to 8 and disagreed with base_loader on grey+alpha, so which libraries happened + # to be installed decided how many input channels a model was fitted with + load_image = function(path, fmt) { + # base_loader() would be the obvious call, but it dispatches on the file extension and so + # cannot read a file the compiled path can. Reading by sniffed format keeps the two in step; + # the 2-D expansion is what base_loader does to a greyscale image. + img <- if(fmt == "png") png::readPNG(path) else jpeg::readJPEG(path) + if(length(dim(img)) == 2L) img <- array(rep(img, 3L), dim = c(dim(img), 3L)) + return(torch::torch_tensor(img, torch::torch_float32())$permute(c(3, 1, 2))) }, load_tiff = function(path) { if(!requireNamespace("tiff", quietly = TRUE)) stop("Package \"tiff\" must be installed to read .tiff files.", call. = FALSE) info <- quiet_tiff(tiff::readTIFF(path, info = TRUE, payload = FALSE)) - # readTIFF returns wrong pixels for separate planes without erroring, so this refuses - # instead of training on them - if(isTRUE(info$planar.config[1] == "separate")) - stop(paste0(path, " stores its bands as separate planes (PLANARCONFIG_SEPARATE), which ", - "tiff::readTIFF() decodes incorrectly. Rewrite the file with pixel interleave."), + # readTIFF returns wrong pixels for separate planes without erroring (AUDIT.md A10), so this + # refuses instead of training on them; the compiled decoder reads that layout correctly + # readTIFF decodes subsampled YCbCr past the end of its scanline buffer: repeated reads of + # the same file return different values, sometimes +/-Inf. Refuse rather than hand that on. + if(isTRUE(grepl("ycbcr", paste(info$color.space, collapse = ""), ignore.case = TRUE))) + stop(paste0(path, " is YCbCr-encoded, which tiff::readTIFF() decodes incorrectly. Install ", + "the libtiff development headers and reinstall cito, or rewrite the file as RGB."), call. = FALSE) - img_list <- quiet_tiff(tiff::readTIFF(path, all = TRUE)) + if(isTRUE(info$sample.format[1] == "int")) + stop(paste0(path, " stores signed samples, which tiff::readTIFF() rescales as if they were ", + "unsigned. Install the libtiff development headers and reinstall cito to read ", + "it with the compiled decoder."), call. = FALSE) + if(isTRUE(info$planar.config[1] == "separate")) + stop(paste0(path, " stores its bands as separate planes (PLANARCONFIG_SEPARATE), which the ", + "fallback TIFF reader decodes incorrectly. Install the libtiff development ", + "headers and reinstall cito to use the compiled decoder, or rewrite the file ", + "with pixel interleave."), call. = FALSE) + # A modern GeoTIFF carries more than the image: overviews at reduced resolution and, in + # essentially every COG, a 1-bit internal mask. readTIFF(all = TRUE) returns them all and + # then fails on the mask, so cito could read neither. Genuine pages share the first + # directory's geometry and bit depth; overviews differ in size, masks in depth. + dirs <- quiet_tiff(tiff::readTIFF(path, info = TRUE, payload = FALSE, all = TRUE)) + keep <- which(dirs$width == dirs$width[1] & dirs$length == dirs$length[1] & + dirs$bits.per.sample == dirs$bits.per.sample[1]) + img_list <- quiet_tiff(tiff::readTIFF(path, all = keep)) + if(!is.list(img_list)) img_list <- list(img_list) img_list <- lapply(img_list, function(img) { img <- torch::torch_tensor(img, dtype = torch::torch_float32()) if(img$ndim == 2) img <- img$unsqueeze(3) diff --git a/cleanup b/cleanup new file mode 100755 index 0000000..3c020d3 --- /dev/null +++ b/cleanup @@ -0,0 +1,2 @@ +#!/bin/sh +rm -f src/Makevars diff --git a/configure b/configure new file mode 100755 index 0000000..f30af18 --- /dev/null +++ b/configure @@ -0,0 +1,67 @@ +#!/bin/sh +# cito: the compiled image decoder is optional. If a library is missing the package still +# builds and falls back to the R readers, so a missing -dev package must never fail the install. +: "${R_HOME=$(R RHOME)}" +PKGCONFIG=$(command -v pkg-config 2>/dev/null) +# R's own compiler, because CRAN's macOS builds use clang with an -arch flag that a bare cc +# does not carry; unquoted on purpose, it is a command plus flags +PROBE_CC=$("${R_HOME}/bin/R" CMD config CC 2>/dev/null) +[ -n "$PROBE_CC" ] || PROBE_CC=${CC:-cc} + +CPPFLAGS="" +LIBS="" + +# compiles and links a probe; a header without a working library is the failure that would +# otherwise only surface when R tries to load the shared object +# The probe must *call* a symbol. A main() that references nothing links against a static +# archive without pulling a single member, so the archive's own unresolved dependencies are +# never seen -- and on macOS, where -undefined dynamic_lookup is in SHLIB_LDFLAGS, that +# produces a cito.so that links and then fails at dlopen, taking the whole package with it. +link_ok () { # $1 header, $2 call, $3 cflags, $4 libs + # jpeglib.h uses FILE but does not include stdio.h itself, so the probe must + echo "#include +#include +#include <$1> +int main(void){ $2; return 0; }" | $PROBE_CC -xc - $3 $4 -o /dev/null 2>/dev/null +} + +try_lib () { # $1 pkg-config name, $2 macro, $3 fallback -l flag, $4 header, $5 probe call + if [ -n "$PKGCONFIG" ] && $PKGCONFIG --exists "$1" 2>/dev/null; then + cf=$($PKGCONFIG --cflags "$1") + for lf in "$($PKGCONFIG --libs "$1")" "$($PKGCONFIG --libs --static "$1")"; do + if link_ok "$4" "$5" "$cf" "$lf"; then + CPPFLAGS="$CPPFLAGS $cf -D$2" + LIBS="$LIBS $lf" + echo "cito: found $1 via pkg-config" + return 0 + fi + done + echo "cito: $1 reported by pkg-config but does not link, trying the plain header" + fi + for d in /usr/include /usr/local/include /usr/include/$(uname -m)-linux-gnu \ + /opt/homebrew/include /opt/R/arm64/include /opt/R/x86_64/include; do + [ -f "$d/$4" ] || continue + libdir="${d%/include}/lib" + ldflag="" + [ -d "$libdir" ] && ldflag="-L$libdir" + if link_ok "$4" "$5" "-I$d" "$ldflag $3"; then + CPPFLAGS="$CPPFLAGS -I$d -D$2" + LIBS="$LIBS $ldflag $3" + echo "cito: found $4 in $d" + return 0 + fi + done + echo "cito: $1 not found, that format will use the R fallback reader" + return 1 +} + +try_lib libjpeg HAVE_JPEG -ljpeg jpeglib.h "struct jpeg_error_mgr e; jpeg_std_error(&e)" +try_lib libpng HAVE_PNG -lpng png.h "png_access_version_number()" +try_lib libtiff-4 HAVE_TIFF -ltiff tiffio.h "TIFFGetVersion()" + +esc () { printf '%s' "$1" | sed -e 's/[&|\\]/\\\\&/g'; } +sed -e "s|@CITO_CPPFLAGS@|$(esc "$CPPFLAGS")|" -e "s|@CITO_LIBS@|$(esc "$LIBS")|" \ + src/Makevars.in > src/Makevars +echo "cito: PKG_CPPFLAGS =$CPPFLAGS" +echo "cito: PKG_LIBS =$LIBS" +exit 0 diff --git a/man/cito.Rd b/man/cito.Rd index ea8bd76..e2859be 100644 --- a/man/cito.Rd +++ b/man/cito.Rd @@ -199,7 +199,8 @@ Useful links: Authors: \itemize{ - \item Christian Amesder \email{Christian.Amesoeder@informatik.uni-regensburg.de} + \item Maximilian Pichler \email{maximilian.pichler@biologie.uni-regensburg.de} (\href{https://orcid.org/0000-0003-2252-8327}{ORCID}) + \item Christian Amesöder \email{Christian.Amesoeder@informatik.uni-regensburg.de} } Other contributors: diff --git a/man/cnn.Rd b/man/cnn.Rd index f7d6710..b717850 100644 --- a/man/cnn.Rd +++ b/man/cnn.Rd @@ -21,6 +21,8 @@ cnn( batchsize = NULL, shuffle = TRUE, data_augmentation = NULL, + num_workers = 0L, + num_threads = 4L, epochs = 100, weights = NULL, early_stopping = Inf, @@ -62,6 +64,10 @@ cnn( \item{data_augmentation}{A list of functions used for data augmentation. Elements must be either functions or strings corresponding to inbuilt data augmentation functions. See details for more information.} +\item{num_workers}{(integer) Number of worker processes for the dataloader. Only effective for folder inputs. 0 decodes in the main process.} + +\item{num_threads}{(integer) Number of threads the compiled image decoder uses per batch. Only effective for folder inputs.} + \item{epochs}{Number of epochs to train the model. Default is 100.} \item{weights}{weights or other values (can be also a matrix) that the likelihood has access to} diff --git a/src/Makevars.in b/src/Makevars.in new file mode 100644 index 0000000..14adae1 --- /dev/null +++ b/src/Makevars.in @@ -0,0 +1,3 @@ +PKG_CPPFLAGS = @CITO_CPPFLAGS@ +PKG_CXXFLAGS = $(SHLIB_OPENMP_CXXFLAGS) +PKG_LIBS = @CITO_LIBS@ $(SHLIB_OPENMP_CXXFLAGS) diff --git a/src/Makevars.win b/src/Makevars.win new file mode 100644 index 0000000..15ba952 --- /dev/null +++ b/src/Makevars.win @@ -0,0 +1,6 @@ +# Rtools ships static libraries, so the dependency closure has to come from --static: +# a bare --libs returns -ltiff alone and the link fails on its transitive dependencies +PKG_CONFIG = pkg-config +PKG_CPPFLAGS = -DHAVE_JPEG -DHAVE_PNG -DHAVE_TIFF $(shell $(PKG_CONFIG) --cflags libtiff-4 libjpeg libpng) +PKG_CXXFLAGS = $(SHLIB_OPENMP_CXXFLAGS) +PKG_LIBS = $(shell $(PKG_CONFIG) --libs --static libtiff-4 libjpeg libpng) $(SHLIB_OPENMP_CXXFLAGS) diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp new file mode 100644 index 0000000..695d0e8 --- /dev/null +++ b/src/RcppExports.cpp @@ -0,0 +1,45 @@ +// Generated by using Rcpp::compileAttributes() -> do not edit by hand +// Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393 + +#include + +using namespace Rcpp; + +#ifdef RCPP_USE_GLOBAL_ROSTREAM +Rcpp::Rostream& Rcpp::Rcout = Rcpp::Rcpp_cout_get(); +Rcpp::Rostream& Rcpp::Rcerr = Rcpp::Rcpp_cerr_get(); +#endif + +// cito_decoder_available +Rcpp::LogicalVector cito_decoder_available(); +RcppExport SEXP _cito_cito_decoder_available() { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + rcpp_result_gen = Rcpp::wrap(cito_decoder_available()); + return rcpp_result_gen; +END_RCPP +} +// cito_decode_batch +Rcpp::List cito_decode_batch(Rcpp::CharacterVector paths, int nthreads); +RcppExport SEXP _cito_cito_decode_batch(SEXP pathsSEXP, SEXP nthreadsSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< Rcpp::CharacterVector >::type paths(pathsSEXP); + Rcpp::traits::input_parameter< int >::type nthreads(nthreadsSEXP); + rcpp_result_gen = Rcpp::wrap(cito_decode_batch(paths, nthreads)); + return rcpp_result_gen; +END_RCPP +} + +static const R_CallMethodDef CallEntries[] = { + {"_cito_cito_decoder_available", (DL_FUNC) &_cito_cito_decoder_available, 0}, + {"_cito_cito_decode_batch", (DL_FUNC) &_cito_cito_decode_batch, 2}, + {NULL, NULL, 0} +}; + +RcppExport void R_init_cito(DllInfo *dll) { + R_registerRoutines(dll, NULL, CallEntries, NULL, NULL); + R_useDynamicSymbols(dll, FALSE); +} diff --git a/src/decode.cpp b/src/decode.cpp new file mode 100644 index 0000000..b4ebb12 --- /dev/null +++ b/src/decode.cpp @@ -0,0 +1,526 @@ +// One decoder for every format cito reads. Decodes into a single pre-allocated buffer at the +// file's native bit depth, so no R double array is ever built. Errors are collected per file and +// raised on the R side: nothing in here may call the R API, because it runs under OpenMP. +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +// each decoder is optional: configure defines the macro only when the library was found, +// so a machine without the -dev packages still builds and falls back to the R readers +extern "C" { +#ifdef HAVE_JPEG +#include +#endif +#ifdef HAVE_PNG +#include +#endif +#ifdef HAVE_TIFF +#include +#endif +} + +enum Fmt { FMT_UNKNOWN = 0, FMT_JPEG, FMT_PNG, FMT_TIFF }; + +struct Meta { + int h, w, c, bits, pages, fallback, sformat; // sformat 0 uint, 1 int, 2 float + Fmt fmt; + int status; // 0 ok, non-zero means msg is set + char msg[220]; +}; + +static void set_err(Meta& m, const char* what) { + m.status = 1; + std::snprintf(m.msg, sizeof(m.msg), "%s", what); +} + +// returns FMT_UNKNOWN for a readable file of unknown type; *openable tells the two apart +static Fmt sniff(const char* path, bool* openable) { + FILE* f = std::fopen(path, "rb"); + *openable = (f != NULL); + if(!f) return FMT_UNKNOWN; + unsigned char b[8] = {0}; + size_t n = std::fread(b, 1, 8, f); + std::fclose(f); + if(n < 4) return FMT_UNKNOWN; + if(b[0] == 0xFF && b[1] == 0xD8 && b[2] == 0xFF) return FMT_JPEG; + if(n >= 8 && b[0] == 0x89 && b[1] == 'P' && b[2] == 'N' && b[3] == 'G') return FMT_PNG; + // classic TIFF (42) and BigTIFF (43), both byte orders; libtiff handles either + if((b[0] == 'I' && b[1] == 'I' && (b[2] == 42 || b[2] == 43) && b[3] == 0) || + (b[0] == 'M' && b[1] == 'M' && b[2] == 0 && (b[3] == 42 || b[3] == 43))) return FMT_TIFF; + return FMT_UNKNOWN; +} + +#ifdef HAVE_JPEG +/* ---- libjpeg: default error_exit calls exit(), so every call gets its own setjmp ---- */ +struct jerr_mgr { struct jpeg_error_mgr pub; jmp_buf jb; char msg[JMSG_LENGTH_MAX]; }; +static void jerr_exit(j_common_ptr ci) { + jerr_mgr* e = reinterpret_cast(ci->err); + (*ci->err->format_message)(ci, e->msg); + std::longjmp(e->jb, 1); +} +static void jerr_silent(j_common_ptr, int) {} + +// no non-trivial C++ object may be live across the setjmp, hence the raw FILE* and manual cleanup +static int jpeg_probe_decode(const char* path, Meta& m, unsigned char* dst) { + struct jpeg_decompress_struct ci; + jerr_mgr jerr; + FILE* f = std::fopen(path, "rb"); + if(!f) { set_err(m, "cannot open file"); return 1; } + ci.err = jpeg_std_error(&jerr.pub); + jerr.pub.error_exit = jerr_exit; + jerr.pub.emit_message = jerr_silent; + if(setjmp(jerr.jb)) { + jpeg_destroy_decompress(&ci); + std::fclose(f); + set_err(m, jerr.msg); + return 1; + } + jpeg_create_decompress(&ci); + jpeg_stdio_src(&ci, f); + jpeg_read_header(&ci, TRUE); + // libjpeg cannot convert CMYK or YCCK to RGB, and asking it to fails only at + // jpeg_start_decompress, by which point the batch buffer is already sized from the probe + if(ci.num_components != 1 && ci.num_components != 3) { + const int nc = ci.num_components; + jpeg_destroy_decompress(&ci); + std::fclose(f); + (void)nc; + m.fallback = 1; // CMYK/YCCK: jpeg::readJPEG reads it, so hand the batch over + return 1; + } + // base_loader() expands greyscale to three channels, so the compiled path must too: + // otherwise the first conv layer depends on which libraries were present at install time + ci.out_color_space = JCS_RGB; + if(dst == NULL) { + m.h = ci.image_height; + m.w = ci.image_width; + m.c = 3; + m.bits = 8; + jpeg_destroy_decompress(&ci); + std::fclose(f); + return 0; + } + if((int)ci.image_height != m.h || (int)ci.image_width != m.w) { + jpeg_destroy_decompress(&ci); + std::fclose(f); + set_err(m, "file changed between probe and decode"); + return 1; + } + jpeg_start_decompress(&ci); + const int stride = ci.output_width * ci.output_components; + while(ci.output_scanline < ci.output_height) { + unsigned char* row = dst + (size_t)ci.output_scanline * stride; + jpeg_read_scanlines(&ci, &row, 1); + } + jpeg_finish_decompress(&ci); + jpeg_destroy_decompress(&ci); + std::fclose(f); + return 0; +} + +#endif /* HAVE_JPEG */ + +#ifdef HAVE_PNG +/* ---- libpng ---- */ +#ifndef PNG_SETJMP_SUPPORTED +#error "cito needs a libpng built with setjmp support" +#endif +// libpng's default error function prints to stderr and leaves the caller with no message, +// so the handler records what failed and only then longjmps +struct png_err_ctx { char msg[200]; }; +static void cito_png_error(png_structp p, png_const_charp msg) { + png_err_ctx* c = (png_err_ctx*) png_get_error_ptr(p); + if(c) std::snprintf(c->msg, sizeof(c->msg), "%s", msg); + std::longjmp(png_jmpbuf(p), 1); +} +static void cito_png_warn(png_structp, png_const_charp) {} + +static bool host_is_little_endian() { + const uint16_t one = 1; + unsigned char b; + std::memcpy(&b, &one, 1); + return b == 1; +} + +static int png_probe_decode(const char* path, Meta& m, unsigned char* dst) { + png_err_ctx ctx; + ctx.msg[0] = '\0'; + FILE* f = std::fopen(path, "rb"); + if(!f) { set_err(m, "cannot open file"); return 1; } + png_structp p = png_create_read_struct(PNG_LIBPNG_VER_STRING, &ctx, cito_png_error, cito_png_warn); + if(!p) { std::fclose(f); set_err(m, "png_create_read_struct failed"); return 1; } + png_infop info = png_create_info_struct(p); + if(!info) { png_destroy_read_struct(&p, NULL, NULL); std::fclose(f); set_err(m, "png_create_info_struct failed"); return 1; } + png_bytep* volatile rows = NULL; + if(setjmp(png_jmpbuf(p))) { + if(rows) std::free(rows); + png_destroy_read_struct(&p, &info, NULL); + std::fclose(f); + set_err(m, ctx.msg[0] ? ctx.msg : "libpng error while reading"); + return 1; + } + png_init_io(p, f); + png_read_info(p, info); + png_uint_32 w = 0, h = 0; + int depth = 0, ctype = 0; + png_get_IHDR(p, info, &w, &h, &depth, &ctype, NULL, NULL, NULL); + if(ctype == PNG_COLOR_TYPE_PALETTE) png_set_palette_to_rgb(p); + if(ctype == PNG_COLOR_TYPE_GRAY && depth < 8) png_set_expand_gray_1_2_4_to_8(p); + const bool has_trns = png_get_valid(p, info, PNG_INFO_tRNS) != 0; + if(has_trns) png_set_tRNS_to_alpha(p); + // base_loader() expands a 2-D grey array to three channels but leaves grey+alpha at two. + // tRNS has just turned GRAY into GRAY_ALPHA, so testing the original colour type would + // expand it as well and give four channels where the R reader gives two. + if(ctype == PNG_COLOR_TYPE_GRAY && !has_trns) png_set_gray_to_rgb(p); + // libpng wants the byte-order transform before update_info, not after + if(depth == 16 && host_is_little_endian()) png_set_swap(p); + png_read_update_info(p, info); + const int channels = png_get_channels(p, info); + const int out_depth = (depth == 16) ? 16 : 8; + if(dst == NULL) { + m.h = (int)h; m.w = (int)w; m.c = channels; m.bits = out_depth; + png_destroy_read_struct(&p, &info, NULL); + std::fclose(f); + return 0; + } + if((int)h != m.h || (int)w != m.w || channels != m.c || out_depth != m.bits) { + png_destroy_read_struct(&p, &info, NULL); + std::fclose(f); + set_err(m, "file changed between probe and decode"); + return 1; + } + const size_t rowbytes = png_get_rowbytes(p, info); + rows = (png_bytep*) std::malloc(sizeof(png_bytep) * h); + if(!rows) { png_destroy_read_struct(&p, &info, NULL); std::fclose(f); set_err(m, "out of memory"); return 1; } + for(png_uint_32 r = 0; r < h; r++) rows[r] = dst + (size_t)r * rowbytes; + png_read_image(p, rows); + std::free(rows); + rows = NULL; + png_destroy_read_struct(&p, &info, NULL); + std::fclose(f); + return 0; +} + +#endif /* HAVE_PNG */ + +#ifdef HAVE_TIFF +/* ---- libtiff: the error handler is global, so the message lands in a thread-local slot ---- */ +#if defined(__GNUC__) || defined(__clang__) +#define CITO_TLS __thread +#else +#define CITO_TLS thread_local +#endif +static CITO_TLS char tiff_msg[200]; +static void tiff_err_handler(const char*, const char* fmt, va_list ap) { + std::vsnprintf(tiff_msg, sizeof(tiff_msg), fmt, ap); +} + +struct TiffCloser { void operator()(TIFF* t) const { if(t) TIFFClose(t); } }; + +static int tiff_probe_decode(const char* path, Meta& m, unsigned char* dst) { + tiff_msg[0] = '\0'; + // libtiff reports through return codes rather than longjmp, so RAII is safe here and + // removes six hand-written TIFFClose calls from the error paths + std::unique_ptr tp(TIFFOpen(path, "r")); + TIFF* t = tp.get(); + if(!t) { set_err(m, tiff_msg[0] ? tiff_msg : "cannot open TIFF"); return 1; } + uint32_t w = 0, h = 0; + uint16_t spp = 1, bps = 8, fmtfield = SAMPLEFORMAT_UINT, planar = PLANARCONFIG_CONTIG, photo = PHOTOMETRIC_MINISBLACK; + TIFFGetField(t, TIFFTAG_IMAGEWIDTH, &w); + TIFFGetField(t, TIFFTAG_IMAGELENGTH, &h); + TIFFGetFieldDefaulted(t, TIFFTAG_SAMPLESPERPIXEL, &spp); + TIFFGetFieldDefaulted(t, TIFFTAG_BITSPERSAMPLE, &bps); + TIFFGetFieldDefaulted(t, TIFFTAG_SAMPLEFORMAT, &fmtfield); + TIFFGetFieldDefaulted(t, TIFFTAG_PLANARCONFIG, &planar); + TIFFGetFieldDefaulted(t, TIFFTAG_PHOTOMETRIC, &photo); + const bool is_float = (bps == 32 && fmtfield == SAMPLEFORMAT_IEEEFP); + const bool is_int = (bps == 16 && fmtfield == SAMPLEFORMAT_INT); + const bool is_uint = ((bps == 8 || bps == 16) && fmtfield == SAMPLEFORMAT_UINT); + if(!is_float && !is_int && !is_uint) { m.fallback = 1; return 1; } + m.sformat = is_float ? 2 : (is_int ? 1 : 0); + // only layouts whose scanline is exactly w*spp*bytes are safe for the flat copy below. + // Palette needs a colormap lookup and subsampled YCbCr packs fewer bytes than w*px, which + // would make the memcpy read past the scanline buffer. Both go back to the R reader. + if(photo != PHOTOMETRIC_MINISBLACK && photo != PHOTOMETRIC_RGB && + photo != PHOTOMETRIC_MINISWHITE && photo != PHOTOMETRIC_PALETTE) { m.fallback = 1; return 1; } + // a palette image stores one index per pixel and expands to RGB, which is what readTIFF does + if(photo == PHOTOMETRIC_PALETTE && (spp != 1 || bps > 8)) { m.fallback = 1; return 1; } + // every layout check must happen here, in the probe: the batch is rejected before the buffer + // is allocated, and a check in the decode pass would come too late to stop it + if(planar == PLANARCONFIG_SEPARATE && TIFFIsTiled(t)) { m.fallback = 1; return 1; } + { + const size_t bytes0 = (size_t)bps/8, px0 = bytes0 * spp; + if(TIFFIsTiled(t)) { + uint32_t tw0 = 0, th0 = 0; + TIFFGetField(t, TIFFTAG_TILEWIDTH, &tw0); + TIFFGetField(t, TIFFTAG_TILELENGTH, &th0); + if((size_t)TIFFTileSize(t) != (size_t)tw0*th0*px0) { m.fallback = 1; return 1; } + } else if(planar == PLANARCONFIG_CONTIG && (size_t)TIFFScanlineSize(t) != (size_t)w*px0) { + m.fallback = 1; + return 1; + } else if(planar == PLANARCONFIG_SEPARATE && (size_t)TIFFScanlineSize(t) != (size_t)w*bytes0) { + m.fallback = 1; + return 1; + } + } + if(dst == NULL) { + m.h = (int)h; m.w = (int)w; m.bits = (int)bps; + m.c = (photo == PHOTOMETRIC_PALETTE) ? 3 : (int)spp; + // A cloud-optimised GeoTIFF stores its overviews as extra directories flagged + // FILETYPE_REDUCEDIMAGE. Counting those as pages would refuse every COG, so only + // full-resolution directories count; load_tiff() stacks genuine pages and this path cannot. + { + int pages = 0; + tdir_t cur = TIFFCurrentDirectory(t); + TIFFSetDirectory(t, 0); + do { + uint32_t sub = 0; + TIFFGetFieldDefaulted(t, TIFFTAG_SUBFILETYPE, &sub); + // overviews carry FILETYPE_REDUCEDIMAGE and internal masks FILETYPE_MASK; neither is + // a page, and a mask directory is present in essentially every modern COG + if(!(sub & (FILETYPE_REDUCEDIMAGE | FILETYPE_MASK))) pages++; + } while(TIFFReadDirectory(t)); + TIFFSetDirectory(t, cur); + m.pages = pages; + } + return 0; + } + const bool palette = (photo == PHOTOMETRIC_PALETTE); + const int want_c = palette ? 3 : (int)spp; + if((int)h != m.h || (int)w != m.w || want_c != m.c || (int)bps != m.bits) { + set_err(m, "file changed between probe and decode"); + return 1; + } + // the palette expansion needs somewhere to put the indices before overwriting dst with RGB + std::vector idx; + unsigned char* target = dst; + if(palette) { + idx.resize((size_t)h*w); + target = idx.data(); + } + const size_t bytes = (size_t)bps/8; + const size_t px = bytes * spp; + if(planar == PLANARCONFIG_SEPARATE) { + // one plane per sample: read each separately and interleave into the HWC buffer + const tmsize_t ss = TIFFScanlineSize(t); + std::vector buf(ss); + for(uint16_t sm = 0; sm < spp; sm++) { + for(uint32_t r = 0; r < h; r++) { + if(TIFFReadScanline(t, buf.data(), r, sm) < 0) { set_err(m, tiff_msg[0] ? tiff_msg : "TIFFReadScanline failed"); return 1; } + unsigned char* dstrow = target + (size_t)r*w*px + (size_t)sm*bytes; + for(uint32_t col = 0; col < w; col++) + std::memcpy(dstrow + (size_t)col*px, buf.data() + (size_t)col*bytes, bytes); + } + } + return 0; + } + if(TIFFIsTiled(t)) { + uint32_t tw = 0, th = 0; + TIFFGetField(t, TIFFTAG_TILEWIDTH, &tw); + TIFFGetField(t, TIFFTAG_TILELENGTH, &th); + const tmsize_t tsz = TIFFTileSize(t); + std::vector tile(tsz); + for(uint32_t y = 0; y < h; y += th) { + for(uint32_t x = 0; x < w; x += tw) { + if(TIFFReadTile(t, tile.data(), x, y, 0, 0) < 0) { set_err(m, tiff_msg[0] ? tiff_msg : "TIFFReadTile failed"); return 1; } + const uint32_t ny = (y + th > h) ? h - y : th; + const uint32_t nx = (x + tw > w) ? w - x : tw; + for(uint32_t r = 0; r < ny; r++) + std::memcpy(target + ((size_t)(y + r)*w + x)*px, tile.data() + (size_t)r*tw*px, (size_t)nx*px); + } + } + } else { + const tmsize_t ss = TIFFScanlineSize(t); + std::vector buf(ss); + for(uint32_t r = 0; r < h; r++) { + if(TIFFReadScanline(t, buf.data(), r, 0) < 0) { set_err(m, tiff_msg[0] ? tiff_msg : "TIFFReadScanline failed"); return 1; } + std::memcpy(target + (size_t)r*w*px, buf.data(), (size_t)w*px); + } + } + if(palette) { + uint16_t *rmap = NULL, *gmap = NULL, *bmap = NULL; + if(!TIFFGetField(t, TIFFTAG_COLORMAP, &rmap, &gmap, &bmap)) { + set_err(m, "palette image without a colormap"); + return 1; + } + // libtiff hands back 16-bit colormap entries; readTIFF scales them to 8 bits the same way + for(size_t i = 0, n = (size_t)h*w; i < n; i++) { + const unsigned v = idx[i]; + dst[i*3 + 0] = (unsigned char)(rmap[v] >> 8); + dst[i*3 + 1] = (unsigned char)(gmap[v] >> 8); + dst[i*3 + 2] = (unsigned char)(bmap[v] >> 8); + } + } + // MINISWHITE is not inverted here: both GDAL and tiff::readTIFF return the stored samples and + // treat the photometric as a display hint, and for training data the stored values are what + // is wanted. Inverting would make this path the only one that disagrees. + return 0; +} + +#endif /* HAVE_TIFF */ + +static int dispatch(const char* path, Meta& m, unsigned char* dst) { + switch(m.fmt) { +#ifdef HAVE_JPEG + case FMT_JPEG: return jpeg_probe_decode(path, m, dst); +#endif +#ifdef HAVE_PNG + case FMT_PNG: return png_probe_decode(path, m, dst); +#endif +#ifdef HAVE_TIFF + case FMT_TIFF: return tiff_probe_decode(path, m, dst); +#endif + case FMT_UNKNOWN: set_err(m, "unrecognised file format (not JPEG, PNG or TIFF)"); return 1; + default: m.fallback = 1; return 1; // built without that library: the R reader takes it + } +} + +// [[Rcpp::export]] +Rcpp::LogicalVector cito_decoder_available() { + Rcpp::LogicalVector v = Rcpp::LogicalVector::create( +#ifdef HAVE_JPEG + true, +#else + false, +#endif +#ifdef HAVE_PNG + true, +#else + false, +#endif +#ifdef HAVE_TIFF + true +#else + false +#endif + ); + v.names() = Rcpp::CharacterVector::create("jpeg", "png", "tiff"); + return v; +} + +// [[Rcpp::export]] +Rcpp::List cito_decode_batch(Rcpp::CharacterVector paths, int nthreads = 1) { + if(nthreads < 1) nthreads = 1; + const int n = paths.size(); + std::vector p(n); + for(int i = 0; i < n; i++) p[i] = Rcpp::as(paths[i]); + + // libtiff's handlers are process-global, so whatever tiff/terra installed has to be put back. + // An R allocation can longjmp past a C++ destructor, so the guard must never span one: it is + // scoped tightly around the two phases that actually call libtiff. +#ifdef HAVE_TIFF + struct TiffHandlerGuard { + TIFFErrorHandler e, w; + TiffHandlerGuard() : e(TIFFSetErrorHandler(tiff_err_handler)), w(TIFFSetWarningHandler(NULL)) {} + ~TiffHandlerGuard() { TIFFSetErrorHandler(e); TIFFSetWarningHandler(w); } + }; +#define CITO_TIFF_GUARD TiffHandlerGuard cito_tiff_guard_ +#else +#define CITO_TIFF_GUARD do {} while(0) +#endif + + std::vector meta(n); + { + CITO_TIFF_GUARD; + for(int i = 0; i < n; i++) { + std::memset(&meta[i], 0, sizeof(Meta)); + bool openable = false; + meta[i].fmt = sniff(p[i].c_str(), &openable); + if(!openable) set_err(meta[i], "cannot open file"); + else dispatch(p[i].c_str(), meta[i], NULL); + } + } + + // anything the compiled path cannot read safely goes back to the R reader for the whole batch + for(int i = 0; i < n; i++) { + if(meta[i].fallback || (!meta[i].status && meta[i].pages > 1)) + return Rcpp::List::create(Rcpp::_["ok"] = false, Rcpp::_["fallback"] = true); + } + + Rcpp::CharacterVector errs(n); + int first_ok = -1; + for(int i = 0; i < n; i++) { + errs[i] = meta[i].status ? meta[i].msg : ""; + if(!meta[i].status && first_ok < 0) first_ok = i; + } + if(first_ok < 0) + return Rcpp::List::create(Rcpp::_["ok"] = false, Rcpp::_["fallback"] = false, + Rcpp::_["errors"] = errs); + + // the reference geometry is the modal one, so one odd file cannot mark every other as mismatched + std::map, int> counts; + for(int i = 0; i < n; i++) + if(!meta[i].status) counts[{meta[i].h, meta[i].w, meta[i].c, meta[i].bits*8 + meta[i].sformat}]++; + // seeded from the first readable file so that a tie goes to it rather than to whichever + // geometry happens to sort first + std::array modal{meta[first_ok].h, meta[first_ok].w, meta[first_ok].c, + meta[first_ok].bits*8 + meta[first_ok].sformat}; + int best_count = counts[modal]; + for(const auto& kv : counts) + if(kv.second > best_count) { best_count = kv.second; modal = kv.first; } + int best = first_ok; + for(int i = 0; i < n; i++) + if(!meta[i].status && + std::array{meta[i].h, meta[i].w, meta[i].c, + meta[i].bits*8 + meta[i].sformat} == modal) { best = i; break; } + const Meta& ref = meta[best]; + for(int i = 0; i < n; i++) { + if(meta[i].status) continue; + if(meta[i].h != ref.h || meta[i].w != ref.w || meta[i].c != ref.c) + std::snprintf(meta[i].msg, sizeof(meta[i].msg), + "geometry %dx%dx%d does not match the batch (%dx%dx%d)", + meta[i].h, meta[i].w, meta[i].c, ref.h, ref.w, ref.c); + else if(meta[i].bits != ref.bits || meta[i].sformat != ref.sformat) { + // Both are reachable at 16 bits: signed elevation and unsigned reflectance in one folder + // is ordinary. Decoding them under one interpretation rescales half the batch silently, + // and which half depends on the shuffle, so the same file would differ between epochs. + return Rcpp::List::create(Rcpp::_["ok"] = false, Rcpp::_["fallback"] = true); + } + else continue; + meta[i].status = 1; + errs[i] = meta[i].msg; + } + + const size_t per = (size_t)ref.h * ref.w * ref.c * (ref.bits/8); + Rcpp::RawVector out(per * n); + unsigned char* base = (unsigned char*) RAW(out); + + { + CITO_TIFF_GUARD; + #ifdef _OPENMP + #pragma omp parallel for num_threads(nthreads) schedule(dynamic) + #endif + for(int i = 0; i < n; i++) { + if(meta[i].status) continue; + // an exception leaving an OpenMP structured block is std::terminate, so a bad_alloc from + // a hostile tile geometry has to be caught here + try { + dispatch(p[i].c_str(), meta[i], base + (size_t)i * per); + } catch(const std::exception& e) { + set_err(meta[i], e.what()); + } catch(...) { + set_err(meta[i], "unknown error while decoding"); + } + } + } + for(int i = 0; i < n; i++) { + // fallback can only be decided in the probe; reaching it here means the file changed + if(meta[i].fallback && !meta[i].status) set_err(meta[i], "file changed between probe and decode"); + if(meta[i].status) errs[i] = meta[i].msg[0] ? meta[i].msg : "decoding failed"; + } + + return Rcpp::List::create( + Rcpp::_["ok"] = true, + Rcpp::_["fallback"] = false, + Rcpp::_["data"] = out, + Rcpp::_["dim"] = Rcpp::IntegerVector::create(n, ref.h, ref.w, ref.c), + Rcpp::_["bits"] = ref.bits, + Rcpp::_["sformat"] = ref.sformat, + Rcpp::_["errors"] = errs); +} diff --git a/tests/testthat/fixtures/cog_u16.tif b/tests/testthat/fixtures/cog_u16.tif new file mode 100644 index 0000000..99d5ae1 Binary files /dev/null and b/tests/testthat/fixtures/cog_u16.tif differ diff --git a/tests/testthat/fixtures/dem_int16.tif b/tests/testthat/fixtures/dem_int16.tif new file mode 100644 index 0000000..761dd92 Binary files /dev/null and b/tests/testthat/fixtures/dem_int16.tif differ diff --git a/tests/testthat/fixtures/float32.tif b/tests/testthat/fixtures/float32.tif new file mode 100644 index 0000000..2b5646b Binary files /dev/null and b/tests/testthat/fixtures/float32.tif differ diff --git a/tests/testthat/fixtures/masked.tif b/tests/testthat/fixtures/masked.tif new file mode 100644 index 0000000..38c7217 Binary files /dev/null and b/tests/testthat/fixtures/masked.tif differ diff --git a/tests/testthat/fixtures/miniswhite.tif b/tests/testthat/fixtures/miniswhite.tif new file mode 100644 index 0000000..70f68e2 Binary files /dev/null and b/tests/testthat/fixtures/miniswhite.tif differ diff --git a/tests/testthat/fixtures/multipage.tif b/tests/testthat/fixtures/multipage.tif new file mode 100644 index 0000000..87dd54b Binary files /dev/null and b/tests/testthat/fixtures/multipage.tif differ diff --git a/tests/testthat/fixtures/palette.tif b/tests/testthat/fixtures/palette.tif new file mode 100644 index 0000000..e3c20be Binary files /dev/null and b/tests/testthat/fixtures/palette.tif differ diff --git a/tests/testthat/fixtures/separate_u16.tif b/tests/testthat/fixtures/separate_u16.tif new file mode 100644 index 0000000..b6105da Binary files /dev/null and b/tests/testthat/fixtures/separate_u16.tif differ diff --git a/tests/testthat/fixtures/ycbcr.tif b/tests/testthat/fixtures/ycbcr.tif new file mode 100644 index 0000000..a239eae Binary files /dev/null and b/tests/testthat/fixtures/ycbcr.tif differ diff --git a/tests/testthat/test-decode.R b/tests/testthat/test-decode.R new file mode 100644 index 0000000..35edfcf --- /dev/null +++ b/tests/testthat/test-decode.R @@ -0,0 +1,241 @@ +library(cito) +source("utils.R") + +# On a CRAN Linux build all three libraries are present, so the compiled decoder is the path +# that actually runs. Nothing else in the suite asserts that it agrees with the R readers it +# replaced, and a disagreement there changes the number of input channels a model is fitted with. + +GREY_TRNS_B64 = "iVBORw0KGgoAAAANSUhEUgAAAAIAAAACCAAAAABX3VL4AAAAAnRSTlMAAHaTzTgAAAAOSURBVHicY+ASYRDRAAAA2ABbxl+YxwAAAABJRU5ErkJggg==" + +# no base64 dependency in the package, and the fixture is 94 bytes +base64enc_decode = function(x) { + alpha = strsplit("ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789+/", "")[[1]] + chars = strsplit(sub("=+$", "", x), "")[[1]] + bits = unlist(lapply(match(chars, alpha) - 1L, function(v) as.integer(intToBits(v))[6:1])) + n = length(bits) %/% 8L + as.raw(vapply(seq_len(n), function(i) sum(bits[((i-1L)*8L+1L):(i*8L)] * 2L^(7:0)), numeric(1))) +} + +make_fixtures = function(dir) { + set.seed(42) + rgb = array(runif(24*20*3), dim = c(24, 20, 3)) + grey = matrix(runif(24*20), 24, 20) + writeables = list( + "rgb.jpg" = function(f) jpeg::writeJPEG(rgb, f, quality = 1), + "grey.jpg" = function(f) jpeg::writeJPEG(grey, f, quality = 1), + "rgb.png" = function(f) png::writePNG(rgb, f), + "grey.png" = function(f) png::writePNG(grey, f), + "greyalpha.png" = function(f) png::writePNG(array(c(grey, matrix(1, 24, 20)), dim = c(24, 20, 2)), f), + "rgba.png" = function(f) png::writePNG(array(c(rgb, matrix(1, 24, 20)), dim = c(24, 20, 4)), f), + "rgb8.tif" = function(f) tiff::writeTIFF(rgb, f, bits.per.sample = 8L), + "rgb16.tif" = function(f) tiff::writeTIFF(rgb, f, bits.per.sample = 16L), + "grey16.tif" = function(f) tiff::writeTIFF(grey, f, bits.per.sample = 16L) + ) + for(nm in names(writeables)) writeables[[nm]](file.path(dir, nm)) + # a greyscale PNG carrying a tRNS chunk decodes to grey+alpha, and must not then be + # expanded to RGB. png::writePNG cannot write tRNS, so the file is stored as bytes. + writeBin(base64enc_decode(GREY_TRNS_B64), file.path(dir, "grey_trns.png")) + sort(list.files(dir, full.names = TRUE)) +} + +reference = function(path) { + if(grepl("\\.tiff?$", path)) { + a = suppressWarnings(tiff::readTIFF(path)) + if(length(dim(a)) == 2) dim(a) = c(dim(a), 1L) + } else { + a = torchvision::base_loader(path) + } + torch::torch_tensor(a)$permute(c(3, 1, 2)) +} + +testthat::test_that("the compiled decoder matches the R readers", { + skip_if_no_torch() + testthat::skip_if(!any(cito:::cito_decoder_available()), "built without any image library") + + dir = file.path(tempdir(), "cito-decode") + dir.create(dir, showWarnings = FALSE) + on.exit(unlink(dir, recursive = TRUE)) + files = make_fixtures(dir) + + for(f in files) { + res = cito:::cito_decode_batch(f, 1L) + if(!isTRUE(res$ok)) next # format not compiled in, the R reader has it + dt = if(res$bits == 8L) "uint8" else "int16" + t = torch::torch_tensor_from_buffer(res$data, shape = res$dim, dtype = dt)$to(dtype = torch::torch_float32()) + if(res$bits == 16L) t = torch::torch_where(t < 0, t + 65536, t) + t = t$permute(c(1, 4, 2, 3))/(2^res$bits - 1) + ref = reference(f) + testthat::expect_equal(dim(t)[-1], dim(ref), info = basename(f)) + testthat::expect_lt(as.numeric(torch::torch_max(torch::torch_abs(t[1, , , ] - ref))), 1e-5) + } +}) + +testthat::test_that("the compiled path and the R fallback agree exactly", { + skip_if_no_torch() + testthat::skip_if(!any(cito:::cito_decoder_available()), "built without any image library") + + dir = file.path(tempdir(), "cito-decode-equiv") + dir.create(dir, showWarnings = FALSE) + on.exit(unlink(dir, recursive = TRUE)) + files = make_fixtures(dir) + # a file with no extension: the compiled path sniffs magic bytes, so the R path must too + file.copy(file.path(dir, "rgb.png"), file.path(dir, "noextension")) + files = sort(list.files(dir, full.names = TRUE)) + + y = torch::torch_tensor(matrix(1)) + read_with = function(path, compiled) { + ds = cito:::cito_dataset(path, y, data_augmentation = NULL) + if(!compiled) ds$decoder = c(jpeg = FALSE, png = FALSE, tiff = FALSE) + ds$.getbatch(1L)[[1]] + } + for(f in files) { + a = read_with(f, TRUE) + b = read_with(f, FALSE) + # a difference here means the fitted architecture depends on which libraries were installed + testthat::expect_identical(dim(a), dim(b), info = basename(f)) + testthat::expect_lt(as.numeric(torch::torch_max(torch::torch_abs(a - b))), 1e-6) + } +}) + +testthat::test_that("the decoder reports per-file errors and refuses layouts it cannot read", { + skip_if_no_torch() + testthat::skip_if(!any(cito:::cito_decoder_available()), "built without any image library") + + dir = file.path(tempdir(), "cito-decode-err") + dir.create(dir, showWarnings = FALSE) + on.exit(unlink(dir, recursive = TRUE)) + writeLines("not an image", file.path(dir, "notes.txt")) + writeBin(as.raw(c(0xFF, 0xD8, 0xFF, 0xE0, 1, 2, 3, 4)), file.path(dir, "truncated.jpg")) + + testthat::expect_false(cito:::cito_decode_batch(file.path(dir, "missing.jpg"), 1L)$ok) + testthat::expect_false(cito:::cito_decode_batch(file.path(dir, "notes.txt"), 1L)$ok) + testthat::expect_false(cito:::cito_decode_batch(file.path(dir, "truncated.jpg"), 1L)$ok) + + # a float TIFF is a layout the compiled path hands back rather than guessing at + tiff::writeTIFF(array(runif(8*8*3), dim = c(8, 8, 3)), file.path(dir, "f.tif"), bits.per.sample = 32L) + testthat::expect_false(cito:::cito_decode_batch(file.path(dir, "f.tif"), 1L)$ok) +}) + +testthat::test_that("cnn() reads a folder whichever decoder is used", { + skip_if_no_torch() + dir = file.path(tempdir(), "cito-decode-cnn") + dir.create(dir, showWarnings = FALSE) + on.exit(unlink(dir, recursive = TRUE)) + set.seed(1) + for(i in 1:8) jpeg::writeJPEG(array(runif(16*16*3), dim = c(16, 16, 3)), + file.path(dir, sprintf("i%02d.JPG", i)), quality = 1) + arch = create_architecture(conv(2), maxPool(), linear(2)) + testthat::expect_error({ + m = cnn(X = dir, Y = runif(8), architecture = arch, loss = "mse", epochs = 1L, + batchsize = 4L, plot = FALSE, verbose = FALSE) + .p = predict(m, newdata = dir) + }, NA) +}) + + +# The layouts remote-sensing data actually arrives in. These are files rather than generated +# fixtures because writing them needs GDAL, which cito does not depend on. Between them they +# cover the cases that cost the most if they regress: a cloud-optimised GeoTIFF, whose overviews +# are extra directories and must not be mistaken for pages; float32, which is what per-pixel +# embedding products such as TESSERA distribute; signed integers, which elevation models use; +# band-interleaved planes; a genuine multi-page file; and subsampled YCbCr, whose scanlines are +# shorter than width times samples and once caused a heap over-read. +eo_fixture = function(name) testthat::test_path("fixtures", name) + +testthat::test_that("the remote-sensing layouts load, and both paths agree where both can read", { + skip_if_no_torch() + testthat::skip_if(!any(cito:::cito_decoder_available()), "built without any image library") + + y = torch::torch_tensor(matrix(1)) + read_with = function(path, compiled) { + ds = cito:::cito_dataset(path, y, data_augmentation = NULL) + if(!compiled) ds$decoder = c(jpeg = FALSE, png = FALSE, tiff = FALSE) + ds$.getbatch(1L)[[1]] + } + + # shape and sample format the compiled path must produce + expected = list(cog_u16.tif = c(4L, 32L, 32L), float32.tif = c(4L, 32L, 32L), + dem_int16.tif = c(1L, 32L, 32L), separate_u16.tif = c(4L, 32L, 32L)) + for(nm in names(expected)) { + t = read_with(eo_fixture(nm), TRUE) + testthat::expect_identical(dim(t)[-1], expected[[nm]], info = nm) + testthat::expect_true(all(is.finite(as.numeric(t$reshape(-1)))), info = nm) + } + + # a COG must not be mistaken for a multi-page file, and both paths must read it the same + for(nm in c("cog_u16.tif", "float32.tif", "multipage.tif", "masked.tif", + "palette.tif", "miniswhite.tif")) { + a = read_with(eo_fixture(nm), TRUE) + b = read_with(eo_fixture(nm), FALSE) + testthat::expect_identical(dim(a), dim(b), info = nm) + testthat::expect_lt(as.numeric(torch::torch_max(torch::torch_abs(a - b))), 1e-6) + } + + # float32 and signed samples are physical quantities and must not be rescaled to [0, 1] + f32 = read_with(eo_fixture("float32.tif"), TRUE) + testthat::expect_lt(as.numeric(f32$min()), -1) + dem = read_with(eo_fixture("dem_int16.tif"), TRUE) + testthat::expect_lt(as.numeric(dem$min()), -100) + testthat::expect_gt(as.numeric(dem$max()), 1000) + + # an internal mask is a 1-bit directory flagged FILETYPE_MASK. Counting it as a page would + # refuse essentially every modern COG, and readTIFF(all = TRUE) chokes on the 1-bit depth. + testthat::expect_identical(dim(read_with(eo_fixture("masked.tif"), TRUE))[-1], c(1L, 32L, 32L)) + + # a palette image expands to RGB through its colormap, and MINISWHITE returns the stored + # samples rather than inverting them, because that is what GDAL and readTIFF both do + testthat::expect_identical(dim(read_with(eo_fixture("palette.tif"), TRUE))[-1], c(3L, 32L, 32L)) + testthat::expect_identical(dim(read_with(eo_fixture("miniswhite.tif"), TRUE))[-1], c(1L, 32L, 32L)) + + # subsampled YCbCr is handed back rather than copied past the end of its scanline buffer + testthat::expect_true(isTRUE(cito:::cito_decode_batch(eo_fixture("ycbcr.tif"), 1L)$fallback)) + # a genuine multi-page file goes to the R reader, which stacks the pages + testthat::expect_true(isTRUE(cito:::cito_decode_batch(eo_fixture("multipage.tif"), 1L)$fallback)) +}) + + +testthat::test_that("a float32 batch survives garbage collection", { + skip_if_no_torch() + testthat::skip_if(!any(cito:::cito_decoder_available()), "built without any image library") + + # torch_tensor_from_buffer() borrows the raw vector, and $to(dtype) copies only when the + # dtype changes -- for a float32 file it does not, so without an explicit clone the tensor + # ends up aliasing memory R has collected, and the batch fills with NaN and +/-Inf. + f = eo_fixture("float32.tif") + ds = cito:::cito_dataset(rep(f, 4L), torch::torch_tensor(matrix(runif(4))), data_augmentation = NULL) + b = ds$.getbatch(1:4)[[1]] + before = as.numeric(torch::torch_max(torch::torch_abs(b))) + invisible(lapply(1:100, function(i) runif(50000))) + gc() + x = as.array(b) + testthat::expect_false(any(is.nan(x))) + testthat::expect_true(all(is.finite(x))) + testthat::expect_equal(max(abs(x)), before, tolerance = 1e-6) +}) + +testthat::test_that("signed and unsigned samples are never mixed in one batch", { + skip_if_no_torch() + testthat::skip_if(!any(cito:::cito_decoder_available()), "built without any image library") + + # a signed elevation model and an unsigned reflectance band in one folder is ordinary, and + # reading them under one interpretation rescales half the batch. Which half would depend on + # the shuffle, so the same file would differ between epochs. + # same geometry and bit depth, differing only in sample format, so that the sample-format + # check is what decides and not the geometry check + u16 = file.path(tempdir(), "refl_u16.tif") + on.exit(unlink(u16)) + tiff::writeTIFF(matrix(seq(0, 1, length.out = 32*32), 32, 32), u16, bits.per.sample = 16L) + res = cito:::cito_decode_batch(c(eo_fixture("dem_int16.tif"), u16), 1L) + testthat::expect_false(isTRUE(res$ok)) + testthat::expect_true(isTRUE(res$fallback)) +}) + +testthat::test_that("layouts the R reader gets wrong are refused, not handed on", { + skip_if_no_torch() + y = torch::torch_tensor(matrix(1)) + # readTIFF returns different pixels on repeated reads of a subsampled YCbCr file, sometimes + # +/-Inf. The compiled path refuses it; the fallback must refuse it too rather than inherit it. + ds = cito:::cito_dataset(eo_fixture("ycbcr.tif"), y, data_augmentation = NULL) + ds$decoder = c(jpeg = FALSE, png = FALSE, tiff = FALSE) + testthat::expect_error(ds$.getbatch(1L), "YCbCr") +})