landscapemetrics: extents do not overlap error-CodePudding

I am using the landscapemetrics package to calculate landscape metrics (e.g., total area, total edge length, etc.) from an input raster. I am using the sample_lsm function to perform these calculations within circular buffers centered on GPS points. Here is my code:

library(raster)
library(terra)
library(landscapemetrics)
library(landscapetools)

landscape <- raster("~/landscape.tif")
show(landscape)

class      : RasterLayer 
dimensions : 8940, 3863, 34535220  (nrow, ncol, ncell)
resolution : 30, 30  (x, y)
extent     : 1728090, 1843980, 1968180, 2236380  (xmin, xmax, ymin, ymax)
crs        :  proj=aea  lat_0=23  lon_0=-96  lat_1=29.5  lat_2=45.5  x_0=0  y_0=0  datum=NAD83  units=m  no_defs 
source     : landscape.tif 
names      : landscape 
values     : 0, 255  (min, max)
attributes :
        ID      COUNT Red Green Blue Opacity NLCD.Land.Cover.Class
 from:   0 7853863229   0     0    0       0          Unclassified
  to : 255          0 255   255  255     255 

points <- read.csv("~/points.csv")
points <- data.matrix(points)
apply(points, 2, range)

     longitude latitude
[1,] -75.14847 39.51212
[2,] -74.13275 41.28174

x <- sample_lsm(landscape, 
                y = points, 
                plot_id = NULL, 
                shape = "circle", 
                size = 10000, 
                what = "lsm_c_te",  
                classes_max = NULL,
                verbose = FALSE)

I am receiving the following error:

Error in .local(x, y, ...) : extents do not overlap

My interpretation of this is that my input raster (landscape) and the GPS points I provided (points) are misaligned, but I don't know why. Both of these objects were created in GIS (CRS is ESRI:102039) and loaded into R Studio. Can someone help explain how to fix this?

CodePudding user response：

It looks like your points have the longitude/latitude coordinate reference system (crs), but your raster has a different crs: proj=aea lat_0=23 lon_0=-96 lat_1=29.5 lat_2=45.5 x_0=0 y_0=0 datum=NAD83 units=m no_defs. You need to transform the point data to the coordinate reference system of the raster (not the other way around). You can do that like this:

library(terra)
landscape <- rast("~/landscape.tif")
points <- read.csv("~/points.csv")
v <- vect(points, c("longitude", "latitude"), crs=" proj=longlat")
pv <- project(v, crs(landscape))

Check to see that the points now overlap with the raster

plot(landscape); points(pv)

If you need the coordinates from SpatVector pv, you can get them like this:

newpoints <- crds(pv)