Faster cc_sea(): switched terra:extract to terra:relate

R/cc_outl.R

@@ -199,6 +199,13 @@ cc_outl <- function(x,
                            area = terra::expanse(ref))
         area <- area[!is.na(area$area), ]
         area <- area[!is.na(area$country), ]
+
+        ## naturalearth uses -99 for missing values?
+        area <- area[area$country != "-99", ]
+
+        ## CN-TW (Taiwan) doesn't map to a country in the GBIF database:
+        ## naturalearth dataset uses CN-TW for Taiwan, GBIF uses TW.
+        area[area$country == "CN-TW", "country"] <- "TW"
 
         # get number of records in GBIF per country as proxy for sampling intensity
         nrec <- vapply(area$country, 

R/cc_sea.R

@@ -137,7 +137,9 @@ cc_sea <- function(x,
     ## have to make sure they are the same and in this case they already are
     ## -----
     ## point-in-polygon test
-    ext_dat <- terra::extract(ref, pts)
+    ##ext_dat <- terra::extract(ref, pts)
+    ext_dat <- terra::relate(pts, ref, "intersects", pairs = TRUE,
+                             na.rm = FALSE)
     out <- !is.na(ext_dat[!duplicated(ext_dat[, 1]), 2])
     out <- data.frame(terra::geom(pts)[, c("x", "y"), drop = FALSE], out)
     colnames(out)[1:2] <- c(lon, lat)
@@ -154,7 +156,9 @@ cc_sea <- function(x,
 
     # select relevant columns
     #point in polygon test
-    ext_dat <- terra::extract(ref, pts)
+    ext_dat <- terra::relate(pts, ref, "intersects", pairs = TRUE,
+                             na.rm = FALSE) 
+    ##ext_dat <- terra::extract(ref, pts)
     out <- !is.na(ext_dat[!duplicated(ext_dat[, 1]), 2])
   }