Stars: 'st_warp' doesn't work for multiple files with proxy

[kadyb]

Looks like st_warp() doesn't work with multiple files loaded as proxy. However, it works for multibands raster with proxy.

library("stars")

# multiband raster
file = system.file("tif/L7_ETMs.tif", package = "stars")
ras = read_stars(file, proxy = TRUE)

dest = st_as_stars(st_bbox(ras), dx = 285)
dest = c(dest, dest, dest, dest, dest, dest)
dest = merge(dest)

st_warp(ras, dest, method = "bilinear", use_gdal = TRUE)

# multiple files
file = system.file("tif/olinda_dem_utm25s.tif", package = "stars")
ras = read_stars(c(file, file), proxy = TRUE) # 'proxy = FALSE' works fine
ras = st_redimension(ras)

dest = st_as_stars(st_bbox(ras), dx = 285)
dest = c(dest, dest)
dest = merge(dest)

st_warp(ras, dest, method = "bilinear", use_gdal = TRUE)
#> Error in CPL_gdal_warper(source, destination, as.integer(resampling_method(options)),  : 
#>   warper: source and destination must have the same number of bands

[edzer]

Yes, ras needs to be a single, two-band file for this to work. Is that problematic?

[kadyb]

It's not a big problem, but I think it will be an improvement from a user perspective. Now, as a workaround, I can create one temporary raster with several bands and it works.

[edzer]

Thanks; that was easier than I thought initially.

[kadyb]

Great. I tested and can confirm that it works. Thanks!

BTW: I noticed that there are significant differences between executions.

system.time(st_warp(ras, dest, method = "bilinear", use_gdal = TRUE))
# times [s]: 64.6, 57.3, 46.9, 59.1, 57.4, 46.8, 55.1, 66.8

[kadyb]

Can you confirm - did you fix the performance in the last commit or the problem was on my side?

Now I get:

system.time(st_warp(ras, dest, method = "bilinear", use_gdal = TRUE))
# times [s]: 31.1, 31.0, 31.1, 31.1, 31.1

[edzer]

Maybe it was https://github.com/r-spatial/stars/commit/912bff030ec2207f460510828b73cb81b777f4be !

[kadyb]

This is also pretty interesting. Downsampling using proxy is faster than processing rasters loaded into memory. I used stars latest version from GitHub. I had 64 GB of RAM. Code and data.

# proxy = TRUE
system.time(st_warp(ras, dest, method = "average", use_gdal = TRUE))
#> 35.56 34.68 35.30 35.38 35.14

Maximum memory peak was ~6 GB.

# proxy = FALSE
system.time(st_warp(ras, dest, method = "average", use_gdal = TRUE))
#> 53.15 53.10 54.39 54.41 54.35

Maximum memory peak was ~28 GB and it looked like the memory was not released after the last iteration.

[edzer]

proxy objects may only access overviews.