R Under development (unstable) (2023-08-09 r84924 ucrt) -- "Unsuffered Consequences" Copyright (C) 2023 The R Foundation for Statistical Computing Platform: x86_64-w64-mingw32/x64 R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > # Create 'stars' object > set.seed(1331) > suppressPackageStartupMessages(library(stars)) > volcano = rbind(volcano, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA) # add NA rows > d = st_dimensions(x = 1:ncol(volcano), y = 1:nrow(volcano)) > (r = st_as_stars(t(volcano))) stars object with 2 dimensions and 1 attribute attribute(s): Min. 1st Qu. Median Mean 3rd Qu. Max. NA's A1 94 108 124 130.1879 150 195 732 dimension(s): from to offset delta point x/y X1 1 61 0 1 FALSE [x] X2 1 99 0 1 FALSE [y] > r = st_set_dimensions(r, 1, offset = 0, delta = 1) > r = st_set_dimensions(r, 2, offset = nrow(volcano), delta = -1) > > # Create points > pnt = st_sample(st_as_sfc(st_bbox(r)), 100) > pnt = st_as_sf(pnt) > > # Extract - 'st_join' > x = st_join(pnt, st_as_sf(r)) > > # Extract - 'st_extract' > y = st_extract(r, pnt) > > # check there are NA's: > any(is.na(x)) [1] TRUE > # Compare > all.equal(x$A1, y[[1]]) [1] TRUE > > ################################################ > if (FALSE) { + + ## tic: segfaults + # check equal results with stars_proxy: + #x = st_extract(stars:::st_as_stars_proxy(r), pnt) + #all.equal(x$A1, y[[1]]) + #all.equal(x, y) + + r = c(r, 2*r, 10*r) + x = st_join(pnt, st_as_sf(r)) + y = st_as_sf(st_extract(r, pnt)) + all.equal(x, y) + + ## tic: segfaults + #x = st_extract(stars:::st_as_stars_proxy(merge(r)), pnt) + #all.equal(st_as_sf(x), y) + + tif = system.file("tif/L7_ETMs.tif", package = "stars") + xp = read_stars(tif, proxy = TRUE) + xm = read_stars(tif, proxy = FALSE) + pts = st_sample(st_as_sfc(st_bbox(xp)), 10) + pts = c(pts, st_as_sfc("POINT(0 0)"), pts) + em = st_extract(xm, pts) + if (utils::packageVersion("sf") >= "0.9-7") { + ep = st_extract(xp, pts) + print(all.equal(ep, em, check.attributes = TRUE)) + } + + # two-attribute objects: + library(stars) + tif = system.file("tif/L7_ETMs.tif", package = "stars") + xp = read_stars(c(tif, tif), proxy = TRUE) + xm = read_stars(c(tif, tif), proxy = FALSE) + pts = st_sample(st_as_sfc(st_bbox(xp)), 10) + pts = c(pts, st_as_sfc("POINT(0 0)"), pts) + em = st_extract(xm, pts) + if (utils::packageVersion("sf") >= "0.9-7") { + ep = st_extract(xp, pts) + print(all.equal(ep, em, check.attributes = TRUE)) + } + + # single-attribute, single raster objects: + tif1 = paste0(tempfile(), ".tif") + write_stars(xm[1,,,1], "x.tif") + xp = read_stars("x.tif", proxy = TRUE) + xm = read_stars("x.tif", proxy = FALSE) + em = st_extract(xm, pts) + if (utils::packageVersion("sf") >= "0.9-7") { + ep = st_extract(xp, pts) + print(all.equal(ep, em, check.attributes = TRUE)) + } + + # multiple-file attributes: + x = c( + "avhrr-only-v2.19810901.nc", + "avhrr-only-v2.19810902.nc", + "avhrr-only-v2.19810903.nc", + "avhrr-only-v2.19810904.nc", + "avhrr-only-v2.19810905.nc", + "avhrr-only-v2.19810906.nc", + "avhrr-only-v2.19810907.nc", + "avhrr-only-v2.19810908.nc", + "avhrr-only-v2.19810909.nc" + ) + file_list = system.file(paste0("netcdf/", x), package = "starsdata") + if (!identical(file_list, "")) { + y = read_stars(file_list, quiet = TRUE) + print(y) + st_crs(y) = "OGC:CRS84" + pts = st_sample(st_as_sfc(st_bbox(y)), 10) + em = st_extract(y, pts) + + (y = read_stars(file_list, quiet = TRUE, proxy = TRUE)) + print(y) + st_crs(y) = "OGC:CRS84" + if (utils::packageVersion("sf") >= "0.9-7") { + ep = st_extract(y, pts) + print(all.equal(em, ep)) + } + } + + # nearest & bilinear comparison: + if (utils::packageVersion("sf") >= "0.9-7") { + set.seed(12331) + s = st_as_stars(matrix(rnorm(16), 4)) + pts = st_sample(st_as_sfc(st_bbox(s)), 10000, type = "regular") + s1 = st_extract(s, pts, bilinear = FALSE) + s2 = st_extract(s, pts, bilinear = TRUE) + s1$s2 = s2[[1]] + names(s1)[c(1,3)] = c("nearest", "bilinear") + print(s1[sample(10000, 5),]) + } + } > > proc.time() user system elapsed 0.85 0.15 1.01