setara dengan titik paket sp dalam poligon / overing menggunakan sf

16

Saya memigrasi kode dari paket sp ke paket sf yang lebih baru. Kode saya sebelumnya, saya memiliki poligon SpatialDataFrame (censimentoMap) dan SpatialPointDataFrame (indirizzi.sp) dan saya mendapatkan id sel poligon ("Cell110") untuk setiap titik yang diletakkan dengan instruksi di bawah ini:

points.data <- over(indirizzi.sp, censimentoMap[,"Cell110"])

Sebenarnya saya membuat dua objek sf:

shape_sf <- st_read(dsn = shape_dsn) shape_sf <- st_transform(x=shape_sf, crs=crs_string) dan indirizzi_sf = st_as_sf(df, coords = c("lng", "lat"), crs = crs_string)

Dan saya mencari setara sf dari instruksi di atas ... Migth it be:

ids<-sapply(st_intersects(x=indirizzi_sf,y=shshape_sfpeCrif), function(z) if (length(z)==0) NA_integer_ else z[1]) cell_ids <- shape_sf[ids,"Cell110"]

Giorgio Spedicato
sumber

Jawaban:

20

Anda bisa mendapatkan hasil yang sama dengan menggunakan st_join: Pertama buat poligon demo dan beberapa poin dengan sf.

library(sf)
library(magrittr)

poly <- st_as_sfc(c("POLYGON((0 0 , 0 1 , 1 1 , 1 0, 0 0))")) %>% 
  st_sf(ID = "poly1")    

pts <- st_as_sfc(c("POINT(0.5 0.5)",
                   "POINT(0.6 0.6)",
                   "POINT(3 3)")) %>%
  st_sf(ID = paste0("point", 1:3))

sekarang lihat hasilnya menggunakan lebih pada objek sp

over(as(pts, "Spatial"), as(polys, "Spatial"))
>#      ID
># 1 poly1
># 2 poly1
># 3  <NA>

sekarang setara dengan sf st_join

st_join(pts, poly, join = st_intersects)
># Simple feature collection with 3 features and 2 fields
># geometry type:  POINT
># dimension:      XY
># bbox:           xmin: 0.5 ymin: 0.5 xmax: 3 ymax: 3
># epsg (SRID):    NA
># proj4string:    NA
>#     ID.x  ID.y               .
># 1 point1 poly1 POINT (0.5 0.5)
># 2 point2 poly1 POINT (0.6 0.6)
># 3 point3  <NA>     POINT (3 3)

atau untuk hasil yang sama persis

as.data.frame(st_join(pts, poly, join = st_intersects))[2] %>% setNames("ID")

>#    ID
># 1 poly1
># 2 poly1
># 3  <NA>
sebdalgarno
sumber