Issues with st_intersection or st_difference to cut polygons spread outside USA

150 Views Asked by At

I have a following use case to cut/diff the polygons that spread outside the boundary. I am trying to migrate from rgdal package to sf packages as they are expiring. The code below works fine.

WorldCountry <- sf::st_read("data/modified.countries.geo.json")
countries <- c("Canada", "Mexico")
data_Map <- WorldCountry[WorldCountry$name %in% countries, ]
data_Map <- as_Spatial(data_Map)
for (i in seq(length(spa_polys), 3)) {
   spa_polys@proj4string <- data_Map@proj4string
   geo_diff <- gDifference(spa_polys[i], data_Map)
   iso@polygons[[i]] <- geo_diff@polygons[[1]]
   iso@polygons[[i]]@ID = as.character(i)
   rm(geo_diff)
 }

With rgdal version

It generates the image with cut polygons. However, while changing the code to use SF packages, it does not correctly calculates the difference for both st_intersection and st_difference.

WorldCountry <- sf::st_read("data/modified.countries.geo.json")
countries <- c("Canada", "Mexico")
data_Map <- WorldCountry[WorldCountry$name %in% countries, ]
data_Map <- as_Spatial(data_Map)
data_Map_sf <- st_as_sf(data_Map, coords = c("longitude", "latitude"), crs = 4326)

for (i in seq(length(spa_polys), 1)) {
  spa_poly_sf <- st_as_sf(spa_polys[i])
  st_crs(spa_poly_sf) <- st_crs(data_Map_sf)
  if (!st_is_longlat(spa_poly_sf)) {
    spa_poly_sf <- st_transform(spa_poly_sf, st_crs(data_Map_sf))
  }
  geo_cut <- st_intersection(spa_poly_sf, data_Map_sf)
  iso$geometry[[i]] <- geo_cut
  iso$ID[i] <- as.character(i)
  rm(geo_cut)
}

As shown in the figure, the polygons are not cut.

with SF

Any suggestions would be greatly helpful.

I have tried both methods using the intersection and difference. I am expecting it to work like it worked with rgdal.

0

There are 0 best solutions below